[Gmsh] Swiss Cheese in GMSH

Ruth V. Sabariego r.sabariego at ulg.ac.be
Wed Mar 30 17:26:25 CEST 2011



On 28/03/11 23:29, Lorenzo Isella wrote:
> Hello,
> And many thanks for your suggestions!
> I somehow managed to put together the script pasted below, but I have 
> a few questions (I must admit I myself do not understand 100% of what 
> I have written).
> (1) first of all: is it correct? I mean: am I generating a 
> computational domain given by a cube minus 5 spheres?
There are 6 volumes in your geometry: the five spheres and the cube 
minus the 5 spheres.
You can check by meshing and using the clipping tool (menu Tools)
> (2) if it is correct, how can I read the positions of the holes from a 
> text file?
Just add to your file:
Include "positions.txt" ;

Regards,
Ruth

> Many thanks
>
> Lorenzo
>
> Function CheeseHole
>
>   // In the following commands we use the reserved variable name
>   // `newp', which automatically selects a new point number. This
>   // number is chosen as the highest current point number, plus
>   // one. (Note that, analogously to `newp', the variables `newc',
>   // `news', `newv' and `newreg' select the highest number amongst
>   // currently defined curves, surfaces, volumes and `any entities
>   // other than points', respectively.)
>
>   p1 = newp; Point(p1) = {x,  y,  z,  lcar3} ;
>   p2 = newp; Point(p2) = {x+r,y,  z,  lcar3} ;
>   p3 = newp; Point(p3) = {x,  y+r,z,  lcar3} ;
>   p4 = newp; Point(p4) = {x,  y,  z+r,lcar3} ;
>   p5 = newp; Point(p5) = {x-r,y,  z,  lcar3} ;
>   p6 = newp; Point(p6) = {x,  y-r,z,  lcar3} ;
>   p7 = newp; Point(p7) = {x,  y,  z-r,lcar3} ;
>
>   c1 = newreg; Circle(c1) = {p2,p1,p7};
>   c2 = newreg; Circle(c2) = {p7,p1,p5};
>   c3 = newreg; Circle(c3) = {p5,p1,p4};
>   c4 = newreg; Circle(c4) = {p4,p1,p2};
>   c5 = newreg; Circle(c5) = {p2,p1,p3};
>   c6 = newreg; Circle(c6) = {p3,p1,p5};
>   c7 = newreg; Circle(c7) = {p5,p1,p6};
>   c8 = newreg; Circle(c8) = {p6,p1,p2};
>   c9 = newreg; Circle(c9) = {p7,p1,p3};
>   c10 = newreg; Circle(c10) = {p3,p1,p4};
>   c11 = newreg; Circle(c11) = {p4,p1,p6};
>   c12 = newreg; Circle(c12) = {p6,p1,p7};
>
>   // We need non-plane surfaces to define the spherical holes. Here we
>   // use ruled surfaces, which can have 3 or 4 sides:
>
>   l1 = newreg; Line Loop(l1) = {c5,c10,c4};   Ruled Surface(newreg) = 
> {l1};
>   l2 = newreg; Line Loop(l2) = {c9,-c5,c1};   Ruled Surface(newreg) = 
> {l2};
>   l3 = newreg; Line Loop(l3) = {c12,-c8,-c1}; Ruled Surface(newreg) = 
> {l3};
>   l4 = newreg; Line Loop(l4) = {c8,-c4,c11};  Ruled Surface(newreg) = 
> {l4};
>   l5 = newreg; Line Loop(l5) = {-c10,c6,c3};  Ruled Surface(newreg) = 
> {l5};
>   l6 = newreg; Line Loop(l6) = {-c11,-c3,c7}; Ruled Surface(newreg) = 
> {l6};
>   l7 = newreg; Line Loop(l7) = {-c2,-c7,-c12};Ruled Surface(newreg) = 
> {l7};
>   l8 = newreg; Line Loop(l8) = {-c6,-c9,c2};  Ruled Surface(newreg) = 
> {l8};
>
>   // We then store the surface loops identification numbers in a list
>   // for later reference (we will need these to define the final
>   // volume):
>
>   theloops[t] = newreg ;
>
>   Surface Loop(theloops[t]) = {l8+1,l5+1,l1+1,l2+1,l3+1,l7+1,l6+1,l4+1};
>
>   thehole = newreg ;
>   Volume(thehole) = theloops[t] ;
>
> Return
>
> lcar3 = .055;
>
> hloc = 0.1;
> Point(1) = {0, 0, 0, hloc};
> Point(2) = {1, 0, 0, hloc};
> Point(3) = {1, 1, 0, hloc};
> Point(4) = {0, 1, 0, hloc};
>
> Point(5) = {0, 0, 1, hloc};
> Point(6) = {1, 0, 1, hloc};
> Point(7) = {1, 1, 1, hloc};
> Point(8) = {0, 1, 1, hloc};
>
> Line(1) = {1,2};
> Line(2) = {2,3};
> Line(3) = {3,4};
> Line(4) = {4,1};
>
> Line(5) = {5,6};
> Line(6) = {6,7};
> Line(7) = {7,8};
> Line(8) = {8,5};
>
> Line(9)  = {1,5};
> Line(10) = {2,6};
> Line(11) = {3,7};
> Line(12) = {4,8};
>
> // bottom
> Line Loop(21) = {-1,-4,-3,-2};
> Plane Surface(31) = {21} ;
>
> // top
> Line Loop(22) = {5,6,7,8};
> Plane Surface(32) = {22} ;
>
> // left
> Line Loop(23) = {1,10,-5,-9};
> Plane Surface(33) = {23} ;
>
> // right
> Line Loop(24) = {12,-7,-11,3};
> Plane Surface(34) = {24} ;
>
> // front
> Line Loop(25) = {2,11,-6,-10};
> Plane Surface(35) = {25} ;
>
> // back
> Line Loop(26) = {9,-8,-12,4};
> Plane Surface(36) = {26} ;
>
> // We can use a `For' loop to generate five holes in the cube:
>
> x = 0 ; y = 0.75 ; z = 0 ; r = 0.09 ;
>
> For t In {1:5}
>
>   x += 0.166 ;
>   z += 0.166 ;
>
>   // We call the `CheeseHole' function:
>
>   Call CheeseHole ;
>
>   // We define a physical volume for each hole:
>
>   Physical Volume (t) = thehole ;
>
>   // We also print some variables on the terminal (note that, since
>   // all variables are treated internally as floating point numbers,
>   // the format string should only contain valid floating point format
>   // specifiers like `%g', `%f', '%e', etc.):
>
>   Printf("Hole %g (center = {%g,%g,%g}, radius = %g) has number %g!",
>      t, x, y, z, r, thehole) ;
>
> EndFor
>
> theloops[0] = newreg ;
>
>
> Surface Loop(theloops[0]) = {31,32,33,34,35,36};
>
> Volume(100) = {theloops[]};
>
> Physical Volume (10) = 100 ;
>
>
>
> On 28/03/11 21:00, David Colignon wrote:
>> Hi Lorenzo,
>>
>> you can start with:
>>
>> http://geuz.org/gmsh/doc/texinfo/gmsh.html#t5_002egeo
>>
>> and the "CheesHole" function within :-)
>>
>> Regards,
>>
>> Dave
>>
>> -- 
>> David Colignon, Ph.D.
>> Collaborateur Logistique du F.R.S.-FNRS
>> CÉCI - Consortium des Équipements de Calcul Intensif
>> ACE - Applied & Computational Electromagnetics
>> Institut Montefiore B28
>> Université de Liège
>> 4000 Liège - BELGIQUE
>> Tél: +32 (0)4 366 37 32
>> Fax: +32 (0)4 366 29 10
>> WWW: http://hpc.montefiore.ulg.ac.be/
>> Agenda: 
>> http://www.google.com/calendar/embed?src=david.colignon%40gmail.com
>>
>> On 28/03/11 20:46, Lorenzo Isella wrote:
>>> Dear All,
>>> I hope that the title of the email will not be misleading.
>>> I have just started to learn how to use fipy, which in turns relies on
>>> Gmsh for non-trivial meshes.
>>> My problem is the following: I am given a computational domain which
>>> looks like a cube having some spherical cavities in
>>> its interior.
>>> That is why I was referring to Swiss cheese in the email title.
>>> To be more precise the situation is the following: all the spherical
>>> cavities have the same radius R and I am given a
>>> simple text file like
>>>
>>> x_1 y_1 z_1
>>> x_2 y_2 z_2
>>> ..............
>>> x_N y_N z_N
>>>
>>> where
>>> x_i y_i z_i
>>> are the 3D coordinates of the i-th cavity.
>>> Now, what I would like to do is to script Gmsh in such a way to
>>> (1) generate a cube having side [-L/2, L/2] along x, y and z
>>> (2) introduce the N cavities at the right positions
>>> (3) mesh the resulting object being able in particular to refine the
>>> mesh around the cavities (which is where I need
>>> more accuracy).
>>>
>>> Is anything like that doable with GMSH? According to what I read on
>>> the fipy mailing list and saw on the gmsh tutorials,
>>> the answer is a definite yes, but some help is needed.
>>> Even a simple script which just a couple of cavities would help me a 
>>> lot.
>>> Any suggestion is appreciated
>>>
>>> Lorenzo
>>>
>>> PS: the overall shape of the computational domain is not essential for
>>> me, provided it is large enough. In other words:
>>> the cube can be replaced with a large sphere or a cylinder if it eases
>>> the numerical work.
>>>
>>> _______________________________________________
>>> gmsh mailing list
>>> gmsh at geuz.org
>>> http://www.geuz.org/mailman/listinfo/gmsh
>>
>
> _______________________________________________
> gmsh mailing list
> gmsh at geuz.org
> http://www.geuz.org/mailman/listinfo/gmsh
>

-- 
Dr. Ir. Ruth V. Sabariego
University of Liege, Dept. of Electrical Engineering&  Computer Science,
Applied&  Computational Electromagnetics (ACE),
phone: +32-4-3663737 - fax: +32-4-3662910 - http://ace.montefiore.ulg.ac.be/