/[escript]/trunk/finley/test/python/slip_stress_old.py
ViewVC logotype

Diff of /trunk/finley/test/python/slip_stress_old.py

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 887 by gross, Thu Nov 2 07:17:07 2006 UTC revision 892 by gross, Wed Nov 8 06:14:29 2006 UTC
# Line 67  class SlippingFault(SaddlePointProblem): Line 67  class SlippingFault(SaddlePointProblem):
67        def solve_f(self,u,p,tol=1.e-8):        def solve_f(self,u,p,tol=1.e-8):
68           self.__pde_u.setTolerance(tol)           self.__pde_u.setTolerance(tol)
69           self.__pde_u.setValue(y_contact=p)           self.__pde_u.setValue(y_contact=p)
70             print "p:",inf(p),sup(p)
71             print "u:",inf(u),sup(u)
72             self.__pde_u.setValue(y_contact=p)
73           return  self.__pde_u.getSolution()           return  self.__pde_u.getSolution()
74    
75        def solve_g(self,u,tol=1.e-8):        def solve_g(self,u,tol=1.e-8):
76           dp=(self.slip-jump(u))*lam_lmbd/FunctionOnContactZero(self.domain).getX()           dp=-(self.slip-jump(u))*lam_lmbd/FunctionOnContactZero(self.domain).getSize()
77             print dp
78           return  dp           return  dp
79    
80    
81  dom=ReadMesh("meshfault3D.fly")  dom=ReadMesh("meshfault3D.fly",integrationOrder=-1)
82  prop=SlippingFault(dom)  prop=SlippingFault(dom)
83  d=dom.getDim()  d=dom.getDim()
84  x=dom.getX()[d-1]  x=dom.getX()[0]
85    # x=dom.getX()[d-1]
86  mask=whereZero(x-inf(x))*numarray.ones((d,))  mask=whereZero(x-inf(x))*numarray.ones((d,))
87  s=numarray.array([0.,1.,1.])  s=numarray.array([0.,1.,1.])
88  x=FunctionOnContactZero(dom).getX()  x=FunctionOnContactZero(dom).getX()
# Line 93  for i in range(3): Line 98  for i in range(3):
98  u0=Vector(0.,Solution(dom))  u0=Vector(0.,Solution(dom))
99  p0=Vector(1.,FunctionOnContactZero(dom))  p0=Vector(1.,FunctionOnContactZero(dom))
100  prop.initialize(fixed_u_mask=mask,slip=Data(s,FunctionOnContactZero(dom)), density=rho,lmbd=lam_lmbd, mu=lam_mu)  prop.initialize(fixed_u_mask=mask,slip=Data(s,FunctionOnContactZero(dom)), density=rho,lmbd=lam_lmbd, mu=lam_mu)
101  u,p=prop.solve(u0,p0,iter_max=100,tolerance=0.5,accepted_reduction=1.1 )  u,p=prop.solve(u0,p0,iter_max=50,tolerance=0.1,accepted_reduction=0.99)
102  saveVTK("dis.xml",u=u)  saveVTK("dis.xml",u=u)
103  saveVTK("fault.xml",sigma=p,s=jump(u))  saveVTK("fault.xml",sigma=p,s=jump(u))

Legend:
Removed from v.887  
changed lines
  Added in v.892

  ViewVC Help
Powered by ViewVC 1.1.26