p4est  1.1
 All Data Structures Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
Data Structures | Functions
p4est_step3.c File Reference

This 2D example program uses p4est to solve a simple advection problem. More...

#include <p4est_vtk.h>
#include <p4est_bits.h>
#include <p4est_extended.h>
#include <p4est_iterate.h>

Data Structures

struct  step3_data_t
 Per-quadrant data for this example. More...
 
struct  step3_ctx_t
 The example parameters. More...
 

Functions

static double step3_initial_condition (double x[], double du[], step3_ctx_t *ctx)
 Compute the value and derivatives of the initial condition. More...
 
static void step3_get_midpoint (p4est_t *p4est, p4est_topidx_t which_tree, p4est_quadrant_t *q, double xyz[3])
 Get the coordinates of the midpoint of a quadrant. More...
 
static void step3_init_initial_condition (p4est_t *p4est, p4est_topidx_t which_tree, p4est_quadrant_t *q)
 Initialize the initial condition data of a quadrant. More...
 
static double step3_error_sqr_estimate (p4est_quadrant_t *q)
 Estimate the square of the approximation error on a quadrant. More...
 
static int step3_refine_err_estimate (p4est_t *p4est, p4est_topidx_t which_tree, p4est_quadrant_t *q)
 Refine by the L2 error estimate. More...
 
static int step3_coarsen_initial_condition (p4est_t *p4est, p4est_topidx_t which_tree, p4est_quadrant_t *children[])
 Coarsen by the L2 error estimate of the initial condition. More...
 
static int step3_coarsen_err_estimate (p4est_t *p4est, p4est_topidx_t which_tree, p4est_quadrant_t *children[])
 Coarsen by the L2 error estimate of the current state approximation. More...
 
static void step3_replace_quads (p4est_t *p4est, p4est_topidx_t which_tree, int num_outgoing, p4est_quadrant_t *outgoing[], int num_incoming, p4est_quadrant_t *incoming[])
 Initialize the state variables of incoming quadrants from outgoing quadrants. More...
 
static void step3_interpolate_solution (p4est_iter_volume_info_t *info, void *user_data)
 Callback function for interpolating the solution from quadrant midpoints to corners. More...
 
static void step3_write_solution (p4est_t *p4est, int timestep)
 Write the state variable to vtk format, one file per process. More...
 
static void step3_quad_divergence (p4est_iter_volume_info_t *info, void *user_data)
 Approximate the divergence of (vu) on each quadrant. More...
 
static void step3_upwind_flux (p4est_iter_face_info_t *info, void *user_data)
 Approximate the flux across a boundary between quadrants. More...
 
static void step3_timestep_update (p4est_iter_volume_info_t *info, void *user_data)
 Compute the new value of the state from the computed time derivative. More...
 
static void step3_reset_derivatives (p4est_iter_volume_info_t *info, void *user_data)
 Reset the approximate derivatives. More...
 
static void step3_minmod_estimate (p4est_iter_face_info_t *info, void *user_data)
 For two quadrants on either side of a face, estimate the derivative normal to the face. More...
 
static void step3_compute_max (p4est_iter_volume_info_t *info, void *user_data)
 Compute the maximum state value. More...
 
static double step3_get_timestep (p4est_t *p4est)
 Compute the timestep. More...
 
static void step3_timestep (p4est_t *p4est, double time)
 Timestep the advection problem. More...
 
int main (int argc, char **argv)
 The main step 3 program. More...
 

Detailed Description

This 2D example program uses p4est to solve a simple advection problem.

It is numerically very simple, and intended to demonstrate several methods of interacting with the p4est data after it has been refined and partitioned. It demonstrates the construction of ghost layers (see p4est_ghost_t in p4est_ghost.h) and communication of ghost-layer data, and it demonstrates interacting with the quadrants and quadrant boundaries through the p4est_iterate() routine (see p4est_iterate.h).

Function Documentation

int main ( int  argc,
char **  argv 
)

The main step 3 program.

Setup of the example parameters; create the forest, with the state variable stored in the quadrant data; refine, balance, and partition the forest; timestep; clean up, and exit.

static int step3_coarsen_err_estimate ( p4est_t p4est,
p4est_topidx_t  which_tree,
p4est_quadrant_t children[] 
)
static

Coarsen by the L2 error estimate of the current state approximation.

Given the maximum global error, we enforce that each quadrant's portion of the error must not exceed its fraction of the total volume of the domain (which is 1).

This function matches the p4est_coarsen_t prototype that is used by p4est_coarsen() and p4est_coarsen_ext().

Parameters
[in]p4estthe forest
[in]which_treethe tree in the forest containing children
[in]childrena family of quadrants
Returns
1 if children should be coarsened, 0 otherwise.
static int step3_coarsen_initial_condition ( p4est_t p4est,
p4est_topidx_t  which_tree,
p4est_quadrant_t children[] 
)
static

Coarsen by the L2 error estimate of the initial condition.

Given the maximum global error, we enforce that each quadrant's portion of the error must not exceed is fraction of the total volume of the domain (which is 1).

Parameters
[in]p4estthe forest
[in]which_treethe tree in the forest containing children
[in]childrena family of quadrants
Returns
1 if children should be coarsened, 0 otherwise.
static void step3_compute_max ( p4est_iter_volume_info_t info,
void *  user_data 
)
static

Compute the maximum state value.

This function updates the maximum value from the value of a single cell.

This function matches the p4est_iter_volume_t prototype used by p4est_iterate().

Parameters
[in]infothe information about this quadrant that has been populated by p4est_iterate()
[in,out]user_datathe user_data given to p4est_iterate(): in this case, it points to the maximum value that will be updated
static double step3_error_sqr_estimate ( p4est_quadrant_t q)
static

Estimate the square of the approximation error on a quadrant.

We compute our estimate by integrating the difference of a constant approximation at the midpoint and a linear approximation that interpolates at the midpoint.

Parameters
[in]qa quadrant
Returns
the square of the error estimate for the state variables contained in q's data.
static void step3_get_midpoint ( p4est_t p4est,
p4est_topidx_t  which_tree,
p4est_quadrant_t q,
double  xyz[3] 
)
static

Get the coordinates of the midpoint of a quadrant.

Parameters
[in]p4estthe forest
[in]which_treethe tree in the forest containing q
[in]qthe quadrant
[out]xyzthe coordinates of the midpoint of q
static double step3_get_timestep ( p4est_t p4est)
static

Compute the timestep.

Find the smallest quadrant and scale the timestep based on that length and the advection velocity.

Parameters
[in]p4estthe forest
Returns
the timestep.
static void step3_init_initial_condition ( p4est_t p4est,
p4est_topidx_t  which_tree,
p4est_quadrant_t q 
)
static

Initialize the initial condition data of a quadrant.

This function matches the p4est_init_t prototype that is used by p4est_new(), p4est_refine(), p4est_coarsen(), and p4est_balance().

Parameters
[in]p4estthe forest
[in]which_treethe tree in the forest containing q
[in,out]qthe quadrant whose data gets initialized
static double step3_initial_condition ( double  x[],
double  du[],
step3_ctx_t ctx 
)
static

Compute the value and derivatives of the initial condition.

Parameters
[in]xthe coordinates
[out]duthe derivative at x
[in]ctxthe example parameters
Returns
the initial condition at x
static void step3_interpolate_solution ( p4est_iter_volume_info_t info,
void *  user_data 
)
static

Callback function for interpolating the solution from quadrant midpoints to corners.

The function p4est_iterate() takes as an argument a p4est_iter_volume_t callback function, which it executes at every local quadrant (see p4est_iterate.h). This function matches the p4est_iter_volume_t prototype.

In this example, we use the callback function to interpolate the state variable to the corners, and write those corners into an array so that they can be written out.

Parameters
[in]infothe information about this quadrant that has been populated by p4est_iterate()
[in,out]user_datathe user_data that was given as an argument to p4est_iterate: in this case, it points to the array of corner values that we want to write. The values for the corner of the quadrant described by info are written during the execution of the callback.
static void step3_minmod_estimate ( p4est_iter_face_info_t info,
void *  user_data 
)
static

For two quadrants on either side of a face, estimate the derivative normal to the face.

This function matches the p4est_iter_face_t prototype used by p4est_iterate().

Parameters
[in]infothe information about this quadrant that has been populated by p4est_iterate()
[in]user_datathe user_data given to p4est_iterate(): in this case, it points to the ghost_data array, which contains the step3_data_t data for all of the ghost cells, which was populated by p4est_ghost_exchange_data()
static void step3_quad_divergence ( p4est_iter_volume_info_t info,
void *  user_data 
)
static

