在double_scalars中遇到Runge-Kutta方法溢出

问题描述

我不是编程专家，但是我试图通过尝试边界条件的错误来对runge-kutta方法进行建模。对我而言，这很好，但是如果尝试扩展时域（即通过增加xStop向向量添加更多项），则该程序会突然爆炸并溢出。

你能告诉我为什么它这么不稳定吗？

import numpy as np
from matplotlib import pyplot as plt

def runge_kutta4(F,x,y,xStop,h):

    def rk4(F,h):
        K0 = h*F(x,y)
        K1 = h*F(x + h/2.0,y + K0/2.0)
        K2 = h*F(x + h/2.0,y + K1/2.0)
        K3 = h*F(x + h,y + K2)
        return (K0 + 2.0*K1 + 2.0*K2 + K3)/6.0
    X = []
    Y = []
    X.append(x)
    Y.append(y)
    while x < xStop:
        h = min(h,xStop - x)
        y = y + rk4(F,h)
        x = x + h
        X.append(x)
        Y.append(y)
    return np.array(X),np.array(Y)

def secant(f,x1,x2,tol,mx):
    for i in range(mx):
        xnew = x2-(x2-x1)/(f(x2)-f(x1)) * f(x2)
        if abs(xnew-x2) < tol: break
        else:
            x1 = x2
            x2 = xnew
    else:
        print('cagasteZ')
    return xnew


def initCond(u): # Initial values of [y,y’,y"];
# use 'u' if unkNown
    return np.array([0.0,0.0,u])

beta = 7.5 #CONDICION DE BORDE DEL OTRO EXTREMO

def res(u): # Boundary condition residual--see Eq. (8.3)
    X,Y = runge_kutta4(F,xStart,initCond(u),h)
    y = Y[len(Y) - 1]
    residuo = y[0] - beta   #ESTE ES EL VALOR DEL OTRO EXTREMO
    return residuo

def F(x,y): # First-order differential equations
    F = np.zeros(3)
    F[0] = y[1]
    F[1] = y[2]
    F[2] = y[1]**2-y[0]*y[2]-1 
    return F

xStart = 0.0 # Start of integration
xStop = 7.5 # End of integration
u1 = 1 # 1st trial value of unkNown init. cond.
u2 = 1.5 # 2nd trial value of unkNown init. cond.
h = 0.01 # initial step size
u = secant(res,u1,u2,0.001,1000)#
X,h)
plt.figure()
plt.plot(X,Y[:,0],'-',X,1],2],'-')
plt.xlabel('$\eta$')
plt.legend(('f','df/d$\eta$','d$^2$f/d$\eta^2$'),loc = 2)
plt.grid(True)
plt.show()

解决方法

暂无找到可以解决该程序问题的有效方法，小编努力寻找整理中！

如果你已经找到好的解决方法，欢迎将解决方案带上本链接一起发送给小编。

小编邮箱:dio#foxmail.com (将#修改为@）