システムとモデリング

SysML, Matlab/Simulink,他モデリング全般について。たまにプロジェクトマネジメント。

Modelicaプログラミング練習

こんばんは!Otepipiです! 今回はプログラミング言語Modelicaを使用してプログラミングの練習をしていきたいと思います。

Modelicaとは

Modelicaはオブジェクト指向プログラミング言語ですが、その中でも物理モデリングに特化しています。 言語仕様は下記にあります。

https://www.modelica.org/documents/ModelicaSpec34.pdf

大きな特長としては非因果モデルの記述が可能という点です。 非因果モデリングは、「出力」「入力」の区別をしない記述方法で、例えばJuliaやPython等で「3x+5=14+6α」といった式からxを求めたい場合

x = 3 + 2*α

とxについて式変形して求めることになりますが、Modelica言語ですと

3x + 5 = 14 + 6a

上記の式のままでxを求めることができます。

プログラミング練習

下記参考書から微分方程式の例題を説いてみます。

微分方程式で数学モデルを作ろう

微分方程式で数学モデルを作ろう

捕食ー被捕食者モデルを解いてみますが、これは以前に以下の記事でJulia言語で解いています。 otepipi.hatenablog.com

例題

$$捕食者-非捕食者の相互作用の単純なモデルは下記で与えられる。$$ $$\frac{dx}{dt}=x+y$$

$$\frac{dy}{dt}=-x+y$$

$$x………捕食者の個体数$$ $$y………非捕食者の個体数$$

$$初期値を下記のように仮定する。$$ $$x(0)=y(0)=1000$$

$$未来のxとyを決定せよ。$$

Modelicaでの実装

model TEST
  Real x (start = 1000);
  Real y (start = 1000);
 equation
 der (x) = x + y ;
 der(y) = -x + y ;
end TEST ;

さらにこのプログラムをコンパイルして実行するには下記の命令が必要です。

simulate(TEST, stopTime=1)

最後に結果をプロットします。

plot({x,y})

結果は下図のようになります。 f:id:Otepipi:20180221194825p:plain

参考:Juliaでの実装

参考に、Julia言語での実装を下に示します。

using ODE
using Plots
function f(t, v)
     (x, y) = v
     dx_dt=x + y
     dy_dt=-x + y
    
    [dx_dt ; dy_dt]
    end;

v_in = [1000.0; 1000.0];
time = 0.0:0.1:1;
x = map(v -> v[1], v);
y = map(v -> v[2], v);
plot(t,y,seriestype=:scatter,label="y")
plot!(t,x,seriestype=:scatter,label="x")

この時のプロットは下記のようになり、Modelicaの結果とよく一致しています。 f:id:Otepipi:20180221195501p:plain

今回は非因果モデリングの実例を示すことができなかったので 次回はModelicaを使った非因果モデリングについて解説したいと思います。