メッシュ生成

前置きとして, GetFEM 用語 に関する簡単な説明をします.

GetFEM はファイル getfem/bgeot_mesh_structure.hgetfem/getfem_mesh.h で定義されるメッシュを格納する固有の構造を持っています.

このオブジェクトは,異なる次元の要素を混合しても,任意の次元に任意の要素を格納できます.

複雑なジオメトリをメッシュ化するには, GetFEM には(非常に)実験的なメッシングプロシージャしかありません.しかし,どんなフォーマットからでもメッシュを簡単に読み込むことができます(いくつかの手順は, getfem/getfem_import.h にあり,いくつかのパブリックドメインメッシュジェネレータからメッシュをロードします).

構造体 getfem::mesh は,境界や要素の集合など,メッシュの領域に関する記述を含むこともできます. これは,凸と凸面のコンテナ getfem::mesh_region を介して処理されます.メッシュの表現については, [remacle2003] を参照してください.

メッシュに要素を追加する

変数 mymesh が次のように定義されているとします.

getfem::mesh mymesh;

このメッシュに新しい要素を挿入するには,点のリストから追加する方法と既存の点のインデックスのリストから追加する2つの方法があります.

メッシュ上に新しい点を入力するには,次のメソッドを使用してください.

i = mymesh.add_point(pt);

ここで, ptbgeot::base_node 型です. インデックス i は,メッシュ上のこの点のインデックスです. すでにメッシュ内に点が存在する場合,新しい点は挿入されず,既存の点のインデックスが返されます. メッシュには,その点の次元である主次元があります. 同じメッシュ内に異なる次元の点を持つことはできません.

新しい要素をメッシュに追加する最も基本的な関数は次の通りです.

j = mymesh.add_convex(pgt, it);

これは, bgeot::pgeometric_trans 型の pgt (基本的にはbg_gt型のインスタンスへのポインタ)を持つテンプレート関数であり,既に存在する点のインデックスのリスト上のイテレータです. たとえば,3次元メッシュに新しい3角形を追加する必要がある場合,最初に3つの点のインデックスを持つ配列を定義する必要があります

std::vector<bgeot::size_type> ind(3);
ind[0] = mymesh.add_point(bgeot::base_node(0.0, 0.0, 0.0));
ind[1] = mymesh.add_point(bgeot::base_node(0.0, 1.0, 0.0));
ind[2] = mymesh.add_point(bgeot::base_node(0.0, 0.0, 1.0));

要素の追加は次のように行います.

mymesh.add_convex(bgeot::simplex_geotrans(2,1), ind.begin());

ここで bgeot::simplex_geotrans(N,1); は次元Nの単体のための通常の線形幾何変換を示します.

単体のために,以下のようにより特殊化された関数が存在します.

mymesh.add_simplex(2, ind.begin());

次の関数を使って点のリストを直接与えることも可能です.

mymesh.add_convex_by_points(pgt, itp);

ここで itp は点の配列のイテレータです. 例えば,次のような場合,

std::vector<bgeot::base_node> pts(3);
pts[0] = bgeot::base_node(0.0, 0.0, 0.0);
pts[1] = bgeot::base_node(0.0, 1.0, 0.0);
pts[2] = bgeot::base_node(0.0, 0.0, 1.0);
mymesh.add_convex_by_points(bgeot::simplex_geotrans(2,1), pts.begin());

このように使用することも可能です.

mymesh.add_simplex_by_points(2, pts.begin());

単体以外の要素についても, mymesh.add_convex_by_pointsmymesh.add_convex で適切な幾何学的変換と共に使用することが可能です.

  • bgeot::parallelepiped_geotrans(N, 1) は次元 N の平行6面体( N = 2 の場合は4辺形, N = 3 の場合は6面体,...)

  • bgeot::prism_geotrans(N, 1) は次元 N のプリズムの通常の変換を表しています(通常のプリズムは N = 3 です.一般化されたプリズムは,次元 N-1 の単体とセグメント)

