segfault while using netgen optimizer caused by `Surface In Volume`
Hi all,
I have to following setup (see below). Netgen
optimizer gives me the best mesh quality for my scenarios and I would like to keep using it. However, in some case, I need to have a surface submerged into a volume in order to setup refinement along it. Unfortunately, I get segfaults once the run time reaches netgen
. I noticed that everything works if I comment Surface{6} In Volume{32};
line. I don't know what is the exact reason. I tried different versions of gmsh and netgen installing with Spack. So, Im sure that I didn't mess up with linking wrong libraries.
I am sure that you (developers) have much better experience with it. Could you, please, check this setup? You can downscale the problem to generate a coarse mesh? I don't know but maybe I am missing something or there is a hack. Or maybe there is a version of gmsh (or netgen) which can easily handle this case.
Thanks in advance.
lc = 6e3;
lc_fault = 100;
faultlength = 20e3;
faultwidth = 25e3 ;
region = 200e3 ;
depth = 80e3;
dip = 40*Pi/180;
stk = 122*Pi/180;
top_depth = 5e3+0.3*lc_fault;
x1=0.5*faultlength*Sin(stk);
y1=0.5*faultlength*Cos(stk);
Point(1) = {0.5*faultlength*Cos(stk), -0.5*faultlength*Sin(stk), -5e3, lc_fault};
Point(4) = {-0.5*faultlength*Cos(stk), 0.5*faultlength*Sin(stk), -5e3, lc_fault};
Point(5) = {-0.5*faultlength*Cos(stk)-faultwidth*Sin(dip)*Sin(stk), 0.5*faultlength*Sin(stk)-faultwidth*Cos(stk)*Sin(dip), -(5e3+1*faultwidth)*Cos(dip), lc_fault};
Point(6) = {0.5*faultlength*Cos(stk)-faultwidth*Sin(dip)*Sin(stk), -0.5*faultlength*Sin(stk)-faultwidth*Cos(stk)*Sin(dip), -(5e3+1*faultwidth)*Cos(dip), lc_fault};
Point(7) = {0.5*region, 0.5*region, -depth, lc};
Point(8) = {0.5*region, -0.5*region, -depth, lc};
Point(9) = {-0.5*region, -0.5*region, -depth, lc};
Point(10) = {-0.5*region, 0.5*region, -depth, lc};
Point(11) = {-0.5*region, 0.5*region, 0, lc};
Point(13) = {0.5*region, 0.5*region, 0, lc};
Point(14) = {0.5*region, -0.5*region, 0, lc};
Point(15) = {-0.5*region, -0.5*region, 0, lc};
Line(1) = {4, 5};
Line(2) = {5, 6};
Line(3) = {6, 1};
Line(4) = {1, 4};
Line(7) = {7, 8};
Line(8) = {8, 9};
Line(9) = {9, 10};
Line(10) = {10, 7};
Line(11) = {7, 13};
Line(12) = {13, 14};
Line(13) = {14, 15};
Line(14) = {15, 9};
Line(15) = {15, 11};
Line(16) = {11, 10};
Line(17) = {11, 13};
Line(18) = {8, 14};
Line Loop(6) = {4, 1, 2, 3};
Ruled Surface(6) = {6};
Line Loop(20) = {8, 9, 10, 7};
Plane Surface(20) = {20};
Line Loop(22) = {9, -16, -15, 14};
Plane Surface(22) = {22};
Line Loop(24) = {10, 11, -17, 16};
Plane Surface(24) = {24};
Line Loop(26) = {7, 18, -12, -11};
Plane Surface(26) = {26};
Line Loop(28) = {8, -14, -13, -18};
Plane Surface(28) = {28};
Line Loop(30) = {13, 15, 17, 12};
Plane Surface(30) = {30};
//Line{4} In Suface{30};
Surface Loop(32) = {26, 20, 28, 22, 24, 30};
Volume(32) = {32};
Surface{6} In Volume{32}; // <<<--- THIS LINE
Field[1] = Box;
Field[1].VIn = 4 * lc_fault;
Field[1].VOut = 0.5*lc;
Field[1].XMax = 20e3;
Field[1].XMin = -20e3;
Field[1].YMax = 25e3;
Field[1].YMin = -15e3;
Field[1].ZMax = 0;
Field[1].ZMin = -30e3;
Field[1].Thickness = 45e3;
Field[2] = Distance;
Field[2].FacesList = {6};
Field[2].FieldX = -1;
Field[2].FieldY = -1;
Field[2].FieldZ = -1;
Field[2].NNodesByEdge = 20;
Field[3] = MathEval;
Field[3].F = Sprintf("0.3*F2 +(F1/5000)^2 + %g", lc_fault);
Field[4] = Min;
Field[4].FieldsList = {1, 3};
Background Field = 4;
Physical Surface(101) = {30};
Physical Surface(105) = {20, 22, 24, 26, 28};
// Physical Surface(103) = {6};
Physical Volume(32) = {32};
Mesh.MshFileVersion = 2.2;