구조화 된 메시와 거친 및 임의의 행렬이있는 3D 포함

Nov 16 2020

예를 들어 구조화 된 3D 메쉬 (포함)를 정의하는 간단한 방법이 있는지 궁금합니다.

코스와 구조화되지 않은 매트릭스로 둘러 쌉니다. 어느 정도 독립적으로 다듬는 것이 가능해야합니다 (물론 포함이 다듬어지면 외부 행렬도 직접 포함 행렬 경계에서 더 미세해질 것입니다).

내 시도는 항상 다음과 같은 매트릭스뿐만 아니라 포함에 대해 구조화되지 않은 메쉬를 생성합니다.

마지막 메쉬를 생성하기 위해 다음 코드를 사용했습니다.

xI = 200; yI = 200; zI = 20;
InclusionRegion = 
  Region[Hexahedron[{{-xI/2, -yI/2, -zI/2}, {xI/2, -yI/2, -zI/2}, {xI/
       2, yI/2, -zI/2}, {-xI/2, yI/2, -zI/2}
     , {-xI/2, -yI/2, zI/2}, {xI/2, -yI/2, zI/2}, {xI/2, yI/2, 
      zI/2}, {-xI/2, yI/2, zI/2}}], Axes -> True];

xM = xI*2; yM = yI*2; zM = zI*2;
MatrixRegion = 
  Region[Hexahedron[{{-xM/2, -yM/2, -zM/2}, {xM/2, -yM/2, -zM/2}, {xM/
       2, yM/2, -zM/2}, {-xM/2, yM/2, -zM/2}
     , {-xM/2, -yM/2, zM/2}, {xM/2, -yM/2, zM/2}, {xM/2, yM/2, 
      zM/2}, {-xM/2, yM/2, zM/2}}], Axes -> True];

mesh = ToElementMesh[
   DiscretizeGraphics[
    RegionDifference[MatrixRegion, InclusionRegion]]
   , "RegionMarker" -> {{{0., 0., 0.}, 1, 10000}, {{xM/2, yM/2, zM/2},
       2, 1000}}
   , MaxCellMeasure -> {"Volume" -> 10000}
   , "MeshOrder" -> 1];

도움이나 제안에 감사드립니다. 미리 감사드립니다

최대

답변

12 TimLaska Nov 17 2020 at 06:20

이 답변은 @ user21을 확장하여 X, Y 및 Z 방향을 따라 포함의 다른 메쉬 밀도를 포함합니다.

현재 메셔 (버전 12.1.1)는 등방성 메시를 생성하는 것을 좋아한다는 점에 유의해야합니다. 0과 각 방향의 요소 수 범위의 매개 변수화 된 (I, J, K) 구조화 된 메시를 생성하여 다양한 메시 밀도를 달성 할 수 있습니다. 그런 다음 I, J, K 공간에서 사용자 조정 된 좌표로 좌표를 다시 조정할 수 있습니다.

먼저 등방성 구조화 된 메시를 만들어 보겠습니다.

nx = 10; ny = 40; nz = 5;
isoMesh =
   ToElementMesh[Cuboid[{0, 0, 0}, {nx, ny, nz}], 
     "MeshOrder" -> 1, MaxCellMeasure -> 1,
     "RegionMarker" -> {{{nx, ny, nz}/2, 1}}, 
     "MeshElementType" -> TetrahedronElement];
isoMesh["Wireframe"]

둘째, I, J, K 공간에서 사용자 크기 조정 좌표로 크기 조정 변환 함수를 만들어 보겠습니다.

scaledToUser = 
  RescalingTransform[{{0, nx}, {0, ny}, {0, nz}}, {{-xI/2, 
     xI/2}, {-yI/2, yI/2}, {-zI/2, zI/2}}];

이제 다음과 같이 좌표를 다시 스케일링하여 내부 메시를 만들 수 있습니다.

innerMesh = 
  ToElementMesh[
   "Coordinates" -> scaledToUser /@ isoMesh["Coordinates"], 
   "MeshElements" -> isoMesh["MeshElements"]];
innerMesh["Wireframe"]

이제의 새로운 정의로 @ user21의 워크 플로를 따르기 만하면 innermeshX, Y, Z 방향을 따라 다양한 메시 밀도를 얻을 수 있습니다.

