如何在Matlab中的bvp4c解算器中包含分段定義的函數

Question

我在使用分段定義函數的bvp4c時遇到了問題。 我測試了代碼，並且當分段定義的函數不變時它工作正常。 問題是，在分段定義函數不恆定的區域，我得到了錯誤的結果（我知道肯定）。

有關如何處理此問題的任何想法或建議？

感謝

function bvp4 
    xlow=0; 
    xhigh=0.30; 
    solinit=bvpinit(linspace(xlow,xhigh,1000),[0 0]); 
    sol = bvp4c(@bvp4ode,@bvp4bc,solinit); 
    xint=[xlow:0.0001:xhigh]; 
    Sxint=deval(sol,xint); 
    Sxint1=abs(sqrt(Sxint)); 
    xint=[xlow:0.0001:xhigh]; 
    plot(xint,Sxint1(1,:),'r') 

    function dydx = bvp4ode(x,y) 
    So=0.00125; 
    s=1.5; 
    dydx = [y(2);  
     ((G(x)+125*f(x)*y(1)*(1+1/s^2)^0.5-1000*9.81*So*H(x))/(1000*0.5*l(x)*(f(x)/8)^0.5)-y(2)*2*(-2/3*x+0.071+2/3*0.08)*(-2/3)*b(x))/H(x)/H(x)]; 


    function res = bvp4bc(ya,yb) 
    res = [ya(1);  yb(1)]; 



    function fval = f(x) 
if  (x >= 0) && (x <= 0.08) 
    fval = 0.0187; 
elseif (x > 0.08) && (x <= 0.17) 
    fval = 0.0298; 
elseif (x > 0.17) && (x <= 0.3) 
    fval= 0.0408; 
end 


function Gval = G(xint) 
if  (xint >= 0) && (xint <= 0.08) 
    Gval = 0.1306; 
elseif (xint > 0.08) && (xint <= 0.17) 
    Gval = 0.1306; 
elseif (xint > 0.17) && (xint <= 0.3) 
    Gval = -0.0337; 
end 

function Hval = H(xint) 
if  (xint >= 0) && (xint < 0.08) 
    Hval = 0.071; 
elseif (xint >= 0.08) && (xint <= 0.17) 
    Hval = -2/3*xint+(0.071+2/3*0.08); 
elseif (xint >0.17) && (xint <= 0.3) 
    Hval = 0.011; 
end 

function bval = b(xint) 
if  (xint >= 0) && (xint < 0.08) 
    bval = 0; 
elseif (xint >= 0.08) && (xint <= 0.17) 
    bval = 1; 
elseif (xint > 0.17) && (xint <= 0.3) 
    bval= 0; 
end 


function lval = l(xint) 
if  (xint >= 0) && (xint <= 0.08) 
    lval = 0.067; 
elseif (xint > 0.08) && (xint <= 0.17) 
    lval = 0.134; 
elseif (xint > 0.17) && (xint <= 0.3) 
    lval= 1.165; 
end 

Answer 1

不可驚喜你突然跳躍脆弱的求解。

任何訂單p求解器預計至少有倍的ODE函數可以連續微分，以合理地調整網格點的網格。本地步長。任何偏差都會導致在奇異點附近出現過度且可能振盪的自適應，導致計算時間過長或步長下溢。

我看到兩種可能性來解決這個問題，使用事件（如果支持BVP）來切換模型/ ODE函數或使用提供的multipoint mechanism將積分間隔分成參數函數恆定的部分。然後，您也可以使用簡單的數組作爲參數，而不是使用許多分支的函數。