/[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 1000 - (hide annotations)
Wed Feb 28 00:34:42 2007 UTC (12 years, 5 months ago) by gross
File MIME type: text/x-python
File size: 9487 byte(s)
coulomb stress failure added
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 tag_map: a tagmap
68     @type tag_map: L{TagMap}
69     @ivar mass_rate: rate of change mass in a tagged region
70     @type mass_rate: C{dict}
71     @ivar density_rate: density in each region
72     @type density_rate: C{Data}
73     """
74     def __init__(self,**kwargs):
75     """
76     """
77     super(DensityChange,self).__init__(**kwargs)
78     self.declareParameter(domain=None,
79     tag_map=None,
80     mass_rate={},
81     density_rate=None)
82    
83     def doInitialization(self):
84     """
85     initialize time integration
86     """
87     self.__volumes={}
88     for i in self.tag_map.getNames():
89     c=Scalar(0.,Function(self.domain))
90     self.tag_map.insert(c,**{i : 1.})
91     c=integrate(c)
92     if c>0: self.__volumes[i] = c
93    
94     def doStepPreprocessing(self, dt):
95     self.density_rate=Scalar(0.,Function(self.domain))
96     d={}
97     for i in self.__volumes:
98 gross 998 if self.mass_rate.has_key(i): d[i]=-self.mass_rate[i]/self.__volumes[i]
99 gross 997 self.tag_map.insert(self.density_rate,**d)
100    
101     class LinearElasticStressChange(Model):
102     """
103     calculates the stress according to an initial garvity field and a consecutive
104     chenge of density based on linar elastic material model. It is assumed that
105     the lame coefficients don't change over time, the direction of gravity is
106     pointing into the x_{dim} direction.
107    
108     @note: add stress changes due to tectonic loading and slip
109    
110     @ivar domain: mining region
111     @type domain: L{Domian}
112     @ivar tag_map: a tagmap
113     @type tag_map: L{TagMap}
114     @ivar displacement: displacement field
115     @type displacement: L{Vector} or C{None}
116     @ivar stress: displacement field
117     @type stress: L{Vector} or C{None}
118     @ivar density: initial density distribution
119     @type density: C{float} or L{Scalar}
120     @ivar density_rate: density rate by tag (may be changed of time)
121     @type density_rate: C{dict}
122     @ivar lame_lambda: elasticity coefficient lambda (assumed to be constant over time)
123     @type lame_lambda: C{float} or L{Scalar}
124     @ivar lame_mu: elasticity coefficient mu (assumed to be constant over time)
125     @type lame_mu: C{float} or L{Scalar}
126     @ivar location_of_fixed_displacement: mask of locations and component with zero displacements
127     @type location_of_fixed_displacement: L{Vector} or C{None}
128     """
129     def __init__(self,**kwargs):
130     """
131     """
132     super(LinearElasticStressChange,self).__init__(**kwargs)
133     self.declareParameter(domain=None,
134     tag_map=None,
135     displacement=None, \
136     stress=None, \
137     density=1., \
138     lame_lambda=2., \
139     lame_mu=1., \
140     location_of_fixed_displacement=None, \
141     density_rate=None)
142    
143     def doInitialization(self):
144     """
145     initialize time integration
146     """
147     if self.stress == None: self.stress=Tensor(0,Function(self.domain))
148     if self.displacement==None: self.displacement=Scalar(0,Solution(self.domain))
149     self.__pde=LinearPDE(self.domain)
150     self.__pde.setSymmetryOn()
151     A =Tensor4(0.,Function(self.domain))
152     for i in range(self.domain.getDim()):
153     for j in range(self.domain.getDim()):
154     A[i,j,j,i] += self.lame_mu
155     A[i,j,i,j] += self.lame_mu
156     A[i,i,j,j] += self.lame_lambda
157     self.__pde.setValue(A=A,q=self.location_of_fixed_displacement)
158    
159     self.__stress=self.stress
160     self.__displacement=self.displacement
161     self.__first=True
162    
163     def doInitialStep(self):
164     """
165     """
166 gross 998 d=self.domain.getDim()
167 gross 997 self.__pde.setValue(Y=-kronecker(Function(self.domain))[d-1]*9.81*self.density)
168     ddisp=self.__pde.getSolution()
169     deps=symmetric(grad(ddisp))
170     self.stress=self.stress+self.lame_mu*deps+self.lame_lambda*trace(deps)*kronecker(Function(self.domain))
171     self.displacement=self.displacement+ddisp
172     self.__first=False
173    
174     def terminateInitialIteration(self):
175     return not self.__first
176    
177     def terminateIteration(self):
178     return not self.__first
179    
180 gross 998 def doInitialPostprocessing(self):
181 gross 997 self.__stress=self.stress
182     self.__displacement=self.displacement
183 gross 998 self.trace("initial displacement range : %s: %s"%(inf(self.displacement[2]),sup(self.displacement[2])))
184     self.trace("initial stress range : %s: %s"%(inf(self.stress),sup(self.stress)))
185 gross 997
186     def doStepPreprocessing(self,dt):
187     self.stress=self.__stress
188     self.displacement=self.__displacement
189     self.__first=True
190    
191 gross 998 def doStep(self,dt):
192 gross 997 """
193     """
194 gross 998 d=self.domain.getDim()
195 gross 997 if not self.density_rate == None:
196     self.__pde.setValue(Y=-kronecker(Function(self.domain))[d-1]*9.81*self.density_rate)
197     ddisp=self.__pde.getSolution()
198     deps=symmetric(grad(ddisp))
199     self.stress=self.stress+dt*(self.lame_mu*deps+self.lame_lambda*trace(deps)*kronecker(Function(self.domain)))
200     self.displacement=self.displacement+dt*ddisp
201     self.__first=False
202    
203     def doStepPostprocessing(self, dt):
204     self.__stress=self.stress
205     self.__displacement=self.displacement
206 gross 998 self.trace("displacement range : %s: %s"%(inf(self.displacement[2]),sup(self.displacement[2])))
207     self.trace("stress range : %s: %s"%(inf(self.stress),sup(self.stress)))
208 gross 1000
209     class CoulombFailureStress(ParameterSet):
210     """
211     calculates the Coulomb failure stress an planes of a given orientation.
212    
213     @ivar domain: mining region
214     @type domain: L{Domian}
215     @ivar tag_map: a tagmap
216     @type tag_map: L{TagMap}
217     @ivar displacement: displacement field
218     @type displacement: L{Vector} or C{None}
219     @ivar stress: displacement field
220     @type stress: L{Vector} or C{None}
221     @ivar density: initial density distribution
222     @type density: C{float} or L{Scalar}
223     @ivar density_rate: density rate by tag (may be changed of time)
224     @type density_rate: C{dict}
225     @ivar lame_lambda: elasticity coefficient lambda (assumed to be constant over time)
226     @type lame_lambda: C{float} or L{Scalar}
227     @ivar lame_mu: elasticity coefficient mu (assumed to be constant over time)
228     @type lame_mu: C{float} or L{Scalar}
229     @ivar location_of_fixed_displacement: mask of locations and component with zero displacements
230     @type location_of_fixed_displacement: L{Vector} or C{None}
231     """
232     def __init__(self,**kwargs):
233     """
234     """
235     super(CoulombFailureStress,self).__init__(**kwargs)
236     self.declareParameter(stress=numarray.zeros((3,3)),
237     friction_coefficient=0.,
238     normal=numarray.array([1.,0.,0.]))
239     def cfs(self):
240     """
241     returns Coulomb failure stress
242     """
243     sn=matrixmult(self.stress,self.normal)
244     nsn=inner(self.normal,sn)
245     nssn=inner(sn,sn)
246     return (sqrt(nssn-nsn**2)-self.friction_coefficient*nsn)/length(self.normal)
247    

  ViewVC Help
Powered by ViewVC 1.1.26