Revision 1384 - (show annotations)
Fri Jan 11 02:29:38 2008 UTC (13 years ago) by phornby
File MIME type: application/x-tex
File size: 14007 byte(s)
Make a temp copy of the trunk before checking in the windows changes


 1 \chapter{The Module \pycad} \label{PYCAD CHAP} 2 3 4 \section{Introduction} 5 6 \pycad provides a simple way to build a mesh for your finite element 7 simulation. You begin by building what we call a {\it Design} using 8 primitive geometric objects, and then to go on to build a mesh from 9 the {\it Design}. The final step of generating the mesh from a {\it 10 Design} uses freely available mesh generation software, such as \gmshextern. 11 12 A {\it Design} is built by defining points, which are used to specify 13 the corners of geometric objects and the vertices of curves. Using 14 points you construct more interesting objects such as lines, 15 rectangles, and arcs. By adding many of these objects into what we 16 call a {\it Design}, you can build meshes for arbitrarily complex 2-D 17 and 3-D structures. 18 19 The example included below shows how to use {\it pycad} to create a 2-D mesh 20 in the shape of a trapezoid with a cutout area. 21 22 \begin{python} 23 from esys.pycad import * 24 from esys.pycad.gmsh import Design 25 from esys.finley import MakeDomain 26 27 # A trapezoid 28 p0=Point(0.0, 0.0, 0.0) 29 p1=Point(1.0, 0.0, 0.0) 30 p2=Point(1.0, 0.5, 0.0) 31 p3=Point(0.0, 1.0, 0.0) 32 l01=Line(p0, p1) 33 l12=Line(p1, p2) 34 l23=Line(p2, p3) 35 l30=Line(p3, p0) 36 c=CurveLoop(l01, l12, l23, l30) 37 38 # A small triangular cutout 39 x0=Point(0.1, 0.1, 0.0) 40 x1=Point(0.5, 0.1, 0.0) 41 x2=Point(0.5, 0.2, 0.0) 42 x01=Line(x0, x1) 43 x12=Line(x1, x2) 44 x20=Line(x2, x0) 45 cutout=CurveLoop(x01, x12, x20) 46 47 # Create the surface with cutout 48 s=PlaneSurface(c, holes=[cutout]) 49 50 # Create a Design which can make the mesh 51 d=Design(dim=2, element_size=0.05) 52 53 # Add the trapezoid with cutout 54 d.addItems(s) 55 56 # Create the geometry, mesh and Escript domain 57 d.setScriptFileName("trapezoid.geo") 58 d.setMeshFileName("trapezoid.msh") 59 domain=MakeDomain(d, integrationOrder=-1, reducedIntegrationOrder=-1, optimizeLabeling=True) 60 61 # Create a file that can be read back in to python with mesh=ReadMesh(fileName) 62 domain.write("trapezoid.fly") 63 \end{python} 64 65 This example is included with the software in 66 \code{pycad/examples/trapezoid.py}. If you have gmsh installed you can 67 run the example and view the geometry and mesh with: 68 69 \begin{python} 70 python trapezoid.py 71 gmsh trapezoid.geo 72 gmsh trapezoid.msh 73 \end{python} 74 75 A \code{CurveLoop} is used to connect several lines into a single curve. 76 It is used in the example above to create the trapezoidal outline for the grid 77 and also for the triangular cutout area. 78 You can use any number of lines when creating a \code{CurveLoop}, but 79 the end of one line must be identical to the start of the next. 80 81 Sometimes you might see us write \code{-c} where \code{c} is a 82 \code{CurveLoop}. This is the reverse curve of the curve \code{c}. 83 It is identical to the original except that its points are traversed 84 in the opposite order. This may make it easier to connect two curves 85 in a \code{CurveLoop}. 86 87 The example python script above calls both 88 \code{d.setScriptFileName()} and \code{d.setMeshFileName()}. You need 89 only call these if you wish to save the gmsh geometry and mesh files. 90 91 Note that the underlying mesh generation software will not accept all 92 the geometries you can create with {\it pycad}. For example, {\it 93 pycad} will happily allow you to create a 2-D {\it Design} that is a 94 closed loop with some additional points or lines lying outside of the 95 enclosed area, but gmsh will fail to create a mesh for it. 96 97 98 99 100 101 102 \section{\pycad Classes} 103 \declaremodule{extension}{esys.pycad} 104 \modulesynopsis{Python geometry description and meshing interface} 105 106 \subsection{Primitives} 107 108 Some of the most commonly-used objects in {\it pycad} are listed here. For a more complete 109 list see the full API documentation. 110 111 \begin{classdesc}{Point}{x1, x2, x3} 112 Create a point with from coordinates. 113 \end{classdesc} 114 115 \begin{classdesc}{Line}{point1, point2} 116 Create a line with between starting and ending points. 117 \end{classdesc} 118 119 \begin{classdesc}{Curve}{point1, point2, ...} 120 Create a \code{Curve}, which is simply a list of points. 121 \end{classdesc} 122 123 \begin{classdesc}{Spline}{curve} 124 Interpret a \code{Curve} using a spline. 125 \end{classdesc} 126 127 \begin{classdesc}{BSpline}{curve} 128 Interpret a \code{Curve} using a b-spline. 129 \end{classdesc} 130 131 \begin{classdesc}{BezierCurve}{curve} 132 Interpret a \code{Curve} using a Bezier curve. 133 \end{classdesc} 134 135 \begin{classdesc}{CurveLoop}{list} 136 Create a closed \code{Curve} connecting the lines and/or points given in the \code{list}. 137 \end{classdesc} 138 139 \begin{classdesc}{Arc}{center_point, start_point, end_point} 140 Create an arc by specifying a center for a circle and start and end points. An arc may subtend an angle of at most $\pi$ radians. 141 \end{classdesc} 142 143 \begin{classdesc}{PlaneSurface}{loop, \optional{holes=[list]}} 144 Create a surface for a 2-D mesh, which may have one or more holes. 145 \end{classdesc} 146 147 \begin{classdesc}{RuledSurface}{list} 148 Create a surface that can be interpolated using transfinite interpolation. 149 \end{classdesc} 150 151 \begin{classdesc}{SurfaceLoop}{list} 152 Create a loop of 2D primitives, which defines the shell of a volume. 153 \end{classdesc} 154 155 \begin{classdesc}{Volume}{loop, \optional{holes=[list]}} 156 Create a volume for a 3-D mesh given a SurfaceLoop, which may have one or more holes. 157 \end{classdesc} 158 159 \begin{classdesc}{PropertySet}{list} 160 Create a PropertySet given a list of 1-D, 2-D or 3-D items. See the section on Properties below for more information. 161 \end{classdesc} 162 163 %============================================================================================================ 164 \subsection{Transformations} 165 166 Sometimes it's convenient to create an object and then make copies at 167 different orientations and in different sizes. Transformations are 168 used to move geometrical objects in the 3-dimensional space and to 169 resize them. 170 171 \begin{classdesc}{Translation}{\optional{b=[0,0,0]}} 172 defines a translation $x \to x+b$. \var{b} can be any object that can be converted 173 into a \numarray object of shape $(3,)$. 174 \end{classdesc} 175 176 \begin{classdesc}{Rotatation}{\optional{axis=[1,1,1], \optional{ point = [0,0,0], \optional{angle=0*RAD} } } } 177 defines a rotation by \var{angle} around axis through point \var{point} and direction \var{axis}. 178 \var{axis} and \var{point} can be any object that can be converted 179 into a \numarray object of shape $(3,)$. 180 \var{axis} does not have to be normalized but must have positive length. The right hand rule~\cite{RIGHTHANDRULE} 181 applies. 182 \end{classdesc} 183 184 185 \begin{classdesc}{Dilation}{\optional{factor=1., \optional{center=[0,0,0]}}} 186 defines a dilation by the expansion/contraction \var{factor} with 187 \var{center} as the dilation center. 188 \var{center} can be any object that can be converted 189 into a \numarray object of shape $(3,)$. 190 \end{classdesc} 191 192 \begin{classdesc}{Reflection}{\optional{normal=[1,1,1], \optional{offset=0}}} 193 defines a reflection on a plane defined in normal form $n^t x = d$ 194 where $n$ is the surface normal \var{normal} and $d$ is the plane \var{offset}. 195 \var{normal} can be any object that can be converted 196 into a \numarray object of shape $(3,)$. 197 \var{normal} does not have to be normalized but must have positive length. 198 \end{classdesc} 199 200 \begin{datadesc}{DEG} 201 A constant to convert from degrees to an internal angle representation in radians. For instance use \code{90*DEG} for $90$ degrees. 202 \end{datadesc} 203 204 \subsection{Properties} 205 206 If you are building a larger geometry you may find it convenient to 207 create it in smaller pieces and then assemble them into the whole. 208 Property sets make this easy, and they allow you to name the smaller 209 pieces for convenience. 210 211 Property sets are used to bundle a set of geometrical objects in a 212 group. The group is identified by a name. Typically a property set 213 is used to mark subregions with share the same material properties or 214 to mark portions of the boundary. For efficiency, the \Design class 215 object assigns a integer to each of its property sets, a so-called tag 216 \index{tag}. The appropriate tag is attached to the elements at 217 generation time. 218 219 See the file \code{pycad/examples/quad.py} for an example using a {\it PropertySet}. 220 221 222 \begin{classdesc}{PropertySet}{name,*items} 223 defines a group geometrical objects which can be accessed through a \var{name} 224 The objects in the tuple \var{items} mast all be \ManifoldOneD, \ManifoldTwoD or \ManifoldThreeD objects. 225 \end{classdesc} 226 227 228 \begin{methoddesc}[PropertySet]{getManifoldClass}{} 229 returns the manifold class \ManifoldOneD, \ManifoldTwoD or \ManifoldThreeD expected from the items 230 in the property set. 231 \end{methoddesc} 232 233 \begin{methoddesc}[PropertySet]{getDim}{} 234 returns the spatial dimension of the items 235 in the property set. 236 \end{methoddesc} 237 238 \begin{methoddesc}[PropertySet]{getName}{} 239 returns the name of the set 240 \end{methoddesc} 241 242 \begin{methoddesc}[PropertySet]{setName}{name} 243 sets the name. This name should be unique within a \Design. 244 \end{methoddesc} 245 246 \begin{methoddesc}[PropertySet]{addItem}{*items} 247 adds a tuple of items. They need to be objects of class \ManifoldOneD, \ManifoldTwoD or \ManifoldThreeD. 248 \end{methoddesc} 249 250 \begin{methoddesc}[PropertySet]{getItems}{} 251 returns the list of items 252 \end{methoddesc} 253 254 \begin{methoddesc}[PropertySet]{clearItems}{} 255 clears the list of items 256 \end{methoddesc} 257 258 \begin{methoddesc}[PropertySet]{getTag}{} 259 returns the tag used for this property set 260 \end{methoddesc} 261 262 \section{Interface to the mesh generation software} 263 \declaremodule{extension}{esys.pycad.gmsh} 264 \modulesynopsis{Python geometry description and meshing interface} 265 266 The class and methods described here provide an interface to the mesh 267 generation software, which is currently gmsh. This interface could be 268 adopted to triangle or another mesh generation package if this is 269 deemed to be desireable in the future. 270 271 \begin{classdesc}{Design}{ 272 \optional{dim=3, \optional{element_size=1., \optional{order=1, \optional{keep_files=False}}}}} 273 The \class{Design} describes the geometry defined by primitives to be meshed. 274 The \var{dim} specifies the spatial dimension. The argument \var{element_size} defines the global 275 element size which is multiplied by the local scale to set the element size at each \Point. 276 The argument \var{order} defines the element order to be used. If \var{keep_files} is set to 277 \True temporary files a kept otherwise they are removed when the instance of the class is deleted. 278 \end{classdesc} 279 280 281 \begin{methoddesc}[Design]{setDim}{\optional{dim=3}} 282 sets the spatial dimension which needs to be $1$, $2$ or $3$. 283 \end{methoddesc} 284 285 \begin{methoddesc}[Design]{getDim}{} 286 returns the spatial dimension. 287 \end{methoddesc} 288 289 \begin{methoddesc}[Design]{setElementOrder}{\optional{order=1}} 290 sets the element order which needs to be $1$ or $2$. 291 \end{methoddesc} 292 293 \begin{methoddesc}[Design]{getElementOrder}{} 294 returns the element order. 295 \end{methoddesc} 296 297 298 \begin{methoddesc}[Design]{setElementSize}{\optional{element_size=1}} 299 set the global element size. The local element size at a point is defined as 300 the global element size multipied by the local scale. The element size must be positive. 301 \end{methoddesc} 302 303 304 \begin{methoddesc}[Design]{getElementSize}{} 305 returns the global element size. 306 \end{methoddesc} 307 308 \begin{memberdesc}[Design]{DELAUNAY} 309 the gmsh Delauny triangulator. 310 \end{memberdesc} 311 312 \begin{memberdesc}[Design]{TETGEN} 313 the TetGen~\cite{TETGEN} triangulator. 314 \end{memberdesc} 315 316 \begin{memberdesc}[Design]{TETGEN} 317 the NETGEN~\cite{NETGEN} triangulator. 318 \end{memberdesc} 319 320 \begin{methoddesc}[Design]{setKeepFilesOn}{} 321 work files are kept at the end of the generation. 322 \end{methoddesc} 323 324 \begin{methoddesc}[Design]{setKeepFilesOff}{} 325 work files are deleted at the end of the generation. 326 \end{methoddesc} 327 328 \begin{methoddesc}[Design]{keepFiles}{} 329 returns \True if work files are kept. Otherwise \False is returned. 330 \end{methoddesc} 331 332 \begin{methoddesc}[Design]{setScriptFileName}{\optional{name=None}} 333 set the filename for the gmsh input script. if no name is given a name with extension "geo" is generated. 334 \end{methoddesc} 335 336 \begin{methoddesc}[Design]{getScriptFileName}{} 337 returns the name of the file for the gmsh script. 338 \end{methoddesc} 339 340 341 \begin{methoddesc}[Design]{setMeshFileName}{\optional{name=None}} 342 sets the name for the gmsh mesh file. if no name is given a name with extension "msh" is generated. 343 \end{methoddesc} 344 345 \begin{methoddesc}[Design]{getMeshFileName}{} 346 returns the name of the file for the gmsh msh 347 \end{methoddesc} 348 349 350 \begin{methoddesc}[Design]{addItems}{*items} 351 adds the tuple of var{items}. An item can be any primitive or a \class{PropertySet}. 352 \warning{If a \PropertySet is added as an item added object that are not 353 part of a \PropertySet are not considered in the messing. 354 } 355 356 \end{methoddesc} 357 358 \begin{methoddesc}[Design]{getItems}{} 359 returns a list of the items 360 \end{methoddesc} 361 362 \begin{methoddesc}[Design]{clearItems}{} 363 resets the items in design 364 \end{methoddesc} 365 366 \begin{methoddesc}[Design]{getMeshHandler}{} 367 returns a handle to the mesh. The call of this method generates the mesh from the geometry and 368 returns a mechnism to access the mesh data. In the current implementation this 369 method returns a file name for a gmsh file containing the mesh data. 370 \end{methoddesc} 371 372 \begin{methoddesc}[Design]{getScriptString}{} 373 returns the gmsh script to generate the mesh as a string. 374 \end{methoddesc} 375 376 \begin{methoddesc}[Design]{getCommandString}{} 377 returns the gmsh command used to generate the mesh as string. 378 \end{methoddesc} 379 380 \begin{methoddesc}[Design]{setOptions}{\optional{algorithm=None, \optional{ optimize_quality=True,\optional{ smoothing=1}}}} 381 sets options for the mesh generator. \var{algorithm} sets the algorithm to be used. 382 The algorithm needs to be \var{Design.DELAUNAY} 383 \var{Design.TETGEN} 384 or \var{Design.NETGEN}. By default \var{Design.DELAUNAY} is used. \var{optimize_quality}=\True invokes an optimization of the mesh quality. \var{smoothing} sets the number of smoothing steps to be applied to the mesh. 385 \end{methoddesc} 386 387 \begin{methoddesc}[Design]{getTagMap}{} 388 returns a \class{TagMap} to map the name \class{PropertySet} in the class to tag numbers generated by gmsh. 389 \end{methoddesc}