境界要素法入門 part2 [境界要素法]
境界要素法入門 part ２
［変分法-3］
∬dxdyの積分領域Rは固定で、(1)は「定積分」でした。
変分法の常套手段として、(1)を部分積分し変分δφの微分を消しますが、その目的は後述します。(1)の場合は、次のグリーンの定理を使います。グリーンの定理とは、ガウスの発散定理の部分積分版です。
∇は勾配で、∬dxdyは固定領域Rで行います。
はRの境界Cの外法線線素ベクトルで、・は内積です。∫・dsはC上の線積分、Δはラプラシアン
です。
(2)でf＝φ，g＝δφとおけば、確かに、
なので、(1)は、
となります。
δφは、(3)を満たす解関数φ(x，y)からのテスト関数のずれを表す、微小だけれど立派な2変数関数でした。可能なテスト関数とは、任意の偏微分可能な関数です。よってδφは、微小な任意の偏微分可能な関数を表してる事になります。ここで面倒なので、偏微分可能なら連続だろうからδφの範囲を、任意の微小な連続関数とします(^^;)。
最初に(3)の1項目です。このまま続行すると最後には、非圧縮性渦なし完全流体の支配方程式を表す偏微分方程式にたどり着くのですが、偏微分方程式の解は一般に、解析領域Rの境界C上の境界条件を与えて初めて、解が一意に決まります。
どんな境界条件かは知りませんが、今はとにかくどんな境界条件にも適用可能なφの条件が知りたい訳ですから、φはとにかく境界条件を満たすとするのが妥当です。すなわちRの境界C上で、ずれδφ＝0です。そして(3)の1項目はC上の線積分でした。δφ＝0から、それは0とおくのが妥当です。
ここで(4)が「不定積分」だったら明らかに、Δφ＝0を導けます。しかし(4)は「定積分」でした。
そこでもう一回、δφの範囲を思い出します。δφは任意の微小な連続関数でした。任意の微小な連続関数って、どうやって作ったら良いでしょう？。
最も乱暴には、「任意の連続関数fを勝手に一つ持ってきて、それに任意に勝手で微小な定数εをかければOKよね？」という方法があります（← 間違ってはいない(^^;)）。
f(x，y)を任意の2変数連続関数として、δφ＝εf(x，y)とします。
でもεは定数で積分とは無関係です。従って、任意の連続関数fに対して、
でなければならない事になります。
ここで［変分学の基本定理］をご紹介します。2次元でも3次元でも同じなので、1次元で紹介します。
［変分学の基本定理］
f(x)とg(x)を（閉区間[a,b]の）連続関数とする。任意のf(x)について、
である条件は、g(x)＝0。
この定理は約半世紀前には本当に、［変分学の基本定理］と言われてたようです。しかし［基本定理］という仰々しい名前に対して証明が余りにこっ恥ずかしいからなのか、いつの間にか明確に語られなくなりました（← 自分は非常に問題な気がします）。実際うるさい事さえ言わなければ上記定理は、馬鹿みたいに簡単に証明できます。
fは任意だから、f＝gの場合がある。この時は、
非負の連続関数の積分が（0以上の値を持つグラフの面積が）0になるなら、g(x)＝0は明らか。逆にg(x)＝0なら、与式は明らか。
・・・という訳で、(5)よりfが任意の連続関数なので（Δφもきっと連続だろうから）、
Δφ＝0
となり、非圧縮性渦なし完全流体の支配方程式が得られます。
もうおわかりだと思いますが、「変分法の常套手段として、(1)を部分積分し変分δφの微分を消す」目的は、［変分学の基本定理］を使いたいがためなんですよ。
［変分学の基本定理］の証明を厳密にやろうと思ったら、「非負の連続関数の積分が0になるなら、g(x)＝0」をちゃんとやる必要があり、正式には測度論まで遡らないと駄目です。それが無理な注文である事はわかりますが、だったら「こんなの明らかよね？」って言えば良いじゃないですか。そうせずに部分積分の目的も明確に言わないくらいなら。いかに「こっ恥ずかし」かろうと・・・。
[変分法-4］
・・・こうして数々の数学的いい加減さを闇に葬って来ました(^^;)。中でも目玉級に問題なのは、テスト関数の番号付けパラメータ、αの存在でしょう。
じつはこれ、関数空間上の距離に相当します。関数空間とは乱暴に言うと、2変数の連続関数全体の集合なんかですが、その要素である関数の間に距離を定義できる事がわかっています。距離の定義も色々ありますけれど、どのような距離を取ろうと本質的には同等です。特に2つの関数f，gについての距離をα(f，g)として、
α(f，g)＝0 ⇔ f＝g
はどんな距離においても成り立ちます。
関数空間に距離αを入れて距離空間とし、その上でJ(α)について通常の微積をやるというのが、変分法の正式なやり方です。これは実数上の微積分だって、じつは同じです。xとyを実数とすれば実数全体の集合には最初から、距離α＝|x－y|が備え付けになってるので普段は意識しないだけです。
普通は関数空間までやる余裕はなく、そこは省略されます。それでほとんど変分の説明もないままに、ほとんど突然に「こんなの明らかだよね？」と、
などが、ドド～ン！と現れます(^^;)。
また［変分学の基本定理］が明確に紹介される事もなく、従って(1)を、
の形に変形する目的も暗中模索のまま、いつのまにか訳もわからずに、
へ達する事になります。普通の神経では耐えられませんな(^^;)。
（執筆 ddt³）
境界要素法入門 part 1 [境界要素法]
境界要素法入門 part 1
変分法１
自分は有限要素法（FEM）の前に境界要素法（BEM）をやっていたという、変な奴です(^^;)。
ちなみにFEMは、Finite Element Methodと書いて、ファイナイト・エレメント・メソッドと読みます（フィニットと読んじゃだめ）。BEMは、Boundory Element Methodと書いて、バウンダリー・エレメント・メソッドと読みます。
BEMにもFEMにも、色んな定式化がありますが、やっぱり変分原理の成り立つケースが、一番綺麗に行きます。変分原理とは何かと言うと、普通は最小エネルギー原理と等価になります。では最小エネルギー原理とは何かと言うと、ポテンシャルエネルギーの定義から明らかだったりするんです。
バネ定数kのバネが天井に鉛直に固定され、先端に質量mの質点がぶら下がってる状態を想像します。重力をmg（gは重力加速度でg=9.0m/s²）として、バネと質点のポテンシャルエネルギーは、次のように計算しますよね？。釣り合い方程式は、
なので（xはバネの伸び）、上式を伸びxで積分して、
ここでE₀は積分定数です。上記積分は簡単にできて、
です。E₀は例えば、x＝0の地点を基準点に選べば0にでき、質点の位置エネルギーの減少分がバネの弾性エネルギーとして貯まった、と読めるようになります。
次に、
という関数を考え、E(x)の最小点を探してみます。とうぜんE(x)をxで微分してdE/dx＝0としますよね？。
その結果は、
となり、釣り合い方程式に戻れます。バネはx＝mg/k伸びて釣り合う、です。
でもこれって、当たり前じゃないですか？。だってE(x)の関数形は、kx－mg＝0を積分して作ったんですから。本質的には変分原理とはこれだけです。
釣り合い方程式（系の支配方程式） ⇔ エネルギーの最小化（変分原理）
一般の場合は今のような考えを、3次元の弾性体や運動系にも適用できるように拡張しただけのものです。運動系の場合、支配方程式は運動方程式です。厳めしく変分原理などと言うと難しそうに聞こえますが、エネルギーの定義まで遡って考えれば、「当然じゃん」という事が案外多いです。
※ 運動系のラグラジアンもよぉ～く考えると、仮想静止系のエネルギーだったりします(^^)。
そういう訳で、系を解析する方法は2つある事になります。支配方程式（釣り合い方程式や運動方程式）を直接扱うか、エネルギーの最小化（変分原理）という形で陰的に解を求めるか。もちろん変分原理が不便だったら誰も使いません。変分原理による数値解法の方が、精度が良かったりコンピュータとの相性が良かったりするケースがけっこうあります。もちろん支配方程式を直接扱う方が、わかりやすいですが。
［変分法-2］
まず具体的な例を紹介します。非圧縮性渦なし完全流体のケースです。このケースでは、速度ポテンシャルφ(x，y)を計算できれば良い事になります（面倒なので2次元にします）。
流体の任意点の速度ベクトルvは
で算定でき、およそ流体力学で知りたい事の全ては、流体各点の速度だからです。このとき最小化するべきエネルギーとして用いるのは、もちろん流体の全運動エネルギーです。
ρは流体密度で、上記はちゃんと運動エネルギーの形をしてますよね？。質点の場合との違いは質量mが、密度ρ（単位体積当たりの質量）になってるところです。積分領域R内の微小面積dxdyがちっちゃいちっちゃい微小な要素の質量はρdxdyです。質量を密度で表したので、∫dxdyで積分しますよ～というだけです。積分領域Rは、解析したい流体の全体積（面積）になります。
ここで前回のxに相当するものはφです。φが系の挙動を決めるからです。そうするとエネルギーJの最小化とは、Jを最小にするような2変数関数φ(x，y)を求めよ、という話になります。この一見解くのが不可能そうな問題は、ゴールドスタインの古典力学（吉岡書店）の方法に従うと割と見通し良く解決できます。もちろん数学的には厳密ではありませんが(^^;)、読んだ中では一番わかりやすかったです。
まずベタに考えます。可能な全てのφ(x，y)に総当たりして(1)を計算し、Jが最小になった時のφを選べば良い訳ですよね？。しかし可能な全てのφ(x，y)とは、偏微分可能な任意の2変数関数全部の事です。それらは、ほぼ全ての連続関数くらいたくさん（無数に）あり総当りはとても無理ですけれど、可能な全てのφ(x，y)に番号付けする事は、数学的理想化としてはありです。例えばφ(x，y，0)，φ(x，y，1)，φ(x，y，2)，・・・みたいに(^^)。
でも全ての連続関数くらいたくさんあるので、自然数n＝0，1，2，・・・では番号の数が不足かも知れません。そこで実数αで番号付けする事を考えます。φ(x，y)をφ(x，y，α)とみなすという話ですが、これをテスト関数と呼びます。そうすると(1)は、
です。
ここで(2)の入出力関係を整理します。流体密度ρは一定（定数）です。∫dxdyの積分領域Rは固定（定数みたいなもの）です。よって(2)右辺の変数は未知関数φ(x，y，α)だけですが、φの変化は番号付けパラメータαが担います。従って積分値Jは、実数αのみの関数J(α)です（← OKですか？(^^;)）。
次にパラメータαは、人間の方で勝手に設定するものです。その設定はかなり勝手にできます。例えば(2)のJ(α)の最小は、α＝0で起きると決められます。つまり解φ(x，y)を、φ(x，y)＝φ(x，y，0)で表すという事です（任意に決められるから、本当はα≠0でも良い）。
重要なのは次の想定です。
せっかくαで番号付けしたんだから、|α|が小さいほどテスト関数φ(x，y，α)は解φ(x，y)＝φ(x，y，0)に近い、として良いだろう。
・・・そして(^^;)、
α→0の時のφ(x，y，α)のφ(x，y)への近づき方は、J(α)がα＝0で微分可能になるくらい滑らかに設定できるはずだ。
・・・と(^^;)。
以上を認めれば、やる事は一つです。「J(α)の最小を知りたいなら、J(α)をαで微分して0とおけ」。
という訳で、「あれぇ～？。不可能そうな話が、普通の微積になっちゃった(^^)」です。
(3)はいわばJ(α)の全微分を取った形です。全体が0になるのはdJ/dα＝0だからです。そうすると×dαは余計に見えますが、これを残しておかないと以後の話が続きません。そこは普通の微積とちょっとは違います(^^;)。それに(3)では、まだ絵に描いた餅です。
(3)の形に対応する(2)の右辺を想像します。(3)はαによるJ(0)からのJ(α)の値のずれを表し、dαは0からのαの微小なずれです。よって微小なdαに対応する、(2)右辺のずれが具体的にわかれば良い訳です。
でδφを定義し、変分と呼びます。つまりδφはαのずれに対応するφのずれで、δφは立派な一つの2変数関数です。
いまα＝dαは微小なので、先の想定からδφも微小です。従って(2)右辺の積分のずれは、δφ(x，y，α)を各点(x，y)ごとで考えて、
※ δφが微小という条件しか使っていません(^^)。
になります。
けっきょく最後にはδφを各点(x，y)ごとに考えるんだから、別にδφをdφと書いても本当は良かったんですよ。しかしそれではφの全微分と混同しそうなので「関数のずれを表してるんだよぉ～」と、慣習的にδφが使われます。それにδでないと、(5)の右辺が妙な感じになるし(^^;)。
意味から（定義(4)に従えば）、
は明らかなので、
を得ます。∫dxdyの積分領域Rは固定で、(6)は「定積分」です。
ここで具体的計算には役に立たないけれど注意すべきは、(5)の積分が0になる根拠こそ、(3)だという事実です。この事実は定式化上は重要なんですが、具体的計算には役立たずなので、ここで(3)はお役御免です(^^;)。
（執筆 ddt³）
【ネムネコ註】
記号「×」は、外積の記号「×」ではなく、普通の掛け算なので混同しないように！！
（４）式の変分δφは、目的とする関数（微分方程式の解）φ(x,y)とテスト関数φ(x,y,α)のズレを表していると考えればよい。
理想的には、δφ=0、つまり、φ(x,y)=φ(x,y,α)となる関数を探し出せればよいが、すべての関数を調べ尽くすことはできない。そこで、このズレδφを最小とする、φ(x,y)をよく表す近似関数φ(x,y,α)を探しだしているのだと考えると分かり良いと思う。
