MFEM  v4.4.0
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
nonlininteg.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2022, Lawrence Livermore National Security, LLC. Produced
2 // at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3 // LICENSE and NOTICE for details. LLNL-CODE-806117.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability visit https://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the BSD-3 license. We welcome feedback and contributions, see file
10 // CONTRIBUTING.md for details.
11 
12 #ifndef MFEM_NONLININTEG
13 #define MFEM_NONLININTEG
14 
15 #include "../config/config.hpp"
16 #include "fe.hpp"
17 #include "coefficient.hpp"
18 #include "fespace.hpp"
19 #include "ceed/operator.hpp"
20 
21 namespace mfem
22 {
23 
24 /** @brief This class is used to express the local action of a general nonlinear
25  finite element operator. In addition it may provide the capability to
26  assemble the local gradient operator and to compute the local energy. */
28 {
29 protected:
31 
32  // CEED extension
33  ceed::Operator* ceedOp;
34 
36 
38  : IntRule(ir), ceedOp(NULL) { }
39 
40 public:
41  /** @brief Prescribe a fixed IntegrationRule to use (when @a ir != NULL) or
42  let the integrator choose (when @a ir == NULL). */
43  virtual void SetIntRule(const IntegrationRule *ir) { IntRule = ir; }
44 
45  /// Prescribe a fixed IntegrationRule to use.
46  void SetIntegrationRule(const IntegrationRule &ir) { SetIntRule(&ir); }
47 
48  /// Set the memory type used for GeometricFactors and other large allocations
49  /// in PA extensions.
50  void SetPAMemoryType(MemoryType mt) { pa_mt = mt; }
51 
52  /// Get the integration rule of the integrator (possibly NULL).
53  const IntegrationRule *GetIntegrationRule() const { return IntRule; }
54 
55  /// Perform the local action of the NonlinearFormIntegrator
56  virtual void AssembleElementVector(const FiniteElement &el,
58  const Vector &elfun, Vector &elvect);
59 
60  /// @brief Perform the local action of the NonlinearFormIntegrator resulting
61  /// from a face integral term.
62  virtual void AssembleFaceVector(const FiniteElement &el1,
63  const FiniteElement &el2,
65  const Vector &elfun, Vector &elvect);
66 
67  /// Assemble the local gradient matrix
68  virtual void AssembleElementGrad(const FiniteElement &el,
70  const Vector &elfun, DenseMatrix &elmat);
71 
72  /// @brief Assemble the local action of the gradient of the
73  /// NonlinearFormIntegrator resulting from a face integral term.
74  virtual void AssembleFaceGrad(const FiniteElement &el1,
75  const FiniteElement &el2,
77  const Vector &elfun, DenseMatrix &elmat);
78 
79  /// Compute the local energy
80  virtual double GetElementEnergy(const FiniteElement &el,
82  const Vector &elfun);
83 
84  /// Method defining partial assembly.
85  /** The result of the partial assembly is stored internally so that it can be
86  used later in the methods AddMultPA(). */
87  virtual void AssemblePA(const FiniteElementSpace &fes);
88 
89  /** The result of the partial assembly is stored internally so that it can be
90  used later in the methods AddMultPA().
91  Used with BilinearFormIntegrators that have different spaces. */
92  virtual void AssemblePA(const FiniteElementSpace &trial_fes,
93  const FiniteElementSpace &test_fes);
94 
95  /** @brief Prepare the integrator for partial assembly (PA) gradient
96  evaluations on the given FE space @a fes at the state @a x. */
97  /** The result of the partial assembly is stored internally so that it can be
98  used later in the methods AddMultGradPA() and AssembleGradDiagonalPA().
99  The state Vector @a x is an E-vector. */
100  virtual void AssembleGradPA(const Vector &x, const FiniteElementSpace &fes);
101 
102  /// Compute the local (to the MPI rank) energy with partial assembly.
103  /** Here the state @a x is an E-vector. This method can be called only after
104  the method AssemblePA() has been called. */
105  virtual double GetLocalStateEnergyPA(const Vector &x) const;
106 
107  /// Method for partially assembled action.
108  /** Perform the action of integrator on the input @a x and add the result to
109  the output @a y. Both @a x and @a y are E-vectors, i.e. they represent
110  the element-wise discontinuous version of the FE space.
111 
112  This method can be called only after the method AssemblePA() has been
113  called. */
114  virtual void AddMultPA(const Vector &x, Vector &y) const;
115 
116  /// Method for partially assembled gradient action.
117  /** All arguments are E-vectors. This method can be called only after the
118  method AssembleGradPA() has been called.
119 
120  @param[in] x The gradient Operator is applied to the Vector @a x.
121  @param[in,out] y The result Vector: @f$ y += G x @f$. */
122  virtual void AddMultGradPA(const Vector &x, Vector &y) const;
123 
124  /// Method for computing the diagonal of the gradient with partial assembly.
125  /** The result Vector @a diag is an E-Vector. This method can be called only
126  after the method AssembleGradPA() has been called.
127 
128  @param[in,out] diag The result Vector: @f$ diag += diag(G) @f$. */
129  virtual void AssembleGradDiagonalPA(Vector &diag) const;
130 
131  /// Indicates whether this integrator can use a Ceed backend.
132  virtual bool SupportsCeed() const { return false; }
133 
134  /// Method defining fully unassembled operator.
135  virtual void AssembleMF(const FiniteElementSpace &fes);
136 
137  /** Perform the action of integrator on the input @a x and add the result to
138  the output @a y. Both @a x and @a y are E-vectors, i.e. they represent
139  the element-wise discontinuous version of the FE space.
140 
141  This method can be called only after the method AssembleMF() has been
142  called. */
143  virtual void AddMultMF(const Vector &x, Vector &y) const;
144 
145  ceed::Operator& GetCeedOp() { return *ceedOp; }
146 
148  {
149  delete ceedOp;
150  }
151 };
152 
153 /** The abstract base class BlockNonlinearFormIntegrator is
154  a generalization of the NonlinearFormIntegrator class suitable
155  for block state vectors. */
157 {
158 public:
159  /// Compute the local energy
160  virtual double GetElementEnergy(const Array<const FiniteElement *>&el,
162  const Array<const Vector *>&elfun);
163 
164  /// Perform the local action of the BlockNonlinearFormIntegrator
167  const Array<const Vector *> &elfun,
168  const Array<Vector *> &elvec);
169 
170  virtual void AssembleFaceVector(const Array<const FiniteElement *> &el1,
171  const Array<const FiniteElement *> &el2,
173  const Array<const Vector *> &elfun,
174  const Array<Vector *> &elvect);
175 
176  /// Assemble the local gradient matrix
177  virtual void AssembleElementGrad(const Array<const FiniteElement*> &el,
179  const Array<const Vector *> &elfun,
180  const Array2D<DenseMatrix *> &elmats);
181 
182  virtual void AssembleFaceGrad(const Array<const FiniteElement *>&el1,
185  const Array<const Vector *> &elfun,
186  const Array2D<DenseMatrix *> &elmats);
187 
189 };
190 
191 
192 /// Abstract class for hyperelastic models
194 {
195 protected:
196  ElementTransformation *Ttr; /**< Reference-element to target-element
197  transformation. */
198 
199 public:
200  HyperelasticModel() : Ttr(NULL) { }
201  virtual ~HyperelasticModel() { }
202 
203  /// A reference-element to target-element transformation that can be used to
204  /// evaluate Coefficient%s.
205  /** @note It is assumed that Ttr_.SetIntPoint() is already called for the
206  point of interest. */
207  void SetTransformation(ElementTransformation &Ttr_) { Ttr = &Ttr_; }
208 
209  /** @brief Evaluate the strain energy density function, W = W(Jpt).
210  @param[in] Jpt Represents the target->physical transformation
211  Jacobian matrix. */
212  virtual double EvalW(const DenseMatrix &Jpt) const = 0;
213 
214  /** @brief Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
215  @param[in] Jpt Represents the target->physical transformation
216  Jacobian matrix.
217  @param[out] P The evaluated 1st Piola-Kirchhoff stress tensor. */
218  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const = 0;
219 
220  /** @brief Evaluate the derivative of the 1st Piola-Kirchhoff stress tensor
221  and assemble its contribution to the local gradient matrix 'A'.
222  @param[in] Jpt Represents the target->physical transformation
223  Jacobian matrix.
224  @param[in] DS Gradient of the basis matrix (dof x dim).
225  @param[in] weight Quadrature weight coefficient for the point.
226  @param[in,out] A Local gradient matrix where the contribution from this
227  point will be added.
228 
229  Computes weight * d(dW_dxi)_d(xj) at the current point, for all i and j,
230  where x1 ... xn are the FE dofs. This function is usually defined using
231  the matrix invariants and their derivatives.
232  */
233  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
234  const double weight, DenseMatrix &A) const = 0;
235 };
236 
237 
238 /** Inverse-harmonic hyperelastic model with a strain energy density function
239  given by the formula: W(J) = (1/2) det(J) Tr((J J^t)^{-1}) where J is the
240  deformation gradient. */
242 {
243 protected:
244  mutable DenseMatrix Z, S; // dim x dim
245  mutable DenseMatrix G, C; // dof x dim
246 
247 public:
248  virtual double EvalW(const DenseMatrix &J) const;
249 
250  virtual void EvalP(const DenseMatrix &J, DenseMatrix &P) const;
251 
252  virtual void AssembleH(const DenseMatrix &J, const DenseMatrix &DS,
253  const double weight, DenseMatrix &A) const;
254 };
255 
256 
257 /** Neo-Hookean hyperelastic model with a strain energy density function given
258  by the formula: \f$(\mu/2)(\bar{I}_1 - dim) + (K/2)(det(J)/g - 1)^2\f$ where
259  J is the deformation gradient and \f$\bar{I}_1 = (det(J))^{-2/dim} Tr(J
260  J^t)\f$. The parameters \f$\mu\f$ and K are the shear and bulk moduli,
261  respectively, and g is a reference volumetric scaling. */
263 {
264 protected:
265  mutable double mu, K, g;
268 
269  mutable DenseMatrix Z; // dim x dim
270  mutable DenseMatrix G, C; // dof x dim
271 
272  inline void EvalCoeffs() const;
273 
274 public:
275  NeoHookeanModel(double mu_, double K_, double g_ = 1.0)
276  : mu(mu_), K(K_), g(g_), have_coeffs(false) { c_mu = c_K = c_g = NULL; }
277 
279  : mu(0.0), K(0.0), g(1.0), c_mu(&mu_), c_K(&K_), c_g(g_),
280  have_coeffs(true) { }
281 
282  virtual double EvalW(const DenseMatrix &J) const;
283 
284  virtual void EvalP(const DenseMatrix &J, DenseMatrix &P) const;
285 
286  virtual void AssembleH(const DenseMatrix &J, const DenseMatrix &DS,
287  const double weight, DenseMatrix &A) const;
288 };
289 
290 
291 /** Hyperelastic integrator for any given HyperelasticModel.
292 
293  Represents @f$ \int W(Jpt) dx @f$ over a target zone, where W is the
294  @a model's strain energy density function, and Jpt is the Jacobian of the
295  target->physical coordinates transformation. The target configuration is
296  given by the current mesh at the time of the evaluation of the integrator.
297 */
299 {
300 private:
301  HyperelasticModel *model;
302 
303  // Jrt: the Jacobian of the target-to-reference-element transformation.
304  // Jpr: the Jacobian of the reference-to-physical-element transformation.
305  // Jpt: the Jacobian of the target-to-physical-element transformation.
306  // P: represents dW_d(Jtp) (dim x dim).
307  // DSh: gradients of reference shape functions (dof x dim).
308  // DS: gradients of the shape functions in the target (stress-free)
309  // configuration (dof x dim).
310  // PMatI: coordinates of the deformed configuration (dof x dim).
311  // PMatO: reshaped view into the local element contribution to the operator
312  // output - the result of AssembleElementVector() (dof x dim).
313  DenseMatrix DSh, DS, Jrt, Jpr, Jpt, P, PMatI, PMatO;
314 
315 public:
316  /** @param[in] m HyperelasticModel that will be integrated. */
318 
319  /** @brief Computes the integral of W(Jacobian(Trt)) over a target zone
320  @param[in] el Type of FiniteElement.
321  @param[in] Ttr Represents ref->target coordinates transformation.
322  @param[in] elfun Physical coordinates of the zone. */
323  virtual double GetElementEnergy(const FiniteElement &el,
325  const Vector &elfun);
326 
327  virtual void AssembleElementVector(const FiniteElement &el,
329  const Vector &elfun, Vector &elvect);
330 
331  virtual void AssembleElementGrad(const FiniteElement &el,
333  const Vector &elfun, DenseMatrix &elmat);
334 };
335 
336 /** Hyperelastic incompressible Neo-Hookean integrator with the PK1 stress
337  \f$P = \mu F - p F^{-T}\f$ where \f$\mu\f$ is the shear modulus,
338  \f$p\f$ is the pressure, and \f$F\f$ is the deformation gradient */
340 {
341 private:
342  Coefficient *c_mu;
343  DenseMatrix DSh_u, DS_u, J0i, J, J1, Finv, P, F, FinvT;
344  DenseMatrix PMatI_u, PMatO_u, PMatI_p, PMatO_p, Z, G, C;
345  Vector Sh_p;
346 
347 public:
349 
350  virtual double GetElementEnergy(const Array<const FiniteElement *>&el,
352  const Array<const Vector *> &elfun);
353 
354  /// Perform the local action of the NonlinearFormIntegrator
357  const Array<const Vector *> &elfun,
358  const Array<Vector *> &elvec);
359 
360  /// Assemble the local gradient matrix
361  virtual void AssembleElementGrad(const Array<const FiniteElement*> &el,
363  const Array<const Vector *> &elfun,
364  const Array2D<DenseMatrix *> &elmats);
365 };
366 
367 
369 {
370 private:
371  Coefficient *Q{};
372  DenseMatrix dshape, dshapex, EF, gradEF, ELV, elmat_comp;
373  Vector shape;
374  // PA extension
375  Vector pa_data;
376  const DofToQuad *maps; ///< Not owned
377  const GeometricFactors *geom; ///< Not owned
378  int dim, ne, nq;
379 
380 public:
382 
383  VectorConvectionNLFIntegrator() = default;
384 
385  static const IntegrationRule &GetRule(const FiniteElement &fe,
387 
388  virtual void AssembleElementVector(const FiniteElement &el,
390  const Vector &elfun,
391  Vector &elvect);
392 
393  virtual void AssembleElementGrad(const FiniteElement &el,
395  const Vector &elfun,
396  DenseMatrix &elmat);
397 
399 
400  virtual void AssemblePA(const FiniteElementSpace &fes);
401 
402  virtual void AssembleMF(const FiniteElementSpace &fes);
403 
404  virtual void AddMultPA(const Vector &x, Vector &y) const;
405 
406  virtual void AddMultMF(const Vector &x, Vector &y) const;
407 };
408 
409 
410 /** This class is used to assemble the convective form of the nonlinear term
411  arising in the Navier-Stokes equations \f$(u \cdot \nabla v, w )\f$ */
414 {
415 private:
416  Coefficient *Q{};
417  DenseMatrix dshape, dshapex, EF, gradEF, ELV, elmat_comp;
418  Vector shape;
419 
420 public:
422 
424 
425  virtual void AssembleElementGrad(const FiniteElement &el,
427  const Vector &elfun,
428  DenseMatrix &elmat);
429 };
430 
431 
432 /** This class is used to assemble the skew-symmetric form of the nonlinear term
433  arising in the Navier-Stokes equations
434  \f$.5*(u \cdot \nabla v, w ) - .5*(u \cdot \nabla w, v )\f$ */
437 {
438 private:
439  Coefficient *Q{};
440  DenseMatrix dshape, dshapex, EF, gradEF, ELV, elmat_comp;
441  Vector shape;
442 
443 public:
445 
447 
448  virtual void AssembleElementGrad(const FiniteElement &el,
450  const Vector &elfun,
451  DenseMatrix &elmat);
452 };
453 
454 }
455 
456 #endif
Abstract class for all finite elements.
Definition: fe_base.hpp:235
virtual void AddMultPA(const Vector &x, Vector &y) const
Method for partially assembled action.
Definition: nonlininteg.cpp:45
virtual void EvalP(const DenseMatrix &J, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
virtual double EvalW(const DenseMatrix &J) const
Evaluate the strain energy density function, W = W(Jpt).
void trans(const Vector &u, Vector &x)
Definition: ex27.cpp:412
virtual void AssembleElementGrad(const Array< const FiniteElement * > &el, ElementTransformation &Tr, const Array< const Vector * > &elfun, const Array2D< DenseMatrix * > &elmats)
Assemble the local gradient matrix.
Class for an integration rule - an Array of IntegrationPoint.
Definition: intrules.hpp:90
virtual void AssembleFaceVector(const Array< const FiniteElement * > &el1, const Array< const FiniteElement * > &el2, FaceElementTransformations &Tr, const Array< const Vector * > &elfun, const Array< Vector * > &elvect)
ceed::Operator & GetCeedOp()
virtual void AssembleGradDiagonalPA(Vector &diag) const
Method for computing the diagonal of the gradient with partial assembly.
Definition: nonlininteg.cpp:57
virtual void AssembleElementGrad(const FiniteElement &el, ElementTransformation &trans, const Vector &elfun, DenseMatrix &elmat)
Assemble the local gradient matrix.
virtual void AssembleElementGrad(const FiniteElement &el, ElementTransformation &trans, const Vector &elfun, DenseMatrix &elmat)
Assemble the local gradient matrix.
void SetIntegrationRule(const IntegrationRule &ir)
Prescribe a fixed IntegrationRule to use.
Definition: nonlininteg.hpp:46
virtual void AssembleH(const DenseMatrix &J, const DenseMatrix &DS, const double weight, DenseMatrix &A) const
Evaluate the derivative of the 1st Piola-Kirchhoff stress tensor and assemble its contribution to the...
virtual void AssembleMF(const FiniteElementSpace &fes)
Method defining fully unassembled operator.
A specialized ElementTransformation class representing a face and its two neighboring elements...
Definition: eltrans.hpp:480
Data type dense matrix using column-major storage.
Definition: densemat.hpp:23
virtual double EvalW(const DenseMatrix &Jpt) const =0
Evaluate the strain energy density function, W = W(Jpt).
virtual void SetIntRule(const IntegrationRule *ir)
Prescribe a fixed IntegrationRule to use (when ir != NULL) or let the integrator choose (when ir == N...
Definition: nonlininteg.hpp:43
Structure for storing mesh geometric factors: coordinates, Jacobians, and determinants of the Jacobia...
Definition: mesh.hpp:1814
virtual void AssembleElementVector(const FiniteElement &el, ElementTransformation &trans, const Vector &elfun, Vector &elvect)
Perform the local action of the NonlinearFormIntegrator.
virtual void EvalP(const DenseMatrix &J, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
virtual void AssembleFaceGrad(const Array< const FiniteElement * > &el1, const Array< const FiniteElement * > &el2, FaceElementTransformations &Tr, const Array< const Vector * > &elfun, const Array2D< DenseMatrix * > &elmats)
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const =0
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
virtual void AssembleElementGrad(const FiniteElement &el, ElementTransformation &trans, const Vector &elfun, DenseMatrix &elmat)
Assemble the local gradient matrix.
virtual void AssembleMF(const FiniteElementSpace &fes)
Method defining fully unassembled operator.
Definition: nonlininteg.cpp:63
virtual void AssemblePA(const FiniteElementSpace &fes)
Method defining partial assembly.
virtual void AssembleElementGrad(const Array< const FiniteElement * > &el, ElementTransformation &Tr, const Array< const Vector * > &elfun, const Array2D< DenseMatrix * > &elmats)
Assemble the local gradient matrix.
const IntegrationRule * GetIntegrationRule() const
Get the integration rule of the integrator (possibly NULL).
Definition: nonlininteg.hpp:53
void EvalCoeffs() const
virtual void AssemblePA(const FiniteElementSpace &fes)
Method defining partial assembly.
Definition: nonlininteg.cpp:25
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 integr...
Definition: nonlininteg.cpp:99
IncompressibleNeoHookeanIntegrator(Coefficient &mu_)
virtual double GetElementEnergy(const Array< const FiniteElement * > &el, ElementTransformation &Tr, const Array< const Vector * > &elfun)
Compute the local energy.
virtual void AddMultGradPA(const Vector &x, Vector &y) const
Method for partially assembled gradient action.
Definition: nonlininteg.cpp:51
virtual double GetLocalStateEnergyPA(const Vector &x) const
Compute the local (to the MPI rank) energy with partial assembly.
Definition: nonlininteg.cpp:18
void SetPAMemoryType(MemoryType mt)
Definition: nonlininteg.hpp:50
virtual void AssembleElementGrad(const FiniteElement &el, ElementTransformation &Ttr, const Vector &elfun, DenseMatrix &elmat)
Assemble the local gradient matrix.
virtual void AssembleElementVector(const FiniteElement &el, ElementTransformation &Tr, const Vector &elfun, Vector &elvect)
Perform the local action of the NonlinearFormIntegrator.
Definition: nonlininteg.cpp:75
virtual double GetElementEnergy(const FiniteElement &el, ElementTransformation &Ttr, const Vector &elfun)
Computes the integral of W(Jacobian(Trt)) over a target zone.
Dynamic 2D array using row-major layout.
Definition: array.hpp:351
virtual void AssembleElementGrad(const FiniteElement &el, ElementTransformation &Tr, const Vector &elfun, DenseMatrix &elmat)
Assemble the local gradient matrix.
Definition: nonlininteg.cpp:91
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Definition: fespace.hpp:88
void SetTransformation(ElementTransformation &Ttr_)
Base class Coefficients that optionally depend on space and time. These are used by the BilinearFormI...
Definition: coefficient.hpp:39
MemoryType
Memory types supported by MFEM.
Definition: mem_manager.hpp:31
virtual void AddMultMF(const Vector &x, Vector &y) const
Definition: nonlininteg.cpp:69
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...
Definition: nonlininteg.cpp:83
Structure representing the matrices/tensors needed to evaluate (in reference space) the values...
Definition: fe_base.hpp:136
ElementTransformation * Ttr
const IntegrationRule * IntRule
Definition: nonlininteg.hpp:30
VectorConvectionNLFIntegrator(Coefficient &q)
virtual void AssembleH(const DenseMatrix &J, const DenseMatrix &DS, const double weight, DenseMatrix &A) const
Evaluate the derivative of the 1st Piola-Kirchhoff stress tensor and assemble its contribution to the...
HyperelasticNLFIntegrator(HyperelasticModel *m)
virtual double EvalW(const DenseMatrix &J) const
Evaluate the strain energy density function, W = W(Jpt).
This class is used to express the local action of a general nonlinear finite element operator...
Definition: nonlininteg.hpp:27
virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS, const double weight, DenseMatrix &A) const =0
Evaluate the derivative of the 1st Piola-Kirchhoff stress tensor and assemble its contribution to the...
virtual void AssembleElementVector(const Array< const FiniteElement * > &el, ElementTransformation &Tr, const Array< const Vector * > &elfun, const Array< Vector * > &elvec)
Perform the local action of the NonlinearFormIntegrator.
virtual void AddMultPA(const Vector &x, Vector &y) const
Method for partially assembled action.
NeoHookeanModel(Coefficient &mu_, Coefficient &K_, Coefficient *g_=NULL)
NeoHookeanModel(double mu_, double K_, double g_=1.0)
Abstract class for hyperelastic models.
virtual void AssembleElementVector(const FiniteElement &el, ElementTransformation &Ttr, const Vector &elfun, Vector &elvect)
Perform the local action of the NonlinearFormIntegrator.
virtual double GetElementEnergy(const Array< const FiniteElement * > &el, ElementTransformation &Tr, const Array< const Vector * > &elfun)
Compute the local energy.
Vector data type.
Definition: vector.hpp:60
virtual bool SupportsCeed() const
Indicates whether this integrator can use a Ceed backend.
NonlinearFormIntegrator(const IntegrationRule *ir=NULL)
Definition: nonlininteg.hpp:37
virtual double GetElementEnergy(const FiniteElement &el, ElementTransformation &Tr, const Vector &elfun)
Compute the local energy.
virtual void AddMultMF(const Vector &x, Vector &y) const
static const IntegrationRule & GetRule(const FiniteElement &fe, ElementTransformation &T)
virtual void AssembleElementVector(const Array< const FiniteElement * > &el, ElementTransformation &Tr, const Array< const Vector * > &elfun, const Array< Vector * > &elvec)
Perform the local action of the BlockNonlinearFormIntegrator.
virtual void AssembleGradPA(const Vector &x, const FiniteElementSpace &fes)
Prepare the integrator for partial assembly (PA) gradient evaluations on the given FE space fes at th...
Definition: nonlininteg.cpp:38