|April 25, 2006, 02:15||
Neumann pressure BC and velocity field
Good time of the day.
I try to create a simple 2-D transient CFD solver (it's one-phase now, fluid is fully compressible, laminar flow). There is a problem with pressure BC.
I use an 1-D channel (constant section) test task to validate the model. The inlet BC is Dirichlet for velocity (mass flux and velocity is given on inlet face) and Neumann for pressure (pressure at inlet face supposed to be equal to pressure in the center of inlet control volume). The outlet BC is Dirichlet for pressure and Neumann for velocity (there are no problems with outlet).
I implemented a monophase SIMPLE algorithm with colocated grid and Rhie-Chow interpolation for suppresing the checkerboard effect.
If it matters, the test parameters are:
Initial velocities are zero, initial pressure is atmospheric at all cells. Temperature is always constant (850 *C, I like this number ).
Inlet velocity: 1 m/s
Inlet mass flux: 0.325 kg/s
Outlet pressure: atmospheric (101300 Pa)
Fluid density: Ideal gas law
Fluid molar mass: 30 g/mol
Fluid viscosity: 0.00002 Pa*s
Channel section: 1 m2
Timestep: 0.005 s
Number of cells: 5 (just a test)
After calculating the first timestep, the pressure field is correct (commercial FlowVision CFD program is used for compare, this software is certificated, www.fv-tech.com, free demo available for unlimited time but max number of cells for demo version is 15000). Excess pressure at inlet: ~120 Pa, at outlet: ~0 Pa (i. e. atmospheric).
The velocity field after first timestep is correct except inlet control volume. Velocity in this CV is about two times lower than it must be (~0.4 m/s instead of ~0.85 m/s). Other CV velocities are seem to be OK (compared with FlowVision). I think this problem is caused by following. Pressure at cell faces for velocity (momentum) equation is calculated by linear interpolation of pressures at neighbour CVs of the face. If we suppose correct pressure field we may not solve for pressure but must get correct velocities after solving momentum equation. In momentum equation, the more pressure difference at cell faces is, the more velocity after the timestep will be (channel section is constant, IMHO convective terms may be neglected due to high pressure differences in this example). So, at all CVs except inlet CV, the pressure difference at faces will be about two times higher than in inlet CV, because for inlet CV pressure on the inlet face equals the pressure at the center of control volume. I tried to illustrate this on a picture below, "=" is pressure value, "O" is CV center, dots - for formatting text only:
Inlet cell pressure (Neuman BC):
Other cell pressures:
It can be seen, that the face pressure difference for inlet cell is about two times higher than for other cells (pressure curve is almost linear).
But I'm sure that the pressure BCs are correct, isn't it? So where is my mistake?
Sorry for bad english, my native language is russian...
[BTW, are there any tags like [b] or in this forum?]
Thanx for your help.