1950 年代に構造解析から誕生した「有限要素法」は偏微分方程式の標準的な数値解法の 一つとして定着しましたが、
ところが、 FreeFEM という見事なアイデアで、この状況はずいぶん改善されたと思います。 このプログラムは、フランスの UPMC (当時 Pierre et Marie Curie 大学)の スタッフによって開発された もので、Gfemと呼ばれるPascal風の言語により、
1) 境界形状 2) 境界条件 3) 偏微分方程式を記述するだけで、メッシュの自動生成から数値解を求める作業のすべてをやって くれますから、何とも気軽にいろいろな問題を解くことができます。
例えば、長径 2、短径 1 の楕円膜全体に垂直方向の力 1 を加えたときの変位を求 める問題は、
border(1,0,2*pi,40) { /* 境界番号 1、パラメータは 0 から 2*pi まで */ x:=2*cos(t); /* 頂点数 40 でパラメトリックに楕円を描く */ y:=sin(t); }; buildmesh(800); /* メッシュ生成 (最大頂点数 800) */ solve(v) { /* 偏微分方程式の指定と求解 */ onbdy(1) v=0; /* 境界条件 */ pde(v) -laplace(v) = 1 /* ポアソン方程式 */ }; plot(u); /* 等ポテンシャル線の作図 */といった極めて簡単な記述で済んでしまいます。
ご興味を持たれた方は、ぜひ前記のURLをお尋ねになってください。C++ で 書かれたソースはLinuxやFreeBSDでも簡単にコンパイルできますし、 Mac や Windows ならバイナリも用意されています。私は Windows を使いませんが、 Makefile を追加する程度で、Windows95 の Visual C++ 4.0 でもコンパイルでき ました。しかし、Unix のほうが圧倒的に使いやすいです。
使い方は、付属のドキュメントと実例をもとに、少し実験するとわかりますが、 最近、出版された、フランス語からの英訳
Brigitte Lucquin, Oliver Pironneau,- Introduction to Scientific Computing (John Wiley & Sons) ISBN 0-471-97266-5の1章にも、かなり多くのサンプルが載っています。下記は、その中の1つ、 「1.8.2 Lame's system」のプログラムと出力例ですが、中央に割け目のある 方形断面ビームの上下の面に引き離す方向の力を加えたときの歪み解析です。 縦方向の面は自由に動けます。対称性を利用して、右上半分の 1/4 を解析して います。
/* Lame's system */ border(1,0,1,20) { x:=1-(t-1)^2; y:=0 }; border(2,1,2,20) { x:=(t-1)^2+1; y:=0 }; border(3,0,1,5) { x:=2; y:=t }; border(4,0,2,10) { x:=2-t; y:=1 }; border(5,0,1,5) { x:=0; y:=1-t }; buildmesh(500); E:=2.15; sigma:=0.29; lambda:=E/(2*(1+sigma)); nu:=E*sigma/((1+sigma)*(1-2*sigma)); mu:=nu+sigma; solve(u,v) { onbdy(3,5) u=0; onbdy(2) v=0; onbdy(4) dnu(v)=1; pde(u) -laplace(u)*mu-dxx(u)*nu-dxy(v)*nu=0; pde(v) -laplace(v)*mu-dyx(u)*nu-dyy(v)*nu=0; }; plotmesh; x=x+0.1*u; y=y+0.1*v; plotmesh;
なお、このすばらしいアイデアを日本でも気軽に使えるようにしたいと思って、 日本語のドキュメント、WEB 版のデモ、気づいた範囲でのパッチを用意してみました。
なお、開発元では、Freefem++、 FreeFEM++、 FreeFEM3Dといった機能強化した新しいバージョンの開発が進んでいて、 http://www.freefem.org/ff++/ で Unix 用のソース、ドキュメント、教材が公開されておます。
当初の FreeFEM++ では complex が使えないとか、固有値問題が扱えないとか、 いろいろ不便なところがありましたが、急速に改良が進み、格段の進化を遂げました。
一時は FreeBSD でも ports が用意されたのですが、残念ながら maintainer が見付からないようで、すぐになくなってしまいました。
以下、FreeBSD に於けるコンパイルで苦労される方向けに、 FreeBSD-11.2 で私自身がコンパイルしたときのメモを付けておきます。
freefem++ が Linux で開発されていることを念頭に置いて compile します。 FreeBSD の compiler は gcc から clang に変わっていますが、 Linux は gcc で、freefem++ も gcc, g++, gfortran を使おうとしますから、 注意が必要で、下記の手順でも FreeFem++ の描画機能の一部を受け持つ ports/libGLU が ライブラリの不整合で Segmentation falut になって使えません。 その対策として、最後に gcc, g++ で compile した libGLU.a を static link した ffglut を作っています。
なお、lang/gcc(gfortran), shells/bash, ftp/wget, lang/perl15x, graphics/freeglut は前もって ports から install しておいてください。
freefem++-3.61-1.tar.gz を展開したディレクトリで下記の手順でコンパイルします。
0) freefem++-3.61-1.tar.gz を展開し src/femlib/CheckPtr.cpp にパッチをあてる --- src/femlib/CheckPtr.cpp.orig +++ src/femlib/CheckPtr.cpp @@ -27,6 +27,8 @@ */ #if __APPLE__ #include+#elif __FreeBSD__ +#include #else #include #endif 1) 環境を Linux に合わせる # ln -s /usr/local/bin/wget /usr/bin/wget # ln -s /usr/local/bin/perl /usr/bin/perl # ln -s /usr/local/bin/bash /bin/bash # mv /usr/bin/make /usr/bin/make.bsd # ln -s /usr/local/bin/gmake /usr/bin/make 2) freefem++-3.61-1 を compile する $ cd freefem++-3.61-1 $ ./configure --enable-download -without-mpi $ make 3) gcc, g++ で compile した libGLU.a を作る freefem++-3.61-1.tar.gz を展開したディレクトリに glu-9.0.0.tar.gz を展開する $ cd freefem++-3.61-1/glu-9.0.0 $ ./configure CC=gcc CXX=g++ PREFIX=`pwd` $ make prefix=`pwd` install 4) ffglut を作る (ftp://ftp.freedesktop.org/pub/mesa/glu/glu-9.0.0.tar.bz2) $ cd freefem++-3.61-1/src/nw $ GLU=../../glu-9.0.0/lib $ g++ -g -DNDEBUG -O3 -DBAMG_LONG_LONG -DNCHECKPTR -fPIC -rdynamic -o ffglut ffglut.o gggg.o ffthreads.o fem.o Mesh3dn.o Mesh2dn.o Mesh1dn.o GQuadTree.o FQuadTree.o Drawing.o mshptg.o ffapi.o ../libMesh/libMesh.a $GLU/libGLU.a -lglut -lGL -lpthread -lblas -ldl -lm -lrt -lgfortran -L/usr/local/lib -lm -ldl -lz -lsz -lhdf5 -lhdf5_hl 5) install する # cd freefem++-3.61-1 # make install # make clean # 動作確認後 6) 環境を BSD に戻す # rm /usr/bin/wget /usr/bin/perl /bin/bash # rm /usr/bin/make # ln -s /usr/local/bin/make.bsd /usr/bin/make
この他、flang が使える場合は gcc/gfortran を使わずにコンパイルすることができます。
なお、FreeFEM++では、解くべき問題を記述する微分方程式として、 導関数の階数が低いWeak Form(弱形式)だけを使うようになりましたので、 偏微分方程式の記述が以前の FreeFEM と違っていることに注意してください。
Weak Form については、境界条件の記述を含めて、 全面改定されたマニュアル(freefem++doc.pdf)にも、 豊富な例題と共に解説がありますし、前期の書籍にも書かれていますが、 日本語で読めるものとしては、下記の記述も簡潔で良いと思います。
菊地文雄,- 有限要素法概説 (サイエンス社) ISBN4-7819-0308-8
ここでは、 円柱導体の表皮効果を計算する例を1つだけつけておきます。
2010-02-01 追記
FreeFem++ も実用レベルになりましたし、FreeBSD でも ports/math/freefem++ に入って きましたので、この解説も、もう要らないと思ったのですが、FreeBSD-8.0 になった ら、 ports/math/freefem++ がなくなって(注1)、ports/math/freefem (FreeFEM-3.5.8) だけ だけになってしまいましたし、偏微分方程式が Strong form (強形式) で記述できると、 数値計算に馴染のない方にはわかりやすいと思いますので、freefem 3.5 の参考資料と して、当面残すことにしました。
私が自社の問題に適用した際に、必要と感じた機能を追加した FreeFem-3.4 のソースコ ードは技術者向け教材 有限要素法によるケーブル(電線)の電気設計 に入れてあります。
平林浩一, 1998-08-26, 2006-02-25, 2010-03-01, 2010-08-04, 2011-03-22, 2011-11-09, 2012-06-11, 2013-04-12, 2018-08-22, 2019-11-30