我已经在python中使用bvp求解器解决了以下bvp问题。
d4y/dx4= 0.00033*V/(0.000001-y)^(2) , y(0)=y'(0)=y(1)=y'(1)=0
在上面的eqn中,“ V”是使用for循环进行了更改的参数。有趣的是,对于V> Vo,上述微分方程的解应该是不稳定的。对于V> Vo,bvp求解器仍然会产生一些错误的值。一旦出现这种不稳定性,如何使求解器停止计算?
对于归一化方程式(y
和V
的小数位数]
y'''' = c/(1-y)**2
我得到c=70.099305
的临界值。对于非常小的c
,解决方案同样很小,接近y(t)=c/24*(t*(1-t))**2
。对于c
接近临界值,网格细化集中在y=0.4
附近的最大值。
c=70.099305
def f(t,u): return [u[1],u[2],u[3],c/(1-u[0])**2]
def bc(u0,u1): return [u0[0], u0[1], u1[0], u1[1]]
t = np.linspace(0,1,5);
u = np.zeros([4,len(t)])
res = solve_bvp(f,bc,t,u, tol=1e-4, max_nodes=5000)
print(res.message)
%matplotlib inline
if res.success:
plt.figure(figsize=(10,5))
t = np.linspace(0,1,502)
plt.plot(t, c/24*(t*(1-t))**2,c='y', lw=3)
plt.plot(t,res.sol(t)[0],'b')
plt.plot(res.x,res.y[0],'xr')
plt.grid(); plt.show()
蓝色-数值解,黄色-小c
的近似值]]