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

Contents of /trunk/finley/test/python/debugIntegrate.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 147 - (show annotations)
Fri Aug 12 01:45:47 2005 UTC (15 years ago) by jgs
Original Path: trunk/esys2/finley/test/python/debugIntegrate.py
File MIME type: text/x-python
File size: 1861 byte(s)
erge of development branch dev-02 back to main trunk on 2005-08-12

1
2 from escript.escript import *
3 from escript.linearPDEs import *
4 from finley import finley
5
6 def makeParallelFaultSystem2D(xne=1,yne=1,xl=1.,yfaults=[0.,1.],integrationOrder=0,p1=0,p2=0):
7 """
8 generates a 2D rectangular domain and a system of parallel faults at yfaults
9
10 input data :
11 xne - the number of elements over the horizontal edge of the rectangular mesh
12 xl - length of the horizontal edge of the domain
13 yfault[y0,...yi,...y1] where y0 and y1 are the extremity of the block and yi the fault position in the (Oy) direction
14 output data :
15 a mesh with horizontal discontinuities at yi positions
16 """
17
18 #yne=xne*1./(xl+1.e-15)
19
20 mshs=[]
21 for i in range(len(yfaults)-1):
22 yl=yfaults[i+1]-yfaults[i]
23 msh=finley.Rectangle(xne,yne,1,xl,yl,p1,p2,useElementsOnFace=True)
24 msh1=msh
25 print "mesh done"
26 msh1.setX(msh.getX()+[0.,yfaults[i]])
27 print "Xnew"
28 mshs.append(msh1)
29 print "apres append"
30 #results=finley.joinFaces(mshs)
31 mesh_joint=finley.JoinFaces(mshs)
32 print "apres faces jointes"
33 return mesh_joint
34
35 def ContractProd(M1,M2):
36 msh=M1.getDomain()
37 ContractProd=inner(M1,M2)
38 return ContractProd
39
40 def GetStrainEnergy(stressStateInit,uStateInit,stress,u):
41 strainStateInit=getStrain(uStateInit)
42 strain=getStrain(u)
43 E_Strain0=ContractProd(stressStateInit,strainStateInit).integrate()[0]
44 E_StrainCurrent=ContractProd(stress,strain).integrate()[0]
45 E_StrainE=-E_Strain0+E_StrainCurrent
46 return E_StrainE
47
48 def main():
49 xne=10
50 yne=20
51 xl=10
52 yl=10
53 yfaults=[0.,5.,10.]
54 numberfaults=1
55
56 mesh=makeParallelFaultSystem2D(xne,yne,xl,yfaults,0,0,0)
57
58 t=Tensor(10.,Function(mesh))
59 t1=Tensor(20.,Function(mesh))
60
61 contract_product=ContractProd(t,t1)
62 integral=integrate(contract_product)
63
64 print "integral",integral
65
66 s=Scalar(10.,Function(mesh))
67 print integrate(s)
68
69 print "done"
70
71 main()

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26