先用欧拉法求得一个初步的近似值,称为预报值,然后用它替代梯形法右端的yi+1再直接计算fi+1,得到校正值yi+1,这样就可以得到预报校正公式。
{
y
p
=
y
n
+
h
f
(
x
n
,
y
n
)
y
c
=
y
n
+
h
f
(
x
n
+
1
,
y
p
)
y
n
+
1
=
1
2
(
y
p
+
y
c
)
\left\{\begin{array}{l}{y_{p}=y_{n}+h f\left(x_{n}, y_{n}\right)} \\ {y_{c}=y_{n}+h f\left(x_{n+1}, y_{p}\right)} \\ {y_{n+1}=\frac{1}{2}\left(y_{p}+y_{c}\right)}\end{array}\right.
⎩⎨⎧yp=yn+hf(xn,yn)yc=yn+hf(xn+1,yp)yn+1=21(yp+yc)
#include <iostream>
#include <iomanip>
#include <vector>
#include <math.h>
using namespace std;
typedef vector<pair<double, double>> Data;
void disp(Data);
double f(double x, double y)
{
// return y - 2 * x / y; //0 1 0.1 0 1
return (y * y + y) / x; // 1,-2,
}
Data Corrected_Euler(double, double, double, double, double);
int main()
{
double x, y0, a, b, h;
cin >> x >> y0 >> h >> a >> b;
Data data = Corrected_Euler(x, y0, h, 1, 3);
disp(data);
system("pause");
return 0;
}
Data Corrected_Euler(double x, double y, double h, double a, double b)
{
Data data;
double yp, yc;
double n = (b - a) / h;
// cout << n<<endl;
for (int i = 0; i <= n; i++)
{
yp = y + h * f(x, y);
yc = y + h * f(x + h, yp);
y = (yp + yc) / 2;
data.push_back(make_pair(x, y));
x += h;
// cout << ++cnt << " " << setprecision(7) << y << endl;
}
return data;
}
void disp(Data d)
{
for (auto ele : d)
cout << setw(4) << ele.first <<" "<< setw(11) << setprecision(8) << ele.second << endl;
}