(57)【特許請求の範囲】
【請求項１】
金属内部の電磁界を計算せずに電磁波の過渡的な挙動を数値解析する電磁波解析装置において、
金属の厚さをｄ、電気伝導率をσとした場合に、角周波数が０に近づくとき金属領域の表面インピーダンスが１／ｄσに近づくような表面インピーダンス境界条件に基づいて、有限差分時間領域法により電磁波の過渡的な挙動を計算する解析手段、
を有することを特徴とする電磁波解析装置。
【請求項２】
前記解析手段は、虚数単位をｊ、金属の透磁率をμ、金属の電気伝導率をσ、金属の厚さをｄ、角周波数をωとした場合に、表面インピーダンスＺ（ω）が、
【数１】
0003636586.tif 000024
で定義された表面インピーダンス境界条件に基づいて、有限差分時間領域法により電磁波の過渡的な挙動を計算することを特徴とする請求項１記載の電磁波解析装置。
【請求項３】
表面インピーダンス境界条件を考慮した電磁波解析を行うべき領域を指定する初期設定手段をさらに有し、
前記解析手段は、前記初期設定手段で指定された領域の電磁波の過渡的な挙動を計算することを特徴とする請求項１記載の電磁波解析装置。
【請求項４】
金属内部の電磁界を計算せずに電磁波の過渡的な挙動を数値解析する電磁波解析プログラムを記録したコンピュータ読み取り可能な記録媒体において、
金属の厚さをｄ、電気伝導率をσとした場合に、角周波数が０に近づくとき金属領域の表面インピーダンスが１／ｄσに近づくような表面インピーダンス境界条件に基づいて、有限差分時間領域法により電磁波の過渡的な挙動を計算する解析手段、
としてコンピュータを機能させることを特徴とする電磁波解析プログラムを記録したコンピュータ読み取り可能な記録媒体。
【発明の詳細な説明】
【０００１】
【発明の属する技術分野】
本発明は電磁波の過渡的な挙動を数値解析する電磁波解析装置及び電磁波解析プログラムを記録したコンピュータ読み取り可能な記録媒体に関し、特に金属領域を有限の電気伝導率を有する媒体として扱う電磁波解析装置及びそのような電磁波解析をコンピュータに行わせるための電磁波解析プログラムを記録したコンピュータ読み取り可能な記録媒体に関する。
【０００２】
【従来の技術】
電子計算機を使った数値計算によって電磁波の過渡的な挙動を解析する手法の１つに有限差分時間領域( ＦＤＴＤ：Finite Difference Time Domain ）法がある。この手法は、マックスウェル方程式を時間と空間とについて差分法で解く手法である。ＦＤＴＤ法は、適用範囲の広さ等から多くの場面で利用されている。
【０００３】
ＦＤＴＤ法を用いて電磁波解析を行う際、解析領域に金属が存在しているときの計算方法としては、以下のような手法がある。
第１の手法は、金属を完全導体、すなわち電気伝導率が無限大であるものとみなして計算を行う手法である（これを「第１の従来例」とする）。この手法では、金属領域の表面に境界条件として電気壁境界条件を設定し「Ｅtan ＝０」とする。ここで、「Ｅtan 」は、金属領域表面上の電界の接線成分を表す。この手法は最も一般的な手法であり、金属内部の電磁界を計算する必要がないため便利である。ただし、この手法を用いた場合、導体損失を考慮できないため、導体損失が無視できないような場合には対応できない。
【０００４】
そこで、金属を完全導体とみなせない場合には、金属領域の電磁界を考慮して計算をする第２の手法で解析できる。
第２の手法は、金属領域に相当する格子（メッシュ）を他の領域の格子と同様に扱い、その金属の電気伝導率を設定する手法である（これを「第２の従来例」とする）。この手法によれば、金属内部の電磁界も計算されるため、導体損失を考慮した解析が可能である。
【０００５】
なお、上記の第１の従来例と第２の従来例との詳細は、「北村他、信学会論文誌  C-I Vol.J76-C-I No5 pp173-180,１９９３年５月」に記載されている。
ところが、第２の従来例の手法を用いた場合、電磁界は金属内部に行くにしたがって指数関数的に急激に変化するので、金属領域の空間離散間隔を非常に小さくしなければならない。しかも、空間離散間隔を小さくした場合には同時に時間離散間隔をも小さくする必要が生ずる。その結果、計算に必要なメモリ容量と計算時間の著しい増加を招くという問題点がある。
【０００６】
したがって、限られたハードウェア資源を用いた電磁波解析を短時間で行うには、上記第２の従来例に代わる別の計算手法が必要となる。そこで、次のような第３の手法が考えられている。
【０００７】
第３の手法は、金属領域表面にＢｅｇｇｓらの示した表面インピーダンス境界条件を設定する手法である（これを「第３の従来例」とする）。金属の厚さ方向をｘ軸とした場合、ｘ軸に垂直な方向の電磁波Ｅy の定常解析（周波数領域）の式は、次のように表される。
【０００８】
【数２】
0003636586.tif 000002
【０００９】
ここで、ωは角周波数、μは金属の透磁率、σは電気伝導率、ｊは虚数単位である。この式（２）に基づいてＦＤＴＤ法による過渡解析（時間領域）を行えば、電磁波の過渡的な挙動が解析できる。この手法によれば、金属内部の電磁界を計算する必要がなく、かつ導体損失を考慮することができる。この手法の詳細は、「J.H.Beggs et al.,IEEE Trans. Antennas ＆Propag.,AP40,No.1,pp.49-118,jpn. 1992 」に記載されている。
【００１０】
【発明が解決しようとする課題】
しかし、第３の従来例の手法を用いた場合、Ｂｅｇｇｓらの表面インピーダンス境界条件は金属の厚さが表皮厚δ（δ＝（２／μσω）1/2 ）よりも十分に厚いことが前提である。そのため、金属の厚さが表皮厚に比べて十分でない場合、例えば低周波数の電磁界成分に対しては適用できないという問題点がある。
【００１１】
本発明はこのような点に鑑みてなされたものであり、金属の導体損失を考慮した電磁波の過渡的な挙動の計算を、金属の厚さ及び角周波数に依存せずに高速に行うことができる電磁波解析装置を提供することを目的とする。
【００１２】
また、本発明の他の目的は、金属の導体損失を考慮した電磁波の過渡的な挙動のコンピュータによる計算を、金属の厚さ及び角周波数に依存せずに高速に行わせることができる電磁波解析プログラムを記録したコンピュータ読み取り可能な記録媒体を提供することである。
【００１３】
【課題を解決するための手段】
本発明では上記課題を解決するために、金属内部の電磁界を計算せずに電磁波の過渡的な挙動を数値解析する電磁波解析装置において、金属の厚さをｄ、電気伝導率をσとした場合に、角周波数が０に近づくとき金属領域の表面インピーダンスが１／δσに近づくような表面インピーダンス境界条件に基づいて、有限差分時間領域法により電磁波の過渡的な挙動を計算する解析手段、を有することを特徴とする電磁波解析装置が提供される。
【００１４】
この電磁波解析装置によれば、角周波数が０に近づくとき金属領域の表面インピーダンスが１／ｄσに近づくような表面インピーダンス境界条件が定められているため、表皮厚δ（δ＝（２／（μσω））1/2 ）が金属の厚さｄよりも大きな値であっても、電磁波の過渡的な挙動を解析できるとともに、金属内部の電磁界を計算しないため、空間離散間隔や空間離散間隔を特別小さくする必要がない。
【００１５】
また、上記課題を解決するために、金属内部の電磁界を計算せずに電磁波の過渡的な挙動を数値解析する電磁波解析プログラムを記録したコンピュータ読み取り可能な記録媒体において、金属の厚さをｄ、電気伝導率をσとした場合に、角周波数が０に近づくとき金属領域の表面インピーダンスが１／ｄσに近づくような表面インピーダンス境界条件に基づいて、有限差分時間領域法により電磁波の過渡的な挙動を計算する解析手段、としてコンピュータを機能させることを特徴とする電磁波解析プログラムを記録したコンピュータ読み取り可能な記録媒体が提供される。
【００１６】
このような記録媒体に記録された電磁波解析プログラムをコンピュータに実行させれば、上記本発明に係る電磁波解析装置に必要な機能がコンピュータにより実現される。
【００１７】
【発明の実施の形態】
以下、本発明の実施の形態を図面を参照して説明する。
図１は、本発明の原理図である。初期設定手段１は、金属の厚さｄ、メッシュ数、メッシュサイズ、各セルの誘電率ε、透磁率μ、電気伝導率σ等の初期設定を行う。また、初期設定手段１は、表面インピーダンス境界条件を考慮した電磁波解析を行うべき領域も指定する。解析手段２は、初期設定手段１で指定された領域の電磁波解析を、初期設定に基づいて行う。すなわち、表皮厚δが無限大（角周波数ωが０）に近づくとき、金属領域の表面インピーダンスＺ（ω）が１／ｄσとなるような表面インピーダンス境界条件を満たした方程式を有限差分時間領域法で解析することにより、電磁波の過渡的な挙動を計算する。そして、解析手段２は、電界の挙動を示す式を解析結果３として出力する。
【００１８】
これにより、金属の導体損失を考慮した電磁波の過渡的な挙動の計算を、金属の厚さに依存せずに高速に行うことができる。
以下に、解析手段２が計算すべき方程式の導出方法について説明する。
【００１９】
図２は、格子状に分割された領域内に置かれた金属１０を示す図である。ここでは、金属１０の厚さ方向にｘ軸をとり、ｘ軸に垂直な方向にｙ軸とｚ軸とをとっている。そして、この金属１０が置かれている３次元空間は、メッシュで区切られている。また、金属１０は厚さがｄであり、十分に広い金属膜であるものとする。金属１０に接した領域（セル）１１が、表面インピーダンス境界条件を考慮して解析を行うべき領域である。
【００２０】
まず、金属１０の断面を流れる電流の総和Ｉは次式で与えられる。
【００２１】
【数３】
0003636586.tif 000003
【００２２】
ここで、Ｉ0 は金属表面の電流密度、ｊは虚数単位、δは表皮厚、μは金属の透磁率、σは金属１０の電気伝導率、ωは角周波数である。すると、表面インピーダンスＺs （ω）は、次式のようになる。
【００２３】
【数４】
0003636586.tif 000004
【００２４】
ここで、Ｒs0は金属の直流でのシート抵抗で、ＲS0＝１／ｄσである。この式（４）は、膜厚が表皮厚程度あるいはそれ以下の金属膜に対しても適用可能な表面インピーダンス境界条件を示している。
【００２５】
式（４）をＦＤＴＤ法に適用するには、表面インピーダンスＺs （ω）を次のように変形する。
【００２６】
【数５】
0003636586.tif 000005
【００２７】
ここで、ａ＝ｄ（μσ）1/2 である。
ところで、金属１０を流れる電流は磁界で表すことができる。したがって、電磁波のｙ軸方向の成分Ｅy は、表面インピーダンスＺs （ω）と磁界のｚ軸方向成分Ｈz （ω）との積で表すことができ、次の式で表される。
【００２８】
【数６】
0003636586.tif 000006
【００２９】
この式（６）により、電磁波の定常解析（周波数領域の解析）を行うことができる。
電磁波の過渡解析（時間領域の解析）を行うには、式（５）におけるＺs ’（ω）に対してフーリエ逆変換を施す。すると、Ｚs ’（ｔ）が以下のように決定される。
【００３０】
【数７】
0003636586.tif 000007
【００３１】
ここで、Ｆ-1はフーリエ逆変換演算を、Θ3 は楕円テータ関数を表す。時間領域における表面インピーダンスと電磁界との関係は、
【００３２】
【数８】
0003636586.tif 000008
【００３３】
であるから、ＦＤＴＤ法においては、式（７）を用いて以下のように定式化できる。
【００３４】
【数９】
0003636586.tif 000009
【００３５】
ここで、τ＝Δｔαである。このようにして次の式
【００３６】
【数１０】
0003636586.tif 000010
【００３７】
が求められる。この式（１０）が、解析手段２による電磁波の過渡的な挙動の計算に用いられる方程式である。
なお、式（１０）のままでは計算開始時刻から毎時刻の磁界成分を全て記憶しておく必要が生じ、計算をする上で不便である。そこで、式（１０）中の積分を、
【００３８】
【数１１】
0003636586.tif 000011
【００３９】
とし、式（１１）を以下のようにＰｒｏｎｙ展開する。
【００４０】
【数１２】
0003636586.tif 000012
【００４１】
すると、式（１０）は、
【００４２】
【数１３】
0003636586.tif 000013
【００４３】
となる。この式（１３）を
【００４４】
【数１４】
0003636586.tif 000014
【００４５】
とおくと、次の式が得られる。
【００４６】
【数１５】
0003636586.tif 000015
【００４７】
これにより逐次的な計算が可能となる。すなわち、式（１５）によれば、現在より１ステップおよび２ステップ前の時刻の磁界成分を記憶しておけばよく、計算開始時刻から毎時刻の磁界成分を全て記憶しておく必要がなくなる。したがって、逐次的に計算をすることができ、コンピュータで計算する際のメモリ容量を低減することができる。
【００４８】
次に、本発明の電磁波解析装置における電磁波解析の処理手順を説明する。
図３は、電磁波解析の処理手順を示すフローチャートである。なお、以下の処理の中で、ステップＳ１のみが図１に示した初期設定手段１が行う処理であり、他の処理は解析手段２が行う。
〔Ｓ１〕初期設定を行う。初期設定の内容は、メッシュ数、メッシュサイズ、各セルの誘電率ε、透磁率μ、電気伝導率σ、及び解析時間Ｔmax 等である。また、時間Ｔの初期値としてＴ＝０としている。
【００４９】
このとき、各セルが境界条件に依存するか否か、及び境界条件に依存する場合には表面インピーダンス境界条件に依存するのか否かを設定する。すなわち、金属の周囲のセルは表面インピーダンス境界条件に依存するものとし、解析領域の端部のセルなどは吸収境界条件や周期境界条件などに依存するものとする。それ以外のセルは、境界条件に依存しないものとする。
〔Ｓ２〕表面インピーダンス境界条件のための初期設定を行う。具体的には、以下の処理を行う。
【００５０】
まず、表面インピーダンス境界条件を適用する面を設定し、面ごとに電気伝導率σと金属の厚さｄとを設定する。
次に、予め与えられているＮの値を設定し、面ごとに、
【００５１】
【数１６】
0003636586.tif 000016
【００５２】
の近似式を十分に満たす係数ｐk 、ｑk を算出する。そして、算出したｐk 、ｑk をｐ（ｋ）、ｑ（ｋ）（ｋ＝１〜Ｎ）のように配列として記憶する。
〔Ｓ３〕境界条件に依存しないセルの電界の計算を行う。
〔Ｓ４〕表面インピーダンス境界条件以外の境界条件に依存したセルの電界の計算を行う。
〔Ｓ５〕表面インピーダンス境界条件に依存したセルの電界の計算を行う。具体的には、次の式
【００５３】
【数１７】
0003636586.tif 000017
【００５４】
の計算をする。ここで、ｋ＝１〜Ｎであり、Ｈｙは表面インピーダンス境界条件に依存する電界Ｅとねじれの位置にある磁界である。そして、
【００５５】
【数１８】
0003636586.tif 000018
【００５６】
により表面インピーダンス境界条件に依存する電界Ｅを計算する。
〔Ｓ６〕時間ＴにΔＴ／２を加算して、磁界の計算を行う。
〔Ｓ７〕境界条件に依存した磁界の計算を行う。
〔Ｓ８〕時間Ｔが解析時間Ｔmax より大きいか否かを判断する。時間Ｔが解析時間Ｔmax を超えていれば、ステップＳ９に進み、超えていなければ時間ＴにΔＴ／２を加算してステップＳ３に進む。これにより、時間Ｔが解析時間Ｔmax に達するまで、電磁波の過渡解析が繰り返し行われる。
〔Ｓ９〕解析時間Ｔmax における電磁波の解析が終了したら、計算結果として得られた電界成分（３軸方向）と磁界成分（３軸方向）の時系列データを出力する。
【００５７】
なお、上記フローチャートのステップＳ３，Ｓ４，Ｓ６，Ｓ７で行われる計算は、表面インピーダンス境界条件を用いずにＦＤＴＤ法による電磁波解析を行う従来の手法でも行われていた計算であるため、計算内容の詳細については省略する。
【００５８】
以上のようにして、金属の導体損失を考慮した電磁波の過渡的な挙動の計算を、金属の厚さに依存せずに高速に行うことができる。すなわち、第２の従来例で説明した手法のようにメッシュを細かくする必要が無いため、第２の従来例に比べて計算に必要なデータ量が少なくですむ。その結果、少ないメモリで高速に演算可能である。また、第３の従来例と異なり、表皮厚に比べて金属の厚さが薄い場合においても電磁波の過渡的解析が可能である。すなわち、低周波の場合においても解析できる。
【００５９】
【実施例】
以下に、マイクロストリップ線路を伝播する電磁波を、本発明の電磁波解析装置により解析した場合の実施例を説明する。
【００６０】
図４は、マイクロストリップ線路を示す図である。このマイクロストリップ線路は、ＳｉＯＮからなる基板２１の両面に、Ａｕからなる薄膜の導体板２２と導体板２３とが形成されている。これらの周囲の領域３１，３２は、空気で満たされている。また、図中の外側の枠が吸収境界４１を示している。
【００６１】
基板２１の厚さは３．０μｍ、導体板２２の厚さは１．０μｍ、導体板２３の厚さは１．０μｍである。導体板２２の上の面から吸収境界４１までの距離は３．０μｍであり、導体板２２の幅は５．０μｍである。導体板２３の下面から吸収境界４１までの距離は２．０μｍである。また、Ａｕの電気伝導率σは、σ＝４．１６７×１０7 （Ｓ／ｍ）であり、ＳｉＯＮの比誘電率εr は、εr ＝５．０である。
【００６２】
ところで、周波数ｆが２４ＧＨｚの場合の表皮厚δは、δ≒０．５μｍ（δ＝１／（πｆμσ）1/2 ）である。従って、２４ＧＨｚよりも低周波の電磁波の解析を行う場合には、導体板２２，２３の厚さが表皮厚に比べて十分に厚いとはいえなくなる。そこで、図示したマイクロストリップ線路の導体損失と周波数との関係を、本発明の電磁波解析装置によって行った場合と、従来の手法を用いて行った場合とを比較した。以下の実施例では、直流から２０ＧＨｚまでの間の周波数帯に関して解析を行った。
【００６３】
なお、導体損失は以下のようにして求めることができる。
図５は、図４のマイクロストリップ線路の側面を示す図である。
導体損失loss( ω，Ｌ）は、所定の距離Ｌだけ離れた位置における２つの導体板間の電位差Ｖ（Ｚ＝０）、Ｖ（Ｚ＝Ｌ）から以下のようにして導き出せる。
【００６４】
【数１９】
0003636586.tif 000019
【００６５】
【数２０】
0003636586.tif 000020
【００６６】
ここで、Ｆはフーリエ変換、Ｒe は複素数の実部をとることを意味する。
図６は、導体損失の解析結果を示す図である。この図は、横軸に周波数（ＧＨｚ）を、縦軸に導体損失（ｄＢ）をとっている。
【００６７】
図中、実線５１で示しているのが測定値である。点線５２で示しているのが本発明による解析結果である。一点破線５３で示しているのが第２の従来例（金属領域に相当する格子を他の領域の格子と同様に扱う手法）による解析結果である。破線５４で示しているのが第３の従来例（金属領域表面にＢｅｇｇｓらの示した表面インピーダンス境界条件を設定する手法）である。
【００６８】
なお、第１の従来例では、金属を完全導体とみなすことが前提となっているため、導体損失を求めても常に「０」の値しか得られない。
この図から分かるように本発明により解析した場合と、第２の従来例により解析した場合には、実際の測定値に非常に近い値が得られた。ところが、第３の従来例により解析した場合には、周波数が低くなるほど誤差が大きくなる。
【００６９】
ここで、本発明による解析と第２の従来例による解析との間の、必要となるメモリ容量と計算の所用時間とを比較してみた。ここで、必要なメモリ容量や計算の所用時間には、メッシュの数が大きく影響する。
【００７０】
図７は、本発明と第２の従来例との間のメッシュの違いを比較する図である。（Ａ）は本発明のメッシュを示す図であり、（Ｂ）は第２の従来例のメッシュを示す図である。
【００７１】
本発明を用いた場合には、ｘ軸方向のメッシュ間隔の最小値ｘmin は０．５μｍ、ｙ軸方向のメッシュ間隔の最小値ｙmin は１．０μｍ、ｚ軸方向のメッシュ間隔の最小値ｚmin は４．０μｍである。
【００７２】
ここで、メッシュ間隔と時間離散間隔とは、以下のような条件を満たしていることが必要である。
【００７３】
【数２１】
0003636586.tif 000021
【００７４】
ここでｃは真空中の光速である。したがって、本発明を用いた場合の時間離散間隔ｄｔは、
【００７５】
【数２２】
0003636586.tif 000022
【００７６】
である。また、メッシュの数は「ｎｘ×ｎｙ×ｎｚ＝１２×１３×ｎｚ」となる。
一方、第２の従来例を用いた場合には、金属２２，２３内部の表面に近い部分ではメッシュ間隔を小さくしなければならない。そのため、ｘ軸方向のメッシュ間隔の最小値ｘmin とｙ軸方向のメッシュ間隔の最小値ｙmin とは、δ／４である。すると、ｘmin は０．１２５μｍ、ｙmin は０．１２５μｍ、ｚmin は４．０μｍとなる。このときの時間離散間隔は、
【００７７】
【数２３】
0003636586.tif 000023
【００７８】
である。また、メッシュの数は「ｎｘ×ｎｙ×ｎｚ＝２４×３６×ｎｚ」となる。
したがって、第２の従来例を用いた場合には、本発明を用いた場合と比べて約５．５倍の主記憶容量と約２７．７倍の計算時間が必要となる。
【００７９】
以上の実施例から、本発明の電磁波解析装置を用いれば、金属の膜厚に関係なく、低周波の電磁界の解析を少ないメモリ容量で高速に実行できることがわかる。
【００８０】
なお、上記の処理機能は、コンピュータによって実現することができる。その場合、電磁波解析装置が有すべき機能の処理内容は、コンピュータで読み取り可能な記録媒体に記録されたプログラムに記述されており、このプログラムをコンピュータで実行することにより、上記処理がコンピュータで実現される。コンピュータで読み取り可能な記録媒体としては、磁気記録装置や半導体メモリ等がある。市場を流通させる場合には、ＣＤ−ＲＯＭ(Compact Disk Read Only Memory) やフロッピーディスク等の可搬型記録媒体にプログラムを格納して流通させたり、ネットワークを介して接続されたコンピュータの記憶装置に格納しておき、ネットワークを通じて他のコンピュータに転送することもできる。コンピュータで実行する際には、コンピュータ内のハードディスク装置等にプログラムを格納しておき、メインメモリにロードして実行する。
【００８１】
【発明の効果】
以上説明したように本発明では、角周波数が０に近づくとき金属領域の表面インピーダンスが１／ｄσに近づくような表面インピーダンス境界条件が定められているため、表皮厚δが金属の厚さｄよりも大きな値であっても、電磁波の過渡的な挙動を解析できるとともに、金属内部の電磁界を計算しないため、空間離散間隔や空間離散間隔を特別小さくする必要がない。その結果、金属の導体損失を考慮した電磁波の過渡的な挙動の計算を、金属の厚さに依存せずに高速に行うことが可能となる。
【００８２】
また、本発明の電磁波解析プログラムを記録したコンピュータ読み取り可能な記録媒体では、記録された電磁波解析プログラムをコンピュータに実行させれば、角周波数が０に近づくとき金属領域の表面インピーダンスが１／ｄσに近づくような表面インピーダンス境界条件を定めた電磁波解析処理をコンピュータに行わせることができる。その結果、金属の導体損失を考慮した電磁波の過渡的な挙動の、金属の厚さに依存しない高速計算をコンピュータに行わせることができる。
【図面の簡単な説明】
【図１】本発明の原理図である。
【図２】格子状に分割された領域内に置かれた金属を示す図である。
【図３】電磁波解析の処理手順を示すフローチャートである。
【図４】マイクロストリップ線路を示す図である。
【図５】図４のマイクロストリップ線路の側面図である。
【図６】導体損失の解析結果を示す図である。
【図７】本発明と第２の従来例との間のメッシュの違いを比較する図である。（Ａ）は本発明のメッシュを示す図であり、（Ｂ）は第２の従来例のメッシュを示す図である。
【符号の説明】
１  初期設定手段
２  解析手段
３  解析結果
１０  金属
１１  セル
【図１】
0003636586.tif 000025
【図２】
0003636586.tif 000026
【図３】
0003636586.tif 000027
【図４】
0003636586.tif 000028
【図５】
0003636586.tif 000029
【図６】
0003636586.tif 000030
【図７】
0003636586.tif 000031