/[escript]/trunk/modellib/test/python/run_convection.py
ViewVC logotype

Contents of /trunk/modellib/test/python/run_convection.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 6651 - (show annotations)
Wed Feb 7 02:12:08 2018 UTC (2 years, 5 months ago) by jfenwick
File MIME type: text/x-python
File size: 5466 byte(s)
Make everyone sad by touching all the files

Copyright dates update

1 ##############################################################################
2 #
3 # Copyright (c) 2003-2018 by The University of Queensland
4 # http://www.uq.edu.au
5 #
6 # Primary Business: Queensland, Australia
7 # Licensed under the Apache License, version 2.0
8 # http://www.apache.org/licenses/LICENSE-2.0
9 #
10 # Development until 2012 by Earth Systems Science Computational Center (ESSCC)
11 # Development 2012-2013 by School of Earth Sciences
12 # Development from 2014 by Centre for Geoscience Computing (GeoComp)
13 #
14 ##############################################################################
15 from __future__ import print_function, division
16
17 __copyright__="""Copyright (c) 2003-2018 by The University of Queensland
18 http://www.uq.edu.au
19 Primary Business: Queensland, Australia"""
20 __license__="""Licensed under the Apache License, version 2.0
21 http://www.apache.org/licenses/LICENSE-2.0"""
22 __url__="https://launchpad.net/escript-finley"
23
24 #
25 #
26 # a very simple convection model in model frame:
27 #
28
29 import os
30 import sys
31 if sys.version_info >= (3,0):
32 from io import StringIO
33 else:
34 #specifically to avoid non-unicode default strings
35 #being passed to a python3 style StringIO that expects unicode
36 from StringIO import StringIO
37
38 import esys.escriptcore.utestselect as unittest
39 from esys.escriptcore.testing import *
40 from esys.escript import hasFeature
41 from esys.escript.modelframe import Link,Simulation
42 from esys.modellib.input import Sequencer,InterpolateOverBox,GaussianProfile,LinearCombination
43 from esys.modellib.flow import SteadyIncompressibleFlow
44 from esys.modellib.temperature import TemperatureAdvection
45 from esys.modellib.materials import SimpleEarthModel,GravityForce
46 from esys.modellib.visualization import WriteVTK
47
48 try:
49 import esys.dudley
50 HAVE_DUDLEY = True
51 except ImportError:
52 HAVE_DUDLEY = False
53
54 try:
55 import esys.finley
56 from esys.modellib.geometry import RectangularDomain, ScalarConstrainerOverBox,VectorConstrainerOverBox
57 HAVE_FINLEY = True
58 except ImportError:
59 HAVE_FINLEY = False
60
61 # TODO: once Amesos2 can deal with block matrices uncomment
62 HAVE_DIRECT = hasFeature("PASO_DIRECT") #or hasFeature('trilinos')
63
64 try:
65 WORKDIR=os.environ['MODELLIB_WORKDIR']
66 except KeyError as e:
67 WORKDIR='.'
68
69 #if this func is inside the test case, the Link() calls blow up
70 def run(dom, stream):
71 temp_val=InterpolateOverBox()
72 temp_val.domain=Link(dom,"domain")
73 temp_val.value_left_bottom_front=1.
74 temp_val.value_right_bottom_front=1.
75 temp_val.value_left_top_front=0.
76 temp_val.value_right_top_front=0.
77 temp_val.value_left_bottom_back=1.
78 temp_val.value_right_bottom_back=1.
79 temp_val.value_left_top_back=0.
80 temp_val.value_right_top_back=0.
81
82 temp_constraints=ScalarConstrainerOverBox()
83 temp_constraints.domain=Link(dom)
84 temp_constraints.top=1
85 temp_constraints.bottom=1
86
87 vel_constraints=VectorConstrainerOverBox()
88 vel_constraints.domain=Link(dom)
89 vel_constraints.left=[1,0,0]
90 vel_constraints.right=[1,0,0]
91 vel_constraints.top=[0,1,0]
92 vel_constraints.bottom=[0,1,0]
93 vel_constraints.front=[0,0,1]
94 vel_constraints.back=[0,0,1]
95
96 mat=SimpleEarthModel()
97 mat.density0=1.
98 mat.viscocity0=1.
99 mat.rayleigh_number=10000.
100 mat.alpha=0.001
101
102 temp=TemperatureAdvection(debug=True)
103 temp.domain=Link(dom)
104 temp.density=Link(mat,"density0")
105 temp.heat_capacity=Link(mat,"heat_capacity")
106 temp.location_fixed_temperature=Link(temp_constraints,"location_of_constraint")
107 temp.fixed_temperature=Link(temp_val,"out")
108 temp.safety_factor=0.01
109 mat.temperature=Link(temp,"temperature")
110 grav=GravityForce()
111 grav.domain=Link(dom,"domain")
112 grav.direction=[0.,-1.,0.]
113 grav.density=Link(mat,"density")
114 grav.gravity=Link(mat,"gravity")
115
116
117 vel=SteadyIncompressibleFlow(debug=True)
118 vel.domain=Link(dom)
119 vel.internal_force=Link(grav,"gravity_force")
120 vel.viscosity=Link(mat,"viscosity")
121 vel.location_prescribed_velocity=Link(vel_constraints,"location_of_constraint")
122 vel.rel_tol=1.e-6
123 temp.velocity=Link(vel,"velocity")
124
125 sq=Sequencer()
126 sq.t_end=0.001
127
128 vis=WriteVTK()
129 vis.t=Link(sq)
130 vis.data0=Link(temp,"temperature")
131 vis.data1=Link(vel,"velocity")
132 vis.dt=0.0001
133 vis.filename=os.path.join(WORKDIR,"temp.vtu")
134
135 per=GaussianProfile()
136 per.domain=Link(dom)
137 per.x_c=[0.5,0.5,0.5]
138 per.A=0.0001
139 per.width=0.01
140 per.r=0
141
142 lc=LinearCombination()
143 lc.f0=1.
144 lc.v0=Link(per,"out")
145 lc.f1=1.
146 lc.v1=Link(temp_val,"out")
147 temp.temperature=Link(lc,"out")
148
149 s=Simulation([sq,vel_constraints, temp_constraints,Simulation([vel],debug=True),temp,vis],debug=True)
150 s.writeXML(stream)
151 s.run()
152
153 @unittest.skipUnless(HAVE_DIRECT, "Direct solver not available")
154 class Test_Convection(unittest.TestCase):
155 def setUp(self):
156 import sys
157 self.old = sys.stdout
158 sys.stdout = StringIO()
159
160 def tearDown(self):
161 import sys
162 sys.stdout = self.old
163
164 @unittest.skipUnless(HAVE_FINLEY, "Finley module not available")
165 def test_order2(self):
166 dom=RectangularDomain()
167 dom.order=2
168 run(dom, sys.stdout)
169
170 @unittest.skipUnless(HAVE_DUDLEY and HAVE_FINLEY, "Dudley module not available")
171 def test_order1(self):
172 dom=RectangularDomain(esys.dudley)
173 dom.order=1
174 run(dom, sys.stdout)
175
176 if __name__ == '__main__':
177 run_tests(__name__, exit_on_failure=True)

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26