Coefficients.
Data Types | |
| type | coef_t |
| Coefficients defined on a given (mesh, \( X_h \)) tuple. Arrays use indices (i,j,k,e): element e, local coordinate (i,j,k). More... | |
Functions/Subroutines | |
| subroutine | coef_init_empty (this, xh, msh) |
| Initialize empty coefs for a space and a mesh. | |
| subroutine | coef_init_all (this, gs_h, scope) |
| Initialize coefficients. | |
| subroutine | coef_free (this) |
| Deallocate coefficients. | |
| subroutine | coef_release_scratch (this) |
| Release the coefficients that COEF_OPERATOR does not retain. | |
| subroutine | coef_generate_dxyzdrst (c) |
| subroutine | coef_generate_geo (c) |
| Generate geometric data for the given mesh. | |
| subroutine | coef_metric_condition (this) |
| Compute the largest condition number of the metric tensor over the mesh. | |
| subroutine | coef_generate_geo_compressed (c) |
| Compute processor-local compressed versions of mappings Gij. | |
| subroutine | coef_generate_mass (c) |
| Generate mass matrix B for the given mesh and space. | |
| subroutine | coef_require_facets (this, who) |
| Abort unless this coef holds the facet areas and normals. | |
| pure real(kind=rp) function, dimension(3) | coef_get_normal (this, i, j, k, e, facet) |
| Facet normal at a point. | |
| pure real(kind=rp) function | coef_get_area (this, i, j, k, e, facet) |
| Facet area at a point. | |
| subroutine | coef_generate_area_and_normal (coef) |
| Generate facet area and surface normals. | |
| subroutine | coef_generate_cyclic_bc (this) |
| subroutine | coef_recompute_metrics (this) |
| Recompute and update geometric factors (ALE) | |
| subroutine | coef_enable_lagged_mass (this) |
| Enable separate memory for lagged B matrices if needed. For eg. when mesh moves. | |
| subroutine | coef_update_lagged_mass (this) |
| Update history: Blaglag = Blag, Blag = B. | |
Variables | |
| real(kind=rp), parameter, public | neko_metric_cond_sp = 1.0e4_rp |
| Largest metric condition number for which single precision storage of the geometric factors \( G_{ij} \) is considered safe. | |
| integer, parameter, public | coef_full = 0 |
| Retain every coefficient. The default, and the only scope that supports recompute_metrics(), generate_cyclic_bc(), get_area(), get_normal() and any consumer of jac, jacinv, Binv or the derivative arrays. | |
| integer, parameter, public | coef_operator = 1 |
| Retain only what applying a discrete operator needs: \( G_{ij} \), h1, h2, B and mult. | |
Function/Subroutine Documentation
◆ coef_enable_lagged_mass()
◆ coef_free()
◆ coef_generate_area_and_normal()
|
private |
◆ coef_generate_cyclic_bc()
|
private |
◆ coef_generate_dxyzdrst()
|
private |
◆ coef_generate_geo()
|
private |
◆ coef_generate_geo_compressed()
|
private |
◆ coef_generate_mass()
|
private |
◆ coef_get_area()
|
private |
- Note
- Pure, see the note on coef_get_normal().
◆ coef_get_normal()
|
private |
- Note
- Pure, so it cannot reject a coef_t that lacks the normals. Consumers assert that themselves at setup, see coef_require_facets().
◆ coef_init_all()
|
private |
- Parameters
-
gs_h Gather-scatter handle carrying the dofmap to build on. scope Which coefficients to retain, COEF_FULL (default) or COEF_OPERATOR. See the scope parameters for what each one keeps.
- Todo:
- Be clever and try to avoid allocating zeroed geom. factors
Default values correspond to no rotation
Definition at line 300 of file coef.f90.
◆ coef_init_empty()
◆ coef_metric_condition()
|
private |
The stored geometric factors are \( G = w_3 J M \), where \( M \) is the metric matrix at a quadrature point. Both \( w_3 \) and the Jacobian are scalar multipliers, so \( \kappa(G) = \kappa(M) \) and the condition number of the stored \( G_{ij} \) is a pure property of the element geometry.
It matters when \( G_{ij} \) is stored in a lower precision than it was computed in. Rounding each entry perturbs it by the unit roundoff \( \epsilon \), hence the eigenvalues by up to \( \epsilon \lambda_{max} \), so the rounded tensor stays positive definite – and the Helmholtz operator therefore stays a valid SPD discretisation – only while \( \kappa(G) \ll 1/\epsilon \): about \( 1.7\times 10^{7} \) in single precision against \( 9\times 10^{15} \) in double. Below that limit the rounding is a \( \sim\epsilon \) perturbation of the diffusion coefficient and the solution moves by the same order; past it the local operator can turn indefinite and a Krylov solve breaks down rather than degrading.
\( \kappa \) grows as the square of the element aspect ratio, and skew compounds it, so a wall resolved boundary layer can come within a factor of ten of the single precision limit while remaining ten orders clear of the double precision one.
- Note
- On device builds the host copy of \( G_{ij} \) is only refreshed at initialization, so this refreshes it itself when called later. It is not called from recompute_metrics(); a moving mesh that deforms significantly should invoke it explicitly.
Definition at line 1288 of file coef.f90.
◆ coef_recompute_metrics()
|
private |
◆ coef_release_scratch()
The derivative arrays, the Jacobian and its inverse are scratch for \( G_{ij} \) and B, and are dead once those exist. Binv, the facet areas and the normals are never read when applying an operator. Called at the end of init; a no-op under COEF_FULL.
- Note
- area and the normals are not allocated at all under COEF_OPERATOR. They are listed here so the set of released coefficients is stated in one place.
Definition at line 807 of file coef.f90.
◆ coef_require_facets()
|
private |
For consumers of the facet metrics to call once at setup, so that being handed a COEF_OPERATOR coef fails at construction with a message naming the consumer, rather than dereferencing a released array mid solve. get_area() and get_normal() are pure and cannot report it themselves, and some consumers read nx, ny and nz directly and never go through them at all.
- Parameters
-
who Name of the consumer, used in the error message.
◆ coef_update_lagged_mass()
Variable Documentation
◆ coef_full
◆ coef_operator
The derivative arrays, jac and jacinv are scratch for \( G_{ij} \) and B, and Binv, area and the facet normals are never read by an operator, so all of them are released once the metrics are built. The work feeding only those is skipped as well: the facet metrics, the gather-scatter behind Binv, the volume reduction and the metric condition estimate.
Intended for multigrid levels, which read nothing else – see hsmg_init() and phmg_init(). The geometry cannot be rebuilt without the derivative arrays, so this scope is incompatible with a moving mesh.
◆ neko_metric_cond_sp
The hard limit is \( 1/\epsilon_{sp} \approx 1.7\times 10^{7} \), past which rounding swamps the smallest eigenvalue of the metric and the element operator can lose positive definiteness. This keeps three orders of margin, which accepts isotropic, graded and channel interior meshes and rejects wall resolved ones. See coef_metric_condition().
