test/python/convection.py

# $Id$

#
#
#    a very simple convection model in model frame:
#

from esys.escript.modelframe import Link,Simulation
from esys.modellib.geometry import RectangularDomain,ScalarConstrainer,VectorConstrainer
from esys.modellib.input import Sequencer,InterpolateOverBox,GaussianProfile,LinearCombination
from esys.modellib.flow import SteadyIncompressibleFlow
from esys.modellib.temperature import TemperatureAdvection
from esys.modellib.materials import SimpleEarthModel,GravityForce
from esys.modellib.visualization import WriteVTK

dom=RectangularDomain()
dom.order=2

temp_constraints=ScalarConstrainer()
temp_constraints.domain=Link(dom)
temp_constraints.top=1
temp_constraints.bottom=1

vel_constraints=VectorConstrainer()
vel_constraints.domain=Link(dom)
vel_constraints.left=[1,0,0]
vel_constraints.right=[1,0,0]
vel_constraints.top=[0,1,0]
vel_constraints.bottom=[0,1,0]
vel_constraints.front=[0,0,1]
vel_constraints.back=[0,0,1]


temp_val=InterpolateOverBox()
temp_val.domain=Link(dom,"domain")
temp_val.right_top_back=Link(dom,"l")
temp_val.value_left_bottom_front=1.
temp_val.value_right_bottom_front=1.
temp_val.value_left_top_front=0.
temp_val.value_right_top_front=0.
temp_val.value_left_bottom_back=1.
temp_val.value_right_bottom_back=1.
temp_val.value_left_top_back=0.
temp_val.value_right_top_back=0.

mat=SimpleEarthModel()
mat.density0=1.
mat.viscocity0=1.
mat.rayleigh_number=10000.
mat.alpha=0.001

temp=TemperatureAdvection(debug=True)
temp.domain=Link(dom)
temp.density=Link(mat,"density0")
temp.heat_capacity=Link(mat,"heat_capacity")
temp.location_fixed_temperature=Link(temp_constraints,"location_of_constraint")
temp.fixed_temperature=Link(temp_val,"out")
temp.safety_factor=0.01
mat.temperature=Link(temp,"temperature")


grav=GravityForce()
grav.domain=Link(dom,"domain")
grav.direction=[0.,-1.,0.]
grav.density=Link(mat,"density")
grav.gravity=Link(mat,"gravity")


vel=SteadyIncompressibleFlow(debug=True)
vel.domain=Link(dom)
vel.internal_force=Link(grav,"gravity_force")
vel.viscosity=Link(mat,"viscosity")
vel.location_prescribed_velocity=Link(vel_constraints,"location_of_constraint")
vel.rel_tol=1.e-6
temp.velocity=Link(vel,"velocity")

sq=Sequencer()
sq.t_end=0.005

vis=WriteVTK()
vis.t=Link(sq)
vis.scalar=Link(temp,"temperature")
vis.vector=Link(vel,"velocity")
vis.stride=5

per=GaussianProfile()
per.domain=Link(dom)
per.x_c=[0.5,0.5,0.5]
per.A=0.0001
per.width=0.01
per.r=0

lc=LinearCombination()
lc.f0=1.
lc.v0=Link(per,"out")
lc.f1=1.
lc.v1=Link(temp_val,"out")
temp.temperature=Link(lc,"out")

# s=Simulation([dom,sq,temp_constraints,vel_constraints,temp_val,per,lc,mat,grav,Simulation([vel],debug=True),temp,vis],debug=True)
s=Simulation([dom,sq,Simulation([vel],debug=True),temp,vis],debug=True)
s.writeXML()
s.run()
1	jgs	147	# $Id$
2
3			#
4			#
5			# a very simple convection model in model frame:
6			#
7
8	jgs	149	from esys.escript.modelframe import Link,Simulation
9			from esys.modellib.geometry import RectangularDomain,ScalarConstrainer,VectorConstrainer
10			from esys.modellib.input import Sequencer,InterpolateOverBox,GaussianProfile,LinearCombination
11			from esys.modellib.flow import SteadyIncompressibleFlow
12			from esys.modellib.temperature import TemperatureAdvection
13			from esys.modellib.materials import SimpleEarthModel,GravityForce
14			from esys.modellib.visualization import WriteVTK
15	jgs	147
16			dom=RectangularDomain()
17			dom.order=2
18
19			temp_constraints=ScalarConstrainer()
20			temp_constraints.domain=Link(dom)
21			temp_constraints.top=1
22			temp_constraints.bottom=1
23
24			vel_constraints=VectorConstrainer()
25			vel_constraints.domain=Link(dom)
26			vel_constraints.left=[1,0,0]
27			vel_constraints.right=[1,0,0]
28			vel_constraints.top=[0,1,0]
29			vel_constraints.bottom=[0,1,0]
30			vel_constraints.front=[0,0,1]
31			vel_constraints.back=[0,0,1]
32
33
34			temp_val=InterpolateOverBox()
35			temp_val.domain=Link(dom,"domain")
36			temp_val.right_top_back=Link(dom,"l")
37			temp_val.value_left_bottom_front=1.
38			temp_val.value_right_bottom_front=1.
39			temp_val.value_left_top_front=0.
40			temp_val.value_right_top_front=0.
41			temp_val.value_left_bottom_back=1.
42			temp_val.value_right_bottom_back=1.
43			temp_val.value_left_top_back=0.
44			temp_val.value_right_top_back=0.
45
46			mat=SimpleEarthModel()
47			mat.density0=1.
48			mat.viscocity0=1.
49			mat.rayleigh_number=10000.
50			mat.alpha=0.001
51
52			temp=TemperatureAdvection(debug=True)
53			temp.domain=Link(dom)
54			temp.density=Link(mat,"density0")
55			temp.heat_capacity=Link(mat,"heat_capacity")
56			temp.location_fixed_temperature=Link(temp_constraints,"location_of_constraint")
57			temp.fixed_temperature=Link(temp_val,"out")
58			temp.safety_factor=0.01
59			mat.temperature=Link(temp,"temperature")
60
61
62			grav=GravityForce()
63			grav.domain=Link(dom,"domain")
64			grav.direction=[0.,-1.,0.]
65			grav.density=Link(mat,"density")
66			grav.gravity=Link(mat,"gravity")
67
68
69			vel=SteadyIncompressibleFlow(debug=True)
70			vel.domain=Link(dom)
71			vel.internal_force=Link(grav,"gravity_force")
72			vel.viscosity=Link(mat,"viscosity")
73			vel.location_prescribed_velocity=Link(vel_constraints,"location_of_constraint")
74			vel.rel_tol=1.e-6
75			temp.velocity=Link(vel,"velocity")
76
77			sq=Sequencer()
78			sq.t_end=0.005
79
80			vis=WriteVTK()
81			vis.t=Link(sq)
82			vis.scalar=Link(temp,"temperature")
83			vis.vector=Link(vel,"velocity")
84			vis.stride=5
85
86	jgs	148	per=GaussianProfile()
87	jgs	147	per.domain=Link(dom)
88			per.x_c=[0.5,0.5,0.5]
89			per.A=0.0001
90			per.width=0.01
91			per.r=0
92
93			lc=LinearCombination()
94			lc.f0=1.
95			lc.v0=Link(per,"out")
96			lc.f1=1.
97			lc.v1=Link(temp_val,"out")
98			temp.temperature=Link(lc,"out")
99
100	jgs	148	# s=Simulation([dom,sq,temp_constraints,vel_constraints,temp_val,per,lc,mat,grav,Simulation([vel],debug=True),temp,vis],debug=True)
101			s=Simulation([dom,sq,Simulation([vel],debug=True),temp,vis],debug=True)
102	jgs	147	s.writeXML()
103			s.run()
Name	Value
svn:eol-style	native
svn:keywords	Author Date Id Revision