BVP4c가 알 수없는 경계를 해결합니다.

Question

4 개의 시스템을 해결하기 위해 bvp4c를 사용하려고합니다. 문제는 경계 중 하나를 알 수 없다는 것입니다.

bvp4c에서 처리 할 수 ​​있습니까? 내 코드에서 L은 내가 알지 못하는 부분이다.

아래에 오류 메시지가 표시됩니다.

function mat4bvp 
L = 8; 
solinit = bvpinit(linspace(0,L,100),@mat4init); 
sol = bvp4c(@mat4ode,@mat4bc,solinit); 
sint = linspace(0,L); 
Sxint = deval(sol,sint); 
end 

% ------------------------------------------------------------ 
function dtdpdxdy = mat4ode(s,y,L) 
Lambda = 0.3536; 
dtdpdxdy = [y(2) 
      -sin(y(1)) + Lambda*(L-s)*cos(y(1)) 
      cos(y(1)) 
      sin(y(1))]; 
end 
% ------------------------------------------------------------ 
function res = mat4bc(ya,yb,L) 
res = [ ya(1) 
     ya(2) 
     ya(3) 
     ya(4) 
     yb(1)]; 
end 
% ------------------------------------------------------------ 
function yinit = mat4init(s) 
yinit = [ cos(s) 
    0 
    0 
    0 
     ]; 
end 

불행히도 다음과 같은 오류 메시지가 표시됩니다.

>> mat4bvp 
Not enough input arguments. 

Error in mat4bvp>mat4ode (line 13) 
      -sin(y(1)) + Lambda*(L-s)*cos(y(1)) 

Error in bvparguments (line 105) 
    testODE = ode(x1,y1,odeExtras{:}); 

Error in bvp4c (line 130) 
    bvparguments(solver_name,ode,bc,solinit,options,varargin); 

Error in mat4bvp (line 4) 
sol = bvp4c(@mat4ode,@mat4bc,solinit); 

Answer 1

가변 끝점을 고정 된 끝점으로 변환하는 트릭은 시간 척도를 변경하는 것입니다. x'(t)=f(t,x(t))는 미분 방정식 인 경우 1에 0에서 t=L*s, s를 설정하고 사용하는 것이 y(s)=x(L*s)

y'(s)=L*x'(L*s)=L*f(L*s,y(s)) 

다음 트릭 대한 관련 미분 방정식을 연산은 미분 방정식의 일부에 의해에 전역 변수를 변형하는 상수 함수로 계산. 따라서 새로운 시스템은

[ y'(s), L'(s) ] = [ L(s)*f(L(s)*s,y(s)), 0 ] 

이고 추가 값은 L입니다. 이

from math import sin, cos 
import numpy as np 
from scipy.integrate import solve_bvp, odeint 
import matplotlib.pyplot as plt 

# The original function with the interval length as parameter 

def fun0(t, y, L): 
    Lambda = 0.3536; 
    #print t,y,L 
    return np.array([ y[1], -np.sin(y[0]) + Lambda*(L-t)*np.cos(y[0]), np.cos(y[0]), np.sin(y[0]) ]); 

# Wrapper function to apply both tricks to transform variable interval length to a fixed interval. 

def fun1(s,y): 
    L = y[-1]; 
    dydt = np.zeros_like(y); 
    dydt[:-1] = L*fun0(L*s, y[:-1], L); 
    return dydt; 

# Implement evaluation of the boundary condition residuals: 

def bc(ya, yb): 
    return [ ya[0],ya[1], ya[2], ya[3], yb[0] ]; 

# Define the initial mesh with 5 nodes: 

x = np.linspace(0, 1, 3) 

# This problem has multiple solutions. Try two initial guesses. 

L_a=8 
L_b=9 

y_a = odeint(lambda y,t: fun1(t,y), [0,0,0,0,L_a], x) 
y_b = odeint(lambda y,t: fun1(t,y), [0,0,0,0,L_b], x) 

# Now we are ready to run the solver. 

res_a = solve_bvp(fun1, bc, x, np.array(y_a.transpose())) 
res_b = solve_bvp(fun1, bc, x, np.array(y_b.transpose())) 

L_a = res_a.sol(0)[-1] 
L_b = res_b.sol(0)[-1] 
print "L_a=%.8f, L_b=%.8f" % (L_a,L_b) 

# Plot the two found solutions. The solution are in a spline form, use this to produce a smooth plot. 

x_plot = np.linspace(0, 1, 100) 
y_plot_a = res_a.sol(x_plot)[0] 
y_plot_b = res_b.sol(x_plot)[0] 

plt.plot(L_a*x_plot, y_plot_a, label='y_a, L=%.8f'%L_a) 
plt.plot(L_b*x_plot, y_plot_b, label='y_b, L=%.8f'%L_b) 
plt.legend() 
plt.xlabel("t") 
plt.ylabel("y") 
plt.show()