pdeCoefficientsToDouble

シンボリック PDE 係数の double 形式への変換

R2021a 以降

構文

coeffs = pdeCoefficientsToDouble(symCoeffs)

coeffs = pdeCoefficientsToDouble(symCoeffs,u)

説明

coeffs = pdeCoefficientsToDouble(symCoeffs) は、構造体 symCoeffs のシンボリック オブジェクトを倍精度数または関数ハンドルに変換します。出力は構造体 coeffs であり、Partial Differential Equation Toolbox™ で specifyCoefficients を呼び出すことにより、PDE モデルの係数を定義するために使用できます。

構造体 coeffs には、次の形式の PDE 系の係数 m、d、c、a、および f が含まれています。

$$m\frac{{\partial }^{2}u}{\partial t^2}+d\frac{\partial u}{\partial t}-\nabla ·\left(c\otimes \nabla u\right)+au=f$$

これは Partial Differential Equation Toolbox で解くことができます。詳細については、Equations You Can Solve Using Partial Differential Equation Toolbox (Partial Differential Equation Toolbox)を参照してください。

coeffs = pdeCoefficientsToDouble(symCoeffs,u) は、構造体 symCoeffs のシンボリック オブジェクトを、従属変数 u をもつ PDE 系の係数用に、倍精度数または関数ハンドルに変換します。

例

シンボリック PDE からの係数の抽出

変数 u(x,y) のラプラシアンを表すシンボリック PDE を作成します。

syms u(x,y) f
pdeeq = laplacian(u,[x y]) == f

pdeeq(x, y) = 
$$\frac{{\partial }^2}{\partial x^2}\mathrm{ }u\left(x,y\right)+\frac{{\partial }^2}{\partial y^2}\mathrm{ }u\left(x,y\right)=f$$

PDE の係数をシンボリック式として抽出し、それらの値を表示します。

symCoeffs = pdeCoefficients(pdeeq,u,'Symbolic',true);
structfun(@disp,symCoeffs)

$$0$$

$$0$$

$$\left(\begin{array}{cc}-1& 0\\ 0& -1\end{array}\right)$$

$$f$$

$$0$$

pdeCoefficients は、シンボリック PDE を以下の形式のスカラー PDE 方程式に変換します。

$$m\frac{{\partial }^{2}u}{\partial t^2}+d\frac{\partial u}{\partial t}-\nabla \cdot (c\nabla u)+au=f$$,

また、係数 a、c、m、d、および f を抽出し、構造体 symCoeffs に格納します。

f の値を選択します。symCoeffs はシンボリック オブジェクトであるため、pdeCoefficientsToDouble を使用して係数を double データ型に変換します。double データ型の係数は、Partial Differential Equation Toolbox の specifyCoefficients 関数の有効な入力です。

symCoeffs = subs(symCoeffs,f,-3);
coeffs = pdeCoefficientsToDouble(symCoeffs)

coeffs = struct with fields:
    a: 0
    c: [4×1 double]
    m: 0
    d: 0
    f: -3

複数の PDE からなる方程式

2 つの 2 階 PDE からなる方程式を解きます。この PDE 系を解くには、pdeCoefficients を使用して PDE の係数をシンボリックに抽出し、pdeCoefficientsToDouble を使用してその係数を倍精度数に変換し、specifyCoefficients を使用して PDE モデルの係数を指定します。

この PDE 系は、固定された構造板に一様な圧力負荷をかけたときの構造板のたわみを表します。従属変数 $$u_1$$ および $$u_2$$ をもつ PDE 系は、次の式で与えられます。

$$-{\nabla }^2{u}_1+u_2=0$$,

$$-D{\nabla }^2{u}_2=p$$,

ここで、$$D$$ は板のたわみ剛性で、次の式で与えられます。

$$D=\frac{E{h}^3}{12(1-{\nu }^2)}$$,

また、$$E$$ は弾性率、$$\nu$$ はポアソン比、$$h$$ は板の厚さ、$$u_1$$ は板の横方向のたわみ、$$p$$ は圧力負荷です。

2 つの連立方程式の PDE モデルを作成します。

model = createpde(2);

矩形形状を作成します。矩形の側面の長さを指定します。次に、この幾何形状を PDE モデルに組み込みます。

