/[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 932 by gross, Fri Jan 19 09:27:15 2007 UTC revision 2248 by gross, Thu Feb 5 23:24:04 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    
40  class Design(design.Design):  class Design(design.Design):
41      """      """
42      design fo gmsh      Design for gmsh.
43      """      """
44      ISO="iso"      DELAUNAY="iso"
     TRI="tri"  
     ANISO="aniso"  
45      NETGEN="netgen"      NETGEN="netgen"
46      TETGEN="tetgen"      TETGEN="tetgen"
47    
48      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):
49         """         """
50         initializes the gmsh design         Initializes the gmsh design.
51    
52         @param dim: patial dimension         @param dim: spatial dimension
53         @element_size: global element size         @param element_size: global element size
54         @order: element order         @param order: element order
55         @keep_files: flag to keep work files.         @param keep_files: flag to keep work files
56         """         """
57         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)
58         self.setScriptFileName()         self.setScriptFileName()
59         self.setMeshFileName()         self.setMeshFileName()
60         self.setOptions()         self.setOptions()
61    
62      def setScriptFileName(self,name=None):      def setScriptFileName(self,name=None):
63         """         """
64         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
65           with extension I{geo} is generated.
66         """         """
67         if name == None:         if name == None:
68             self.__scriptname=tempfile.mkstemp(suffix=".geo")[1]             tmp_f_id=tempfile.mkstemp(suffix=".geo")
69               self.__scriptname=tmp_f_id[1]
70               os.close(tmp_f_id[0])
71         else:         else:
72             self.__scriptname=name             self.__scriptname=name
73             self.setKeepFilesOn()             self.setKeepFilesOn()
74    
75      def getScriptFileName(self):      def getScriptFileName(self):
76         """         """
77         returns the name of the file for the gmsh script         Returns the name of the gmsh script file.
78         """         """
79         return self.__scriptname         return self.__scriptname
80    
81      def setMeshFileName(self, name=None):      def setMeshFileName(self, name=None):
82         """         """
83         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
84           extension I{msh} is generated.
85         """         """
86         if name == None:         if name == None:
87             self.__mshname=tempfile.mkstemp(suffix=".msh")[1]             tmp_f_id=tempfile.mkstemp(suffix=".msh")
88               self.__mshname=tmp_f_id[1]
89               os.close(tmp_f_id[0])
90         else:         else:
91             self.__mshname=name             self.__mshname=name
92             self.setKeepFilesOn()             self.setKeepFilesOn()
93    
94      def getMeshFileName(self):      def getMeshFileName(self):
95         """         """
96         returns the name of the file for the gmsh msh         Returns the name of the gmsh mesh file.
97         """         """
98         return self.__mshname         return self.__mshname
99      def setOptions(self,algorithm=None,optimize_quality=True,smoothing=3):  
100        def setOptions(self,algorithm=None,optimize_quality=True,smoothing=1, curvature_based_element_size=False):
101          """          """
102          sets options for the mesh generator          Sets options for the mesh generator.
103          """          """
104          if algorithm==None: algorithm=self.ISO          if curvature_based_element_size:
105                  print "information: gmsh does not support curvature based element size anymore. Option ignored."
106            if algorithm==None: algorithm=self.DELAUNAY
107          self.__algo=algorithm          self.__algo=algorithm
108          self.__optimize_quality=optimize_quality          self.__optimize_quality=optimize_quality
109          self.__smoothing=smoothing          self.__smoothing=smoothing
110    
111      def __del__(self):      def __del__(self):
112          """          """
113          clean up          Cleans up.
114          """          """
115          if not self.keepFiles():          if not self.keepFiles() :
116                 os.unlink(self.getScriptFileName())              os.unlink(self.getScriptFileName())
117                 os.unlink(self.getMeshFileName())              os.unlink(self.getMeshFileName())
118      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  
119      def getCommandString(self):      def getCommandString(self):
120          """          """
121          returns the gmsh comand          Returns the gmsh command line.
122          """          """
123          if self.__optimize_quality:          if self.__optimize_quality:
124                opt="-optimize "                opt="-optimize "
125          else:          else:
126                opt=""                opt=""
127          exe="gmsh -%s -algo %s -smooth %s %s -v 0 -order %s -o %s %s"%(self.getDim(),  
128                                                                         self.__algo,          exe="gmsh -format msh -%s -algo %s -smooth %s %s-v 0 -order %s -o %s %s" % (
129                                                                         self.__smoothing,                  self.getDim(), self.__algo, self.__smoothing, opt,
130                                                                         opt,                  self.getElementOrder(), self.getMeshFileName(),
131                                                                         self.getElementOrder(),                  self.getScriptFileName())
                                                                        self.getMeshFileName(),  
                                                                        self.getScriptFileName())  
