Notebook

Modelica.Media を使って室温の変化を調べる¶

Modelica.Media ライブラリを使用して、圧力を一定に保った部屋を加熱するモデルです。部屋には通気孔があり、膨張した空気が流出します。これは 1. Modelicaのクラスの概要の ClassExample9 として作成したもので、Modelica.Mediaライブラリの BaseProperties モデルの使い方を示す例題です。

部屋の容積を $V = 22.0$ m3 (4畳半程度)、圧力を $p_{amb} = 101325$ Pa (1気圧)、初期温度を $T_{amb} = 293.15$ K、加熱量を $Q_{flow} = 100$ Wとします。

空気の密度を $\rho$ 、比内部エネルギーを $u$ とすると、室内の空気の質量 $M$ と内部エネルギー $U$ は、

$M = \rho V$ と

$U = M u$ となります。質量保存式とエネルギー保存式は、

$\frac{dM}{dt} = m_{flow}$ と

$\frac{dU}{dt} = Q_{flow} + h \cdot m_{flow}$ となります。圧力は、

$p = p_{amb}$ となります。

これらの方程式を実装したモデルを作成します。 $m_{flow}$ は質量流量で流出する場合は負になります。 $h$ は比エンタルピです。

In [1]:

%%writefile ClassExample9_1.mo
package ClassExample9_1
  import Modelica.Media;
  import SI = Modelica.SIunits;
  
  model RoomA
    package Medium = Media.Air.DryAirNasa;
    parameter SI.Temperature T_amb = 293.15;
    parameter SI.Pressure p_amb = 101325;
    parameter SI.Volume V = 22.0;
    parameter SI.HeatFlowRate Q_flow = 100;
    
    Medium.BaseProperties medium;
    SI.MassFlowRate m_flow(start = 0.0);
    SI.Mass M;
    SI.Energy U;
  equation
    M = medium.d * V;
    U = medium.u * M;
    der(M) = m_flow;
    der(U) = Q_flow + medium.h * m_flow;
    medium.p = p_amb;
  initial equation
     medium.T = T_amb;
  end RoomA;

end ClassExample9_1;

Overwriting ClassExample9_1.mo

モデルをコンパイルして FMU を作成します。

In [2]:

from pymodelica import compile_fmu
fmu1 = compile_fmu('ClassExample9_1.RoomA', 'ClassExample9_1.mo')

FMU をロードして、シミュレーションを実行します。

In [3]:

from pyfmi import load_fmu
model1 = load_fmu(fmu1)
res1 = model1.simulate(final_time = 2000)

Final Run Statistics: --- 

 Number of steps                                 : 3
 Number of function evaluations                  : 9
 Number of Jacobian evaluations                  : 1
 Number of function eval. due to Jacobian eval.  : 1
 Number of error test failures                   : 0
 Number of nonlinear iterations                  : 5
 Number of nonlinear convergence failures        : 0
 Number of state function evaluations            : 4

Solver options:

 Solver                   : CVode
 Linear multistep method  : BDF
 Nonlinear solver         : Newton
 Linear solver type       : DENSE
 Maximal order            : 5
 Tolerances (absolute)    : 0.0005
 Tolerances (relative)    : 0.0001

Simulation interval    : 0.0 - 2000.0 seconds.
Elapsed simulation time: 0.000372886657715 seconds.

シミュレーション結果の部屋の温度と流出する空気の質量流量をプロットします。流出流量を正にするためにマイナス符号を付けています。

In [4]:

%matplotlib notebook
t1 = res1['time']
T1 = res1['medium.T']
m_flow1 = res1['m_flow']

import matplotlib.pyplot as plt
plt.figure(1)
plt.subplot(2,1,1)
plt.xlim(0,2000)
plt.plot(t1,T1)
plt.legend(['Temperature [K]'], loc='lower right')

plt.subplot(2,1,2)
plt.xlim(0,2000)
plt.plot(t1,-m_flow1)
plt.legend(['Mass Flow [kg/s]'])
plt.show()

RoomA の熱や流体の出入り口をコネクタで表現したモデル RoomB を作成します。熱の出入りを HeatPort_a で表し、通気孔に流体の出入りを FluidPort_b で表します。そして、RoomB と熱源を表す FixedHeatFlow および圧力境界を表す Boundary_pT を接続したテストモデル RoomBTest を作成します。

In [5]:

%%writefile ClassExample9_2.mo
package ClassExample9_2
  import Modelica.Media;
  import SI = Modelica.SIunits;

  model RoomB
    replaceable package Medium = Media.Air.DryAirNasa;
    parameter SI.Volume V = 22.0;
    Modelica.Fluid.Interfaces.FluidPort_b port_b(redeclare package Medium = Medium);
    Modelica.Thermal.HeatTransfer.Interfaces.HeatPort_a port_a;
    Medium.BaseProperties medium;
    SI.Mass M;
    SI.Energy U;
  equation
    M = medium.d * V;
    U = medium.u * M;
    der(M) = port_b.m_flow;
    der(U) = port_a.Q_flow + actualStream(port_b.h_outflow) * port_b.m_flow;
    port_b.p = medium.p;
    port_b.h_outflow = medium.h;
    port_a.T = medium.T;
    initial equation
    medium.T = 293.15;
  end RoomB;
  
  model RoomBTest
    replaceable package Medium = Media.Air.DryAirNasa;
    RoomB roomB1(redeclare package Medium = Medium);
    Modelica.Thermal.HeatTransfer.Sources.FixedHeatFlow fixedHeatFlow1(Q_flow = 100);
    Modelica.Fluid.Sources.Boundary_pT boundary(
      redeclare package Medium = Medium,
      T = 293.15, nPorts = 1, p = 101325);
  equation
    connect(roomB1.port_b, boundary.ports[1]);
    connect(fixedHeatFlow1.port, roomB1.port_a);
  end RoomBTest;
  
end ClassExample9_2;

Overwriting ClassExample9_2.mo

In [6]:

fmu2 = compile_fmu('ClassExample9_2.RoomBTest', 'ClassExample9_2.mo')

In [7]:

from pyfmi import load_fmu
model2 = load_fmu(fmu2)
res2 = model2.simulate(final_time = 2000)

Final Run Statistics: --- 

 Number of steps                                 : 3
 Number of function evaluations                  : 9
 Number of Jacobian evaluations                  : 1
 Number of function eval. due to Jacobian eval.  : 1
 Number of error test failures                   : 0
 Number of nonlinear iterations                  : 5
 Number of nonlinear convergence failures        : 0
 Number of state function evaluations            : 4

Solver options:

 Solver                   : CVode
 Linear multistep method  : BDF
 Nonlinear solver         : Newton
 Linear solver type       : DENSE
 Maximal order            : 5
 Tolerances (absolute)    : 0.0005
 Tolerances (relative)    : 0.0001

Simulation interval    : 0.0 - 2000.0 seconds.
Elapsed simulation time: 0.000428915023804 seconds.

In [8]:

t2 = res2['time']
T2 = res2['roomB1.medium.T']
m_flow2 = res2['boundary.ports[1].m_flow']

import matplotlib.pyplot as plt
plt.figure(2)
plt.subplot(2,1,1)
plt.xlim(0,2000)
plt.plot(t2,T2)
plt.legend(['Temperature [K]'], loc='lower right')

plt.subplot(2,1,2)
plt.xlim(0,2000)
plt.plot(t2,m_flow2)
plt.legend(['Mass Flow [kg/s]'])

plt.show()

In [ ]: