modellib/py_src/temperature.py

# $Id$


from esys.escript import *
from esys.modelframe import Model,IterationDivergenceError
from esys.linearPDEs import AdvectivePDE,LinearPDE

class TemperatureDiffusion(Model):
       """ """

       def __init__(self,debug=False):
           Model.__init__(self,debug=debug)
           self.declareParameter(domain=None, \
                                 tend=0., \
                                 dt=0.1, \
                                 temperature=1., \
                                 density=1., \
                                 c_p=1., \
                                 thermal_permabilty=1., \
                                 reference_temperature=0., \
                                 radiation_coefficient=0., \
                                 thermal_source=0., \
                                 location_fixed_temperature=Data(), \
                                 iterate=True, \
                                 tol=1.e-8, \
                                 implicit=True)
           self.iter=0

       def doInitialization(self,t):
           self.tn=t
           self.pde=LinearPDE(self.domain)
           self.pde.setSymmetryOn()

       def doIterationInitialization(self,dt):
            self.iter=0
            self.T_last=self.temperature
            self.diff=1.e400

       def doIterationStep(self,dt):
          T=self.temperature
          diff=self.diff
          dim=self.pde.getDim()
          self.iter+=1
          rhocp=self.density*self.c_p
          self.pde.setValue(A=self.thermal_permabilty*kronecker(dim), \
                              D=rhocp/dt, \
                              Y=self.thermal_source+rhocp/dt*self.T_last, \
                              d=self.radiation_coefficient, \
                              y=self.radiation_coefficient*self.reference_temperature, \
                              q=self.location_fixed_temperature, \
                              r=self.T_last)
          if isinstance(self,TemperatureAdvection): self.pde.setValue(C=self.velocity[:dim]*rhocp)
          self.pde.setTolerance(self.tol*1.e-2)
          self.temperature=self.pde.getSolution()
          self.diff=Lsup(T-self.temperature)
          if diff<=self.diff:
              raise IterationDivergenceError,"no improvement in the temperature iteration"

       def terminate(self):
          if not self.implicit:
              return True
          elif self.iter<1:
              return False
          else:
             return self.diff<self.tol*Lsup(self.temperature)

       def doIterationFinalization(self,dt):
          self.tn+=dt

       def getSafeTimeStepSize(self,dt):
           return self.dt

       def finalize(self):
            return self.tn>=self.tend

class TemperatureAdvection(Model):
       """ """

       def __init__(self,debug=False):
           Model.__init__(self,debug=debug)
           self.declareParameter(velocity=numarray.zeros([3]))

       def doInitialization(self,t):
           self.tn=t
           self.pde=AdvectivePDE(self.domain)

       def getSafeTimeStepSize(self,dt):
           v=Lsup(self.velocity)
           if v>0.:
               return min(self.dt,Lsup(self.pde.getDomain().getSize())/v)
           else:
               return self.dt

       
if __name__=="__main__":
   from esys.modelframe import Link,Simulation,ExplicitSimulation
   from esys.visualization import WriteVTK
   from esys.materials import MaterialTable
   from esys.geometry import RectangularDomain,ScalarConstrainer
   from esys.input import InterpolatedTimeProfile,GausseanProfile

   dom=RectangularDomain()
   constraints=ScalarConstrainer()
   constraints.domain=Link(dom)
   constraints.top=1
   constraints.bottom=1
   
   mt=MaterialTable()
   
   pf=InterpolatedTimeProfile()
   pf.t=[0.,0.25,0.5,0.75]
   pf.values=[0.,1.,1.,0.]
   
   q=GausseanProfile()
   q.domain=Link(dom)
   q.width=0.05
   q.x_c=numarray.array([0.5,0.5,0.5])
   q.r=0.01
   q.A=Link(pf,"out")
   
   tt=TemperatureDiffusion()
   tt.domain=Link(dom)
   tt.tend=1.
   tt.dt=0.1
   tt.temperature=0.
   tt.density=Link(mt)
   tt.c_p=Link(mt)
   tt.thermal_permabilty=Link(mt)
   tt.reference_temperature=0.
   tt.radiation_coefficient=Link(mt)
   tt.thermal_source=Link(q,"out")
   tt.location_fixed_temperature=Link(constraints,"location_of_constraint")
   tt.implicit=True
   
   vis=WriteVTK()
   vis.scalar=Link(tt,"temperature")
   
   s=ExplicitSimulation([dom,constraints,pf,q,Simulation([mt,tt],debug=True),vis],debug=True)
   # s=Simulation([dom,constraints,pf,q,Simulation([mt,tt]),vis])
   s.writeXML()
   s.run()
1	jgs	127	# $Id$
2
3
4			from esys.escript import *
5			from esys.modelframe import Model,IterationDivergenceError
6			from esys.linearPDEs import AdvectivePDE,LinearPDE
7
8			class TemperatureDiffusion(Model):
9			""" """
10
11			def __init__(self,debug=False):
12			Model.__init__(self,debug=debug)
13			self.declareParameter(domain=None, \
14			tend=0., \
15			dt=0.1, \
16			temperature=1., \
17			density=1., \
18			c_p=1., \
19			thermal_permabilty=1., \
20			reference_temperature=0., \
21			radiation_coefficient=0., \
22			thermal_source=0., \
23			location_fixed_temperature=Data(), \
24			iterate=True, \
25			tol=1.e-8, \
26			implicit=True)
27			self.iter=0
28
29			def doInitialization(self,t):
30			self.tn=t
31			self.pde=LinearPDE(self.domain)
32			self.pde.setSymmetryOn()
33
34			def doIterationInitialization(self,dt):
35			self.iter=0
36			self.T_last=self.temperature
37			self.diff=1.e400
38
39			def doIterationStep(self,dt):
40			T=self.temperature
41			diff=self.diff
42			dim=self.pde.getDim()
43			self.iter+=1
44			rhocp=self.density*self.c_p
45			self.pde.setValue(A=self.thermal_permabilty*kronecker(dim), \
46			D=rhocp/dt, \
47			Y=self.thermal_source+rhocp/dt*self.T_last, \
48			d=self.radiation_coefficient, \
49			y=self.radiation_coefficient*self.reference_temperature, \
50			q=self.location_fixed_temperature, \
51			r=self.T_last)
52			if isinstance(self,TemperatureAdvection): self.pde.setValue(C=self.velocity[:dim]*rhocp)
53			self.pde.setTolerance(self.tol*1.e-2)
54			self.temperature=self.pde.getSolution()
55			self.diff=Lsup(T-self.temperature)
56			if diff<=self.diff:
57			raise IterationDivergenceError,"no improvement in the temperature iteration"
58
59			def terminate(self):
60			if not self.implicit:
61			return True
62			elif self.iter<1:
63			return False
64			else:
65			return self.diff<self.tol*Lsup(self.temperature)
66
67			def doIterationFinalization(self,dt):
68			self.tn+=dt
69
70			def getSafeTimeStepSize(self,dt):
71			return self.dt
72
73			def finalize(self):
74			return self.tn>=self.tend
75
76			class TemperatureAdvection(Model):
77			""" """
78
79			def __init__(self,debug=False):
80			Model.__init__(self,debug=debug)
81			self.declareParameter(velocity=numarray.zeros([3]))
82
83			def doInitialization(self,t):
84			self.tn=t
85			self.pde=AdvectivePDE(self.domain)
86
87			def getSafeTimeStepSize(self,dt):
88			v=Lsup(self.velocity)
89			if v>0.:
90			return min(self.dt,Lsup(self.pde.getDomain().getSize())/v)
91			else:
92			return self.dt
93
94
95
96			if __name__=="__main__":
97			from esys.modelframe import Link,Simulation,ExplicitSimulation
98			from esys.visualization import WriteVTK
99			from esys.materials import MaterialTable
100			from esys.geometry import RectangularDomain,ScalarConstrainer
101			from esys.input import InterpolatedTimeProfile,GausseanProfile
102
103			dom=RectangularDomain()
104			constraints=ScalarConstrainer()
105			constraints.domain=Link(dom)
106			constraints.top=1
107			constraints.bottom=1
108
109			mt=MaterialTable()
110
111			pf=InterpolatedTimeProfile()
112			pf.t=[0.,0.25,0.5,0.75]
113			pf.values=[0.,1.,1.,0.]
114
115			q=GausseanProfile()
116			q.domain=Link(dom)
117			q.width=0.05
118			q.x_c=numarray.array([0.5,0.5,0.5])
119			q.r=0.01
120			q.A=Link(pf,"out")
121
122			tt=TemperatureDiffusion()
123			tt.domain=Link(dom)
124			tt.tend=1.
125			tt.dt=0.1
126			tt.temperature=0.
127			tt.density=Link(mt)
128			tt.c_p=Link(mt)
129			tt.thermal_permabilty=Link(mt)
130			tt.reference_temperature=0.
131			tt.radiation_coefficient=Link(mt)
132			tt.thermal_source=Link(q,"out")
133			tt.location_fixed_temperature=Link(constraints,"location_of_constraint")
134			tt.implicit=True
135
136			vis=WriteVTK()
137			vis.scalar=Link(tt,"temperature")
138
139			s=ExplicitSimulation([dom,constraints,pf,q,Simulation([mt,tt],debug=True),vis],debug=True)
140			# s=Simulation([dom,constraints,pf,q,Simulation([mt,tt]),vis])
141			s.writeXML()
142			s.run()
Name	Value
svn:eol-style	native
svn:keywords	Author Date Id Revision