/[escript]/trunk/pycad/py_src/gmsh.py
ViewVC logotype

Diff of /trunk/pycad/py_src/gmsh.py

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1005 by ksteube, Fri Mar 2 06:50:48 2007 UTC revision 2319 by gross, Thu Mar 19 03:19:41 2009 UTC
# Line 1  Line 1 
1  # $Id:$  
2    ########################################################
3    #
4    # Copyright (c) 2003-2008 by University of Queensland
5    # Earth Systems Science Computational Center (ESSCC)
6    # http://www.uq.edu.au/esscc
7    #
8    # Primary Business: Queensland, Australia
9    # Licensed under the Open Software License version 3.0
10    # http://www.opensource.org/licenses/osl-3.0.php
11    #
12    ########################################################
13    
14    __copyright__="""Copyright (c) 2003-2008 by University of Queensland
15    Earth Systems Science Computational Center (ESSCC)
16    http://www.uq.edu.au/esscc
17    Primary Business: Queensland, Australia"""
18    __license__="""Licensed under the Open Software License version 3.0
19    http://www.opensource.org/licenses/osl-3.0.php"""
20    __url__="http://www.uq.edu.au/esscc/escript-finley"
21    
22  """  """
23  mesh generation using gmsh  mesh generation using gmsh
# Line 11  mesh generation using gmsh Line 30  mesh generation using gmsh
30  @var __date__: date of the version  @var __date__: date of the version
31  """  """
32    
   
33  __author__="Lutz Gross, l.gross@uq.edu.au"  __author__="Lutz Gross, l.gross@uq.edu.au"
 __copyright__="""  Copyright (c) 2007 by ACcESS MNRF  
                     http://www.access.edu.au  
                 Primary Business: Queensland, Australia"""  
 __license__="""Licensed under the Open Software License version 3.0  
              http://www.opensource.org/licenses/osl-3.0.php"""  
 __url__="http://www.iservo.edu.au/esys/escript"  
 __version__="$Revision:$"  
 __date__="$Date:$"  
34    
35  import design  import design
36  import tempfile  import tempfile
37  import os  import os
38    from primitives import Point, Spline, BezierCurve, BSpline, Line, Arc, CurveLoop, RuledSurface, PlaneSurface, SurfaceLoop, Volume, PropertySet, Ellipse
39    from esys.escript import getMPIWorldMax, getMPIRankWorld
40    
41  class Design(design.Design):  class Design(design.Design):
42      """      """
43      design fo gmsh      Design for gmsh.
44      """      """
45      DELAUNAY="iso"      DELAUNAY="iso"
46      NETGEN="netgen"      NETGEN="netgen"
47      TETGEN="tetgen"      TETGEN="tetgen"
48    
49      def __init__(self,dim=3,element_size=1.,order=1,keep_files=False):      def __init__(self,dim=3,element_size=1.,order=1,keep_files=False):
50         """         """
51         initializes the gmsh design         Initializes the gmsh design.
52    
53         @param dim: patial dimension         @param dim: spatial dimension
54         @param element_size: global element size         @param element_size: global element size
55         @param order: element order         @param order: element order
56         @param keep_files: flag to keep work files.         @param keep_files: flag to keep work files
57         """         """
58         design.Design.__init__(self,dim=dim,element_size=element_size,order=order,keep_files=keep_files)         design.Design.__init__(self,dim=dim,element_size=element_size,order=order,keep_files=keep_files)
59         self.setScriptFileName()         self.setScriptFileName()
60         self.setMeshFileName()         self.setMeshFileName()
61         self.setOptions()         self.setOptions()
62    
63      def setScriptFileName(self,name=None):      def setScriptFileName(self,name=None):
64         """         """
65         set the filename for the gmsh input script. if no name is given a name with extension geo is generated.         Sets the filename for the gmsh input script. If no name is given a name
66           with extension I{geo} is generated.
67         """         """
68         if name == None:         if name == None:
69             self.__scriptname=tempfile.mkstemp(suffix=".geo")[1]             tmp_f_id=tempfile.mkstemp(suffix=".geo")
70               self.__scriptname=tmp_f_id[1]
71               os.close(tmp_f_id[0])
72         else:         else:
73             self.__scriptname=name             self.__scriptname=name
74             self.setKeepFilesOn()             self.setKeepFilesOn()
75    
76      def getScriptFileName(self):      def getScriptFileName(self):
77         """         """
78         returns the name of the file for the gmsh script         Returns the name of the gmsh script file.
79         """         """
80         return self.__scriptname         return self.__scriptname
81    
82      def setMeshFileName(self, name=None):      def setMeshFileName(self, name=None):
83         """         """
84         sets the name for the gmsh mesh file. if no name is given a name with extension msh is generated.         Sets the name for the gmsh mesh file. If no name is given a name with
85           extension I{msh} is generated.
86         """         """
87         if name == None:         if name == None:
88             self.__mshname=tempfile.mkstemp(suffix=".msh")[1]             tmp_f_id=tempfile.mkstemp(suffix=".msh")
89               self.__mshname=tmp_f_id[1]
90               os.close(tmp_f_id[0])
91         else:         else:
92             self.__mshname=name             self.__mshname=name
93             self.setKeepFilesOn()             self.setKeepFilesOn()
94    
95      def getMeshFileName(self):      def getMeshFileName(self):
96         """         """
97         returns the name of the file for the gmsh msh         Returns the name of the gmsh mesh file.
98         """         """
99         return self.__mshname         return self.__mshname
100    
101      def setOptions(self,algorithm=None,optimize_quality=True,smoothing=1, curvature_based_element_size=False):      def setOptions(self,algorithm=None,optimize_quality=True,smoothing=1, curvature_based_element_size=False):
102          """          """
103          sets options for the mesh generator          Sets options for the mesh generator.
104          """          """
105            if curvature_based_element_size:
106                  print "information: gmsh does not support curvature based element size anymore. Option ignored."
107          if algorithm==None: algorithm=self.DELAUNAY          if algorithm==None: algorithm=self.DELAUNAY
         self.__curvature_based_element_size=curvature_based_element_size  