特殊化された関数も存在します

mymesh.add_parallelepiped(N, it);
mymesh.add_parallelepiped_by_points(N, itp);
mymesh.add_prism(N, it);
mymesh.add_prism_by_points(N, itp);

点の配列における点の順序は,単体要素にとっては重要ではありません(あなたが単体の向きを気にする場合を除いては).他の要素については, 通常の1次要素の頂点数 (1次要素) で示される頂点の順序を守ることが重要です.

../_images/getfemuserelem.png

通常の1次要素の頂点数

より高次の変換を含む汎用的な規則は,頂点数が対応するLagrange有限要素法( 付録A.有限要素法リスト を参照)の1つに従うことであることに注意してください.

メッシュから要素を削除する

メッシュから要素を削除するには,単純に以下を使用します

mymesh.sup_convex(i);

ここで, i は要素のインデックスです.

シンプルな構造化メッシュ

平行6面体の領域の場合, getfem/getfem_regular_meshes.h で定義されている3つの関数から単体,平行6面体,またはプリズム要素を持つ構造化メッシュを得ることができます.

最も単純な関数は

void regular_unit_mesh(mesh& m, std::vector<size_type> nsubdiv,
                       bgeot::pgeometric_trans pgt, bool noised = false);

メッシュ m をシンプル/平行6面体/プリズムの規則的なメッシュで塗りつぶします( pgt の値に応じて). 各方向のセル数は nsubdiv で与えられます. 次の例では,単位正方形に2次3角形のメッシュを作成します(メッシュは後でスケーリングして変換できます)

std::vector<getfem::size_type> nsubdiv(2);
nsubdiv[0] = 10; nsubdiv[1] = 20;
regular_unit_mesh(m, nsubdiv, bgeot::simplex_geotrans(2,2));

もっと特殊な規則的なメッシュ関数も利用できます

getfem::parallelepiped_regular_simplex_mesh(mymesh, N, org, ivect, iref);
getfem::parallelepiped_regular_prism_mesh(mymesh, N, org, ivect, iref);
getfem::parallelepiped_regular_pyramid_mesh(mymesh, N, org, ivect, iref);
getfem::parallelepiped_regular_mesh(mymesh, N, org, ivect, iref);

ここで mymesh は構造化されたメッシュが構築されるメッシュ変数で, N は次元です(単体の場合は4,プリズムの場合は5,平行6面体の場合は無制限です). orgbgeot::base_node と設定し,メッシュの起点を表します. ivect は平行四辺形の領域を構築するための N ベクトルの配列上のイテレータです. iref は, 各方向の分割数を表す N の自然数の配列です.

例えば, \(10\times~10\times~10\) の単位立方体の4面体を持つメッシュを作成するには,次のように記述します

getfem::mesh mymesh;
bgeot::base_node org(0.0, 0.0, 0.0);
std::vector<bgeot::base_small_vector> vect(3);
vect[0] = bgeot::base_small_vector(0.1, 0.0, 0.0);
vect[1] = bgeot::base_small_vector(0.0, 0.1, 0.0);
vect[2] = bgeot::base_small_vector(0.0, 0.0, 0.1);
std::vector<int> ref(3);
ref[0] = ref[1] = ref[2] = 10;
getfem::parallelepiped_regular_simplex_mesh(mymesh, 3, org, vect.begin(), ref.begin());

注釈

base_nodebase_small_vector はどちらも小さなベクトルクラス(16個以上の要素を格納できません)で,幾何用の点や幾何用のベクトルを記述するために使われます. それらのメモリフットプリントは std::vector よりも小さくなっています.

メッシュ領域

メッシュオブジェクトには,多くの getfem::mesh_region オブジェクト( getfem/getfem_mesh_region.h で宣言)を含めることができます. これらのオブジェクトは,凸と凸面のセットのコンテナです. それらは,境界を定義するため,または並列ソルバーのメッシュの区画などに使用されます.

