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

Revision 3975 - (show annotations)
Thu Sep 20 01:54:06 2012 UTC (7 years, 1 month ago) by caltinay
File MIME type: text/x-python
File size: 2226 byte(s)
```Merged symbolic branch into trunk. Curious what daniel and spartacus have to
say...

```
 1 ######################################################## 2 # 3 # Copyright (c) 2008-2012 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 newx=interpolate(coordinates + displacement, ContinuousFunction(mesh)) 80 mesh.setX(newx) 81 82 time += dt 83 84 vel_mag = length(velocity) 85 86 #save velocity and pressure output 87 saveVTK("vel.%2.2i.vtu"%(t),vel=vel_mag,vec=velocity,pressure=pressure) 88 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