Hi, <div><br></div><div>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. <div><br></div><div>Very strange!</div></div><div><br></div><div>Regards,</div><div>Yingjie</div><div><br></div><div>=====================Code starts from here (also in attachment)========================</div><div><div><br></div><div>// define geometry-specific parameters</div><div>mm = 1.e-3;</div><div>gap = 1*mm; //Magnet-Magnet gap [m], gap > 0</div><div><br></div><div>DefineConstant[</div><div>  cub = {50*mm, Name "Parameters/1Magnet bottom size [m]"}</div><div>  hite = {10*mm, Name "Parameters/2Magnet hieght [m]"}</div><div>  fdep = {60*mm, Name "Parameters/3Frame depth [m]"}</div><div>  inf = {100*mm, Name "Parameters/4Air box distance [m]"}</div><div>  lc1 = {TotalMemory <= 2048 ? 30*mm : 10*mm, Name "Parameters/5Mesh size on magnets [m]"}</div><div>  lc2 = {TotalMemory <= 2048 ? 90*mm : 20*mm, Name "Parameters/6Mesh size at infinity [m]"}</div><div>];</div><div><br></div><div>// change global Gmsh options</div><div>Mesh.Optimize = 1; // optimize quality of tetrahedra</div><div>Mesh.VolumeEdges = 0; // hide volume edges</div><div>Geometry.ExactExtrusion = 0; // to allow rotation of extruded shapes</div><div>Solver.AutoMesh = 2; // always remesh if necessary (don't reuse mesh on disk)</div><div><br></div><div>//parameters used by this macro</div><div>// cx, cy, cz: center of sphere</div><div>// rad: radius of sphere</div><div>// slc: mesh size</div><div><br></div><div>Macro MakeCircleZ //facing z direction</div><div>  cp = newp; Point(cp) = {cx, cy, cz, slc}; //center</div><div>  dd[] = {rad, 0, -rad, 0};</div><div>  For i In {0:3}</div><div>    spoint[i] = newp;</div><div><span class="Apple-tab-span" style="white-space:pre">    </span>Point(spoint[i]) = {cx+dd[i], cy+dd[3-i], cz, slc};</div><div><span class="Apple-tab-span" style="white-space:pre">  </span>If (i>0)</div><div><span class="Apple-tab-span" style="white-space:pre">  </span>   quarc[i-1] = newreg; //quarter arcs</div><div><span class="Apple-tab-span" style="white-space:pre">  </span>   Circle(quarc[i-1]) = {spoint[i-1], cp, spoint[i]};</div><div><span class="Apple-tab-span" style="white-space:pre">   </span>EndIf</div><div>  EndFor</div><div>  //closing arc</div><div>  quarc[3] = newreg; Circle(quarc[3]) = {spoint[3], cp, spoint[0]};</div><div>  zloop = newreg; Line Loop(zloop) = {quarc[]};</div><div>  zcircle = newreg; Plane Surface(zcircle) = {zloop};</div><div>Return</div><div><br></div><div>cx = 0; cy=0; cz = -gap-hite;</div><div>rad = cub; slc = lc1;</div><div>Call MakeCircleZ;</div><div>mag1[] = Extrude {0, 0, hite} { Surface{zcircle}; };</div><div>jfl1 = zloop; //joint face loop 1</div><div><br></div><div>cx = 0; cy=0; cz = gap+hite;</div><div>rad = cub; slc = lc1;</div><div>Call MakeCircleZ;</div><div>mag2[] = Extrude {0, 0, -hite} { Surface{zcircle}; };</div><div>jfl2 = zloop; //joint face 2</div><div><br></div><div>//Delete {Volume{mag1[1]};}</div><div>//Delete {Volume{mag2[1]};}</div><div>Physical Volume("Magnets") = {mag1[1], mag2[1]};</div><div><br></div><div>//create steel frame around the magnet</div><div>p1 = newp; Point(p1) = {-2*cub, -fdep, -hite-gap, lc1};</div><div>p2 = newp; Point(p2) = { 2*cub, -fdep, -hite-gap, lc1};</div><div>p3 = newp; Point(p3) = { 2*cub, -fdep,  hite+gap, lc1};</div><div>p4 = newp; Point(p4) = {-2*cub, -fdep,  hite+gap, lc1};</div><div>l1 = newl; Line(l1) = {p1,p2}; l2 = newl; Line(l2) = {p2,p3};</div><div>l3 = newl; Line(l3) = {p3,p4}; l4 = newl; Line(l4) = {p4,p1};</div><div><br></div><div>q1 = newp; Point(q1) = {-2*cub, fdep, -hite-gap, lc1};</div><div>q2 = newp; Point(q2) = { 2*cub, fdep, -hite-gap, lc1};</div><div>q3 = newp; Point(q3) = { 2*cub, fdep,  hite+gap, lc1};</div><div>q4 = newp; Point(q4) = {-2*cub, fdep,  hite+gap, lc1};</div><div>e1 = newl; Line(e1) = {q1,q2}; e2 = newl; Line(e2) = {q2,q3};</div><div>e3 = newl; Line(e3) = {q3,q4}; e4 = newl; Line(e4) = {q4,q1};</div><div><br></div><div>d1 = newl; Line(d1) = {p1,q1}; d2 = newl; Line(d2) = {p2,q2};</div><div>d3 = newl; Line(d3) = {p3,q3}; d4 = newl; Line(d4) = {p4,q4};</div><div><br></div><div>ll1 = newll; Line Loop(ll1) = {l1,l2,l3,l4};</div><div><br></div><div>lex1 = newll; Line Loop(lex1) = {l1, d2, -e1, -d1};</div><div>fex1 = news; Plane Surface(fex1) = {lex1}; </div><div>hof1 = news; Plane Surface(hof1) = {lex1, jfl1}; //holed face 1</div><div><br></div><div>lex2 = newll; Line Loop(lex2) = {l2, d3, -e2, -d2};</div><div>fex2 = news; Plane Surface(fex2) = {lex2}; </div><div><br></div><div>lex3 = newll; Line Loop(lex3) = {l3, d4, -e3, -d3};</div><div>fex3 = news; Plane Surface(fex3) = {lex3}; </div><div>hof2 = news; Plane Surface(hof2) = {lex3, jfl2}; //holed face 2</div><div><br></div><div>lex4 = newll; Line Loop(lex4) = {l4, d1, -e4, -d4};</div><div>fex4 = news; Plane Surface(fex4) = {lex4}; </div><div><br></div><div>ll1d = newll; Line Loop(ll1d) = {e1,e2,e3,e4};</div><div><br></div><div>hite2 = hite + gap + 1.5*cub;</div><div>p1 = newp; Point(p1) = {-3.5*cub, -fdep, -hite2, lc1};</div><div>p2 = newp; Point(p2) = { 3.5*cub, -fdep, -hite2, lc1};</div><div>p3 = newp; Point(p3) = { 3.5*cub, -fdep,  hite2, lc1};</div><div>p4 = newp; Point(p4) = {-3.5*cub, -fdep,  hite2, lc1};</div><div>l1 = newl; Line(l1) = {p1,p2}; l2 = newl; Line(l2) = {p2,p3};</div><div>l3 = newl; Line(l3) = {p3,p4}; l4 = newl; Line(l4) = {p4,p1};</div><div><br></div><div>q1 = newp; Point(q1) = {-3.5*cub, fdep, -hite2, lc1};</div><div>q2 = newp; Point(q2) = { 3.5*cub, fdep, -hite2, lc1};</div><div>q3 = newp; Point(q3) = { 3.5*cub, fdep,  hite2, lc1};</div><div>q4 = newp; Point(q4) = {-3.5*cub, fdep,  hite2, lc1};</div><div>e1 = newl; Line(e1) = {q1,q2}; e2 = newl; Line(e2) = {q2,q3};</div><div>e3 = newl; Line(e3) = {q3,q4}; e4 = newl; Line(e4) = {q4,q1};</div><div><br></div><div>d1 = newl; Line(d1) = {p1,q1}; d2 = newl; Line(d2) = {p2,q2};</div><div>d3 = newl; Line(d3) = {p3,q3}; d4 = newl; Line(d4) = {p4,q4};</div><div><br></div><div>ll2 = newll; Line Loop(ll2) = {l1,l2,l3,l4};</div><div><br></div><div>let1 = newll; Line Loop(let1) = {l1, d2, -e1, -d1};</div><div>fet1 = news; Plane Surface(fet1) = {let1}; </div><div><br></div><div>let2 = newll; Line Loop(let2) = {l2, d3, -e2, -d2};</div><div>fet2 = news; Plane Surface(fet2) = {let2}; </div><div><br></div><div>let3 = newll; Line Loop(let3) = {l3, d4, -e3, -d3};</div><div>fet3 = news; Plane Surface(fet3) = {let3}; </div><div><br></div><div>let4 = newll; Line Loop(let4) = {l4, d1, -e4, -d4};</div><div>fet4 = news; Plane Surface(fet4) = {let4}; </div><div><br></div><div>ll2d = newll; Line Loop(ll2d) = {e1,e2,e3,e4};</div><div><br></div><div>s1 = news; Plane Surface(s1) = {ll2, ll1};</div><div>s1d = news; Plane Surface(s1d) = {ll2d, ll1d};</div><div><br></div><div>//the frame, comment out the next five lines to delete it</div><div>slframe = newsl; </div><div>Surface Loop(slframe)={ s1, fex1, fex2, fex3, fex4,</div><div>                        s1d, fet1, fet2, fet3, fet4};</div><div>frame = newv; Volume(frame) = {slframe};</div><div>Physical Volume("Frame") = {frame};</div><div><br></div><div>//now create the surface enclosing both mags and the frame</div><div>solids = newsl; </div><div>Surface Loop(solids) = { s1d, fet1, fet2, fet3, fet4, s1,</div><div>    hof1, mag1[2], mag1[3], mag1[4], mag1[5], mag1[0], fex2,</div><div>    hof2, mag2[2], mag2[3], mag2[4], mag2[5], mag2[0], fex4};</div><div><br></div><div>//parameters used by this macro</div><div>// cx, cy, cz: center of sphere</div><div>// rad: radius of sphere</div><div>// slc: mesh size of sphere</div><div><br></div><div>Macro MakeSphereSurface</div><div>  cp = newp; Point(cp) = {cx, cy, cz, slc}; //center</div><div>  np = newp; Point(np) = {cx, cy, cz+rad, slc}; //north</div><div>  sp = newp; Point(sp) = {cx, cy, cz-rad, slc}; //south</div><div>  dd[] = {rad, 0, -rad, 0};</div><div>  For i In {0:3}</div><div>    ap = newp;  pe[i] = ap; //points on the equator</div><div><span class="Apple-tab-span" style="white-space:pre">       </span>Point(ap) = {cx+dd[i], cy+dd[3-i], cz, slc};</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>na[i] = newreg; Circle(na[i]) = {ap, cp, np}; //north arc</div><div><span class="Apple-tab-span" style="white-space:pre">    </span>sa[i] = newreg; Circle(sa[i]) = {ap, cp, sp}; //south arc</div><div><span class="Apple-tab-span" style="white-space:pre">    </span>If (i>0)</div><div><span class="Apple-tab-span" style="white-space:pre">  </span>   eq[i-1] = newreg; //equator</div><div><span class="Apple-tab-span" style="white-space:pre">  </span>   Circle(eq[i-1]) = {pe[i-1], cp, ap};</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>EndIf</div><div>  EndFor</div><div>  //closing arc of the equator</div><div>  eq[3] = newreg; Circle(eq[3]) = {ap, cp, pe[0]};</div><div>  //now the non-plane surfaces: sloop[]</div><div>  For i In {0:3}</div><div>    lloop = newreg; Line Loop(lloop)={eq[i], -na[i], na[(i+1)%4]};</div><div><span class="Apple-tab-span" style="white-space:pre">  </span>sloop[2*i] = newreg; Ruled Surface(sloop[2*i]) = {lloop};</div><div>    lloop = newreg; Line Loop(lloop)={eq[i], -sa[i], sa[(i+1)%4]};</div><div>    sloop[2*i+1] = newreg; Ruled Surface(sloop[2*i+1]) = {lloop};</div><div>  EndFor</div><div>  thesurf = newreg;</div><div>  Surface Loop(thesurf) = {sloop[]};</div><div>Return</div><div><br></div><div>// create air box around magnets</div><div>BoundingBox; // recompute model bounding box</div><div>cx = (General.MinX + General.MaxX) / 2;</div><div>cy = (General.MinY + General.MaxY) / 2;</div><div>cz = (General.MinZ + General.MaxZ) / 2;</div><div>lx = 2*inf + General.MaxX - General.MinX;</div><div>ly = 2*inf + General.MaxY - General.MinZ;</div><div>lz = 2*inf + General.MaxZ - General.MinZ;</div><div>rad = (lx + ly + lz)/3;</div><div>slc = lc2;</div><div><br></div><div>Call MakeSphereSurface ;</div><div><br></div><div>vair = newv; Volume(vair) = {thesurf, solids};</div><div>Physical Volume("Air") = {vair}; // air</div><div>Physical Surface("Infty") = sloop[];</div><div><br></div><div>// adjust value here for correct merge result</div><div>//Geometry.Tolerance = 1e-6; </div><div>//Coherence Mesh;</div></div>