mymesh.region(30).add(2);   // adds convex 2 into region 30
mymesh.region(30).add(3);   // adds convex 3 into region 30
mymesh.region(30).add(4,3); // adds face 3 of convex 4 into region 30
mymesh.region(30).sup(3);   // Removes convex 3 from region 30
mymesh.sup_convex(4);       // Removes convex 4 from both the mesh and all the regions
for (getfem::mr_visitor i(mymesh.region(30)); !i.finished(); ++i) {
  cout << "convex: " << i.cv() << " face:" << i.f() << endl;
}

getfem::mesh オブジェクトのメソッド

リストは網羅的ではありません.

getfem::mesh::dim()

メッシュの主次元.

getfem::mesh::points_index()

メッシュの有効な点のすべてのインデックスを表す dal::bit_vector オブジェクトを返します(下記参照).

getfem::mesh::points()

各インデックスの点を与えます(bgeot::base_node).

getfem::mesh::convex_index()

メッシュの有効な要素のすべてのインデックスを表す dal::bit_vector オブジェクトを返します(下記参照).

bgeot::mesh_structure::structure_of_convex(i)

インデックス i の要素の構造の記述を与えます.この関数は bgeot::pconvex_structure を返します.

bgeot::convex_structure::nb_faces()

bgeot::pconvex_structure の面の数.

bgeot::convex_structure::nb_points()

bgeot::pconvex_structure の頂点の数.

bgeot::convex_structure::dim()

bgeot::pconvex_structure の内在次元.

bgeot::convex_structure::nb_points_of_face(f)

bgeot::pconvex_structure のローカルインデックス f の面の頂点の数.

bgeot::convex_structure::ind_points_of_face(f)

return a container with the local indexes of all vertices of the face of local index f of bgeot::pconvex_structure. For instance mesh.structure_of_convex(i)->ind_points_of_face(f)[0] is the local index of the first vertex.

bgeot::convex_structure::face_structure(f)

bgeot::pconvex_structure のローカルインデックス f の構造( bgeot::pconvex_structure )を与えます.

getfem::mesh::ind_points_of_convex(i)

bgeot::pconvex_structure の頂点のグローバルインデックスを持つコンテナを提供します.

getfem::mesh::points_of_convex(i)

bgeot::pconvex_structure の頂点を持つコンテナを返します. これは bgeot::base_node の配列です.

getfem::mesh::convex_to_point(ipt)

グローバルインデックス ipt の点にアタッチされたすべての要素のインデックスを持つコンテナを返します.

getfem::mesh::neighbors_of_convex(ic, f)

要素 ic 以外の要素 ic のローカルインデックス f に共通面を持つ mesh 内のすべての要素のインデックスを持つコンテナを返します.

getfem::mesh::neighbor_of_convex(ic, f)

要素 ic 以外の要素 ic のローカルインデックスfの共通面を持つ mesh の最初の要素のインデックスを返します.見つからなければ size_type(-1) を返します.

getfem::mesh::is_convex_having_neighbor(ic, f)

要素 ic がローカルインデックス f の面に対して隣接要素を持っているかどうかを返します.

getfem::mesh::clear()

メッシュからすべての要素と点を削除します.

getfem::mesh::optimize_structure()

構造体をコンパクトにします(番号付けに穴がないように点と凸の番号を付け直します).

getfem::mesh::trans_of_convex(i)

インデックス i の要素( bgeot::pgeometric_trans )の幾何学的変換を返します. 幾何学的変換の詳細については, プロジェクトの説明 を参照してください.

getfem::mesh::normal_of_face_of_convex(ic, f, pt)

ローカル座標(参照要素内の座標)の点でローカルインデックス f の面で要素の外向き法線を表す bgeot::base_small_vector を返します. 幾何変換が線形の場合,点 pt は影響を与えません. これは単位法線ではなく,結果として生じるベクトルのノルムは,基準要素の面の表面と実要素の面の表面との間の比です.

