![]() |
MFEM v4.7.0
Finite element discretization library
|
A TMOP integrator class based on any given TMOP_QualityMetric and TargetConstructor. More...
#include <tmop.hpp>
Public Member Functions | |
TMOP_Integrator (TMOP_QualityMetric *m, TargetConstructor *tc, TMOP_QualityMetric *hm) | |
TMOP_Integrator (TMOP_QualityMetric *m, TargetConstructor *tc) | |
~TMOP_Integrator () | |
void | ReleasePADeviceMemory (bool copy_to_host=true) |
void | SetIntegrationRules (IntegrationRules &irules, int order) |
Prescribe a set of integration rules; relevant for mixed meshes. | |
void | IntegrateOverTarget (bool integ_over_target_) |
void | SetCoefficient (Coefficient &w1) |
Sets a scaling Coefficient for the quality metric term of the integrator. | |
void | EnableLimiting (const GridFunction &n0, const GridFunction &dist, Coefficient &w0, TMOP_LimiterFunction *lfunc=NULL) |
Limiting of the mesh displacements (general version). | |
void | EnableLimiting (const GridFunction &n0, Coefficient &w0, TMOP_LimiterFunction *lfunc=NULL) |
Adds a limiting term to the integrator with limiting distance function (dist in the general version of the method) equal to 1. | |
void | EnableAdaptiveLimiting (const GridFunction &z0, Coefficient &coeff, AdaptivityEvaluator &ae) |
Restriction of the node positions to certain regions. | |
void | EnableAdaptiveLimiting (const ParGridFunction &z0, Coefficient &coeff, AdaptivityEvaluator &ae) |
Parallel support for adaptive limiting. | |
void | EnableSurfaceFitting (const GridFunction &s0, const Array< bool > &smarker, Coefficient &coeff, AdaptivityEvaluator &ae) |
Fitting of certain DOFs to the zero level set of a function. | |
void | EnableSurfaceFitting (const ParGridFunction &s0, const Array< bool > &smarker, Coefficient &coeff, AdaptivityEvaluator &ae) |
Parallel support for surface fitting to the zero level set of a function. | |
void | EnableSurfaceFittingFromSource (const ParGridFunction &s_bg, ParGridFunction &s0, const Array< bool > &smarker, Coefficient &coeff, AdaptivityEvaluator &ae, const ParGridFunction &s_bg_grad, ParGridFunction &s0_grad, AdaptivityEvaluator &age, const ParGridFunction &s_bg_hess, ParGridFunction &s0_hess, AdaptivityEvaluator &ahe) |
Fitting of certain DOFs in the current mesh to the zero level set of a function defined on another (finer) source mesh. | |
void | EnableSurfaceFitting (const GridFunction &pos, const Array< bool > &smarker, Coefficient &coeff) |
Fitting of certain DOFs to given positions in physical space. | |
void | GetSurfaceFittingErrors (const Vector &pos, real_t &err_avg, real_t &err_max) |
bool | IsSurfaceFittingEnabled () |
void | SetLimitingNodes (const GridFunction &n0) |
Update the original/reference nodes used for limiting. | |
virtual real_t | GetElementEnergy (const FiniteElement &el, ElementTransformation &T, const Vector &elfun) |
Computes the integral of W(Jacobian(Trt)) over a target zone. | |
virtual real_t | GetRefinementElementEnergy (const FiniteElement &el, ElementTransformation &T, const Vector &elfun, const IntegrationRule &irule) |
Computes the mean of the energies of the given element's children. | |
virtual real_t | GetDerefinementElementEnergy (const FiniteElement &el, ElementTransformation &T, const Vector &elfun) |
virtual void | AssembleElementVector (const FiniteElement &el, ElementTransformation &T, const Vector &elfun, Vector &elvect) |
Perform the local action of the NonlinearFormIntegrator. | |
virtual void | AssembleElementGrad (const FiniteElement &el, ElementTransformation &T, const Vector &elfun, DenseMatrix &elmat) |
Assemble the local gradient matrix. | |
TMOP_QualityMetric & | GetAMRQualityMetric () |
void | UpdateAfterMeshTopologyChange () |
void | ParUpdateAfterMeshTopologyChange () |
virtual void | AssemblePA (const FiniteElementSpace &) |
Method defining partial assembly. | |
virtual void | AssembleGradPA (const Vector &, const FiniteElementSpace &) |
Prepare the integrator for partial assembly (PA) gradient evaluations on the given FE space fes at the state x. | |
virtual real_t | GetLocalStateEnergyPA (const Vector &) const |
Compute the local (to the MPI rank) energy with partial assembly. | |
virtual void | AddMultPA (const Vector &, Vector &) const |
Method for partially assembled action. | |
virtual void | AddMultGradPA (const Vector &, Vector &) const |
Method for partially assembled gradient action. | |
virtual void | AssembleGradDiagonalPA (Vector &) const |
Method for computing the diagonal of the gradient with partial assembly. | |
DiscreteAdaptTC * | GetDiscreteAdaptTC () const |
void | EnableNormalization (const GridFunction &x) |
Computes the normalization factors of the metric and limiting integrals using the mesh position given by x. | |
void | ParEnableNormalization (const ParGridFunction &x) |
void | EnableFiniteDifferences (const GridFunction &x) |
Enables FD-based approximation and computes dx. | |
void | EnableFiniteDifferences (const ParGridFunction &x) |
void | SetFDhScale (real_t dxscale_) |
bool | GetFDFlag () const |
real_t | GetFDh () const |
void | SetExactActionFlag (bool flag_) |
Flag to control if exact action of Integration is effected. | |
void | UpdateSurfaceFittingWeight (real_t factor) |
Update the surface fitting weight as surf_fit_coeff *= factor;. | |
real_t | GetSurfaceFittingWeight () |
Get the surface fitting weight. | |
void | ComputeUntangleMetricQuantiles (const Vector &x, const FiniteElementSpace &fes) |
virtual void | AssemblePA (const FiniteElementSpace &trial_fes, const FiniteElementSpace &test_fes) |
![]() | |
virtual void | SetIntRule (const IntegrationRule *ir) |
Prescribe a fixed IntegrationRule to use (when ir != NULL) or let the integrator choose (when ir == NULL). | |
void | SetIntegrationMode (Mode m) |
void | SetNURBSPatchIntRule (NURBSMeshRules *pr) |
For patchwise integration, SetNURBSPatchIntRule must be called. | |
bool | HasNURBSPatchIntRule () const |
bool | Patchwise () const |
void | SetIntegrationRule (const IntegrationRule &ir) |
Prescribe a fixed IntegrationRule to use. | |
void | SetPAMemoryType (MemoryType mt) |
const IntegrationRule * | GetIntegrationRule () const |
Get the integration rule of the integrator (possibly NULL). | |
virtual void | AssembleFaceVector (const FiniteElement &el1, const FiniteElement &el2, FaceElementTransformations &Tr, const Vector &elfun, Vector &elvect) |
Perform the local action of the NonlinearFormIntegrator resulting from a face integral term. | |
virtual void | AssembleFaceGrad (const FiniteElement &el1, const FiniteElement &el2, FaceElementTransformations &Tr, const Vector &elfun, DenseMatrix &elmat) |
Assemble the local action of the gradient of the NonlinearFormIntegrator resulting from a face integral term. | |
virtual bool | SupportsCeed () const |
Indicates whether this integrator can use a Ceed backend. | |
virtual void | AssembleMF (const FiniteElementSpace &fes) |
Method defining fully unassembled operator. | |
virtual void | AddMultMF (const Vector &x, Vector &y) const |
ceed::Operator & | GetCeedOp () |
virtual | ~NonlinearFormIntegrator () |
Friends | |
class | TMOPNewtonSolver |
class | TMOPComboIntegrator |
Additional Inherited Members | |
![]() | |
enum | Mode { ELEMENTWISE = 0 , PATCHWISE = 1 , PATCHWISE_REDUCED = 2 } |
A TMOP integrator class based on any given TMOP_QualityMetric and TargetConstructor.
Represents
|
inline |
[in] | m | TMOP_QualityMetric for r-adaptivity (not owned). |
[in] | tc | Target-matrix construction algorithm to use (not owned). |
[in] | hm | TMOP_QualityMetric for h-adaptivity (not owned). |
|
inline |
|
inlineprotected |
Method for partially assembled gradient action.
All arguments are E-vectors. This method can be called only after the method AssembleGradPA() has been called.
Reimplemented from mfem::NonlinearFormIntegrator.
Definition at line 337 of file tmop_pa.cpp.
Definition at line 114 of file tmop_pa_h2m.cpp.
Definition at line 117 of file tmop_pa_h3m.cpp.
Definition at line 94 of file tmop_pa_h2m_c0.cpp.
Definition at line 98 of file tmop_pa_h3m_c0.cpp.
Method for partially assembled action.
Perform the action of integrator on the input x and add the result to the output y. Both x and y are E-vectors, i.e. they represent the element-wise discontinuous version of the FE space.
This method can be called only after the method AssemblePA() has been called.
Reimplemented from mfem::NonlinearFormIntegrator.
Definition at line 315 of file tmop_pa.cpp.
Definition at line 203 of file tmop_pa_p2.cpp.
Definition at line 239 of file tmop_pa_p3.cpp.
Definition at line 152 of file tmop_pa_p2_c0.cpp.
Definition at line 160 of file tmop_pa_p3_c0.cpp.
|
protected |
Definition at line 182 of file tmop_pa_h2d.cpp.
|
protected |
Definition at line 252 of file tmop_pa_h3d.cpp.
|
protected |
Definition at line 84 of file tmop_pa_h2d_c0.cpp.
|
protected |
Definition at line 111 of file tmop_pa_h3d_c0.cpp.
|
virtual |
Assemble the local gradient matrix.
Reimplemented from mfem::NonlinearFormIntegrator.
|
protected |
|
protected |
|
virtual |
Perform the local action of the NonlinearFormIntegrator.
Reimplemented from mfem::NonlinearFormIntegrator.
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
virtual |
Method for computing the diagonal of the gradient with partial assembly.
The result Vector diag is an E-Vector. This method can be called only after the method AssembleGradPA() has been called.
[in,out] | diag | The result Vector: |
Reimplemented from mfem::NonlinearFormIntegrator.
Definition at line 290 of file tmop_pa.cpp.
|
virtual |
Prepare the integrator for partial assembly (PA) gradient evaluations on the given FE space fes at the state x.
The result of the partial assembly is stored internally so that it can be used later in the methods AddMultGradPA() and AssembleGradDiagonalPA(). The state Vector x is an E-vector.
Reimplemented from mfem::NonlinearFormIntegrator.
Definition at line 23 of file tmop_pa.cpp.
|
protected |
Definition at line 349 of file tmop_pa_h2s.cpp.
|
protected |
Definition at line 441 of file tmop_pa_h3s.cpp.
|
protected |
Definition at line 150 of file tmop_pa_h2s_c0.cpp.
|
protected |
Definition at line 169 of file tmop_pa_h3s_c0.cpp.
|
virtual |
Method defining partial assembly.
The result of the partial assembly is stored internally so that it can be used later in the methods AddMultPA().
Reimplemented from mfem::NonlinearFormIntegrator.
Definition at line 215 of file tmop_pa.cpp.
|
virtual |
The result of the partial assembly is stored internally so that it can be used later in the methods AddMultPA(). Used with BilinearFormIntegrators that have different spaces.
Reimplemented from mfem::NonlinearFormIntegrator.
Definition at line 114 of file nonlininteg.cpp.
|
protected |
Definition at line 51 of file tmop_pa.cpp.
|
protected |
Definition at line 167 of file tmop_pa.cpp.
|
protected |
|
protected |
|
protected |
|
protected |
void mfem::TMOP_Integrator::ComputeUntangleMetricQuantiles | ( | const Vector & | x, |
const FiniteElementSpace & | fes ) |
Computes quantiles needed for UntangleMetrics. Note that in parallel, the ParFiniteElementSpace must be passed as argument for consistency across MPI ranks.
|
protected |
|
inlineprotected |
void mfem::TMOP_Integrator::EnableAdaptiveLimiting | ( | const GridFunction & | z0, |
Coefficient & | coeff, | ||
AdaptivityEvaluator & | ae ) |
Restriction of the node positions to certain regions.
Adds the term
[in] | z0 | Function z0 that controls the adaptive limiting. |
[in] | coeff | Coefficient c for the above integral. |
[in] | ae | AdaptivityEvaluator to compute z(x) from z0(x0). |
void mfem::TMOP_Integrator::EnableAdaptiveLimiting | ( | const ParGridFunction & | z0, |
Coefficient & | coeff, | ||
AdaptivityEvaluator & | ae ) |
void mfem::TMOP_Integrator::EnableFiniteDifferences | ( | const GridFunction & | x | ) |
void mfem::TMOP_Integrator::EnableFiniteDifferences | ( | const ParGridFunction & | x | ) |
void mfem::TMOP_Integrator::EnableLimiting | ( | const GridFunction & | n0, |
Coefficient & | w0, | ||
TMOP_LimiterFunction * | lfunc = NULL ) |
void mfem::TMOP_Integrator::EnableLimiting | ( | const GridFunction & | n0, |
const GridFunction & | dist, | ||
Coefficient & | w0, | ||
TMOP_LimiterFunction * | lfunc = NULL ) |
Limiting of the mesh displacements (general version).
Adds the term
[in] | n0 | Original mesh node coordinates (x0 above). |
[in] | dist | Allowed displacement in physical space (d above). |
[in] | w0 | Coefficient scaling the limiting integral. |
[in] | lfunc | TMOP_LimiterFunction defining the function f. If NULL, a TMOP_QuadraticLimiter will be used. The TMOP_Integrator assumes ownership of this pointer. |
void mfem::TMOP_Integrator::EnableNormalization | ( | const GridFunction & | x | ) |
void mfem::TMOP_Integrator::EnableSurfaceFitting | ( | const GridFunction & | pos, |
const Array< bool > & | smarker, | ||
Coefficient & | coeff ) |
Fitting of certain DOFs to given positions in physical space.
Having a set S of marked nodes (or DOFs) and their target positions in physical space x_t, we move these nodes to the target positions during the optimization process. This function adds to the TMOP functional the term
[in] | pos | The desired positions for the mesh nodes. |
[in] | smarker | Indicates which DOFs will be aligned. |
[in] | coeff | Coefficient c for the above integral. |
void mfem::TMOP_Integrator::EnableSurfaceFitting | ( | const GridFunction & | s0, |
const Array< bool > & | smarker, | ||
Coefficient & | coeff, | ||
AdaptivityEvaluator & | ae ) |
Fitting of certain DOFs to the zero level set of a function.
Having a level set function s0(x0) on the starting mesh, and a set of marked nodes (or DOFs), we move these nodes to the zero level set of s0. If s(x) is the image of s0(x0) on the current mesh, this function adds to the TMOP functional the term
[in] | s0 | The level set function on the initial mesh. |
[in] | smarker | Indicates which DOFs will be aligned. |
[in] | coeff | Coefficient c for the above integral. |
[in] | ae | AdaptivityEvaluator to compute s(x) from s0(x0). |
void mfem::TMOP_Integrator::EnableSurfaceFitting | ( | const ParGridFunction & | s0, |
const Array< bool > & | smarker, | ||
Coefficient & | coeff, | ||
AdaptivityEvaluator & | ae ) |
void mfem::TMOP_Integrator::EnableSurfaceFittingFromSource | ( | const ParGridFunction & | s_bg, |
ParGridFunction & | s0, | ||
const Array< bool > & | smarker, | ||
Coefficient & | coeff, | ||
AdaptivityEvaluator & | ae, | ||
const ParGridFunction & | s_bg_grad, | ||
ParGridFunction & | s0_grad, | ||
AdaptivityEvaluator & | age, | ||
const ParGridFunction & | s_bg_hess, | ||
ParGridFunction & | s0_hess, | ||
AdaptivityEvaluator & | ahe ) |
Fitting of certain DOFs in the current mesh to the zero level set of a function defined on another (finer) source mesh.
Having a level set function s_bg(x_bg) on a source/background mesh, a set of marked nodes (or DOFs) in the current mesh, we move the marked nodes to the zero level set of s_bg. This functionality is used for surface fitting and tangential relaxation.
[in] | s_bg | The level set function on the background mesh. |
[in] | s0 | The level set function (automatically) interpolated on the initial mesh. |
[in] | smarker | Marker for aligned DOFs in the current mesh. |
[in] | coeff | Coefficient c for the fitting penalty term. |
[in] | ae | Interpolates s(x) from s_bg(x_bg). |
[in] | s_bg_grad | Gradient of s_bg on the background mesh. |
[in] | s0_grad | Gradient of s0 on the initial mesh. |
[in] | age | Interpolates s_grad(x) from s_bg_grad(x_bg). |
[in] | s_bg_hess | Hessian of s(x) on the background mesh. |
[in] | s0_hess | Hessian of s0 on the initial mesh. |
[in] | ahe | Interpolates s_hess(x) from s_bg_hess(x_bg). See the pmesh-fitting miniapp for details on usage. |
|
inlineprotected |
|
inline |
|
virtual |
|
inline |
|
virtual |
Computes the integral of W(Jacobian(Trt)) over a target zone.
[in] | el | Type of FiniteElement. |
[in] | T | Mesh element transformation. |
[in] | elfun | Physical coordinates of the zone. |
Reimplemented from mfem::NonlinearFormIntegrator.
|
protected |
Compute the local (to the MPI rank) energy with partial assembly.
Here the state x is an E-vector. This method can be called only after the method AssemblePA() has been called.
Reimplemented from mfem::NonlinearFormIntegrator.
Definition at line 362 of file tmop_pa.cpp.
Definition at line 171 of file tmop_pa_w2.cpp.
Definition at line 184 of file tmop_pa_w3.cpp.
Definition at line 131 of file tmop_pa_w2_c0.cpp.
Definition at line 143 of file tmop_pa_w3_c0.cpp.
|
virtual |
Computes the mean of the energies of the given element's children.
In addition to the inputs for GetElementEnergy, this function requires an IntegrationRule to be specified that will give the decomposition of the given element based on the refinement type being considered.
real_t mfem::TMOP_Integrator::GetSurfaceFittingWeight | ( | ) |
|
inlineprotected |
|
inline |
|
inline |
void mfem::TMOP_Integrator::ParEnableNormalization | ( | const ParGridFunction & | x | ) |
void mfem::TMOP_Integrator::ParUpdateAfterMeshTopologyChange | ( | ) |
void mfem::TMOP_Integrator::ReleasePADeviceMemory | ( | bool | copy_to_host = true | ) |
|
inline |
Sets a scaling Coefficient for the quality metric term of the integrator.
With this addition, the integrator becomes
Note that the Coefficient is evaluated in the physical configuration and not in the target configuration which may be undefined.
|
inline |
|
inline |
|
inline |
Prescribe a set of integration rules; relevant for mixed meshes.
This function has priority over SetIntRule(), if both are called.
|
inline |
|
protected |
void mfem::TMOP_Integrator::UpdateAfterMeshTopologyChange | ( | ) |
|
protected |
Definition at line 179 of file tmop_pa.cpp.
void mfem::TMOP_Integrator::UpdateSurfaceFittingWeight | ( | real_t | factor | ) |
|
friend |
|
friend |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
const FiniteElementSpace* mfem::TMOP_Integrator::fes |
const GeometricFactors* mfem::TMOP_Integrator::geom |
|
protected |
|
protected |
|
protected |
const IntegrationRule* mfem::TMOP_Integrator::ir |
|
protected |
|
protected |
|
protected |
|
mutable |
|
protected |
|
protected |
|
protected |
|
protected |
const DofToQuad* mfem::TMOP_Integrator::maps_lim = nullptr |
|
protected |
|
protected |
|
protected |
struct { ... } mfem::TMOP_Integrator::PA |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |