Gerris
Functions

timestep.c File Reference

Timestepping. More...

Functions

void gfs_reset_gradients (GfsDomain *domain, guint dimension, GfsVariable **g)
void gfs_scale_gradients (GfsDomain *domain, guint dimension, GfsVariable **g)
void gfs_correct_normal_velocities (GfsDomain *domain, guint dimension, GfsVariable *p, GfsVariable **g, gdouble dt)
void gfs_velocity_face_sources (GfsDomain *domain, GfsVariable **u, gdouble dt, GfsFunction *alpha, GfsVariable **g)
void gfs_update_gradients (GfsDomain *domain, GfsVariable *p, GfsFunction *alpha, GfsVariable **g)
void gfs_mac_projection (GfsDomain *domain, GfsMultilevelParams *par, gdouble dt, GfsVariable *p, GfsFunction *alpha, GfsVariable **g, void(*divergence_hook)(GfsDomain *domain, gdouble dt, GfsVariable *div))
void gfs_correct_centered_velocities (GfsDomain *domain, guint dimension, GfsVariable **g, gdouble dt)
void gfs_approximate_projection (GfsDomain *domain, GfsMultilevelParams *par, gdouble dt, GfsVariable *p, GfsFunction *alpha, GfsVariable *res, GfsVariable **g, void(*divergence_hook)(GfsDomain *domain, gdouble dt, GfsVariable *div))
void gfs_predicted_face_velocities (GfsDomain *domain, guint d, GfsAdvectionParams *par)
void gfs_diffusion (GfsDomain *domain, GfsMultilevelParams *par, GfsVariable *v, GfsVariable *rhs, GfsVariable *rhoc, GfsVariable *metric)
void gfs_add_sinking_velocity (GfsDomain *domain, GfsAdvectionParams *par)
void gfs_remove_sinking_velocity (GfsDomain *domain, GfsAdvectionParams *par)
void gfs_centered_velocity_advection_diffusion (GfsDomain *domain, guint dimension, GfsAdvectionParams *par, GfsVariable **gmac, GfsVariable **g, GfsFunction *alpha)
void gfs_tracer_advection_diffusion (GfsDomain *domain, GfsAdvectionParams *par, GfsFunction *alpha)

Detailed Description

Timestepping.


Function Documentation

void gfs_add_sinking_velocity ( GfsDomain *  domain,
GfsAdvectionParams *  par 
)
Parameters:
domaina #GfsDomain.
parthe advection parameters.

Adds the sinking velocity to the MAC velocity field of domain.

Here is the call graph for this function:

void gfs_approximate_projection ( GfsDomain *  domain,
GfsMultilevelParams *  par,
gdouble  dt,
GfsVariable *  p,
GfsFunction *  alpha,
GfsVariable *  res,
GfsVariable **  g,
void(*)(GfsDomain *domain, gdouble dt, GfsVariable *div)  divergence_hook 
)
Parameters:
domaina #GfsDomain.
parthe projection control parameters.
dtthe timestep.
pthe pressure.
alphathe Poisson equation gradient weight.
resthe residual or NULL.
gwhere to store the pressure gradient.
divergence_hooka hook function or NULL.

Corrects the centered velocity field on the leaf level of domain using an approximate projection. The resulting centered velocity field is approximately divergence free. The (potential) pressure field is also obtained as a by-product.

The residual field of the par projection parameters is set to the norm of the residual (on the MAC grid) after the projection. The niter field of the par projection parameters is set to the number of iterations performed to solve the Poisson equation. The other projection parameters are not modified.

The Poisson equation for the pressure is first solved on a MAC grid where the MAC velocities are obtained from the centered velocities by simple averaging. The resulting pressure gradients (defined on the faces) are then averaged down on the center of the cells to correct the centered velocity.

Here is the call graph for this function:

void gfs_centered_velocity_advection_diffusion ( GfsDomain *  domain,
guint  dimension,
GfsAdvectionParams *  par,
GfsVariable **  gmac,
GfsVariable **  g,
GfsFunction *  alpha 
)
Parameters:
domaina #GfsDomain.
dimensionthe number of dimensions (2 or 3).
parthe advection parameters.
gmacthe MAC pressure gradient.
gthe pressure gradient.
alphathe inverse of density or NULL.

Advects the (centered) velocity field using the current face-centered (MAC) velocity field and par->flux to compute the velocity flux through the faces of each cell.

For each component of the velocity, before calling the par->flux function the face values are first defined (at time t + dt/2) and can then be used within the par->flux function.

"Small" cut cells are treated using a cell-merging approach to avoid any restrictive CFL stability condition.

Here is the call graph for this function:

void gfs_correct_centered_velocities ( GfsDomain *  domain,
guint  dimension,
GfsVariable **  g,
gdouble  dt 
)
Parameters:
domaina #GfsDomain.
dimensionthe number of dimensions (2 or 3).
gthe pressure gradient.
dtthe timestep.

Corrects the velocity field of domain using the pressure gradient stored in g[].

The g[] variables are freed by this function.

Here is the call graph for this function:

Here is the caller graph for this function:

