(57)【要約】調整環境内で気流を生成する空調システムの動作を制御するシステムは、調整環境内の点のセットにおいて気流の測定値をもたらすセンサのセットと、気流の以前の状態を気流の現在の状態に遷移させる第１項と遷移を支援する第２項との組み合わせを含む気流動態のモデルを記憶するメモリと、気流の以前の状態を、気流動態のモデルに従って時間的に前方に遷移させて、気流の現在の状態における観測誤差を低減することによって、調整環境内の気流の現在の状態を推定する推定器と、気流の現在の状態に基づいて空調システムを制御するコントローラとを備える。第１項は、気流の物理モデルの有限次元空間上への投影を含む。気流の物理モデルは、調整環境の物理パラメータを含み、投影は、第１項において調整環境の物理パラメータを維持する。第２項は、物理パラメータの有界不確実性の範囲、負の利得、及び、点のセットにおける気流の測定値と気流の物理モデルによる点のセットにおける気流の推定値との間の観測誤差の関数を含む。
2021508120.tif 000002
【特許請求の範囲】
【請求項１】
調整環境内で気流を生成する空調システムの動作を制御するシステムであって、
前記調整環境内の点のセットにおいて前記気流の測定値をもたらすセンサのセットと、
前記気流の以前の状態を前記気流の現在の状態に遷移させる第１項と、前記遷移を支援する第２項との組み合わせを含む気流動態のモデルを記憶するメモリであって、前記第１項は、前記気流の物理モデルの有限次元空間上への投影を含み、前記気流の物理モデルは、前記調整環境の物理パラメータを含み、前記投影は、前記第１項において前記調整環境の前記物理パラメータを維持し、少なくとも１つの物理パラメータの値は、有界不確実性の範囲を規定する前記有界不確実性を有して不確実であり、前記第２項は、前記有界不確実性の範囲、負の利得、及び、前記点のセットにおける前記気流の測定値と前記気流の物理モデルによる前記点のセットにおける前記気流の推定値との間の観測誤差の関数を含む、メモリと、
前記気流の以前の状態を、気流動態のモデルに従って時間的に前方に遷移させて、前記気流の現在の状態における前記観測誤差を低減することによって、前記調整環境内の前記気流の現在の状態を推定する推定器と、
前記気流の現在の状態に基づいて前記空調システムを制御するコントローラと、
を備える、システム。
【請求項２】
前記推定器は、前記調整環境内の前記気流の状態を、複数の制御ステップについて反復的に更新する、請求項１に記載のシステム。
【請求項３】
前記推定器及び前記コントローラは、プロセッサを用いて実施される、請求項１に記載のシステム。
【請求項４】
前記第２項は、前記物理パラメータの有界不確実性の前記範囲を、前記調整環境における前記観測誤差のエネルギー関数の導関数に関連付ける、請求項１に記載のシステム。
【請求項５】
前記第２項は、前記エネルギー関数の前記導関数を強制的に負定値にする、請求項４に記載のシステム。
【請求項６】
前記エネルギー関数は、前記観測誤差の二乗の積分のリアプノフ関数である、請求項５に記載のシステム。
【請求項７】
前記気流動態のモデルの前記第１項は、前記調整環境内の前記気流の時空間連続値の偏微分方程式（ＰＤＥ）の前記物理パラメータを有する常微分方程式（ＯＤＥ）を含み、これにより、前記第２項は、前記気流の物理モデルにおける前記ＯＤＥを変更する、請求項１に記載のシステム。
【請求項８】
前記物理パラメータは、前記調整環境内の空気の粘度、前記調整環境の幾何学的形状、前記調整環境内の物体の数及びタイプのうちの１つ又は組み合わせを含む、請求項７に記載のシステム。
【請求項９】
前記ＰＤＥは、ブシネスク方程式によって定義され、前記ＯＤＥの次元は、前記ブシネスク方程式に従った前記ＰＤＥの次元未満である、請求項７に記載のシステム。
【請求項１０】
前記ＰＤＥと前記ＯＤＥとを結び付ける投影演算子は、前記ＯＤＥについての基底関数の関数である、請求項７に記載のシステム。
【請求項１１】
前記基底関数は、固有直交分解（ＰＯＤ）基底関数及び動的モード分解（ＤＭＤ）基底関数のうちの一方又は組み合わせを含む、請求項１０に記載のシステム。
【請求項１２】
気流動態のモデルの利得を反復的に更新して前記観測誤差を低減するように構成されたプロセッサを更に備える、請求項１に記載のシステム。
【請求項１３】
前記プロセッサは、前記利得を求めるのにデータ駆動型最適化を用いる、請求項１２に記載のシステム。
【請求項１４】
前記データ駆動型最適化は、極値探索最適化、ガウス過程信頼上限最適化、及び疑似焼きなまし最適化のうちの１つ又は組み合わせである、請求項１３に記載のシステム。
【請求項１５】
前記コントローラは、
熱交換器を通して空気流量を調整するファンを有する蒸発器と、
前記熱交換器を通して前記空気流量を調整するファンを有する凝縮器と、
前記制御するシステムを通して冷媒を圧縮及び圧送するための速度を有する圧縮器と、
蒸気圧縮システムの高圧部分と低圧部分との間で調整可能な圧力降下を提供する膨張弁と、
を備える前記空調システムの少なくとも１つのコンポーネントに対する制御入力を決定する、請求項１に記載のシステム。
【請求項１６】
調整環境内で気流を生成する空調システムの動作を制御する方法であって、前記方法は、前記方法を実施する記憶された命令に結合されるプロセッサを用い、前記命令は、前記プロセッサによって実行されると、前記方法のステップを実行し、前記ステップは、
前記調整環境内の点のセットにおいて前記気流を測定するステップと、
前記気流の以前の状態を前記気流の現在の状態に遷移させる第１項と、前記遷移を支援する第２項との組み合わせを含む気流動態のモデルを索出するステップであって、前記第１項は、前記気流の物理モデルの有限次元空間上への投影を含み、前記気流の物理モデルは、前記調整環境の物理パラメータを含み、前記投影は、前記第１項において前記調整環境の前記物理パラメータを維持し、少なくとも１つの物理パラメータの値は、有界不確実性の範囲を規定する前記有界不確実性を有して不確実であり、前記第２項は、前記有界不確実性の前記範囲、負の利得、及び、前記点のセットにおける前記気流の測定値と前記気流の物理モデルによる前記点のセットにおける前記気流の推定値との間の観測誤差との関数を含む、ステップと、
前記気流の以前の状態を、気流動態のモデルに従って時間的に前方に遷移させて、前記気流の現在の状態における前記観測誤差を低減することによって、前記調整環境内の前記気流の現在の状態を推定するステップと、
前記気流の現在の状態に基づいて前記空調システムを制御するステップと、
を含む、方法。
【請求項１７】
前記第２項は、前記物理パラメータの有界不確実性の前記範囲を、前記調整環境における前記観測誤差のエネルギー関数の導関数に関連付け、前記第２項は、前記エネルギー関数の前記導関数を強制的に負定値にし、前記エネルギー関数は、前記観測誤差の二乗の積分のリアプノフ関数である、請求項１６に記載の方法。
【請求項１８】
前記気流動態のモデルの前記第１項は、前記調整環境内の前記気流の時空間連続値の偏微分方程式（ＰＤＥ）の前記物理パラメータを有する常微分方程式（ＯＤＥ）を含み、これにより、前記第２項は、前記気流の物理モデルにおける前記ＯＤＥを変更し、前記ＰＤＥは、ブシネスク方程式によって定義され、前記ＯＤＥの次元は、前記ブシネスク方程式に従った前記ＰＤＥの次元未満である、請求項１６に記載の方法。
【請求項１９】
前記物理パラメータは、前記調整環境内の空気の粘度、前記調整環境の幾何学的形状、前記調整環境内の物体の数及びタイプのうちの１つ又は組み合わせを含む、請求項１６に記載の方法。
【請求項２０】
方法を実行するためにプロセッサによって実行可能なプログラムが具現化された非一時的コンピュータ可読記憶媒体であって、前記方法は、
調整環境内の点のセットにおいて気流を測定することと、
前記気流の以前の状態を前記気流の現在の状態に遷移させる第１項と、前記遷移を支援する第２項との組み合わせを含む気流動態のモデルを索出することであって、前記第１項は、前記気流の物理モデルの有限次元空間上への投影を含み、前記気流の物理モデルは、前記調整環境の物理パラメータを含み、前記投影は、前記第１項において前記調整環境の前記物理パラメータを維持し、少なくとも１つの物理パラメータの値は、有界不確実性の範囲を規定する前記有界不確実性を有して不確実であり、前記第２項は、前記有界不確実性の前記範囲、負の利得、及び、前記点のセットにおける前記気流の測定値と前記気流の物理モデルによる前記点のセットにおける前記気流の推定値との間の観測誤差の関数を含む、索出することと、
前記気流の以前の状態を、気流動態のモデルに従って時間的に前方に遷移させて、前記気流の現在の状態における前記観測誤差を低減することによって、前記調整環境内の前記気流の現在の状態を推定することと、
前記気流の現在の状態に基づいて空調システムを制御することと、
を含む、非一時的コンピュータ可読記憶媒体。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は包括的には空調システムに関し、より詳細には、空調システムの動作を制御することに関する。
【背景技術】
【０００２】
暖房、換気及び空調（ＨＶＡＣ）システム等の空調システムは、産業用及び住宅用として広く使用されている。例えば、空調システムは、特定の熱量を送達する蒸気圧縮サイクルのために、変速圧縮器、可変位置弁、及び変速ファンのうちの１つ又は組み合わせを含むことができる。特定の熱量を送達する蒸気圧縮システムへのコマンド入力は多くの場合に一通りではなく、コンポーネントの種々の組み合わせによって、消費するエネルギー量が異なる。それゆえ、エネルギーを最小化するとともに、効率を最大化する入力の組み合わせを用いて蒸気圧縮システムを動作させることが望ましい。
【０００３】
これまで、エネルギー効率を最大化する方法は、空調システムの物理的過程の数学モデルを使用することに依存してきた。それらのモデルに基づく方法は、蒸気圧縮システムのコンポーネントのコマンド入力がシステムの熱力学的挙動及び消費されるエネルギーに及ぼす影響を記述しようと試みている。それらの方法では、モデルを用いて、熱負荷要件を満たし、かつエネルギーを最小化する入力の組み合わせを予測してきた。
【０００４】
ＨＶＡＣシステムの動作は、調整環境の或るエリアから別のエリアへの空気の移動を規定する調整環境内の気流を変化させる。しかしながら、気流の物理モデルは、無限次元であり、リアルタイム制御応用において用いるには複雑すぎるものである。加えて、気流のモデルは、空調システムの動作中に変化する場合もある。例えば、特許文献１を参照されたい。
【先行技術文献】
【特許文献】
【０００５】
【特許文献１】米国特許出願公開第２０１６／０２５８６４４号
【発明の概要】
【発明が解決しようとする課題】
【０００６】
そのため、システムの動作中に生成された気流のリアルタイムの知識を用いて、空調システムを制御するシステム及び方法が必要とされている。
【課題を解決するための手段】
【０００７】
いくつかの実施形態の目的は、動作の熱負荷要件が満たされ、システムの性能が改善されるような、空調システムの動作を制御するシステム及び方法を提供することである。１つの実施形態の更なる目的は、空調システムを制御する際に気流のモデルを使用するようなシステム及び方法を提供することである。別の実施形態の更なる目的は、空調システムの動作中に空調システムの性能を経時的に改善するようなシステム及び方法を提供することである。
【０００８】
いくつかの実施形態は、空調システムは経時的に変動する、という理解に基づいている。蒸気圧縮システムの動作を或る時点で正確に記述するモデルが、例えば、徐々に冷媒が漏れること、又は熱交換器への腐食の蓄積が原因でシステムが変化するので、その後の時点で正確ではなくなる場合がある。さらに、空調システムのモデルは、部屋の大きさ等の設備特有の特性を無視することが多く、それによって、空調システムのモデルは、システムの動作から逸脱する。
【０００９】
いくつかの実施形態は、空調システムの動作中に生成される気流の知識を用いて、システムの動作を最適化することができる、という認識に基づいている。しかしながら、気流のモデルは、空調システムの動作中に変化する可能性もある。したがって、システムの動作中に生成される気流のリアルタイム知識を用いて空調システムを制御するシステム及び方法が必要とされている。
【００１０】
いくつかの実施形態は、モデルが気流動態を正確に表し、調整環境内の変化に適応できるように、気流動態のモデルを決定し、更新する問題に対処する。残念ながら、気流の物理モデル、すなわち気流動態の物理特性を捕捉するモデルは、無限次元であり、リアルタイム制御応用において用いるには複雑すぎる。したがって、いくつかの実施形態は、リアルタイムの更新及び制御に適した低次元のモデルを使用する。例えば、いくつかの実施形態において、モデルは、ブシネスク方程式による気流の物理モデルより少ない数のパラメータを有する次数低減モデルを含み、ブシネスク方程式は偏微分方程式（ＰＤＥ）であり、次数低減モデルは常微分方程式（ＯＤＥ）である。
【００１１】
そのために、いくつかの実施形態は、例えば、無限次元の複雑な物理モデルを、より低い次元、例えば、有限次元のより単純なモデルに縮約することを目指すモデル縮約問題に対処する。これに関連して、モデル縮約問題は、複雑なモデルによるシステムの性能測定値と、より単純なモデルによるシステムの性能測定値との間の誤差が最小化されるように、より単純なモデルの全てのパラメータ又は係数を決定することを意味する。例えば、気流モデル化推定及び制御の用途において、気流を表すＰＤＥモデルを、より次元が低く、より複雑でないＯＤＥモデルに変換する必要がある。そのような縮約によって、空調システム等の気流システムのリアルタイム制御が可能になる。
【００１２】
いくつかの実施形態は、ＰＤＥモデルは、２つのタイプのパラメータ、すなわち、内部パラメータ及び外部パラメータを有する、という理解に基づいている。内部パラメータとは、モデルの外部では何も意味を持たないモデルの内部利得を指す。外部パラメータとは、気流に影響を与える世界の物理パラメータであり、気流動態のモデルからは独立して存在している。このような物理パラメータの例としては、調整環境内の空気の粘度、調整環境の幾何学的形状、調整環境内の物体、例えば人物等の熱源の数及びタイプのうちの１つ又は組み合わせが挙げられる。
【００１３】
いくつかの実施形態は、モデル複雑度の低減、例えば、ＰＤＥモデルのＯＤＥモデルへの変換は、複雑度低減モデルにおいてＰＤＥモデルの物理パラメータを維持するべきである、という認識に基づいている。このようにして、複雑度低減モデルは、ＰＤＥモデルをより正確に記述することができる。しかしながら、ＰＤＥモデルの物理パラメータは、多くの場合に曖昧であり、すなわち不確実性を含む。例えば、空気の粘度は変動する可能性があり、及び／又は、異なる数の人物が、異なる時点に調整環境内に存在する可能性がある。いくつかの実施形態は、物理パラメータの値が、本明細書において有界不確実性（bounded uncertainty）の範囲と呼ばれる、有界不確実性を有して不確実である、という認識に基づいている。しかしながら、有界不確実性の範囲内で、物理パラメータの値は未知であり、変動する可能性がある。
【００１４】
いくつかの実施形態は、調整環境における観測誤差のエネルギー関数に、有界不確実性の範囲内の物理パラメータの任意の値を減少させることによって、気流動態のモデルの物理パラメータについての知識の欠如を補償することができる、という認識に基づいている。そのため、いくつかの実施形態は、エネルギー関数及び時間に対するその導関数を解析するのにリアプノフ（Lyapunov）手法を用いる。いくつかの実施形態は、このような手法が、気流の状態推定に対する物理パラメータの不確実性の影響を緩和することができる、という理解に基づいている。
【００１５】
リアプノフ関数の導関数の解析中に、いくつかの実施形態は、不確実性の実際の値に関係なく、リアプノフ関数の導関数を常に負にする項を発見した。この項は、不確実性の実際の値に関係なく、リアプノフ関数の導関数を常に負にするために、有界不確実性の範囲及び負の利得の積を含み、これによって、さらに、不確実性の値に関係なく、エネルギー関数を時間反復にわたって減少させ、ひいては、有界不確実性の範囲内の物理パラメータの任意の値について、観測誤差を時間反復にわたってゼロまで減少させる。
【００１６】
いくつかの実施形態は、物理パラメータの有界不確実性の範囲をエネルギー関数の導関数に関連付けて、気流の状態の以前の状態から現在の状態への遷移の収束を支援する項によって、気流動態のモデルを変更する。そのために、いくつかの実施形態では、気流動態のモデルは、気流の以前の状態をその現在の状態に遷移させる第１項と、その遷移を支援する第２項とを含む。具体的には、第１項は、気流の物理モデルの有限次元空間上への投影を含み、気流の物理モデルは、調整環境の物理パラメータを含み、投影は、第１項において調整環境の物理パラメータを維持する。第２項は、有界不確実性の範囲、負の利得、及び、点のセットにおける気流の測定値と気流のモデルによる点のセットにおける気流の推定値との間の観測誤差の関数を含む。
【００１７】
このようにして、気流動態のモデルにおいて第１項と第２項とを組み合わせることによって、その精度を維持しながらモデルの複雑度を低減することが可能になる。種々の実施形態において、気流の状態は反復的に特定される。第２項は、観測誤差のエネルギー関数の導関数を強制的に負定値にするように設計されるので、気流の状態の反復更新は、複数の制御ステップにわたって観測誤差を強制的に低減させる。
【００１８】
定義
実施形態を説明する際に、全体を通じて（上記を含む）以下の定義が適用可能である。
【００１９】
「制御システム」又は「コントローラ」は、他のデバイス又はシステムの挙動を管理、命令、誘導又は規制するデバイス又はデバイスのセットを指す。制御システムはソフトウェア又はハードウェアのいずれかによって実現することができ、１つ又はいくつかのモジュールを含むことができる。フィードバックループを含む、制御システムは、マイクロプロセッサを用いて実現することができる。制御システムは、埋込システムとすることができる。
【００２０】
「空調システム」、又は加熱、換気及び空調（ＨＶＡＣ）システムは、熱力学、流体力学及び／又は熱伝達の原理に基づいて、蒸気圧縮サイクルを用いて、システムのコンポーネントを通して冷媒を移動させるシステムを指す。空調システムは、建物の居住者に外気のみを供給するシステムから、建物の温度を制御するだけのシステム、温度及び湿度を制御するシステムまでの非常に広範なシステムのセットに及ぶ。
【００２１】
「空調システムのコンポーネント」は、制御システムによって制御可能な動作を有するシステムの任意のコンポーネントを指す。それらのコンポーネントは、限定はしないが、冷媒を圧縮し、システムの中に圧送するための可変速度を有する圧縮器と、システムの高圧部分と低圧部分との間で調整可能な圧力降下を与えるための膨張弁と、それぞれ熱交換器を通る空気流量を調整するための変速ファンを組み込む、蒸発熱交換器及び凝縮熱交換器とを含む。
【００２２】
「蒸発器」は、熱交換器を通り抜ける冷媒が熱交換器の長さにわたって蒸発し、その結果、熱交換器の出口における冷媒の比エンタルピーが熱交換器の入口における冷媒の比エンタルピーより高くなり、冷媒が全体として液体から気体に変化する、蒸気圧縮システム内の熱交換器を指す。空調システム内に１つ以上の蒸発器が存在する場合がある。
【００２３】
「凝縮器」は、熱交換器を通り抜ける冷媒が熱交換器の長さにわたって凝縮し、その結果、熱交換器の出口における冷媒の比エンタルピーが熱交換器の入口における冷媒の比エンタルピーより低くなり、冷媒が全体として気体から液体に変化する、蒸気圧縮システム内の熱交換器を指す。空調システム内に１つ以上の凝縮器が存在する場合がある。
【００２４】
「制御信号のセット」又は「制御入力」は、蒸気圧縮システムのコンポーネントの動作を制御する特定の入力値を指す。制御信号のセットは、限定はしないが、圧縮器の速度値、膨張弁の位置、蒸発器内のファンの速度、及び凝縮器内のファンの速度を含む。
【００２５】
「設定点」は、空調システム等のシステムが動作の結果として到達し、保持することを目標とする目標値を指す。設定点という用語は、特定の制御信号のセット、並びに熱力学パラメータ及び環境パラメータの任意の特定の値に適用される。
【００２６】
「中央処理ユニット（ＣＰＵ）」又は「プロセッサ」は、ソフトウェア命令を読み出して実行するコンピュータ又はコンピュータのコンポーネントを指す。
【００２７】
「モジュール」又は「ユニット」は、タスク又はタスクの一部を実行するコンピュータにおける基本コンポーネントを指す。「モジュール」又は「ユニット」は、ソフトウェアによってもハードウェアによっても実施することができる。
【図面の簡単な説明】
【００２８】
【図１Ａ】１つの実施形態による、調整環境内で気流を生成する空調システム１００の動作を制御するシステム９９のブロック図である。
【図１Ｂ】いくつかの実施形態に従って調整される部屋１６０を空調する一例の概略図である。
【図１Ｃ】いくつかの実施形態による制御システム９９のブロック図である。
【図２Ａ】いくつかの実施形態による、気流動態のモデル１１０の概略図である。
【図２Ｂ】いくつかの実施形態による、気流動態の物理モデルを、有限次元空間上に投影することの概略図である。
【図３】いくつかの実施形態による、ＰＤＥモデル３７５とＯＤＥモデル３３７との間の関係を示す概略図である。
【図４】いくつかの実施形態による、無限次元モデル３７５を、簡略化された有限次元モデル３３７内に投影する投影演算子３８７を選択することの概略図である。
【図５】気流動態のモデルの物理パラメータの不確実性を補償するためにいくつかの実施形態によって用いられるエネルギーベースの方法のブロック図である。
【図６】１つの実施形態による、空調システムを制御する反復方法６００のブロック図である。
【図７】いくつかの実施形態による、気流動態のモデルの少なくとも１つの利得を自動調節する方法のブロック図である。
【発明を実施するための形態】
【００２９】
図１Ａは、１つの実施形態による、調整環境内で気流を生成する空調システム１００の動作を制御するシステム９９のブロック図を示している。空調システム１００は、蒸発器ファン１１４、凝縮器ファン１１３、膨張弁１１１、及び圧縮器１１２等のコンポーネントのうちの１つ又は組み合わせを備えることができる。システム９９は、例えば、サーモスタットからの設定点１１５、及びセンサ１３０からの読み値を受信し、コンポーネントの動作を制御するための制御信号のセットを出力する役割を担うコントローラ１２０を用いて、空調システム１００を制御する。コントローラ１２０は、制御信号のセットを、対応するコンポーネントのための特定の制御入力のセットに変換するために制御デバイスのセットに動作可能に接続される。例えば、コントローラは、圧縮器制御デバイス１２２、膨張弁制御デバイス１２１、蒸発器ファン制御デバイス１２４及び凝縮器ファン制御デバイス１２３に接続される。
【００３０】
コントローラは、気流の現在の状態１５０を推定するように構成された推定器１４０に動作可能に接続される。種々の実施形態において、コントローラ１２０は、気流の現在の状態１５０を用いて空調システム１００を制御する。気流の状態の例としては、空調システムの動作中に調整される空気の部屋における流れの値と温度の値とのうちの一方又は組み合わせが挙げられる。このようにして、コントローラは、所与の熱負荷について設定点値が達成されるように、空調システムの動作を制御する。例えば、監視用コントローラは、空調システムの少なくとも１つのコンポーネントについて少なくとも１つの制御入力を決定及び／又は更新して、気流の状態１５０を用いて求められる性能指標を最適化する。
【００３１】
気流の現在の状態１５０を特定するために、推定器１４０は、気流動態のモデル１１０を用いる。モデル１１０は、気流の現在の状態における観測誤差を低減するために、気流動態のモデルに従って気流の以前の状態を時間的に前方に遷移させることによって、推定器１４０が気流の状態１５０を特定することができる。種々の実施態様において、観測誤差とは、点のセットにおける気流の測定値と、気流のモデルによる点のセットにおける気流の対応する推定値との間の誤差である。このようにして、観測誤差は、調整環境内で気流の状態を測定するセンサの数を削減することができる部分観測誤差である。
【００３２】
図１Ｂは、いくつかの実施形態に従って調整される部屋１６０を空調する一例の概略図を示している。この例では、部屋１６０は、ドア１６１と、少なくとも１つの窓１６５とを有する。部屋の温度及び気流は、換気ユニット１０１を通して、システム１００等の空調システムによって制御される。部屋内の所与の点における気流の速度を測定する少なくとも１つの気流センサ１３１、及び室温を測定する少なくとも１つの温度センサ１３５等の、センサ１３０のセットが室内に配置される。他のタイプの環境、例えば、複数のＨＶＡＣユニットを備える部屋、又は複数の部屋を有する住宅を考えることができる。
【００３３】
コントローラ１２０は、ハードウェアにおいて実現することもできるし、プロセッサ、例えば、マイクロプロセッサにおいて実行されるソフトウェアプログラムとしても実現することができ、プロセッサは、固定又は可変の制御周期サンプリング間隔において、センサ１３０からの測定値、所望の設定点１１５、気流の状態１５０を受信し、この情報を用いて、入力、例えば、空調システム１００を動作させる制御信号を決定する。
【００３４】
推定器１４０は、ハードウェアにおいて実現することもできるし、コントローラ１２０と同じ又は異なるプロセッサのいずれかであるプロセッサにおいて実行されるソフトウェアプログラムとしても実現することができ、プロセッサは、固定又は可変の制御周期サンプリング間隔において、センサ１３０からの測定値を受信し、新たな測定値及び以前の測定値並びに気流の以前の状態を用いて、気流の現在の状態１５０を特定する。
【００３５】
図１Ｃは、いくつかの実施形態による制御システム９９のブロック図を示している。システム９９は、コンピュータ、例えば、気流動態のモデル１１０を記憶するメモリ１９５に接続された単一の中央処理ユニット（ＣＰＵ）又は複数のＣＰＵプロセッサ１９０の形態におけるコンピュータを備える。制御システム９９は、空調システム１００の動作からのフィードバック信号１０３を用いる。いくつかの実施態様では、フィードバック信号は、センサ１３０の測定値によって表される。制御システム９９は、空調システム１００を動作させる制御信号１０４を決定する。
【００３６】
いくつかの実施形態は、気流の物理モデル、すなわち気流動態の物理特性を捕捉するモデルは、無限次元であり、リアルタイム制御応用において用いるには複雑すぎる、という認識に基づいている。したがって、いくつかの実施形態は、リアルタイム更新及び制御に適した低次元すなわち有限次元のモデル１１０を用いる。例えば、いくつかの実施形態において、モデルは、ブシネスク方程式による気流の物理モデルより少ない数のパラメータを有する次数低減モデルを含み、ブシネスク方程式は偏微分方程式（ＰＤＥ）であり、次数低減モデルは常微分方程式（ＯＤＥ）である。
【００３７】
そのために、いくつかの実施形態は、例えば、無限次元の複雑な物理モデルを、より低い次元、例えば、有限次元のより単純なモデルに縮約することを目指すモデル縮約問題に対処する。これに関連して、モデル縮約問題は、複雑なモデルによるシステムの性能測定値と、より単純なモデルによるシステムの性能測定値との間の誤差が最小化されるように、より単純なモデルの全てのパラメータ又は係数を決定することを意味する。例えば、気流モデル化推定及び制御の用途において、気流を表すＰＤＥモデルを、次元がより低く、より複雑でないＯＤＥモデルに変換する必要がある。そのような縮約によって、空調システム等の気流システムのリアルタイム制御が可能になる。
【００３８】
図２Ａは、いくつかの実施形態による、気流動態のモデル１１０の概略図を示している。それらの実施形態では、モデル１１０は、気流の以前の状態を気流の現在の状態に遷移させる第１項２１０と、この遷移を支援する第２項２２０との組み合わせを含む。いくつかの実施形態では、第１項は、気流動態の物理モデルを有限次元空間上に投影することを含む。気流動態の物理モデルは、気流の以前の状態を、気流の現在の状態へと時間的に前方に遷移させるので、第１項２１０も、そのような遷移を支援する。加えて、気流の物理モデルは、調整環境の物理パラメータを含むので、種々の実施形態において、この投影は、第１項において調整環境の物理パラメータを維持する。
【００３９】
図２Ｂは、いくつかの実施形態による、気流動態の物理モデルを、有限次元空間上に投影することの概略図を示している。いくつかの実施形態は、気流動態の物理モデル２３０は、２つのタイプのパラメータ、すなわち、内部パラメータ２３３及び外部パラメータ２３５を有する、という理解に基づいている。内部パラメータ２３３とは、モデルの外部では何の意味も持たないモデルの内部利得を指す。外部パラメータ２３５とは、気流に影響を与える世界の物理パラメータであり、気流動態のモデルからは独立して存在している。このような物理パラメータの例としては、調整環境内の空気の粘度、調整環境の幾何学的形状、調整環境内の物体、例えば人物等の熱源の数及びタイプうちの１つ又は組み合わせが挙げられる。
【００４０】
いくつかの実施形態は、モデル複雑度の低減、例えば、気流動態の次数低減モデル２４０をもたらす、物理モデル２３０の有限空間への投影２３７は、外部パラメータ２３５を維持するべきであり、一方、内部パラメータ２４３は内部パラメータ２３３とは異なるものとすることができる、という認識に基づいている。このようにして、複雑度低減モデル２４０は、物理モデル２３０をより正確に表すことができる。
【００４１】
しかしながら、物理モデルの物理パラメータは、多くの場合に曖昧であり、すなわち不確実性を含む。例えば、空気の粘度は変動する可能性があり、及び／又は、異なる数の人物が、異なる時点に調整環境内に存在する可能性がある。いくつかの実施形態は、物理パラメータの値が、本明細書において有界不確実性の範囲と呼ばれる、有界不確実性を有して不確実である、という認識に基づいている。しかしながら、有界不確実性の範囲内では、物理パラメータの値は未知であり、変動する可能性がある。
【００４２】
そのために、気流動態のモデル１１０の第１項２１０は、有界不確実性の範囲を定義する有界不確実性を有して不確実である少なくとも１つの物理パラメータ２３５の値を有する次数低減モデル２４０を含む。例えば、いくつかの実施態様では、物理パラメータ２３５は、調整環境内の空気の粘度である。調整環境又は類似の環境内の空気の粘度の統計解析は、有界不確実性の範囲に、望ましい統計的信頼度を与えることができる。しかしながら、この範囲内で、粘度の値は変動する可能性がある。
【００４３】
いくつかの実施形態は、調整環境における観測誤差のエネルギー関数に、有界不確実性の範囲内の物理パラメータの任意の値について減少することを強制することによって、気流動態のモデルの物理パラメータについての知識の欠如を補償することができる、という理解に基づいている。そのために、いくつかの実施形態は、エネルギー関数及び時間に対するその導関数を解析するのにリアプノフベースの方法を用いる。いくつかの実施形態は、このような手法によって、気流の状態推定に対する物理パラメータの不確実性の影響を緩和することが可能になる、という理解に基づいている。
【００４４】
リアプノフ関数の導関数の解析中に、いくつかの実施形態は、不確実性の実際の値に関係なく、リアプノフ関数の導関数を常に負にする項を発見した。この項は、不確実性の実際の値に関係なく、リアプノフ関数の導関数を常に負にするために、有界不確実性の範囲、負の利得、及び、点のセットにおける気流の測定値と気流のモデルによる点のセットにおける気流の推定値との間の観測誤差の積を含み、これによって、さらに、不確実性の値に関係なく、エネルギー関数を時間反復にわたって減少させ、ひいては、有界不確実性の範囲内の物理パラメータの任意の値について、観測誤差を時間反復にわたってゼロまで減少させる。
【００４５】
そのために、モデル１１０の第２項２２０は、有界不確実性の範囲、負の利得、及び点のセットにおける気流の測定値と、気流のモデルによる点のセットにおける気流の推定値との間の観測誤差の関数を含む。このようにして、気流動態のモデルにおいて第１項と第２項とを組み合わせることによって、精度を維持しながらモデルの複雑度を低減することが可能になる。
【００４６】
種々の実施形態において、気流の状態は反復的に特定される。第２項は、観測誤差のエネルギー関数の導関数を強制的に負定値にするように設計されるので、気流の状態の反復更新は、複数の制御ステップにわたって観測誤差を強制的に低減させる。
【００４７】
いくつかの実施形態では、気流動態の物理モデル２３０は、気流の状態を時間的に前方に遷移させる偏微分方程式（ＰＤＥ）によって記述することができる。対照的に、気流動態の次数低減モデル２４０と、特にその第１項２１０とは、気流の状態を時間的に前方に遷移させる常微分方程式（ＯＤＥ）である。
【００４８】
図３は、いくつかの実施形態による、ＰＤＥモデル３７５とＯＤＥモデル３３７との間の関係を示す概略図を示している。いくつかの実施形態では、ＰＤＥモデル３７５は、無限次元であるのに対し、ＯＤＥモデル３３７は、有限次元である。いくつかの実施形態は、投影演算子３８７を用いて低次元空間、すなわち有限次元空間上にＰＤＥモデルを投影して、ＯＤＥモデル３３７に到達する。このようにして、気流動態のＯＤＥモデルにより、気流動態を空調システムのリアルタイム制御に用いることができる。
【００４９】
いくつかの実施態様では、ＰＤＥモデル３７５は、ブシネスク方程式によって定義される。ブシネスク方程式は、
【数１】
2021508120.tif 000003
によって与えられる。
【００５０】
この方程式では、Ｔは温度スカラー変数であり、
【数２】
2021508120.tif 000004
は３次元における速度ベクトルである。ここで、μは粘度であるとともにレイノルド数の逆数である。ｋは熱拡散係数である。ｐは圧力スカラー変数を表す。演算子Δ及び∇は、
【数３】
2021508120.tif 000005
と定義される。
【００５１】
ブシネスク方程式は、部屋内の気流と温度との間の結合を記述するが、無限次元である。加えて、ブシネスク方程式は、ＰＤＥの境界条件が経時的に変化する、すなわち、ＰＤＥに関与する係数のうちの１つが変化するという不確実条件を有する。こうした数学的変化は、部屋内のいくつかの実際の変化を反映している。例えば、部屋内の窓又はドアの開閉は、ブシネスクＰＤＥの境界条件を変化させる。同様に、経日変化及び季節的変化等の天候変化は、部屋内の温度と部屋外の温度との差に影響を与える。これによって、さらに、ＰＤＥ係数のうちのいくつかが影響を受け、例えば、レイノルド数が影響を受ける可能性がある。
【００５２】
図４は、いくつかの実施形態による、無限次元モデル３７５を、簡略化された有限次元モデル３３７に投影する投影演算子３８７を選択することの概略図を示している。種々の実施形態において、投影演算子は、ＰＤＥモデルの物理パラメータを維持する必要がある。そのために、いくつかの実施形態は、次数低減モデル（ＲＯＭ）３３７と呼ばれることもある有限次元モデルのために、投影演算子３８７ベースの基底関数４６０を選択する。基底関数は、直接数値シミュレーション（ＤＮＳ：direct numerical simulation）と呼ばれる大量の数値解によって得られるＰＤＥの厳密解のスナップショットから得られる。
【００５３】
異なる実施形態は、投影演算子３８７を求めるのに異なる種類の基底関数４６０を用いる。例えば、１つの実施形態は、固有直交分解（ＰＯＤ）基底関数４４０を用いる。別の実施形態は、動的モード分解（ＤＭＤ）基底関数４５０を用いる。投影演算子の数学的定式化の一例が以下に提供される。
【００５４】
投影演算子３８７が求められた後、いくつかの実施形態は、演算子３８７を用いて、無限次元のＰＤＥモデル３７５を、ＰＤＥモデルの物理パラメータを維持しながら、有限次元のＯＤＥモデル３３７に変換する。
【００５５】
図５は、気流動態のモデルの物理パラメータの不確実性を補償するためにいくつかの実施形態によって用いられるエネルギーベースの方法のブロック図を示している。例えば、１つの実施形態は、投影３８７によって求められたＯＤＥモデル３３７の所与のエネルギー値を評価する（５２０）ために、リアプノフ関数５４０を用いる。評価５２０は、モデル３３７とともに用いられると、物理パラメータの不確実性に対して耐性がある補正されたモデル５３０をもたらす補正項５１０を導出するのに用いられる。
【００５６】
例えば、いくつかの実施形態では、エネルギー関数及び時間に対するその導関数を解析する。例えば、こうした実施形態では、エネルギー関数の導関数を強制的に負定値にする補正項を設計する。いくつかの実施形態は、このような手法によって、気流の状態の推定に対する物理パラメータの不確実性の影響を緩和することができる、という理解に基づいている。
【００５７】
いくつかの実施態様では、エネルギー関数は、観測誤差の二乗の積分のリアプノフ関数である。これにより、リアプノフ関数の導関数の解析により、不確実性の実際の値に関係なく、補正項５１０はリアプノフ関数の導関数を常に負にすることができることが実証される。
【００５８】
そのため、種々の実施形態において、この補正項は、有界不確実性の範囲、負の利得、及び、点のセットにおける気流の測定値と気流のモデルによる点のセットにおける気流の推定値との間の観測誤差の積を含む。１つの実施形態では、気流の状態は、複数の制御ステップについて調整環境内で反復的に更新され、観測誤差は徐々に減少していく。
【００５９】
図６は、１つの実施形態による、空調システムを制御する反復方法６００のブロック図を示している。この実施形態は、補正されたＯＤＥモデル５３０を時間的に前に進めて、調整環境内の点／場所のセットにおける気流の状態を推定する。この点のセットは、センサ１３０による気流の状態の測定値６３７に対応し、測定値６３７と推定値６７５との間の差を、本明細書においては観測誤差６１０と呼ぶ。
【００６０】
具体的には、方法６００は反復的に実行されるので、この方法は、気流の以前の状態、すなわち、以前の制御ステップについて特定された状態６７５を、気流動態のモデル５３０に従って時間的に前方に遷移させることによって、調整環境内の気流の現在の状態６７５を特定する。以前の状態は、気流の現在の状態における観測誤差６１０を低減させる等のために遷移される。
【００６１】
１つの実施形態は、測定値６３７に対応する点のセットのみについての観測誤差を最小化しながら、調整環境全体、すなわち部屋６８５の全部の場所における気流の現在の状態６７５を推定する。１つの実施態様は、計算を簡略化するためにデータフィッティング技法を用いる。さらに又は代替的に、１つの実施形態は、部屋６８５の全ての場所において気流の全ての値を推定するのに逆投影演算子６９５を用いる。元々ＰＤＥの解をより低次元空間に投影するために用いられた投影が、今度は、低次元空間から、ＰＤＥの解、すなわち部屋全体を包含した、高次元空間における気流を復元するように反転される。
【００６２】
調整環境全体の気流の状態６８５が特定されると、気流の状態を用いて、調整環境内の何らかの所望の温度及び快適レベルを達成するように空調システムの動作が制御される（６９０）。例えば、気流の状態を、所望の快適環境について取得された気流状態の所望の設定点と比較することができ、その後、所望の気流設定点と気流の状態との間の誤差を空調ユニットに送信して、部屋内の温度又は換気を調節することができる。
【００６３】
いくつかの実施形態は、気流の補正されたモデルのいくつかの利得も、状態推定の精度を改善するために更新することができる、という理解に基づいている。例えば、補正項５１０に起因して、気流動態のＯＤＥモデル５３０は、物理パラメータの不確実性を補償するために負の利得を有することができる。しかしながら、負の利得の最適値を、オフラインで明示的に求めるのは困難になり得る。
【００６４】
例えば、いくつかの実施形態は、気流動態のモデルの少なくとも１つの利得を反復的に更新して観測誤差を低減するように構成される。この実施形態は、利得を求めるのにデータ駆動型最適化を用いることができる。データ駆動型最適化の例として、極値探索最適化、ガウス過程信頼上限（Gaussian process upper-confidence bound）最適化、及び疑似焼きなまし最適化のうちの１つ又は組み合わせが挙げられる。
【００６５】
図７は、いくつかの実施形態による、気流動態のモデルの利得を自動調節する方法のブロック図を示している。図７の方法は、ＯＤＥモデル５３０によって推定された気流の値を用い、その後、推定値を気流の測定値７３７と比較する。これらの値の間の観測誤差７６０を用いて、性能コスト関数７１０が評価される。コスト関数７１０は、安定性制約７２５を満たしながら、補正項５１０の利得の最良のあり得る値を探索する最適化７２０において用いられる。次に、利得の更新された値７３０が、ＯＤＥモデル５３０を更新するのに用いられる。
【００６６】
最適化７２０は、任意のデータ駆動型最適化７４４とすることができる。例えば、１つの実施形態は、極値探索最適化７４１を用いる。別の実施形態は、ガウス過程信頼上限（ＧＰ−ＵＣＢ）最適化７４２を用いる。
【００６７】
気流モデルの第１項の例
いくつかの実施形態は、以下の形式、
【数４】
2021508120.tif 000006
の非線形ＰＤＥ３７５システムの場合の状態推定問題を検討し、式中、ｚ０∈Ｄ（Ａ）⊂Ｈであり、Ａは、Ｃ０半群をヒルベルト空間Ｈ上に生成する線形演算子であり、Ｂ：Ｒｍ→Ｈは入力演算子であり、Ｃ：Ｄ（Ａ）→Ｒｐは測定演算子であり、ｈは、以下を満たすより高次の項を含む。
【００６８】
関数ｈ：Ｄ（Ａ）×Ｒｍ→［Ｄ（Ａ）］は、ｈ（０，０）＝０と、対（ｚ，ｕ）ごとに、
【数５】
2021508120.tif 000007
及び
【数６】
2021508120.tif 000008
を満たす全ての
【数７】
2021508120.tif 000009
について、
【数８】
2021508120.tif 000010
となるような正の定数εｚ、εｕ、Ｌｚ及びＬｕが存在するという局所的リプシッツ仮定とを満たす。
【００６９】
低次元部分空間
【数９】
2021508120.tif 000011
について、いくつかの実施形態は、以下の構造を有する推定器のＯＤＥモデル３３７を検討する。
【数１０】
2021508120.tif 000012
ただし、
【数１１】
2021508120.tif 000013
であり、
【数１２】
2021508120.tif 000014
、    
【数１３】
2021508120.tif 000015
、    
【数１４】
2021508120.tif 000016
及び
【数１５】
2021508120.tif 000017
が求められる。
【数１６】
2021508120.tif 000018
についての可能な選択肢は、Ａの支配固有関数のセット（モード近似）、又は（１）のシミュレーションの集合の固有直交分解（ＰＯＤ）を実行し、トランケートすることによって取得される基底関数のセット（ＰＯＤ近似）によって張られた空間とすることができる。
【００７０】
Ｈから
【数１７】
2021508120.tif 000019
への直交投影子３８７を、
【数１８】
2021508120.tif 000020
とする。したがって、状態推定における誤差を
【数１９】
2021508120.tif 000021
と定義する。
【００７１】
【数２０】
2021508120.tif 000022
が、ｒ個の支配ＰＯＤ基底関数４４０のスパンであり、かつ、ＴＰＯＤ３８７が、固有軌道ｚについての対応する投影であるとき、ＴＰＯＤは、全ての投影Ｔにわたる投影誤差
【数２１】
2021508120.tif 000023
を、次元ｒのＨの部分空間に最小化する。したがって、ＴＰＯＤの特定の選択は自然なもの４４０であるが、他の基底関数４６０を用いて投影を設計することができる。例えば、ＤＭＤ基底関数を用いることができる（４５０）。
【００７２】
異なる実施態様は、推定器（２）におけるＢｃ及びＧの異なる値を選択する。１つの実施形態は、収束を保証するために、全ての
【数２２】
2021508120.tif 000024
及びｕ∈Ｒｍについて以下の関係、すなわち
【数２３】
2021508120.tif 000025
及び
【数２４】
2021508120.tif 000026
を構築する。
【００７３】
気流モデルの第２項の例
第２項を導出する（２３７）ために、いくつかの実施形態は、補正項２２０をＯＤＥモデル３３７に加算することによって、前節で設計された公称推定器をロバストにする（robustify）ために制御理論を用いる。
【００７４】
ＰＤＥシステム（１）３７５は、ｈに対する不確実性を以下のように認めるものとする。
【数２５】
2021508120.tif 000027
ただし、Δｈ：Ｈ→Ｈが以下の仮定を満たすことを条件とする。
【００７５】
一方、不確実性は以下を満たす。
【数２６】
2021508120.tif 000028
【００７６】
推定器（２）の動態によれば、推定器の収束は、ベクトル場Ｇの設計に依拠する。上記で提示した公称設計をロバストにするために、いくつかの実施形態は、リアプノフ再設計手法を用い、不確実項Δｈを補償するために追加の項をΔＧ２２０に加算する。
【００７７】
ロバストなＯＤＥ推定器５３０は、
【数２７】
2021508120.tif 000029
と記述される。ただし、
【数２８】
2021508120.tif 000030
は、ｅの指数安定性に対するΔｈの影響を補償するように設計される。
【００７８】
この事例では、関連した誤差動態は、
【数２９】
2021508120.tif 000031
と記述され、ただし、ΔＧ２２０は、ｋ＜０と、
【数３０】
2021508120.tif 000032
を満たす任意の
【数３１】
2021508120.tif 000033
とについて、
【数３２】
2021508120.tif 000034
によって与えられる。
【００７９】
気流モデルの学習ベース調節の例
上記で提示した受動的ロバスト化２２０は、漸近的性能を保証する。しかしながら、いくつかの実施形態は、推定器の過渡的性能を改善することも対象としている。この理由から、いくつかの実施形態は、内部パラメータ、例えば推定器フィードバック利得ｋ７３０を学習する、図７に関連して示した能動学習によって前の受動的ロバスト推定器を補足することによって、この推定器を改善する。
【００８０】
そのために、いくつかの実施形態は、受動的ロバストＯＤＥモデル５３０を能動学習と統合して、モデル５３０の性能を改善する。実際、オンライン調節から利益を得ることができる１つのパラメータは、ロバストｏｄｅモデル利得ｋ７３０である。
【００８１】
ロバストｏｄｅモデル利得７３０のこの最適値を見つけるために、いくつかの実施形態は、ロバストｏｄｅモデルがシステム状態を推定しながら、利得をオンラインで自動調節するのにデータ駆動型最適化７２０を用いる。フィードバック利得は、
【数３３】
2021508120.tif 000035
と記述される。式中、ｋｎｏｍｉｎａｌは推定利得の公称値を表しており、δｋは推定器の過渡的性能を改善するために利得の必要な調整である。
【００８２】
学習コスト関数７１０は、
【数３４】
2021508120.tif 000036
であり、式中、Ｔ＞０であり、
【数３５】
2021508120.tif 000037
は、推定器（７）、（９）の解であり、ｙは測定された出力である。さらに、解析目的で、Ｑについて以下の仮定が必要である。すなわち、（１２）におけるコスト関数Ｑ（δｋ）は、δｋ＝δｋ＊において極小値を有する。
【００８３】
１つの実施形態は、δｋを調節するために、以下の時変振幅ベース極値探索アルゴリズム７４１を用い、
【数３６】
2021508120.tif 000038
式中、δｋ＞０、ωｋ＞０、εｋ＞０は調節利得である。
【００８４】
別の実施形態は、ロバストｏｄｅモデル利得７３０を最適化するのにガウス過程信頼上限アルゴリズム７４２を用いる。更に他の実施形態は、他のデータ駆動型最適化アルゴリズム７４４、例えば、強化学習アルゴリズム、ニューラルネットワーク、ディープニューラルネットワーク等を用いる。
【００８５】
本発明の上記の実施形態は数多くの方法のいずれかにおいて実現することができる。例えば、それらの実施形態は、ハードウェア、ソフトウェア又はその組み合わせを用いて実現することができる。ソフトウェアにおいて実現されるとき、そのソフトウェアコードは、単一のコンピュータ内に設けられるにしても、複数のコンピュータ間に分散されるにしても、任意の適切なプロセッサ、又はプロセッサの集合体において実行することができる。そのようなプロセッサは集積回路として実現することができ、集積回路コンポーネント内に１つ以上のプロセッサが含まれる。しかしながら、プロセッサは、任意の適切な構成の回路を用いて実現することができる。
【００８６】
また、本発明の実施形態を、方法としても具現化することができ、その一例が提供された。その方法の一部として実行される動作は、任意の適切な方法において順序化することができる。したがって、例示的な実施形態において順次の動作として示される場合であっても、例示されるのとは異なる順序において動作が実行される実施形態を構成することもでき、異なる順序は、いくつかの動作を同時に実行することを含むことができる。
【００８７】
請求項要素を変更するために特許請求の範囲において「第１の」、「第２の」等の序数の用語を使用することは、それだけで、或る請求項要素が別の請求項要素よりも優先度が高いこと、優位であること、若しくは上位にあることを、又は方法の動作が実行される時間的な順序を暗示するのではなく、請求項要素を区別するために、或る特定の名称を有する１つの請求項要素を（序数用語を使用しなければ）同じ名称を有する別の要素から区別するラベルとして単に使用される。
【図１Ａ】
2021508120.tif 000039
【図１Ｂ】
2021508120.tif 000040
【図１Ｃ】
2021508120.tif 000041
【図２Ａ】
2021508120.tif 000042
【図２Ｂ】
2021508120.tif 000043
【図３】
2021508120.tif 000044
【図４】
2021508120.tif 000045
【図５】
2021508120.tif 000046
【図６】
2021508120.tif 000047
2021508120.tif 000048
【図７】
2021508120.tif 000049
【国際調査報告】
2021508120.tif 000050
2021508120.tif 000051