(57)【要約】  風向及び風速の所定期間における時系列計測値を取得し、前記時系列計測値をもとに、風向及び風速の範囲を、各々ｍ、ｎ個の分割範囲に分割し、風速の前記分割範囲毎に風速代表値を設定し、それぞれの前記分割範囲に含まれる前記時系列計測値の前記所定期間での頻度を求めることによってｍ×ｎ行列の風向・風速頻度分布を作成し、前記風速代表値及び前記風向・風速頻度分布を風向・風速情報とする風向・風速情報入力工程と；粉塵の発塵源の情報を入力する発塵源情報入力工程と；前記風向・風速情報と、前記発塵源の情報と、を用いて、任意の座標点における粉塵濃度ｃを計算する粉塵濃度計算工程と；前記粉塵濃度ｃに基づいて、任意の降下地点における降下粉塵量を計算する降下粉塵量計算工程と；を備える降下粉塵量推定方法。
2010001925.tif 000002
【特許請求の範囲】
【請求項１】
風向及び風速の所定期間における時系列計測値を取得し、前記時系列計測値をもとに、風向及び風速の範囲を、各々ｍ、ｎ個の分割範囲に分割し、風速の前記分割範囲毎に風速代表値を設定し、それぞれの前記分割範囲に含まれる前記時系列計測値の前記所定期間での頻度を求めることによってｍ×ｎ行列の風向・風速頻度分布を作成し、前記風速代表値及び前記風向・風速頻度分布を風向・風速情報とする風向・風速情報入力工程と；
粉塵の発塵源の情報を入力する発塵源情報入力工程と；
前記風向・風速情報と、前記発塵源の情報と、を用いて、任意の座標点における粉塵濃度ｃを計算する粉塵濃度計算工程と；
前記粉塵濃度ｃに基づいて、任意の降下地点における降下粉塵量を計算する降下粉塵量計算工程と；
を備えることを特徴とする降下粉塵量推定方法。
【請求項２】
前記粉塵濃度計算工程は、前記粉塵の地表での反射率βを更に用いて前記粉塵濃度ｃを計算することを特徴とする請求項１に記載の降下粉塵量推定方法。
【請求項３】
前記発塵源情報入力工程では、前記発塵源の３次元空間上での座標と、前記発塵源の発塵強度Ｑと、前記粉塵の粒子径と、前記粉塵の密度と、前記発塵源からの有効発塵高さＨｅと、が入力されることを特徴とする請求項２に記載の降下粉塵量推定方法。
【請求項４】
前記粉塵濃度計算工程は、３次元座標軸上の点（ｘ，ｙ，ｚ）における粉塵濃度ｃ（ｘ，ｙ，ｚ）を（１）式を用いて計算し；
2010001925.tif 000019
ｕは風速であり、Ｋｙ、Ｋｚは、各々、ｙ軸、ｚ軸方向の乱流拡散係数であり、ｗは粒子の終末沈降速度であり、前記発塵源のｘ座標は０であり、前記発塵源のｙ座標はｙ０である；
ことを特徴とする請求項３に記載の降下粉塵量推定方法。
【請求項５】
前記降下粉塵量計算工程は、（２）式又は（３）式を用いて降下粉塵量を計算し；
2010001925.tif 000020
2010001925.tif 000021
ｃ（ｘ，ｙ，ｚ）は、３次元座標軸上の点（ｘ，ｙ，ｚ）における粉塵濃度であり；
Ｃは降下粉塵量であり、
Ｋzはｚ軸方向の乱流拡散係数である；
ことを特徴とする請求項２に記載の降下粉塵量推定方法。
【請求項６】
風向及び風速の所定期間における時系列計測値を取得し、前記時系列計測値をもとに、風向及び風速の範囲を、各々ｍ、ｎ個の分割範囲に分割し、風速の前記分割範囲毎に風速代表値を設定し、それぞれの前記分割範囲に含まれる前記時系列計測値の前記所定期間での頻度を求めることによってｍ×ｎ行列の風向・風速頻度分布を作成し、前記風速代表値及び前記風向・風速頻度分布を風向・風速情報とする風向・風速情報入力部と；
粉塵の発塵源の情報が入力される発塵源情報入力部と；
前記風向・風速情報と、前記発塵源の情報と、を用いて、任意の座標点における粉塵濃度ｃを計算する粉塵濃度計算部と；
前記粉塵濃度ｃに基づいて、任意の降下地点における降下粉塵量を計算する降下粉塵量計算部と；
を備えることを特徴とする降下粉塵量推定装置。
【請求項７】
風向及び風速の所定期間における時系列計測値を取得し、前記時系列計測値をもとに、風向及び風速の範囲を、各々ｍ、ｎ個の分割範囲に分割し、風速の前記分割範囲毎に風速代表値を設定し、それぞれの前記分割範囲に含まれる前記時系列計測値の前記所定期間での頻度を求めることによってｍ×ｎ行列の風向・風速頻度分布を作成し、前記風速代表値及び前記風向・風速頻度分布を風向・風速情報とする風向・風速情報入力手順と；
粉塵の発塵源の情報を入力する発塵源情報入力手順と；
前記風向・風速情報と、前記発塵源の情報と、を用いて、任意の座標点における粉塵濃度ｃを計算する粉塵濃度計算手順と；
前記粉塵濃度ｃに基づいて、任意の降下地点における降下粉塵量を計算する降下粉塵量計算手順と；
を備えることを特徴とする降下粉塵量推定プログラム。
【請求項８】
請求項７に記載の降下粉塵量推定プログラムを記録したコンピューター読み取り可能な記録媒体。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、大気中の粉塵の移流・拡散挙動の計算結果に基づいた降下粉塵量の推定方法、装置、及びコンピュータプログラムに関する。
本願は、２００８年７月１日に、日本に出願された特願２００８−１７２５０１号に基づき優先権を主張し、その内容をここに援用する。
【背景技術】
【０００２】
工場の煙突等の発塵源（粉塵・煤塵の発生源）から放出される粉塵粒子の風による飛散挙動を計算し、降下粉塵量を推定することは、降下粉塵が地域住民へ及ぼす影響を評価する際に極めて重要となる。【０００３】
非特許文献１では、以下の方法が提案されている：降下粉塵量は、発塵源と降下点の１次元距離に依存すると仮定する；発塵源と降下点間の距離、粉塵の発生速度、粉塵の粒度分布と粉塵の密度、有効発塵高さ、及び、風向と風速との頻度分布を入力値として、下の（４）式の計算を行なう；（４）式の計算結果に風向・風速の頻度分布を乗じて、降下粉塵量を推定する。
また、実験式または経験式に基づき降下粉塵量を推定する装置及び方法が提案されている（特許文献１参照）。【０００４】
2010001925.tif 000003
【０００５】
Ｑは発塵強度（粉塵の発生速度）、Ｈｅは有効発塵高さ、ｗは粒子径で決定される終末沈降速度、ｘは発塵源と降下点の風下距離、ｕは風速、Γはガンマ関数、Ｃ（ｘ）は降下粉塵量を表す。
一方、飛散物質が搬送気体の密度と同一とみなして良い場合には、３次元のガウシアンプルームモデルが提案されている。この場合、下の（５）式により、任意の空間位置（ｘ，ｙ，ｚ）における飛散物質濃度を計算することができる（非特許文献２参照）。【０００６】
2010001925.tif 000004
【０００７】
ｘは発塵源から風下方向の座標軸上（ｘ軸）上の任意の点、ｙはｘ軸と水平平面上で垂直に交わる座標軸（ｙ軸）上の任意の点、ｚはｘ軸とｙ軸で形成される水平平面に対し鉛直方向の座標軸（ｚ軸）上の任意の点、ｃ（ｘ，ｙ，ｚ）は、任意の座標点（ｘ，ｙ，ｚ）における粉塵濃度、ｕは風速、Ｑは発塵強度、σｙ、σｚは、各々、ｙ軸、ｚ軸方向の有視煙の広がり、Ｈｅは有効発塵高さ、ｙ０は発塵源の存在するｙ軸上の座標点を示す。また、発塵源の存在するｘ軸上の座標点は０としている。
【先行技術文献】
【特許文献】
【０００８】
【特許文献１】日本国特開２００５−１２８６９１号公報
【非特許文献】
【０００９】
【非特許文献１】C.H.Bosanquet et al., Proc Inst Mech Engrs, Vol.162, p355 (1950)
【非特許文献２】風の気象学，竹内清秀，東京大学出版
【非特許文献３】化学工学便覧  改訂四版
【非特許文献４】Briggs G. A. Plume rise, U.S. AEC (1969)
【発明の概要】
【発明が解決しようとする課題】
【００１０】
しかしながら、非特許文献１に開示されている推定方法では、降下粉塵量が、発塵源と降下点間の風下方向の１次元距離に依存すると仮定している。このため、図１に示すように、ひとつの降下点に対し、発塵源と降下点の距離、発塵強度、粉塵の粒度分布と粉塵の密度、及び有効発塵高さが全て等しい発塵源１と発塵源２が存在する場合は、発塵源１と発塵源２の寄与は同じと推定される。しかし、実際は、粉塵粒子は３次元方向に乱流拡散する。このため、降下点からの距離が近い発塵源２の寄与が、発塵源１より大きくなるという事実を記述できないという問題があった。また、（４）式はBosanquetが１９５０年に発表した実験式であるため、数式の精度は、実験式を導出したときの実験の環境に依存し、一般性に欠けるという問題があった。
一方、非特許文献２に開示されている推定方法では、飛散物質が搬送気体の密度と同一であるという仮定が必要であった。粉塵粒子の場合は、この仮定が満足されないため、この推定方法は、粉塵粒子の挙動を精度良く推定することができない。【００１１】
本発明は前述のような点に鑑みてなされたものであり、従来の実験式と比較して、実現象の原理原則に忠実かつ正確な推定を行うことができる方法を提供することを目的とする。
【課題を解決するための手段】
【００１２】
上記の課題を解決するために、本発明の各態様は、以下の手段を採用した。
（１）本発明は、風向及び風速の所定期間における時系列計測値を取得し、前記時系列計測値をもとに、風向及び風速の範囲を、各々ｍ、ｎ個の分割範囲に分割し、風速の前記分割範囲毎に風速代表値を設定し、それぞれの前記分割範囲に含まれる前記時系列計測値の前記所定期間での頻度を求めることによってｍ×ｎ行列の風向・風速頻度分布を作成し、前記風速代表値及び前記風向・風速頻度分布を風向・風速情報とする風向・風速情報入力工程と；粉塵の発塵源の情報を入力する発塵源情報入力工程と；前記風向・風速情報と、前記発塵源の情報と、を用いて、任意の座標点における粉塵濃度ｃを計算する粉塵濃度計算工程と；前記粉塵濃度ｃに基づいて、任意の降下地点における降下粉塵量を計算する降下粉塵量計算工程と；を備える降下粉塵量推定方法である。
（２）前記（１）に記載の降下粉塵量推定方法において、前記粉塵濃度計算工程は、前記粉塵の地表での反射率βを更に用いて粉塵濃度ｃを計算してもよい。
（３）前記（２）に記載の降下粉塵量推定方法において、前記発塵源情報入力工程では、前記発塵源の３次元空間上での座標と、前記発塵源の発塵強度Ｑと、前記粉塵の粒子径と、前記粉塵の密度と、前記発塵源からの有効発塵高さＨｅと、が入力されてもよい。
（４）前記（３）に記載の降下粉塵量推定方法において、前記粉塵濃度計算工程は、３次元座標軸上の点（ｘ，ｙ，ｚ）における粉塵濃度ｃ（ｘ，ｙ，ｚ）を（１）式を用いて計算し；
2010001925.tif 000005
ｕは風速であり、Ｋｙ、Ｋｚは、各々、ｙ軸、ｚ軸方向の乱流拡散係数であり、ｗは粒子の終末沈降速度であり、前記発塵源のｘ座標は０であり、前記発塵源のｙ座標はｙ０であってもよい。
（５）前記（２）に記載の降下粉塵量推定方法において、前記降下粉塵量計算工程は、（２）式又は（３）式を用いて降下粉塵量を計算し；
2010001925.tif 000006
ｃ（ｘ，ｙ，ｚ）は、３次元座標軸上の点（ｘ，ｙ，ｚ）における粉塵濃度であり；Ｃは降下粉塵量であり、Ｋzはｚ軸方向の乱流拡散係数であってもよい。
（６）また、本発明は、風向及び風速の所定期間における時系列計測値を取得し、前記時系列計測値をもとに、風向及び風速の範囲を、各々ｍ、ｎ個の分割範囲に分割し、風速の前記分割範囲毎に風速代表値を設定し、それぞれの前記分割範囲に含まれる前記時系列計測値の前記所定期間での頻度を求めることによってｍ×ｎ行列の風向・風速頻度分布を作成し、前記風速代表値及び前記風向・風速頻度分布を風向・風速情報とする風向・風速情報入力部と；粉塵の発塵源の情報が入力される発塵源情報入力部と；前記風向・風速情報と、前記発塵源の情報と、を用いて、任意の座標点における粉塵濃度ｃを計算する粉塵濃度計算部と；前記粉塵濃度ｃに基づいて、任意の降下地点における降下粉塵量を計算する降下粉塵量計算部と；を備える降下粉塵量推定装置である。
（７）また、本発明は、風向及び風速の所定期間における時系列計測値を取得し、前記時系列計測値をもとに、風向及び風速の範囲を、各々ｍ、ｎ個の分割範囲に分割し、風速の前記分割範囲毎に風速代表値を設定し、それぞれの前記分割範囲に含まれる前記時系列計測値の前記所定期間での頻度を求めることによってｍ×ｎ行列の風向・風速頻度分布を作成し、前記風速代表値及び前記風向・風速頻度分布を風向・風速情報とする風向・風速情報入力手順と；粉塵の発塵源の情報を入力する発塵源情報入力手順と；前記風向・風速情報と、前記発塵源の情報と、を用いて、任意の座標点における粉塵濃度ｃを計算する粉塵濃度計算手順と；前記粉塵濃度ｃに基づいて、任意の降下地点における降下粉塵量を計算する降下粉塵量計算手順と；をコンピューターに実行させる降下粉塵量推定プログラムである。
（８）また、本発明は、上記（７）に記載のプログラムを記録したコンピューター読み取り可能な記録媒体である。
【発明の効果】
【００１３】
上記（１）の発明によれば、従来の実験式と比較して、実現象の原理原則に忠実かつ正確に粉塵挙動の推定ができる。これにより、発塵源から発生した粉塵が、降下粉塵量として市街地にどれだけ影響を及ぼすかの定量的な評価ができる。これとともに、集塵機等の発塵抑制対策設備の最適な規模を見積もる際の設備設計指標が得られる。特に、領域内の風向・風速が一律でない場合でも、実測した風向・風速に従って正確に粉塵挙動の推定ができる。また、風向・風速を特定の数の分割範囲に分割して計算を行なうため、計算量の節約ができる。このため、複雑な系についても計算が実現可能であり、有意な推定結果を得られる。
上記（２）の発明によれば、粉塵の粒径等の条件に拠って、前記粉塵の地表での反射及び沈着の効果が無視できない場合でも、正確に粉塵挙動の推定ができる。
上記（３）の発明によれば、発塵源の詳細な特性が測定できる場合、この特性に従って、正確に粉塵挙動の推定ができる。
上記（４）の発明によれば、３次元的な乱流拡散係モデルを用い、測定した風速、乱流拡散係数、粉塵粒子の終末沈降速度などの特性に従って、正確に粉塵挙動の推定ができる。
上記（５）の発明によれば、粉塵の粒径が大きい場合や、粉塵の密度が高い場合など、重力による沈降が無視できない場合でも、正確に粉塵挙動の推定ができる。
【図面の簡単な説明】
【００１４】
【図１】降下点と２つの発塵源の位置関係を説明するための図である。
【図２】発塵源からの粉塵飛散挙動を説明するための図である。
【図３】粉塵濃度および降下煤塵量を計算するための第１のフローチャートである。
【図４】粉塵濃度および降下煤塵量を計算するための第２のフローチャートである。
【図５Ａ】Pasquill-Giffordによる煙の広がりを決定するための図である。
【図５Ｂ】Pasquill-Giffordによる煙の広がりを決定するための図である。
【図６】本発明の実施形態に係る方法で計算した粉塵飛散の濃度分布の一例を示す図である。
【図７】ガウシアンプルームモデルで計算した粉塵飛散の濃度分布の一例を示す図である。
【図８】本発明の実施形態に係る方法で計算した降下煤塵量分布の一例を示す図である。
【図９】非特許文献１のBosanquet式で計算した降下煤塵濃度分布の一例を示す図である。
【図１０】降下煤塵量のデポジットゲージによる実測値の一例を示す図である。
【図１１】本発明の実施形態に係る装置構成の概要を示す図である。
【図１２】粉塵濃度および降下煤塵量推定装置として機能するコンピュータシステムのハードウエア構成例を示す図である。
【発明を実施するための形態】
【００１５】
本発明は、発塵源から降下点への風による３次元方向の移流・拡散による飛散挙動を、粉塵の物質収支にもとづき、理論的に求め、降下粉塵量を推定する。
以下、添付図面を参照して、本発明の好適な実施形態について説明する。
本実施形態では、図２に示すように、発塵源の一例として、工場の煙突からの粉塵粒子の飛散を考える。粉塵粒子は、Ｈｅの有効発塵高さから風下方向に設定したｘ軸を中心軸として飛散する。発塵源はｘ＝０の位置に存在するとして、鉛直方向をｚ軸に、ｘ軸と交差する水平面の軸をｙ軸とする。
【００１６】
発塵源から風によって飛散する粉塵粒子の濃度をｃ（ｘ，ｙ，ｚ）とする。濃度ｃ（ｘ，ｙ，ｚ）について、物質収支は、下の（６）式で記述できる。
【００１７】
2010001925.tif 000007
【００１８】
（６）式において、ｕは風速を表し、上述したようにｘ軸は風向の中心に設定し、風はｘ軸方向のみに吹くと仮定する。
（６）式の左辺第１項は、粉塵粒子が風の移流によって飛散する量を示す。（６）式の右辺第１項は、粉塵粒子が乱流拡散によってｙ軸方向に飛散する量を示す。（６）式の右辺第２項は、粉塵粒子が乱流拡散によってｚ軸方向に飛散する量を示す。（６）式の右辺第３項は、粉塵粒子が重力によって沈降する量を示す。ここで、粉塵粒子が乱流拡散でｘ軸方向に飛散する量は、風の移流によりｘ軸方向に飛散する量と比較して小さいため、無視できると仮定した。Ｋｙ、Ｋｚは、各々、ｙ軸方向とｚ軸方向の乱流拡散係数を表す。
【００１９】
下の（７）式は、発塵源における粉塵粒子の境界条件の設定に関する式である。Ｑ・δ（ｙ−ｙ０，ｚ−Ｈｅ）はｘ＝０、ｙ＝ｙｏ、ｚ＝Ｈｅの位置から発塵強度Ｑで粉塵粒子が発生していることを示す。δはデルタ関数である。ここで、デルタ関数とは、（８）式の条件を満たす実超関数である。ここで、発塵強度Ｑは、ローボリュームサンプラー等の手段により求めることができる。
【００２０】
2010001925.tif 000008
【００２１】
2010001925.tif 000009
【００２２】
下の（９）式は、地表における粉塵粒子の境界条件の設定に関する式である。βは粉塵粒子の地表での反射率を示し、β＝０で完全沈着、β＝１で完全反射を示す。βは粉塵粒子の大きさに依存する係数である。通常の計算では、地表面に到達し沈着する粉塵粒子が計算対象になることが殆どであるため、β＝０で近似しても問題ない。正確にβの値を設定したい場合は、粉塵の発生速度、有効発塵高さおよび粉塵の粒子径分布が既知の煙突等の発塵源からの降下煤塵量を、デポジットゲージ等の手段で計測する。そして、各粒子径毎に、降下煤塵の計算値が降下煤塵の実測値に対応するように、βに適切な値を設定する。
【００２３】
2010001925.tif 000010
【００２４】
（６）式を、（７）式と（９）式の境界条件のもとに解析的に解くと、（１）式が導出できる。
【００２５】
2010001925.tif 000011
【００２６】
降下点における降下煤塵量Ｃ（ｘ，ｙ）は、地表ｚ＝０における任意の平面座標点（ｘ，ｙ）での粉塵濃度ｃ（ｘ，ｙ，ｚ＝０）の（９）式で算出した計算値に、粉塵粒子の終末沈降速度ｗを乗じた以下の（２）式により計算できる。
【００２７】
2010001925.tif 000012
【００２８】
空気中で粉塵が重力落下するとき、粉塵の重力と浮力がバランスする結果、時間がたつと速度が一定になる。これを終末沈降速度といい、（１０）式で計算できる。
【００２９】
2010001925.tif 000013
ここで、ρＳは粉塵粒子の密度、ρａは大気の密度、ｄｋは流刑頻度ｋに相当する粉塵粒子の直径、μは大気の粘性係数を示す。
【００３０】
また、降下点における降下煤塵量Ｃ（ｘ，ｙ）は、地表ｚ＝０における任意の平面座標点（ｘ，ｙ）での粉塵濃度ｃ（ｘ，ｙ，ｚ＝０）のｚ微分値を（９）式をもとに計算し、Ｋｚと（１−β）を乗じた（３）式でも計算できる。
【００３１】
2010001925.tif 000014
ここでＫｙはｚ軸方向の乱流拡散係数、βは粒子の地表での反射率を示す。
【００３２】
以下、（１）式及び（２）式を使用して、降下煤塵量を計算する手順の一例について、図３および図４のフローチャートを参照しながら説明する。  先ず、風速、風力情報の入力ステップについて説明する。最初のステップＳ０１において、風向、風速の頻度分布を算出する。一定期間の風向、風速の計測値をもとに、風向、風速データを風向についてｍ個、風速についてｎ個のカテゴリーに分類する。風速ｉ、風向ｊのカテゴリに属するデータの頻度ｂｉｊを算出し、ｎ行ｍ列の行列を作成する。このとき、この行列の各要素ｂｉｊの総和が１００％になるように調節する。
通常、風向データの分類は、１方位以上の何方位でも計算でき、入手できる方位データに鑑みて設定すればよい。好ましくは気象庁のアメダスデータに対応させ、１６方位を採用し、ｍ＝１６とする。
風速データの分類は、風速の実測値を直接使用して計算できる。好ましくは、設定期間の最小風速と最大風速との間の区間を、（微風、弱風、通常風、強風）、又は、（微風、弱風、通常風、強風、超強風）等の分類に４または５分割し、ｎ＝４またはｎ＝５とすればよい。ｎが６以上でも計算できるが、場合によっては計算が煩雑になる。分割した風速の範囲のそれぞれにおいて、風速の代表値を選定する。代表値の選定に当たっては、各分割区間に属する風速データのうち、発生頻度の最も大きい風速値を採用してもよいし、分割区間に属する風速データの平均値を採用してもよい。
また、風向計と風速計は、周囲にある建物等の障害物の影響を受けない位置に設置することが好ましい。気象庁が近傍にある場合は、アメダス等の気象庁観測データを使用してもよい。
また、風向と風速を計測する期間の設定は、降下煤塵量の計算期間と同一に設定する。そして、設定した期間において、風速、風向の頻度分布の統計量を決定できるように、風速、風向のサンプリング周期を決定する。例えば、ｍ＝１６、ｎ＝４の場合において、１ヶ月間の期間を設定する場合は、１時間周期の風速、風向データを採取すれば十分である。
【００３３】
ステップＳ０２において、降下煤塵の初期値をゼロに設定する。
【００３４】
次に、発塵情報の入力ステップについて説明する。ステップＳ０３において、発塵源の情報を入力する。具体的には、発塵源のｘ−ｙ座標（ｘ＝０，ｙ＝ｙ０），煙突高さＨ，発塵強度Ｑ，粉塵粒子径ｄ、粉塵粒子密度ρ、排ガス風量Ｗ，排ガス風速Ｖ等を入力する。
【００３５】
次に、降下地点情報の入力ステップについて説明する。ステップＳ０４において、降下地点の情報を入力する。具体的には、降下煤塵量を推定する降下地点の座標（ｘ，ｙ）を入力する。
なお、上述した、風向、風速情報の入力ステップ、発塵情報の入力ステップ、降下地点情報の入力ステップは、入力の順番は問わず、入れ替わっても構わない。
【００３６】
次に、図４に示す粉塵濃度計算ステップＳ１００に入る。本ステップでは、上記各ステップで入力された情報を用いて、降下点における粉塵濃度を計算する。
【００３７】
ステップＳ０５において、風向・風速の各頻度分布ｂｉｊにおける風速ｕｉ（代表風速値）を、以降のステップでの計算に使用する風速ｕと定義する。風速ｕは、後述する（１１）式を使った有効発塵高さの計算と、（１）式を使った降下点における粉塵濃度計算に使用される。
ステップＳ０６において、粉塵が工場煙突等、周囲の空気より高温で排出される場合は、ステップＳ０３で入力した煙突高さＨ、排ガス風量Ｗ、排ガス風速Ｖをもとに、有効発塵高さＨｅを、例えば（１１）式により計算する。
【００３８】
2010001925.tif 000015
ここで、Ｈは煙突高さ、Ｑｅ＿ｇａｓは排ガス流量、Ｔ１は大気温度、ΔＴは排ガス温度とＴ１との差、ｇは重力加速度、ｄθ／ｄｚは大気の温度勾配（通常０．０３を用いる）、ｗは粒子の終末沈降速度、ｕは風速を表す。
【００３９】
一方、粉塵が、風により飛散する場合は、カメラ等の装置で発塵状況を記録し、粉塵が最大濃度である高さを色の濃さから識別判定し、発塵高さを決定する。この発塵高さを、有効発塵高さＨｅとして設定する。
【００４０】
粉塵粒子径については、ローボリュームサンプラーで採取した粉塵粒子サンプルを、粒子径分布測定装置に投入し、粒度分布を計測する。計測した粒度分布をもとに、最大頻度径等の定義を使い、代表粒子径を決定する。
また、粒子径の計測範囲を分割し、分割した各々の範囲における代表粒子径を、粉塵粒子径として用いてもよい。この場合は、分割した粒子径の範囲に存在する粉塵粒子の重量比をあらかじめ求めておく。そして、重量の総和が１００％になるように、各平均粒子径毎に頻度分布を計算する。そして、各平均粒子径毎の降下煤塵量に、前記の各平均粒子径毎に計算した頻度分布を乗じて、降下煤塵量の総量を求める。粒子径分布測定装置は、沈降法、光透過式等の方法がある。（例えば、非特許文献４、参照）
ステップＳ０７において、ステップＳ０３で入力した粉塵粒子の粒子径ｄ、密度ρをもとに、終末沈降速度ｗを、前述の（１０）式により計算する。
【００４１】
ステップＳ０８において、乱流拡散係数Ｋｙ、Ｋｚを決定する。ＫｙおよびＫｚは、図５Ａ，５Ｂに示すPasquill-Giffordが米国の草原でトレーサー実験を行い観測した煙の広がりσｙ、σｚの実験値をもとに、下記の（１２）式に基づき換算した値を設定する。（例えば、非特許文献４を参照）
【００４２】
2010001925.tif 000016
【００４３】
ここで、図５Ａ，５Ｂの記号Ａ、Ｂ、Ｃ、Ｄ、Ｅ、Ｆは大気安定度を表し、Ａが非常に不安定、Ｂが不安定、Ｃがやや不安定、Ｄが中立、Ｅが安定、Ｆが非常に安定な状態に相当する。大気安定度は、乱流拡散係数Ｋｙ、Ｋｚの大きさに関係し、図５Ａ，５Ｂから大気安定度が安定になるにつれ、乱流拡散係数Ｋｙ、Ｋｚは小さくなる。ステップＳ０９において、反射率βを入力する。
【００４４】
ステップＳ０１０において、前述の（１）式により、ｚ＝０における降下点の座標位置（ｘ，ｙ）における粉塵粒子の濃度ｃ（ｘ，ｙ，ｚ＝０）を計算する。
【００４５】
次に、降下煤塵量計算ステップＳ２００に入る。ステップＳ０１１において、ステップＳ０１０で計算した粉塵粒子の濃度ｃ（ｘ，ｙ，ｚ＝０）とステップＳ０７で計算した粉塵粒子の終末沈降速度ｗを使って、前述の（２）式により、降下煤塵量を計算する。この降下煤塵量を、ΔＣ（ｘ，ｙ）と定義する。
【００４６】
ステップＳ０１２では、降下煤塵量Ｃ（ｘ，ｙ）に、ステップＳ０１１で計算した各風向・風速の頻度分布ｂｉｊに相当する降下煤塵量の計算値ΔＣ（ｘ，ｙ）を加算する。
【００４７】
ステップＳ０１３では、ｉ＜ｎ、ｊ＜ｍの条件を満たしていれば、ステップＳ０１４のループを実行し、ステップＳ０５からＳ０１２の計算を繰り返すことで、降下煤塵量Ｃ（ｘ，ｙ）を求めることができる。ステップ０１４のループは、まず、ｊの値を１つずつ更新し、ｊ＝ｍに到達するまで、ステップＳ０５からステップＳ０１２の計算を繰りかえし、ｊ＝ｍになれば、ｉの値を同様に、１つずつ更新し、ｉ＝ｎに到達するまで、ステップＳ０５からステップＳ０１２の計算を繰りかえすことを表す。
【００４８】
Ｓ０１３において、ｉ＝ｎおよびｊ＝ｍになれば、Ｅｎｄに移行し、計算を終了する。
【００４９】
上述の手法は、粉塵が風により飛散する際の、粉塵の挙動を、物理現象の原理原則に沿うように記述したものであり、煤塵などの工場から排気される粉塵、自動車の排気ガス中に含まれる粉塵、花粉症の原因となる花粉、黄砂、および砂漠の砂等の飛散現象に適用できる。
上述の手法によって得られた粉塵の挙動の推定結果に基づき、例えば、煙突等の発塵源の数や設置位置、発塵強度等を適切に設計した上で発塵源を設置することができる。また、集塵機等の発塵抑制対策設備を最適に配置することができる。これによって、粉塵の飛散による環境への影響等を適切に制御することができる。
【実施例】
【００５０】
本発明による粉塵濃度の計算結果を図６に示す。計算の前提条件は、有効発塵高さＨｅ＝５０ｍ、粒子径６０μｍの粉塵粒子の終末沈降速度ｗ＝０．１ｍ／ｓｅｃ、発塵強度Ｑ＝１０ｇ／ｓ、風速ｕ＝６ｍ／ｓｅｃ、反射率β＝０．０であり、ｙ＝０のｘ軸−ｚ軸断面の粉塵粒子の濃度分布を示す。解析対象となる粉塵の粒子径は、３０μｍから９０μｍの範囲に分布しており、その平均値である６０μｍを採用した。ｘを発塵源から風下方向の距離（ｍ）とすると、乱流拡散係数Ｋｙ，Ｋｚはそれぞれ、１００４．７／ｘ（ｍ２／ｓｅｃ），３９０９．６／ｘ（ｍ２／ｓｅｃ）である。【００５１】
比較例として、図７に、ガウシアンプルームモデルとして知られている（１３）式により計算した、上記と同一条件での計算結果を示す。（１３）式では粒子の重力沈降の効果が考慮されていないため、図６の場合と異なり、粉塵濃度のコンター線の下方減衰効果が観察できない。また、粒子は完全反射することを仮定しているので、地表における粉塵粒子の濃度勾配が図６の結果と異なっているのが判る。【００５２】
2010001925.tif 000017
【００５３】
図８から図１０は、地表における降下煤塵量分布のコンター図であり、それぞれ、本実施例の方法、Bosanquet式の方法、デポジットゲージによる実測値に基づいて描画されている。縦軸はｘ軸、横軸はｙ軸である。図中の濃度が高いほど降下煤塵量が多いことを示す。
図８は、本手法を使用して計算した降下煤塵量の計算結果を示す。風向は１６方位を採用し、ｍ＝１６とした。風速は、１年間の最小風速０．１ｍ／ｓ、最大風速９ｍ／ｓの間を、微風、弱風、通常風、強風の分類に対応させるべく４分割し、ｎ＝４とした。具体的には、当該領域での過去の風速の発生頻度をもとに、各々、０．１〜１．５ｍ／ｓ、１．５〜３ｍ／ｓ、３〜６ｍ／ｓ、６〜９ｍ／ｓの区間に分割した。各区間の代表風速として、１ｍ／ｓ、２ｍ／ｓ、４．５ｍ／ｓ、７ｍ／ｓを計算に採用した。
風向、風速データは、気象庁のアメダスデータを使用し、観測データ間隔１時間の１ヶ月間のデータから風向、風速の頻度分布を算出した。表１に、計算に使用した風向、風速の頻度分布を示す。【００５４】
【表１】
2010001925.tif 000018
【００５５】
対象とする粒子径が３０μｍ以上と比較的大きいため、粒子の地面での反射率としては、粒子の完全沈着条件であるβ＝０を採用した。
図９は、比較例として、文献１にある（１）式で計算した降下煤塵量の計算結果を示す。
図１０は、降下煤塵量の実測値をもとに描いたコンター図である。当該領域をほぼ均等に８０分割し、各々の分割領域の中心に１台のデポジットゲージを設置した。総計８０台のデポジットゲージに１０日間堆積した煤塵の重量を計測して、降下煤塵量を求めた。図８は、図９と比較し、実測値である図１０とコンター形状がより良く一致しているのが判る。【００５６】
図１１に、本発明による粉塵濃度および降下煤塵量推定装置１００の概要を示す。装置１００は、粉塵濃度計算手段１４と降下煤塵量計算手段１５とを備え、所定の降下地点において地表に降下し沈着する粉塵量を降下煤塵量として推定する。粉塵濃度計算手段１４には、風向・風速情報入力手段１１と、発塵源情報入力手段１２と、降下地点情報入力手段１３とからの情報が入力され、さらに粉塵地表反射率βが入力される。【００５７】
風向・風速情報入力手段１１は、風向及び風速の所定期間における時系列計測値をもとに、風向と風速の範囲を、各々ｍ、ｎ個の分割範囲に分割する。また、分割したｎ個の風速の分割範囲毎に代表値を設定する。更に、それぞれの風速・風速の分割範囲に含まれる時系列計測値の数が、前記所定期間の全計測値数に占める頻度を求め、ｍ×ｎの頻度分布行列を設定する。さらに、求めた風向、風速の頻度分布の値を風向、風速情報として、粉塵濃度計算手段１４に入力する。【００５８】
発塵源情報入力手段１２に、発塵源の３次元空間上での座標を指定し、発塵源の発塵強度（粉塵の発生速度）、前記粉塵の粒子径、及び、前記粉塵の密度を入力する。また、計算または実測により算出した前記発塵源からの有効発塵高さを粉塵濃度計算手段１４に入力する。【００５９】
降下地点情報入力手段１３は、地表における粉塵濃度を評価するための降下地点の２次元空間における座標を指定し、この座標を粉塵濃度計算手段１４に入力する。【００６０】
粉塵濃度計算手段１４は、風向・風速情報入力手段１１、発塵源情報入力手段１２、および、降下地点情報入力手段１３により入力された全入力情報と、更に入力される粉塵の地表での反射率β等の入力値をもとに、発塵源から風によって周囲に移流および拡散される粉塵について、３次元空間における任意の座標点（ｘ，ｙ，ｚ）における粉塵濃度ｃ（ｘ，ｙ，ｚ）を、（１）式により計算する。【００６１】
降下煤塵量計算手段１５では、前記粉塵濃度計算手段１４において計算した任意の座標点のうちｚ＝０の降下地点において地表に降下し、沈着する粉塵量を計算する。このとき、Ｃを降下煤塵量として（２）式を用いてもよい。又は、Ｋｙをｚ軸方向の乱流拡散係数とし、βを粒子の地表での反射率として（３）式を用いてもよい。【００６２】
図１２には、粉塵濃度および降下煤塵量推定装置１００として機能するコンピュータシステムのハードウェア構成例を示す。降下煤塵量推定装置１００は、ＣＰＵ２０と、入力装置２１と、表示装置２２と、記録装置２３とを含み、各部はバス２４を介して接続される。【００６３】
入力装置２１には、降下煤塵量の推定に必要な風向・風速情報、発塵情報、降下地点情報、反射率βが入力される。ＣＰＵ２０は、入力装置２１に入力された情報をもとに、有効発塵高さＨｅ、粉塵粒子の終末沈降速度、乱流拡散係数Ｋｙ，Ｋｚを決定し、降下地点における降下煤塵量を計算する。【００６４】
表示装置２２は、ＣＰＵ２０で計算した複数の降下地点における降下煤塵量の計算値をもとに、降下煤塵量のコンター図を表示する。【００６５】
記録装置２３は、入力装置２１に入力された全情報およびＣＰＵ２０で計算した各降下地点における降下煤塵量を記録する。記録装置２３はＲＯＭ、ＲＡＭ、ＨＤ等により構成されている。記録装置２３には、前述した降下煤塵量推定装置１００としての動作を制御するコンピュータプログラムも格納される。【００６６】
ＣＰＵ２０がこのコンピュータプログラムを実行することによって降下煤塵量推定装置１００の機能、または処理を実現する。また、記録装置２３にデータベースが格納される。【００６７】
なお、本発明の降下煤塵量推定装置は、複数の機器から構成されるシステムに適用しても、一つの機器からなる装置に適用してもよい。【００６８】
また、本発明の目的は、前述した機能を実現するコンピュータプログラムをシステム或いは装置に供給し、そのシステム或いは装置のコンピュータ（ＣＰＵ若しくはＭＰＵ）が実行することによっても達成され、この場合、コンピュータプログラム自体が本発明を構成することになる。以上、本発明を種々の実施形態と共に説明したが、本発明はこれらの実施形態にのみ限定されるものではなく、本発明の範囲内で変更等が可能である。
【産業上の利用可能性】
【００６９】
本発明によれば、従来の実験式と比較して、実現象の原理原則に忠実かつ正確に粉塵挙動の推定ができる。これにより、発塵源から発生した粉塵が、降下煤塵量として市街地にどれだけ影響を及ぼすかの定量的な評価ができる。これとともに、集塵機等の発塵抑制対策設備の最適な規模を見積もる際の設備設計指標が得られる。
【符号の説明】
【００７０】
２０    ＣＰＵ
２１    入力装置
２２    表示装置
２３    記録装置
１００    降下煤塵量推定装置
【図１】
2010001925.tif 000022
【図２】
2010001925.tif 000023
【図３】
2010001925.tif 000024
【図４】
2010001925.tif 000025
【図５Ａ】
2010001925.tif 000026
【図５Ｂ】
2010001925.tif 000027
【図１１】
2010001925.tif 000028
【図１２】
2010001925.tif 000029
【図６】
2010001925.tif 000030
【図７】
2010001925.tif 000031
【図８】
2010001925.tif 000032
【図９】
2010001925.tif 000033
【図１０】
2010001925.tif 000034
【手続補正書】
【提出日】平成21年12月1日(2009.12.1)
【手続補正１】
【補正対象書類名】特許請求の範囲
【補正対象項目名】全文
【補正方法】変更
【補正の内容】
【特許請求の範囲】
【請求項１】
風向及び風速の所定期間における時系列計測値を取得し、前記時系列計測値をもとに、風向及び風速の範囲を、各々ｍ、ｎ個の分割範囲に分割し、風速の前記分割範囲毎に風速代表値を設定し、それぞれの前記分割範囲に含まれる前記時系列計測値の前記所定期間での頻度を求めることによってｍ×ｎ行列の風向・風速頻度分布を作成し、前記風速代表値及び前記風向・風速頻度分布を風向・風速情報とする風向・風速情報入力工程と；
粉塵の発塵源の情報を入力する発塵源情報入力工程と；
前記風向・風速情報と、前記発塵源の情報と、を用いて、任意の座標点における粉塵濃度ｃを計算する粉塵濃度計算工程と；
前記粉塵濃度ｃに基づいて、任意の降下地点における降下粉塵量を計算する降下粉塵量計算工程と；
を備え、
前記粉塵濃度計算工程は、前記粉塵の地表での反射率βを更に用いて前記粉塵濃度ｃを計算する、
ことを特徴とする降下粉塵量推定方法。
【請求項２】
（削除）
【請求項３】
前記発塵源情報入力工程では、前記発塵源の３次元空間上での座標と、前記発塵源の発塵強度Ｑと、前記粉塵の粒子径と、前記粉塵の密度と、前記発塵源からの有効発塵高さＨｅと、が入力されることを特徴とする請求項１に記載の降下粉塵量推定方法。
【請求項４】
前記粉塵濃度計算工程は、３次元座標軸上の点（ｘ，ｙ，ｚ）における粉塵濃度ｃ（ｘ，ｙ，ｚ）を（１）式を用いて計算し；
2010001925.tif 000035
ｕは風速であり、Ｋｙ、Ｋｚは、各々、ｙ軸、ｚ軸方向の乱流拡散係数であり、ｗは粒子の終末沈降速度であり、前記発塵源のｘ座標は０であり、前記発塵源のｙ座標はｙ０である；
ことを特徴とする請求項３に記載の降下粉塵量推定方法。
【請求項５】
前記降下粉塵量計算工程は、（２）式又は（３）式を用いて降下粉塵量を計算し；
2010001925.tif 000036
ｃ（ｘ，ｙ，ｚ）は、３次元座標軸上の点（ｘ，ｙ，ｚ）における粉塵濃度であり；
Ｃは降下粉塵量であり、
Ｋｚはｚ軸方向の乱流拡散係数である；
ことを特徴とする請求項１に記載の降下粉塵量推定方法。
【請求項６】
風向及び風速の所定期間における時系列計測値を取得し、前記時系列計測値をもとに、風向及び風速の範囲を、各々ｍ、ｎ個の分割範囲に分割し、風速の前記分割範囲毎に風速代表値を設定し、それぞれの前記分割範囲に含まれる前記時系列計測値の前記所定期間での頻度を求めることによってｍ×ｎ行列の風向・風速頻度分布を作成し、前記風速代表値及び前記風向・風速頻度分布を風向・風速情報とする風向・風速情報入力部と；
粉塵の発塵源の情報が入力される発塵源情報入力部と；
前記風向・風速情報と、前記発塵源の情報と、を用いて、任意の座標点における粉塵濃度ｃを計算する粉塵濃度計算部と；
前記粉塵濃度ｃに基づいて、任意の降下地点における降下粉塵量を計算する降下粉塵量計算部と；
を備え、
前記粉塵濃度計算部は、前記粉塵の地表での反射率βを更に用いて前記粉塵濃度ｃを計算する、
ことを特徴とする降下粉塵量推定装置。
【請求項７】
風向及び風速の所定期間における時系列計測値を取得し、前記時系列計測値をもとに、風向及び風速の範囲を、各々ｍ、ｎ個の分割範囲に分割し、風速の前記分割範囲毎に風速代表値を設定し、それぞれの前記分割範囲に含まれる前記時系列計測値の前記所定期間での頻度を求めることによってｍ×ｎ行列の風向・風速頻度分布を作成し、前記風速代表値及び前記風向・風速頻度分布を風向・風速情報とする風向・風速情報入力手順と；
粉塵の発塵源の情報を入力する発塵源情報入力手順と；
前記風向・風速情報と、前記発塵源の情報と、を用いて、任意の座標点における粉塵濃度ｃを計算する粉塵濃度計算手順と；
前記粉塵濃度ｃに基づいて、任意の降下地点における降下粉塵量を計算する降下粉塵量計算手順と；
を備え、
前記粉塵濃度計算手順は、前記粉塵の地表での反射率βを更に用いて前記粉塵濃度ｃを計算する、
ことを特徴とする降下粉塵量推定プログラム。
【請求項８】
請求項７に記載の降下粉塵量推定プログラムを記録したコンピューター読み取り可能な記録媒体。
【国際調査報告】
2010001925.tif 000037
2010001925.tif 000038
2010001925.tif 000039