Contents of /trunk/doc/examples/usersguide/fluid.py

Revision 5288 - (show annotations)
Tue Dec 2 23:18:40 2014 UTC (4 years, 9 months ago) by sshaw
File MIME type: text/x-python
File size: 2864 byte(s)
```fixing tests for cases where required domains not built
```
 1 from __future__ import division 2 from __future__ import print_function 3 ############################################################################## 4 # 5 # Copyright (c) 2008-2014 by University of Queensland 6 7 # 8 # Primary Business: Queensland, Australia 9 # Licensed under the Open Software License version 3.0 10 11 # 12 # Development until 2012 by Earth Systems Science Computational Center (ESSCC) 13 # Development 2012-2013 by School of Earth Sciences 14 # Development from 2014 by Centre for Geoscience Computing (GeoComp) 15 # 16 ############################################################################## 17 18 ######## August 2008 ######## 19 ########## Leon Graham ########## 20 ## Newtonian fluid using StokesProblemCartesian class## 21 22 from esys.escript import * 23 try: 24 from esys.finley import Rectangle 25 HAVE_FINLEY = True 26 except ImportError: 27 HAVE_FINLEY = False 28 from esys.escript.linearPDEs import LinearPDE 29 from esys.escript.models import StokesProblemCartesian 30 from esys.weipa import saveVTK 31 if not HAVE_FINLEY: 32 print("Finley module not available") 33 else: 34 #physical constants 35 eta=1.0 36 rho=100.0 37 g=10.0 38 39 #solver settings 40 tolerance=1.0e-4 41 max_iter=200 42 t_end=50 43 t=0.0 44 time=0 45 verbose='TRUE' 46 useUzawa='TRUE' 47 48 #define mesh 49 H=2.0 50 L=1.0 51 W=1.0 52 mesh = esys.finley.Rectangle(l0=L, l1=H, order=-1, n0=20, n1=20, useElementsOnFace=0) # use linear macro elements for pressure 53 coordinates = mesh.getX() 54 55 #gravitational force 56 Y=Vector(0.0, Function(mesh)) 57 Y[1]=-rho*g 58 59 #element spacing 60 h=Lsup(mesh.getSize()) 61 62 #boundary conditions for slip at base 63 boundary_cond=whereZero(coordinates[1])*[0.0,1.0]+whereZero(coordinates[0])*[1.0,0.0] 64 65 #velocity and pressure vectors 66 velocity=Vector(0.0, Solution(mesh)) 67 pressure=Scalar(0.0, ReducedSolution(mesh)) 68 69 #Stokes Cartesian 70 solution=StokesProblemCartesian(mesh) 71 solution.setTolerance(tolerance) 72 73 while t <= t_end: 74 75 print(" ----- Time step = %s -----"%( t )) 76 print("Time = %s seconds"%( time )) 77 78 solution.initialize(fixed_u_mask=boundary_cond,eta=eta,f=Y) 79 velocity,pressure=solution.solve(velocity,pressure,max_iter=max_iter,verbose=verbose,usePCG=True) 80 81 print("Max velocity =", Lsup(velocity), "m/s") 82 83 #Courant condition 84 dt=0.4*h/(Lsup(velocity)) 85 print("dt", dt) 86 87 #displace the mesh 88 displacement = velocity * dt 89 coordinates = mesh.getX() 90 newx=interpolate(coordinates + displacement, ContinuousFunction(mesh)) 91 mesh.setX(newx) 92 93 time += dt 94 95 vel_mag = length(velocity) 96 97 #save velocity and pressure output 98 saveVTK("vel.%2.2i.vtu"%(t),vel=vel_mag,vec=velocity,pressure=pressure) 99 t = t+1.0

Properties

Name Value
svn:mergeinfo /branches/lapack2681/doc/examples/usersguide/fluid.py:2682-2741 /branches/pasowrap/doc/examples/usersguide/fluid.py:3661-3674 /branches/py3_attempt2/doc/examples/usersguide/fluid.py:3871-3891 /branches/restext/doc/examples/usersguide/fluid.py:2610-2624 /branches/ripleygmg_from_3668/doc/examples/usersguide/fluid.py:3669-3791 /branches/symbolic_from_3470/doc/examples/usersguide/fluid.py:3471-3974 /branches/symbolic_from_3470/ripley/test/python/doc/examples/usersguide/fluid.py:3517-3974 /trunk/ripley/test/python/doc/examples/usersguide/fluid.py:3480-3515