Задача
Нехай надана залежність y=f(x) у вигляді таблиці, яка містить n значень. Крім того, нехай задана функція F(x), в яку входить m+1 параметрів - апроксимуюча функція:
y = F(x; a0, a1, … , am) . (1)
Задачею апроксимації є отримання таких значень параметрів апроксимуючої функції, які забезпечують найкращу згоду між апроксимуючою функцією і наданою таблицею f(x).
Апроксимуюча функція звичайно має сенс математичної моделі деякого елементу системи, а сама задача апроксимації називається у цьому випадку задачею ідентифікації параметрів моделі.
У якості критерію згоди найчастіше використовують суму квадратів відхилень:
. (2)
У такому разі говорять, що задача розв’язується методом найменших квадратів (метод НК), хоча більш правильно було б говорити про розв’язування задачі апроксимації з використанням критерію найменших квадратів.
Алгоритм
Таким чином, задача апроксимації у загальному випадку зводиться до задачі багатовимірної оптимізації. Головними методами для розв’язування такої задачі є:
- метод покоординатного спуску;
- метод градієнта (найшвидшого спуску);
- метод Ньютона (параболічної екстраполяції);
- метод вкладених алгоритмів.
Неприємною особливістю задачі оптимізації, до якої зводиться задача апроксимації, може бути, яружність цільової функції (2), яка часто спостерігається для нелінійних моделей.
Нехай тепер маємо лінійну модель:
.
і таблицю залежності y = f(x) . З наведених вище відомостей виходить, що параметри a0, a1 можна отримати шляхом розв’язування системи двох лінійних рівнянь:
g00a0 + g01a1 = b0
g10a0 + g11a1 = b1
де
Нижче наведено функцію, яка виконує лінійну апроксимацію табличних даних.
// Приклад 1
void linappr(int N, double* x, double* y, double& a0, double& a1)
{ double D, g01=0, g11=0, b0=0, b1=0;
int i;
for (i=0;i<N;i++) { g01+=x[i]; g11+=sqr(x[i]);
b0+=y[i]; b1+=x[i]*y[i];
}
D=N*g11-g01*g01;
a0=(b0*g11-b1*g01)/D;
a1=(N*b1-g01*b0)/D;
}
13.10. Розв’язування нелінійних рівнянь
Більшість реальних об’єктів та процесів моделюються рівняннями, які містять у собі різного роду нелінійності. Наприклад, при дослідженні пружинного маятника, коефіцієнт пружності пружини можна вважати постійним, але це буде істинно тільки на невеликому інтервалі розтягань пружини. Для знаходження точного рішення змінність коефіцієнту пружності треба ураховувати. Урахування нелінійностей звичайно призводить до утворення складних математичних виразів, для яких не можливо записати розв’язок задачі у явній аналітичній формі. У таких випадках нас виручають чисельні методи.
Задача розв’язування рівняння виду F(x)=0
Надано рівняння F(x)=0 і інтервал [a,b], який містить у собі шуканий розв’язок. Обчислити значення x, яке задовольняє рівнянню з наданою точністю eps. Надана точність означає, що знайдене x задовольняє умові: F(x) £ eps .
Метод бісекції
Алгоритм передбачає послідовне наближення до точного розв’язку. Вважається, що для обраного інтервалу [a,b] виконується умова: F(a)×F(b) £ 0 .
1. Обчислити точку x=(b+a)/2.
2. Якщо |f(x)|<eps алгоритм завершити. Значення x - шуканий розв’язок.
Якщо ні, перейти до наступного п.
3. Якщо F(a)×F(x) < 0 , виконати присвоєння b = x .
Якщо ні, то виконати присвоєння a = x .
4. Перейти до п. 1.
Алгоритм реалізовано у вигляді функції Bisec, що наведена нижче.
void Bisec(funx F, double a, double b, double eps,
double& x)
{ for (;;) { x = (a+b)/2.;
if (b-a<eps) break;
if (F(a)*F(x)<0) b=x; else a=x;
}
}
Першим параметром функції є параметр F типу funx. Це є покажчик на ту функцію, яка подає наше рівняння F(x) = 0 . У зверненні до функції Bisec першим параметром пишеться просто ім’я функції, яка є лівою частиною рівняння. Сам тип funx визначений у модулі syst.h таким чином:
typedef double (*funx)(double);
Давайте подивимось, як буде виглядати програма в цілому, якщо розв’язується таке рівняння:
.
Розв’язок шукається на інтервалі [a,b] = [-1, 2] з точністю 1×10-6.
// Приклад 1
#include <syst.h>