108          self.__algo=algorithm          self.__algo=algorithm
109          self.__optimize_quality=optimize_quality          self.__optimize_quality=optimize_quality
110          self.__smoothing=smoothing          self.__smoothing=smoothing
111    
112      def __del__(self):      def __del__(self):
113          """          """
114          clean up          Cleans up.
115          """          """
116          if not self.keepFiles():          if not self.keepFiles() :
117                 os.unlink(self.getScriptFileName())              os.unlink(self.getScriptFileName())
118                 os.unlink(self.getMeshFileName())              os.unlink(self.getMeshFileName())
119      def getScriptString(self):  
         """  
         returns the gmsh script to generate the mesh  
         """  
         prim=self.getAllPrimitives()  
         out="// generated by esys.pycad\n"  
         for p in prim:  
            out+=p.getGmshCommand(self.getElementSize())+"\n"  
         return out  
120      def getCommandString(self):      def getCommandString(self):
121          """          """
122          returns the gmsh comand          Returns the gmsh command line.
123          """          """
124          if self.__optimize_quality:          if self.__optimize_quality:
125                opt="-optimize "                opt="-optimize "
126          else:          else:
127                opt=""                opt=""
128          if self.__curvature_based_element_size:  
129                clcurv="-clcurv "          exe="gmsh -format msh -%s -algo %s -smooth %s %s-v 0 -order %s -o %s %s" % (
130          else:                  self.getDim(), self.__algo, self.__smoothing, opt,
131                clcurv=""                  self.getElementOrder(), self.getMeshFileName(),
132                    self.getScriptFileName())
         exe="gmsh -%s -algo %s %s-smooth %s %s-v 0 -order %s -o %s %s"%(self.getDim(),  
                                                                        self.__algo,  
                                                                        clcurv,  
                                                                        self.__smoothing,  
                                                                        opt,  
                                                                        self.getElementOrder(),  
                                                                        self.getMeshFileName(),  
                                                                        self.getScriptFileName())  
133          return exe          return exe
134    
135      def getMeshHandler(self):      def getMeshHandler(self):
136          """          """
137          returns a handle to a mesh meshing the design. In the current implementation          Returns a handle to a mesh meshing the design. In the current
138          a mesh file name in gmsh format is returned.          implementation a mesh file name in gmsh format is returned.
139          """          """
         open(self.getScriptFileName(),"w").write(self.getScriptString())  
