Implements a discontinuous Galerkin solver for (single-phase) time-of-flight using reordering. More...
#include <TofDiscGalReorder.hpp>
Public Member Functions | |
TofDiscGalReorder (const UnstructuredGrid &grid, const ParameterGroup ¶m) | |
Construct solver. More... | |
void | solveTof (const double *darcyflux, const double *porevolume, const double *source, std::vector< double > &tof_coeff) |
Solve for time-of-flight. More... | |
void | solveTofTracer (const double *darcyflux, const double *porevolume, const double *source, const SparseTable< int > &tracerheads, std::vector< double > &tof_coeff, std::vector< double > &tracer_coeff) |
Solve for time-of-flight and a number of tracers. More... | |
Additional Inherited Members | |
![]() | |
void | reorderAndTransport (const UnstructuredGrid &grid, const double *darcyflux) |
const std::vector< int > & | sequence () const |
const std::vector< int > & | components () const |
Implements a discontinuous Galerkin solver for (single-phase) time-of-flight using reordering.
The equation solved is:
in which is the fluid velocity,
is time-of-flight and
is the porosity. This is a boundary value problem, and
is specified to be zero on all inflow boundaries. The user may specify the polynomial degree of the basis function space used, but only degrees 0 and 1 are supported so far.
Opm::TofDiscGalReorder::TofDiscGalReorder | ( | const UnstructuredGrid & | grid, |
const ParameterGroup & | param | ||
) |
Construct solver.
[in] | grid | A 2d or 3d grid. |
[in] | param | Parameters for the solver. The following parameters are accepted (defaults):
|
void Opm::TofDiscGalReorder::solveTof | ( | const double * | darcyflux, |
const double * | porevolume, | ||
const double * | source, | ||
std::vector< double > & | tof_coeff | ||
) |
Solve for time-of-flight.
[in] | darcyflux | Array of signed face fluxes. |
[in] | porevolume | Array of pore volumes. |
[in] | source | Source term. Sign convention is: (+) inflow flux, (-) outflow flux. |
[out] | tof_coeff | Array of time-of-flight solution coefficients. The values are ordered by cell, meaning that the K coefficients corresponding to the first cell come before the K coefficients corresponding to the second cell etc. K depends on degree and grid dimension. |
void Opm::TofDiscGalReorder::solveTofTracer | ( | const double * | darcyflux, |
const double * | porevolume, | ||
const double * | source, | ||
const SparseTable< int > & | tracerheads, | ||
std::vector< double > & | tof_coeff, | ||
std::vector< double > & | tracer_coeff | ||
) |
Solve for time-of-flight and a number of tracers.
[in] | darcyflux | Array of signed face fluxes. |
[in] | porevolume | Array of pore volumes. |
[in] | source | Source term. Sign convention is: (+) inflow flux, (-) outflow flux. |
[in] | tracerheads | Table containing one row per tracer, and each row contains the source cells for that tracer. |
[out] | tof_coeff | Array of time-of-flight solution coefficients. The values are ordered by cell, meaning that the K coefficients corresponding to the first cell comes before the K coefficients corresponding to the second cell etc. K depends on degree and grid dimension. |
[out] | tracer_coeff | Array of tracer solution coefficients. N*K per cell, where N is equal to tracerheads.size(). All K coefs for a tracer are consecutive, and all tracers' coefs for a cell come before those for the next cell. |