modellib/py_src/geometry.py

# $Id$

__copyright__="""  Copyright (c) 2006 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"""


from esys.escript import *
from esys.escript.modelframe import Model,ParameterSet
from esys.pycad import TagMap
from esys import finley

class FinleyReader(ParameterSet):
       """
       reads finley mesh file.

       @ivar source: mesh file in finley or gmsh format
       @type source: C{DataSource}
       @ivar intergrationOrder: integration order, default -1 (in).
       @type intergrationOrder: C{int}
       @ivar reducedIntegrationOrder: reduced integration order, default -1 (in).
       @type reducedIntegrationOrder: C{int}
       @ivar optimizeLabeling: switches on optimization of the labeling of the nodes
       @type optimizeLabeling: C{bool}
       """
       def __init__(self,**kwargs):
          """
          initializes the object
          """
          super(FinleyReader,self).__init__(**kwargs)
          self.declareParameter(source="none",
                                region_tag_map_source=None,
                                surface_tag_map_source=None,
                                optimizeLabeling=True,
                                reducedIntegrationOrder=-1,
                                integrationOrder=-1)
          self.__domain=None
          self.__region_tag_map=None
          self.__surface_tag_map=None


       def domain(self):
          """
          returns the domain

          @return: the domain
          @rtype: L{Domain}
          """
          if self.__domain == None:
             if  self.source.fileformat == "fly":
                self.__domain=finley.ReadMesh(self.source.getLocalFileName(),self.integrationOrder) 
             elif self.source.fileformat == "gmsh":
                self.__domain=finley.ReadGmsh(self.source.getLocalFileName(),self.integrationOrder,self.reducedIntegrationOrder, self.optimizeLabeling) 
             else:
                raise TypeError("unknown mesh file format %s."%self.source.fileformat)
             self.trace("mesh read from %s in %s format."%(self.source.getLocalFileName(), self.source.fileformat))           
          return self.__domain

       def region_tag_map(self):
          """
          returns the map from regional tag names to tag integers used in the mesh

          @return: the tag map
          @rtype: L{TagMap}
          """
          if self.__region_tag_map == None:
               self.__region_tag_map = TagMap()
               if  self.region_tag_map_source != None:
                   self.__region_tag_map.fillFromXML(open(self.region_tag_map_source.getLocalFileName()))
               self.trace("region tag map read from %s in %s format."%(self.region_tag_map_source.getLocalFileName(), self.region_tag_map_source.fileformat))           
          return self.__region_tag_map

       def surface_tag_map(self):
          """
          returns the map from surface tag names to tag integers used in the mesh

          @return: the tag map
          @rtype: L{TagMap}
          """
          if self.__surface_tag_map == None:
               self.__surface_tag_map = TagMap()
               if  self.surface_tag_map_source != None:
                   self.__surface_tag_map.fillFromXML(open(self.surface_tag_map_source.getLocalFileName()))
               self.trace("surface tag map read from %s in %s format."%(self.surface_tag_map_source.getLocalFileName(), self.surface_tag_map_source.fileformat))           
          return self.__surface_tag_map
          
                       
class RectangularDomain(ParameterSet):
       """
       Generates a mesh over a rectangular domain finley.

       @ivar dim: spatial dimension, default =2 (in).
       @type dim: spatial dimension
       @ivar l: spatial lengths, default [1.,1.,1.] (in).
       @type l: C{list} of C{floats}s
       @ivar n: number of elements, default [10,10,10] (in).
       @type n: C{list} of C{int}s
       @ivar order: element order, default 1 (in).
       @type order: C{int}
       @ivar periodic: flags for periodicity, default [False,False,False] (in).
       @type periodic: C{list} of C{bool}s
       @ivar intergrationOrder: integration order, default -1 (in).
       @type intergrationOrder: C{int}
       """
       def __init__(self,**kwargs):
           """
           initializes the object
           """
           super(RectangularDomain,self).__init__(**kwargs)
           self.declareParameter(dim=2,\
                                 l=[1.,1.,1.],\
                                 n=[10,10,10], \
                 order=1,\
                                 periodic=[False,False,False],
                                 integrationOrder=-1)
           self.__domain=None

       def domain(self):
           """
           returns the domain

           @return: the domain
           @rtype: L{Domain}
           """
           if self.__domain==None:
              if self.dim==2:
                   self.__domain=finley.Rectangle(n0=self.n[0],\
                                                n1=self.n[1],\
                                                l0=self.l[0],\
                                                l1=self.l[1],\
                                                order=self.order, \
                                                periodic0=self.periodic[0], \
                                                periodic1=self.periodic[1], \
                                                integrationOrder=self.integrationOrder)
              else:
                   self.__domain=finley.Brick(n0=self.n[0],\
                                            n1=self.n[1],\
                                            n2=self.n[2],\
                                            l0=self.l[0],\
                                            l1=self.l[1],\
                                            l2=self.l[2],\
                                            order=self.order, \
                                            periodic0=self.periodic[0], \
                                            periodic1=self.periodic[1], \
                                            periodic2=self.periodic[2], \
                                            integrationOrder=self.integrationOrder)
           return self.__domain

