Diffusion solver

Name

Diffusion solver -- 

Synopsis


#include <gfs.h>


void        gfs_diffusion_rhs               (GfsDomain *domain,
                                             GfsVariable *v,
                                             GfsVariable *rhs,
                                             GfsVariable *dia,
                                             gdouble beta);
void        gfs_diffusion_residual          (GfsDomain *domain,
                                             GfsVariable *u,
                                             GfsVariable *rhs,
                                             GfsVariable *dia,
                                             GfsVariable *res);
void        gfs_diffusion_coefficients      (GfsDomain *domain,
                                             GfsSourceDiffusion *d,
                                             gdouble dt,
                                             GfsVariable *dia,
                                             GfsFunction *alpha,
                                             gdouble beta);
void        gfs_diffusion_cycle             (GfsDomain *domain,
                                             guint levelmin,
                                             guint depth,
                                             guint nrelax,
                                             GfsVariable *u,
                                             GfsVariable *rhs,
                                             GfsVariable *dia,
                                             GfsVariable *res);

Description

Details

gfs_diffusion_rhs ()

void        gfs_diffusion_rhs               (GfsDomain *domain,
                                             GfsVariable *v,
                                             GfsVariable *rhs,
                                             GfsVariable *dia,
                                             gdouble beta);

Adds to the rhs variable of cell the right-hand side of the diffusion equation for variable v.

The diffusion coefficients must have been already set using gfs_diffusion_coefficients().

domain :

a GfsDomain.

v :

a GfsVariable.

rhs :

a GfsVariable.

dia :

the diagonal weight.

beta :

the implicitness parameter (0.5 Crank-Nicholson, 1. backward Euler).


gfs_diffusion_residual ()

void        gfs_diffusion_residual          (GfsDomain *domain,
                                             GfsVariable *u,
                                             GfsVariable *rhs,
                                             GfsVariable *dia,
                                             GfsVariable *res);

Sets the res variable of each leaf cell of domain to the residual of the diffusion equation for v.

The diffusion coefficients must have been set using gfs_diffusion_coefficients() and the right-hand side using gfs_diffusion_rhs().

domain :

a GfsDomain.

u :

the variable to use as left-hand side.

rhs :

the right-hand side.

dia :

the diagonal weight.

res :

the residual.


gfs_diffusion_coefficients ()

void        gfs_diffusion_coefficients      (GfsDomain *domain,
                                             GfsSourceDiffusion *d,
                                             gdouble dt,
                                             GfsVariable *dia,
                                             GfsFunction *alpha,
                                             gdouble beta);

Initializes the face coefficients for the diffusion equation.

domain :

a GfsDomain.

d :

a GfsSourceDiffusion.

dt :

the time-step.

dia :

where to store the diagonal weight.

alpha :

the inverse of density or NULL.

beta :

the implicitness parameter (0.5 Crank-Nicholson, 1. backward Euler).


gfs_diffusion_cycle ()

void        gfs_diffusion_cycle             (GfsDomain *domain,
                                             guint levelmin,
                                             guint depth,
                                             guint nrelax,
                                             GfsVariable *u,
                                             GfsVariable *rhs,
                                             GfsVariable *dia,
                                             GfsVariable *res);

Apply one multigrid iteration to the diffusion equation for u.

The initial value of res on the leaves of root must be set to the residual of the diffusion equation using gfs_diffusion_residual().

The diffusion coefficients must be set using gfs_diffusion_coefficients().

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

domain :

the domain on which to solve the diffusion 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 right-hand side.

dia :

the diagonal weight.

res :

the residual.