/[escript]/trunk/doc/examples/cookbook/heatrefraction_mesher002.py
ViewVC logotype

Contents of /trunk/doc/examples/cookbook/heatrefraction_mesher002.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2658 - (show annotations)
Thu Sep 10 02:58:44 2009 UTC (9 years, 7 months ago) by ahallam
File MIME type: text/x-python
File size: 3897 byte(s)
Updates to all files scripts to support MPI testing proceedure. Updates to cookbook, new section on functino spaces/domains (needs work). Finalising first 3 chapters for editing.
1
2 ########################################################
3 #
4 # Copyright (c) 2003-2009 by University of Queensland
5 # 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 __copyright__="""Copyright (c) 2003-2009 by University of Queensland
15 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 Author: Antony Hallam antony.hallam@uqconnect.edu.au
24 """
25 ############################################################FILE HEADER
26 # heatrefraction_mesher002.py
27 # Create a 3 block fault and overburden style model in 2d using pycad
28 # meshing tools.
29
30 #######################################################EXTERNAL MODULES
31 from esys.pycad import *
32 from esys.pycad.gmsh import Design
33 from esys.finley import MakeDomain
34 import os
35 import numpy as np
36 from math import *
37 from cblib import needdirs, getLoopCoords
38
39 #################################################ESTABLISHING VARIABLES
40 # where to put output files
41 save_path= os.path.join("data","heatrefrac002")
42 needdirs([save_path])
43
44 ################################################ESTABLISHING PARAMETERS
45 # Overall Domain
46 p0=Point(0.0, 0.0, 0.0)
47 p1=Point(0.0, -6000.0, 0.0)
48 p2=Point(5000.0, -6000.0, 0.0)
49 p3=Point(5000.0, 0.0, 0.0)
50
51 ####################################################DOMAIN CONSTRUCTION
52 l01=Line(p0, p1)
53 l12=Line(p1, p2)
54 l23=Line(p2, p3)
55 l30=Line(p3, p0)
56
57 c=CurveLoop(l01, l12, l23, l30)
58
59 # Generate Material Boundary
60 p4=Point(0.0, -2400.0, 0.0)
61 p5=Point(2000.0, -2400.0, 0.0)
62 p6=Point(3000.0, -6000.0, 0.0)
63 p7=Point(5000.0, -2400.0, 0.0)
64
65 # Create TOP BLOCK
66 tbl1=Line(p0,p4)
67 tbl2=Line(p4,p5)
68 tbl3=Line(p5,p7)
69 tbl4=Line(p7,p3)
70 tblockloop = CurveLoop(tbl1,tbl2,tbl3,tbl4,l30)
71 tblock = PlaneSurface(tblockloop)
72 # python plotting polygon
73 tpg = loopcoords(tblockloop)
74 np.savetxt(os.path.join(save_path,"toppg"),tpg,delimiter=" ")
75
76 # Create BOTTOM BLOCK LEFT
77 bbll1=Line(p4,p1)
78 bbll2=Line(p1,p6)
79 bbll3=Line(p6,p5)
80 bbll4=-tbl2
81 bblockloopl = CurveLoop(bbll1,bbll2,bbll3,bbll4)
82 bblockl = PlaneSurface(bblockloopl)
83
84 #clockwise check
85 #bblockloopl2=CurveLoop(mysp,Line(x2,p2),Line(p2,p1),Line(p1,x1))
86 bpg = loopcoords(bblockloopl)
87 np.savetxt(os.path.join(save_path,"botpgl"),bpg,delimiter=" ")
88
89 # Create BOTTOM BLOCK RIGHT
90 bbrl1=Line(p6,p2)
91 bbrl2=Line(p2,p7)
92 bbrl3=-tbl3
93 bbrl4=-bbll3
94 bblockloopr = CurveLoop(bbrl1,bbrl2,bbrl3,bbrl4)
95 bblockr = PlaneSurface(bblockloopr)
96
97 #clockwise check
98 #bblockloopr2=CurveLoop(mysp,Line(x2,p2),Line(p2,p1),Line(p1,x1))
99 bpg = loopcoords(bblockloopr)
100 np.savetxt(os.path.join(save_path,"botpgr"),bpg,delimiter=" ")
101
102 #############################################EXPORTING MESH FOR ESCRIPT
103 # Create a Design which can make the mesh
104 d=Design(dim=2, element_size=200)
105 # Add the subdomains and flux boundaries.
106 d.addItems(PropertySet("top",tblock),\
107 PropertySet("bottomleft",bblockl),\
108 PropertySet("bottomright",bblockr),\
109 PropertySet("linebottom",bbll21, bbrl1))
110
111 # Create the geometry, mesh and Escript domain
112 d.setScriptFileName(os.path.join(save_path,\
113 "heatrefraction_mesh003.geo"))
114
115 d.setMeshFileName(os.path.join(save_path,"heatrefraction_mesh003.msh"))
116 domain=MakeDomain(d, integrationOrder=-1,\
117 reducedIntegrationOrder=-1,\
118 optimizeLabeling=True)
119 # Create a file that can be read back in to python with
120 # mesh=ReadMesh(fileName)
121 domain.write(os.path.join(save_path,"heatrefraction_mesh003.fly"))
122

  ViewVC Help
Powered by ViewVC 1.1.26