Господа, возникло множество вопросов по методу Эйлера и реализации численного интегрирования системы дифференциальных уравнений. Удалось с горем пополам накопать исходный код для решения данной задачи, но к сожалению не все понятно в ней, где разобрался поставил коменты. Но все равно остается много неясностей для меня. Например за что отвечает функция f0, f1 и т.д. Было бы здорово если смогли бы указать на ошибки если они есть, на более правильное решение, ну и прокоментировать места которые вызывают вопросы. Спасибо!
//---------------------------------------------------------------------------
#pragma hdrstop
#include <stdio.h>
#include <math.h>
//---------------------------------------------------------------------------
#pragma argsused
double f0(double t);
double f1(double t, double y[]);
double f2(double t, double y[]);
double f3(double t, double y[]);
void eiler( double t, double y[], int n, double tn, double tk, int m, double delt);
int main(int argc, char* argv[])
{
int n; //порядок системы дифференциальных уравнений
double delt; //шаг интегрирования
double tn; //начальное значение интервала интегрирования;//
double tk; //конечное значение интервала интегрирования
int m ; //шаг интегрирования
double y[3];
double ys[3];
double t;
int j;
int i;
printf("Input n="); //Порядок интегрирования
scanf("%d" , &n);
if(n > 3)
{
return-1;
}
printf("Input delt=");
scanf("%lf" , &delt);
printf("Input tn=");
scanf("%lf" , &tn);
printf("Input tk=");
scanf("%lf" , &tk);
m = (tk-tn)/delt;
//Считываем начальные данные
for ( j= 0; j <n; j++)
{
printf("Input y^(%d",j);
printf(") =");
scanf("%lf" , &y[j]);
}
// решение
t = tn;
for ( j= 0; j <=m; j++)
{
eiler(t, y, n, tn,tk, m, delt);
t = t + delt;
printf("x= %lf", t);
for (i = 0; i <n ; i++)
{
printf(" y^%d", i);
printf("= %lf", y[i]);
}
printf("\n");
}
scanf("%s" , &n);
return 0;
}
//---------------------------------------------------------------------------
double f0(double t)
{
return(5+t*t);
}
//---------------------------------------------------------------------------
double f1(double t,double y[])
{
return(5+t*t-y[0]);
}
//---------------------------------------------------------------------------
double f2(double t, double y[])
{
return(5+t*t-3*y[1]-y[0]);
}
//---------------------------------------------------------------------------
double f3(double t, double y[])
{
return(5+t*t-2*y[2]-3*y[1]-y[0]);
}
void eiler( double t, double y[], int n, double tn, double tk, int m, double delt)
{
switch (n )
{
case 3:
{
//Порядок 3
y[2] = y[2] + delt * f3(t, y);
y[1] = y[1] + delt * y[2];
y[0] = y[0] + delt * y[1];
break;
}
case 2:
{
//Порядок 2
y[1] = y[1] + delt * f2(t, y);
y[0] = y[0] + delt * y[1];
break;
}
case 1:
{
//Порядок 1
y[0] = y[0] + delt * f1(t, y);
break;
}
default:
{
//Порядок непонятен выводим 0
y[0] = f0(t);
break;
}
}
}