/////////////////////////////////////
// Отже, обчислені 4 перших значення. Цього достатньо для початку методу
// Адамса для кроку h.
// B y0 ... y3 лежать 4 значення функцій (_НЕ_ПРОІЗВОДНИХ!!!).
// A в q0 ... q3 лежать значення _проізводних_ цих функцій, помножених на h
// q0 .. q3, а також y0 .. y3 являють собою черги з 4 елементами
again:// Обчислюємо нове значення функції Yi (Це Y [i +1])
for (j = 0; j
{//Всі прирощення
dq2 = q3 [j] - q2 [j]; dq1 = q2 [j] - q1 [j]; dq0 = q1 [j] - q0 [j];
d2q1 = dq2 - dq1; d2q0 = dq1 - dq0;
d3q0 = d2q1 - d2q0;
// нове значення функції (у ya поки що)
ya [j] = y3 [j] + (q3 [j] + (dq2/2.) + (5. * d2q1/12.) + (3. * D3q0/8.)); p>// Зрушуємо всі масиви на 1 вперед і додаємо в чергу нове
// значення функції
y0 [j] = y1 [j]; y1 [j] = y2 [j]; y2 [j] = y3 [j]; y3 [j] = ya [j];
// Просто зсуваємо q, нічого поки що ні додаючи
q0 [j] = q1 [j]; q1 [j] = q2 [j]; q2 [j] = q3 [j];
}
// У чергу в якості q3 ложим нове значення
f (y3, q3, xi);// q3 = f (xi, y3);
for (j = 0; j
// Чергове значення функції обчислено. Следующии крок
xi + = h;
// Продовжити інтегрування?
if (xi
// Якщо перший раз тут, то прорахувати ще раз з кроком h/2
if (flag == 0)
flag = 1;// Порівнювати вже буде з чим
else
{
// Не перший раз - оцінити похибку
// Зараз в y3 - значення тільки що обчисленої функції,
// а в y2 - занчение функції, обчисленої з кроком h * 2
// по відношенню до поточного
for (j = 0; j
{eps2 = fabs (((y3 [j] - y2 [j])/y2 [j]));
if (eps2> eps) break;// Якщо похибка занадто великі
}
if (j == n)// Якщо все ОК
{//Копіюємо результат
for (j = 0; j
free (k1);// Звільняємо пам'ять
return;// Повертаємося в main
}
}
// З якихось причин виходу з функції не сталося -
// тоді зменшуємо крок у 2 рази і повторюємо
// все, починаючи з методу Рунге-Кутта
h/= 2.;// Зменшити крок
goto start;// Повторити розрахунок спочатку, з новими параметрами
}
int main ()
{
double y [3], xs, xe;
int i;
y [0] = 1.; y [1] = 0.1; y [2] = 0.;// Початкові умови
xs = .0; xe = .1;// Початок інтегрування
printf ("x =% 5.3lg, y (% 4.2lg) =% 10.3lg n", xs, xs, y [0]);
for (i = 0; i <20; i + +)
{
Adams (func, y, 3, xs, xe, 10, 1.e-3);
xs + = 0.1; xe + = 0.1;
printf ("x =% 5.3lg, y (% 4.2lg) = % 10.3lg n ", xs, xs, y [0]);
}
return 0;
}
Результат рішення задачі 17 на ЕОМ
Для роботи програму необхідно скомпілювати в моделі не нижче SMALL. Використовувався компілятор Micro $ oft C 6.00 з однойменного пакета. Після запуску програма виводить наступне:
Програма чисельного інтегрування системи диференціальних
рівнянь екстраполяціонним методом Адамса
Розробник: студент гр. ПС-146
Нечаєв Леонід Володимирович
17.03.2004
Диференціальне рівняння має вигляд y'' '+ 2y'' + 3y' + y = x ^ 2 + 5
Отже, залежність y [x]:
x = 0, y (0) = 1
x = 0.1, y (0.1) = 1.01
x = 0.2, y (0.2) = 1.02
x = 0.3, y (0.3) = 1.04
x = 0.4, y (0.4) = 1.07
x = 0.5, y (0.5) = 1.11
x = 0.6, y (0.6) = 1.16
x = 0.7, y (0.7) = 1.22
x = 0.8, y (0.8) = 1.28
x = 0.9, y (0.9) = 1.37
x = +1, y (1) = 1.46
x = 1.1, y (1.1) = 1.56
x = 1.2, y (1.2) = 1.67
x = 1.3, y (1.3) = 1.79
x = 1.4, y (1.4) = 1.92
x = 1.5, y (1.5) = 2.06
x = 1.6, y (1.6) = 2.21
x = 1.7, y (1.7) = 2.36
x = 1.8, y (1.8) = 2.52
x = 1.9, y (1.9) = 2.69
x = 2, y (2) = 2.86
Висновок:
Перевіряємо рішення у програмі Mathematica 4.2. Результати, отримані з точністю до 2 знаків після комою не відрізняються від отриманих. Задача вирішена вірно. <В
Програма для вирішення завдання 30.
Умова задачі 30.
Розробити програму апроксимації функції методом найменших квадратів для моделі за таблицею результатів експерименту
X 1
X 2
Y
1
1
0
-1
-1
...
