CFD Online Discussion Forums

CFD Online Discussion Forums (https://www.cfd-online.com/Forums/)
-   OpenFOAM Running, Solving & CFD (https://www.cfd-online.com/Forums/openfoam-solving/)
-   -   Convergence moving mesh (https://www.cfd-online.com/Forums/openfoam-solving/59322-convergence-moving-mesh.html)

lr103476 February 13, 2006 07:43

I have created a moving mesh s
 
I have created a moving mesh solver by combining icoFoamAutoMotion and interFoam. The 2d motion of the cylinder is specified by movingFlapFvMesh, xampl=1.0, xfreq=0.1, yampl=1.0, yfreq=0.1 and thetaampl=0.0. The inlet velocity it 1 m/s and corresponds to Re=150.

The timestep is 0.005. Up to time=1.225 nice covergence is obtained, but at time=1.225 the solution begins to diverge, as can be seen from the Courant numbers shown below. Furthermore the number of iterations of solving motionUx and motionUy is increasing every iteration. Therefore I think that the motionSolver causes this problem, but I have no clue why. Any suggestions?

Thanks.

Greetings, Frank

==================================================
Output
==================================================
Mean and max Courant Numbers = 0.023978 0.948195
Time = 1.225

Correct mesh motion diffusion field.
ICCG: Solving for motionUx, Initial residual = 0.00074787, Final residual = 9.94298e-06, No Iterations 280
ICCG: Solving for motionUy, Initial residual = 0.000497424, Final residual = 9.96133e-06, No Iterations 187
Correcting 2-D mesh motion ...done
volume continuity errors : sum local = 4.06953e-15, global = -1.90209e-18
BICCG: Solving for Ux, Initial residual = 0.000139335, Final residual = 3.63493e-06, No Iterations 1
BICCG: Solving for Uy, Initial residual = 0.00080448, Final residual = 4.75049e-07, No Iterations 2
ICCG: Solving for p, Initial residual = 0.84611, Final residual = 9.43944e-07, No Iterations 290
time step continuity errors : sum local = 4.11744e-11, global = 1.12633e-14, cumulative = 6.75289e-10
ICCG: Solving for p, Initial residual = 0.859509, Final residual = 9.95797e-07, No Iterations 290
time step continuity errors : sum local = 4.79839e-11, global = 6.3742e-15, cumulative = 6.75296e-10

Wall patch = 3
Wall patch name = cyl_wall
Uav = (0.282238 -0.717762 5.26235e-28)
Aref = 1.44219
Lref = 0.564692
DragCoefficient = 5.82549
pressureDragCoefficient = 4.83947
viscDragCoefficient = 0.986017
LiftCoefficient = 2.19606
MomentCoefficient = 8.31774

ExecutionTime = 9948.97 s


Mean and max Courant Numbers = 0.0239854 1.13739
Time = 1.23

Correct mesh motion diffusion field.
ICCG: Solving for motionUx, Initial residual = 0.000752869, Final residual = 9.99251e-06, No Iterations 281
ICCG: Solving for motionUy, Initial residual = 0.000500919, Final residual = 9.95567e-06, No Iterations 187
Correcting 2-D mesh motion ...done
volume continuity errors : sum local = 4.09456e-15, global = 2.59593e-18
BICCG: Solving for Ux, Initial residual = 0.000173378, Final residual = 4.67296e-06, No Iterations 1
BICCG: Solving for Uy, Initial residual = 0.00112283, Final residual = 7.11956e-07, No Iterations 2
ICCG: Solving for p, Initial residual = 0.890628, Final residual = 9.45922e-07, No Iterations 295
time step continuity errors : sum local = 5.84773e-11, global = -7.53318e-14, cumulative = 6.7522e-10
ICCG: Solving for p, Initial residual = 0.901723, Final residual = 9.32168e-07, No Iterations 296
time step continuity errors : sum local = 6.38908e-11, global = 1.57041e-13, cumulative = 6.75377e-10

Wall patch = 3
Wall patch name = cyl_wall
Uav = (0.284429 -0.715571 5.79729e-29)
Aref = 1.44221
Lref = 0.559905
DragCoefficient = 6.26111
pressureDragCoefficient = 5.26338
viscDragCoefficient = 0.997727
LiftCoefficient = 2.30624
MomentCoefficient = 10.1442

ExecutionTime = 10008.3 s


Mean and max Courant Numbers = 0.0240167 1.89597
Time = 1.235

Correct mesh motion diffusion field.
ICCG: Solving for motionUx, Initial residual = 0.000757996, Final residual = 9.99529e-06, No Iterations 281
ICCG: Solving for motionUy, Initial residual = 0.000504364, Final residual = 9.90425e-06, No Iterations 188
Correcting 2-D mesh motion ...done
volume continuity errors : sum local = 4.10338e-15, global = 1.22539e-18
BICCG: Solving for Ux, Initial residual = 0.000229742, Final residual = 6.51542e-06, No Iterations 1
BICCG: Solving for Uy, Initial residual = 0.00160232, Final residual = 1.13808e-06, No Iterations 2
ICCG: Solving for p, Initial residual = 0.923514, Final residual = 9.85298e-07, No Iterations 299
time step continuity errors : sum local = 8.725e-11, global = -8.65098e-13, cumulative = 6.74512e-10
ICCG: Solving for p, Initial residual = 0.932281, Final residual = 9.66381e-07, No Iterations 301
time step continuity errors : sum local = 9.52416e-11, global = 8.27088e-13, cumulative = 6.75339e-10

Wall patch = 3
Wall patch name = cyl_wall
Uav = (0.286627 -0.713373 -7.71628e-28)
Aref = 1.44222
Lref = 0.55508
DragCoefficient = 7.31698
pressureDragCoefficient = 6.28853
viscDragCoefficient = 1.02845
LiftCoefficient = 2.55426
MomentCoefficient = 14.3745

ExecutionTime = 10067.6 s


Mean and max Courant Numbers = 0.0240893 3.38602
Time = 1.24

Correct mesh motion diffusion field.
ICCG: Solving for motionUx, Initial residual = 0.000763107, Final residual = 9.99048e-06, No Iterations 280
ICCG: Solving for motionUy, Initial residual = 0.000507756, Final residual = 9.93261e-06, No Iterations 188
Correcting 2-D mesh motion ...done
volume continuity errors : sum local = 4.08207e-15, global = -1.90597e-18
BICCG: Solving for Ux, Initial residual = 0.000328665, Final residual = 1.3123e-06, No Iterations 2
BICCG: Solving for Uy, Initial residual = 0.00234532, Final residual = 2.49132e-06, No Iterations 2
ICCG: Solving for p, Initial residual = 0.945383, Final residual = 9.50596e-07, No Iterations 306
time step continuity errors : sum local = 1.20953e-10, global = -1.79658e-13, cumulative = 6.7516e-10
ICCG: Solving for p, Initial residual = 0.953924, Final residual = 9.49527e-07, No Iterations 308
time step continuity errors : sum local = 1.33918e-10, global = 5.25042e-14, cumulative = 6.75212e-10

Wall patch = 3
Wall patch name = cyl_wall
Uav = (0.288832 -0.711168 5.26306e-28)
Aref = 1.44221
Lref = 0.550215
DragCoefficient = 9.85785
pressureDragCoefficient = 8.74477
viscDragCoefficient = 1.11308
LiftCoefficient = 3.1429
MomentCoefficient = 24.2969

ExecutionTime = 10126.9 s


Mean and max Courant Numbers = 0.0242311 6.68236
Time = 1.245

Correct mesh motion diffusion field.
ICCG: Solving for motionUx, Initial residual = 0.000768247, Final residual = 9.8493e-06, No Iterations 282
ICCG: Solving for motionUy, Initial residual = 0.000511277, Final residual = 9.82936e-06, No Iterations 193
Correcting 2-D mesh motion ...done
volume continuity errors : sum local = 4.09053e-15, global = 8.99978e-19
BICCG: Solving for Ux, Initial residual = 0.000517687, Final residual = 2.34858e-06, No Iterations 2
BICCG: Solving for Uy, Initial residual = 0.00348416, Final residual = 5.18035e-06, No Iterations 3
ICCG: Solving for p, Initial residual = 0.955581, Final residual = 9.54336e-07, No Iterations 312
time step continuity errors : sum local = 1.69908e-10, global = 1.25309e-12, cumulative = 6.76465e-10
ICCG: Solving for p, Initial residual = 0.96498, Final residual = 9.80104e-07, No Iterations 314
time step continuity errors : sum local = 1.91259e-10, global = -1.58831e-12, cumulative = 6.74877e-10

Wall patch = 3
Wall patch name = cyl_wall
Uav = (0.291044 -0.708956 7.71984e-28)
Aref = 1.44217
Lref = 0.545312
DragCoefficient = 15.0292
pressureDragCoefficient = 13.489
viscDragCoefficient = 1.54019
LiftCoefficient = 4.10684
MomentCoefficient = 48.2943

ExecutionTime = 10187 s


Mean and max Courant Numbers = 0.0244966 14.8136
Time = 1.25

Correct mesh motion diffusion field.
ICCG: Solving for motionUx, Initial residual = 0.000773288, Final residual = 9.89781e-06, No Iterations 282
ICCG: Solving for motionUy, Initial residual = 0.000514758, Final residual = 9.9967e-06, No Iterations 188
Correcting 2-D mesh motion ...done
volume continuity errors : sum local = 4.04741e-15, global = 8.71871e-19
BICCG: Solving for Ux, Initial residual = 0.000889756, Final residual = 2.69979e-06, No Iterations 3
BICCG: Solving for Uy, Initial residual = 0.00519126, Final residual = 1.8881e-06, No Iterations 5
ICCG: Solving for p, Initial residual = 0.950479, Final residual = 9.4596e-07, No Iterations 318
time step continuity errors : sum local = 2.22139e-10, global = 1.45724e-12, cumulative = 6.76334e-10
ICCG: Solving for p, Initial residual = 0.957963, Final residual = 8.6076e-07, No Iterations 320
time step continuity errors : sum local = 2.21044e-10, global = -1.29554e-12, cumulative = 6.75039e-10

Wall patch = 3
Wall patch name = cyl_wall
Uav = (0.293263 -0.706737 -1.80998e-28)
Aref = 1.44214
Lref = 0.54037
DragCoefficient = 26.3516
pressureDragCoefficient = 23.0254
viscDragCoefficient = 3.32619
LiftCoefficient = 6.17307
MomentCoefficient = 104.607

ExecutionTime = 10246.6 s


Mean and max Courant Numbers = 0.0249365 40.0673
Time = 1.255

Correct mesh motion diffusion field.
ICCG: Solving for motionUx, Initial residual = 0.000778549, Final residual = 9.97842e-06, No Iterations 282
ICCG: Solving for motionUy, Initial residual = 0.000518411, Final residual = 9.89038e-06, No Iterations 193
Correcting 2-D mesh motion ...done
volume continuity errors : sum local = 4.08214e-15, global = 1.5051e-18
BICCG: Solving for Ux, Initial residual = 0.00178422, Final residual = 3.98896e-06, No Iterations 7
BICCG: Solving for Uy, Initial residual = 0.00891935, Final residual = 7.66749e-06, No Iterations 8
ICCG: Solving for p, Initial residual = 0.924408, Final residual = 9.65802e-07, No Iterations 322
time step continuity errors : sum local = 2.86659e-10, global = 4.15199e-13, cumulative = 6.75454e-10
ICCG: Solving for p, Initial residual = 0.940882, Final residual = 9.8166e-07, No Iterations 322
time step continuity errors : sum local = 3.33021e-10, global = -5.05583e-13, cumulative = 6.74948e-10

Wall patch = 3
Wall patch name = cyl_wall
Uav = (0.295489 -0.704511 2.74279e-28)
Aref = 1.44219
Lref = 0.53539
DragCoefficient = 61.6947
pressureDragCoefficient = 52.343
viscDragCoefficient = 9.35167
LiftCoefficient = 18.4061
MomentCoefficient = 201.47

ExecutionTime = 10307 s


Mean and max Courant Numbers = 0.0259063 127.598
Time = 1.26

Correct mesh motion diffusion field.
ICCG: Solving for motionUx, Initial residual = 0.000783877, Final residual = 9.81292e-06, No Iterations 284
ICCG: Solving for motionUy, Initial residual = 0.000521931, Final residual = 9.95579e-06, No Iterations 194
Correcting 2-D mesh motion ...done
volume continuity errors : sum local = 4.0859e-15, global = -2.14512e-18
BICCG: Solving for Ux, Initial residual = 0.00572169, Final residual = 7.68117e-06, No Iterations 12
BICCG: Solving for Uy, Initial residual = 0.0344701, Final residual = 7.918e-06, No Iterations 14
ICCG: Solving for p, Initial residual = 0.863708, Final residual = 9.84167e-07, No Iterations 324
time step continuity errors : sum local = 3.85863e-10, global = -8.52682e-13, cumulative = 6.74096e-10
ICCG: Solving for p, Initial residual = 0.89722, Final residual = 9.61545e-07, No Iterations 323
time step continuity errors : sum local = 4.96108e-10, global = -1.1644e-12, cumulative = 6.72931e-10

Wall patch = 3
Wall patch name = cyl_wall
Uav = (0.297722 -0.702278 4.16906e-28)
Aref = 1.44221
Lref = 0.53037
DragCoefficient = 194.264
pressureDragCoefficient = 178.762
viscDragCoefficient = 15.5021
LiftCoefficient = 68.4224
MomentCoefficient = 491.493

ExecutionTime = 10368 s


Mean and max Courant Numbers = 0.0287663 420.85
Time = 1.265

Correct mesh motion diffusion field.
ICCG: Solving for motionUx, Initial residual = 0.000788983, Final residual = 9.7713e-06, No Iterations 285
ICCG: Solving for motionUy, Initial residual = 0.00052554, Final residual = 9.99632e-06, No Iterations 194
Correcting 2-D mesh motion ...done
volume continuity errors : sum local = 4.07999e-15, global = 4.6044e-19
BICCG: Solving for Ux, Initial residual = 0.0340294, Final residual = 3.50395e-06, No Iterations 20
BICCG: Solving for Uy, Initial residual = 0.290113, Final residual = 1.78799e-06, No Iterations 22
ICCG: Solving for p, Initial residual = 0.691357, Final residual = 9.72633e-07, No Iterations 380
time step continuity errors : sum local = 1.24537e-09, global = 6.79261e-12, cumulative = 6.79724e-10
ICCG: Solving for p, Initial residual = 0.997215, Final residual = 8.99401e-07, No Iterations 393
time step continuity errors : sum local = 1.00642e-07, global = 4.05372e-10, cumulative = 1.0851e-09

Wall patch = 3
Wall patch name = cyl_wall
Uav = (0.299962 -0.700038 2.34445e-28)
Aref = 1.44222
Lref = 0.525311
DragCoefficient = 3812.82
pressureDragCoefficient = 916.412
viscDragCoefficient = 2896.41
LiftCoefficient = 9426.75
MomentCoefficient = -102740

ExecutionTime = 10430 s


Mean and max Courant Numbers = 0.571589 111497
Time = 1.27

Correct mesh motion diffusion field.
ICCG: Solving for motionUx, Initial residual = 0.000794283, Final residual = 9.78199e-06, No Iterations 285
ICCG: Solving for motionUy, Initial residual = 0.000529197, Final residual = 9.90328e-06, No Iterations 196
Correcting 2-D mesh motion ...done
volume continuity errors : sum local = 4.07254e-15, global = 1.28338e-19
BICCG: Solving for Ux, Initial residual = 0.999564, Final residual = 3.97533e-06, No Iterations 164
BICCG: Solving for Uy, Initial residual = 0.993602, Final residual = 7.81896e-06, No Iterations 160
ICCG: Solving for p, Initial residual = 0.45766, Final residual = 6.70541e-07, No Iterations 664
time step continuity errors : sum local = 1.8692e-08, global = -3.35061e-11, cumulative = 1.05159e-09
ICCG: Solving for p, Initial residual = 0.998543, Final residual = 9.70955e-07, No Iterations 750
time step continuity errors : sum local = 1.49225e-05, global = -3.6547e-08, cumulative = -3.54954e-08

Wall patch = 3
Wall patch name = cyl_wall
Uav = (0.302209 -0.697791 -3.63624e-28)
Aref = 1.4422
Lref = 0.520213
DragCoefficient = -2.70433e+07
pressureDragCoefficient = -2.69181e+07
viscDragCoefficient = -125286
LiftCoefficient = -6.50234e+06
MomentCoefficient = -1.03616e+08

ExecutionTime = 10500.2 s


Mean and max Courant Numbers = 426.603 3.37259e+07
Time = 1.275

Correct mesh motion diffusion field.
ICCG: Solving for motionUx, Initial residual = 0.000799647, Final residual = 9.80335e-06, No Iterations 285
ICCG: Solving for motionUy, Initial residual = 0.00053282, Final residual = 9.96991e-06, No Iterations 196
Correcting 2-D mesh motion ...done
volume continuity errors : sum local = 4.08114e-15, global = -1.70846e-18
BICCG: Solving for Ux, Initial residual = 0.999992, Final residual = 3.38725, No Iterations 1001
BICCG: Solving for Uy, Initial residual = 0.999993, Final residual = 1.81754, No Iterations 1001
ICCG: Solving for p, Initial residual = 1, Final residual = 0.423881, No Iterations 5001
time step continuity errors : sum local = 5.22248e+06, global = -2.21042, cumulative = -2.21042
ICCG: Solving for p, Initial residual = 0.999492, Final residual = 10.9301, No Iterations 5001
time step continuity errors : sum local = 4.02462e+11, global = 272.815, cumulative = 270.604

Wall patch = 3
Wall patch name = cyl_wall
Uav = (0.304463 -0.695537 7.76322e-28)
Aref = 1.44217
Lref = 0.515075
DragCoefficient = -4.22157e+16
pressureDragCoefficient = -4.06296e+16
viscDragCoefficient = -1.58615e+15
LiftCoefficient = 1.55157e+16
MomentCoefficient = -5.22663e+17

liu February 13, 2006 12:36

In my experience, interFoam is
 
In my experience, interFoam is not ready for open channel flow, yet.

lr103476 February 13, 2006 12:44

I use the icoFoamAutoMotion eq
 
I use the icoFoamAutoMotion equations and related algorithms, piso loop and stuff. The moving mesh part is copied from interFoam, since that solver works with the movingFvMesh class.

Regards,
Frank

guido_adriaensen February 14, 2006 02:14

Hoi Frank, Je piso-loop con
 
Hoi Frank,

Je piso-loop convergeert niet. Als je kijkt naar de initial residual van p, moet deze minimaal 1 orde lager zijn bij de tweede iteratie in de loop, dit is niet het geval (Hier is die waarde zelfs iets hoger):

ICCG: Solving for p, Initial residual = 0.84611, Final residual = 9.43944e-07, No Iterations 290
time step continuity errors : sum local = 4.11744e-11, global = 1.12633e-14, cumulative = 6.75289e-10
ICCG: Solving for p, Initial residual = 0.859509, Final residual = 9.95797e-07, No Iterations 290

Je kunt proberen om het aantal piso correctors op te voeren of een kleinere tijdsstap nemen. Succes

Guido

lr103476 February 17, 2006 05:46

Ok, more piso correctors is no
 
Ok, more piso correctors is not improving the convergence.

On the other hand, a much smaller timestep is actually leading to convergence at the beginning of the motion. After a certain time (which is larger than the previous case) also this solution is diverges.

When I turn off the mesh.move() function in my solver no problem occurs. So the problem lies in the piso loop when the mesh moves.

Thanks, Frank

lr103476 February 21, 2006 06:08

Could someone please comment o
 
Could someone please comment on my solver, shown above.

When I move the mesh, the piso loop diverges, also when I decrease the timestep or increase the number of correctors.

When I turn off the mesh motion, the flow solver works fine. In addition when I turn off the flow solver instead of the mesh motion, the moving mesh is calculated well.

Do I have to include some extra terms in the piso loop to account for the mesh motion?

Thanks in advance.

Regards, Frank

hjasak February 21, 2006 11:10

The solver looks fine (well, a
 
The solver looks fine (well, at least not obviously wrong). However, your residual behaviour indicates a serious problem and my guess is that you have messed up the boundary conditions.

Note that the wall boundary needs to be specified in terms of absolute velocity, so a moving wall with no flux through it will not have fixedValue (0 0 0). If this is the problem use the movingWallVelocity b.c., which corrects the normal component based on the mesh motion.

As a test, have a look at the fluxes through your wall boundaries - if all is well, they should be exactly zero.

Hrv

lr103476 February 22, 2006 10:21

As far as I know I already use
 
As far as I know I already used the movingWallVelocity b.c.

My boundary conditions in 0/U looks like (cyl_wall is the moving surface):

dimensions [0 1 -1 0 0 0 0];
internalField uniform (1 0 0);
boundaryField
{
far_field
{
type symmetryPlane;
}
outlet
{
type zeroGradient;
}
inlet
{
type fixedValue;
value uniform (1 0 0);
}
cyl_wall
{
type movingWallVelocity;
value uniform (0 0 0);
}
frontAndBackPlanes
{
type empty;
}
}

===========================================
And 0/motionU looks like:

dimensions [0 1 -1 0 0 0 0];
internalField uniform (0 0 0);
boundaryField
{
far_field
{
type symmetryPlane;
}
outlet
{
type fixedValue;
value uniform (0 0 0);
}
inlet
{
type fixedValue;
value uniform (0 0 0);
}
cyl_wall
{
type fixedValue;
value uniform (0 0 0);
}
frontAndBackPlanes
{
type empty;
}
}

===========================================
When I used movingFlepFvMesh the phi of some surface elements after 1 sec looks like:

1/phi:

cyl_wall
{
type calculated;
value nonuniform List<scalar>
280
(
0
1.73472e-18
0
1.73472e-18
-1.73472e-18
0
0
0
0
-1.73472e-18
0
0
0
1.73472e-18
0
0
8.67362e-19
-8.67362e-19
0
-8.67362e-19
.
.
.
.
.
}

What is wrong with this b.c. implementation. I think I mess thing up with U and motionU........The fixedValue of motionU is overruled by movingFlapFvMesh, right?

Regards, Frank

lr103476 February 23, 2006 07:09

Summarising: What b.c. do y
 
Summarising:

What b.c. do you generally have to apply to a moving body?

* motionU: fixedValue on the moving surface (this is being overruled by the movingFvMesh class)

* U: movingWallVelocity (uniform (0,0,0)) on the moving surface (according to dr Jasak this will correct the boundary flux)

Is this correct? or must by movingWallVelocity be applied to both U and motionU?

If these settings are correct, my problem of the diverging piso loop has its origins elsewhere.

Thanks,
Frank

hjasak February 23, 2006 09:29

motionU: you can do whatever y
 
motionU: you can do whatever you like, this is only used to specify the motion. Typically, if you wish to specify the motion of an object, you would use the fixedValue b.c. on its surface (makes sense).

U: movingWallVelocity. I usually stick to this. It will correct the normal component of the velocity to achieve zero flux through the boundary based on the vertex motion (wherever it may come from).

There's not connection between the motionU b.c. and the U b.c. - they are totally independent.

I'm not sure what may be going wrong. Try setting very diffusive discretisation and if that still fails you've messed up the boundary conditions or something even more subtly. Have you tried looking at the mesh motion (motion only, no flow), to see if it behaves the way you think it should?

This should work without any trouble,

Hrv

lr103476 February 24, 2006 07:29

Thanks for your answers, but m
 
Thanks for your answers, but my set-up is still not working. I tried central, upwind and blended schemes for the convective flux, but the solution still diverges after 1 second or so.

When I read my ../1/phi file the values at the boundary are of the order 10e-18. Does this flux need to be exactly zero or is 10e-18 small enough.

I also switched off the solver such that the mesh motion can be visualised. In that case the mesh looks very nice at all time up to t=20. No diverges occurs. The solution divergences only when the piso loop is solved.

Maybe the movingFlapFvMesh coupling to the fixedValue of motionU is not done correctly....

Frank

hjasak February 24, 2006 08:56

I'll have a look if you want.
 
I'll have a look if you want. Could you please pack the code (ready to run, if the mesh is made with blockMesh, please delete unnecessary files, with time-step, discretisation, initial fields etc set up) and the case (clean up the object files and similar garbage + please include everything I need to compile it up).

Hrv

lr103476 February 24, 2006 12:30

Thanks, I just sent you th
 
Thanks,

I just sent you the case and solver files. Both are ready to run/compile.

Frank

hjasak February 27, 2006 17:50

Very easy: it does not work be
 
Very easy: it does not work because you've set the time step to give the Courant number of 3.5, which is jolly high. I have dropped the time step by a factor of 4 and with the Co number of around 1 it works beautifully for the first second or so and then the Co number goes higher again which makes it blow up. Nothing wrong with the solution, just the normal flow development on a newrrt fine mesh.

Then, I have added automatic time-step control and it works beautifully.

Hrv

lr103476 February 28, 2006 07:00

Thanks again for your time :-)
 
Thanks again for your time :-).

In my first case (dt=0.02) also the Co of around 3.5 gives me a good solution for the first second. After that the Co always blows up which is rather independent of the chosen timestep.

I really like a constant time step. Since I need to simulate for 100 seconds or so, this is easy to plan and estimate calculation times. Though, this seems to be impossible in my case.

Frank

hjasak March 4, 2006 15:22

I did the first 5 seconds and
 
I did the first 5 seconds and made a little animation for you (it is ugly but proves the point). The automatic time-step control is well-behaved and there's no problem with the simulation (I am using 1 non-orthogonal corrector). If you want the regular output dumps with automatic time-step control, e.g. every 1 second or similar, use adjustableRunTime writeControl option.

cylinder animation

Hrv

lr103476 March 6, 2006 05:54

Allright, thanks. Is there act
 
Allright, thanks. Is there actually some sort of guideline how many non-orthogonal correctors are needed? Is the accuracy or the efficiency affected by this choice, or not?

Frank

hjasak March 6, 2006 06:12

You have several options. I s
 
You have several options. I suspect the problem you are getting is that the mesh deformation may become quite sensitive to the convergence of the motion solver. The issue is that you've got a really big mesh qith a little bit of motion in the middle. Ideally, you would want to subset a piece of the mesh around the cylinder which is big enough to accomodate the motion but smaller than the complete domain. This would make yor simulation much faster (less work in mesh motion) and probably cause less problems. However, it is a bit more difficult to set up (you will be dealing with 2 meshes, where the motion mesh is a subset of the first one).

If you want to see exactly what is going on, try running the non-rothogonality check from time to time during the mesh motion. You will see that the non-orthogonality changes and this influences the convergence behaviour. You can tackle it in a number of different ways, but basically they fall into the two categories you mention
- a less accurate but more robust discretisation for the laplacian: limited Gauss. You lose accuracy but it does not affect the speed
- non-orthogonal correctors: (typically, 1 is plenty). Increases the cost of simulation but is more accurate than the limited snGrad stuff from above.

Enjoy,

Hrv

lr103476 November 18, 2006 08:43

Hi everybody, In my simulat
 
Hi everybody,

In my simulations concerning moving meshes, I need a very small timestep in OpenFOAM compared to other CFD solvers (e.g. Fluent).

In OpenFOAM I use second-order Crank-Nicholson, so in fact a smaller timestep than Fluent (which is first order in time for moving meshes) should work.

The mesh motion is OK for larger timesteps, but for convergence of the icoDyMFoam solver a much much smaller timestep is necessary. Otherwise the Courant number increases till the solution diverges.

Why is Crank-Nicholson in OpenFOAM not unconditionally stable, as proved in many textbooks (e.g Ferziger, Peric or Hirsch).

Regards, Frank

hjasak November 18, 2006 10:55

Why is Crank-Nicholson in Open
 
Quote:

Why is Crank-Nicholson in OpenFOAM not unconditionally stable, as proved in many textbooks (e.g Ferziger, Peric or Hirsch).
Ferziger, Peric, what page? Hirsch, what page? (cannot find it)

I have shown in my PhD that Crank-Nicolson is not unconditionally unstable (based on numerical diffusion) so what you say is definitely news to me.

Would you clarify please.

Hrv

lr103476 November 18, 2006 11:14

mmm, its been a while that I r
 
mmm, its been a while that I read those books. I always thought that Von Neumann stability analysis shows this unconditional stability of Crank-Nicholson. I guess you're right, p.150 of Ferziger says that the above is true, but in some cases it may become unstable. It also says that the stability limits are problem dependent. (I borrowed Hirsch once, so I cannot review that book now).

In my simulations I need at least about 20000 time steps in one flapping (wing) period. This is 10 times more than in Fluent. Are there ways to tweak OpenFOAM such that a larger timestep can be used at the same accuracy?

What are your ideas about this?

Regards, Frank

ztukovic November 18, 2006 11:45

Hi Frank, Please note that
 
Hi Frank,

Please note that the backward scheme is used for second-order temporal discretisation in Fluent and you are using the Crank-Nicolson scheme in OpenFoam.

Regards, Zeljko

lr103476 November 18, 2006 12:24

Hi Zeljko, Fluent 6.2 is on
 
Hi Zeljko,

Fluent 6.2 is only working with first-order temporal discretisation (Euler) when the dynamic mesh module is used. Therefore using second order Crank-Nicholson in OpenFOAM seems to be more efficient. But, up to now I need very small time-steps with which I am not happy.

Regards, Frank

ztukovic November 18, 2006 12:59

I always use the backward temp
 
I always use the backward temporal scheme with the moving mesh in OpenFOAM.

Zeljko

hjasak November 19, 2006 06:12

Hi Frank, The answer is: go
 
Hi Frank,

The answer is: go back to basics. First, you know that running transient SIMPLE will have no Courant number limit. Start with that, and first order discretisation in time - if you like, you can now have a time step of half an hour.

Of course, with the above the temporal accuracy is a disaster. The next step is to tune the time-step to get the right time resolution for your problem. When you do that, you can switch to second-order temporal scheme (backward it time, like Zeljko says) or Crank-Nicolson.

Being aware that Crank Nicolson is borderline stable, you will now need to monitor the discretisation (dispersion) error to see if it fouls up your solution and again adjust the time-step accordingly.

Keep in mind that you are perfectl welcome to configure the solver to run transient SIMPLE with Crank-Nicolson - I suspect you will need to write some top-level code to do that well.

Finally, it is really important to compare like with like. I know exactly what commercial CFD is doing on cases like this (been there) and comments like " This is 10 times more than in Fluent" when you are comparing first order in a commercial code and full second order in space and time in OpenFOAM are completely meaningless.

Hrv

lr103476 November 19, 2006 06:47

Thanks for your time again! Yo
 
Thanks for your time again! You are completely right that I should change one thing at a time to adjust and tune the OpenFOAM settings for my case.

To solve for the temporal resolution I need about at least 200 timesteps per plunging period. OpenFOAM (with Piso and backward or Crank-Nicholson) needs more than 40.000 timesteps to keep the max Courant number below 2. I use adjustableRunTime.

It seems that I need to adjust icoDyMFoam to use SIMPLE pressure/velocity coupling. Why are all the transient solvers using PISO, by the way? Are you aware of any major problems w.r.t. implementing SIMPLE in the transient codes?

Regards, Frank

hjasak November 19, 2006 08:37

Well, my second PhD supervisor
 
Well, my second PhD supervisor at Imperial College was Raad Issa of the PISO fame. ;-) In any case, PISO is more efficient for cases where the time-step is determined by the need for temporal accuracy and the mesh is good-looking (ie. no small cells in regions of fast flow).

I have written several transient SIMPLE solvers (see development library) and there's never any problem. BTW, the "hot" thing at the moment is non-iterative time advancement for LES, which is shown to be second order and solves the pressure equation only once per time-step.

Hrv

lr103476 November 19, 2006 09:09

Allright. I guess my cases are
 
Allright. I guess my cases are rather extreme, in the sense that large acceleration and velocities occur during plunging (and later rotation). Also the mesh is very nice in terms of skewness and orthogonality, but the cells near the wing's surface are very small.

Therefore, I will give SIMPLE a try since apparentrly the time-step not is dictated by the temporal accuracy. I have your development code, so I am going to try to combine transientSimpleFoam and icoDyMFoam.

Regards and thanks again for your patience,
Frank

gtg627e November 19, 2007 12:28

Dr. Jasak, I read in your p
 
Dr. Jasak,

I read in your post above (from Sunday, November 19, 2006 - 06:37 am) that LES simulations may use non-iterative time advancement.
Can this be done in OpenFoam 1.4.1?
Which solvers and settings would one use to achieve this?

Thank you,

Alessandro

hjasak November 19, 2007 12:49

Hello, Sunday, November 19,
 
Hello,

Sunday, November 19, 2006 - 06:37 am: almost exactly a year ago... You got me confused for a moment.

The reference I refer to came from Fluent News, 2004 and offers two options:

References:

1. R.I. Issa, "Solution of the Implicitly Discretized Fluid Flow Equations by Operator-Splitting," Journal of Computational Physics 62, p. 40-65, 1985.

2. J.B. Perot, "An Analysis of the Fractional-Step Method," Journal of Computational Physics 108, p. 51-58, 1993.

The first one is standard, ie. you will find it in most OpenFOAM solvers (among other things, the guy was my PhD supervisor so I knew about this practically from my first days in CFD). I haven't implemented the second one in OpenFOAM yet, but this should not be too hard. There also some modifications to the settings you can do to allow large time-steps when using PISO, but that's a different topic.

Hope this is useful,

Hrv

gtg627e November 19, 2007 14:09

Dr. Jasak, Thank you for th
 
Dr. Jasak,

Thank you for the references.
I will take a look at them and try to understand their implementation in OF.

appreciate your help,

Alessandro


All times are GMT -4. The time now is 20:50.