/[escript]/trunk/modellib/py_src/geometry.py
ViewVC logotype

Contents of /trunk/modellib/py_src/geometry.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 953 - (show annotations)
Fri Feb 9 08:55:28 2007 UTC (16 years, 1 month ago) by gross
File MIME type: text/x-python
File size: 21956 byte(s)
TagMap is supported by modellib now
1 # $Id$
2
3 __copyright__=""" Copyright (c) 2006 by ACcESS MNRF
4 http://www.access.edu.au
5 Primary Business: Queensland, Australia"""
6 __license__="""Licensed under the Open Software License version 3.0
7 http://www.opensource.org/licenses/osl-3.0.php"""
8
9
10 from esys.escript import *
11 from esys.escript.modelframe import Model,ParameterSet
12 from esys.pycad import TagMap
13 from esys import finley
14
15 class FinleyReader(ParameterSet):
16 """
17 reads finley mesh file.
18
19 @ivar source: mesh file in finley or gmsh format
20 @type source: C{DataSource}
21 @ivar intergrationOrder: integration order, default -1 (in).
22 @type intergrationOrder: C{int}
23 @ivar reducedIntegrationOrder: reduced integration order, default -1 (in).
24 @type reducedIntegrationOrder: C{int}
25 @ivar optimizeLabeling: switches on optimization of the labeling of the nodes
26 @type optimizeLabeling: C{bool}
27 """
28 def __init__(self,**kwargs):
29 """
30 initializes the object
31 """
32 super(FinleyReader,self).__init__(**kwargs)
33 self.declareParameter(source="none",
34 region_tag_map_source=None,
35 surface_tag_map_source=None,
36 optimizeLabeling=True,
37 reducedIntegrationOrder=-1,
38 integrationOrder=-1)
39 self.__domain=None
40 self.__region_tag_map=None
41 self.__surface_tag_map=None
42
43
44 def domain(self):
45 """
46 returns the domain
47
48 @return: the domain
49 @rtype: L{Domain}
50 """
51 if self.__domain == None:
52 if self.source.fileformat == "fly":
53 self.__domain=finley.ReadMesh(self.source.getLocalFileName(),self.integrationOrder)
54 elif self.source.fileformat == "gmsh":
55 self.__domain=finley.ReadGmsh(self.source.getLocalFileName(),self.integrationOrder,self.reducedIntegrationOrder, self.optimizeLabeling)
56 else:
57 raise TypeError("unknown mesh file format %s."%self.source.fileformat)
58 self.trace("mesh read from %s in %s format."%(self.source.getLocalFileName(), self.source.fileformat))
59 return self.__domain
60
61 def region_tag_map(self):
62 """
63 returns the map from regional tag names to tag integers used in the mesh
64
65 @return: the tag map
66 @rtype: L{TagMap}
67 """
68 if self.__region_tag_map == None:
69 self.__region_tag_map = TagMap()
70 if self.region_tag_map_source != None:
71 self.__region_tag_map.fillFromXML(open(self.region_tag_map_source.getLocalFileName()))
72 self.trace("region tag map read from %s in %s format."%(self.region_tag_map_source.getLocalFileName(), self.region_tag_map_source.fileformat))
73 return self.__region_tag_map
74
75 def surface_tag_map(self):
76 """
77 returns the map from surface tag names to tag integers used in the mesh
78
79 @return: the tag map
80 @rtype: L{TagMap}
81 """
82 if self.__surface_tag_map == None:
83 self.__surface_tag_map = TagMap()
84 if self.surface_tag_map_source != None:
85 self.__surface_tag_map.fillFromXML(open(self.surface_tag_map_source.getLocalFileName()))
86 self.trace("surface tag map read from %s in %s format."%(self.surface_tag_map_source.getLocalFileName(), self.surface_tag_map_source.fileformat))
87 return self.__surface_tag_map
88
89
90 class RectangularDomain(ParameterSet):
91 """
92 Generates a mesh over a rectangular domain finley.
93
94 @ivar dim: spatial dimension, default =2 (in).
95 @type dim: spatial dimension
96 @ivar l: spatial lengths, default [1.,1.,1.] (in).
97 @type l: C{list} of C{floats}s
98 @ivar n: number of elements, default [10,10,10] (in).
99 @type n: C{list} of C{int}s
100 @ivar order: element order, default 1 (in).
101 @type order: C{int}
102 @ivar periodic: flags for periodicity, default [False,False,False] (in).
103 @type periodic: C{list} of C{bool}s
104 @ivar intergrationOrder: integration order, default -1 (in).
105 @type intergrationOrder: C{int}
106 """
107 def __init__(self,**kwargs):
108 """
109 initializes the object
110 """
111 super(RectangularDomain,self).__init__(**kwargs)
112 self.declareParameter(dim=2,\
113 l=[1.,1.,1.],\
114 n=[10,10,10], \
115 order=1,\
116 periodic=[False,False,False],
117 integrationOrder=-1)
118 self.__domain=None
119
120 def domain(self):
121 """
122 returns the domain
123
124 @return: the domain
125 @rtype: L{Domain}
126 """
127 if self.__domain==None:
128 if self.dim==2:
129 self.__domain=finley.Rectangle(n0=self.n[0],\
130 n1=self.n[1],\
131 l0=self.l[0],\
132 l1=self.l[1],\
133 order=self.order, \
134 periodic0=self.periodic[0], \
135 periodic1=self.periodic[1], \
136 integrationOrder=self.integrationOrder)
137 else:
138 self.__domain=finley.Brick(n0=self.n[0],\
139 n1=self.n[1],\
140 n2=self.n[2],\
141 l0=self.l[0],\
142 l1=self.l[1],\
143 l2=self.l[2],\
144 order=self.order, \
145 periodic0=self.periodic[0], \
146 periodic1=self.periodic[1], \
147 periodic2=self.periodic[2], \
148 integrationOrder=self.integrationOrder)
149 return self.__domain
150
151 class UpdateGeometry(Model):
152 """
153 applies a displacement field to a domain
154
155 @ivar displacement: displacements applied to the original mesh coordinates (in).
156 @type displacement: L{escript.Vector}
157 @ivar domain: domain
158 @type domain: L{escript.Domain}
159 """
160 def __init__(self,**kwargs):
161 """
162 set-up the object
163 """
164 super(UpdateGeometry, self).__init__(**kwargs)
165 self.declareParameter(domain=None,\
166 displacement=None)
167
168
169 def doInitialization(self):
170 """
171 initialize model
172 """
173 self.__x=self.domain.getX()
174 self.__reset=True
175
176 def doStepPreprocessing(self,dt):
177 """
178 applies the current L{displacement} to mesh nodes if required.
179 """
180 if self.__reset:
181 self.trace("mesh nodes updated.")
182 self.domain.setX(self.__x+self.displacement)
183 self.__reset=False
184
185 def doStep(self,dt):
186 """
187 applies the current L{displacement} to mesh nodes.
188 """
189 self.trace("mesh nodes updated.")
190 self.domain.setX(self.__x+self.displacement)
191 self.__reset=True
192
193 def doStepPostprocessing(self,dt):
194 """
195 marks nodes as beeing updated.
196 """
197 self.__reset=False
198
199 class ConstrainerOverBox(Model):
200 """
201 Creates a characteristic function for the location of constraints
202 for all components of a value and selects the value from an initial value
203 ate these locations.
204
205 In the case that the spatial dimension is two, the arguments front and back are ignored.
206
207 @ivar domain: domain (in).
208 @ivar left: True to set a constraint at the left face of the domain (x[0]=min x[0]), default False (in).
209 @ivar right: True to set a constraint at the left face of the domain (x[0]=max x[0]), default False (in).
210 @ivar top: True to set a constraint at the left face of the domain (x[1]=min x[1]), default False (in).
211 @ivar bottom: True to set a constraint at the left face of the domain (x[1]=max x[1]), default False (in).
212 @ivar front: True to set a constraint at the left face of the domain (x[2]=min x[2]), default False (in).
213 @ivar back: True to set a constraint at the left face of the domain (x[2]=max x[2]), default False (in).
214 @ivar tol: absolute tolerance for "x=max x" condition, default 1.e-8 (in).
215 """
216 def __init__(self,**kwargs):
217 super(ConstrainerOverBox, self).__init__(**kwargs)
218 self.declareParameter(domain=None, \
219 value=None, \
220 left=False, \
221 right=False, \
222 top=False, \
223 bottom=False, \
224 front=False, \
225 back=False, \
226 tol=1.e-8)
227 self.__value_of_constraint = None
228 self.__location_of_constraint=None
229 def location_of_constraint(self):
230 """
231 return the values used to constrain a solution
232
233 @return: the mask marking the locations of the constraints
234 @rtype: L{escript.Scalar}
235 """
236 if self.__location_of_constraint == None: self.__setOutput()
237 return self.__location_of_constraint
238
239 def value_of_constraint(self):
240 """
241 return the values used to constrain a solution
242
243 @return: values to be used at the locations of the constraints. If
244 L{value} is not given C{None} is rerturned.
245 @rtype: L{escript.Scalar}
246 """
247 if self.__location_of_constraint == None: self.__setOutput()
248 return self.__value_of_constraint
249
250 def __setOutput(self):
251 x=self.domain.getX()
252 val=self.value
253 if isinstance(val, int) or isinstance(val, float):
254 shape=()
255 elif isinstance(val, list) or isinstance(val, tuple) :
256 shape=(len(val),)
257 elif isinstance(val, numarray.NumArray):
258 shape=val.shape
259 else:
260 shape=val.getShape()
261 self.__location_of_constraint=Data(0,shape,x.getFunctionSpace())
262 if self.domain.getDim()==3:
263 x0,x1,x2=x[0],x[1],x[2]
264 if self.left: self.__location_of_constraint+=whereZero(x0-inf(x0),self.tol)
265 if self.right: self.__location_of_constraint+=whereZero(x0-sup(x0),self.tol)
266 if self.front: self.__location_of_constraint+=whereZero(x1-inf(x1),self.tol)
267 if self.back: self.__location_of_constraint+=whereZero(x1-sup(x1),self.tol)
268 if self.bottom: self.__location_of_constraint+=whereZero(x2-inf(x2),self.tol)
269 if self.top: self.__location_of_constraint+=whereZero(x2-sup(x2),self.tol)
270 else:
271 x0,x1=x[0],x[1]
272 if self.left: self.__location_of_constraint+=whereZero(x0-inf(x0),self.tol)
273 if self.right: self.__location_of_constraint+=whereZero(x0-sup(x0),self.tol)
274 if self.bottom: self.__location_of_constraint+=whereZero(x1-inf(x1),self.tol)
275 if self.top: self.__location_of_constraint+=whereZero(x1-sup(x1),self.tol)
276 self.__value_of_constraint=self.__location_of_constraint*self.value
277 class ScalarConstrainerOverBox(Model):
278 """
279 Creates a characteristic function for the location of constraints
280 for a scalar value and selects the value from an initial value
281 ate these locations.
282
283 In the case that the spatial dimension is two, the arguments front and back are ignored.
284
285 @ivar domain: domain (in).
286 @ivar left: True to set a constraint at the left face of the domain (x[0]=min x[0]), default False (in).
287 @ivar right: True to set a constraint at the left face of the domain (x[0]=max x[0]), default False (in).
288 @ivar top: True to set a constraint at the left face of the domain (x[1]=min x[1]), default False (in).
289 @ivar bottom: True to set a constraint at the left face of the domain (x[1]=max x[1]), default False (in).
290 @ivar front: True to set a constraint at the left face of the domain (x[2]=min x[2]), default False (in).
291 @ivar back: True to set a constraint at the left face of the domain (x[2]=max x[2]), default False (in).
292 @ivar tol: absolute tolerance for "x=max x" condition, default 1.e-8 (in).
293 """
294 def __init__(self,**kwargs):
295 super(ScalarConstrainerOverBox, self).__init__(**kwargs)
296 self.declareParameter(domain=None, \
297 value=None, \
298 left=False, \
299 right=False, \
300 top=False, \
301 bottom=False, \
302 front=False, \
303 back=False, \
304 tol=1.e-8)
305 self.__value_of_constraint = None
306 self.__location_of_constraint=None
307 def location_of_constraint(self):
308 """
309 return the values used to constrain a solution
310
311 @return: the mask marking the locations of the constraints
312 @rtype: L{escript.Scalar}
313 """
314 if self.__location_of_constraint == None: self.__setOutput()
315 return self.__location_of_constraint
316
317 def value_of_constraint(self):
318 """
319 return the values used to constrain a solution
320
321 @return: values to be used at the locations of the constraints. If
322 L{value} is not given C{None} is rerturned.
323 @rtype: L{escript.Scalar}
324 """
325 if self.__location_of_constraint == None: self.__setOutput()
326 return self.__value_of_constraint
327
328 def __setOutput(self):
329 x=self.domain.getX()
330 self.__location_of_constraint=Scalar(0,x.getFunctionSpace())
331 if self.domain.getDim()==3:
332 x0,x1,x2=x[0],x[1],x[2]
333 if self.left: self.__location_of_constraint+=whereZero(x0-inf(x0),self.tol)
334 if self.right: self.__location_of_constraint+=whereZero(x0-sup(x0),self.tol)
335 if self.front: self.__location_of_constraint+=whereZero(x1-inf(x1),self.tol)
336 if self.back: self.__location_of_constraint+=whereZero(x1-sup(x1),self.tol)
337 if self.bottom: self.__location_of_constraint+=whereZero(x2-inf(x2),self.tol)
338 if self.top: self.__location_of_constraint+=whereZero(x2-sup(x2),self.tol)
339 else:
340 x0,x1=x[0],x[1]
341 if self.left: self.__location_of_constraint+=whereZero(x0-inf(x0),self.tol)
342 if self.right: self.__location_of_constraint+=whereZero(x0-sup(x0),self.tol)
343 if self.bottom: self.__location_of_constraint+=whereZero(x1-inf(x1),self.tol)
344 if self.top: self.__location_of_constraint+=whereZero(x1-sup(x1),self.tol)
345 if self.value:
346 self.__value_of_constraint=self.__location_of_constraint*self.value
347
348 class VectorConstrainerOverBox(Model):
349 """
350 Creates a characteristic function for the location of constraints vector value.
351 In the case that the spatial dimension is two, the arguments front and
352 back as well as the third component of each argument is ignored.
353
354 @ivar domain: domain
355 @ivar left: list of three boolean. left[i]==True sets a constraint for the i-th component at the left face of the domain (x[0]=min x[0]),
356 default [False,False,False] (in).
357 @ivar right: list of three boolean. left[i]==True sets a constraint for the i-th component at the right face of the domain (x[0]=max x[0]),
358 default [False,False,False] (in).
359 @ivar top: list of three boolean. left[i]==True sets a constraint for the i-th component at the top face of the domain (x[1]=min x[1]),
360 default [False,False,False] (in).
361 @ivar bottom: list of three boolean. left[i]==True sets a constraint for the i-th component at the bottom face of the domain (x[1]=min x[1]),
362 default [False,False,False] (in).
363 @ivar front: list of three boolean. left[i]==True sets a constraint for the i-th component at the front face of the domain (x[2]=min x[2]),
364 default [False,False,False] (in).
365 @ivar back: list of three boolean. left[i]==True sets a constraint for the i-th component at the back face of the domain (x[2]=max x[2]),
366 default [False,False,False] (in).
367 @ivar tol: absolute tolerance for "x=max x" condition, default 1.e-8 (in).
368 """
369 def __init__(self, **kwargs):
370 super(VectorConstrainerOverBox, self).__init__(**kwargs)
371 self.declareParameter(domain=None, \
372 value=None, \
373 left=[0,0,0], \
374 right=[0,0,0], \
375 top=[0,0,0], \
376 bottom=[0,0,0], \
377 front=[0,0,0], \
378 back=[0,0,0], \
379 tol=1.e-8)
380 self.__value_of_constraint = None
381 self.__location_of_constraint=None
382
383 def location_of_constraint(self):
384 """
385 return the values used to constrain a solution
386
387 @return: the mask marking the locations of the constraints
388 @rtype: L{escript.Vector}
389 """
390 if self.__location_of_constraint == None: self.__setOutput()
391 return self.__location_of_constraint
392
393 def value_of_constraint(self):
394 """
395 return the values used to constrain a solution
396
397 @return: values to be used at the locations of the constraints. If
398 L{value} is not given C{None} is rerturned.
399 @rtype: L{escript.Vector}
400 """
401 if self.__location_of_constraint == None: self.__setOutput()
402 return self.__value_of_constraint
403
404 def __setOutput(self):
405 x=self.domain.getX()
406 self.__location_of_constraint=Vector(0,x.getFunctionSpace())
407 if self.domain.getDim()==3:
408 x0,x1,x2=x[0],x[1],x[2]
409 left_mask=whereZero(x0-inf(x0),self.tol)
410 if self.left[0]: self.__location_of_constraint+=left_mask*[1.,0.,0.]
411 if self.left[1]: self.__location_of_constraint+=left_mask*[0.,1.,0.]
412 if self.left[2]: self.__location_of_constraint+=left_mask*[0.,0.,1.]
413 right_mask=whereZero(x0-sup(x0),self.tol)
414 if self.right[0]: self.__location_of_constraint+=right_mask*[1.,0.,0.]
415 if self.right[1]: self.__location_of_constraint+=right_mask*[0.,1.,0.]
416 if self.right[2]: self.__location_of_constraint+=right_mask*[0.,0.,1.]
417 front_mask=whereZero(x1-inf(x1),self.tol)
418 if self.front[0]: self.__location_of_constraint+=front_mask*[1.,0.,0.]
419 if self.front[1]: self.__location_of_constraint+=front_mask*[0.,1.,0.]
420 if self.front[2]: self.__location_of_constraint+=front_mask*[0.,0.,1.]
421 back_mask=whereZero(x1-sup(x1),self.tol)
422 if self.back[0]: self.__location_of_constraint+=back_mask*[1.,0.,0.]
423 if self.back[1]: self.__location_of_constraint+=back_mask*[0.,1.,0.]
424 if self.back[2]: self.__location_of_constraint+=back_mask*[0.,0.,1.]
425 bottom_mask=whereZero(x2-inf(x2),self.tol)
426 if self.bottom[0]: self.__location_of_constraint+=bottom_mask*[1.,0.,0.]
427 if self.bottom[1]: self.__location_of_constraint+=bottom_mask*[0.,1.,0.]
428 if self.bottom[2]: self.__location_of_constraint+=bottom_mask*[0.,0.,1.]
429 top_mask=whereZero(x2-sup(x2),self.tol)
430 if self.top[0]: self.__location_of_constraint+=top_mask*[1.,0.,0.]
431 if self.top[1]: self.__location_of_constraint+=top_mask*[0.,1.,0.]
432 if self.top[2]: self.__location_of_constraint+=top_mask*[0.,0.,1.]
433 if self.value:
434 self.__value_of_constraint=self.__location_of_constraint*self.value
435 else:
436 x0,x1=x[0],x[1]
437 left_mask=whereZero(x0-inf(x0),self.tol)
438 if self.left[0]: self.__location_of_constraint+=left_mask*[1.,0.]
439 if self.left[1]: self.__location_of_constraint+=left_mask*[0.,1.]
440 right_mask=whereZero(x0-sup(x0),self.tol)
441 if self.right[0]: self.__location_of_constraint+=right_mask*[1.,0.]
442 if self.right[1]: self.__location_of_constraint+=right_mask*[0.,1.]
443 bottom_mask=whereZero(x1-inf(x1),self.tol)
444 if self.bottom[0]: self.__location_of_constraint+=bottom_mask*[1.,0.]
445 if self.bottom[1]: self.__location_of_constraint+=bottom_mask*[0.,1.]
446 top_mask=whereZero(x1-sup(x1),self.tol)
447 if self.top[0]: self.__location_of_constraint+=top_mask*[1.,0.]
448 if self.top[1]: self.__location_of_constraint+=top_mask*[0.,1.]
449 if self.value:
450 self.__value_of_constraint=self.__location_of_constraint*self.value[:2]
451
452 # vim: expandtab shiftwidth=4:

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26