/[escript]/trunk/doc/cookbook/heatrefraction001.py
ViewVC logotype

Diff of /trunk/doc/cookbook/heatrefraction001.py

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

revision 2596 by ahallam, Wed Aug 5 22:42:36 2009 UTC revision 2597 by ahallam, Thu Aug 6 03:30:09 2009 UTC
# Line 38  from esys.escript.pdetools import * Line 38  from esys.escript.pdetools import *
38  qin=70*Milli*W/(m*m) #our heat source temperature is now zero  qin=70*Milli*W/(m*m) #our heat source temperature is now zero
39  Ti=290.15*K # Kelvin #the starting temperature of our iron bar  Ti=290.15*K # Kelvin #the starting temperature of our iron bar
40    
41    # the folder to gett our outputs from, leave blank "" for script path -
42    # note these depen. are generated from heatrefraction_mesher001.py
43    saved_path = "data/heatrefrac001"
44    
45  ###### 2 BLOCK MODEL #########  ###### 2 BLOCK MODEL #########
46  ## DOMAIN  ## DOMAIN
47  ## Anticline  ## Anticline
48  mymesh = ReadMesh("heatrefraction_mesh001.fly")  mymesh=ReadMesh(os.path.join(saved_path,"heatrefraction_mesh001.fly"))
49  tpg = np.loadtxt("toppg")  tpg = np.loadtxt(os.path.join(saved_path,"toppg"))
50  tpgx = tpg[:,0]  tpgx = tpg[:,0]
51  tpgy = tpg[:,1]  tpgy = tpg[:,1]
52  bpg = np.loadtxt("botpg")  bpg = np.loadtxt(os.path.join(saved_path,"botpg"))
53  bpgx = bpg[:,0]  bpgx = bpg[:,0]
54  bpgy = bpg[:,1]  bpgy = bpg[:,1]
55  ## Syncline  ## Syncline
# Line 73  mypde.setValue(A=kappa*kronecker(mymesh) Line 77  mypde.setValue(A=kappa*kronecker(mymesh)
77  # ... set initial temperature ....  # ... set initial temperature ....
78  x=mymesh.getX()  x=mymesh.getX()
79    
80  qH=qin*whereZero(x[1]-inf(x[1]))  qH=Scalar(0,FunctionOnBoundary(mymesh))
81    qH.setTaggedValue("linebottom",qin)
82  mypde.setValue(q=whereZero(x[1]),r=Ti)  mypde.setValue(q=whereZero(x[1]),r=Ti)
83  mypde.setValue(y=qH)#,r=17*Celsius)  mypde.setValue(y=qH)#,r=17*Celsius)
84    
85  # get steady state solution and export to vtk.  # get steady state solution and export to vtk.
86  T=mypde.getSolution()  T=mypde.getSolution()
87  saveVTK("tempheatrefract.xml",sol=T, q=-kappa*grad(T))  #saveVTK("tempheatrefract.xml",sol=T, q=-kappa*grad(T))
88    
89  # rearrage mymesh to suit solution function space        # rearrage mymesh to suit solution function space      
90  oldspacecoords=mymesh.getX()  oldspacecoords=mymesh.getX()
# Line 131  CKL = pl.fill(tpgx,tpgy,'brown',bpgx,bpg Line 136  CKL = pl.fill(tpgx,tpgy,'brown',bpgx,bpg
136  #~ CK = pl.contourf(xi,yi,ziK,2)  #~ CK = pl.contourf(xi,yi,ziK,2)
137  CS = pl.contour(xi,yi,zi,5,linewidths=0.5,colors='k')  CS = pl.contour(xi,yi,zi,5,linewidths=0.5,colors='k')
138  pl.clabel(CS, inline=1, fontsize=8)  pl.clabel(CS, inline=1, fontsize=8)
139    pl.title("Heat Refraction across a clinal structure.")
140    pl.xlabel("Horizontal Displacement (m)")
141    pl.ylabel("Depth (m)")
142  #~ CB = pl.colorbar(CS, shrink=0.8, extend='both')  #~ CB = pl.colorbar(CS, shrink=0.8, extend='both')
143  pl.savefig("temp.png")  pl.savefig(os.path.join(saved_path,"heatrefraction001_cont.png"))
144    
145  QUIV = pl.quiver(qulocs[:,0],qulocs[:,1],qu[:,0],qu[:,1],angles='xy',color="white")  QUIV=pl.quiver(qulocs[:,0],qulocs[:,1],qu[:,0],qu[:,1],angles='xy',color="white")
146  pl.savefig("temp2.png")  pl.title("Heat Refraction across a clinal structure \n with gradient quivers.")
147    pl.savefig(os.path.join(saved_path,"heatrefraction001_contqu.png"))
148    
149    pl.clf()
150    CKL = pl.fill(tpgx,tpgy,'white',bpgx,bpgy,'white',zorder=-1000)
151    pl.savefig('model.png')

Legend:
Removed from v.2596  
changed lines
  Added in v.2597

  ViewVC Help
Powered by ViewVC 1.1.26