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

Diff of /trunk/doc/examples/cookbook/example09m.py

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

revision 3056 by ahallam, Sun Jul 4 21:52:36 2010 UTC revision 3057 by ahallam, Tue Jul 6 02:10:59 2010 UTC
# Line 73  l56=Line(p5, p6) Line 73  l56=Line(p5, p6)
73  l67=Line(p6, p7)  l67=Line(p6, p7)
74  l74=Line(p7, p4)  l74=Line(p7, p4)
75  l45=Line(p4, p5)  l45=Line(p4, p5)
 l05=Line(p0, p5)  
 l16=Line(p1, p6)  
 l27=Line(p2, p7)  
 l34=Line(p3, p4)  
76    
77  # Join line segments to create domain boundaries and then surfaces  # Join line segments to create domain boundaries and then surfaces
78  ctop=CurveLoop(l01, l12, l23, l30);     stop=PlaneSurface(ctop)  ctop=CurveLoop(l01, l12, l23, l30);     stop=PlaneSurface(ctop)
79  cbot=CurveLoop(-l67, -l56, -l45, -l74); sbot=PlaneSurface(cbot)  cbot=CurveLoop(-l67, -l56, -l45, -l74); sbot=PlaneSurface(cbot)
 cx0 =CurveLoop(l05, l56, -l16, -l01);   sx0=PlaneSurface(cx0)  
 cy0 =CurveLoop(-l30, l34, l45, -l05);   sy0=PlaneSurface(cy0)  
 cxy =CurveLoop(-l34, -l23, l27, l74);   sxy=PlaneSurface(cxy)  
 cyx =CurveLoop(-l12, l16, l67, -l27);   syx=PlaneSurface(cyx)  
   
 # Create the volume.  
 sloop=SurfaceLoop(stop,sbot,sx0,sy0,sxy,syx)  
 model=Volume(sloop)  
   
80    
81  #create interface  #create interface
82  sspl=51 #number of discrete points in spline  sspl=51 #number of discrete points in spline
# Line 137  nsintf_ar=[] Line 124  nsintf_ar=[]
124  for i in range(0,sspl):  for i in range(0,sspl):
125      #create all the points required      #create all the points required
126      point_ar=[Point(xs[i],ys[j],zs[i,j]) for j in range(0,sspl)]      point_ar=[Point(xs[i],ys[j],zs[i,j]) for j in range(0,sspl)]
127    
128        if (i == 0):
129            ip0=point_ar[0]
130            ip3=point_ar[-1]
131        if (i == sspl-1):
132            ip1=point_ar[0]
133            ip2=point_ar[-1]
134    
135      #create the spline and store it      #create the spline and store it
136      spl_ar.append(Spline(*tuple(point_ar)))      spl_ar.append(Spline(*tuple(point_ar)))
137      #create the lines along the x axis and x axis at ymax      #create the lines along the x axis and x axis at ymax
138      if (i < sspl-1):      if (i > 0):
139          #print i,xs[i],ys[0],zs[i,0]          #print i,xs[i],ys[0],zs[i,0]
140          tp0=Point(xs[i],ys[0],zs[i,0])          tp2=point_ar[0]; tp3=point_ar[-1]
141          tp1=Point(xs[i+1],ys[0],zs[i+1,0])          tl0=Line(cp0,tp2); tl1=Line(cp1,tp3)
         tp2=Point(xs[i],ys[-1],zs[i,-1])  
         tp3=Point(xs[i+1],ys[-1],zs[i+1,-1])  
         tl0=Line(tp0,tp1); tl1=Line(tp2,tp3)  
