|
Simple MAC
1.0
A Marker and Cell Method Solver for the Navier Stokes Equations
|
|
Simple MAC
1.0
A Marker and Cell Method Solver for the Navier Stokes Equations
|
The MAC Method solver and it's persistent data. More...
Public Member Functions | |
| subroutine | initzero () |
| Initialize all variables with zero values. | |
| subroutine | returnu (array) |
| Accessor for U-Velocity. | |
| subroutine | returnv (array) |
| Accessor for V-Velocity. | |
| subroutine | returnp (array) |
| Accessor for Pressure. | |
| subroutine | ghostcondition () |
| Applies ghost cell boundary conditions at the walls. | |
| subroutine | lidcondition () |
| Applies a moving lid to the top of the cavity, and a no slip condition at every other wall. | |
| subroutine | calctstep () |
| Calculates are maximum allowable timestep, \( \Delta t \), based on our stability conditions. | |
| subroutine | calcfngn () |
| Calculates intermediate step of \( F_n \) and \( G_n \). Reference Harlow & Welch (1965) for explanation. | |
| subroutine | calcqn () |
| \( Q_n \) is calculated as a source for the RHS of the pressure Poisson equation. | |
| subroutine | calcvel () |
| This subroutine calculates velocity of the field and stores the data in the module. | |
| subroutine | vort (w) |
| Calculates vorticity based on velocity field. | |
| subroutine | poisson () |
| Poisson solver for the pressure field. | |
| subroutine | parpoisson () |
| Parallelized Poisson solver. | |
Public Attributes | |
| double precision, dimension(nx, ny) | u |
| Instantaneous velocity in x-direction, \( u_{i,j}^{n+1} \). | |
| double precision, dimension(nx, ny) | v |
| Instantaneous velocity in y-direction, \( v_{i,j}^{n+1} \). | |
| double precision, dimension(nx, ny) | p |
| Scalar pressure at cell center \( P_{i,j} \). | |
| double precision, dimension(nx, ny) | fn |
| See Harlow & Welch 1965. | |
| double precision, dimension(nx, ny) | gn |
| See Harlow & Welch 1965. | |
| double precision, dimension(nx, ny) | q |
| Pressure source term for Poisson. | |
| subroutine solver::calcfngn | ( | ) |
Calculates intermediate step of \( F_n \) and \( G_n \). Reference Harlow & Welch (1965) for explanation.
These variables are calculated based on the previous time-step's velocity. \( F_n \) is:
\[ F_{i+\frac{1}{2},j}^n = u_{i+\frac{1}{2},j} + \Delta t \left [ \frac{u^n_{i+\frac{3}{2},j} - 2u^n_{i+\frac{1}{2},j} + u^n_{i-\frac{1}{2},j}}{Re(\Delta x)^2} + \frac{u^n_{i+\frac{1}{2},j-1} - 2u^n_{i+\frac{1}{2},j} + u^n_{i+\frac{1}{2},j-1}}{Re(\Delta y)^2} - \frac{(uv)^n_{i+\frac{1}{2},j+\frac{1}{2}} - (uv)^n_{i-\frac{1}{2},j-\frac{1}{2}}}{\Delta y} \right ] \]
\( G_n \) is defined as:
\[ G_{i,j+\frac{1}{2}}^n = v_{i,j+\frac{1}{2}} + \Delta t \left [ \frac{v^n_{i+1,j+\frac{1}{2}} - 2v^n_{i,j+\frac{1}{2}} + v^n_{i-1,j+\frac{1}{2}}}{Re(\Delta x)^2} + \frac{v^n_{i,j+\frac{3}{2}} - 2v^n_{i,j+\frac{1}{2}} + u^n_{i,j-\frac{1}{2}}}{Re(\Delta y)^2} - \frac{(uv)^n_{i+\frac{1}{2},j+\frac{1}{2}} - (uv)^n_{i-\frac{1}{2},j-\frac{1}{2}}}{\Delta y} - \frac{(v^n_{i,j+1})^2 - (v^n_{i,j})^2}{\Delta y} \right ] \]
| subroutine solver::calcqn | ( | ) |
\( Q_n \) is calculated as a source for the RHS of the pressure Poisson equation.
It is calculated based on the values of \( F_n \) and \( Q_n \) calculated in the routine calcfngn().
\[ Q^n_{i,j} = \left [ \frac{P_{i-1,j} -2P_{i+1,j}}{(\Delta x)^2} + \frac{P_{i,j-1} - 2P_{i,j} + P_{i,j+1}}{(\Delta y)^2}\right ]^{n+1} \]
which may be discretized to show:
\[ Q^n_{i,j} = \frac{1}{\Delta t} \left [ \frac{F^n_{i+\frac{1}{2},j} - F^n_{i-\frac{1}{2},j}}{\Delta x} + \frac{G^n_{i,j+\frac{1}{2}} - G^n_{i,j-\frac{1}{2}}}{\Delta y} \right ]. \]
| subroutine solver::calctstep | ( | ) |
Calculates are maximum allowable timestep, \( \Delta t \), based on our stability conditions.
The first stability condition to be evaluated is:
\[ \Delta t = Re \Delta x \Delta y * r \]
Where \( r \) must be less than 0.25. The second stability condition, called the CFL condition, is shown to be:
\[ \Delta t = \frac{4}{Re (|u|+|v|)^2} \]
We take the minimum of these two conditions to ensure stability.
| subroutine solver::calcvel | ( | ) |
This subroutine calculates velocity of the field and stores the data in the module.
Values of the velocity may be accessed in Python from the routines returnu() and returnv(). The velocity is calculated by the equations:
\[ u_{i+\frac{1}{2},j}^{n+1} = F_{i+\frac{1}{2},j}^n - \frac{\Delta t}{\Delta x}(P_{i+1,j}^{n+1} - P_{i,j}^{n+1}) \]
\[ v_{i+,j+\frac{1}{2} }^{n+1} = G_{i,j+\frac{1}{2} }^n - \frac{\Delta t}{\Delta x}(P_{i,j+1}^{n+1} - P_{i,j}^{n+1}) \]
| subroutine solver::ghostcondition | ( | ) |
| subroutine solver::initzero | ( | ) |
| subroutine solver::lidcondition | ( | ) |
| subroutine solver::parpoisson | ( | ) |
| subroutine solver::poisson | ( | ) |
Poisson solver for the pressure field.
The equation for the pressure field is defined to be:
\[ \nabla^2 P = Q_n \]
For two dimensions the equation becomes
\[ \frac{\partial P}{\partial x} + \frac{\partial P}{\partial y} - Q_n = 0 \]
In discretized form, using an explicit Gauss-Seidel method, we see
\[ P^{n+1}_{i,j} = \frac{1}{4} \left [ P^{n+1}_{i-1,j} + P^{n+1}_{i,j-1} + P^n_{i+1,j} + P^n_{i,j+1} - (\Delta x)^2 Q^n_{i,j}\right ] \]
| subroutine solver::returnp | ( | double precision, dimension(nx,ny), intent(out) | array | ) |
| subroutine solver::returnu | ( | double precision, dimension(nx,ny), intent(out) | array | ) |
| subroutine solver::returnv | ( | double precision, dimension(nx,ny), intent(out) | array | ) |
| subroutine solver::vort | ( | double precision, dimension(nx,ny), intent(out) | w | ) |
| double precision, dimension(nx,ny) solver::fn |
| double precision, dimension(nx,ny) solver::gn |
| double precision, dimension(nx,ny) solver::p |
| double precision, dimension(nx,ny) solver::q |
| double precision, dimension(nx,ny) solver::u |
| double precision, dimension(nx,ny) solver::v |
1.8.1.1