class UpdateGeometry(Model):
      """
      applies a displacement field to a domain
      
      @ivar displacement: displacements applied to the original mesh coordinates (in).
      @type displacement: L{escript.Vector}
      @ivar domain: domain
      @type domain: L{escript.Domain}
      """
      def __init__(self,**kwargs):
           """
           set-up the object
           """
           super(UpdateGeometry, self).__init__(**kwargs)
           self.declareParameter(domain=None,\
                                 displacement=None)


      def doInitialization(self):
         """
         initialize model
         """
         self.__x=self.domain.getX()
         self.__reset=True
         
      def doStepPreprocessing(self,dt):
         """
         applies the current L{displacement} to mesh nodes if required.
         """
         if self.__reset:
            self.trace("mesh nodes updated.")
            self.domain.setX(self.__x+self.displacement)
         self.__reset=False

      def doStep(self,dt):
         """
         applies the current L{displacement} to mesh nodes. 
         """
         self.trace("mesh nodes updated.")
         self.domain.setX(self.__x+self.displacement)
         self.__reset=True

      def doStepPostprocessing(self,dt):
         """
         marks nodes as beeing updated.
         """
         self.__reset=False

class ConstrainerOverBox(Model):
      """
      Creates a characteristic function for the location of constraints 
      for all components of a value and selects the value from an initial value 
      ate these locations.

      In the case that the spatial dimension is two, the arguments front and back are ignored.

      @ivar domain: domain (in).
      @ivar left:  True to set a constraint at the left face of the domain (x[0]=min x[0]), default False (in).
      @ivar right: True to set a constraint at the left face of the domain (x[0]=max x[0]), default False (in).
      @ivar top: True to set a constraint at the left face of the domain (x[1]=min x[1]), default False (in).
      @ivar bottom: True to set a constraint at the left face of the domain (x[1]=max x[1]), default False (in).
      @ivar front: True to set a constraint at the left face of the domain (x[2]=min x[2]), default False (in).
      @ivar back: True to set a constraint at the left face of the domain (x[2]=max x[2]), default False (in).
      @ivar tol: absolute tolerance for "x=max x" condition, default 1.e-8 (in).
      """
      def __init__(self,**kwargs):
           super(ConstrainerOverBox, self).__init__(**kwargs)
           self.declareParameter(domain=None, \
                                 value=None,  \
                                 left=False, \
                                 right=False, \
                                 top=False, \
                                 bottom=False, \
                                 front=False, \
                                 back=False, \
                                 tol=1.e-8)
           self.__value_of_constraint = None
           self.__location_of_constraint=None
      def location_of_constraint(self):
          """
          return the values used to constrain a solution

          @return: the mask marking the locations of the constraints
          @rtype: L{escript.Scalar}
          """
          if self.__location_of_constraint == None: self.__setOutput()
          return self.__location_of_constraint
         
      def value_of_constraint(self):
          """
          return the values used to constrain a solution

          @return: values to be used at the locations of the constraints. If 
                  L{value} is not given C{None} is rerturned.
          @rtype: L{escript.Scalar}
          """
          if self.__location_of_constraint == None: self.__setOutput()
          return self.__value_of_constraint
         
      def __setOutput(self):
          if self.__location_of_constraint == None:
             x=self.domain.getX()
             val=self.value
             if isinstance(val, int) or isinstance(val, float):
                shape=()
             elif isinstance(val, list) or isinstance(val, tuple) :
                shape=(len(val),)
             elif isinstance(val, numarray.NumArray):
                 shape=val.shape
             elif val == None:
                  shape=()
             else: 
                 shape=val.getShape()
             self.__location_of_constraint=Data(0,shape,x.getFunctionSpace())
             if self.domain.getDim()==3:
                   x0,x1,x2=x[0],x[1],x[2]
                   if self.left: self.__location_of_constraint+=whereZero(x0-inf(x0),self.tol)
                   if self.right: self.__location_of_constraint+=whereZero(x0-sup(x0),self.tol)
                   if self.front: self.__location_of_constraint+=whereZero(x1-inf(x1),self.tol)
                   if self.back: self.__location_of_constraint+=whereZero(x1-sup(x1),self.tol)
                   if self.bottom: self.__location_of_constraint+=whereZero(x2-inf(x2),self.tol)
                   if self.top: self.__location_of_constraint+=whereZero(x2-sup(x2),self.tol)
             else:
                   x0,x1=x[0],x[1]
                   if self.left: self.__location_of_constraint+=whereZero(x0-inf(x0),self.tol)
                   if self.right: self.__location_of_constraint+=whereZero(x0-sup(x0),self.tol)
                   if self.bottom: self.__location_of_constraint+=whereZero(x1-inf(x1),self.tol)
                   if self.top: self.__location_of_constraint+=whereZero(x1-sup(x1),self.tol)
             if not self.value == None:
                   self.__value_of_constraint=self.__location_of_constraint*self.value
