水中シルトフェンスの変形・移動を可能とする
3
次
元モデリング：流体・シルトフェンス連成計算の
実現
山田 進
∗,∗∗町田 昌彦
∗.∗∗田中 みのり
†関 克己
†有川 太郎
† ∗日本原子力研究開発機構 システム計算科学センター
∗∗日本原子力研究開発機構 安全研究センター
†中央大学 理工学部
概要. シルトフェンスとは，河川や沿岸等で土木工事により発生する汚濁物質の拡散抑 制のため使用されるカーテン状の敷設物であり，一定の流れ場の下では鉛直および水平方 向の変形が発生する．本論文では，流体の運動と連成し，その２つの変形を同時に再現可 能なアルゴリズムを提案する．実際に水路にシルトフェンスを敷設した場合のシミュレー ションを行い，同じ条件での実験結果との比較から，提案アルゴリズムの妥当性を確認 する．Three-dimensional Modeling for Deformation and Movement
of Siltfence in Water: Coupled Simulation of Water Flow and
Siltfence
Susumu Yamada
∗,∗∗Masahiko Machida
∗,∗∗Minori Tanaka
†Katsumi Seki
†Taro Arikawa
†∗
Center for Computational Science & e-Systems, Japan Atomic Energy Agency
∗∗Nuclear Safety Research Center, Japan Atomic Energy Agency
†
Faculty of Science and Engineering, Chuo University
Abstract. Siltfence is a curtain-like structure to prevent the spread of suspended pol-lutants generated by constructions in rivers or coasts, and its shape deforms horizontally and vertically. In this paper, we propose a modeling method to simulate the siltfence de-formed in the two directions by coupling with water flow dynamics. We preform numerical simulations on a siltfence installed in flume using the proposed method, and we confirm its validity by comparing between the simulation and flume experiment.
1.
はじめに
河川及び湖沼や海洋沿岸等では様々な目的で土木工事が行われる．その際，発生する汚 濁物質の拡散を抑制し環境への影響を低減させるため，汚濁防止膜（シルトフェンス）が
水中シルトフェンスの変形・移動を可能とする
3
次
元モデリング：流体・シルトフェンス連成計算の
実現
山田 進
∗,∗∗町田 昌彦
∗.∗∗田中 みのり
†関 克己
†有川 太郎
† ∗日本原子力研究開発機構 システム計算科学センター
∗∗日本原子力研究開発機構 安全研究センター
†中央大学 理工学部
概要. シルトフェンスとは，河川や沿岸等で土木工事により発生する汚濁物質の拡散抑 制のため使用されるカーテン状の敷設物であり，一定の流れ場の下では鉛直および水平方 向の変形が発生する．本論文では，流体の運動と連成し，その２つの変形を同時に再現可 能なアルゴリズムを提案する．実際に水路にシルトフェンスを敷設した場合のシミュレー ションを行い，同じ条件での実験結果との比較から，提案アルゴリズムの妥当性を確認 する．Three-dimensional Modeling for Deformation and Movement
of Siltfence in Water: Coupled Simulation of Water Flow and
Siltfence
Susumu Yamada
∗,∗∗Masahiko Machida
∗,∗∗Minori Tanaka
†Katsumi Seki
†Taro Arikawa
†∗
Center for Computational Science & e-Systems, Japan Atomic Energy Agency
∗∗Nuclear Safety Research Center, Japan Atomic Energy Agency
†
Faculty of Science and Engineering, Chuo University
Abstract. Siltfence is a curtain-like structure to prevent the spread of suspended pol-lutants generated by constructions in rivers or coasts, and its shape deforms horizontally and vertically. In this paper, we propose a modeling method to simulate the siltfence de-formed in the two directions by coupling with water flow dynamics. We preform numerical simulations on a siltfence installed in flume using the proposed method, and we confirm its validity by comparing between the simulation and flume experiment.
1.
はじめに
河川及び湖沼や海洋沿岸等では様々な目的で土木工事が行われる．その際，発生する汚 濁物質の拡散を抑制し環境への影響を低減させるため，汚濁防止膜（シルトフェンス）が 広く用いられてきた[6]．一例として，関西空港建設ではシルトフェンスが大規模に敷設 された他，最近では沖縄・辺野古での土木工事での利用等が挙げられる．本論文では，そ のシルトフェンスが一定の流れ場の下，変形・移動する様子を再現するアルゴリズムを提 案する． シルトフェンスとは膜体（カーテン）状の仕切りであり，その透水性はほとんど無視で きる．したがって，上記のような沿岸等での土木工事の際，Fig. 1のように設置すること で，発生する微細な土粒子などの拡散・流出を防ぐ効果があると考えられている． また，東京電力ホールディングス(株)・福島第一原子力発電所港湾内には，放射性物質 拡散の抑制を目的として設置され[22]，放射性物質の拡散抑制の効果が実際に見られると の報告があった [9]．その後，電力事業者らはこの報告を基に，原子力発電所での事故に よる放射性物質の海洋拡散を抑制するため，シルトフェンスの設置を検討している．原子 力発電所では，幾重にも放射性物質の環境放出を防護するための対策が取られているが， 万が一，原子炉外への揮発性放射性物質の放出が始まった際は，Fig. 2のように放水砲に よる散水によって放射性物質を降下させ，港湾内に設置するシルトフェンスで降下水中の 放射性物質の海洋拡散を抑制するという方策が環境拡散抑制手段の一つとして検討されて いる[4, 5, 9, 10, 17, 21]． 以上，シルトフェンスは様々な状況に応じてその活用が想定されるが，その根拠となる 汚濁物質拡散抑制効果を明確に評価することは難しい．構造上，発生した汚濁水を完全に 閉鎖するものではなく，敷設する現場による違いも大きいことから，その効果を評価する に当たっては科学的評価手法の開発が必要である．本論文では，その科学的評価手法開発 SiltfenceFig. 1. A schematic figure of siltfence installed at a construction site in coastal area. The siltfence can trap a part of the suspended sediments generated by the construction. Water cannon Drain Primary containment vessel Siltfence
Fig. 2. Schematic figure of drop-ping of radioactive materials with a water cannon and suppressing their diffusion by a siltfence.
Float fixed on the surface
Weight Siltfence
Siltfence Float
Fixed on the bottom
Flow direction
Fig. 3. Schematic figure of typical two types of siltfences. The left hand type is called as the standing type, and the right hand one the hanging type. Each dash line indicates the initial position of each type of the siltfence under a constant flow. Here, dcis the effective
siltfence height. を目的として，シルトフェンスの挙動をシミュレーションする方法を提案する． 一般に，シルトフェンスの設置方法には，Fig. 3のように自立型と垂下型があり，自立 型は潜水工事等が必須であるため，敷設が容易な垂下型が用いられることが多い．従っ て，本論文では垂下型を対象にする．一般に，垂下型では，上端に浮体，下端に重錘を設 置し，上端の両端を固定する．上両端のみを固定する際，膜の長さを固定端長（固定位置 を直線で結んだ長さ）より長くして敷設すると（通常の敷設では，固定端長を19.5とす ると膜長は20.0とするのが一般的である[6]），一定の流れがシルトフェンスに垂直方向 に発生している状況下では，Fig. 4のように変形することが知られている．ここで，鉛直 方向の変形は「ふかれ」，水平方向の変形は「たわみ」と呼ぶ． シルトフェンス設置の目的である汚濁物質の拡散抑制効果が上記の変形に依存すること は明らかであり，変形を予測できるモデリングが重要な研究開発課題であることが分か る．これまで，問題を簡単化するため，上端を直線で固定することで「たわみ」が発生し ないモデルに対して，一定の流れ場の下，「ふかれ」の度合いのみを調べる研究が行われ 多くの研究成果が得られてきた．例えば，小田らは，「ふかれ」たシルトフェンスの形状 を曲線と見做し，各点にかかる力を解析することで抗力と流速の関係式を導出し，関係式 に含まれる実験定数を水槽実験等から求めることで，おおよそ普遍的な「ふかれ」と流速 の関係を評価することに成功している[13]．最近では，武田らが実験定数を含め，水槽実 験結果から相関式を求める等の発展がある[20]． しかし，実際の利用では水路等にて上端を直線的に固定するような利用例はなく，シル トフェンスの「ふかれ」と「たわみ」の両者が同時に発生し，その形状はFig. 4のように 3次元的に変形し複雑な様相を示す．本論文の目的は，ある一定の流れ場の下，両者の変 形自由度を同時に取り扱い，その変形を予測可能なシミュレーション方法を提案すること にある． これまでに，シルトフェンスの汚濁物質拡散抑制効果を調べるため，計算機シミュレー
Float fixed on the surface
Weight Siltfence
Siltfence Float
Fixed on the bottom
Flow direction
Fig. 3. Schematic figure of typical two types of siltfences. The left hand type is called as the standing type, and the right hand one the hanging type. Each dash line indicates the initial position of each type of the siltfence under a constant flow. Here, dcis the effective
siltfence height. を目的として，シルトフェンスの挙動をシミュレーションする方法を提案する． 一般に，シルトフェンスの設置方法には，Fig. 3のように自立型と垂下型があり，自立 型は潜水工事等が必須であるため，敷設が容易な垂下型が用いられることが多い．従っ て，本論文では垂下型を対象にする．一般に，垂下型では，上端に浮体，下端に重錘を設 置し，上端の両端を固定する．上両端のみを固定する際，膜の長さを固定端長（固定位置 を直線で結んだ長さ）より長くして敷設すると（通常の敷設では，固定端長を19.5とす ると膜長は20.0とするのが一般的である[6]），一定の流れがシルトフェンスに垂直方向 に発生している状況下では，Fig. 4のように変形することが知られている．ここで，鉛直 方向の変形は「ふかれ」，水平方向の変形は「たわみ」と呼ぶ． シルトフェンス設置の目的である汚濁物質の拡散抑制効果が上記の変形に依存すること は明らかであり，変形を予測できるモデリングが重要な研究開発課題であることが分か る．これまで，問題を簡単化するため，上端を直線で固定することで「たわみ」が発生し ないモデルに対して，一定の流れ場の下，「ふかれ」の度合いのみを調べる研究が行われ 多くの研究成果が得られてきた．例えば，小田らは，「ふかれ」たシルトフェンスの形状 を曲線と見做し，各点にかかる力を解析することで抗力と流速の関係式を導出し，関係式 に含まれる実験定数を水槽実験等から求めることで，おおよそ普遍的な「ふかれ」と流速 の関係を評価することに成功している[13]．最近では，武田らが実験定数を含め，水槽実 験結果から相関式を求める等の発展がある[20]． しかし，実際の利用では水路等にて上端を直線的に固定するような利用例はなく，シル トフェンスの「ふかれ」と「たわみ」の両者が同時に発生し，その形状はFig. 4のように 3次元的に変形し複雑な様相を示す．本論文の目的は，ある一定の流れ場の下，両者の変 形自由度を同時に取り扱い，その変形を予測可能なシミュレーション方法を提案すること にある． これまでに，シルトフェンスの汚濁物質拡散抑制効果を調べるため，計算機シミュレー (a) (b) Flow direction float weight fixed fixed
Fig. 4. Schematic figure of deformation of siltfence installed in a flume. The fence has floats at the top and weights at the bottom, and one of its top corners is fixed to one side of the flume and the other to the opposite side. Water current changes the shape of the siltfence from its initial state (yellow) to its deformed one (green). Here, we call the vertical shift of the bottom layer, that is (a), ”fukare”. In addition, the horizontal shift of the top one, that is (b), is called ”tawami”.
ションによる取り組みが行われてきた．これは，シルトフェンスの計算機シミュレーショ ンが可能となれば，様々な活用場面に即して，シミュレーションを行うことで，抑制効果 について評価可能となる他，抑制効果を向上させるための敷設方法の提案等も可能とな る等，大きな発展が期待できるからである．しかし，これまでの取り扱いでは，直接にシ ミュレーションすることは困難と考えられ，変形するシルトフェンスの代わりに，一旦， 「ふかれ」が発生したとして，アクリル板のような変形しない薄い直方体の障害物を設置 し，アクリル板の深さを，シルトフェンスの「ふかれ」が起こった後の実効深さとして反 映させ，流れ場のシミュレーション[14, 15]を行う等の方法で拡散抑制効果等が評価され てきた．しかし，実際のシルトフェンスの変形形状とは異なり，流れ場も現実を再現して いない等の課題がある．最近では「たわみ」が生じないよう，上端全体を直線として固定 したシルトフェンスを対象として，流体と膜体変形を連成する数値シミュレーションのア ルゴリズムが，金山ら[7, 8]により提案されたが，「ふかれ」と「たわみ」の両者を同時に 取り扱えるものは，著者の知る限り未だ皆無である． そこで，本論文では，金山らによる提案モデルを改良し，「ふかれ」と「たわみ」の２ つの変形を同時に取り扱うことが可能なモデルを提案する．モデルの妥当性を検証するた め，実際に一様な流れのある水路において敷設した垂下型シルトフェンスの変形をシミュ レーションし，同じ条件下での水路実験との比較を行う．結果として，両者の一致は極め て良く，提案したモデリングの妥当性が示される． 本論文の構成は以下の通りである．第2章では，シルトフェンスの変形に対し，従来の 解析的方法を概説した後，金山らにより提案された一様な流れ場中において生じる「ふか れ」を再現する流体との連成シミュレーションモデルを詳述し，その課題を記す．第3章 で金山らのモデルを発展させ，「ふかれ」及び「たわみ」を統一的に取り扱い可能とする モデルを提案する．第4章では，提案モデルを用いた数値シミュレーション結果を示すと
共に，同一条件下での水槽実験結果との比較を行う．第5章はまとめと結論である．
2.
シルトフェンス変形の理論及びシミュレーションモデリ
ングの現状
本章では，これまでに用いられてきたシルトフェンス変形の理論的取り扱いをレビュー し，汚濁物質拡散抑制効果を評価する上での課題を示す他，金山らにより示されたシミュ レーションモデリングを説明し，その改良に向けた方向性について議論する．2.1
流体の基礎方程式および計算方法
本研究では，垂下型シルトフェンス周囲の流体の運動を記述する（鉛直方向への運動成 分が重要）ため，流体の方程式に対し静水圧近似等を用いない．従って，流速(u, v, w)は 直交座標を用いた場合，次の式で記述される． •連続の式 ∂u ∂x + ∂v ∂y + ∂w ∂z =0 (2.1) •運動方程式 ∂u ∂t +u ∂u ∂x +v ∂u ∂y +w ∂u ∂z =− 1 ρ ∂p ∂x + ν ( ∂2u ∂x2 + ∂2u ∂y2 + ∂2u ∂z2 ) +Rx (2.2) ∂v ∂t +u ∂v ∂x +v ∂v ∂y +w ∂v ∂z =− 1 ρ ∂p ∂y + ν ( ∂2v ∂x2 + ∂2v ∂y2 + ∂2v ∂z2 ) +Ry (2.3) ∂w ∂t +u ∂w ∂x +v ∂w ∂y +w ∂w ∂z =− 1 ρ ∂p ∂z + ν ( ∂2w ∂x2 + ∂2w ∂y2 + ∂2w ∂z2 ) +Rz− g (2.4) ここで，pは圧力，νは動粘性係数(=10−6m2/s)，ρは流体の密度(=1000kg/m3)，gは重 力加速度(=9.8m/s2)である．R x, Ry, Rzは渦動粘性成分であり， Rx = ∂ ∂x ( 2AM∂u ∂x ) + ∂ ∂y { AM ( ∂u ∂y + ∂v ∂x )} + ∂ ∂z { AM ( ∂u ∂z + ∂w ∂x )} Ry = ∂ ∂y ( 2AM∂v ∂y ) + ∂ ∂x { AM ( ∂u ∂y + ∂v ∂x )} + ∂ ∂z { AM ( ∂v ∂z + ∂w ∂y )} Rz = ∂ ∂z ( 2AM∂w ∂z ) + ∂ ∂z { AM ( ∂u ∂z + ∂w ∂x )} + ∂ ∂y { AM ( ∂v ∂z + ∂w ∂y )} で与えられ，AM，即ち，渦動粘性係数については，Cs = 0.173としたスマゴリンスキーモデル[18]（Large Eddy Simulation乱流モデル）を用いて
Am = 12C2s(∆x∆y∆z) 2 3 ( ∂u ∂x )2 + ( ∂v ∂y )2 + ( ∂w ∂z )2 + (∂u ∂y + ∂∂vx )2 +(∂∂vz + ∂∂wy)2+(∂∂wx + ∂∂uz)2 2 1 2
共に，同一条件下での水槽実験結果との比較を行う．第5章はまとめと結論である．
2.
シルトフェンス変形の理論及びシミュレーションモデリ
ングの現状
本章では，これまでに用いられてきたシルトフェンス変形の理論的取り扱いをレビュー し，汚濁物質拡散抑制効果を評価する上での課題を示す他，金山らにより示されたシミュ レーションモデリングを説明し，その改良に向けた方向性について議論する．2.1
流体の基礎方程式および計算方法
本研究では，垂下型シルトフェンス周囲の流体の運動を記述する（鉛直方向への運動成 分が重要）ため，流体の方程式に対し静水圧近似等を用いない．従って，流速(u, v, w)は 直交座標を用いた場合，次の式で記述される． •連続の式 ∂u ∂x + ∂v ∂y + ∂w ∂z =0 (2.1) •運動方程式 ∂u ∂t +u ∂u ∂x +v ∂u ∂y +w ∂u ∂z =− 1 ρ ∂p ∂x + ν ( ∂2u ∂x2 + ∂2u ∂y2 + ∂2u ∂z2 ) +Rx (2.2) ∂v ∂t +u ∂v ∂x +v ∂v ∂y +w ∂v ∂z =− 1 ρ ∂p ∂y + ν ( ∂2v ∂x2 + ∂2v ∂y2 + ∂2v ∂z2 ) +Ry (2.3) ∂w ∂t +u ∂w ∂x +v ∂w ∂y +w ∂w ∂z =− 1 ρ ∂p ∂z + ν ( ∂2w ∂x2 + ∂2w ∂y2 + ∂2w ∂z2 ) +Rz− g (2.4) ここで，pは圧力，νは動粘性係数(=10−6m2/s)，ρは流体の密度 (=1000kg/m3)，gは重 力加速度(=9.8m/s2)である．R x, Ry, Rzは渦動粘性成分であり， Rx = ∂ ∂x ( 2AM∂u ∂x ) + ∂ ∂y { AM ( ∂u ∂y + ∂v ∂x )} + ∂ ∂z { AM ( ∂u ∂z + ∂w ∂x )} Ry = ∂ ∂y ( 2AM∂v ∂y ) + ∂ ∂x { AM ( ∂u ∂y + ∂v ∂x )} + ∂ ∂z { AM ( ∂v ∂z + ∂w ∂y )} Rz = ∂ ∂z ( 2AM∂w ∂z ) + ∂ ∂z { AM ( ∂u ∂z + ∂w ∂x )} + ∂ ∂y { AM ( ∂v ∂z + ∂w ∂y )} で与えられ，AM，即ち，渦動粘性係数については，Cs = 0.173としたスマゴリンスキーモデル[18]（Large Eddy Simulation乱流モデル）を用いて
Am = 12C2s(∆x∆y∆z) 2 3 ( ∂u ∂x )2 + ( ∂v ∂y )2 + ( ∂w ∂z )2 + (∂u ∂y + ∂∂vx )2 +(∂∂vz + ∂∂wy)2+(∂∂wx + ∂∂uz)2 2 1 2 で与える．ここで，∆x，∆y，∆zはシミュレーションを実施する際の x, y, z方向の格子間 隔である．以下で実行する本論文のシミュレーションでは，上記方程式を差分化しスタッ ガード格子を用いたSMAC(Simplified Marker and Cell)法[1]で計算する．
2.2
「ふかれ」及び「たわみ」の理論的解析
一様な流れの存在下，垂下型シルトフェンスを敷設した場合，Fig. 4にて記したように 「ふかれ」と「たわみ」の変形が起こる ．「ふかれ」は汚濁水の拡散抑制に対し最も重要 な変形自由度であり，「ふかれ」が大きくなると，水の流れに対抗する効果は殆ど失われ， 抑制効果は消失すると直感的に理解できる．従って，膜体としての抑制効果は「ふかれ」 を考慮した膜体の「有効膜高さ」(Fig. 3のdc)により，おおよそ決まっていることが分か る．小田らは，その重要なシルトフェンスの有効膜高さdcについて，以下の評価式を提 案している[13]． dc= W σ sin ( σd W ) , (2.5) σ =B1 2ρ ( Uh h − dc )2 . (2.6) ここで，dは膜状丈，Uは無限遠方における平均流速，Wは重錘の単位幅水中重量，σは 膜の単位長さ当たりの流体応力，Bは実験定数，hは膜の設置位置での水深である．本式 から，有効膜高さdc は，重錘の重量等のシルトフェンス自体の特徴量と水の流速等の外 部物理量との関係で表現される．つまり，上記の関係式を基にすると，シルトフェンスの 諸元を基に，どれ程の流速で，シルトフェンスがどの程度「ふかれ」るかが予測可能とな る．しかし，理論式中には実験定数Bが含まれており，これについては取り得る値の範囲 について複数の報告[2, 13, 19]がある一方，不定性があることから，理論を改良するか， 数値シミュレーションによるモデリングを提案する等の新たな研究が求められてきた．こ のBについては，最近，武田らにより実験結果を基に実験相関式が導出される等の研究報 告[20]もあり，今後も実験定数 Bの決定も含めた「ふかれ」および「有効膜高さ」は重 要な研究対象になると考えられる． その一方，「たわみ」については，「ふかれ」と比べ研究例は少ないが，椹木らにより理 論的な研究結果が報告されている[16]．対象となる「たわみ」は，流体からシルトフェン ス膜体要素が受ける流体応力とフェンス膜体要素部分にかかる膜体自身の張力との釣り合 いから求められる．この考え方に基づき，椹木らは上端と下端の両端の 4点を固定した 場合の「たわみ」の形状を解析的に求めており，基本的な考え方は「ふかれ」の解析的な 理論と同等である．しかし，実際の垂下型シルトフェンスにおいては，下端は固定され ないため，「たわみ」と「ふかれ」を独立に取り扱うことはできないことが分かる．なお， 垂下型シルトフェンスに対し，両者を独立に取り扱う理論式を用いた解析が行われた報 告[11,12]もあり，「たわみ」の少ない場合等では，近似的な理論として有用性があると考えられるが，現実的な敷設状況に対する評価の際には，両者を同時に解析できるモデリン グが強く求められる．
2.3
金山らのシミュレーションモデリングとその問題点
前節にてシルトフェンス変形の自由度である「ふかれ」および「たわみ」に対する理論 解析について，従来の研究例を紹介したが，シルトフェンスは様々な場面で利用されてお り，利用状況（活用する現場）に合わせ，その変形と汚濁物質の拡散抑制効果を議論すべ きであり，普遍的なシミュレーションモデリングの開発が期待される．実際，これまで に，金山らは2次元膜体モデルの連成計算法を提案した他，それを拡張し「たわみ」のな いシルトフェンスに対する3次元モデル用の流体連成計算方法を報告した[7, 8]．その3 次元モデル用計算方法[8]では，Fig. 5のようにシルトフェンスを長さが変わらない部材 （Fig. 5の実線）とそれらを接合させる節点からなる離散化したモデルで表現する．各節 点はFig. 5に示すように，シルトフェンスの部分領域を代表し，全ての部分領域を組み合 わせると過不足なくシルトフェンスの形状となる．金山らは，シルトフェンスの流れに対 する負荷抵抗力と形状変化の計算方法を以下のように提案した[8]． まず，Fig. 6のように，ある節点Oに対して作用している力を全て考え，それらの力に よりシルトフェンスの変形が静的に支えられていると仮定する．つまり，その節点が代表 する部分領域の重量W，結合している下方の節点C およびDからの力（FcおよびFd）， ２つの上方の節点AおよびBの結合方向の力（Fa およびFb）の他にシルトフェンスの法 線方向の力（流れから受ける力）F に対して， W + Fc+Fd +F + Fa+Fb = 0 (2.7) が成り立つと仮定する．この節点Oにおける力F の反力はシルトフェンスの復元力であ り，流れの計算に反映させる．ここで節点Oにおけるシルトフェンスの法線方向 ⃗f0は， ⃗f0 =(OA × ⃗⃗ OD + ⃗OB × ⃗OA + ⃗OC × ⃗OB + ⃗OD × ⃗OC) で与られる．なお，節点Oが最下段の場合，下方の節点Cおよび Dは存在しないため， その力は考慮せず，重錘の重力が代わりに働くとする．以上より，力の釣り合いの方程式 (2.7) を最下段から上段に向かって作成すると，各節点ごとに未知数３つに対して３つの 方程式が作成され，これらの方程式を計算することで，各節点におけるシルトフェンスの 法線方向に生じる復元力(−F)を計算することができる．この復元力をシルトフェンスに よる流れ場に与える抵抗力として流体の運動方程式(2.2)-(2.4)に加え，シルトフェンスの 影響を考慮した流れ場も計算される．ただし，実際の流れ場の計算では，ある節点におけ る抵抗力をその節点が代表しているシルトフェンスの部分領域に均等に分配する必要があ り，節点が代表している部分領域がどの計算セルにどの程度の割合で存在しているかを求 め，その割合に合わせてその節点の抵抗力を各セルに分配し，流れの運動方程式に反映さ せて，流れ場の計算を実施する．えられるが，現実的な敷設状況に対する評価の際には，両者を同時に解析できるモデリン グが強く求められる．
2.3
金山らのシミュレーションモデリングとその問題点
前節にてシルトフェンス変形の自由度である「ふかれ」および「たわみ」に対する理論 解析について，従来の研究例を紹介したが，シルトフェンスは様々な場面で利用されてお り，利用状況（活用する現場）に合わせ，その変形と汚濁物質の拡散抑制効果を議論すべ きであり，普遍的なシミュレーションモデリングの開発が期待される．実際，これまで に，金山らは2次元膜体モデルの連成計算法を提案した他，それを拡張し「たわみ」のな いシルトフェンスに対する3次元モデル用の流体連成計算方法を報告した[7, 8]．その3 次元モデル用計算方法[8]では，Fig. 5のようにシルトフェンスを長さが変わらない部材 （Fig. 5の実線）とそれらを接合させる節点からなる離散化したモデルで表現する．各節 点はFig. 5に示すように，シルトフェンスの部分領域を代表し，全ての部分領域を組み合 わせると過不足なくシルトフェンスの形状となる．金山らは，シルトフェンスの流れに対 する負荷抵抗力と形状変化の計算方法を以下のように提案した[8]． まず，Fig. 6のように，ある節点Oに対して作用している力を全て考え，それらの力に よりシルトフェンスの変形が静的に支えられていると仮定する．つまり，その節点が代表 する部分領域の重量W，結合している下方の節点C およびDからの力（FcおよびFd）， ２つの上方の節点 AおよびBの結合方向の力（Fa およびFb）の他にシルトフェンスの法 線方向の力（流れから受ける力）Fに対して， W + Fc+Fd +F + Fa+Fb =0 (2.7) が成り立つと仮定する．この節点Oにおける力F の反力はシルトフェンスの復元力であ り，流れの計算に反映させる．ここで節点Oにおけるシルトフェンスの法線方向 ⃗f0は， ⃗f0 = (OA × ⃗⃗ OD + ⃗OB × ⃗OA + ⃗OC × ⃗OB + ⃗OD × ⃗OC) で与られる．なお，節点 Oが最下段の場合，下方の節点C およびDは存在しないため， その力は考慮せず，重錘の重力が代わりに働くとする．以上より，力の釣り合いの方程式 (2.7) を最下段から上段に向かって作成すると，各節点ごとに未知数３つに対して３つの 方程式が作成され，これらの方程式を計算することで，各節点におけるシルトフェンスの 法線方向に生じる復元力(−F)を計算することができる．この復元力をシルトフェンスに よる流れ場に与える抵抗力として流体の運動方程式(2.2)-(2.4)に加え，シルトフェンスの 影響を考慮した流れ場も計算される．ただし，実際の流れ場の計算では，ある節点におけ る抵抗力をその節点が代表しているシルトフェンスの部分領域に均等に分配する必要があ り，節点が代表している部分領域がどの計算セルにどの程度の割合で存在しているかを求 め，その割合に合わせてその節点の抵抗力を各セルに分配し，流れの運動方程式に反映さ せて，流れ場の計算を実施する．Fig. 5. Modeling of siltfence using grid points. Fig. 6. Balance of force at grid point O. 次に，流れによるシルトフェンスの形状変化であるが，上記の流体計算により得られた 流れ場の情報を利用し，各節点を点要素として流れに合わせて移動させるが，節点位置と 流速を計算する位置は異なるため，節点周辺の流速情報を補間することで節点位置での流 速を計算して利用する．また，節点の移動に合わせてその節点が代表している部分領域も 同じように移動させる．このとき，各節点での流速は一様ではないため，移動後のすべて の部分領域を組み合わせても，本来のシルトフェンスの形状とは一致しない．金山らのモ デル化では，斜め方向（および左右端では上下方向）の節点間は長さが変化しない部材で 接続しているため，これらの距離が変化しないことを利用し，補正対象の１つの節点とそ れと結合する位置補正済みの上層の２つの節点で構成される平面上で，節点間の距離が適 切になるように，補正対象の節点の位置を補正する． この金山らの方法は流体シミュレーションをベースとして，シルトフェンスを膜体要素 （節点と節点をつなぐ部材）を用いてモデリングし，流体との連成計算として，その変形 および流動場を予測するものであり，計算順序を以下に示す． 1. SMAC法を用いて流体の連続の式(2.1)および運動方程式(2.2)-(2.4)から流れ場を 計算する． 2. 計算した流れ場に基づいてシルトフェンスの膜体要素である節点を移動させる． 3. 節点間の距離を考慮して節点の位置を補正する． 4. 式(2.7)に従って各節点に生じるシルトフェンスの復元力を計算する． 5. 4.で計算した各節点のシルトフェンスの復元力を節点が代表する部分領域に均等に 分散して，流体の運動方程式(2.2)-(2.4)に付加する． 6. 1に戻る． この連成計算により，流れ場によりシルトフェンスは変形し，その変形の反作用力とし
て抵抗力が計算され，流動場の変化に反映されることから，汚濁物質の移流及び拡散も 併せてシミュレーション可能であり，シルトフェンスの汚濁物質拡散抑制効果をシミュ レーションにより推定できると考えられる．本方法の妥当性を示す結果としては，「たわ み」が生じないように設置した矩形水路シミュレーションにおいて，シルトフェンスの変 形形状（有効膜高さdc）が，小田らにより提案された理論的な評価式(2.5), (2.6)と整合 したことが挙げられる．本論文の目的は，この方法を分析し，さらに，様々な実際の現場 に対して，適用可能とするための改良方法を提案することにある．著者らはその目的のた め，まず金山らの方法を実装し，Fig. 4のような水路に敷設したシルトフェンスのシミュ レーションを行った．その結果，以下の問題があることが分かった．金山らの方法は上端 を直線状に固定した場合に生じる「ふかれ」のみを取り扱う場合には十分だが，水平方向 の節点距離への配慮がないため，最上段の節点を節点間隔を考慮しながら流れに合わせて 移動させることができない．即ち，最上段の「たわみ」を再現することができない．した がって，シルトフェンス全体として「たわみ」を取り扱えないことが分かる．また，条件 によっては，節点で構成される菱形が水平方向に潰れ，水平方向の節点間距離が本来の距 離よりも長くなることを避けることができないことも分かった．実際，後述のシミュレー ション結果にて示すが，「たわみ」を考慮しない場合でも，10%以上も水平方向に長くな る場合が観察される．以上から，上端を直線状に固定した場合でも，下端等で「たわみ」 が生じてしまうようなケース（3次元的変形をするケース）には適用できないことが分か る．以上の課題を克服するためのモデリングを次章以降にて提案する．
3.
「ふかれ」と「たわみ」を同時に取り扱うシミュレーショ
ンモデリング
3.1
最急降下法を用いた節点位置の補正方法
前節の節点位置の補正では，上層の節点との間隔のみを固定し補正を行った．この補正 では，鉛直方向にのみ節点間隔を固定するが，水平方向については拘束せず，伸長・縮小 を妨げない．しかし，実際のシルトフェンスは，両方向に対し，殆ど伸長しない一方，折 れ曲がりやしわが生じることがあり，縮小することは許す必要があると考えられる．これ らのシルトフェンスの現実的挙動を考慮し，鉛直・水平方向共に，本来の間隔より縮小す る自由度を与える一方，伸長については強い拘束を与えるモデル化を行う．そこで，縮ん だ時より伸びた時の復元力が大きい非対称なバネでモデルの節点を接続していると想定 し，水平方向を含む補正前の節点間の距離と本来の距離の差を用いた評価関数を作成し， 最急降下法を用いて評価関数を最小化するように節点を補正する方法を提案する．このよ うな方法を導入することで，節点位置の補正は上記のバネでつながった弾性膜体とする物 理的モデリングの下，自然な形で実施される． まず，評価関数の作成方法を説明する．jは節点iと斜め方向，鉛直方向，または水平て抵抗力が計算され，流動場の変化に反映されることから，汚濁物質の移流及び拡散も 併せてシミュレーション可能であり，シルトフェンスの汚濁物質拡散抑制効果をシミュ レーションにより推定できると考えられる．本方法の妥当性を示す結果としては，「たわ み」が生じないように設置した矩形水路シミュレーションにおいて，シルトフェンスの変 形形状（有効膜高さ dc）が，小田らにより提案された理論的な評価式(2.5), (2.6)と整合 したことが挙げられる．本論文の目的は，この方法を分析し，さらに，様々な実際の現場 に対して，適用可能とするための改良方法を提案することにある．著者らはその目的のた め，まず金山らの方法を実装し，Fig. 4のような水路に敷設したシルトフェンスのシミュ レーションを行った．その結果，以下の問題があることが分かった．金山らの方法は上端 を直線状に固定した場合に生じる「ふかれ」のみを取り扱う場合には十分だが，水平方向 の節点距離への配慮がないため，最上段の節点を節点間隔を考慮しながら流れに合わせて 移動させることができない．即ち，最上段の「たわみ」を再現することができない．した がって，シルトフェンス全体として「たわみ」を取り扱えないことが分かる．また，条件 によっては，節点で構成される菱形が水平方向に潰れ，水平方向の節点間距離が本来の距 離よりも長くなることを避けることができないことも分かった．実際，後述のシミュレー ション結果にて示すが，「たわみ」を考慮しない場合でも，10%以上も水平方向に長くな る場合が観察される．以上から，上端を直線状に固定した場合でも，下端等で「たわみ」 が生じてしまうようなケース（3次元的変形をするケース）には適用できないことが分か る．以上の課題を克服するためのモデリングを次章以降にて提案する．
3.
「ふかれ」と「たわみ」を同時に取り扱うシミュレーショ
ンモデリング
3.1
最急降下法を用いた節点位置の補正方法
前節の節点位置の補正では，上層の節点との間隔のみを固定し補正を行った．この補正 では，鉛直方向にのみ節点間隔を固定するが，水平方向については拘束せず，伸長・縮小 を妨げない．しかし，実際のシルトフェンスは，両方向に対し，殆ど伸長しない一方，折 れ曲がりやしわが生じることがあり，縮小することは許す必要があると考えられる．これ らのシルトフェンスの現実的挙動を考慮し，鉛直・水平方向共に，本来の間隔より縮小す る自由度を与える一方，伸長については強い拘束を与えるモデル化を行う．そこで，縮ん だ時より伸びた時の復元力が大きい非対称なバネでモデルの節点を接続していると想定 し，水平方向を含む補正前の節点間の距離と本来の距離の差を用いた評価関数を作成し， 最急降下法を用いて評価関数を最小化するように節点を補正する方法を提案する．このよ うな方法を導入することで，節点位置の補正は上記のバネでつながった弾性膜体とする物 理的モデリングの下，自然な形で実施される． まず，評価関数の作成方法を説明する．jは節点iと斜め方向，鉛直方向，または水平 方向で結合している節点とする．この際， di, j= √ (xi − xj)2+(yi− yj)2+(zi − zj)2 は節点i, j間の距離であり，この距離と本来の節点間距離Di, jの差の関数として評価関数 を定義すれば良い．また，「たわみ」があるケースでは，シルトフェンスが大きく変形す る可能性があるため，ある程度自由に節点間距離が短くなれるように，本来の節点間距離 Di, jのs(< 1)倍までは評価値を0とし，それより節点間の距離が短くなる場合は，本来の 節点間距離との差に対応した評価値を与えることとする．一方，シルトフェンスが伸長す ることは，構造上想定していないため，節点間隔が本来の間隔よりも長くなる際は，その 差に対応した評価値を与えるが，節点間隔が長くなることは明らかに避ける必要があるた め，短い場合よりも大きいパラメータを乗じることとする．上記のような評価関数を定義 した場合，本来の距離との乖離を補正するには，最急降下法を利用し，乖離が大きいほど 評価関数の変化量が大きくなるようにするのが一般的である．従って，以下のように適切 に設定したパラメータ µ1, µ2 を乗じた距離の差の２乗を評価関数に用いるのが合理的で ある． Fi, j= µ1(di, j− Di, j)2 di, j ≥ Di, j 0 sDi, j ≤ di, j < Di, j µ2(sDi, j− di, j)2 di, j < sDi, j (3.1) ただし，本評価関数を用いても，任意の評価値になる節点の組み合わせは無数にあるた め*1，節点の補正量も評価関数に加えることでこの不定性を回避する．即ち，評価関数と して d =∑ i ∑j Fi, j+ λ(∆x2i + ∆y2i + ∆z2i) (3.2) を採用し，この値を最急降下法により最小化することで，シルトフェンスの形状を維持す る適切な補正が可能になる．ここで，∆xi, ∆yi, ∆zi は節点iのx, y, z方向の移動量，λは係 数パラメータであり，以下で実施する数値実験ではモデルのサイズに合わせて設定する． 次に補正の手順だが，金山らは，Fig. 7のように固定端である最上段から下方向に一段 ごとに補正を行う方法を提案している．ここで，青い節点は基準になる節点，緑の点は位 置を補正する節点，赤線は補正の際に考慮する距離である．一方，本提案では，Fig. 8に 示すように，二段ごとに前項で提案した評価関数を利用した最急降下法を用いて補正を行 う．その際，ある節点の補正には，その節点よりも上段か同じ段の節点の位置情報のみを 参考にする（ある節点は，それより下段の節点位置によって補正されない）．また，この 二段ごとの補正は局所的であるため，シルトフェンスの全体の形が不自然になる可能性が あるため，上記の補正後に全体の節点に対して同様の評価関数を用いた再補正を行う．こ *1一例として，上端を直線的に固定している場合のシルトフェンスに対して，上端を軸として，その周りに 任意の回転をさせても，節点位置は変化するが，評価値は変化しない．(a) Grid points correction in the first step. (b) Grid points correction in the second step.
Fig. 7. Relationship between the corrected grid points (green) and their reference ones (blue) in Kanayama’s correction scheme.
(a) Grid points correction in the first step. (b) Grid points correction in the second step.
Fig. 8. Relationship between the corrected grid points (green) and their reference ones (blue) in the proposed correction scheme.
の補正においても二段ごとの補正と同様に，ある節点を補正する際にはその節点よりも上 段か同じ段の節点の位置情報しか参考にしない．この補正方法を用いることで，「たわみ」 が発生するような場面も含めて様々な条件でシルトフェンスの変形がシミュレーション可 能となった．
(a) Grid points correction in the first step. (b) Grid points correction in the second step.
Fig. 7. Relationship between the corrected grid points (green) and their reference ones (blue) in Kanayama’s correction scheme.
(a) Grid points correction in the first step. (b) Grid points correction in the second step.
Fig. 8. Relationship between the corrected grid points (green) and their reference ones (blue) in the proposed correction scheme.
の補正においても二段ごとの補正と同様に，ある節点を補正する際にはその節点よりも上 段か同じ段の節点の位置情報しか参考にしない．この補正方法を用いることで，「たわみ」 が発生するような場面も含めて様々な条件でシルトフェンスの変形がシミュレーション可 能となった．
Fig. 9. Area corresponding to grid point at top of siltfence. The left figure shows the case when all points on the top layer are fixed. The right figure shows the one when only the ends of the top layer are fixed.
3.2
最上段の節点の取り扱い
シルトフェンスの最上段の節点に対し「たわみ」を許し両端の節点のみを固定する場 合，最上段の節点も流れに合わせて移動させる．実際のシミュレーションでは，最上段の 節点を流速にあわせて移動させ，水平方向に対して，最急降下法を用いて節点位置を補正 する．また，Fig.5に示したように，シルトフェンスの上端が固定されている場合，上端 からの力は上から2段目の節点にかかる力より計算されるが，「たわみ」が発生する場合， 最上段の節点も移動するため，Fig. 9のように最上段の節点に生じる力より計算する．3.3
シルトフェンスから流れ場へ与える影響
シルトフェンスが流動場へ与える力として，2.1節で示したシルトフェンス自身の重み から生じる成分があり，これは最急降下法による補正方法を用いた場合でも同様の方法で 与えることができる．しかし，上記の補正方法により，節点位置を移動させることで，シ ルトフェンスの移動分の量の水も移動させる必要があるため，その分の力も考慮する必要 がある*1．ここでは簡単のため，x軸方向についてのみ考える．まず，ある節点の位置が x0，その位置での流速がvxであるとする．∆t秒後，流れと補正により節点が x1 に移動 したとする．本来の節点の流速はvxのはずであるが，補正によりその流速は(x1 − x0)/∆t になるため，この節点には((x1− x0)/∆t − vx)/∆t(= ax)の加速度が作用していることにな る．また，この節点が担当するシルトフェンスの領域を yz平面に射影した面積をSx と すると，この加速度は体積Sxvx∆t(= Vx)の流体に影響していると考えられる．そのため， 流体の密度をρとするとある節点におけるシルトフェンスの補正によるx方向の力は axVxρ *1金山らの方法では，節点位置の補正を水を通過させない膜面内に制限しているため，補正による節点の移 動によって水は移動しないと考えており，この力を考慮していない[8]．となる．y, z方向の力も同様に計算できる．これらの力をシルトフェンスの自重による抵 抗力F に加えることで，補正による流れ場への影響を正確に考慮することが可能となる．
4.
数値計算
4.1
金山らの方法と比較
今回提案したシルトフェンスの変形計算の有効性を確認するため，Fig. 10に示す節点 間隔を縦横とも1mとし，シルトフェンスを表現している節点の１つを仮に移動させ，金 山らの方法と今回提案した方法でシルトフェンス形状の補正を行い，その補正後の形状 の比較を行う．ここで，提案方法を用いて計算する際の最急降下法は，アルミホのルー ル[3]（パラメータはα =0.00001，β =0.5）に基づいて計算する．また，評価関数(3.1) および (3.2)のパラメータはµ1 = 1, µ2 = 0.05（斜め方向と鉛直方向）,0（横方向）, s = 1, λ =0.1と設定した． 2つの補正方法の違いを示すため，Fig. 10に示したA点及びB点をそれぞれ膜面に対 し法線方向に0.1メートル移動させ，各々の方法で補正を行った際の各方向の節点間の距 離と本来の距離の差（ずれ）をTable 1に示す*1．ここでの横方向(Horizontal)，斜め方向 (Diagonal)，鉛直方向(Vertical)の基準の長さは，各々，1m，√12+0.52(≈ 1.118)m，2m である．この結果から，金山らの方法では補正により生じる節点間の距離のずれは横方向 のみであるが，提案方法においては，ずれをすべての方向に対して適度に生じさせ，本来 の距離よりも長くなることに強い制限をかけているため，伸長する量は金山らの方法の 10分の1程度であることが確認される． 次に，Fig. 11に節点ごとの移動距離を示す．この結果から，金山らの方法では，移動さ せた節点に対して斜め下方向の節点や，真下方向の節点の位置が大きく移動しているが， そこから離れた節点はほとんど変動していないことが確認できる．これは，シルトフェン スが面として連続し，つながっていることを考えると不自然である．一方，提案方法では 金山らの方法と同様の節点以外も移動しており，特に，横方向の長さを考慮していること から，移動させた節点と同じ高さの層の節点も移動していることが確認できる．また，各 節点の移動量は金山の方法よりも小さいことも確認できる．以上から，提案方法は大域的 に節点の配置を考慮して補正を行っていると考えられる．4.2
水路での数値実験（手法及び理論の比較）
本節では，Fig. 12に示す水深hが15m，横幅が22m，長さが200mの水路において， 境界の流速を 0.1m/秒から1.0m/秒まで0.1m/秒刻みで変化させ，金山らの方法と提案方 *1金山らの方法では，斜め方向及び鉛直方向の節点間隔は固定しているため，間隔の違いは数値計算誤差に よるものしかない．そのため，この表では省略する．となる．y, z方向の力も同様に計算できる．これらの力をシルトフェンスの自重による抵 抗力F に加えることで，補正による流れ場への影響を正確に考慮することが可能となる．
4.
数値計算
4.1
金山らの方法と比較
今回提案したシルトフェンスの変形計算の有効性を確認するため，Fig. 10に示す節点 間隔を縦横とも1mとし，シルトフェンスを表現している節点の１つを仮に移動させ，金 山らの方法と今回提案した方法でシルトフェンス形状の補正を行い，その補正後の形状 の比較を行う．ここで，提案方法を用いて計算する際の最急降下法は，アルミホのルー ル[3]（パラメータはα =0.00001，β =0.5）に基づいて計算する．また，評価関数(3.1) および (3.2)のパラメータはµ1 = 1, µ2 = 0.05（斜め方向と鉛直方向）,0（横方向）, s = 1, λ =0.1と設定した． 2つの補正方法の違いを示すため，Fig. 10に示したA点及びB点をそれぞれ膜面に対 し法線方向に0.1メートル移動させ，各々の方法で補正を行った際の各方向の節点間の距 離と本来の距離の差（ずれ）をTable 1に示す*1．ここでの横方向(Horizontal)，斜め方向 (Diagonal)，鉛直方向(Vertical)の基準の長さは，各々，1m，√12+0.52(≈ 1.118)m，2m である．この結果から，金山らの方法では補正により生じる節点間の距離のずれは横方向 のみであるが，提案方法においては，ずれをすべての方向に対して適度に生じさせ，本来 の距離よりも長くなることに強い制限をかけているため，伸長する量は金山らの方法の 10分の1程度であることが確認される． 次に，Fig. 11に節点ごとの移動距離を示す．この結果から，金山らの方法では，移動さ せた節点に対して斜め下方向の節点や，真下方向の節点の位置が大きく移動しているが， そこから離れた節点はほとんど変動していないことが確認できる．これは，シルトフェン スが面として連続し，つながっていることを考えると不自然である．一方，提案方法では 金山らの方法と同様の節点以外も移動しており，特に，横方向の長さを考慮していること から，移動させた節点と同じ高さの層の節点も移動していることが確認できる．また，各 節点の移動量は金山の方法よりも小さいことも確認できる．以上から，提案方法は大域的 に節点の配置を考慮して補正を行っていると考えられる．4.2
水路での数値実験（手法及び理論の比較）
本節では，Fig. 12に示す水深hが15m，横幅が22m，長さが200mの水路において， 境界の流速を 0.1m/秒から1.0m/秒まで0.1m/秒刻みで変化させ，金山らの方法と提案方 *1金山らの方法では，斜め方向及び鉛直方向の節点間隔は固定しているため，間隔の違いは数値計算誤差に よるものしかない．そのため，この表では省略する．Fig. 10. Schematic figure of the target siltfence. We corect the positions of the grid points form the initial state created by shifting point A or B using each correction scheme. Table 1. Deviation from original distance between neighboring grid points using Kanayama’s and proposed correction schemes after moving a grid at point A and B (see Fig. 10).
Interval of deviation (m)
Kanayama’s scheme Proposed scheme
Horizontal Horizontal Diagonal Vertical Point A [-0.0000, 0.0049] [-0.0011, 0.0003] [-0.0040, 0.0007] [-0.0005, 0.0000] Point B [-0.0196, 0.0098] [-0.0041, 0.0009] [-0.0059, 0.0007] [-0.0006, 0.0000] 法によるシルトフェンスの変形・移動計算を行う．流体計算にはSMAC法を利用し，格 子間隔は1m×1m×1mとした．また，水路の壁面・底面での流れはnon-slip条件とした． さらに，シルトフェンスの諸元として，膜丈長 d を10m，重錘の単位幅水中重量W を 1980N/m = (200kg f /m)，シルトフェンス自体の重さは無視し，上流端から7.5mの位置 に最上段を固定して設置し，モデル化のための節点間隔は縦1m，横1mにした．この条 件の下，時間刻みdt を0.01 秒とし，30000ステップ（実時間で300秒）のシミュレー ションを実施した．本シミュレーションにおいて金山の方法は2.3節に記した計算順序で 連成計算を行うが，今回提案した方法では3.の節点の補正の際に，今回提案した補正方法 を用いること，および，5.のシルトフェンスの自重による復元力を運動方程式に加える際 に3.3節で説明した力も一緒に加えることの2点が異なっている． 提案方法によるシルトフェンスの変形計算のため，評価関数(3.1)および (3.2)のパラ
Kanayama’s correction scheme. Proposed correction scheme. (a) Correction for the shift at the point A.
Kanayama’s correction scheme. Proposed correction scheme. (b) Correction for the shift at the point B.
Fig. 11. Amount of correction of grid points by the two correction schemes.
メータは，µ1 =1, µ2 = 0.05（斜め方向と鉛直方向）,0（横方向）, s = 1, λ = 0.1とし，アル ミホのルール（パラメータはα =0.00001，β =0.5）に基づいて計算する．また，どちら の補正方法でも，計算で得られた補正点が水路外になった場合は，各軸に平行に強制移動 させ，水路内に戻るように再補正する． これらの条件の下，シミュレーションにより得られた下層の中央部分の「ふかれ」量と 小田ら[13]により提案された推定式(2.5) および(2.6) による「ふかれ」量の比較をFig. 13に示す．なお，式(2.5) の適用範囲は2π/d ≤ deであり，deがこの範囲よりも小さい 場合は他の先行研究と同様に，式(2.5)のsinの値を1とする．この結果から，提案方法 による結果は，金山らの方法の結果とほぼ一致しており，実験定数をB = 1.5とした場合 の小田らの推定式と整合している．また，提案方法と金山らの方法を用いた際の流速が 0.5m/sの300秒後のシルトフンスと水路中央での流れ場の縦断図をFig. 14に示すが，こ の結果から，2つの方法において，縦断面でのシルトフェンスの「ふかれ」量と流れ場は ほぼ同じになっていることが確認できる． 次に，シルトフェンスの形状の詳細な違いを比較するため，提案方法を用いて計算され たシルトフェンスの形状（下流上方から俯瞰）を，金山らの方法で計算した節点の位置と
Kanayama’s correction scheme. Proposed correction scheme. (a) Correction for the shift at the point A.
Kanayama’s correction scheme. Proposed correction scheme. (b) Correction for the shift at the point B.
Fig. 11. Amount of correction of grid points by the two correction schemes.
メータは，µ1 =1, µ2 = 0.05（斜め方向と鉛直方向）,0（横方向）, s = 1, λ = 0.1とし，アル ミホのルール（パラメータはα =0.00001，β =0.5）に基づいて計算する．また，どちら の補正方法でも，計算で得られた補正点が水路外になった場合は，各軸に平行に強制移動 させ，水路内に戻るように再補正する． これらの条件の下，シミュレーションにより得られた下層の中央部分の「ふかれ」量と 小田ら[13]により提案された推定式(2.5) および(2.6)による「ふかれ」量の比較をFig. 13に示す．なお，式(2.5) の適用範囲は2π/d ≤ deであり，deがこの範囲よりも小さい 場合は他の先行研究と同様に，式(2.5)のsinの値を1とする．この結果から，提案方法 による結果は，金山らの方法の結果とほぼ一致しており，実験定数をB = 1.5とした場合 の小田らの推定式と整合している．また，提案方法と金山らの方法を用いた際の流速が 0.5m/sの300秒後のシルトフンスと水路中央での流れ場の縦断図をFig. 14に示すが，こ の結果から，2つの方法において，縦断面でのシルトフェンスの「ふかれ」量と流れ場は ほぼ同じになっていることが確認できる． 次に，シルトフェンスの形状の詳細な違いを比較するため，提案方法を用いて計算され たシルトフェンスの形状（下流上方から俯瞰）を，金山らの方法で計算した節点の位置と
Fig. 12. Schematic figure of flume. Fig. 13. Comparison of effective siltfence height with Oda’s formula.
0 5 10 15 20 25 30 0 3 6 9 12 15 x-direction (m) Depth (m) あｃ
(a) Kanayama’s scheme.
あｃ 0 5 10 15 20 25 30 0 3 6 9 12 15 x-direction (m) Depth (m) (b) Proposed scheme.
Fig. 14. Flow field and shape of siltfence.
の差（ずれ）をFig. 15に示す．この結果から，シルトフェンス中心付近のずれは小さい が，端に近づくとずれが大きくなることが確認できる．また，流速ごとの各方向の節点間 隔と本来の間隔との差（ずれ）をTable 2に示す．この結果から，金山らの方法では，シ ルトフェンスを表現する節点の斜めおよび鉛直方向の距離が変化しないため，横方向の節 点距離を伸縮させることでその形状を調整するため，基準の距離と比較し，横方向の距離 が最大で10%程度長くなることが確認できる．一方，提案方法を用いた場合，すべての 方向の間隔を伸縮させることで，その形状を調整し再現するため，横方向に対しては最大 でも1%程度しか伸長しないことが確認できる． 以上から，提案方法はシルトフェンスを大きく伸長させることなく適切に変形している ことが確認できる．これは，金山らの方法で見られる横方向の節点間隔の不自然な伸長を 克服できていることを示している．
Fig. 15. Overhead view of siltfence simulated by proposed method. The color of grid point indicates the gap between each point and the corresponding one calculated by Kanayama’s scheme.
Table 2. Deviation from original distance between neighboring grid points using Kanayama’s and proposed correction scheme.
Interval of deviation (m)
velocity Kanayama’s scheme Proposed scheme
(m/sec) Horizontal Horizontal Diagonal Vertical 0.1 [-0.3467, 0.2088] [-0.4166, 0.0054] [-0.0555, 0.0038] [-0.0519, 0.0017] 0.2 [-0.3959, 0.2534] [-0.4392, 0.0074] [-0.0453, 0.0059] [-0.0376, 0.0016] 0.3 [-0.2977, 0.0553] [-0.3337, 0.0056] [-0.0315, 0.0030] [-0.0134, 0.0006] 0.4 [-0.2702, 0.0428] [-0.3702, 0.0045] [-0.0350, 0.0031] [-0.0153, 0.0007] 0.5 [-0.3351, 0.0409] [-0.2788, 0.0050] [-0.0206, 0.0031] [-0.0025, 0.0007] 0.6 [-0.3704, 0.0258] [-0.3549, 0.0046] [-0.0105, 0.0029] [-0.0013, 0.0007] 0.7 [-0.2421, 0.0538] [-0.2653, 0.0029] [-0.0308, 0.0035] [-0.0055, 0.0010] 0.8 [-0.5139, 0.1011] [-0.3360, 0.0054] [-0.0237, 0.0037] [ 0.0001, 0.0010] 0.9 [-0.2267, 0.0145] [-0.2059, 0.0073] [-0.0096, 0.0027] [-0.0019, 0.0014] 1.0 [-0.3253, 0.0437] [-0.1808, 0.0103] [-0.0208, 0.0039] [-0.0012, 0.0011]
4.3
「たわみ」のあるシルトフェンス
Fig. 16のような水路に「たわみ」が生じるシルトフェンスを設置した際の提案方法の有 効性を評価する．なお，比較検証のため，水槽実験も同様の条件で実施する．実験では， 流量が2.5m3/分のポンプを1～4個利用して流れを生成する．この時，水路の断面積を考 慮すると，１つのポンプあたり，0.1302m/秒の流速が得られることがわかる．また，シルFig. 15. Overhead view of siltfence simulated by proposed method. The color of grid point indicates the gap between each point and the corresponding one calculated by Kanayama’s scheme.
Table 2. Deviation from original distance between neighboring grid points using Kanayama’s and proposed correction scheme.
Interval of deviation (m)
velocity Kanayama’s scheme Proposed scheme
(m/sec) Horizontal Horizontal Diagonal Vertical 0.1 [-0.3467, 0.2088] [-0.4166, 0.0054] [-0.0555, 0.0038] [-0.0519, 0.0017] 0.2 [-0.3959, 0.2534] [-0.4392, 0.0074] [-0.0453, 0.0059] [-0.0376, 0.0016] 0.3 [-0.2977, 0.0553] [-0.3337, 0.0056] [-0.0315, 0.0030] [-0.0134, 0.0006] 0.4 [-0.2702, 0.0428] [-0.3702, 0.0045] [-0.0350, 0.0031] [-0.0153, 0.0007] 0.5 [-0.3351, 0.0409] [-0.2788, 0.0050] [-0.0206, 0.0031] [-0.0025, 0.0007] 0.6 [-0.3704, 0.0258] [-0.3549, 0.0046] [-0.0105, 0.0029] [-0.0013, 0.0007] 0.7 [-0.2421, 0.0538] [-0.2653, 0.0029] [-0.0308, 0.0035] [-0.0055, 0.0010] 0.8 [-0.5139, 0.1011] [-0.3360, 0.0054] [-0.0237, 0.0037] [ 0.0001, 0.0010] 0.9 [-0.2267, 0.0145] [-0.2059, 0.0073] [-0.0096, 0.0027] [-0.0019, 0.0014] 1.0 [-0.3253, 0.0437] [-0.1808, 0.0103] [-0.0208, 0.0039] [-0.0012, 0.0011]
4.3
「たわみ」のあるシルトフェンス
Fig. 16のような水路に「たわみ」が生じるシルトフェンスを設置した際の提案方法の有 効性を評価する．なお，比較検証のため，水槽実験も同様の条件で実施する．実験では， 流量が2.5m3/分のポンプを1～4個利用して流れを生成する．この時，水路の断面積を考 慮すると，１つのポンプあたり，0.1302m/秒の流速が得られることがわかる．また，シル トフェンスの最下段に設置した重錘の単位幅水中重量は5.194N(=0.53kgf/m)とし，シル トフェンスの幅は水路幅0.8mに対して0.82mと設定した． 一方，数値シミュレーションでは，1つのポンプによる流速の半分である0.0651m/秒を 基準とした1～8倍の流速（ポンプ0.5台から4台）を上流の境界に一様に加える．流体計 算にはSMAC法を利用し，時間刻み幅を0.001 秒，メッシュサイズはx, y方向は0.05m とし，z方向は 0.025mとする．また，壁面・底面の条件は3.2 節と同様non-slip条件と し，シルトフェンスの節点間隔は横0.05125m，縦0.04mとし，シルトフェンス自体の重 さは無視する．ただし，初期状態での横方向の間隔は 0.05mとし，ステップごとに徐々 に伸ばし，500ステップ(0.5秒経過)後の時点で0.05125mになるようにしている．評価 関数(3.1)および(3.2)のパラメータはµ1 = 100, µ2 =0.1（斜め方向と鉛直方向）, 0（横 方向）, λ = 1と設定する．また，「たわみ」のあるシルトフェンスを対象にしているため， 節点間が若干縮んでも復元力が働かないように s = 0.8とした． 以上の条件で30000ステップ(30秒)計算した際の各方向の節点間の間隔と本来の間隔 の差をTable 3に，「ふかれ」量，水路実験による「ふかれ」量および小田らによる推定式 から導かれる「ふかれ」量をFig. 17に示す．ただし，小田らの推定式において，本来U は無限遠での平均流速であるが，この実験では上流と下流で水路の幅が異なり流速が異な るため，ここでは U として上流の境界流速を採用する．この結果から，提案方法により シルトフェンスが本来の長さ以上に伸びることが適切に制限されていること，B = 0.65と した小田らの推定式と良好な一致が見られる．また，水路実験結果と「ふかれ」量の比較 をすると，ポンプが1台および２台の時（境界の流速が0.1302m/秒および0.2604m/秒の 時）の水路実験での値はシミュレーションや B = 0.65とした推定式から導かれる値より も大くなる（膜の有効深さは小さい）．これは，シミュレーションにおいてはシルトフェ ンス自体の重さ（および，その硬さ）が考慮されておらず，流速が小さい場合，これらの 影響が無視できないため「ふかれ」量が異なるものと考えられる．また，シルトフェンス の形状をFig. 18に，実際の水路実験のスナップショットを Fig. 19に示す．この結果か ら，シルトフェンスの形状は（流速が遅い場合の「ふかれ」量が異なることを別にして）， どちらの場合もシルトフェンスの中央付近が膨らむ形状になっており，おおよそ同じ形状 を示していることが確認できる． 以上の結果から，適切な補正パラメータを設定した提案方法を用いることで，「たわみ」 のあるシルトフェンスでも，「ふかれ」量やその形状を再現することができると結論付け られる．5.
まとめ
本論文では，水路内に設置されたシルトフェンスの振る舞いをシミュレーションするた め，最急降下法を用いたモデル化を行った．金山らの方法では，シルトフェンスを節点群 で表現し，左右端の鉛直方向および斜め方向の節点間の距離を固定しシルトフェンスの移Table 3. Deviation from original distance between neighboring grid points.
Boundary flow Interval of deviation (m)
velocity (m/sec) Horizontal Diagonal Vertical 0.0651 [-0.0019, 0.0000] [-0.0133, 0.0000] [-0.0199, 0.0000] 0.1302 [-0.0021, 0.0000] [-0.0171, 0.0000] [-0.0164,-0.0018] 0.1953 [-0.0037, 0.0001] [-0.0168, 0.0000] [-0.0188, 0.0000] 0.2604 [-0.0092, 0.0001] [-0.0167, 0.0000] [-0.0067, 0.0000] 0.3255 [-0.0048, 0.0001] [-0.0163, 0.0000] [-0.0001, 0.0000] 0.3906 [-0.0004, 0.0002] [-0.0161, 0.0001] [-0.0011, 0.0000] 0.4557 [-0.0015, 0.0002] [-0.0159, 0.0001] [-0.0001, 0.0000] 0.5208 [-0.0012, 0.0002] [-0.0158, 0.0001] [-0.0000, 0.0000]
Fig. 16. Schematic figure of flume. Fig. 17. Comparison of amount of blown out siltfence between numerical simulation and flume experiment. 動を補正している．その結果，横方向の節点間隔が実際の間隔よりも伸長する傾向が表 れ，「たわみ」が生じる場合のシミュレーションには適さない．一方，提案方法は最急降 下法を用いて左右端の鉛直方向および斜め方向だけでなく横方向の節点間隔も考慮する． こうして，金山らの方法よりも多くの計算量を必要とするが，各方向の節点間隔を現実的 な長さに適切に調節することが可能となる．実際の「たわみ」のあるシルトフェンスを水 路に設置した際の実験と比較し，シルトフェンス自体の重さが無視できる十分な流速であ れば，同様の形状を再現できることを確認した．なお，流れ場の計算に用いているSMAC 法では圧力の修正量を計算するため，線形方程式を解くが，この計算コストが大きく，全 体として補正のための計算時間の増加は大きくない．一方，金山らの方法は，適切に未知 数と条件を決定しているため，今回提案した方法のように反復計算を行う必要がなく計算