Approximate the divergence of (vu) on each quadrant.

We use piecewise constant approximations on each quadrant, so the value is always 0.

Like step3_interpolate_solution(), this function matches the p4est_iter_volume_t prototype used by p4est_iterate().

Parameters
[in]infothe information about the quadrant populated by p4est_iterate()
[in]user_datanot used
static int step3_refine_err_estimate ( p4est_t p4est,
p4est_topidx_t  which_tree,
p4est_quadrant_t q 
)
static

Refine by the L2 error estimate.

Given the maximum global error, we enforce that each quadrant's portion of the error must not exceed is fraction of the total volume of the domain (which is 1).

This function matches the p4est_refine_t prototype that is used by p4est_refine() and p4est_refine_ext().

Parameters
[in]p4estthe forest
[in]which_treethe tree in the forest containing q
[in]qthe quadrant
Returns
1 if q should be refined, 0 otherwise.
static void step3_replace_quads ( p4est_t p4est,
p4est_topidx_t  which_tree,
int  num_outgoing,
p4est_quadrant_t outgoing[],
int  num_incoming,
p4est_quadrant_t incoming[] 
)
static

Initialize the state variables of incoming quadrants from outgoing quadrants.

The functions p4est_refine_ext(), p4est_coarsen_ext(), and p4est_balance_ext() take as an argument a p4est_replace_t callback function, which allows one to setup the quadrant data of incoming quadrants from the data of outgoing quadrants, before the outgoing data is destroyed. This function matches the p4est_replace_t prototype.

In this example, we linearly interpolate the state variable of a quadrant that is refined to its children, and we average the midpoints of children that are being coarsened to the parent.

Parameters
[in]p4estthe forest
[in]which_treethe tree in the forest containing children
[in]num_outgoingthe number of quadrants that are being replaced: either 1 if a quadrant is being refined, or P4EST_CHILDREN if a family of children are being coarsened.
[in]outgoingthe outgoing quadrants
[in]num_incomingthe number of quadrants that are being added: either P4EST_CHILDREN if a quadrant is being refined, or 1 if a family of children are being coarsened.
[in,out]incomingquadrants whose data are initialized.
static void step3_reset_derivatives ( p4est_iter_volume_info_t info,
void *  user_data 
)
static

Reset the approximate derivatives.

p4est_iterate() has an invariant to the order of callback execution: the p4est_iter_volume_t callback will be executed on a quadrant before the p4est_iter_face_t callbacks are executed on its faces. This function resets the derivative stored in the quadrant's data before step3_minmod_estimate() updates the derivative based on the face neighbors.

This function matches the p4est_iter_volume_t prototype used by p4est_iterate().

Parameters
[in]infothe information about this quadrant that has been populated by p4est_iterate()
[in]user_datanot used
static void step3_timestep ( p4est_t p4est,
double  time 
)
static

Timestep the advection problem.

Update the state, refine, repartition, and write the solution to file.

Parameters
[in,out]p4estthe forest, whose state is updated
[in]timethe end time
static void step3_timestep_update ( p4est_iter_volume_info_t info,
void *  user_data 
)
static

Compute the new value of the state from the computed time derivative.

We use a simple forward Euler scheme.

The derivative was computed by a p4est_iterate() loop by the callbacks step3_quad_divergence() and step3_upwind_flux(). Now we multiply this by the timestep and add to the current solution.

This function matches the p4est_iter_volume_t prototype used by p4est_iterate().

Parameters
[in]infothe information about this quadrant that has been populated by p4est_iterate()
[in]user_datathe user_data given to p4est_iterate(): in this case, it points to the timestep.
static void step3_upwind_flux ( p4est_iter_face_info_t info,
void *  user_data 
)
static

Approximate the flux across a boundary between quadrants.

We use a very simple upwind numerical flux.

This function matches the p4est_iter_face_t prototype used by p4est_iterate().

Parameters
[in]infothe information about the quadrants on either side of the interface, populated by p4est_iterate()
[in]user_datathe user_data given to p4est_iterate(): in this case, it points to the ghost_data array, which contains the step3_data_t data for all of the ghost cells, which was populated by p4est_ghost_exchange_data()
static void step3_write_solution ( p4est_t p4est,
int  timestep 
)
static

Write the state variable to vtk format, one file per process.

Parameters
[in]p4estthe forest, whose quadrant data contains the state
[in]timestepthe timestep number, used to name the output files