class ScalarConstrainerOverBox(Model):
      """
      Creates a characteristic function for the location of constraints 
      for a scalar value and selects the value from an initial value 
      ate these locations.

      In the case that the spatial dimension is two, the arguments front and back are ignored.

      @ivar domain: domain (in).
      @ivar left:  True to set a constraint at the left face of the domain (x[0]=min x[0]), default False (in).
      @ivar right: True to set a constraint at the left face of the domain (x[0]=max x[0]), default False (in).
      @ivar top: True to set a constraint at the left face of the domain (x[1]=min x[1]), default False (in).
      @ivar bottom: True to set a constraint at the left face of the domain (x[1]=max x[1]), default False (in).
      @ivar front: True to set a constraint at the left face of the domain (x[2]=min x[2]), default False (in).
      @ivar back: True to set a constraint at the left face of the domain (x[2]=max x[2]), default False (in).
      @ivar tol: absolute tolerance for "x=max x" condition, default 1.e-8 (in).
      """
      def __init__(self,**kwargs):
           super(ScalarConstrainerOverBox, self).__init__(**kwargs)
           self.declareParameter(domain=None, \
                                 value=None,  \
                                 left=False, \
                                 right=False, \
                                 top=False, \
                                 bottom=False, \
                                 front=False, \
                                 back=False, \
                                 tol=1.e-8)
           self.__value_of_constraint = None
           self.__location_of_constraint=None
      def location_of_constraint(self):
          """
          return the values used to constrain a solution

          @return: the mask marking the locations of the constraints
          @rtype: L{escript.Scalar}
          """
          if self.__location_of_constraint == None: self.__setOutput()
          return self.__location_of_constraint
         
      def value_of_constraint(self):
          """
          return the values used to constrain a solution

          @return: values to be used at the locations of the constraints. If 
                  L{value} is not given C{None} is rerturned.
          @rtype: L{escript.Scalar}
          """
          if self.__location_of_constraint == None: self.__setOutput()
          return self.__value_of_constraint
         
      def __setOutput(self):
          x=self.domain.getX()
          self.__location_of_constraint=Scalar(0,x.getFunctionSpace())
          if self.domain.getDim()==3:
                x0,x1,x2=x[0],x[1],x[2]
                if self.left: self.__location_of_constraint+=whereZero(x0-inf(x0),self.tol)
                if self.right: self.__location_of_constraint+=whereZero(x0-sup(x0),self.tol)
                if self.front: self.__location_of_constraint+=whereZero(x1-inf(x1),self.tol)
                if self.back: self.__location_of_constraint+=whereZero(x1-sup(x1),self.tol)
                if self.bottom: self.__location_of_constraint+=whereZero(x2-inf(x2),self.tol)
                if self.top: self.__location_of_constraint+=whereZero(x2-sup(x2),self.tol)
          else:
                x0,x1=x[0],x[1]
                if self.left: self.__location_of_constraint+=whereZero(x0-inf(x0),self.tol)
                if self.right: self.__location_of_constraint+=whereZero(x0-sup(x0),self.tol)
                if self.bottom: self.__location_of_constraint+=whereZero(x1-inf(x1),self.tol)
                if self.top: self.__location_of_constraint+=whereZero(x1-sup(x1),self.tol)
          if not self.value == None:
              self.__value_of_constraint=self.__location_of_constraint*self.value

