💻

Julia, ModelingToolkit の練習 - RC回路

2024/07/15に公開

プログラミング言語 Julia の ModelingToolkit の チュートリアルをトレース。
RC回路を非因果モデリングで解く。

Tutorial Acausal Component-Based Modeling

以下のチュートリアルのトレース
https://docs.sciml.ai/ModelingToolkit/stable/tutorials/acausal_components/

RC回路のイメージ

参考

https://ja.wikipedia.org/wiki/RC回路
https://www.jstage.jst.go.jp/article/sicejl/53/4/53_278/_pdf

Juliadでのプログラミング (写経)

部品の定義

使用するパッケージの読み込み

using ModelingToolkit, Plots, DifferentialEquations
using ModelingToolkit: t_nounits as t, D_nounits as D
gr(
    lw = 2, size=(700,500), legend=:outertopright,
    xtickfont = font(10), ytickfont = font(10), markersize=8,
    framestyle=:box
)

電気回路の部品をモデリングする。
回路部品は、入出力端子を持つ。端子は、2つの値、電圧と電流を持つ。
Pinコンポーネントをコネクタとして定義する。
2つのPinが接続されると、キルヒホッフの法則を満たすように定義をする。
キルヒホッフの法則は、電流の和はゼロになり、接続された端子間の電圧は等しくなる。
connect = Flow で、電流の和はゼロ、通常の接続では、変数は等しくなる。

@connector Pin begin
    v(t)
    i(t), [connect = Flow]
end;

次に、接地ノードを定義します。
接地ノードは、定電圧リザーバー(通常は V = 0 とされる)に接続されたピンです。
このコンポーネントを定義するには、Pin サブコンポーネントを持ち、
その Pin の電圧がゼロに等しいことを指定します。

@mtkmodel Ground begin
    @components begin
        g = Pin()
    end
    @equations begin
        g.v ~ 0
    end
end;

次に、OnePort を定義します。
OnePort は、2 つのピンを持つすべての単純な電気部品の抽象化です。
正のピンと負のピンの間の電圧差が部品の電圧になる。
2 つのピン間の電流の合計はゼロでなければならず、部品の電流は正のピンの電流に等しくなります。

@mtkmodel OnePort begin
    @components begin
        p = Pin()
        n = Pin()
    end
    @variables begin
        v(t)
        i(t)
    end
    @equations begin
        v ~ p.v - n.v # 電圧差が、部品の電圧
        0 ~ p.i + n.i # 電流の和はゼロ
        i ~ p.i           # 電流は、正のピンとおなじになる。 
    end
end;

抵抗を定義します。
抵抗は、正と負の 2 つのピンを持ち、オームの法則 (v = i * r) に従うオブジェクトです。
抵抗器の電圧は、2 つのピン間の電圧差として与えられます。
一方、電荷保存則により、入力電流と出力電流は等しくなければならず、これは電流の方向に関係なく電流の合計がゼロでなければならないことを意味します。これにより、抵抗器の方程式が得られます。

OnePort を継承して、抵抗Rをデフォルト値= 1 として追加。

@mtkmodel Resistor begin
    @extend OnePort() # OnePortの継承。oneport.v などの名前空間はつくられない。
    @parameters begin
        R = 1.0 # sets the default resistance
    end
    @equations begin
        v ~ i * R # OnePortに新しい式(オームの法則)を拡張
    end
end;

コンデンサを定義します。
OnePortを継承し、コンデンサ容量を追加、デフォルト値 = 1 として定義。
コンデンサの電圧の時間変化(微分値)は、電流/コンデンサ容量 の関係を追加定義。

@mtkmodel Capacitor begin
    @extend OnePort()
    @parameters begin
        C = 1.0
    end
    @equations begin
        D(v) ~ i/C
    end
end;

次に、定電圧電源項を定義します。
これは、2 つのピンを持つオブジェクトと考えることができ、オブジェクト自体は一定電圧に保たれ、電流を生成します。
これは次のようにモデル化できます。

@mtkmodel ConstantVoltage begin
    @extend OnePort()
    @parameters begin
        V = 1.0
    end
    @equations begin
        V ~ v
    end
end;

部品の接続と回路シミュレーション

定義した部品を使って、以下のRC回路を作りシミュレーションする。

@mtkmodel RCModel begin
    @components begin
        resistor = Resistor(R = 1.0)
        capacitor = Capacitor(C = 1.0)
        source = ConstantVoltage(V = 1.0)
        ground = Ground()
    end
    @equations begin
        connect(source.p, resistor.p)
        connect(resistor.n, capacitor.p)
        connect(capacitor.n, source.n)
        connect(capacitor.n, ground.g)
    end
end;
@mtkbuild rc_model = RCModel(resistor.R = 2.0)
u0 = [rc_model.capacitor.v => 0.0
      rc_model.capacitor.p.i => 0.0]
prob = ODEProblem(rc_model, u0, (0, 10.0))
sol = solve(prob)
plot(sol)

コンデンサーの電圧の変化をシミュレート、グラフ化した。

抵抗にかかる電圧は、以下になる。

plot(sol, idxs = [rc_model.resistor.v])

Discussion