1

3 つの微分方程式系をモデル化しようとしています。これは、アーク長sに対してパラメーター化された液滴モデルです。

方程式は次のとおりです。
dx/ds=cos(シータ)
dz/ds=sin(シータ)
(シータ)/ds=2*b+c*z-sin(シータ)/x

初期条件は、s=0 で x、z、および theta がすべて 0 であることです。d(theta)/ds の特異点を回避するために、s=0で d(theta)/ds=b という条件も設定しています。私はすでにこのコードを書いています:

[s,x]=ode23(@(s,x)drpivp(s,x,p),sspan,x0);

%where p contains two parameters and x0 contains initial angle theta, x, z values.  
%droplet ODE function:  
function drpivp = drpivp(s,x,p);  
%x(1)=theta  
%x(2)=x  
%x(3)=z  
%b is curvature at apex  
%c is capillarity constant  
b=p(1);  
c=p(2);  
drpivp=[2/p(1)+p(2)*x(3)-sin(x(1))/x(2); cos(x(1)); sin(x(1))];

これにより、スパイラルアウトするソリューションが得られます。1 つの液滴プロファイルを作成する代わりに、多数のプロファイルを作成します。もちろん、ここでは式を適切に初期化していません。s=0 でのシータに別の式を使用する方法がわからないからです。

問題は、通常の s=0 ではなく、d(theta)/ds=b という初期条件を含めるにはどうすればよいかということです。これは、matlab の組み込みソルバーを使用して可能ですか?

ありがとう。

4

1 に答える 1

1

これを行うにはいくつかの方法がありますが、最も簡単な方法は、単純に if ステートメントを数式に追加することです。

function drpivp = drpivp(s,x,p);  
%x(1)=theta  
%x(2)=x  
%x(3)=z  
%b is curvature at apex  
%c is capillarity constant  
b=p(1);  
c=p(2); 
if (s == 0)
    drpivp=[b; cos(x(1)); sin(x(1))];
else
    drpivp=[2/p(1)+p(2)*x(3)-sin(x(1))/x(2); cos(x(1)); sin(x(1))];
end
于 2013-04-15T22:53:19.223 に答える