class VectorConstrainerOverBox(Model):
      """
      Creates a characteristic function for the location of constraints vector value.
      In the case that the spatial dimension is two, the arguments front and
      back as well as the third component of each argument is ignored.

      @ivar domain: domain
      @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]),
                       default [False,False,False] (in).
      @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]), 
                default [False,False,False] (in).
      @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]), 
                default [False,False,False] (in).
      @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]), 
                default [False,False,False] (in).
      @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]), 
                default [False,False,False] (in).
      @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]), 
                default [False,False,False] (in).
      @ivar tol: absolute tolerance for "x=max x" condition, default 1.e-8 (in).
      """
      def __init__(self, **kwargs):
           super(VectorConstrainerOverBox, self).__init__(**kwargs)
           self.declareParameter(domain=None, \
                                 value=None,  \
                                 left=[0,0,0],  \
                                 right=[0,0,0],  \
                                 top=[0,0,0],  \
                                 bottom=[0,0,0],  \
                                 front=[0,0,0], \
                                 back=[0,0,0], \
                                 tol=1.e-8)
           self.__value_of_constraint = None
           self.__location_of_constraint=None

      def location_of_constraint(self):
          """
          return the values used to constrain a solution

          @return: the mask marking the locations of the constraints
          @rtype: L{escript.Vector}
          """
          if self.__location_of_constraint == None: self.__setOutput()
          return self.__location_of_constraint
         
      def value_of_constraint(self):
          """
          return the values used to constrain a solution

          @return: values to be used at the locations of the constraints. If 
                  L{value} is not given C{None} is rerturned.
          @rtype: L{escript.Vector}
          """
          if self.__location_of_constraint == None: self.__setOutput()
          return self.__value_of_constraint
         
      def __setOutput(self):
          x=self.domain.getX()
          self.__location_of_constraint=Vector(0,x.getFunctionSpace())
          if self.domain.getDim()==3:
             x0,x1,x2=x[0],x[1],x[2]
             left_mask=whereZero(x0-inf(x0),self.tol)
             if self.left[0]: self.__location_of_constraint+=left_mask*[1.,0.,0.]
             if self.left[1]: self.__location_of_constraint+=left_mask*[0.,1.,0.]
             if self.left[2]: self.__location_of_constraint+=left_mask*[0.,0.,1.]
             right_mask=whereZero(x0-sup(x0),self.tol)
             if self.right[0]: self.__location_of_constraint+=right_mask*[1.,0.,0.]
             if self.right[1]: self.__location_of_constraint+=right_mask*[0.,1.,0.]
             if self.right[2]: self.__location_of_constraint+=right_mask*[0.,0.,1.]
             front_mask=whereZero(x1-inf(x1),self.tol)
             if self.front[0]: self.__location_of_constraint+=front_mask*[1.,0.,0.]
             if self.front[1]: self.__location_of_constraint+=front_mask*[0.,1.,0.]
             if self.front[2]: self.__location_of_constraint+=front_mask*[0.,0.,1.]
             back_mask=whereZero(x1-sup(x1),self.tol)
             if self.back[0]: self.__location_of_constraint+=back_mask*[1.,0.,0.]
             if self.back[1]: self.__location_of_constraint+=back_mask*[0.,1.,0.]
             if self.back[2]: self.__location_of_constraint+=back_mask*[0.,0.,1.]
             bottom_mask=whereZero(x2-inf(x2),self.tol)
             if self.bottom[0]: self.__location_of_constraint+=bottom_mask*[1.,0.,0.]
             if self.bottom[1]: self.__location_of_constraint+=bottom_mask*[0.,1.,0.]
             if self.bottom[2]: self.__location_of_constraint+=bottom_mask*[0.,0.,1.]
             top_mask=whereZero(x2-sup(x2),self.tol)
             if self.top[0]: self.__location_of_constraint+=top_mask*[1.,0.,0.]
             if self.top[1]: self.__location_of_constraint+=top_mask*[0.,1.,0.]
             if self.top[2]: self.__location_of_constraint+=top_mask*[0.,0.,1.]
             if not self.value == None:
                self.__value_of_constraint=self.__location_of_constraint*self.value
          else:
             x0,x1=x[0],x[1]
             left_mask=whereZero(x0-inf(x0),self.tol)
             if self.left[0]: self.__location_of_constraint+=left_mask*[1.,0.]
             if self.left[1]: self.__location_of_constraint+=left_mask*[0.,1.]
             right_mask=whereZero(x0-sup(x0),self.tol)
             if self.right[0]: self.__location_of_constraint+=right_mask*[1.,0.]
             if self.right[1]: self.__location_of_constraint+=right_mask*[0.,1.]
             bottom_mask=whereZero(x1-inf(x1),self.tol)
             if self.bottom[0]: self.__location_of_constraint+=bottom_mask*[1.,0.]
             if self.bottom[1]: self.__location_of_constraint+=bottom_mask*[0.,1.]
             top_mask=whereZero(x1-sup(x1),self.tol)
             if self.top[0]: self.__location_of_constraint+=top_mask*[1.,0.]
             if self.top[1]: self.__location_of_constraint+=top_mask*[0.,1.]
             if not self.value == None:
                self.__value_of_constraint=self.__location_of_constraint*self.value[:2]

