改进的欧拉公式
流程图以及方程式如下(以下代码是按照方程式写就)
由于使用改进的欧拉公式需要已知函数表达式的具体形式,故而这里用了书上的具体的习题
习题如下:
代码如下:
# coding=gbk;
#改进的欧拉方法需要函数的具体表达式
#这里就暂且用书上的题目所给公式了
def compute_fx(x,y):
return y-2*x/y;
#进行数据初始化
#需要输入初值x0 y0 步长h 以及计算的次数N
def RawDataIn(x,y,h1,N1):
temp_list=input("请输入初值x0,y0;步长h以及总共需要计算的次数N: ").split(" ");
x=float(temp_list[0]);
y=float(temp_list[1]);
h1=float(temp_list[2]);
N1=float(temp_list[3]);
return x,y,h1,N1;
def Modified_Euler(x_initial,y_initial,h1,N1):
n=0;
y1=[]; #用来记录矫正之后的y的值
while n<N1:
n+=1;
#先计算预测 再计算校正
#这是预测
x1=x_initial+h1;
y1_pre=y_initial+h1*compute_fx(x_initial,y_initial);
#这是校正
y1_modify=y_initial+(compute_fx(x_initial, y_initial)+compute_fx(x1, y1_pre))*h1/2;
y1.append(y1_modify);
x_initial=x1;
y_initial=y1_modify;
return y1;
x0=y0=h=N=0;
x0,y0,h,N=RawDataIn(x0, y0, h, N);
result=[];
result=Modified_Euler(x0, y0, h, N);
for i in range(1,int(N)+1):
X=x0+h*i;
print("{x:.3f} --> {y:.5f}".format(x=X, y=result[i-1]));
代码运行结果:
书上结果:
遇事不决,可问春风