Creating new discretization schemas for PDE

Method of Lines provides a separation of the spatial and temporal discretization, any combination may be possible to combine (even though it can lead to numerical unstabilities). For this reason, Simflowny separates the spatial discretization from the temporal discretization schema in two documents, Spatial Operator Discretization and PDE discretization Schema respectively. The first one corresponds to spatial discretization. The last document is an integration of an abstract spatial discretization into the time discretization.

Creating new Spatial Operator Discretizations

As stated previously, a spatial operator discretization document is used in Simflowny to define how the PDE spatial derivatives are going to be discretized from continuum.

In order to create an spatial discretization, from the document manager click on the Add button and select "Spatial Operator Discretization". The editor will open with the empty skeleton of a spatial operator discretization.

Head

Head information is commonly used in several documents in Simflowny. See more detail in Head.

Discretizations

The main core of this document is inside the Discretizations element. Equations must be composed from terms having arbitrary-order derivatives. Each derivative can be discretized in two ways, directly or recursively. By directly, we mean that a second order derivative needs a second order discretization. By recursively, we mean that a second order derivative can be discretized with a first order discretization, and the result of this first discretization is also discretized by a first order discretization again. Evidently, the result may be not the same, but it is not wrong in principle. For this reason, it is a good practice to always set a first order discretization, even though higher order discretizations are also set. Doing this, there is assured that any arbitrary order PDE can be solve with the spatial operator discretization.

Having said that, a clarification of the elements to define a Discretization is required. The first element is Derivative Type, in which the user must define the order of the derivatives targeted to this discretization. Derivatives with order equal or higher than two have a particularity. The discretization of crossed derivatives is not always possible or straightforward. To distinguish those cases, the Coordinate Combination element is provided. The filling of this element is made through a string matching the expression i[ijk]*. The value have to start with i, representing a first discrete coordinate. Then any of the three indices must appear to fit a lengh of the derivative order. If using, for example, ij, it is stated that the Derivative applies to the second order crossed derivatives ($\partial_x\partial_y\mbox{ }u$, $\partial_z\partial_x\mbox{ }u$, etc.) . In the picture below there is the discretization for second-order derivatives with same coordinates ($\partial_x\partial_x\mbox{ }u$, $\partial_y\partial_y\mbox{ }u$, etc.).

images/UserGuide/secondOrderDerivativeDiscretization.png

The following property to set into Derivative is the Main Schema. The transformation rule that will perform the spatial discretization of this derivative must be set. In order to know in detail what a transformation rule is, please, check Transformation Rule documentation .

Normally, spatial discretization uses information from neighbour cells. The number of neighbour cells depends on the schema. But, what happens if the neighbour information is not available, for instance, near the domain limit or in a region interphase. In Simflowny, there are three posibilites to deal with this problem:

Using physical boundary conditions. Configured in the problem, the boundary conditions set a simple behaviour of cells in domain limits.
Setting unbalanced spatial discretizations. There is the possibility to specify reduced or unbalanced schemas as the cell is near the limits.
Extrapolation. The values are calculated through extrapolation from valid cells. If any of the two previous options is applied, this is the default behaviour.

The optional tag Boundary Schemas is used to set the second option. To do so, the user must provide the transformation rules and where they must be applied. For instance, for a Main Schema that takes 2 cells by each side we may define the following operator:

images/UserGuide/SpatialOperatorDiscretizationBoundaries.png

For each cell inmediately near the boundary (Distance To Boundary = 1) we set unbalanced schemas Forward Differences 4th Order and Backward Differences 4th Order, which only need information from the interior. For the points at distance 2 there is set a Centered Differences 2nd Order method, which only needs information from 1 near cell in each direction.

Notice that for distance = 2, the same schema is used in the Lower and Upper side of the boundary. Simflowny internally simplifies the code to avoid repeating calls to the same functions.

Creating new PDE Discretization Schemas

In order to create a discretization schema for PDE, from the document manager click on the Add button and select "PDE Discretization schema". The editor will open with an empty discretization schema.

Head

Head information is commonly used in several documents in Simflowny. See more detail in Head.

Parameters (Optional)

PDE Discretization Schemas can be parametrized. The Schema Parameters can be added from the root element. Each Parameter must have a name, a type and, optionaly, a default value. The default value is used to guide other user of this schema on a proper value for the parameter.

Schema

The Schema tag is the main tag in the document. It is organized in Step elements, which can be added as much as desired. For instance, if using a third order Runge-Kutta, three Step elements are set.

Every step can have several elements, most of them similar for mesh or particle discretization, but necessarily set separately.

Mesh/Particle Pre Step (optional)

This section allows the schema author to set a transformation rule with an arbitrary Simml algorithm that will be executed at the beginning of the step. The transformation rule must be set, as well as the input and output variables.

Spatial Discretization

Spatial discretization can be set for conservative terms (Conservative Term Discretization) and non conservative terms (Non Conservative Term Discretization), both optional.

For conservative discretization, the Flux Variables and Discretization tags must be set. Flux Variables tag is used to set the variable wildcards for input and output for the flux formulae. For instance, as seen in the image below, for a flux variable with input $f_p and output Flux$i$f, imagine that the problem has fields A and B in a 2D space, and that both fields have fluxes for both coordinates. The code generated would be equivalent at something like:

