test/python/drucker_prager.py

# $Id$

__copyright__="""  Copyright (c) 2006 by ACcESS MNRF
                    http://www.access.edu.au
                Primary Business: Queensland, Australia"""
__license__="""Licensed under the Open Software License version 3.0
             http://www.opensource.org/licenses/osl-3.0.php"""

# run by scons build/posix/modellib/test/python/drucker_prager.passed from check out cd.

import os
from esys.modellib.mechanics import DruckerPrager
from esys.escript.modelframe import Link,Simulation
from esys.modellib.geometry import RectangularDomain, VectorConstrainer, UpdateGeometry
from esys.modellib.input import Sequencer, InterpolateOverBox
from esys.modellib.visualization import WriteVTK

try:
   WORKDIR=os.environ['MODELLIB_WORKDIR']
except KeyError:
   WORKDIR='.'


debug=True

dom=RectangularDomain(debug)
dom.l=[1.,1.,1.]
dom.n=[50,30,2]
dom.order=2
dom.integrationOrder=2


sq=Sequencer(debug)
sq.t=0
sq.t_end=0.8
sq.dt_max=0.02

iob=InterpolateOverBox(debug)
iob.domain=Link(dom,"domain")
iob.value_left_bottom_front=[1.,0.,0.]
iob.value_right_bottom_front=[0.,0.,0.]
iob.value_left_bottom_back=[1.,0.,0.]
iob.value_right_bottom_back=[0.,0.,0.]
iob.value_left_top_front=[1.,0.,0.]
iob.value_right_top_front=[0.,0.,0.]
iob.value_left_top_back=[1.,0.,0.]
iob.value_right_top_back=[0.,0.,0.]

iob.value_left_bottom_front=[1.,0.]
iob.value_right_bottom_front=[0.,0.]
iob.value_left_bottom_back=[1.,0.]
iob.value_right_bottom_back=[0.,0.]
iob.value_left_top_front=[1.,0.]
iob.value_right_top_front=[0.,0.]
iob.value_left_top_back=[1.,0.]
iob.value_right_top_back=[0.,0.]

m=DruckerPrager(debug)
m.domain=Link(dom,"domain")

cv=VectorConstrainer(debug)
cv.domain=Link(dom,"domain")
cv.value=Link(iob,"out")

cv.left=[True, False, False] 
cv.right= [True, False, False]
cv.bottom= [False, True, False]
cv.top= [False, False, False]
cv.front= [False, False, True]
cv.back= [False, False, False]

m.velocity=Link(iob,"out")
m.prescribed_velocity=Link(cv,"value_of_constraint")
m.location_prescribed_velocity=Link(cv,"location_of_constraint")
m.rel_tol=0.0001

m.expansion_coefficient= 0.
m.bulk_modulus=1000.
m.shear_modulus=1.
m.plastic_stress=0.
m.friction_parameter=0.
m.dilatancy_parameter=0.
m.shear_length=m.shear_modulus*0.75*10000.


ug=UpdateGeometry(debug)
ug.domain=Link(dom,"domain")
ug.displacement=Link(m,"displacement")

vis=WriteVTK()
vis.t=Link(sq)
vis.scalar=Link(m,"plastic_stress")
vis.vector=Link(m,"velocity")
vis.tensor=Link(m,"stress")
vis.dt=0.5
vis.filename=WORKDIR+"/temp.xml"