getfem::mesh::convex_area_estimate(ic)

ic 凸面の面積の推定値を与えます.

getfem::mesh::convex_quality_estimate(ic)

要素 ic の品質の概算を与えます.

getfem::mesh::convex_radius_estimate(ic)

要素 ic の半径の推定値を与えます.

getfem::mesh::region(irg)

getfem::mesh_region を返します. 領域はメッシュに格納され,一連の凸の番号と凸の面を含むことができます.

getfem::mesh::has_region(irg)

インデックス irg が作成されている場合はtrueを返します.

凸/凸面コンテナ getfem::mesh_region のメソッドは次のとおりです:

getfem::mesh_region::add(ic)

その領域にインデックス `` ic`` の凸を追加します.

getfem::mesh_region::add(ic, f)

ic 番目の凸の面番号 f を追加してください.

getfem::mesh_region::sup(ic)
getfem::mesh_region::sup(ic, f)

その領域から凸または凸面を除去します.

getfem::mesh_region::is_in(ic)
getfem::mesh_region::is_in(ic, f)

凸(または凸面)が領域内にある場合はtrueを返します.

getfem::mesh_region::is_only_faces()

領域に凸が含まれていない場合はtrueを返します.

getfem::mesh_region::is_only_convexes()

領域に凸面が含まれていない場合はtrueを返します.

getfem::mesh_region::index()

領域に格納されている(または面が格納されている)凸のリストを含む dal::bit_vector を返します.

getfem::mesh_region に対する反復は, getfem::mr_visitor で行う必要があります

getfem::mesh_region &rg = mymesh.region(2);
for (getfem::mr_visitor i(rg); !i.finished(); ++i) {
  cout << "contains convex " < < i.cv();
  if (i.is_face()) cout  << "face " << i.f() << endl;
}

dal::bit_vector を使います

オブジェクト dal::bit_vectorgetfem/dal_bit_vector.h で宣言されています)は, GetFEM で頻繁に使用される構造体です. これは std::bitsetstd::vector<bool> に非常に近いものですが,非負の整数の集合を表現する追加の機能を持ち,それらを繰り返し処理します.

nndal::bit_vector であると宣言されている場合,2つの命令 nn.add(6) または nn[6] = true は等価であり,整数6がその集合に追加されることを意味します.

同様に nn.sup(6) または nn[6] = false で整数6を集合から削除します. 命令nn.add(6,4)は集合に6,7,8,9を追加します.

dal::bit_vector を反復するにはイテレーターを通常どおり使用することができますが,ほとんどの場合,このオブジェクトは整数の集合を表しているため,集合に含まれる整数を反復したいだけです. これを行う最も簡単な方法は,疑似反復子 dal::bv_visitor を使用することです.

例えば,メッシュ上の点を標準出力に出力するコードを次に示します

for (dal::bv_visitor i(mymesh.points_index()); !i.finished(); ++i)
  cout << "Point of index " << i << " of the mesh: " << mymesh.points()[i] << endl;

面番号

通常の要素の面数を図 通常の要素の面に対しての記数 に示します.

../_images/getfemuserelemf.png

通常の要素の面に対しての記数

凸と点はgf_mオブジェクト内でグローバルに番号が付けられていますが,面のグローバルな番号付けはありません.そのため,与えられた面を参照する唯一の方法は凸面の番号と凸面の局所面番号を与えますことです.

メッシュの保存とロード

GetFEM ファイル形式から

getfem / getfem_mesh.h には,ファイルからメッシュをロードし,メッシュをファイルに書き込むための2つのメソッドが定義されています.

getfem::mesh::write_to_file(const std::string &name)

メッシュをファイルに保存します.

getfem::mesh::read_from_file(const std::string &name)