140          cmd = self.getCommandString()          cmd = self.getCommandString()
141          ret = os.system(cmd) / 256          if getMPIRankWorld() == 0:
142      if ret > 0:              open(self.getScriptFileName(),"w").write(self.getScriptString())
143        raise RuntimeError, "Could not build mesh: %s"%cmd              ret = os.system(cmd) / 256
144      else:          else:
145                ret=0
146            ret=getMPIWorldMax(ret)
147            if ret > 0:
148              raise RuntimeError, "Could not build mesh: %s"%cmd
149            else:
150            return self.getMeshFileName()            return self.getMeshFileName()
151    
152        def getScriptString(self):
153            """
154            Returns the gmsh script to generate the mesh.
155            """
156            h=self.getElementSize()
157            out="// generated by esys.pycad\n"
158            for prim in self.getAllPrimitives():
159               p=prim.getUnderlyingPrimitive()
160               if isinstance(p, Point):
161                   c=p.getCoordinates()
162                   out+="Point(%s) = {%s , %s, %s , %s };\n"%(p.getID(),c[0],c[1],c[2], p.getLocalScale()*h)
163    
164               elif isinstance(p, Spline):
165                   out+="Spline(%s) = {%s};\n"%(p.getID(),self.__mkArgs(p.getControlPoints()))
166    
167               elif isinstance(p, BezierCurve):
168                   out+="Bezier(%s) = {%s};\n"%(p.getID(),self.__mkArgs(p.getControlPoints()))
169    
170               elif isinstance(p, BSpline):
171                   out+="BSpline(%s) = {%s};\n"%(p.getID(),self.__mkArgs(p.getControlPoints()))
172    
173               elif isinstance(p, Line):
174                   out+="Line(%s) = {%s, %s};\n"%(p.getID(),p.getStartPoint().getDirectedID(),p.getEndPoint().getDirectedID())
175    
176               elif isinstance(p, Arc):
177                  out+="Circle(%s) = {%s, %s, %s};\n"%(p.getID(),p.getStartPoint().getDirectedID(),p.getCenterPoint().getDirectedID(),p.getEndPoint().getDirectedID())
178    
179               elif isinstance(p, Ellipse):
180                  out+="Ellipse(%s) = {%s, %s, %s, %s};\n"%(p.getID(),p.getStartPoint().getDirectedID(),p.getCenterPoint().getDirectedID(),p.getPointOnMainAxis().getDirectedID(), p.getEndPoint().getDirectedID())
181    
182               elif isinstance(p, CurveLoop):
183                   out+="Line Loop(%s) = {%s};\n"%(p.getID(),self.__mkArgs(p.getCurves()))
184    
185               elif isinstance(p, RuledSurface):
186                   out+="Ruled Surface(%s) = {%s};\n"%(p.getID(),p.getBoundaryLoop().getDirectedID())
187    
188               elif isinstance(p, PlaneSurface):
189                   line=self.__mkArgs(p.getHoles())
190                   if len(line)>0:
191                     out+="Plane Surface(%s) = {%s, %s};\n"%(p.getID(),p.getBoundaryLoop().getDirectedID(), line)
192                   else:
193                     out+="Plane Surface(%s) = {%s};\n"%(p.getID(),p.getBoundaryLoop().getDirectedID())
194    
195               elif isinstance(p, SurfaceLoop):
196                   out+="Surface Loop(%s) = {%s};\n"%(p.getID(),self.__mkArgs(p.getSurfaces()))
197    
198               elif isinstance(p, Volume):
199                   line=self.__mkArgs(p.getHoles())
200                   if len(line)>0:
201                     out+="Volume(%s) = {%s, %s};\n"%(p.getID(),p.getSurfaceLoop().getDirectedID(), line)
202                   else:
203                     out+="Volume(%s) = {%s};\n"%(p.getID(),p.getSurfaceLoop().getDirectedID())
204    
205               elif isinstance(p, PropertySet):
206                   if p.getNumItems()>0:
207                      dim=p.getDim()
208                      line="Physical "
209                      if dim==0:
210                          line+="Point"
211                      elif dim==1:
212                          line+="Line"
213                      elif dim==2:
214                          line+="Surface"
215                      else:
216                          line+="Volume"
217                      out+=line+"(" + str(p.getID()) + ") = {"+self.__mkArgs(p.getItems())+"};\n"
218    
219               else:
220                   raise TypeError("unable to pass %s object to gmsh."%str(type(p)))
221            return out
222    
223    
224        def __mkArgs(self,args):
225            line=""
226            for i in args:
227                if len(line)>0:
228                    line+=", %s"%i.getDirectedID()
229                else:
230                    line="%s"%i.getDirectedID()
231            return line
232    

Legend:
Removed from v.1005  
changed lines
  Added in v.2319

  ViewVC Help
Powered by ViewVC 1.1.26