s=Simulation([sq,m,vis],debug=True)
s.run()
1	gross	821	# $Id$
2
3			__copyright__=""" Copyright (c) 2006 by ACcESS MNRF
4			http://www.access.edu.au
5			Primary Business: Queensland, Australia"""
6			__license__="""Licensed under the Open Software License version 3.0
7			http://www.opensource.org/licenses/osl-3.0.php"""
8
9			# run by scons build/posix/modellib/test/python/drucker_prager.passed from check out cd.
10
11	gross	829	import os
12	gross	821	from esys.modellib.mechanics import DruckerPrager
13			from esys.escript.modelframe import Link,Simulation
14	gross	823	from esys.modellib.geometry import RectangularDomain, VectorConstrainer, UpdateGeometry
15	gross	821	from esys.modellib.input import Sequencer, InterpolateOverBox
16	gross	829	from esys.modellib.visualization import WriteVTK
17	gross	821
18	gross	829	try:
19			WORKDIR=os.environ['MODELLIB_WORKDIR']
20			except KeyError:
21			WORKDIR='.'
22
23
24	gross	821	debug=True
25
26			dom=RectangularDomain(debug)
27	gross	834	dom.l=[1.,1.,1.]
28	gross	838	dom.n=[50,30,2]
29	gross	834	dom.order=2
30	gross	838	dom.integrationOrder=2
31	gross	821
32
33			sq=Sequencer(debug)
34			sq.t=0
35	gross	834	sq.t_end=0.8
36	gross	838	sq.dt_max=0.02
37	gross	821
38			iob=InterpolateOverBox(debug)
39			iob.domain=Link(dom,"domain")
40			iob.value_left_bottom_front=[1.,0.,0.]
41			iob.value_right_bottom_front=[0.,0.,0.]
42			iob.value_left_bottom_back=[1.,0.,0.]
43			iob.value_right_bottom_back=[0.,0.,0.]
44			iob.value_left_top_front=[1.,0.,0.]
45			iob.value_right_top_front=[0.,0.,0.]
46			iob.value_left_top_back=[1.,0.,0.]
47			iob.value_right_top_back=[0.,0.,0.]
48
49	gross	829	iob.value_left_bottom_front=[1.,0.]
50			iob.value_right_bottom_front=[0.,0.]
51			iob.value_left_bottom_back=[1.,0.]
52			iob.value_right_bottom_back=[0.,0.]
53			iob.value_left_top_front=[1.,0.]
54			iob.value_right_top_front=[0.,0.]
55			iob.value_left_top_back=[1.,0.]
56			iob.value_right_top_back=[0.,0.]
57
58	gross	821	m=DruckerPrager(debug)
59			m.domain=Link(dom,"domain")
60
61			cv=VectorConstrainer(debug)
62			cv.domain=Link(dom,"domain")
63			cv.value=Link(iob,"out")
64	gross	825
65	gross	821	cv.left=[True, False, False]
66	gross	825	cv.right= [True, False, False]
67			cv.bottom= [False, True, False]
68	gross	821	cv.top= [False, False, False]
69	gross	825	cv.front= [False, False, True]
70	gross	821	cv.back= [False, False, False]
71
72	gross	829	m.velocity=Link(iob,"out")
73	gross	821	m.prescribed_velocity=Link(cv,"value_of_constraint")
74			m.location_prescribed_velocity=Link(cv,"location_of_constraint")
75	gross	838	m.rel_tol=0.0001
76	gross	821
77	gross	825	m.expansion_coefficient= 0.
78	gross	834	m.bulk_modulus=1000.
79	gross	825	m.shear_modulus=1.
80			m.plastic_stress=0.
81			m.friction_parameter=0.
82			m.dilatancy_parameter=0.
83	gross	838	m.shear_length=m.shear_modulus0.7510000.
84	gross	821
85	gross	825
86	gross	823	ug=UpdateGeometry(debug)
87			ug.domain=Link(dom,"domain")
88			ug.displacement=Link(m,"displacement")
89
90	gross	829	vis=WriteVTK()
91			vis.t=Link(sq)
92	gross	832	vis.scalar=Link(m,"plastic_stress")
93	gross	829	vis.vector=Link(m,"velocity")
94			vis.tensor=Link(m,"stress")
95	gross	836	vis.dt=0.5
96	gross	829	vis.filename=WORKDIR+"/temp.xml"
97
98			s=Simulation([sq,m,vis],debug=True)
99	gross	821	s.run()