ModelicaのBuildingsライブラリのThermalZones.Detailed.MixedAirで室温/熱負荷計算を学ぶ_その1
概要
OpenModelicaとJModelicaでLBNLのBuildingsライブラリを使う。初心者なのでいろいろ教えて頂けるとありがたい。
前回までは換気回路網計算を見ていたが、今回は熱回路網計算によって自然室温計算をしているモデルBuildings.ThermalZones.Detailed.Examples.MixedAirFreeResponseを見ていく。本当はBuildings.ThermalZones.Detailed.Validation.BESTEST.Cases6xx.Case600FFの方を見たかったのだけれどどうも計算が途中で止まるのでまずはこちらをみていくことにする。理解が深まったらなにか原因がわかるかもしれない。OpenModelicaでは実行できなかったのでコードを読むのはOpenModelica、実行するのはJModelicaとした。
使用バージョン
-OpenModelica1.12.0 →Modelica標準ライブラリ3.2.2
-JModelica3.4 →Modelica標準ライブラリ3.2.2
-Buildings 5.1.0
Buildings.ThermalZonesについて
ThermalZonesは建物の年間熱負荷計算とか室温変動の非定常計算などをするためのパッケージなのだと思う。室温の計算のためにBuildings.FluidやBuildings.HeatTransferなどのモデルを組み合わせてモデルを作ってくれているようだ。下の表のようにまずは大きく二つに分かれていてDetailedが詳細計算で、ReducedOrderの方はVerein Deutscher Ingenieureというドイツ技術者協会の規格のVDI 6007にある室と建物の過渡的な熱応答の計算方法を実装しているらしい。
| 名称 | 説明 |
|---|---|
| Detailed | 室モデルのパッケージ |
| ReducedOrder | VDI 6007に基づく簡易化したモデル |
どちらも気になるのだがDetailedの方をみることにする。Detailedに含まれるパッケージは下の表の通り。
| 名称 | 説明 |
|---|---|
| UsersGuide | ユーザーガイド |
| CFD | CFDにより室内空気を計算するモデル |
| MixedAir | 完全混合空気による室のモデル |
| Constructions | 室モデルで使われる壁体構成のためのモデル |
| FLEXLAB | FLEXLAB（多分LBNLの実験施設）をモデル化するために使うモデル |
| Types | 型定義のパッケージ |
| Examples | モデルの使用方法とテストモデルを示すモデル群 |
| Validation | 検証モデル群 |
| BaseClasses | Buildings.ThermalZones.Detailedのための基底クラスのパッケージ |
メインとなるモデルはCFDとMixedAirの二つで、CFDの方は室内の温度などの分布をFast Fluid Dynamics (FFD) のプログラムで計算してModelicaと連成しているらしくBuildings.ThermalZones.Detailed.Examples.FFD.UsersGuideあたりを読みながら使うのだと思う。面白そうだが 計算時間とか精度とかどんなものなのだろう。いずれ使ってみたい。
今回はMixedAirの方を使った検証モデルについて見ていくことにする。このモデルは室内空気が完全混合しているものとして一つの代表点で室温をあらわしている。よくある熱回路網計算である。ちなみにConstructionsはMixedAirの実装に使われているモデルで、FLEXLABはMixedAirを使ってLBNLの実験施設のFLEXLABのセルを再現しているものだと思われる。
MixedAirの方にはコンピュータプログラムによる建物のエネルギー性能評価の試験方法のアメリカの規格(ANSI/ASHRAE Standard 140-2007)に含まれるBESTESTによる方法に従って検証をしたファイルがある。この規格の最新版は2017なので最新の規格に対応できているかまではわからない。ASHRAEはいくつかの規格を公開しているものの140は公開していない。買うと110ドルするのでどこかで見れないものだろうか。いくつかファイルがあるのだが、Case600FFが単室の自然室温の計算をしているもので、これに空調を加えたものがCase600、そしてCase600を継承して多くのファイルが派生している。Cases6xx系統は軽量な（熱容量の小さい）壁体のモデルでCases9xx系統が重い（熱容量の大きい）壁体のモデルになっている。これにも取り組みたいのだが、まずはExamplesにあるもう少しシンプルなモデルMixedAirFreeResponseをみることにする。
MixedAirFreeResponseの構成
最初にモデル全体の構成を確認する。
- 室のモデル（MixedAir）roo
- 気象データweaDat
- 土壌温度TSoi
- 間仕切conOutを介した隣室温度TBou
- 室内の発熱multiplex3_1とそれにつないだ固定値の3種の発熱
- 窓のブラインド等の開閉制御信号replicatorとそれにつないだ固定値
これらを用いて二日間分の単室の室温を計算している。
図において線でつながっていないアイコンがあるがこれはrecord(変数のまとまり)である。
今回パラメータとして定義しているのは室を構成する素材の数と、recordによる素材の構成の定義である。
- 窓の数 parameter Integer nConExtWin = 1;
- 土壌に接する壁の数 parameter Integer nConBou = 1;
- 隣室に接する壁の数 parameter Integer nSurBou = 1;
の3つをパラメータとしており外壁の数もparameter Integer nConExt = 2;としてよいはずだがこれはなぜかrooの方で定義されている。土壌や隣室は勝手にそう呼んでいるだけで、nConBouは室rooの中でConBouで指定した壁面に接する温度、nSurBouは室rooの空気に直接接する温度をつなぐ面の数のことである。
素材の構成はライブラリに登録されたものをそのまま使うことも、素材を指定して厚みを与えて自分で構成を作ることもできる。素材の物性値を変更したい場合は自分で素材のrecordを作ってやる必要がある。
外壁parameter Buildings.HeatTransfer.Data.OpaqueConstructions.Insulation100Concrete200 matLayExt;
間仕切壁parameter Buildings.HeatTransfer.Data.OpaqueConstructions.Brick120 matLayPar;
屋根（上が外側として指定のはず）
parameter Buildings.HeatTransfer.Data.OpaqueConstructions.Generic matLayRoo(
material={
HeatTransfer.Data.Solids.InsulationBoard(x=0.2),
HeatTransfer.Data.Solids.Concrete(x=0.2)},
final nLay=2);
床（上が外側として指定のはず）
parameter Buildings.HeatTransfer.Data.OpaqueConstructions.Generic matLayFlo(
material={
HeatTransfer.Data.Solids.Concrete(x=0.2),
HeatTransfer.Data.Solids.InsulationBoard(x=0.15),
HeatTransfer.Data.Solids.Concrete(x=0.05)},
final nLay=3);
窓の素材
parameter Buildings.HeatTransfer.Data.GlazingSystems.DoubleClearAir13Clear glaSys(
UFra=2,
shade=Buildings.HeatTransfer.Data.Shades.Gray(),
haveInteriorShade=false,
haveExteriorShade=false);
以上が素材のrecordで、当然これらもアイコンをダブルクリックして変数に値を入力できる。
屋根のmatLayRooの場合以下のような感じになる。両面の長波放射率と日射吸収率はデフォルト値がはいっている。外気側のroughnessは何に使う値かわからない。対流熱伝達率でもかわるのだろうか。
モデルの頭でpackage MediumA = Buildings.Media.Airとしてrooの中で再宣言するのはいつもと同じである。
気象データweaDat
気象データはconnect(weaDat.weaBus, roo.weaBus)としてweaBusをつなぎこんでいる。ファイルを適当に指定して、湿球温度の計算は切っておく。
Buildings.BoundaryConditions.WeatherData.ReaderTMY3 weaDat(
filNam=Modelica.Utilities.Files.loadResource("modelica://Buildings/Resources/weatherdata/USA_IL_Chicago-OHare.Intl.AP.725300_TMY3.mos"),
computeWetBulbTemperature=false)
間仕切conOutと隣室温度TBou
surf_surBouという室rooの空気に直接接続するheatportの配列に多層壁のモデルの配列をconnect(roo.surf_surBou, conOut.port_a)としてつないでいる。多層壁の反対側は15℃固定とした隣室側表面温度TBouにconnect(TBou.port,conOut.port_b)としてつないでいる。多層壁モデルはstateAtSurfaceをtrueにすると表温度を計算するので、今回のように片側の温度を固定している場合はFalseにしてやる必要がある。Aに壁の面積、layersに壁体構成のrecordを入力する。port_aがrecordで上側に定義したものとつながるようになっている。
Buildings.HeatTransfer.Sources.FixedTemperature TBou[nSurBou](each T=288.15);
Buildings.HeatTransfer.Conduction.MultiLayer conOut[nSurBou](
each A=6*4,
each layers=matLayPar,
each steadyStateInitial=true,
each stateAtSurface_a=true,
each stateAtSurface_b=false);
土壌温度TSoi
室rooの壁体conBouに接続するsurf_conBouというheatportの配列に10℃の固定温度TSoiをconnect(TSoi.port, roo.surf_conBou)として与えている。
Buildings.HeatTransfer.Sources.FixedTemperature TSoi[nConBou](each T = 283.15)
室内の発熱multiplex3_1
室内の発熱として放射、対流、潜熱流の順番の3つのheatportの配列qGai_flowにmultiplex3_1から発熱をconnect(multiplex3_1.y, roo.qGai_flow)として与えている。放射だけならheaPorRad、対流だけならheaPorAirというheatportがそれぞれ単独で存在するようだ。潜熱流としてMixedAirに内部発熱を定義したいときは要素３つの配列にしてやる必要があるのだと思う。
Modelica.Blocks.Routing.Multiplex3は3つのinputを順にならべた配列にして１つのoutputにするモデルである。inputが配列でも問題ないが今回は単独の数値3つをまとめているのでoutputは要素3つの配列である。値はすべて0なので実質発熱はない。
Modelica.Blocks.Routing.Multiplex3 multiplex3_1; Modelica.Blocks.Sources.Constant qRadGai_flow(k=0); Modelica.Blocks.Sources.Constant qConGai_flow(k=0); Modelica.Blocks.Sources.Constant qLatGai_flow(k=0);
connect(qRadGai_flow.y, multiplex3_1.u1[1]); connect(qConGai_flow.y, multiplex3_1.u2[1]); connect(qLatGai_flow.y, multiplex3_1.u3[1]);
窓のブラインド等の開閉制御信号replicator
uShaという窓のブラインド等の開閉制御信号を入力する配列にReplicatorで作成した配列をconnect(roo.uSha, replicator.y)としてつないでいる。Replicatorはinputの値がnout個並ぶ配列を出力するモデルであり、複数の窓に一括で同じ制御をしたいときなどに便利だと思う。今回は定数で0をconnect(uSha.y, replicator.u) としてつなげて複製しているので日射遮蔽はないものとしていると思われる。
Modelica.Blocks.Sources.Constant uSha(k=0) Modelica.Blocks.Routing.Replicator replicator(nout=max(1,nConExtWin))
これとは別にuWinという窓の状態(電気のオンオフで透過率が変わるなど？)を表すReal型の入力ポート（たぶん日射をカットする割合を入力する）がある。uShaもuWinも0～1の値を入力することになる。窓周りは詳しく調べておきたいところだがとりあえずこんなものがあるのか程度で留めておく。
壁や窓を含む室モデルroo
このモデルが今回の本題である。設定は以下のようになっている。
Buildings.ThermalZones.Detailed.MixedAir roo(
redeclare package Medium = MediumA,
AFlo=6*4,
hRoo=2.7,
nConExt=2,
datConExt(layers={matLayRoo, matLayExt},
A={6*4, 6*3},
til={Buildings.Types.Tilt.Ceiling, Buildings.Types.Tilt.Wall},
azi={Buildings.Types.Azimuth.S, Buildings.Types.Azimuth.W}),
nConExtWin=nConExtWin,
datConExtWin(
layers={matLayExt},
each A=4*3,
glaSys={glaSys},
each hWin=2,
each wWin=4,
ove(wR={0},wL={0}, gap={0.1}, dep={1}),
each fFra=0.1,
each til=Buildings.Types.Tilt.Wall,
azi={Buildings.Types.Azimuth.S}),
nConPar=1,
datConPar(layers={matLayPar}, each A=10,
each til=Buildings.Types.Tilt.Wall),
nConBou=1,
datConBou(layers={matLayFlo}, each A=6*4,
each til=Buildings.Types.Tilt.Floor,
each stateAtSurface_a = false),
nSurBou=1,
surBou(each A=6*3,
each absIR=0.9,
each absSol=0.9,
each til=Buildings.Types.Tilt.Wall),
energyDynamics=Modelica.Fluid.Types.Dynamics.FixedInitial,
T_start=273.15+22,
lat=0.73268921998722) "Room model"
床面積、天井高さ、各構成要素の数とその設定などが入力されている。各構成要素の数の入力にパラメーターの変数を使っていないところがあるがミスだろうか。ConExtが外壁で気象データとつながる壁体、ConExtwinが窓、ConParが室内壁(壁の両面がこの室に面している)、ConBouが土壌、SurBouが今回の間仕切壁となっている。それぞれの要素を表にすると以下のような使い分けになっている。この名称の頭にnをつけて要素数、datをつけてそれぞれの設定を記述している。
| 名称 | 用途 | 境界条件を与えるポート |
|---|---|---|
| ConExt | 外壁 | weabus |
| ConExtwin | 窓 | weabus |
| ConPar | 室内の壁の蓄熱とか放射の効果を見たいとき | なし(両面とも室内空気に面する) |
| ConBou | 境界の温度を気象データ以外から与えたいとき | surf_conBou |
| SurBou | 壁のモデルをMixedAir以外でつくりたいとき | surf_surBou |
surBouは壁体がモデルの外でモデル化されるので素材の指定はないのだが、これの指定の変数名はdatSurBouでなくsurBouであることには注意が必要かも知れない。長波放射率（室内の壁同士の放射熱交換用）absIR=0.9と日射吸収率absSol=0.9（窓から入ってきた日射の計算用）を設定している。
Buildingsライブラリで使っている定数として自然換気で風向を指定するときにAzimuthについては確認した。今回は壁体の角度tilに使っている角度の定数Buildings.Types.Tiltを確認する。こちらも単位はラジアンで角度の定数が入っていて、壁面の方向が0[rad]が上向き（水平屋根面など）となっている。
constant Modelica.SIunits.Angle Ceiling = 0 ; constant Modelica.SIunits.Angle Floor = Modelica.Constants.pi; constant Modelica.SIunits.Angle Wall = Modelica.Constants.pi/2;
窓の設定は庇や袖壁が設定できるようになっている。hWinが窓高さ、wWinが窓の幅、ove(窓の右端から庇の右端の距離,窓の左端から庇の左端の距離, 窓の上面から庇まで距離,庇の深さ)が庇の設定で、このあたりのパラメータは窓のガラス部分のサイズを基準に入力するはずである。右というのが外側から窓に向かって右側であるらしい。またfFraは窓面全体面積にたいするサッシ面積の割合となっている。ちなみに袖壁はsidFin(窓の上面より上に出る長さ,窓の側面からの距離, 袖壁の深さ)袖壁は窓の最下端以下の部分からはじまり、窓の両側に対象的に配置される。
初期温度は室の空気だけでなくここで定義した壁体すべてにも適用される。緯度は気象データのものと一致している必要があるようなので太陽位置などを計算するのだろう。太陽位置は日射量を直達成分と拡散成分に分解したときと壁面の角度に応じて計算するときで一致していなければならない。
今回は使用していないが、MixedAirにはportsとC_flowはfluidportとCO2などの物質の濃度のポートがある。このポートとAirflowのモデルをつかって換気計算も同時に行うことができる。（計算負荷がどの程度になるかは不明）
OMEdit上での設定も確認してみる。
dummyConとdummyGlaSysというダミーのパラメータをインスタンスにあたえるための設定があるが使いどころがよく分からない。
linearizeRadiationは放射計算を温度の４乗で行うと計算負荷が高くなるので放射伝達率などで近似的に線形で計算する設定だと思われる。
hIntFixedは室内側対流熱伝達率でhExtFixedは外気側対流熱伝達率である。ただしデフォルトではこの値は使われておらず、下のintConModとextConModでFixedを選択することでこの値を使用することになると思われる。
steadyStateWindowは窓を定常計算モデルとするかの設定である。窓は熱容量が小さいので定常計算（熱容量なし）とするのはよくあるが、デフォルトでは熱容量を考慮しておりしかもこれがシミュレーションの高速化に寄与するらしい。
mSenFacの室の家具などの熱容量の設定がここにあるが、今回は特別見込まず空気の熱容量のみになっている。
ホモトピー法を使用するかどうかの設定。
ここは実験的に時間のサンプリングの処理により計算の高速化を図ろうとする設定らしいが将来的になくなるかもしれないらしい。
圧力と温度の初期値の設定。
計算の実行
OpenModelicaで実行しようとすると以下のエラーがでる。
[3]～[5]はこのファイルを開いたときの警告、[1]と[2]は実行時のエラーである。壁体の配列のあたりでなにかよくないことが起こっているのかもしれない。窓のペアガラスの中空層の空気の部分もファイルを開いたときにエラーとなっている。
とりあえずJModelicaで実行してみる。今回は前回JModelicaを使ったときよりバージョンがあがって2.4を使っている。ファイルが見つからないという警告がでなくなった。少し初期値がみつからないと警告してくるがとりあえず計算はできた。
室温roo.air.vol.Tと外気温度roo.weaBus.DryBul、屋根roo.conEXT[1]と間仕切りconOutの表面温度は以下の通り。計算開始時にがくっと温度が変化している。屋根表面の温度が外気より下がるのは夜間放射の影響と思われる。昼間は日射の影響により屋根や間仕切りの室内側表面温度に比べて室温が上昇している。家具の熱容量は考慮していないがRCの壁を使っているので室温変動はこんなものなのだろうか。
まとめ
とりあえずMixedAirの使い方が分かった。あいかわらず計算結果から知りたい値を探し出すのがたいへんである。伝熱の各要素などの理解を深めたい気もするが、次回はとりあえずBESTESTの計算をどうにかうまくできないかと思う。計算が止まる原因の解決が難しそうなら多数室の計算とか空調時の熱負荷の算出に挑戦してもよいかもしれない。
先日LBNLがModelicaやFMI関係のスタッフを募集していた。
BuildingsライブラリやOpen Building Control, Spawn of EnergyPlus それからIBPSA librariesにも関わる事をするようだ。Modelica関係の開発がさらに活発化するかもしれない。Job openings in LBNL’s Simulation Research Group https://t.co/5bqOxLwkKA
— Michael Wetter (@michael_wetter) 2018年10月15日
