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

Revision 3901 - (show annotations)
Wed May 23 08:09:45 2012 UTC (7 years, 4 months ago) by jfenwick
File MIME type: text/x-python
File size: 2175 byte(s)
```Rectangle and Brick use rich face element by default.
This is to address issue #653

Updated unit tests to explicitly not use them where they didn't
specify previoulsy

```
 1 ######################################################## 2 # 3 # Copyright (c) 2008-2010 by University of Queensland 4 # Earth Systems Science Computational Center (ESSCC) 5 6 # 7 # Primary Business: Queensland, Australia 8 # Licensed under the Open Software License version 3.0 9 10 # 11 ######################################################## 12 13 ######## August 2008 ######## 14 ########## Leon Graham ########## 15 ## Newtonian fluid using StokesProblemCartesian class## 16 17 from esys.escript import * 18 import esys.finley 19 from esys.escript.linearPDEs import LinearPDE 20 from esys.escript.models import StokesProblemCartesian 21 from esys.weipa import saveVTK 22 23 #physical constants 24 eta=1.0 25 rho=100.0 26 g=10.0 27 28 #solver settings 29 tolerance=1.0e-4 30 max_iter=200 31 t_end=50 32 t=0.0 33 time=0 34 verbose='TRUE' 35 useUzawa='TRUE' 36 37 #define mesh 38 H=2.0 39 L=1.0 40 W=1.0 41 mesh = esys.finley.Rectangle(l0=L, l1=H, order=-1, n0=20, n1=20, useElementsOnFace=0) # use linear macro elements for pressure 42 coordinates = mesh.getX() 43 44 #gravitational force 45 Y=Vector(0.0, Function(mesh)) 46 Y[1]=-rho*g 47 48 #element spacing 49 h=Lsup(mesh.getSize()) 50 51 #boundary conditions for slip at base 52 boundary_cond=whereZero(coordinates[1])*[0.0,1.0]+whereZero(coordinates[0])*[1.0,0.0] 53 54 #velocity and pressure vectors 55 velocity=Vector(0.0, Solution(mesh)) 56 pressure=Scalar(0.0, ReducedSolution(mesh)) 57 58 #Stokes Cartesian 59 solution=StokesProblemCartesian(mesh) 60 solution.setTolerance(tolerance) 61 62 while t <= t_end: 63 64 print(" ----- Time step = %s -----"%( t )) 65 print("Time = %s seconds"%( time )) 66 67 solution.initialize(fixed_u_mask=boundary_cond,eta=eta,f=Y) 68 velocity,pressure=solution.solve(velocity,pressure,max_iter=max_iter,verbose=verbose,usePCG=True) 69 70 print("Max velocity =", Lsup(velocity), "m/s") 71 72 #Courant condition 73 dt=0.4*h/(Lsup(velocity)) 74 print("dt", dt) 75 76 #displace the mesh 77 displacement = velocity * dt 78 coordinates = mesh.getX() 79 mesh.setX(coordinates + displacement) 80 81 time += dt 82 83 vel_mag = length(velocity) 84 85 #save velocity and pressure output 86 saveVTK("vel.%2.2i.vtu"%(t),vel=vel_mag,vec=velocity,pressure=pressure) 87 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