escript/py_src/pdetools.py

# $Id$

"""
Provides some tools related to PDEs. 

Currently includes:
    - Projector - to project a discontinuous
"""

import escript
import linearPDEs
import numarray
import util

class Projector:
  """
  The Projector is a factory which projects a discontiuous function onto a
  continuous function on the a given domain.
  """
  def __init__(self, domain, reduce = True, fast=True):
    """
    Create a continuous function space projector for a domain.

    @param domain: Domain of the projection.
    @param reduce: Flag to reduce projection order (default is True)
    @param fast: Flag to use a fast method based on matrix lumping (default is true)
    """
    self.__pde = linearPDEs.LinearPDE(domain)
    if fast:
      self.__pde.setSolverMethod(linearPDEs.LinearPDE.LUMPING)
    self.__pde.setSymmetryOn()
    self.__pde.setReducedOrderTo(reduce)
    self.__pde.setValue(D = 1.)
    return

  def __del__(self):
    return

  def __call__(self, input_data):
    """
    Projects input_data onto a continuous function

    @param input_data: The input_data to be projected.
    """
    out=escript.Data(0.,input_data.getShape(),what=escript.ContinuousFunction(self.__pde.getDomain()))
    if input_data.getRank()==0:
        self.__pde.setValue(Y = input_data)
        out=self.__pde.getSolution()
    elif input_data.getRank()==1:
        for i0 in range(input_data.getShape()[0]):
           self.__pde.setValue(Y = input_data[i0])
           out[i0]=self.__pde.getSolution()
    elif input_data.getRank()==2:
        for i0 in range(input_data.getShape()[0]):
           for i1 in range(input_data.getShape()[1]):
              self.__pde.setValue(Y = input_data[i0,i1])
              out[i0,i1]=self.__pde.getSolution()
    elif input_data.getRank()==3:
        for i0 in range(input_data.getShape()[0]):
           for i1 in range(input_data.getShape()[1]):
              for i2 in range(input_data.getShape()[2]):
                 self.__pde.setValue(Y = input_data[i0,i1,i2])
                 out[i0,i1,i2]=self.__pde.getSolution()
    else:
        for i0 in range(input_data.getShape()[0]):
           for i1 in range(input_data.getShape()[1]):
              for i2 in range(input_data.getShape()[2]):
                 for i3 in range(input_data.getShape()[3]):
                    self.__pde.setValue(Y = input_data[i0,i1,i2,i3])
                    out[i0,i1,i2,i3]=self.__pde.getSolution()
    return out


class Locator:
     """
     Locator provides access to the values of data objects at a given
     spatial coordinate x.  
     
     In fact, a Locator object finds the sample in the set of samples of a
     given function space or domain where which is closest to the given
     point x. 
     """

     def __init__(self,where,x=numarray.zeros((3,))):
       """
       Initializes a Locator to access values in Data objects on the Doamin 
       or FunctionSpace where for the sample point which
       closest to the given point x.
       """
       if isinstance(where,escript.FunctionSpace):
          self.__function_space=where
       else:
          self.__function_space=escript.ContinuousFunction(where)
       self.__id=util.length(x[:self.__function_space.getDim()]-self.__function_space.getX()).mindp()

     def __str__(self):
       """
       Returns the coordinates of the Locator as a string.
       """
       return "<Locator %s>"%str(self.getX())

     def getFunctionSpace(self):
        """
    Returns the function space of the Locator.
    """
        return self.__function_space

     def getId(self):
        """
    Returns the identifier of the location.
    """
        return self.__id

     def getX(self):
        """
    Returns the exact coordinates of the Locator.
    """
        return self(self.getFunctionSpace().getX())

     def __call__(self,data):
        """
    Returns the value of data at the Locator of a Data object otherwise 
    the object is returned.
    """
        return self.getValue(data)

     def getValue(self,data):
        """
    Returns the value of data at the Locator if data is a Data object 
    otherwise the object is returned.
    """
        if isinstance(data,escript.Data):
           if data.getFunctionSpace()==self.getFunctionSpace():
             out=data.convertToNumArrayFromDPNo(self.getId()[0],self.getId()[1])
           else:
             out=data.interpolate(self.getFunctionSpace()).convertToNumArrayFromDPNo(self.getId()[0],self.getId()[1])
           if data.getRank()==0:
              return out[0]
           else:
              return out
        else:
           return data

