just to sum up:
each single part can be meshed successfully,
but as soon as you put all parts in one file,
gmsh will report error of self intersection.
could this be a bug?
-----Original Messages-----
*From:*"蓝颖杰" <y...@phbs.pku.edu.cn>
*Sent Time:*2017-02-15 13:00:14 (Wednesday)
*To:* g...@geuz.org
*Cc:*
*Subject:* meshing trouble
Hi,
I have fought with this issue for a long time, but still can't resolve it,
any help would be appreciated. I have solid objects in the air
(with spherical out surface) that needed to be excluded. Unfortunately the
solid objects (a frame and two magnets) touch each other. So I
made a surface loop to include all the solids to define the air body. Gmsh
can't mesh the model, complaining about intersecting surfaces. If
you delete the frame (see code comment on how to delete), everything else
being the same, Gmsh can mesh it without any problems.
Very strange!
Regards,
Yingjie
=====================Code starts from here (also in
attachment)========================
// define geometry-specific parameters
mm = 1.e-3;
gap = 1*mm; //Magnet-Magnet gap [m], gap > 0
DefineConstant[
cub = {50*mm, Name "Parameters/1Magnet bottom size [m]"}
hite = {10*mm, Name "Parameters/2Magnet hieght [m]"}
fdep = {60*mm, Name "Parameters/3Frame depth [m]"}
inf = {100*mm, Name "Parameters/4Air box distance [m]"}
lc1 = {TotalMemory <= 2048 ? 30*mm : 10*mm, Name "Parameters/5Mesh size on
magnets [m]"}
lc2 = {TotalMemory <= 2048 ? 90*mm : 20*mm, Name "Parameters/6Mesh size at
infinity [m]"}
];
// change global Gmsh options
Mesh.Optimize = 1; // optimize quality of tetrahedra
Mesh.VolumeEdges = 0; // hide volume edges
Geometry.ExactExtrusion = 0; // to allow rotation of extruded shapes
Solver.AutoMesh = 2; // always remesh if necessary (don't reuse mesh on
disk)
//parameters used by this macro
// cx, cy, cz: center of sphere
// rad: radius of sphere
// slc: mesh size
Macro MakeCircleZ //facing z direction
cp = newp; Point(cp) = {cx, cy, cz, slc}; //center
dd[] = {rad, 0, -rad, 0};
For i In {0:3}
spoint[i] = newp;
Point(spoint[i]) = {cx+dd[i], cy+dd[3-i], cz, slc};
If (i>0)
quarc[i-1] = newreg; //quarter arcs
Circle(quarc[i-1]) = {spoint[i-1], cp, spoint[i]};
EndIf
EndFor
//closing arc
quarc[3] = newreg; Circle(quarc[3]) = {spoint[3], cp, spoint[0]};
zloop = newreg; Line Loop(zloop) = {quarc[]};
zcircle = newreg; Plane Surface(zcircle) = {zloop};
Return
cx = 0; cy=0; cz = -gap-hite;
rad = cub; slc = lc1;
Call MakeCircleZ;
mag1[] = Extrude {0, 0, hite} { Surface{zcircle}; };
jfl1 = zloop; //joint face loop 1
cx = 0; cy=0; cz = gap+hite;
rad = cub; slc = lc1;
Call MakeCircleZ;
mag2[] = Extrude {0, 0, -hite} { Surface{zcircle}; };
jfl2 = zloop; //joint face 2
//Delete {Volume{mag1[1]};}
//Delete {Volume{mag2[1]};}
Physical Volume("Magnets") = {mag1[1], mag2[1]};
//create steel frame around the magnet
p1 = newp; Point(p1) = {-2*cub, -fdep, -hite-gap, lc1};
p2 = newp; Point(p2) = { 2*cub, -fdep, -hite-gap, lc1};
p3 = newp; Point(p3) = { 2*cub, -fdep, hite+gap, lc1};
p4 = newp; Point(p4) = {-2*cub, -fdep, hite+gap, lc1};
l1 = newl; Line(l1) = {p1,p2}; l2 = newl; Line(l2) = {p2,p3};
l3 = newl; Line(l3) = {p3,p4}; l4 = newl; Line(l4) = {p4,p1};
q1 = newp; Point(q1) = {-2*cub, fdep, -hite-gap, lc1};
q2 = newp; Point(q2) = { 2*cub, fdep, -hite-gap, lc1};
q3 = newp; Point(q3) = { 2*cub, fdep, hite+gap, lc1};
q4 = newp; Point(q4) = {-2*cub, fdep, hite+gap, lc1};
e1 = newl; Line(e1) = {q1,q2}; e2 = newl; Line(e2) = {q2,q3};
e3 = newl; Line(e3) = {q3,q4}; e4 = newl; Line(e4) = {q4,q1};
d1 = newl; Line(d1) = {p1,q1}; d2 = newl; Line(d2) = {p2,q2};
d3 = newl; Line(d3) = {p3,q3}; d4 = newl; Line(d4) = {p4,q4};
ll1 = newll; Line Loop(ll1) = {l1,l2,l3,l4};
lex1 = newll; Line Loop(lex1) = {l1, d2, -e1, -d1};
fex1 = news; Plane Surface(fex1) = {lex1};
hof1 = news; Plane Surface(hof1) = {lex1, jfl1}; //holed face 1
lex2 = newll; Line Loop(lex2) = {l2, d3, -e2, -d2};
fex2 = news; Plane Surface(fex2) = {lex2};
lex3 = newll; Line Loop(lex3) = {l3, d4, -e3, -d3};
fex3 = news; Plane Surface(fex3) = {lex3};
hof2 = news; Plane Surface(hof2) = {lex3, jfl2}; //holed face 2
lex4 = newll; Line Loop(lex4) = {l4, d1, -e4, -d4};
fex4 = news; Plane Surface(fex4) = {lex4};
ll1d = newll; Line Loop(ll1d) = {e1,e2,e3,e4};
hite2 = hite + gap + 1.5*cub;
p1 = newp; Point(p1) = {-3.5*cub, -fdep, -hite2, lc1};
p2 = newp; Point(p2) = { 3.5*cub, -fdep, -hite2, lc1};
p3 = newp; Point(p3) = { 3.5*cub, -fdep, hite2, lc1};
p4 = newp; Point(p4) = {-3.5*cub, -fdep, hite2, lc1};
l1 = newl; Line(l1) = {p1,p2}; l2 = newl; Line(l2) = {p2,p3};
l3 = newl; Line(l3) = {p3,p4}; l4 = newl; Line(l4) = {p4,p1};
q1 = newp; Point(q1) = {-3.5*cub, fdep, -hite2, lc1};
q2 = newp; Point(q2) = { 3.5*cub, fdep, -hite2, lc1};
q3 = newp; Point(q3) = { 3.5*cub, fdep, hite2, lc1};
q4 = newp; Point(q4) = {-3.5*cub, fdep, hite2, lc1};
e1 = newl; Line(e1) = {q1,q2}; e2 = newl; Line(e2) = {q2,q3};
e3 = newl; Line(e3) = {q3,q4}; e4 = newl; Line(e4) = {q4,q1};
d1 = newl; Line(d1) = {p1,q1}; d2 = newl; Line(d2) = {p2,q2};
d3 = newl; Line(d3) = {p3,q3}; d4 = newl; Line(d4) = {p4,q4};
ll2 = newll; Line Loop(ll2) = {l1,l2,l3,l4};
let1 = newll; Line Loop(let1) = {l1, d2, -e1, -d1};
fet1 = news; Plane Surface(fet1) = {let1};
let2 = newll; Line Loop(let2) = {l2, d3, -e2, -d2};
fet2 = news; Plane Surface(fet2) = {let2};
let3 = newll; Line Loop(let3) = {l3, d4, -e3, -d3};
fet3 = news; Plane Surface(fet3) = {let3};
let4 = newll; Line Loop(let4) = {l4, d1, -e4, -d4};
fet4 = news; Plane Surface(fet4) = {let4};
ll2d = newll; Line Loop(ll2d) = {e1,e2,e3,e4};
s1 = news; Plane Surface(s1) = {ll2, ll1};
s1d = news; Plane Surface(s1d) = {ll2d, ll1d};
//the frame, comment out the next five lines to delete it
slframe = newsl;
Surface Loop(slframe)={ s1, fex1, fex2, fex3, fex4,
s1d, fet1, fet2, fet3, fet4};
frame = newv; Volume(frame) = {slframe};
Physical Volume("Frame") = {frame};
//now create the surface enclosing both mags and the frame
solids = newsl;
Surface Loop(solids) = { s1d, fet1, fet2, fet3, fet4, s1,
hof1, mag1[2], mag1[3], mag1[4], mag1[5], mag1[0], fex2,
hof2, mag2[2], mag2[3], mag2[4], mag2[5], mag2[0], fex4};
//parameters used by this macro
// cx, cy, cz: center of sphere
// rad: radius of sphere
// slc: mesh size of sphere
Macro MakeSphereSurface
cp = newp; Point(cp) = {cx, cy, cz, slc}; //center
np = newp; Point(np) = {cx, cy, cz+rad, slc}; //north
sp = newp; Point(sp) = {cx, cy, cz-rad, slc}; //south
dd[] = {rad, 0, -rad, 0};
For i In {0:3}
ap = newp; pe[i] = ap; //points on the equator
Point(ap) = {cx+dd[i], cy+dd[3-i], cz, slc};
na[i] = newreg; Circle(na[i]) = {ap, cp, np}; //north arc
sa[i] = newreg; Circle(sa[i]) = {ap, cp, sp}; //south arc
If (i>0)
eq[i-1] = newreg; //equator
Circle(eq[i-1]) = {pe[i-1], cp, ap};
EndIf
EndFor
//closing arc of the equator
eq[3] = newreg; Circle(eq[3]) = {ap, cp, pe[0]};
//now the non-plane surfaces: sloop[]
For i In {0:3}
lloop = newreg; Line Loop(lloop)={eq[i], -na[i], na[(i+1)%4]};
sloop[2*i] = newreg; Ruled Surface(sloop[2*i]) = {lloop};
lloop = newreg; Line Loop(lloop)={eq[i], -sa[i], sa[(i+1)%4]};
sloop[2*i+1] = newreg; Ruled Surface(sloop[2*i+1]) = {lloop};
EndFor
thesurf = newreg;
Surface Loop(thesurf) = {sloop[]};
Return
// create air box around magnets
BoundingBox; // recompute model bounding box
cx = (General.MinX + General.MaxX) / 2;
cy = (General.MinY + General.MaxY) / 2;
cz = (General.MinZ + General.MaxZ) / 2;
lx = 2*inf + General.MaxX - General.MinX;
ly = 2*inf + General.MaxY - General.MinZ;
lz = 2*inf + General.MaxZ - General.MinZ;
rad = (lx + ly + lz)/3;
slc = lc2;
Call MakeSphereSurface ;
vair = newv; Volume(vair) = {thesurf, solids};
Physical Volume("Air") = {vair}; // air
Physical Surface("Infty") = sloop[];
// adjust value here for correct merge result
//Geometry.Tolerance = 1e-6;
//Coherence Mesh;
_______________________________________________
gmsh mailing list
gmsh@onelab.info
http://onelab.info/mailman/listinfo/gmsh
