/[escript]/trunk/modellib/py_src/crustal/setups.py
ViewVC logotype

Annotation of /trunk/modellib/py_src/crustal/setups.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1044 - (hide annotations)
Mon Mar 19 07:29:31 2007 UTC (12 years, 2 months ago) by gross
File MIME type: text/x-python
File size: 9155 byte(s)
clear name tagging is supported now.
1 gross 987 """
2     set up of mining activities in modelframe
3    
4     @var __author__: name of author
5     @var __licence__: licence agreement
6     @var __url__: url entry point on documentation
7     @var __version__: version
8     @var __date__: date of the version
9     """
10    
11     __copyright__=""" Copyright (c) 2006 by ACcESS MNRF
12     http://www.access.edu.au
13     Primary Business: Queensland, Australia"""
14     __license__="""Licensed under the Open Software License version 3.0
15     http://www.opensource.org/licenses/osl-3.0.php"""
16     __author__="Lutz Gross, l.gross@uq.edu.au"
17     __url__="http://www.iservo.edu.au/esys/escript"
18     __version__="$Revision$"
19     __date__="$Date$"
20    
21 gross 1000 from esys.escript.modelframe import Model, ParameterSet
22 gross 997 from esys.escript import *
23     from esys.escript.linearPDEs import LinearPDE
24 gross 987 from mines import parse
25 gross 1000 import numarray
26 gross 987
27     class MiningHistory(Model):
28     """
29     manages the history of mines. It mandels the time steping according to the available
30     data and a dictionary of mass changes per year for all the mines.
31    
32     @ivar history: mine site history file.
33     @type history: C{DataSource}
34     @ivar mass_changes: current mass change per year.
35     @type mass_changes: C{DataSource}
36    
37     """
38     def __init__(self,**kwargs):
39     """
40     """
41     super(MiningHistory,self).__init__(**kwargs)
42     self.declareParameter(t=0.,
43     history=None,
44     mass_changes=None)
45    
46     def doInitialization(self):
47     """
48     initialize time integration
49     """
50     self.__minesite=parse(open(self.history.getLocalFileName(),'r'))
51     self.mass_changes=self.__minesite.getMassChanges(self.t)
52    
53     def doStepPreprocessing(self, dt):
54     self.mass_changes=self.__minesite.getMassChanges(self.t)
55    
56     def getSafeTimeStepSize(self, dt):
57     return self.__minesite.getNextTimeMarker(self.t)-self.t
58    
59    
60 gross 997 class DensityChange(Model):
61     """
62     translates local mass extraction into density changes.
63     "local" refers a taged region.
64 gross 987
65 gross 997 @ivar domain: mining region
66     @type domain: L{Domian}
67     @ivar mass_rate: rate of change mass in a tagged region
68     @type mass_rate: C{dict}
69     @ivar density_rate: density in each region
70     @type density_rate: C{Data}
71     """
72     def __init__(self,**kwargs):
73     """
74     """
75     super(DensityChange,self).__init__(**kwargs)
76     self.declareParameter(domain=None,
77     mass_rate={},
78     density_rate=None)
79    
80     def doInitialization(self):
81     """
82     initialize time integration
83     """
84     self.__volumes={}
85 gross 1044 for i in getTagNames(self.domain):
86     c=integrate(insertTaggedValues(Function(self.domain),**{i : 1.}))
87 gross 997 if c>0: self.__volumes[i] = c
88    
89     def doStepPreprocessing(self, dt):
90     d={}
91     for i in self.__volumes:
92 gross 998 if self.mass_rate.has_key(i): d[i]=-self.mass_rate[i]/self.__volumes[i]
93 gross 1044 self.density_rate=inserTaggedValues(Scalar(0.,Function(self.domain)),**d)
94 gross 997
95     class LinearElasticStressChange(Model):
96     """
97     calculates the stress according to an initial garvity field and a consecutive
98     chenge of density based on linar elastic material model. It is assumed that
99     the lame coefficients don't change over time, the direction of gravity is
100     pointing into the x_{dim} direction.
101    
102     @note: add stress changes due to tectonic loading and slip
103    
104     @ivar domain: mining region
105     @type domain: L{Domian}
106     @ivar displacement: displacement field
107     @type displacement: L{Vector} or C{None}
108     @ivar stress: displacement field
109     @type stress: L{Vector} or C{None}
110     @ivar density: initial density distribution
111     @type density: C{float} or L{Scalar}
112     @ivar density_rate: density rate by tag (may be changed of time)
113     @type density_rate: C{dict}
114     @ivar lame_lambda: elasticity coefficient lambda (assumed to be constant over time)
115     @type lame_lambda: C{float} or L{Scalar}
116     @ivar lame_mu: elasticity coefficient mu (assumed to be constant over time)
117     @type lame_mu: C{float} or L{Scalar}
118     @ivar location_of_fixed_displacement: mask of locations and component with zero displacements
119     @type location_of_fixed_displacement: L{Vector} or C{None}
120     """
121     def __init__(self,**kwargs):
122     """
123     """
124     super(LinearElasticStressChange,self).__init__(**kwargs)
125     self.declareParameter(domain=None,
126     displacement=None, \
127     stress=None, \
128     density=1., \
129     lame_lambda=2., \
130     lame_mu=1., \
131     location_of_fixed_displacement=None, \
132     density_rate=None)
133    
134     def doInitialization(self):
135     """
136     initialize time integration
137     """
138     if self.stress == None: self.stress=Tensor(0,Function(self.domain))
139     if self.displacement==None: self.displacement=Scalar(0,Solution(self.domain))
140     self.__pde=LinearPDE(self.domain)
141     self.__pde.setSymmetryOn()
142     A =Tensor4(0.,Function(self.domain))
143     for i in range(self.domain.getDim()):
144     for j in range(self.domain.getDim()):
145     A[i,j,j,i] += self.lame_mu
146     A[i,j,i,j] += self.lame_mu
147     A[i,i,j,j] += self.lame_lambda
148     self.__pde.setValue(A=A,q=self.location_of_fixed_displacement)
149    
150     self.__stress=self.stress
151     self.__displacement=self.displacement
152     self.__first=True
153    
154     def doInitialStep(self):
155     """
156     """
157 gross 998 d=self.domain.getDim()
158 gross 997 self.__pde.setValue(Y=-kronecker(Function(self.domain))[d-1]*9.81*self.density)
159     ddisp=self.__pde.getSolution()
160     deps=symmetric(grad(ddisp))
161     self.stress=self.stress+self.lame_mu*deps+self.lame_lambda*trace(deps)*kronecker(Function(self.domain))
162     self.displacement=self.displacement+ddisp
163     self.__first=False
164    
165     def terminateInitialIteration(self):
166     return not self.__first
167    
168     def terminateIteration(self):
169     return not self.__first
170    
171 gross 998 def doInitialPostprocessing(self):
172 gross 997 self.__stress=self.stress
173     self.__displacement=self.displacement
174 gross 998 self.trace("initial displacement range : %s: %s"%(inf(self.displacement[2]),sup(self.displacement[2])))
175     self.trace("initial stress range : %s: %s"%(inf(self.stress),sup(self.stress)))
176 gross 997
177     def doStepPreprocessing(self,dt):
178     self.stress=self.__stress
179     self.displacement=self.__displacement
180     self.__first=True
181    
182 gross 998 def doStep(self,dt):
183 gross 997 """
184     """
185 gross 998 d=self.domain.getDim()
186 gross 997 if not self.density_rate == None:
187     self.__pde.setValue(Y=-kronecker(Function(self.domain))[d-1]*9.81*self.density_rate)
188     ddisp=self.__pde.getSolution()
189     deps=symmetric(grad(ddisp))
190     self.stress=self.stress+dt*(self.lame_mu*deps+self.lame_lambda*trace(deps)*kronecker(Function(self.domain)))
191     self.displacement=self.displacement+dt*ddisp
192     self.__first=False
193    
194     def doStepPostprocessing(self, dt):
195     self.__stress=self.stress
196     self.__displacement=self.displacement
197 gross 998 self.trace("displacement range : %s: %s"%(inf(self.displacement[2]),sup(self.displacement[2])))
198     self.trace("stress range : %s: %s"%(inf(self.stress),sup(self.stress)))
199 gross 1000
200     class CoulombFailureStress(ParameterSet):
201     """
202     calculates the Coulomb failure stress an planes of a given orientation.
203    
204     @ivar domain: mining region
205     @type domain: L{Domian}
206     @ivar displacement: displacement field
207     @type displacement: L{Vector} or C{None}
208     @ivar stress: displacement field
209     @type stress: L{Vector} or C{None}
210     @ivar density: initial density distribution
211     @type density: C{float} or L{Scalar}
212     @ivar density_rate: density rate by tag (may be changed of time)
213     @type density_rate: C{dict}
214     @ivar lame_lambda: elasticity coefficient lambda (assumed to be constant over time)
215     @type lame_lambda: C{float} or L{Scalar}
216     @ivar lame_mu: elasticity coefficient mu (assumed to be constant over time)
217     @type lame_mu: C{float} or L{Scalar}
218     @ivar location_of_fixed_displacement: mask of locations and component with zero displacements
219     @type location_of_fixed_displacement: L{Vector} or C{None}
220     """
221     def __init__(self,**kwargs):
222     """
223     """
224     super(CoulombFailureStress,self).__init__(**kwargs)
225     self.declareParameter(stress=numarray.zeros((3,3)),
226     friction_coefficient=0.,
227     normal=numarray.array([1.,0.,0.]))
228     def cfs(self):
229     """
230     returns Coulomb failure stress
231     """
232     sn=matrixmult(self.stress,self.normal)
233     nsn=inner(self.normal,sn)
234     nssn=inner(sn,sn)
235     return (sqrt(nssn-nsn**2)-self.friction_coefficient*nsn)/length(self.normal)
236    

  ViewVC Help
Powered by ViewVC 1.1.26