ファイルからメッシュを読み込みます.

メッシュを読み込んで情報を抽出する方法の例を以下に示します

#include <getfem/getfem_mesh.h>

getfem::mesh mymesh;

int main(int argc, char *argv[]) {
  try {
    // read the mesh from the file name given by the first argument
    mymesh.read_from_file(std::string(argv[1]));

    // List all the convexes
    dal::bit_vector nn = mymesh.convex_index();
    bgeot::size_type i;
    for (i << nn; i != bgeot::size_type(-1); i << nn) {
      cout << "Convex of index " << i << endl;
      bgeot::pconvex_structure cvs = mymesh.structure_of_convex(i);
      cout << "Number of vertices: " << cvs->nb_points() << endl;
      cout << "Number of faces: " << cvs->nb_faces() << endl;
      for (bgeot::short_type f = 0; f < cvs->nb_faces(); ++f) {
        cout << "face " << f << " has " << cvs->nb_points_of_face(f);
        cout << " vertices with local indexes: ";
        for (bgeot::size_type k = 0; k < cvs->nb_points_of_face(f); ++k)
          cout << cvs->ind_points_of_face(f)[k] << " ";
        cout << " and global indexes: ";
        for (bgeot::size_type k = 0; k < cvs->nb_points_of_face(f); ++k)
          cout << mymesh.ind_points_of_convex(i)[cvs->ind_points_of_face(f)[k]] << " ";
      }
    }
  } GMM_STANDARD_CATCH_ERROR; // catches standard errors
}

メッシュをインポートする

getfem/getfem_import.h は次の関数を提供します:

void import_mesh(const std::string& fmtfilename, mesh& m);

ここで,文字列 fmtfilename には,ファイル形式("gid", "gmsh", "cdb", "noboite", "am_fmt", "emc2_mesh", または"structured")の記述子と,その後にコロンとファイル名が含まれていなければなりません(書式記述子がない場合,そのファイルはネイティブのgetfemメッシュであり, mesh::read_from_file() メソッドが使用されているものとみなされます). 例えば:

getfem::mesh m;
getfem::import_mesh("gid:../tests/meshes/tripod.GiD.msh",m);

あるいは,関数

void import_mesh(const std::string& filename, const std::string& fmt,
                 mesh& m);

前述の書式指定子の1つである文字列 fmt と同等の方法で使用できます.

"gid"フォーマット指定子は, GiD によって生成されるメッシュ用であり, gmsh は,オープンソースメッシュ生成プログラム Gmsh によって生成されるメッシュ用です. "cdb" 書式指定子はCDWRITEコマンドでブロック形式で出力された ANSYS モデルからメッシュを読み込むためのものです. 現在, ANSYS 要素タイプ42,45,73,82,87,89,90,92,95,162,182,183,185,186,187,および191をインポートすることができますが,これらの要素にリンクされた有限要素技術は含まず,そのジオメトリのみが含まれます. "noboite"形式はTetMesh-GHS3D用で, "am_fmt"と "emc2_mesh"は "EMC2"(ただし2Dのみ)で構築されたファイル用です.

structured 形式は通常のメッシュの短い仕様です.その場合の fmtfilename の残りの部分はファイル名ではなく,次の形式の文字列です

getfem::import_mesh("structured:GT='GT_PK(2,1)';"
                    "NSUBDIV=[5,5];"
                    "ORG=[0,0];"
                    "SIZES=[1,1];"
                    "NOISED=0", m);

ここで GT は幾何変換の名前, NSUBDIV は各座標(デフォルト値2)の細分数のベクトル, ORG はメッシュの原点(デフォルト値 [0,0,...] ), SIZES は 各方向のサイズのベクトル(デフォルト値 [1, 1, ...]NOISED = 1 の場合は,メッシュ内部の節点がランダムに振られます(デフォルト値 NOISED = 0 ). GT 以外のパラメータはオプションです.