Diff of /trunk/escript/py_src/heat.py

revision 2158 by caltinay, Mon Dec 15 07:17:47 2008 UTC revision 2169 by caltinay, Wed Dec 17 03:08:58 2008 UTC
# Line 37  import util Line 37  import util
37  from linearPDEs import TransportPDE  from linearPDEs import TransportPDE
38
39  class TemperatureCartesian(TransportPDE):  class TemperatureCartesian(TransportPDE):
40        """      """
41        solves      Represents and solves the temperature advection-diffusion problem
42
43            rhocp(T_{,t} + v_i T_{,i} - ( k T_{,i})_i = Q      M{rhocp(T_{,t} + v_i T_{,i} - ( k T_{,i})_i = Q}
44
45                   k T_{,i}*n_i=surface_flux      M{k T_{,i}*n_i=surface_flux} and M{T_{,t} = 0} where C{given_T_mask}>0.
46
47                   T_{,t} = 0  where given_T_mask>0      If surface_flux is not given 0 is assumed.
48
49        if surface_flux is not given 0 is assumed.      Typical usage::
50
51        typical usage:          sp = TemperatureCartesian(domain)
52            sp.setTolerance(1.e-4)
53              sp=TemperatureCartesian(domain)          t = 0
54              sp.setTolerance(1.e-4)          T = ...
55              t=0          sp.setValues(rhocp=...,  v=..., k=..., given_T_mask=...)
56              T=...          sp.setInitialTemperature(T)
57              sp.setValues(rhocp = ..,  v=.., k=.., given_T_mask=..)          while t < t_end:
58              sp.setInitialTemperature(T)              sp.setValue(Q=...)
59              while t < t_end:              T = sp.getTemperature(dt)
60                  sp.setValue(Q= ...)              t += dt
61                  T=sp.getTemperature(dt)      """
62                  t+=dt      def __init__(self,domain,theta=0.5,**kwargs):
"""
def __init__(self,domain,theta=0.5,**kwargs):
63          """          """
65
66          @param domain: domain of the problem          @param domain: domain of the problem
67          @param theta: method control: theta=1. backward Euler          @param theta: method control:
68                                        theta=0.5: Crank-Nicholson scheme              - theta=1.0: backward Euler
69                                        theta=0.: forward Euler (not recommendable)              - theta=0.5: Crank-Nicholson scheme
70                - theta=0.0: forward Euler (not recommended)
71          """          """
72          TransportPDE.__init__(self,domain,numEquations=1,theta=theta,**kwargs)          TransportPDE.__init__(self,domain,numEquations=1,theta=theta,**kwargs)
73          self.setReducedOrderOn()          self.setReducedOrderOn()
74          self.__rhocp=None          self.__rhocp=None
75          self.__v=None          self.__v=None
76
77        def setInitialTemperature(self,T):      def setInitialTemperature(self,T):
78            """
79            Same as L{setInitialSolution}.
80            """
81          self.setInitialSolution(T)          self.setInitialSolution(T)
82
84             if rhocp!=None:          if rhocp!=None:
85                 self.__rhocp=rhocp              self.__rhocp=rhocp
86             if v!=None:          if v!=None:
87                 self.__v=v              self.__v=v
88             if rhocp!=None:          if rhocp!=None:
89                 super(TemperatureCartesian,self).setValue(M=self.__rhocp)              super(TemperatureCartesian,self).setValue(M=self.__rhocp)
90             if (rhocp!=None or v!=None) and self.__rhocp!=None and self.__v!=None:          if (rhocp!=None or v!=None) and self.__rhocp!=None and self.__v!=None:
91                 super(TemperatureCartesian,self).setValue(C=-self.__rhocp*self.__v)              super(TemperatureCartesian,self).setValue(C=-self.__rhocp*self.__v)
92             if k!=None:          if k!=None:
93                 super(TemperatureCartesian,self).setValue(A=-k*util.kronecker(self.getDomain()))              super(TemperatureCartesian,self).setValue(A=-k*util.kronecker(self.getDomain()))
94             if Q!=None:          if Q!=None:
95                 super(TemperatureCartesian,self).setValue(Y=Q)              super(TemperatureCartesian,self).setValue(Y=Q)
96             if surface_flux!=None:          if surface_flux!=None:
97                 super(TemperatureCartesian,self).setValue(y=surface_flux)              super(TemperatureCartesian,self).setValue(y=surface_flux)