142          linex0_ar.append(tl0); linexy_ar.append(tl1)          linex0_ar.append(tl0); linexy_ar.append(tl1)
143        cp0=point_ar[0]; cp1=point_ar[-1]
144                    
145  for i in range(0,sspl):  for i in range(0,sspl):
146      #create the mini surfaces via curveloops and then ruledsurfaces      #create the mini surfaces via curveloops and then ruledsurfaces
# Line 161  for i in range(0,sspl): Line 154  for i in range(0,sspl):
154  sintf=SurfaceLoop(*tuple(sintf_ar))  sintf=SurfaceLoop(*tuple(sintf_ar))
155    
156  # for each side  # for each side
157  ip0=Point(xs[0],ys[0],zs[0,0])  #ip0=Point(xs[0],ys[0],zs[0,0])
158  ip1=Point(xs[-1],ys[0],zs[-1,0])  #ip1=Point(xs[-1],ys[0],zs[-1,0])
159  ip2=Point(xs[-1],ys[-1],zs[-1,-1])  #ip2=Point(xs[-1],ys[-1],zs[-1,-1])
160  ip3=Point(xs[0],ys[-1],zs[0,-1])  #ip3=Point(xs[0],ys[-1],zs[0,-1])
161    
162  linte_ar=[]; #lines for vertical edges  linte_ar=[]; #lines for vertical edges
163  linhe_ar=[]; #lines for horizontal edges  linhe_ar=[]; #lines for horizontal edges
# Line 183  linhe_ar.append(Line(ip2,ip3)) Line 176  linhe_ar.append(Line(ip2,ip3))
176  linhe_ar.append(Line(ip3,ip0))  linhe_ar.append(Line(ip3,ip0))
177    
178  cintfa_ar=[]; cintfb_ar=[] #curveloops for above and below interface on sides  cintfa_ar=[]; cintfb_ar=[] #curveloops for above and below interface on sides
179  cintfa_ar.append(CurveLoop(linte_ar[0],linhe_ar[0],-linte_ar[2],-l01))  cintfa_ar.append(CurveLoop(-linte_ar[2],-l01,linte_ar[0],*tuple(linex0_ar)))
180  cintfa_ar.append(CurveLoop(linte_ar[2],linhe_ar[1],-linte_ar[4],-l12))  cintfa_ar.append(CurveLoop(linte_ar[2],spl_ar[-1],-linte_ar[4],-l12))
181  cintfa_ar.append(CurveLoop(linte_ar[4],linhe_ar[2],-linte_ar[6],-l23))  cintfa_ar.append(CurveLoop(-linte_ar[6],-l23,linte_ar[4],*tuple(linexy_ar)))
182  cintfa_ar.append(CurveLoop(linte_ar[6],linhe_ar[3],-linte_ar[0],-l30))  cintfa_ar.append(CurveLoop(linte_ar[6],-spl_ar[0],-linte_ar[0],-l30))
183    
184  cintfb_ar.append(CurveLoop(linte_ar[1],l56,-linte_ar[3],-linhe_ar[0]))  #cintfb_ar.append(CurveLoop(linte_ar[1],l56,-linte_ar[3],-linhe_ar[0]))
185  cintfb_ar.append(CurveLoop(linte_ar[3],l67,-linte_ar[5],-linhe_ar[1]))  #cintfb_ar.append(CurveLoop(linte_ar[3],l67,-linte_ar[5],-linhe_ar[1]))
186  cintfb_ar.append(CurveLoop(linte_ar[5],l74,-linte_ar[7],-linhe_ar[2]))  #cintfb_ar.append(CurveLoop(linte_ar[5],l74,-linte_ar[7],-linhe_ar[2]))
187  cintfb_ar.append(CurveLoop(linte_ar[7],l45,-linte_ar[1],-linhe_ar[3]))  #cintfb_ar.append(CurveLoop(linte_ar[7],l45,-linte_ar[1],-linhe_ar[3]))
188    
189    cintfb_ar.append(-CurveLoop(linte_ar[3],-l56,-linte_ar[1],*tuple(linex0_ar)))
190    cintfb_ar.append(-CurveLoop(linte_ar[5],-l67,-linte_ar[3],-spl_ar[-1]))
191    cintfb_ar.append(CurveLoop(linte_ar[5],l74,-linte_ar[7],*tuple(linexy_ar)))
192    cintfb_ar.append(CurveLoop(linte_ar[7],l45,-linte_ar[1],spl_ar[0]))
193    
194    
195    
196  sintfa_ar=[PlaneSurface(cintfa_ar[i]) for i in range(0,4)]  sintfa_ar=[PlaneSurface(cintfa_ar[i]) for i in range(0,4)]
197  sintfb_ar=[PlaneSurface(cintfb_ar[i]) for i in range(0,4)]  sintfb_ar=[PlaneSurface(cintfb_ar[i]) for i in range(0,4)]
198    
199  vintfa=Volume(SurfaceLoop(stop,*tuple(sintfa_ar+nsintf_ar)))  vintfa=Volume(SurfaceLoop(stop,*tuple(sintfa_ar+sintf_ar)))
200  vintfb=Volume(SurfaceLoop(sbot,*tuple(sintfb_ar+sintf_ar)))  vintfb=Volume(SurfaceLoop(sbot,*tuple(sintfb_ar+sintf_ar)))
201    
202    # Create the volume.
203    #sloop=SurfaceLoop(stop,sbot,*tuple(sintfa_ar+sintfb_ar))
204    #model=Volume(sloop)
205    
206    
207  #############################################EXPORTING MESH FOR ESCRIPT  #############################################EXPORTING MESH FOR ESCRIPT
208  # Create a Design which can make the mesh  # Create a Design which can make the mesh
209  d=Design(dim=3, element_size=200.0*m)  d=Design(dim=3, element_size=200.0*m)
210  # Add the subdomains and flux boundaries.  # Add the subdomains and flux boundaries.
211  d.addItems(model,stop,sbot,sx0,sy0,sxy,syx)  #d.addItems(stop,sbot)
212  #d.addItems(PropertySet('intf',sintf))  #d.addItems(PropertySet('intf',sintf))
213    
214  d.addItems(PropertySet('vintfa',vintfa),PropertySet('vintfb',vintfb))  #d.addItems(PropertySet('vintfa',vintfa))
215    #d.addItems(PropertySet('vintfb',vintfb))
216    #d.addItems(PropertySet('l1',linte_ar[0]))
217    d.addItems(*tuple(sintfb_ar))
218    d.addItems(sbot)
219    d.addItems(*tuple(sintf_ar))
220    
221  d.setScriptFileName(os.path.join(save_path,"example09m.geo"))  d.setScriptFileName(os.path.join(save_path,"example09m.geo"))
222    

Legend:
Removed from v.3056  
changed lines
  Added in v.3057

  ViewVC Help
Powered by ViewVC 1.1.26