void gfs_correct_normal_velocities ( GfsDomain *  domain,
guint  dimension,
GfsVariable *  p,
GfsVariable **  g,
gdouble  dt 
)
Parameters:
domaina #GfsDomain.
dimensionthe number of dimensions (2 or 3).
pthe pressure field.
gwhere to store the pressure gradient or NULL.
dtthe timestep.

Corrects the normal velocity field of domain using p and and dt.

Assumes that the Poisson weighting coefficients have already been computed using gfs_poisson_coefficients().

Also fills the g variables (if not NULL) with the centered gradient of p.

Here is the call graph for this function:

Here is the caller graph for this function:

void gfs_diffusion ( GfsDomain *  domain,
GfsMultilevelParams *  par,
GfsVariable *  v,
GfsVariable *  rhs,
GfsVariable *  rhoc,
GfsVariable *  metric 
)
Parameters:
domaina #GfsDomain.
parthe multilevel parameters.
va #GfsVariable.
rhsthe right-hand side.
rhocthe mass.
metricthe metric term.

Solves a diffusion equation for variable v using a Crank-Nicholson scheme with multilevel relaxations.

Diffusion coefficients must have been set using gfs_diffusion_coefficients() and a right-hand side defined using gfs_diffusion_rhs().

Here is the call graph for this function:

void gfs_mac_projection ( GfsDomain *  domain,
GfsMultilevelParams *  par,
gdouble  dt,
GfsVariable *  p,
GfsFunction *  alpha,
GfsVariable **  g,
void(*)(GfsDomain *domain, gdouble dt, GfsVariable *div)  divergence_hook 
)
Parameters:
domaina #GfsDomain.
parthe projection control parameters.
dtthe timestep.
pthe pressure.
alphathe Poisson equation gradient weight.
gwhere to store the pressure gradient.
divergence_hooka hook function or NULL.

Corrects the face-centered velocity field (MAC field) on the leaf level of domain using an exact (MAC) projection. The resulting face-centered velocity field is (almost) exactly divergence free. The (potential) pressure field is also obtained as a by-product as well as its gradient at the center of the leaf cells of the domain. The gradient is stored in newly-allocated g[] variables and is obtained by simple averaging from the face values to the center. The newly-allocated g[] variables should be freed when not needed anymore.

The residual field of the par projection parameters is set to the norm of the residual after the projection. The niter field of the par projection parameters is set to the number of iterations performed to solve the Poisson equation. The other projection parameters are not modified.

Here is the call graph for this function:

void gfs_predicted_face_velocities ( GfsDomain *  domain,
guint  d,
GfsAdvectionParams *  par 
)
Parameters:
domaina #GfsDomain.
dthe number of dimensions (2 or 3).
parthe advection parameters.

Fills the face (MAC) normal velocities of each leaf cell of domain with the predicted values at time t + dt/2 using a godunov type advection scheme.

Here is the call graph for this function:

void gfs_remove_sinking_velocity ( GfsDomain *  domain,
GfsAdvectionParams *  par 
)
Parameters:
domaina #GfsDomain.
parthe advection parameters.

Removes the sinking velocity from the MAC velocity field of domain.

Here is the call graph for this function:

void gfs_reset_gradients ( GfsDomain *  domain,
guint  dimension,
GfsVariable **  g 
)
Parameters:
domaina #GfsDomain.
gan array of dimension #GfsVariables.
dimensionthe number of dimension.

Resets the gradient vector g.

Here is the caller graph for this function:

void gfs_scale_gradients ( GfsDomain *  domain,
guint  dimension,
GfsVariable **  g 
)
Parameters:
domaina #GfsDomain.
dimensionthe number of dimensions.
gthe components of the gradient.

Scales the gradient accumulated in g (typically using gfs_correct_normal_velocities()).

Here is the call graph for this function:

Here is the caller graph for this function:

void gfs_tracer_advection_diffusion ( GfsDomain *  domain,
GfsAdvectionParams *  par,
GfsFunction *  alpha 
)
Parameters:
domaina #GfsDomain.
parthe advection parameters.
alphathe specific volume or NULL.

Advects the v field of par using the current face-centered (MAC) velocity field.

Here is the call graph for this function:

Here is the caller graph for this function:

void gfs_update_gradients ( GfsDomain *  domain,
GfsVariable *  p,
GfsFunction *  alpha,
GfsVariable **  g 
)
Parameters:
domaina #GfsDomain.
pthe pressure.
alphathe Poisson equation gradient weight.
gwhere to store the pressure gradient.

Updates the gradients in g using p and alpha.

Here is the call graph for this function:

void gfs_velocity_face_sources ( GfsDomain *  domain,
GfsVariable **  u,
gdouble  dt,
GfsFunction *  alpha,
GfsVariable **  g 
)
Parameters:
domaina #GfsDomain.
uthe velocity vector.
dtthe timestep.
alphathe specific volume.
gthe gradient vector.

Add source terms on the velocity component of each cell faces.

Here is the call graph for this function:

Here is the caller graph for this function:

 All Data Structures Files Functions Variables