# Annotation of /trunk/finley/test/python/linearElastic.py

Revision 617 - (hide annotations)
Wed Mar 22 02:58:17 2006 UTC (15 years, 1 month ago) by elspeth
File MIME type: text/x-python
File size: 1655 byte(s)
```More copyright.

```
 1 gross 447 """ \$Id\$ """ 2 elspeth 617 __copyright__=""" Copyright (c) 2006 by ACcESS MNRF 3 4 Primary Business: Queensland, Australia""" 5 __license__="""Licensed under the Open Software License version 3.0 6 7 gross 447 from esys.escript import * 8 from esys.escript.linearPDEs import LinearPDE 9 import esys.finley as finley 10 11 press0=1. 12 lamb=1. 13 nu=0.3 14 15 # this sets the hook tensor: 16 def setHookTensor(w,l,n): 17 C=Tensor4(0.,w) 18 for i in range(w.getDim()): 19 for j in range(w.getDim()): 20 C[i,i,j,j]+=l 21 C[j,i,j,i]+=n 22 C[j,i,i,j]+=n 23 24 return C 25 26 # generate mesh: here 10x20 mesh of order 2 27 domain=finley.Rectangle(10,20,1,l0=0.5,l1=1.0) 28 # get handel to nodes and elements: 29 e=Function(domain) 30 fe=FunctionOnBoundary(domain) 31 n=ContinuousFunction(domain) 32 # 33 # set a mask msk of type vector which is one for nodes and components set be a constraint: 34 # 35 msk=whereZero(n.getX()[0])*[1.,1.] 36 # 37 # set the normal stress components on face elements. 38 # faces tagged with 21 get the normal stress [0,-press0]. 39 # 40 # now the pressure is set to zero for x0 coordinates bigger then 0.1 41 press=whereNegative(fe.getX()[0]-0.1)*200000.*[1.,0.] 42 # assemble the linear system: 43 mypde=LinearPDE(domain) 44 mypde.setValue(A=setHookTensor(e,lamb,nu),y=press,q=msk,r=[0,0]) 45 mypde.setSymmetryOn() 46 # solve for the displacements: 47 u_d=mypde.getSolution() 48 # get the gradient and calculate the stress: 49 g=grad(u_d) 50 stress=lamb*trace(g)*kronecker(domain)+nu*(g+transpose(g)) 51 # write the hydrostatic pressure: 52 saveVTK("result.xml",displacement=u_d,pressure=trace(stress)/domain.getDim())

Name Value
svn:executable *