TALK=T;RUN( 1, 1)
DISPLAY
This case describes how the 'moving grid' method (see MOFOR
entry in Encyclopedia, POLIS) is applied to simulation
of a sphere falling through air. The 'moving grid' method is
implemented by In-Form formulae written in the Q1 file.
A sphere falls through air experiences a force in direction
opposite to its motion. Terminal velocity is achieved when
the drag force is equal in magnitude but opposite in direction
to the gravity force propelling the sphere.
In this case, the quadratic drag equation is used to calculate
the drag force as follows.
Fd = - 0.5 * rho * U**2.0 * A * Cd
where Fd is the force of drag
rho is the density of air
U is the speed of the sphere
A is the area of the projection of the sphere on a plane
perpendicular to the direction of motion.
Cd is the drag coefficient
The velocity as a function of time for the sphere can be derived
as follows:
U(t)= A1 * (EXP(A2*t)-1)/(EXP(A2*t)+1)
where A1 = (2*M*G/rho*A*Cd)**0.5)
A2 = ((2*rho*A*Cd*G)/M)**0.5
where G is the gravity,9.81
M is the mass of the sphere
The following terminal velcoty is asymptotically approached:
Utm = (2*M*G/rho*A*Cd)**0.5
The accelaration is given by the derivative of the velocity
function with respect to time,
a = 4*A1**2.0*A3*EXP(A2*t)/(EXP(A2*t)+1)**2.0
where A3 is 0.5*rho*A*Cd/M
The fluid used is air; the radius of the sphere is 1.0 and
the drag coefficient, Cd used is 1.0
This file contains a marco of commands which cause the
Viewer when the Macro button is pressed to display the
animation automatically. Note that the probe value of the
last contour indicates the terminal velocity.
ENDDIS
VRV USE
DOMAIN ON
* Setting object visibility and painting status
OBJECT SHOW TYPE BLOCKAGE
OBJECT PAINT TYPE BLOCKAGE OFF
OBJECT WIREFRAME TYPE BLOCKAGE OFF
OBJECT SHOW TYPE OUTLET
OBJECT PAINT TYPE OUTLET OFF
OBJECT WIREFRAME TYPE OUTLET OFF
VARIABLE Pressure; CON ON
VECTOR ON
ANIMATE 1 20 1
VARIABLE Velocity;vector on
PAUSE
ENDUSE
************************************************************
Q1 created by VDI menu, Version 2007, Date 23/11/07
CPVNAM=VDI;SPPNAM=Core
************************************************************
Echo DISPLAY / USE settings
************************************************************
IRUNN = 1 ;LIBREF = 805
************************************************************
Group 1. Run Title
TEXT(Falling sphere with air resistance )
************************************************************
Group 2. Transience
STEADY=F
* Set overall time and no. of steps
RSET(U,0.000000E+00,4.000000E+00,80)
* Modify regions
************************************************************
Groups 3, 4, 5 Grid Information
* Overall number of cells, RSET(M,NX,NY,NZ,tolerance)
RSET(M,18,18,31)
************************************************************
Group 6. Body-Fitted coordinates
************************************************************
Group 7. Variables: STOREd,SOLVEd,NAMEd
ONEPHS = T
* Non-default variable names
NAME(149) =VLSQ ; NAME(150) =PRPS
* Solved variables list
SOLVE(P1 ,U1 ,V1 ,W1 )
* Stored variables list
STORE(PRPS,VLSQ)
* Additional solver options
SOLUTN(P1 ,Y,Y,Y,N,N,Y)
SOLUTN(V1 ,Y,Y,Y,N,N,Y)
SOLUTN(W1 ,Y,Y,Y,N,N,Y)
************************************************************
Group 8. Terms & Devices
************************************************************
Group 9. Properties
RHO1 = 1.000000E+00
ENUL = 1.000000E-05
CP1 = 1.005000E+03
ENUT = 0.000000E+00
************************************************************
Group 10.Inter-Phase Transfer Processes
************************************************************
Group 11.Initialise Var/Porosity Fields
FIINIT(PRPS) = -1.000000E+00
No PATCHes used for this Group
INIADD = F
************************************************************
Group 12. Convection and diffusion adjustments
No PATCHes used for this Group
************************************************************
Group 13. Boundary & Special Sources
No PATCHes used for this Group
EGWF = T
************************************************************
Echo InForm settings for Group 13
SAVE13BEGIN
real(aa1,aa2,aa3,aradi,vol1,area1,acd)
Radius of the sphere
aradi=1.0
Drag coefficient, Cd
acd=1.0
the cross area of the sphere
area1=3.14159*aradi**2.0
the volume of the sphere
vol1=4.0*3.14159*aradi**3/3.0
A1=(2*9.81*Vol/(area*cd))**0.5
aa1=(2.0*9.81*vol1/area1/acd)**0.5
A2=(2*9.81*A*Cd)/Vol)**0.5
aa2=(2.0*9.81*area1*acd/vol1)**0.5
A3=0.5*A*Cd/Vol
aa3=0.5*area1*acd/vol1
Uin = aa1*(exp(aa2*tim)-1)/(exp(aa2*tim)+1)
patch(in,low,1,nx,1,ny,1,1,1,lstep)
(source of p1 at in is :aa1:*(exp(:aa2:*tim)-1)/(exp(:aa2:*tim)+1)*$
rho1)
(source of w1 at in is :aa1:*(exp(:aa2:*tim)-1)/(exp(:aa2:*tim)+1) $
with onlyms)
patch(acel,phasem,1,nx,1,ny,1,nz,1,lstep)
(source of w1 at acel is 4*9.81*exp(:aa2:*(tim-0.025))/(exp(:aa2:*($
tim-0.025))+1)^2.0)
Print out in the result the accelaration
(stored var E1 is 4*9.81*exp(:aa2:*(tim-0.025))/(exp(:aa2:*(tim-0.0$
25))+1)^2.0)
SAVE13END
************************************************************
Group 14. Downstream Pressure For PARAB
************************************************************
Group 15. Terminate Sweeps
LSWEEP = 20
RESFAC = 1.000000E-03
************************************************************
Group 16. Terminate Iterations
LITER (P1 ) = 250
************************************************************
Group 17. Relaxation
RELAX(P1 ,LINRLX, 1.000000E+00)
************************************************************
Group 18. Limits
VARMAX(U1 ) = 1.000000E+06 ;VARMIN(U1 ) =-1.000000E+06
VARMAX(V1 ) = 1.000000E+06 ;VARMIN(V1 ) =-1.000000E+06
VARMAX(W1 ) = 1.000000E+06 ;VARMIN(W1 ) =-1.000000E+06
************************************************************
Group 19. EARTH Calls To GROUND Station
USEGRD = T ;USEGRX = T
NAMSAT =MOSG
ASAP = T
PARSOL = T
CALFOR = T
IDISPB = 1 ;IDISPC = 100
RG( 1) = 2.000000E+00
SPEDAT(SET,GXMONI,TRANSIENT,L,F)
************************************************************
Group 20. Preliminary Printout
ECHO = T
DISTIL = T ;NULLPR = F
NDST = 0
DSTTOL = 1.000000E-02
************************************************************
Group 21. Print-out of Variables
OUTPUT(VLSQ,N,N,Y,N,N,N)
OUTPUT(PRPS,N,N,Y,N,N,N)
************************************************************
Group 22. Monitor Print-Out
IXMON = 9 ;IYMON = 9 ;IZMON = 1
NPRMON = 100000
NPRMNT = 1
TSTSWP = -1
************************************************************
Group 23.Field Print-Out & Plot Control
NPRINT = 100000
NTPRIN = -1 ;ISTPRF = 1 ;ISTPRL = 100000
YZPR = T
ISWPRF = 1 ;ISWPRL = 100000
No PATCHes used for this Group
************************************************************
Group 24. Dumps For Restarts
IDISPA = 2 ;IDISPB = 1 ;IDISPC = 100
CSG1 ='M'
EX(P1)= 1.545E+00;EX(U1)= 1.283E-01;EX(V1)= 1.283E-01
EX(W1)= 5.154E+00;EX(E1)= 9.243E-06
GVIEW(P,-1.000000E+00,0.000000E+00,0.000000E+00)
GVIEW(UP,0.000000E+00,1.000000E+00,0.000000E+00)
> DOM, SIZE, 4.000000E+00, 4.000000E+00, 1.000000E+01
> DOM, MONIT, 1.777778E+00, 1.777778E+00, 1.538462E-01
> DOM, SCALE, 1.000000E+00, 1.000000E+00, 1.000000E+00
> GRID, RSET_X_1, 5, 1.000000E+00
> GRID, RSET_X_2, 9, 1.000000E+00
> GRID, RSET_X_3, 4, 1.000000E+00
> GRID, RSET_Y_1, 5, 1.000000E+00
> GRID, RSET_Y_2, 9, 1.000000E+00
> GRID, RSET_Y_3, 4, 1.000000E+00
> GRID, RSET_Z_1, 13, 1.000000E+00
> GRID, RSET_Z_2, 6, 1.000000E+00
> GRID, RSET_Z_3, 12, 1.000000E+00
> DOM, MOMCEN, 0.000000E+00, 0.000000E+00, 0.000000E+00
> OBJ, NAME, B2
> OBJ, POSITION, 0.000000E+00, 0.000000E+00, 1.000000E+01
> OBJ, SIZE, 4.000000E+00, 4.000000E+00, 0.000000E+00
> OBJ, ROT-MODE, OLD
> OBJ, GEOMETRY, cube12
> OBJ, ROTATION24, 1
> OBJ, TYPE, OUTLET
> OBJ, PRESSURE, 0.000000E+00
> OBJ, TEMPERATURE, 0.000000E+00
> OBJ, COEFFICIENT, 1.000000E+03
> OBJ, TIME_LIMITS, ALWAYS_ACTIVE
> OBJ, NAME, B3
> OBJ, POSITION, 1.000000E+00, 1.000000E+00, 4.000000E+00
> OBJ, SIZE, 2.000000E+00, 2.000000E+00, 2.000000E+00
> OBJ, GEOMETRY, sphere
> OBJ, ROTATION24, 1
> OBJ, TYPE, BLOCKAGE
> OBJ, MATERIAL, 198,Solid with smooth-wall friction
> OBJ, TIME_LIMITS, ALWAYS_ACTIVE
STOP