/[escript]/trunk/doc/examples/usersguide/slip.py
ViewVC logotype

Annotation of /trunk/doc/examples/usersguide/slip.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2667 - (hide annotations)
Thu Sep 17 01:49:11 2009 UTC (11 years ago) by jfenwick
File MIME type: text/x-python
File size: 2275 byte(s)
Renamed the main cookbook tex file to match our convention.
Replaced doc/cookbook/figures/heatrefraction002contqu.pdf with
a version which is actually pdf. However it needs to be regnerated since
it it sideways.

The examples have had their copyright notices fixed (dates were too early).
sb2.py has been removed since it uses pyvisi.

scons will now build the cookbook as parts of a docs build.
Also in reposnse to :
scons cookbook_pdf


1 gross 2654
2     ########################################################
3     #
4 jfenwick 2667 # Copyright (c) 2009 by University of Queensland
5 gross 2654 # Earth Systems Science Computational Center (ESSCC)
6     # http://www.uq.edu.au/esscc
7     #
8     # Primary Business: Queensland, Australia
9     # Licensed under the Open Software License version 3.0
10     # http://www.opensource.org/licenses/osl-3.0.php
11     #
12     ########################################################
13    
14 jfenwick 2667 __copyright__="""Copyright (c) 2009 by University of Queensland
15 gross 2654 Earth Systems Science Computational Center (ESSCC)
16     http://www.uq.edu.au/esscc
17     Primary Business: Queensland, Australia"""
18     __license__="""Licensed under the Open Software License version 3.0
19     http://www.opensource.org/licenses/osl-3.0.php"""
20     __url__="https://launchpad.net/escript-finley"
21    
22    
23     from esys.escript import *
24     from esys.escript.linearPDEs import LinearPDE
25     from esys.escript.models import FaultSystem
26     from esys.finley import Rectangle
27 gross 2663 from esys.escript.unitsSI import DEG
28 gross 2654 #... set some parameters ...
29     lam=1.
30     mu=1
31     slip_max=1.
32    
33     mydomain = Rectangle(l0=1.,l1=1.,n0=16, n1=16) # n1 need to be multiple of 4!!!
34     # .. create the fault system
35     fs=FaultSystem(dim=2)
36 gross 2663 fs.addFault(V0=[0.5,0.25], strikes=90*DEG, ls=0.5, tag=1)
37 gross 2654 # ... create a slip distribution on the fault:
38     p, m=fs.getParametrization(mydomain.getX(),tag=1)
39     p0,p1= fs.getW0Range(tag=1)
40     s=m*(p-p0)*(p1-p)/((p1-p0)/2)**2*slip_max*[0.,1.]
41     # ... calculate stress according to slip:
42     D=symmetric(grad(s))
43     chi, d=fs.getSideAndDistance(D.getFunctionSpace().getX(),tag=1)
44     sigma_s=(mu*D+lam*trace(D)*kronecker(mydomain))*chi
45     #... open symmetric PDE ...
46     mypde=LinearPDE(mydomain)
47     mypde.setSymmetryOn()
48     #... set coefficients ...
49     C=Tensor4(0.,Function(mydomain))
50     for i in range(mydomain.getDim()):
51     for j in range(mydomain.getDim()):
52     C[i,i,j,j]+=lam
53     C[j,i,j,i]+=mu
54     C[j,i,i,j]+=mu
55     # ... fix displacement in normal direction
56     x=mydomain.getX()
57     msk=whereZero(x[0])*[1.,0.] + whereZero(x[0]-1.)*[1.,0.] \
58     +whereZero(x[1])*[0.,1.] + whereZero(x[1]-1.)*[0.,1.]
59     mypde.setValue(A=C,X=-0.5*sigma_s,q=msk)
60     #... solve pde ...
61     mypde.getSolverOptions().setVerbosityOn()
62     v=mypde.getSolution()
63     # .. write the displacement to file:
64     D=symmetric(grad(v))
65     sigma=(mu*D+lam*trace(D)*kronecker(mydomain))+0.5*sigma_s
66     saveVTK("slip.vtu",disp=v+0.5*chi*s, stress= sigma)

  ViewVC Help
Powered by ViewVC 1.1.26