\begin{equation*} FluxiA = F^x_A(A^n) FluxjA = F^y_A(A^n) FluxiB = F^x_B(B^n) FluxjB = F^y_B(B^n) \end{equation*}

where $F^x_A$ is the flux formulae for field A at coordinate x.

As can be seen from the image above, outputs of Flux Variables can (and use to) be used as inputs for the Discretization tag.

In Discretization tag, the transformation rule and its inputs must be specified.

Regarding non-conservative discretization. A variable wildcard must be set in Input Variable tag. This wildcard variable used to be the output from the last step (or variable at previous time for the first step). Then, the derivative input variables are defined for non-conservative spatial discretization.

Notice that for particle discretization, there is also needed normalization information for both conservative and non-conservative terms. The normalization information to fill is similar the explained previously. That normalization is needed for stability purposes, and consists on a calculation over the same particle neighbours that the discretization. The normalization formula depend strongly on the discretization formula used.

The last element to set in the Spatial Discretization is the Output variable. The value will represent the wildcard variable that is going to be used as input in the Time Integration elements. Defining this variable allows the user to use a variable out of the scope of the current Step if desired. The variable may use Simml infix tags. The document SimML: language reference details the use of the Simml tags.

images/UserGuide/schemaOutputVariable.png

From the picture above, it can be seen the output variable RHS_$f being the output of the spatial discretization and being used as one of the inputs of the Runge-Kutta Time Integration. $f represents the wildcard of a field, so RHS_$f stands for variables with a prefix RHS_. It could generate RHS_alpha, RHS_phi, etc. Anyway, for further detail go to SimML: language reference.

Notice that the user does not choose the specific spatial discretization. It is set in the discretization policy at discretization time. So the user can define a discretization schema independently from the spatial discretization, what is the objective of the Method of Lines.

Mesh/Particle Dissipation (Optional)

A Step can have an optional element Dissipation. Adding this tag will add a dissipation to the output variable from spatial discretization. The transformation rule which performs the dissipation must be set.

Particle Kernel Discretization (Optional)

This section is used to set kernel and kernel gradient functions for particle discretizations. When using $k or $kg simml tags within the discretization schema, they will be replaced by the kernel and kernel gradient calls.

Particle Interpolation Methods (Optional)

Particle methods use interpolation to get averages for values from neighbour particles. This section allows to set the interpolation and normalization methods to do that.

Particle Advective Term (Optional)

Evolution equation reference frame could be either Eulerian or Lagrangian. When discretizing using particles, Eulerian evolution equations need advective terms. Simflowny will automatically add the following advective term to the eulerian field:

\begin{equation*} v_i * \partial_i (u) \end{equation*}

being $v_i$ the velocity and $u$ the eulerian field.

In order to discretize this terms, the schema must specify the derivative and normalization functions and their respective input variables.

Time Discretization

The Time Discretization element is used to set the time discretization functions. Currently, the user must introduce the explicit discretization, and optionally an implicit discretization (to solve models with stiff terms) The explicit discretization is set through referencing the proper transformation rule for the Step. The input variables is the next requirement. Here, the user must be skilled enough to use variables already available at that point. For instance, output variables from spatial discretization from current and previous Step, output varaibles from Time Discretization from previous steps, schema parameters, and the initial variable.

The initial variable is a wildcard variable which is the initial input variable at the begining of every evolution step during the simulation. For the first evolution step, those variables contain the values from initialization, from following evolution steps they contain the values from the output of the last Step in the discretization schema. The Math element representing this variable is f_p (or f_p($cc) to have the value at the current cell).

Regarding the implicit discretization, the user must specify a name for the stiff right hand side variable and the value to the stiff coefficient, which is the diagonal coefficient of the current implicit step.

If particles with velocity are to be discretized, Position Input Variables tag must be defined. The particle positions are treated as evolution fields, so they are discretized using the same time discretization transformation rule.

The last element to set in a Time Discretization is the Output Variable. Similarly to the one in Spatial Discretization, the variable can be used in further steps. The last Step must always have the value f_n($cc), which stands for $u{n+1}(i)$.

Mesh/Particle Post Step (Optional)

Similar to Pre Step, but the code will be executed at the end of the current step.

Timestep Factor

Every step must indicate the time in which the step is at finalization. This is specially important when using the Berger Oliger algorithm in multilevel simulations. The factor value is related to a complete time increment (dt).

Compact Support Ratio (Optional)

Schemas using particles must set this element. The tag defines the compact support ratio for the interpolation kernel method. That compact support ratio multiplied to the influence radius (known at runtime) will indicate the maximum region of interest of a particle.

Time interpolation (Optional)

In multilevel simulations, Simflowny accepts Berger-Oliger AMR algorithm. Different time integration schemes may need specific time interpolations. The Time Interpolation element can be set to configure the transformation rule where the time interpolation algorithm is set.

If no time interpolation is provided, there is a simple linear default algorithm.

Wiki

Simflowny / UserGuide_Schema

Creating new discretization schemas for PDE