len = 10.0;
gdm = [3 4 0 len len 0 0 0 len len]';
g = decsg(gdm,'S1',('S1')');
geometryFromEdges(model,g);

方程式の物理パラメーターの値を指定します。外部の圧力 $$p$$ は、任意の値を取ることができるシンボリック変数 pres とします。

E = 1.0e6;
h_thick = 0.1;
nu = 0.3;
D = E*h_thick^3/(12*(1 - nu^2));
syms pres

PDE 系をシンボリック方程式として宣言します。PDE の係数を抽出し、その係数をシンボリック形式で返します。

syms u1(x,y) u2(x,y)
pdeeq = [-laplacian(u1) + u2; -D*laplacian(u2) - pres];
symCoeffs = pdeCoefficients(pdeeq,[u1 u2],'Symbolic',true)

symCoeffs = struct with fields:
    m: 0
    a: [2×2 sym]
    c: [4×4 sym]
    f: [2×1 sym]
    d: 0

係数 m、a、c、f、および d を表示します。

structfun(@disp,symCoeffs)

$$0$$

$$\left(\begin{array}{cc}0& 1\\ 0& 0\end{array}\right)$$

$$\left(\begin{array}{cccc}1& 0& 0& 0\\ 0& 1& 0& 0\\ 0& 0& \frac{25000}{273}& 0\\ 0& 0& 0& \frac{25000}{273}\end{array}\right)$$

$$\left(\begin{array}{c}0\\ \mathrm{pres}\end{array}\right)$$

$$0$$

関数 subs を使用して pres の値を代入します。subs の出力はシンボリック オブジェクトであるため、Partial Differential Equation Toolbox への有効な入力となるように、係数を pdeCoefficientsToDouble 関数を使用して double データ型に変換します。

symCoeffs = subs(symCoeffs,pres,2);
coeffs = pdeCoefficientsToDouble(symCoeffs)

coeffs = struct with fields:
    a: [4×1 double]
    c: [16×1 double]
    m: 0
    d: 0
    f: [2×1 double]

PDE モデルの PDE 係数を指定します。

specifyCoefficients(model,'m',coeffs.m,'d',coeffs.d, ...
    'c',coeffs.c,'a',coeffs.a,'f',coeffs.f);

バネの剛性を指定します。4 つの辺すべてで分布するバネを定義し、境界条件を指定します。

k = 1e7;
bOuter = applyBoundaryCondition(model,'neumann','Edge',(1:4), ...
    'g',[0 0],'q',[0 0; k 0]);

幾何形状のメッシュ サイズを指定し、PDE モデルのメッシュを生成します。

hmax = len/20;
generateMesh(model,'Hmax',hmax);

モデルを解きます。

res = solvepde(model);

ノード位置における解にアクセスします。

u = res.NodalSolution;

板の横方向のたわみをプロットします。

numNodes = size(model.Mesh.Nodes,2);
figure;
pdeplot(model,'XYData',u(1:numNodes),'contour','on')
title 'Transverse Deflection'

板の中央における横方向のたわみを求めます。

wMax = min(u(1:numNodes,1))

wMax = 
-0.2763

この結果を、解析的に計算した、板の中央における横方向のたわみと比較します。

pres = 2;
wMax = -.0138*pres*len^4/(E*h_thick^3)

wMax = 
-0.2760

PDE 系の従属変数の指定

R2023a 以降

4 つの従属変数 $\mathit{A}$、$\mathit{E}$、$\mathit{P}$、および $\mathit{T}$ をもつ 4 つの方程式からなる PDE 系を作成します。

syms A(x,y,z,t) E(x,y,z,t) P(x,y,z,t) T(x,y,z,t)
eqn1 = diff(A,t) == -A*E*(-exp(T));
eqn2 = diff(P,t) == -P*E*(-exp(T));
eqn3 = diff(E,t) == -A*E*(-exp(T)) - P*E*(-exp(T));
eqn4 = diff(T,t) == laplacian(T,[x y z]) + A*E*(-exp(T)) + P*E*(-exp(T));
pdeeq = [eqn1; eqn2; eqn3; eqn4]

pdeeq(x, y, z, t) = 
$$\begin{array}{l}\left(\begin{array}{c}\frac{\partial }{\partial t}\mathrm{ }A\left(x,y,z,t\right)={\sigma }_2\\ \frac{\partial }{\partial t}\mathrm{ }P\left(x,y,z,t\right)={\sigma }_1\\ \frac{\partial }{\partial t}\mathrm{ }\text{E}\left(x,y,z,t\right)={\sigma }_2+{\sigma }_1\\ \frac{\partial }{\partial t}\mathrm{ }T\left(x,y,z,t\right)=\frac{{\partial }^2}{\partial x^2}\mathrm{ }T\left(x,y,z,t\right)+\frac{{\partial }^2}{\partial y^2}\mathrm{ }T\left(x,y,z,t\right)+\frac{{\partial }^2}{\partial z^2}\mathrm{ }T\left(x,y,z,t\right)-{\sigma }_2-{\sigma }_1\end{array}\right)\\ \\ \mathrm{where}\\ \\ \mathrm{ }{\sigma }_1={\mathrm{e}}^{T\left(x,y,z,t\right)}\hspace{0.17em}\text{E}\left(x,y,z,t\right)\hspace{0.17em}P\left(x,y,z,t\right)\\ \\ \mathrm{ }{\sigma }_2={\mathrm{e}}^{T\left(x,y,z,t\right)}\hspace{0.17em}A\left(x,y,z,t\right)\hspace{0.17em}\text{E}\left(x,y,z,t\right)\end{array}$$

PDE 系の係数をシンボリック式として抽出します。従属変数を表す変数 u を指定します。

u = [A E P T];
symCoeffs = pdeCoefficients(pdeeq,u,'Symbolic',true)

symCoeffs = struct with fields:
    m: [4×4 sym]
    a: [4×4 sym]
    c: [12×12 sym]
    f: [4×1 sym]
    d: [4×4 sym]

係数 m、a、c、f、および d を表示します。

structfun(@disp,symCoeffs)

$$\left(\begin{array}{cccc}0& 0& 0& 0\\ 0& 0& 0& 0\\ 0& 0& 0& 0\\ 0& 0& 0& 0\end{array}\right)$$

$$\begin{array}{l}\left(\begin{array}{cccc}-{\sigma }_2& 0& 0& 0\\ 0& {\sigma }_1& 0& 0\\ -{\sigma }_2& {\sigma }_1& 0& 0\\ {\sigma }_2& {\mathrm{e}}^{T\left(x,y,z,t\right)}\hspace{0.17em}P\left(x,y,z,t\right)& 0& 0\end{array}\right)\\ \\ \mathrm{where}\\ \\ \mathrm{ }{\sigma }_1=-{\mathrm{e}}^{T\left(x,y,z,t\right)}\hspace{0.17em}P\left(x,y,z,t\right)\\ \\ \mathrm{ }{\sigma }_2={\mathrm{e}}^{T\left(x,y,z,t\right)}\hspace{0.17em}\text{E}\left(x,y,z,t\right)\end{array}$$

$$\left(\begin{array}{cccccccccccc}0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0& 0& 0& 0& 1& 0& 0\\ 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 1& 0\\ 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 0& 1\end{array}\right)$$

$$\left(\begin{array}{c}0\\ 0\\ 0\\ 0\end{array}\right)$$

$$\left(\begin{array}{cccc}1& 0& 0& 0\\ 0& 0& 1& 0\\ 0& 1& 0& 0\\ 0& 0& 0& 1\end{array}\right)$$

Partial Differential Equation Toolbox の関数 specifyCoefficients に対して有効な入力となるように、係数を double データ型に変換します。symCoeffs の係数 a は従属変数を含んでおり、定数ではないため、2 番目の入力引数なしで pdeCoefficientsToDouble(symCoeffs) を呼び出そうとすると、エラーが返される可能性があります。代わりに、pdeCoefficientsToDouble を使用するときに、従属変数 u として 2 番目の引数を指定します。

coeffs = pdeCoefficientsToDouble(symCoeffs,u)

coeffs = struct with fields:
    a: @makeCoefficient/coefficientFunction
    c: [144×1 double]
    m: 0
    d: [16×1 double]
    f: 0

入力引数

symCoeffs — シンボリック形式の PDE 係数

シンボリック式の構造体

シンボリック形式の PDE 係数。シンボリック式の構造体として指定します。

u — PDE の従属変数

シンボリック関数

R2023a 以降

PDE の従属変数。シンボリック関数として指定します。引数 u には、2 次元または 3 次元の定常変数または時間依存変数が含まれていなければなりません。

例: syms E(x,y,z,t) B(x,y,z,t); u = [E B];

出力引数

coeffs — double で演算を行う PDE の係数

倍精度数と関数ハンドルの構造体

double で演算を行う PDE の係数。関数 specifyCoefficients が必要とする倍精度数と関数ハンドルの構造体として返されます。この構造体のフィールドは、a、c、m、d、および f です。Partial Differential Equation Toolbox が必要とする形式で係数を解釈する方法については、以下を参照してください。

• c Coefficient for specifyCoefficients (Partial Differential Equation Toolbox)

• m, d, or a Coefficient for specifyCoefficients (Partial Differential Equation Toolbox)

• f Coefficient for specifyCoefficients (Partial Differential Equation Toolbox)

pdeCoefficientsToDouble が関数ハンドルとして係数を返す際、関数ハンドルは入力引数として location と state の 2 つの構造体を取り、倍精度出力を返します。関数ハンドルは、コマンド ウィンドウで @makeCoefficient/coefficientFunction として表示されます。coeffs.f の関数ハンドルの式をシンボリック形式で表示するには、coeffs.f('show') を使用します。

生成されたすべての係数を specifyCoefficients で使用できない場合もあります。たとえば、m が非ゼロの場合、係数 d は特別な行列形式を取らなければなりません。詳細については、d Coefficient When m Is Nonzero (Partial Differential Equation Toolbox)を参照してください。

バージョン履歴

R2021a で導入

R2023a: PDE 係数を変換する際の従属変数の指定

PDE 系のシンボリック係数を倍精度数または関数ハンドルに変換する際、PDE 系の従属変数 u を pdeCoefficientsToDouble の 2 番目の入力引数として指定することができます。例については、PDE 系の従属変数の指定を参照してください。

参考

syms | diff | laplacian | pdeCoefficients | specifyCoefficients (Partial Differential Equation Toolbox)

トピック

• 非線形熱伝導に関する偏微分方程式の求解