страница 1 |
|||||||||||||||||||||||||||||||||||||||||||
Похожие работы
|
Решение уравнения в точках 0,25; 0,5; 0,75; Найти точное решение оду - страница №1/1
![]() Вычислить с заданной точностью (0.001) решение уравнения в точках 0,25; 0,5; 0,75; 1. Найти точное решение ОДУ. Сравнить приближенное решение с точным. Задача Коши, согласно заданию лабораторной работы, имеет вид ![]() ![]() ![]() Найдем точное решение задачи Коши. Поскольку мы имеем линейное неоднородное уравнение второго порядка, найдем сначала общее решение линейного однородного уравнения, решая характеристическое уравнение ![]() ![]() ![]() Получим, что частное решение неоднородного уравнения имеет вид Теперь приступим к численному решению задачи Коши. Численное решение методом Рунге-Кутта в нашем случае облегчается тем, что фактически требуется задание только одной функции, поскольку ![]() ![]() Для решения задачи воспользуемся языком программирования C++, а функцию для решения методом Рунге-Кутта сделаем такой, что она будет не просто вычислять новые значения y и z, но и записывать их в старые переменные, что позволит легко использовать ее в циклах. Для обобщенности решения будем передавать ей в качестве параметра указатель на указанную выше функцию ![]() Вот как выглядит код этой функции для вычисления очередного шага методом Рунге-Кутта: double RK(double& x, double& y, double& ys, double h, double (*f)(double,double,double)) { double k1,k2,k3,k4; double m1,m2,m3,m4; k1 = h*f(x,y,ys); m1 = h*(ys); k2 = h*f(x+h/2.,y+m1/2.,m2 = ys+k1/2.); m2*=h; k3 = h*f(x+h/2.,y+m2/2., m3 = ys+k2/2.); m3*=h; k4 = h*f(x+=h, y+m3, m4 = ys+k3); m4*=h; ys+=(k1+2.*(k2+k3)+k4)/6.; return y+=(m1+2.*(m2+m3)+m4)/6.; }; Полностью автоматизируем решение нашей задачи, для чего создадим два массива — для решений в интересующих нас точках с шагом h и с шагом h/2, и организуем вычисления в цикле, со все уменьшающимся шагом. После каждого вычисления будем проверять, достигнута ли требуемая точность, и если нет — то цикл будет повторен с уменьшенным значением (перед этим следует переписать значения из массива для половинного шага в первый массив). Таким образом, расчетная часть программы имеет вид: double Y(double x, double y, double ys) { return 3+ys; } int main(int argc, const char * argv[]) { double eps = 0.0001; double yh[4], zh[4], yh2[4], zh2[4]; double h = 0.25; // Начальное значение шага int N = 4; // Количество шагов int n = 1; // Записывать каждое n-е значение // Заполняем таблицу для шага h double x = 0, y = 6, z = 2; for(int i = 1; i <= N; ++i) { RK(x,y,z,h,Y); if (i%n == 0) { yh[i/n-1] = y; zh[i/n-1] = z; } } // Цикл, пока не достигнем точности. for(;;) { // Уменьшаем шаг N *= 2; n *= 2; h /= 2.0; x = 0; y = 6; z = 2; for(int i = 1; i <= N; ++i) { RK(x,y,z,h,Y); if (i%n == 0) { yh2[i/n-1] = y; zh2[i/n-1] = z; } } // Проверяем точность bool ok = true; for(int i = 0; i < 4; ++i) { if (fabs(yh[i] - yh2[i])/15 > eps) { ok = false; break; } if (fabs(zh[i] - zh2[i])/15 > eps) { ok = false; break; } } if (ok == false) { // Точности не хватает! // Переписываем массив и считаем с новым шагом for(int i = 0; i < 4; ++i) { yh[i] = yh2[i]; zh[i] = zh2[i]; } continue; } break; } // Выводим таблицу значений для данного шага x = 0; y = 6; z = 2;
5*exp(x)-3,5*exp(x)-3-z); for(int i = 0; i < N; ++i) { RK(x,y,z,h,Y); printf(fmt,x,y,z, 5*exp(x)+1-3*x,5*exp(x)+1-3*x-y, 5*exp(x)-3,5*exp(x)-3-z); } }
x y y' y точное Погрешность y' точное Погрешность 0.000000 6.000000 2.000000 6.000000 0 2.000000 0 0.125000 6.290741 2.665741 6.290742 1.3e-006 2.665742 1.3e-006 0.250000 6.670124 3.420124 6.670127 2.94e-006 3.420127 2.94e-006 0.375000 7.149952 4.274952 7.149957 5e-006 4.274957 5e-006 0.500000 7.743599 5.243599 7.743606 7.56e-006 5.243606 7.56e-006 0.625000 8.466219 6.341219 8.466230 1.07e-005 6.341230 1.07e-005 0.750000 9.334986 7.584986 9.335000 1.46e-005 7.585000 1.46e-005 0.875000 10.369357 8.994357 10.369376 1.92e-005 8.994376 1.92e-005 1.000000 11.591384 10.591384 11.591409 2.49e-005 10.591409 2.49e-005 Как видно из расчета, при шаге 0.125 погрешность расчета существенно меньше указанной в задании 0.0001. |
ещё >> |