Local solver for the nonlinear eigenvalue problem induced by the diffusion equation.
More...
|
| DiffusionEigs (const constructor_value_type &tau=1.) |
| Constructor for local solver. More...
|
|
template<typename hyEdgeT , typename SmallMatInT , typename SmallMatOutT > |
SmallMatOutT & | trace_to_flux (const SmallMatInT &lambda_values_in, SmallMatOutT &lambda_values_out, hyEdgeT &hyper_edge, const lSol_float_t eig=0.) const |
| Evaluate local contribution to matrix–vector multiplication. More...
|
|
template<class hyEdgeT > |
std::array< unsigned int, 2 *hyEdge_dimT > | node_types (hyEdgeT &hyper_edge) const |
| Fill array with 1 if the node is of Dirichlet type and 0 elsewhen. More...
|
|
template<class hyEdgeT , typename SmallMatT > |
SmallMatT & | make_initial (SmallMatT &lambda_values, hyEdgeT &hyper_edge, const lSol_float_t time=0.) const |
| L2 project given analytical functio to skeletal space. More...
|
|
template<class hyEdgeT , typename SmallMatInT , typename SmallMatOutT > |
SmallMatOutT & | jacobian_of_trace_to_flux (const SmallMatInT &lambda_values_in, SmallMatOutT &lambda_values_out, const lSol_float_t eig, const SmallMatInT &lambda_vals, const lSol_float_t eig_val, hyEdgeT &hyper_edge) const |
| Evaluate local contribution to matrix–vector multiplication. More...
|
|
template<class hyEdgeT > |
std::array< lSol_float_t, 1U > | errors (const std::array< std::array< lSol_float_t, n_shape_bdr_ >, 2 *hyEdge_dimT > &lambda_values, hyEdgeT &hy_edge, const lSol_float_t eig=0.) const |
| Local squared contribution to the L2 error. More...
|
|
template<typename abscissa_float_t , std::size_t abscissas_sizeT, class input_array_t , class hyEdgeT > |
std::array< std::array< lSol_float_t, Hypercube< hyEdge_dimT >::pow(abscissas_sizeT)>, system_dim > | bulk_values (const std::array< abscissa_float_t, abscissas_sizeT > &abscissas, const input_array_t &lambda_values, hyEdgeT &hyper_edge, const lSol_float_t time=0.) const |
| Evaluate local local reconstruction at tensorial products of abscissas. More...
|
|
|
template<typename hyEdgeT > |
SmallSquareMat< n_loc_dofs_, lSol_float_t > | assemble_loc_matrix (const lSol_float_t tau, hyEdgeT &hyper_edge, const lSol_float_t eig) const |
| Assemble local matrix for the local solver. More...
|
|
template<typename hyEdgeT , typename SmallMatT > |
SmallVec< n_loc_dofs_, lSol_float_t > | assemble_rhs_from_lambda (const SmallMatT &lambda_values, hyEdgeT &hyper_edge) const |
| Assemble local right-hand for the local solver (from skeletal). More...
|
|
template<typename hyEdgeT , typename SmallMatT > |
std::array< lSol_float_t, n_loc_dofs_ > | solve_local_problem (const SmallMatT &lambda_values, hyEdgeT &hyper_edge, const lSol_float_t eig) const |
| Solve local problem (with right-hand side from skeletal). More...
|
|
template<typename hyEdgeT > |
std::array< std::array< lSol_float_t, n_shape_bdr_ >, 2 *hyEdge_dimT > | primal_at_boundary (const std::array< lSol_float_t, n_loc_dofs_ > &coeffs, hyEdgeT &hyper_edge) const |
| Evaluate primal variable at boundary. More...
|
|
template<typename hyEdgeT > |
std::array< std::array< lSol_float_t, n_shape_bdr_ >, 2 *hyEdge_dimT > | dual_at_boundary (const std::array< lSol_float_t,(hyEdge_dimT+1) *n_shape_fct_ > &coeffs, hyEdgeT &hyper_edge) const |
| Evaluate dual variable at boundary. More...
|
|
template<unsigned int hyEdge_dimT, unsigned int poly_deg, unsigned int quad_deg, template< unsigned int, typename > typename parametersT = DiffusionParametersDefault, typename lSol_float_t = double>
class LocalSolver::DiffusionEigs< hyEdge_dimT, poly_deg, quad_deg, parametersT, lSol_float_t >
Local solver for the nonlinear eigenvalue problem induced by the diffusion equation.
- Template Parameters
-
hyEdge_dimT | Dimension of a hyperedge, i.e., 1 is for PDEs defined on graphs, 2 is for PDEs defined on surfaces, and 3 is for PDEs defined on volumes. |
space_dimT | The dimension of the surrounding space. |
poly_deg | The polynomial degree of test and trial functions. |
quad_deg | The order of the quadrature rule. |
parametersT | Struct depending on templates space_dimTP and lSol_float_TP that contains static parameter functions. Defaults to above functions included in DiffusionParametersDefault . |
lSol_float_t | The floating point type calculations are executed in. Defaults to double. |
- Authors
- Guido Kanschat, Heidelberg University, 2019–2020.
-
Andreas Rupp, Heidelberg University, 2019–2020.
template<unsigned int hyEdge_dimT, unsigned int poly_deg, unsigned int quad_deg, template< unsigned int, typename > typename parametersT = DiffusionParametersDefault, typename lSol_float_t = double>
template<typename hyEdgeT >
Assemble local matrix for the local solver.
- Template Parameters
-
hyEdgeT | The geometry type / typename of the considered hyEdge's geometry. |
- Parameters
-
tau | Penalty parameter. |
hyper_edge | The geometry of the considered hyperedge (of typename GeomT). |
eig | Eigenvalue for which the matrix is assembled. |
- Return values
-
loc_mat | Matrix of the local solver. |
template<unsigned int hyEdge_dimT, unsigned int poly_deg, unsigned int quad_deg, template< unsigned int, typename > typename parametersT = DiffusionParametersDefault, typename lSol_float_t = double>
template<typename hyEdgeT , typename SmallMatT >
Assemble local right-hand for the local solver (from skeletal).
The right hand side needs the values of the global degrees of freedom. Note that we basically follow the lines of
B. Cockburn, J. Gopalakrishnan, and R. Lazarov. Unified hybridization of discontinuous Galerkin, mixed, and continuous Galerkin methods for second order elliptic problems. SIAM Journal on Numerical Analysis, 47(2):1319–1365, 2009
and discriminate between local solvers with respect to the skeletal variable and with respect to the global right-hand side. This assembles the local right-hand side with respect to the skeletal variable.
- Template Parameters
-
hyEdgeT | The geometry type / typename of the considered hyEdge's geometry. |
SmallMatT | The data type of the lambda_values . |
- Parameters
-
lambda_values | Global degrees of freedom associated to the hyperedge. |
hyper_edge | The geometry of the considered hyperedge (of typename GeomT). |
- Return values
-
loc_rhs | Local right hand side of the locasl solver. |
template<unsigned int hyEdge_dimT, unsigned int poly_deg, unsigned int quad_deg, template< unsigned int, typename > typename parametersT = DiffusionParametersDefault, typename lSol_float_t = double>
template<typename abscissa_float_t , std::size_t abscissas_sizeT, class input_array_t , class hyEdgeT >
std::array<std::array<lSol_float_t, Hypercube<hyEdge_dimT>::pow(abscissas_sizeT)>, system_dim> LocalSolver::DiffusionEigs< hyEdge_dimT, poly_deg, quad_deg, parametersT, lSol_float_t >::bulk_values |
( |
const std::array< abscissa_float_t, abscissas_sizeT > & |
abscissas, |
|
|
const input_array_t & |
lambda_values, |
|
|
hyEdgeT & |
hyper_edge, |
|
|
const lSol_float_t |
time = 0. |
|
) |
| const |
Evaluate local local reconstruction at tensorial products of abscissas.
- Template Parameters
-
absc_float_t | Floating type for the abscissa values. |
abscissas_sizeT | Size of the array of array of abscissas. |
input_array_t | Type of input array. |
hyEdgeT | The geometry type / typename of the considered hyEdge's geometry. |
- Parameters
-
abscissas | Abscissas of the supporting points. |
lambda_values | The values of the skeletal variable's coefficients. |
hyper_edge | The geometry of the considered hyperedge (of typename GeomT). |
time | Time. |
- Return values
-
func_values | Function values at tensorial points. |
template<unsigned int hyEdge_dimT, unsigned int poly_deg, unsigned int quad_deg, template< unsigned int, typename > typename parametersT = DiffusionParametersDefault, typename lSol_float_t = double>
template<typename hyEdgeT >
std::array<std::array<lSol_float_t, n_shape_bdr_>, 2 * hyEdge_dimT> LocalSolver::DiffusionEigs< hyEdge_dimT, poly_deg, quad_deg, parametersT, lSol_float_t >::dual_at_boundary |
( |
const std::array< lSol_float_t,(hyEdge_dimT+1) *n_shape_fct_ > & |
coeffs, |
|
|
hyEdgeT & |
hyper_edge |
|
) |
| const |
|
inlineprivate |
Evaluate dual variable at boundary.
Function to evaluate dual variable of the solution. This function is needed to calculate the local numerical fluxes.
- Template Parameters
-
hyEdgeT | The geometry type / typename of the considered hyEdge's geometry. |
- Parameters
-
coeffs | Coefficients of the local solution. |
hyper_edge | The geometry of the considered hyperedge (of typename GeomT). |
- Return values
-
bdr_coeffs | Coefficients of respective (dim-1) dimensional function at boundaries. |
template<unsigned int hyEdge_dimT, unsigned int poly_deg, unsigned int quad_deg, template< unsigned int, typename > typename parametersT = DiffusionParametersDefault, typename lSol_float_t = double>
template<class hyEdgeT >
std::array<lSol_float_t, 1U> LocalSolver::DiffusionEigs< hyEdge_dimT, poly_deg, quad_deg, parametersT, lSol_float_t >::errors |
( |
const std::array< std::array< lSol_float_t, n_shape_bdr_ >, 2 *hyEdge_dimT > & |
lambda_values, |
|
|
hyEdgeT & |
hy_edge, |
|
|
const lSol_float_t |
eig = 0. |
|
) |
| const |
|
inline |
Local squared contribution to the L2 error.
- Template Parameters
-
hyEdgeT | The geometry type / typename of the considered hyEdge's geometry. |
- Parameters
-
lambda_values | The values of the skeletal variable's coefficients. |
hy_edge | The geometry of the considered hyperedge (of typename GeomT). |
eig | Approximated eigenvalue. |
- Return values
-
vec_b | Local part of vector b. |
template<unsigned int hyEdge_dimT, unsigned int poly_deg, unsigned int quad_deg, template< unsigned int, typename > typename parametersT = DiffusionParametersDefault, typename lSol_float_t = double>
static constexpr unsigned int LocalSolver::DiffusionEigs< hyEdge_dimT, poly_deg, quad_deg, parametersT, lSol_float_t >::hyEdge_dim |
( |
| ) |
|
|
inlinestaticconstexpr |
Dimension of hyper edge type that this object solves on.
template<unsigned int hyEdge_dimT, unsigned int poly_deg, unsigned int quad_deg, template< unsigned int, typename > typename parametersT = DiffusionParametersDefault, typename lSol_float_t = double>
template<typename parameters >
static constexpr bool LocalSolver::DiffusionEigs< hyEdge_dimT, poly_deg, quad_deg, parametersT, lSol_float_t >::is_dirichlet |
( |
const unsigned int |
node_type | ) |
|
|
inlinestaticconstexprprivate |
Find out whether a node is of Dirichlet type.
template<unsigned int hyEdge_dimT, unsigned int poly_deg, unsigned int quad_deg, template< unsigned int, typename > typename parametersT = DiffusionParametersDefault, typename lSol_float_t = double>
template<class hyEdgeT , typename SmallMatInT , typename SmallMatOutT >
SmallMatOutT& LocalSolver::DiffusionEigs< hyEdge_dimT, poly_deg, quad_deg, parametersT, lSol_float_t >::jacobian_of_trace_to_flux |
( |
const SmallMatInT & |
lambda_values_in, |
|
|
SmallMatOutT & |
lambda_values_out, |
|
|
const lSol_float_t |
eig, |
|
|
const SmallMatInT & |
lambda_vals, |
|
|
const lSol_float_t |
eig_val, |
|
|
hyEdgeT & |
hyper_edge |
|
) |
| const |
|
inline |
Evaluate local contribution to matrix–vector multiplication.
This corresponds to the evaluation of the Jacobian evaluated at lambda_vals
, eig_val
in the direction lambda_values
, eig
.
- Template Parameters
-
hyEdgeT | The geometry type / typename of the considered hyEdge's geometry. |
SmallMatInT | Data type of lambda_values_in . |
SmallMatOutT | Data type of lambda_values_out |
- Parameters
-
lambda_values_in | Local part of vector x. |
lambda_values_out | Local part that will be added to A * x. |
eig | Eigenvalue in whose direction Jacobian is evaluated. |
lambda_vals | Lambda at which Jacobian is evaluated. |
eig_val | Eigenvalue at which Jacobian is evaluated. |
hyper_edge | The geometry of the considered hyperedge (of typename hyEdgeT). |
- Return values
-
vecAx | Local part of vector A * x. |
template<unsigned int hyEdge_dimT, unsigned int poly_deg, unsigned int quad_deg, template< unsigned int, typename > typename parametersT = DiffusionParametersDefault, typename lSol_float_t = double>
template<class hyEdgeT , typename SmallMatT >
SmallMatT& LocalSolver::DiffusionEigs< hyEdge_dimT, poly_deg, quad_deg, parametersT, lSol_float_t >::make_initial |
( |
SmallMatT & |
lambda_values, |
|
|
hyEdgeT & |
hyper_edge, |
|
|
const lSol_float_t |
time = 0. |
|
) |
| const |
|
inline |
L2 project given analytical functio to skeletal space.
- Template Parameters
-
hyEdgeT | The geometry type / typename of the considered hyEdge's geometry. |
SmallMatT | The data tyepe of lambda_values . |
- Parameters
-
lambda_values | Local part of vector x. (Redundant) |
hyper_edge | The geometry of the considered hyperedge (of typename hyEdgeT). |
time | Time. (Redundant) |
- Return values
-
bdr_values | Coefficients of L2 projection. |
template<unsigned int hyEdge_dimT, unsigned int poly_deg, unsigned int quad_deg, template< unsigned int, typename > typename parametersT = DiffusionParametersDefault, typename lSol_float_t = double>
static constexpr unsigned int LocalSolver::DiffusionEigs< hyEdge_dimT, poly_deg, quad_deg, parametersT, lSol_float_t >::n_glob_dofs_per_node |
( |
| ) |
|
|
inlinestaticconstexpr |
template<unsigned int hyEdge_dimT, unsigned int poly_deg, unsigned int quad_deg, template< unsigned int, typename > typename parametersT = DiffusionParametersDefault, typename lSol_float_t = double>
static constexpr unsigned int LocalSolver::DiffusionEigs< hyEdge_dimT, poly_deg, quad_deg, parametersT, lSol_float_t >::node_system_dimension |
( |
| ) |
|
|
inlinestaticconstexpr |
Dimension of of the solution evaluated with respect to a hypernode.
template<unsigned int hyEdge_dimT, unsigned int poly_deg, unsigned int quad_deg, template< unsigned int, typename > typename parametersT = DiffusionParametersDefault, typename lSol_float_t = double>
template<class hyEdgeT >
std::array<unsigned int, 2 * hyEdge_dimT> LocalSolver::DiffusionEigs< hyEdge_dimT, poly_deg, quad_deg, parametersT, lSol_float_t >::node_types |
( |
hyEdgeT & |
hyper_edge | ) |
const |
|
inline |
Fill array with 1 if the node is of Dirichlet type and 0 elsewhen.
- Template Parameters
-
hyEdgeT | The geometry type / typename of the considered hyEdge's geometry. |
- Parameters
-
hyper_edge | The geometry of the considered hyperedge (of typename hyEdgeT). |
- Return values
-
result | Array encoding edge types. |
template<unsigned int hyEdge_dimT, unsigned int poly_deg, unsigned int quad_deg, template< unsigned int, typename > typename parametersT = DiffusionParametersDefault, typename lSol_float_t = double>
template<typename hyEdgeT >
std::array<std::array<lSol_float_t, n_shape_bdr_>, 2 * hyEdge_dimT> LocalSolver::DiffusionEigs< hyEdge_dimT, poly_deg, quad_deg, parametersT, lSol_float_t >::primal_at_boundary |
( |
const std::array< lSol_float_t, n_loc_dofs_ > & |
coeffs, |
|
|
hyEdgeT & |
hyper_edge |
|
) |
| const |
|
inlineprivate |
Evaluate primal variable at boundary.
Function to evaluate primal variable of the solution. This function is needed to calculate the local numerical fluxes.
- Template Parameters
-
hyEdgeT | The geometry type / typename of the considered hyEdge's geometry. |
- Parameters
-
coeffs | Coefficients of the local solution. |
hyper_edge | The geometry of the considered hyperedge (of typename GeomT). |
- Return values
-
bdr_coeffs | Coefficients of respective (dim-1) dimensional function at boundaries. |
template<unsigned int hyEdge_dimT, unsigned int poly_deg, unsigned int quad_deg, template< unsigned int, typename > typename parametersT = DiffusionParametersDefault, typename lSol_float_t = double>
template<typename hyEdgeT , typename SmallMatT >
std::array<lSol_float_t, n_loc_dofs_> LocalSolver::DiffusionEigs< hyEdge_dimT, poly_deg, quad_deg, parametersT, lSol_float_t >::solve_local_problem |
( |
const SmallMatT & |
lambda_values, |
|
|
hyEdgeT & |
hyper_edge, |
|
|
const lSol_float_t |
eig |
|
) |
| const |
|
inlineprivate |
Solve local problem (with right-hand side from skeletal).
- Template Parameters
-
hyEdgeT | The geometry type / typename of the considered hyEdge's geometry. |
SmallMatT | The data type of the lambda_values . |
- Parameters
-
lambda_values | Global degrees of freedom associated to the hyperedge. |
hyper_edge | The geometry of the considered hyperedge (of typename GeomT). |
eig | Eigenvalue for which the local problem is solved. |
- Return values
-
loc_sol | Solution of the local problem. |
template<unsigned int hyEdge_dimT, unsigned int poly_deg, unsigned int quad_deg, template< unsigned int, typename > typename parametersT = DiffusionParametersDefault, typename lSol_float_t = double>
static constexpr unsigned int LocalSolver::DiffusionEigs< hyEdge_dimT, poly_deg, quad_deg, parametersT, lSol_float_t >::system_dimension |
( |
| ) |
|
|
inlinestaticconstexpr |
Dimension of of the solution evaluated with respect to a hyperedge.
template<unsigned int hyEdge_dimT, unsigned int poly_deg, unsigned int quad_deg, template< unsigned int, typename > typename parametersT = DiffusionParametersDefault, typename lSol_float_t = double>
template<typename hyEdgeT , typename SmallMatInT , typename SmallMatOutT >
SmallMatOutT& LocalSolver::DiffusionEigs< hyEdge_dimT, poly_deg, quad_deg, parametersT, lSol_float_t >::trace_to_flux |
( |
const SmallMatInT & |
lambda_values_in, |
|
|
SmallMatOutT & |
lambda_values_out, |
|
|
hyEdgeT & |
hyper_edge, |
|
|
const lSol_float_t |
eig = 0. |
|
) |
| const |
|
inline |
Evaluate local contribution to matrix–vector multiplication.
This corresponds to an evaluation of the residual in the nonlinear context
- Template Parameters
-
hyEdgeT | The geometry type / typename of the considered hyEdge's geometry. |
SmallMatInT | Data type of lambda_values_in . |
SmallMatOutT | Data type of lambda_values_out |
- Parameters
-
lambda_values_in | Local part of vector x. |
lambda_values_out | Local part that will be added to A * x. |
hyper_edge | The geometry of the considered hyperedge (of typename hyEdgeT). |
eig | Eigenvalue. |
- Return values
-
vecAx | Local part of vector A * x. |
template<unsigned int hyEdge_dimT, unsigned int poly_deg, unsigned int quad_deg, template< unsigned int, typename > typename parametersT = DiffusionParametersDefault, typename lSol_float_t = double>
Dimension of of the solution evaluated with respect to a hypernode.
This allows to the use of this quantity as template parameter in member functions.
template<unsigned int hyEdge_dimT, unsigned int poly_deg, unsigned int quad_deg, template< unsigned int, typename > typename parametersT = DiffusionParametersDefault, typename lSol_float_t = double>
Dimension of of the solution evaluated with respect to a hypernode.
This allows to the use of this quantity as template parameter in member functions.
template<unsigned int hyEdge_dimT, unsigned int poly_deg, unsigned int quad_deg, template< unsigned int, typename > typename parametersT = DiffusionParametersDefault, typename lSol_float_t = double>
(Globally constant) penalty parameter for HDG scheme.