Poisson solver

Name

Poisson solver -- 

Synopsis


#include <gfs.h>


void        gfs_relax                       (GfsDomain *domain,
                                             gint max_depth,
                                             GfsVariable *u,
                                             GfsVariable *rhs);
void        gfs_residual                    (GfsDomain *domain,
                                             FttTraverseFlags flags,
                                             gint max_depth,
                                             GfsVariable *u,
                                             GfsVariable *rhs,
                                             GfsVariable *res);
void        gfs_poisson_coefficients        (GfsDomain *domain,
                                             GfsVariable *c,
                                             gdouble rho);
void        gfs_poisson                     (GfsDomain *domain,
                                             guint levelmin,
                                             guint depth,
                                             guint nrelax,
                                             GfsVariable *u,
                                             GfsVariable *rhs);

void        gfs_correct_normal_velocities   (const FttCellFace *face,
                                             const GfsAdvectionParams *par);
void        gfs_correct_centered_velocity   (const FttCellFace *face,
                                             const GfsAdvectionParams *par);
struct      GfsProjectionParams;
void        gfs_projection_params_init      (GfsProjectionParams *par);
void        gfs_projection_params_read      (GfsProjectionParams *par,
                                             GtsFile *fp);
void        gfs_projection_params_write     (GfsProjectionParams *par,
                                             FILE *fp);
void        gfs_mac_projection              (GfsDomain *domain,
                                             GfsProjectionParams *par,
                                             GfsAdvectionParams *apar);
void        gfs_approximate_projection      (GfsDomain *domain,
                                             GfsProjectionParams *par,
                                             GfsAdvectionParams *apar);
GfsNorm     gfs_domain_norm_residual        (GfsDomain *domain,
                                             FttTraverseFlags flags,
                                             gint max_depth);

Description

Details

gfs_relax ()

void        gfs_relax                       (GfsDomain *domain,
                                             gint max_depth,
                                             GfsVariable *u,
                                             GfsVariable *rhs);

Apply one pass of a Jacobi relaxation to all the leaf cells of domain with a level inferior or equal to max_depth and to all the cells at level max_depth. The relaxation should converge (if the right-hand-side rhs verifies the solvability conditions) toward the solution of a Poisson equation for u at the maximum depth.

domain : the domain to relax.
max_depth : the maximum depth of the domain to relax.
u : the variable to use as left-hand side.
rhs : the variable to use as right-hand side.


gfs_residual ()

void        gfs_residual                    (GfsDomain *domain,
                                             FttTraverseFlags flags,
                                             gint max_depth,
                                             GfsVariable *u,
                                             GfsVariable *rhs,
                                             GfsVariable *res);

For each cell of domain, computes the sum of the residual over the volume of the cell for a Poisson equation with u as left-hand-side and rhs as right-hand-side. Stores the result in res.

domain : a domain.
flags : which types of cells are to be visited.
max_depth : maximum depth of the traversal.
u : the variable to use as left-hand side.
rhs : the variable to use as right-hand side.
res : the variable to use to store the residual.


gfs_poisson_coefficients ()

void        gfs_poisson_coefficients        (GfsDomain *domain,
                                             GfsVariable *c,
                                             gdouble rho);

Initializes the face coefficients for the Poisson equation.

domain : a GfsDomain.
c : the volume fraction.
rho : the relative density.


gfs_poisson ()

void        gfs_poisson                     (GfsDomain *domain,
                                             guint levelmin,
                                             guint depth,
                                             guint nrelax,
                                             GfsVariable *u,
                                             GfsVariable *rhs);

Apply one multigrid iteration to the Poisson equation defined by u and rhs.

The initial value of GFS_RES on the leaves of root must be set to the residual of the Poisson equation (using gfs_residual()).

The face coefficients must be set using gfs_poisson_coefficients().

The values of u on the leaf cells are updated as well as the values of GFS_RES (i.e. the cell tree is ready for another iteration).

domain : the domain on which to solve the Poisson equation.
levelmin : the top level of the multigrid hierarchy.
depth : the total depth of the domain.
nrelax : the number of relaxations to apply at each level.
u : the variable to use as left-hand side.
rhs : the variable to use as right-hand side.


gfs_correct_normal_velocities ()

void        gfs_correct_normal_velocities   (const FttCellFace *face,
                                             const GfsAdvectionParams *par);

Corrects the velocity normal to face with the gradient of the pressure in the corresponding direction.

Also adds half this gradient to the corresponding GFS_G variable of the cells on either sides of the face. If this function is called for each face of a domain, the GFS_G variables in each cell will then contain the (interpolated) value of the pressure gradient.

This functions assumes that the face coefficients for the Poisson problem have been initialized using gfs_poisson_coefficients().

face : a FttCellFace.
par : the advection parameters.


gfs_correct_centered_velocity ()

void        gfs_correct_centered_velocity   (const FttCellFace *face,
                                             const GfsAdvectionParams *par);

(Half) Corrects the centered velocity in each cell on either side of face with the pressure gradient. For the correction to be complete, this function has to be called for each face of the domain.

face : a FttCellFace.
par : the advection parameters.


struct GfsProjectionParams

struct GfsProjectionParams {
  gdouble tolerance;
  guint nrelax;
  guint minlevel;
  guint nitermax;

  guint niter;
  GfsNorm residual_before, residual;
};


gfs_projection_params_init ()

void        gfs_projection_params_init      (GfsProjectionParams *par);

par : 


gfs_projection_params_read ()

void        gfs_projection_params_read      (GfsProjectionParams *par,
                                             GtsFile *fp);

par : 
fp : 


gfs_projection_params_write ()

void        gfs_projection_params_write     (GfsProjectionParams *par,
                                             FILE *fp);

Writes in fp a text representation of the projection parameters par.

par : the projection parameters.
fp : a file pointer.


gfs_mac_projection ()

void        gfs_mac_projection              (GfsDomain *domain,
                                             GfsProjectionParams *par,
                                             GfsAdvectionParams *apar);

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 the GFS_G variables and is obtained by simple averaging from the face values to the center).

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.

domain : a GfsDomain.
par : the projection control parameters.
apar : the advection parameters.


gfs_approximate_projection ()

void        gfs_approximate_projection      (GfsDomain *domain,
                                             GfsProjectionParams *par,
                                             GfsAdvectionParams *apar);

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.

domain : a GfsDomain.
par : the projection control parameters.
apar : the advection parameters.


gfs_domain_norm_residual ()

GfsNorm     gfs_domain_norm_residual        (GfsDomain *domain,
                                             FttTraverseFlags flags,
                                             gint max_depth);

Traverses the domain defined by domain using gfs_domain_cell_traverse() and gathers norm statistics about the volume weighted relative residual (i.e. the sum of the residual over the volume defined by each cell divided by the total volume of the cell).

domain : the domain to obtain the norm from.
flags : which types of cells are to be visited.
max_depth : maximum depth of the traversal.
Returns : a GfsNorm containing the norm statistics about the volume weighted relative residual.