132          return exe          return exe
133    
134      def getMeshHandler(self):      def getMeshHandler(self):
135          """          """
136          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
137          a mesh file name in gmsh format is returned.          implementation a mesh file name in gmsh format is returned.
138          """          """
139          open(self.getScriptFileName(),"w").write(self.getScriptString())          f = open(self.getScriptFileName(),"w")
140          os.system(self.getCommandString())          f.write(self.getScriptString())
141          return self.getMeshFileName()          f.close()
142  # to finley          cmd = self.getCommandString()
143  def MakeDomain(design,integrationOrder=-1,optimizeLabeling=True, richFaceElements=False):          ret = os.system(cmd) / 256
144          mshname=design.getMeshHandler()          if ret > 0:
145          print "MakeDomain :",mshname            raise RuntimeError, "Could not build mesh: %s"%cmd
146          1/0          else:
147          convertGmshToFinley(open(mshname,"r"),open(filename,"w"),dim=self.getDim())            return self.getMeshFileName()
148          if not self.keepTmpFiles(): os.unlink(mshname)  
149          return exe      def getScriptString(self):
150          line=gmsh_file.readline().split()          """
151          while len(line)>0:          Returns the gmsh script to generate the mesh.
152             print line          """
153             line=gmsh_file.readline().split()          h=self.getElementSize()
154            out="// generated by esys.pycad\n"
155            for prim in self.getAllPrimitives():
156               p=prim.getUnderlyingPrimitive()
157               if isinstance(p, Point):
158                   c=p.getCoordinates()
159                   out+="Point(%s) = {%s , %s, %s , %s };\n"%(p.getID(),c[0],c[1],c[2], p.getLocalScale()*h)
160    
161               elif isinstance(p, Spline):
162                   out+="Spline(%s) = {%s};\n"%(p.getID(),self.__mkArgs(p.getControlPoints()))
163    
164               elif isinstance(p, BezierCurve):
165                   out+="Bezier(%s) = {%s};\n"%(p.getID(),self.__mkArgs(p.getControlPoints()))
166    
167               elif isinstance(p, BSpline):
168                   out+="BSpline(%s) = {%s};\n"%(p.getID(),self.__mkArgs(p.getControlPoints()))
169    
170               elif isinstance(p, Line):
171                   out+="Line(%s) = {%s, %s};\n"%(p.getID(),p.getStartPoint().getDirectedID(),p.getEndPoint().getDirectedID())
172    
173               elif isinstance(p, Arc):
174                  out+="Circle(%s) = {%s, %s, %s};\n"%(p.getID(),p.getStartPoint().getDirectedID(),p.getCenterPoint().getDirectedID(),p.getEndPoint().getDirectedID())
175    
176               elif isinstance(p, Ellipse):
177                  out+="Ellipse(%s) = {%s, %s, %s, %s};\n"%(p.getID(),p.getStartPoint().getDirectedID(),p.getCenterPoint().getDirectedID(),p.getPointOnMainAxis().getDirectedID(), p.getEndPoint().getDirectedID())
178    
179               elif isinstance(p, CurveLoop):
180                   out+="Line Loop(%s) = {%s};\n"%(p.getID(),self.__mkArgs(p.getCurves()))
181    
182               elif isinstance(p, RuledSurface):
183                   out+="Ruled Surface(%s) = {%s};\n"%(p.getID(),p.getBoundaryLoop().getDirectedID())
184    
185               elif isinstance(p, PlaneSurface):
186                   line=self.__mkArgs(p.getHoles())
187                   if len(line)>0:
188                     out+="Plane Surface(%s) = {%s, %s};\n"%(p.getID(),p.getBoundaryLoop().getDirectedID(), line)
189                   else:
190                     out+="Plane Surface(%s) = {%s};\n"%(p.getID(),p.getBoundaryLoop().getDirectedID())
191    
192               elif isinstance(p, SurfaceLoop):
193                   out+="Surface Loop(%s) = {%s};\n"%(p.getID(),self.__mkArgs(p.getSurfaces()))
194    
195               elif isinstance(p, Volume):
196                   line=self.__mkArgs(p.getHoles())
197                   if len(line)>0:
198                     out+="Volume(%s) = {%s, %s};\n"%(p.getID(),p.getSurfaceLoop().getDirectedID(), line)
199                   else:
200                     out+="Volume(%s) = {%s};\n"%(p.getID(),p.getSurfaceLoop().getDirectedID())
201    
202               elif isinstance(p, PropertySet):
203                   if p.getNumItems()>0:
204                      dim=p.getDim()
205                      line="Physical "
206                      if dim==0:
207                          line+="Point"
208                      elif dim==1:
209                          line+="Line"
210                      elif dim==2:
211                          line+="Surface"
212                      else:
213                          line+="Volume"
214                      out+=line+"(" + str(p.getID()) + ") = {"+self.__mkArgs(p.getItems())+"};\n"
215    
216               else:
217                   raise TypeError("unable to pass %s object to gmsh."%str(type(p)))
218            return out
219    
220    
221        def __mkArgs(self,args):
222            line=""
223            for i in args:
224                if len(line)>0:
225                    line+=", %s"%i.getDirectedID()
226                else:
227                    line="%s"%i.getDirectedID()
228            return line
229    

Legend:
Removed from v.932  
changed lines
  Added in v.2248

  ViewVC Help
Powered by ViewVC 1.1.26