# vim: expandtab shiftwidth=4:
1	# $Id$
2
3	"""
4	Provides some tools related to PDEs.
5
6	Currently includes:
7	- Projector - to project a discontinuous
8	"""
9
10	import escript
11	import linearPDEs
12	import numarray
13	import util
14
15	class Projector:
16	"""
17	The Projector is a factory which projects a discontiuous function onto a
18	continuous function on the a given domain.
19	"""
20	def __init__(self, domain, reduce = True, fast=True):
21	"""
22	Create a continuous function space projector for a domain.
23
24	@param domain: Domain of the projection.
25	@param reduce: Flag to reduce projection order (default is True)
26	@param fast: Flag to use a fast method based on matrix lumping (default is true)
27	"""
28	self.__pde = linearPDEs.LinearPDE(domain)
29	if fast:
30	self.__pde.setSolverMethod(linearPDEs.LinearPDE.LUMPING)
31	self.__pde.setSymmetryOn()
32	self.__pde.setReducedOrderTo(reduce)
33	self.__pde.setValue(D = 1.)
34	return
35
36	def __del__(self):
37	return
38
39	def __call__(self, input_data):
40	"""
41	Projects input_data onto a continuous function
42
43	@param input_data: The input_data to be projected.
44	"""
45	out=escript.Data(0.,input_data.getShape(),what=escript.ContinuousFunction(self.__pde.getDomain()))
46	if input_data.getRank()==0:
47	self.__pde.setValue(Y = input_data)
48	out=self.__pde.getSolution()
49	elif input_data.getRank()==1:
50	for i0 in range(input_data.getShape()[0]):
51	self.__pde.setValue(Y = input_data[i0])
52	out[i0]=self.__pde.getSolution()
53	elif input_data.getRank()==2:
54	for i0 in range(input_data.getShape()[0]):
55	for i1 in range(input_data.getShape()[1]):
56	self.__pde.setValue(Y = input_data[i0,i1])
57	out[i0,i1]=self.__pde.getSolution()
58	elif input_data.getRank()==3:
59	for i0 in range(input_data.getShape()[0]):
60	for i1 in range(input_data.getShape()[1]):
61	for i2 in range(input_data.getShape()[2]):
62	self.__pde.setValue(Y = input_data[i0,i1,i2])
63	out[i0,i1,i2]=self.__pde.getSolution()
64	else:
65	for i0 in range(input_data.getShape()[0]):
66	for i1 in range(input_data.getShape()[1]):
67	for i2 in range(input_data.getShape()[2]):
68	for i3 in range(input_data.getShape()[3]):
69	self.__pde.setValue(Y = input_data[i0,i1,i2,i3])
70	out[i0,i1,i2,i3]=self.__pde.getSolution()
71	return out
72
73
74	class Locator:
75	"""
76	Locator provides access to the values of data objects at a given
77	spatial coordinate x.
78
79	In fact, a Locator object finds the sample in the set of samples of a
80	given function space or domain where which is closest to the given
81	point x.
82	"""
83
84	def __init__(self,where,x=numarray.zeros((3,))):
85	"""
86	Initializes a Locator to access values in Data objects on the Doamin
87	or FunctionSpace where for the sample point which
88	closest to the given point x.
89	"""
90	if isinstance(where,escript.FunctionSpace):
91	self.__function_space=where
92	else:
93	self.__function_space=escript.ContinuousFunction(where)
94	self.__id=util.length(x[:self.__function_space.getDim()]-self.__function_space.getX()).mindp()
95
96	def __str__(self):
97	"""
98	Returns the coordinates of the Locator as a string.
99	"""
100	return "<Locator %s>"%str(self.getX())
101
102	def getFunctionSpace(self):
103	"""
104	Returns the function space of the Locator.
105	"""
106	return self.__function_space
107
108	def getId(self):
109	"""
110	Returns the identifier of the location.
111	"""
112	return self.__id
113
114	def getX(self):
115	"""
116	Returns the exact coordinates of the Locator.
117	"""
118	return self(self.getFunctionSpace().getX())
119
120	def __call__(self,data):
121	"""
122	Returns the value of data at the Locator of a Data object otherwise
123	the object is returned.
124	"""
125	return self.getValue(data)
126
127	def getValue(self,data):
128	"""
129	Returns the value of data at the Locator if data is a Data object
130	otherwise the object is returned.
131	"""
132	if isinstance(data,escript.Data):
133	if data.getFunctionSpace()==self.getFunctionSpace():
134	out=data.convertToNumArrayFromDPNo(self.getId()[0],self.getId()[1])
135	else:
136	out=data.interpolate(self.getFunctionSpace()).convertToNumArrayFromDPNo(self.getId()[0],self.getId()[1])
137	if data.getRank()==0:
138	return out[0]
139	else:
140	return out
141	else:
142	return data
143
144	# vim: expandtab shiftwidth=4:
Name	Value
svn:eol-style	native
svn:keywords	Author Date Id Revision