finalMesh[
  "Wireframe"["MeshElement" -> "MeshElements", 
    "MeshElementStyle" -> (Directive[FaceForm[#], 
              EdgeForm[]] &  /@ {Orange, Blue}), 
    PlotRange -> {All, All, {-zM, zI/2}}]]
finalMesh[
  "Wireframe"["MeshElement" -> "MeshElements", 
    "MeshElementStyle" -> (Directive[FaceForm[#], 
              EdgeForm[]] &  /@ {Orange, Blue}), 
    PlotRange -> {All, {0, yI/2}, {-zM, zI/2}}]]

구조화 된 육각 메시 접근 방식

주석에서 언급했듯이 포함을 위해 구조화 된 헥스 메시를 사용하려면 현재 버전의 Mathematica 가 3D에서 피라미드 및 쐐기 유형 요소를 지원하지 않기 때문에 전체 메시를 통해 전파하고 싶을 것입니다 .

해결하려는 물리학의 특성에 따라 종종 인터페이스 영역에 급격한 기울기가있을 수 있습니다. 이 경우 도메인으로 기하 급수적으로 증가하는 인터페이스에 미세 요소 레이어가있는 경계 레이어 메시 (또는 이방성 메시)를 사용하면 솔루션에 도움이 될 수 있습니다. 이러한 유형의 메시는 요소 수 측면에서 매우 경제적 일 수 있습니다.

워크 플로우

도우미 기능

먼저 이방성 메시를 만들기 위해 몇 가지 도우미 함수를 정의합니다.

(*Import required FEM package*)
Needs["NDSolve`FEM`"];
(* Define Some Helper Functions For Structured Quad Mesh*)
pointsToMesh[data_] :=
  MeshRegion[Transpose[{data}], 
   Line@Table[{i, i + 1}, {i, Length[data] - 1}]];
unitMeshGrowth[n_, r_] := 
 Table[(r^(j/(-1 + n)) - 1.)/(r - 1.), {j, 0, n - 1}]
meshGrowth[x0_, xf_, n_, r_] := (xf - x0) unitMeshGrowth[n, r] + x0
firstElmHeight[x0_, xf_, n_, r_] := 
 Abs@First@Differences@meshGrowth[x0, xf, n, r]
lastElmHeight[x0_, xf_, n_, r_] := 
 Abs@Last@Differences@meshGrowth[x0, xf, n, r]
findGrowthRate[x0_, xf_, n_, fElm_] := 
 Quiet@Abs@
   FindRoot[firstElmHeight[x0, xf, n, r] - fElm, {r, 1.0001, 100000}, 
     Method -> "Brent"][[1, 2]]
meshGrowthByElm[x0_, xf_, n_, fElm_] := 
 N@Sort@Chop@meshGrowth[x0, xf, n, findGrowthRate[x0, xf, n, fElm]]
meshGrowthByElm0[len_, n_, fElm_] := meshGrowthByElm[0, len, n, fElm]
flipSegment[l_] := (#1 - #2) & @@ {First[#], #} &@Reverse[l];
leftSegmentGrowth[len_, n_, fElm_] := meshGrowthByElm0[len, n, fElm]
rightSegmentGrowth[len_, n_, fElm_] := Module[{seg},
  seg = leftSegmentGrowth[len, n, fElm];
  flipSegment[seg]
  ]
reflectRight[pts_] := With[{rt = ReflectionTransform[{1}, {Last@pts}]},
  Union[pts, Flatten[rt /@ Partition[pts, 1]]]]
reflectLeft[pts_] := 
 With[{rt = ReflectionTransform[{-1}, {First@pts}]},
  Union[pts, Flatten[rt /@ Partition[pts, 1]]]]
extendMesh[mesh_, newmesh_] := Union[mesh, Max@mesh + newmesh]

RegionProduct텐서 곱 메시를 사용하여 텐서 곱 그리드 생성

이제 다음 워크 플로에 표시된대로 수평, 수직 및 깊이 방향으로 서로 다른 미세 조정 전략이있는 여러 세그먼트를 붙일 수 있습니다.

(*Define parameters*)
(*Lengths*)
h = 100;(*Horizontal*)
v = 10;(*Vertical*)
d = h;(*Depth*)
(*Number of elements per segment*)
nh = 10;
nv = 10;
nd = 10;
(*Association for Clearer Region Assignment*)
reg = <|"main" -> 1, "incl" -> 2|>;
(*Create mesh segments*)
(*Horizontal segments*)
(* left segment *)
(*First element is 1/50th of seg length*)
sh = rightSegmentGrowth[h, nh, h/50];
Print["Horizontal and depth segments"]
rh = pointsToMesh@(reflectRight@reflectRight[sh] - 2 h)
(*Vertical segment*)
sv = rightSegmentGrowth[v, nv, v/50];
Print["Vertical segment"]
rv = pointsToMesh@(reflectRight@reflectRight[sv] - 2 v)
(*Create tensor product grid with RegionProduct*)
rp = RegionProduct[rh, rv, rh];
(*Show the mesh*)
Print["Highlighted RegionProduct mesh"]
HighlightMesh[rp, Style[1, Orange]]

영역 마커를 사용하여 MeshRegion을 ElementMesh로 변환

(*Extract Coords from RegionProduct*)
crd = MeshCoordinates[rp];
(*grab hexa element incidents RegionProduct mesh*)
inc = Delete[0] /@ MeshCells[rp, 3];
mesh = ToElementMesh["Coordinates" -> crd, 
   "MeshElements" -> {HexahedronElement[inc]}];
(*Extract bmesh*)
bmesh = ToBoundaryMesh[mesh];
(*Inclusion RegionMember Function*)
Ω3Dinclusion = Cuboid[{-h, -v, -h}, {h, v, h}];
rmf = RegionMember[Ω3Dinclusion];
regmarkerfn = If[rmf[#], reg["main"], reg["incl"]] &;
(*Get mean coordinate of each hexa for region marker assignment*)
mean = Mean /@ GetElementCoordinates[mesh["Coordinates"], #] & /@ 
    ElementIncidents[mesh["MeshElements"]] // First;
regmarkers = regmarkerfn /@ mean;
(*Create and view element mesh*)
Print["Converted Hexa Element Mesh Cutaway Drawing"]
mesh = ToElementMesh["Coordinates" -> mesh["Coordinates"], 
   "MeshElements" -> {HexahedronElement[inc, regmarkers]}];
mesh[
  "Wireframe"["MeshElement" -> "MeshElements", 
    "MeshElementStyle" -> (Directive[Opacity[0.5], FaceForm[#](*, 
              EdgeForm[]*)] &  /@ {Blue, Orange}), 
  ViewPoint -> {-1.5, 0.8, -3}, ViewVertical -> {0, 1, 0}, 
    PlotRange -> {{0, 2 h}, {0, 2 v}, {0, 2 h}}]]

완전히 구조화 된 헥스 메시를 사용하여 인터페이스에서 매우 세밀하게 다듬어 진 상당히 경제적 인 메시 (46656 헥스 요소)를 만들었습니다.

15 user21 Nov 16 2020 at 20:18

Acoustic Cloak Model 의 PDEModel 컬렉션에 비슷한 예가 있습니다 . 다음은 3D 버전입니다.

일부 설정 :

Needs["NDSolve`FEM`"]
xI = 200; yI = 200; zI = 20;
xM = xI*2; yM = yI*2; zM = zI*2;

내부 메쉬를 만드는 것으로 시작합니다.

innerMesh = 
 ToElementMesh[Cuboid[{-xI/2, -yI/2, -zI/2}, {xI/2, yI/2, zI/2}], 
  "MeshOrder" -> 1, MaxCellMeasure -> 10000, 
  "RegionMarker" -> {{{0., 0., 0.}, 1}}, 
  "MeshElementType" -> TetrahedronElement]

innerMesh["Wireframe"]

마커가 있는지 확인하십시오.

innerMesh["MeshElementMarkerUnion"]
{1} 

다음으로 외부 모양에 대한 경계 메쉬를 만듭니다.

bmesh1 = ToBoundaryMesh[
  Cuboid[{-xM/2, -yM/2, -zM/2}, {xM/2, yM/2, zM/2}]]

내부 메쉬에서 경계 메쉬를 추출합니다.

bmesh2 = ToBoundaryMesh[innerMesh]

FEMAddOn 을 사용하면 다음과 같이 결합 할 수 있습니다.

ResourceFunction["FEMAddOnsInstall"][]

Needs["FEMAddOns`"]
bmesh = BoundaryElementMeshJoin[bmesh1, bmesh2]

bmesh["Wireframe"]

이제 핵심 포인트입니다. 전체 외부 메쉬를 생성 할 때 경계에 새 노드가 삽입되지 않도록합니다. "SteinerPoints"-> False를 설정하면됩니다.

outerMesh = ToElementMesh[bmesh,
  "SteinerPoints" -> False,
  "RegionHoles" -> {{0, 0, 0}},
  "RegionMarker" -> {{{xM/2, yM/2, zM/2}, 2, 1000}}, 
  MaxCellMeasure -> {"Volume" -> 10000}, "MeshOrder" -> 1]

이제 내부 재질 영역에 정렬되는 내부 및 외부 메쉬가 있으므로 최종 전체 메쉬를 만들 수 있습니다.

innerCoordinates = innerMesh["Coordinates"];
outerCoordinates = outerMesh["Coordinates"];
finalMesh = 
 ToElementMesh[
  "Coordinates" -> Join[outerCoordinates, innerCoordinates], 
  "MeshElements" -> 
   Flatten[{outerMesh["MeshElements"], 
     MapThread[
      TetrahedronElement, {ElementIncidents[
         innerMesh["MeshElements"]] + Length[outerCoordinates], 
       ElementMarkers[innerMesh["MeshElements"]]}]}]]

마커가 있는지 확인하십시오.

finalMesh["MeshElementMarkerUnion"]
{1,2}

그리고 시각화 :

finalMesh[
 "Wireframe"["MeshElement" -> "MeshElements", 
  "MeshElementStyle" -> (Directive[FaceForm[#], 
       EdgeForm[]] & /@ {Orange, Blue}), 
  PlotRange -> {All, All, {-zM, zI/2}}]]

2 차 메시를 생성하려면 다음과 같이 할 수 있습니다.

MeshOrderAlteration[finalMesh, 2]