未知パラメーターをもつ BVP を解く
この例では bvp4c
を使用して未知パラメーターをもつ境界値問題を解く方法を説明します。
Mathieu 方程式は、区間 で次により定義されます。
.
パラメーターが の場合、境界条件は次のようになります。
,
.
ただし、これは を定数の倍数までしか決定しないため、次の特殊解を指定する 3 番目の条件が必要となります。
.
MATLAB® でこの方程式系を解くには、方程式、境界条件、および初期推定をコード化してから、境界値問題ソルバー bvp4c
を呼び出す必要があります。(ここで行ったように) 必要な関数をファイルの最後にローカル関数として含めることも、あるいは個別の名前付きファイルとして MATLAB パスのディレクトリに保存することもできます。
方程式のコード化
方程式をコード化する関数を作成します。この関数にはシグネチャ dydx = mat4ode(x,y,lambda)
がなければなりません。ここでは以下のようになります。
x
は独立変数です。y
は従属変数です。lambda
は固有値を表す未知パラメーターです。
代入 および を使用して、Mathieu 方程式を 1 次系として記述できます。
,
.
対応する関数は次のようになります。
function dydx = mat4ode(x,y,lambda) % equation being solved dydx = [y(2) -(lambda - 2*q*cos(2*x))*y(1)]; end
メモ: 関数はすべて例の最後にローカル関数として含まれます。
境界条件のコード化
次に、境界点で境界条件の残差値を返す関数を記述します。この関数にはシグネチャ res = mat4bc(ya,yb,lambda)
がなければなりません。ここでは以下のようになります。
ya
は、区間 の始点における境界条件の値です。yb
は、区間 の終点における境界条件の値です。lambda
は固有値を表す未知パラメーターです。
この問題には区間 に 3 つの境界条件があります。残差値を計算するには、境界条件を の形式にしなければなりません。この形式では、境界条件は次となります。
,
,
.
対応する関数は次のようになります。
function res = mat4bc(ya,yb,lambda) % boundary conditions res = [ya(2) yb(2) ya(1)-1]; end
初期推定の作成
最後に、解の初期推定を作成します。両方の解要素 および に対する初期推定と、未知パラメーター を指定しなければなりません。計算できるのは初期推定に近い固有値および固有関数のみです。計算済みの固有関数が 4 番目の固有値に対応する可能性を高めるには、正確な定性的動作をもつ初期推定を選択する必要があります。
この問題の場合、余弦関数は 3 つの境界条件を満たすため、適切な初期推定を行います。 と の推定を返す関数を使用して、 の初期推定をコード化します。
function yinit = mat4init(x) % initial guess function yinit = [cos(4*x) -4*sin(4*x)]; end
区間 の 10 点のメッシュ、初期推定関数、 の値に対する 15 の推定を使用して bvpinit
を呼び出します。
lambda = 15; solinit = bvpinit(linspace(0,pi,10),@mat4init,lambda);
方程式の求解
ODE 関数、境界条件関数、初期推定を使用して bvp4c
を呼び出します。
sol = bvp4c(@mat4ode, @mat4bc, solinit);
パラメーターの値
bvp4c
を使って求められる、未知パラメーター の値を表示します。この値は Mathieu 方程式の 4 番目の固有値 () です。
fprintf('Fourth eigenvalue is approximately %7.3f.\n',... sol.parameters)
Fourth eigenvalue is approximately 17.097.
解のプロット
deval
を使用して、区間 の 100 点で bvp4c
によって計算される解を評価します。
xint = linspace(0,pi); Sxint = deval(sol,xint);
両方の解要素をプロットします。プロットは、4 番目の固有値 に関連付けられている固有関数 (およびその導関数) を示しています。
plot(xint,Sxint) axis([0 pi -4 4]) title('Eigenfunction of Mathieu''s Equation.') xlabel('x') ylabel('y') legend('y','y''')
ローカル関数
ここでは、BVP ソルバー bvp4c
が解を計算するために呼び出すローカル補助関数を紹介しています。あるいは、これらの関数を独自のファイルとして MATLAB パスのディレクトリに保存することもできます。
function dydx = mat4ode(x,y,lambda) % equation being solved q = 5; dydx = [y(2) -(lambda - 2*q*cos(2*x))*y(1)]; end %------------------------------------------- function res = mat4bc(ya,yb,lambda) % boundary conditions res = [ya(2) yb(2) ya(1)-1]; end %------------------------------------------- function yinit = mat4init(x) % initial guess function yinit = [cos(4*x) -4*sin(4*x)]; end %-------------------------------------------