# vim: expandtab shiftwidth=4:
1	jgs	127	# $Id$
2
3	elspeth	628	__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	jgs	127
9	elspeth	628
10	jgs	149	from esys.escript import *
11			from esys.escript.modelframe import Model,ParameterSet
12	gross	953	from esys.pycad import TagMap
13	jgs	149	from esys import finley
14	jgs	127
15	gross	823	class FinleyReader(ParameterSet):
16	jgs	149	"""
17	gross	819	reads finley mesh file.
18	jgs	147
19	gross	938	@ivar source: mesh file in finley or gmsh format
20			@type source: C{DataSource}
21	gross	823	@ivar intergrationOrder: integration order, default -1 (in).
22			@type intergrationOrder: C{int}
23	gross	938	@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	gross	398	"""
28	gross	918	def __init__(self,**kwargs):
29	gross	823	"""
30			initializes the object
31			"""
32	gross	918	super(FinleyReader,self).__init__(**kwargs)
33	gross	823	self.declareParameter(source="none",
34	gross	953	region_tag_map_source=None,
35			surface_tag_map_source=None,
36	gross	938	optimizeLabeling=True,
37			reducedIntegrationOrder=-1,
38			integrationOrder=-1)
39	gross	823	self.__domain=None
40	gross	953	self.__region_tag_map=None
41			self.__surface_tag_map=None
42	gross	398
43	gross	953
44	gross	823	def domain(self):
45			"""
46			returns the domain
47
48			@return: the domain
49			@rtype: L{Domain}
50			"""
51	gross	911	if self.__domain == None:
52	gross	938	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	gross	950	self.trace("mesh read from %s in %s format."%(self.source.getLocalFileName(), self.source.fileformat))
59	gross	823	return self.__domain
60	gross	953
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	gross	823
89	gross	398
90	gross	823	class RectangularDomain(ParameterSet):
91	gross	398	"""
92			Generates a mesh over a rectangular domain finley.
93
94	gross	819	@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	gross	823	@ivar intergrationOrder: integration order, default -1 (in).
105			@type intergrationOrder: C{int}
106	jgs	147	"""
107	gross	918	def __init__(self,**kwargs):
108	gross	819	"""
109			initializes the object
110			"""
111	gross	918	super(RectangularDomain,self).__init__(**kwargs)
112	jgs	147	self.declareParameter(dim=2,\
113	jgs	127	l=[1.,1.,1.],\
114			n=[10,10,10], \
115			order=1,\
116	gross	823	periodic=[False,False,False],
117			integrationOrder=-1)
118			self.__domain=None
119	jgs	127
120	gross	823	def domain(self):
121	gross	819	"""
122	gross	823	returns the domain
123
124			@return: the domain
125			@rtype: L{Domain}
126	gross	819	"""
127	gross	911	if self.__domain==None:
128	gross	823	if self.dim==2:
129	gross	911	self.__domain=finley.Rectangle(n0=self.n[0],\
130	gross	823	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	jgs	147
151	gross	823	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	gross	918	def __init__(self,**kwargs):
161	gross	823	"""
162			set-up the object
163			"""
164	gross	918	super(UpdateGeometry, self).__init__(**kwargs)
165	gross	823	self.declareParameter(domain=None,\
166			displacement=None)
167	jgs	147
168	gross	823
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	gross	953	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	gross	993	if self.__location_of_constraint == None:
252			x=self.domain.getX()
253			val=self.value
254			if isinstance(val, int) or isinstance(val, float):
255			shape=()
256			elif isinstance(val, list) or isinstance(val, tuple) :
257			shape=(len(val),)
258			elif isinstance(val, numarray.NumArray):
259			shape=val.shape
260			elif val == None:
261			shape=()
262			else:
263			shape=val.getShape()
264			self.__location_of_constraint=Data(0,shape,x.getFunctionSpace())
265			if self.domain.getDim()==3:
266			x0,x1,x2=x[0],x[1],x[2]
267			if self.left: self.__location_of_constraint+=whereZero(x0-inf(x0),self.tol)
268			if self.right: self.__location_of_constraint+=whereZero(x0-sup(x0),self.tol)
269			if self.front: self.__location_of_constraint+=whereZero(x1-inf(x1),self.tol)
270			if self.back: self.__location_of_constraint+=whereZero(x1-sup(x1),self.tol)
271			if self.bottom: self.__location_of_constraint+=whereZero(x2-inf(x2),self.tol)
272			if self.top: self.__location_of_constraint+=whereZero(x2-sup(x2),self.tol)
273			else:
274			x0,x1=x[0],x[1]
275			if self.left: self.__location_of_constraint+=whereZero(x0-inf(x0),self.tol)
276			if self.right: self.__location_of_constraint+=whereZero(x0-sup(x0),self.tol)
277			if self.bottom: self.__location_of_constraint+=whereZero(x1-inf(x1),self.tol)
278			if self.top: self.__location_of_constraint+=whereZero(x1-sup(x1),self.tol)
279			if not self.value == None:
280			self.__value_of_constraint=self.__location_of_constraint*self.value
281	gross	911	class ScalarConstrainerOverBox(Model):
282	gross	823	"""
283			Creates a characteristic function for the location of constraints
284			for a scalar value and selects the value from an initial value
285			ate these locations.
286	jgs	127
287	gross	823	In the case that the spatial dimension is two, the arguments front and back are ignored.
288	jgs	127
289	gross	823	@ivar domain: domain (in).
290			@ivar left: True to set a constraint at the left face of the domain (x[0]=min x[0]), default False (in).
291			@ivar right: True to set a constraint at the left face of the domain (x[0]=max x[0]), default False (in).
292			@ivar top: True to set a constraint at the left face of the domain (x[1]=min x[1]), default False (in).
293			@ivar bottom: True to set a constraint at the left face of the domain (x[1]=max x[1]), default False (in).
294			@ivar front: True to set a constraint at the left face of the domain (x[2]=min x[2]), default False (in).
295			@ivar back: True to set a constraint at the left face of the domain (x[2]=max x[2]), default False (in).
296			@ivar tol: absolute tolerance for "x=max x" condition, default 1.e-8 (in).
297			"""
298	gross	918	def __init__(self,**kwargs):
299			super(ScalarConstrainerOverBox, self).__init__(**kwargs)
300	jgs	127	self.declareParameter(domain=None, \
301	gross	823	value=None, \
302	jgs	147	left=False, \
303			right=False, \
304			top=False, \
305			bottom=False, \
306			front=False, \
307	gross	819	back=False, \
308	gross	823	tol=1.e-8)
309			self.__value_of_constraint = None
310			self.__location_of_constraint=None
311			def location_of_constraint(self):
312	jgs	149	"""
313	gross	823	return the values used to constrain a solution
314
315			@return: the mask marking the locations of the constraints
316			@rtype: L{escript.Scalar}
317	jgs	149	"""
318	gross	908	if self.__location_of_constraint == None: self.__setOutput()
319	gross	823	return self.__location_of_constraint
320
321			def value_of_constraint(self):
322			"""
323			return the values used to constrain a solution
324
325			@return: values to be used at the locations of the constraints. If
326			L{value} is not given C{None} is rerturned.
327			@rtype: L{escript.Scalar}
328			"""
329	gross	908	if self.__location_of_constraint == None: self.__setOutput()
330	gross	823	return self.__value_of_constraint
331
332			def __setOutput(self):
333	gross	819	x=self.domain.getX()
334	gross	823	self.__location_of_constraint=Scalar(0,x.getFunctionSpace())
335	gross	819	if self.domain.getDim()==3:
336			x0,x1,x2=x[0],x[1],x[2]
337	gross	823	if self.left: self.__location_of_constraint+=whereZero(x0-inf(x0),self.tol)
338			if self.right: self.__location_of_constraint+=whereZero(x0-sup(x0),self.tol)
339			if self.front: self.__location_of_constraint+=whereZero(x1-inf(x1),self.tol)
340			if self.back: self.__location_of_constraint+=whereZero(x1-sup(x1),self.tol)
341			if self.bottom: self.__location_of_constraint+=whereZero(x2-inf(x2),self.tol)
342			if self.top: self.__location_of_constraint+=whereZero(x2-sup(x2),self.tol)
343	gross	819	else:
344			x0,x1=x[0],x[1]
345	gross	823	if self.left: self.__location_of_constraint+=whereZero(x0-inf(x0),self.tol)
346			if self.right: self.__location_of_constraint+=whereZero(x0-sup(x0),self.tol)
347			if self.bottom: self.__location_of_constraint+=whereZero(x1-inf(x1),self.tol)
348			if self.top: self.__location_of_constraint+=whereZero(x1-sup(x1),self.tol)
349	gross	993	if not self.value == None:
350	gross	823	self.__value_of_constraint=self.__location_of_constraint*self.value
351	jgs	147
352	gross	911	class VectorConstrainerOverBox(Model):
353	jgs	149	"""
354	gross	819	Creates a characteristic function for the location of constraints vector value.
355			In the case that the spatial dimension is two, the arguments front and
356			back as well as the third component of each argument is ignored.
357	jgs	127
358	gross	819	@ivar domain: domain
359			@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]),
360			default [False,False,False] (in).
361			@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]),
362			default [False,False,False] (in).
363			@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]),
364			default [False,False,False] (in).
365			@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]),
366			default [False,False,False] (in).
367			@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]),
368			default [False,False,False] (in).
369			@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]),
370			default [False,False,False] (in).
371			@ivar tol: absolute tolerance for "x=max x" condition, default 1.e-8 (in).
372	jgs	127	"""
373	gross	918	def __init__(self, **kwargs):
374			super(VectorConstrainerOverBox, self).__init__(**kwargs)
375	jgs	127	self.declareParameter(domain=None, \
376	gross	823	value=None, \
377	jgs	147	left=[0,0,0], \
378			right=[0,0,0], \
379			top=[0,0,0], \
380			bottom=[0,0,0], \
381	gross	819	front=[0,0,0], \
382	gross	821	back=[0,0,0], \
383	gross	823	tol=1.e-8)
384			self.__value_of_constraint = None
385			self.__location_of_constraint=None
386
387			def location_of_constraint(self):
388	jgs	149	"""
389	gross	823	return the values used to constrain a solution
390
391			@return: the mask marking the locations of the constraints
392			@rtype: L{escript.Vector}
393	jgs	149	"""
394	gross	908	if self.__location_of_constraint == None: self.__setOutput()
395	gross	823	return self.__location_of_constraint
396
397			def value_of_constraint(self):
398			"""
399			return the values used to constrain a solution
400	jgs	149
401	gross	823	@return: values to be used at the locations of the constraints. If
402			L{value} is not given C{None} is rerturned.
403			@rtype: L{escript.Vector}
404			"""
405	gross	908	if self.__location_of_constraint == None: self.__setOutput()
406	gross	823	return self.__value_of_constraint
407
408			def __setOutput(self):
409			x=self.domain.getX()
410			self.__location_of_constraint=Vector(0,x.getFunctionSpace())
411			if self.domain.getDim()==3:
412			x0,x1,x2=x[0],x[1],x[2]
413			left_mask=whereZero(x0-inf(x0),self.tol)
414			if self.left[0]: self.__location_of_constraint+=left_mask*[1.,0.,0.]
415			if self.left[1]: self.__location_of_constraint+=left_mask*[0.,1.,0.]
416			if self.left[2]: self.__location_of_constraint+=left_mask*[0.,0.,1.]
417			right_mask=whereZero(x0-sup(x0),self.tol)
418			if self.right[0]: self.__location_of_constraint+=right_mask*[1.,0.,0.]
419			if self.right[1]: self.__location_of_constraint+=right_mask*[0.,1.,0.]
420			if self.right[2]: self.__location_of_constraint+=right_mask*[0.,0.,1.]
421			front_mask=whereZero(x1-inf(x1),self.tol)
422			if self.front[0]: self.__location_of_constraint+=front_mask*[1.,0.,0.]
423			if self.front[1]: self.__location_of_constraint+=front_mask*[0.,1.,0.]
424			if self.front[2]: self.__location_of_constraint+=front_mask*[0.,0.,1.]
425			back_mask=whereZero(x1-sup(x1),self.tol)
426			if self.back[0]: self.__location_of_constraint+=back_mask*[1.,0.,0.]
427			if self.back[1]: self.__location_of_constraint+=back_mask*[0.,1.,0.]
428			if self.back[2]: self.__location_of_constraint+=back_mask*[0.,0.,1.]
429			bottom_mask=whereZero(x2-inf(x2),self.tol)
430			if self.bottom[0]: self.__location_of_constraint+=bottom_mask*[1.,0.,0.]
431			if self.bottom[1]: self.__location_of_constraint+=bottom_mask*[0.,1.,0.]
432			if self.bottom[2]: self.__location_of_constraint+=bottom_mask*[0.,0.,1.]
433			top_mask=whereZero(x2-sup(x2),self.tol)
434			if self.top[0]: self.__location_of_constraint+=top_mask*[1.,0.,0.]
435			if self.top[1]: self.__location_of_constraint+=top_mask*[0.,1.,0.]
436			if self.top[2]: self.__location_of_constraint+=top_mask*[0.,0.,1.]
437	gross	993	if not self.value == None:
438	gross	823	self.__value_of_constraint=self.__location_of_constraint*self.value
439			else:
440			x0,x1=x[0],x[1]
441			left_mask=whereZero(x0-inf(x0),self.tol)
442			if self.left[0]: self.__location_of_constraint+=left_mask*[1.,0.]
443			if self.left[1]: self.__location_of_constraint+=left_mask*[0.,1.]
444			right_mask=whereZero(x0-sup(x0),self.tol)
445			if self.right[0]: self.__location_of_constraint+=right_mask*[1.,0.]
446			if self.right[1]: self.__location_of_constraint+=right_mask*[0.,1.]
447			bottom_mask=whereZero(x1-inf(x1),self.tol)
448			if self.bottom[0]: self.__location_of_constraint+=bottom_mask*[1.,0.]
449			if self.bottom[1]: self.__location_of_constraint+=bottom_mask*[0.,1.]
450			top_mask=whereZero(x1-sup(x1),self.tol)
451			if self.top[0]: self.__location_of_constraint+=top_mask*[1.,0.]
452			if self.top[1]: self.__location_of_constraint+=top_mask*[0.,1.]
453	gross	993	if not self.value == None:
454	gross	823	self.__value_of_constraint=self.__location_of_constraint*self.value[:2]
455
456	jgs	149	# vim: expandtab shiftwidth=4:
Name	Value
svn:eol-style	native
svn:keywords	Author Date Id Revision