MFEM  v4.6.0
Finite element discretization library
nonlininteg.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2023, 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"
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 public:
30  enum Mode
31  {
32  ELEMENTWISE = 0, /**< Element-wise integration (default) */
33  PATCHWISE = 1, /**< Patch-wise integration (NURBS meshes) */
34  PATCHWISE_REDUCED = 2, /**< Patch-wise integration (NURBS meshes) with
35  reduced integration rules. */
36  };
37 
38 protected:
40 
41  Mode integrationMode = Mode::ELEMENTWISE;
42 
43  // Prescribed integration rules (not reduced approximate rules).
45 
46  // CEED extension
48 
50 
52  : IntRule(ir), ceedOp(NULL) { }
53 
54 public:
55  /** @brief Prescribe a fixed IntegrationRule to use (when @a ir != NULL) or
56  let the integrator choose (when @a ir == NULL). */
57  virtual void SetIntRule(const IntegrationRule *ir) { IntRule = ir; }
58 
60 
61  /// For patchwise integration, SetNURBSPatchIntRule must be called.
63  bool HasNURBSPatchIntRule() const { return patchRules != nullptr; }
64 
65  bool Patchwise() const { return integrationMode != Mode::ELEMENTWISE; }
66 
67  /// Prescribe a fixed IntegrationRule to use.
68  void SetIntegrationRule(const IntegrationRule &ir) { SetIntRule(&ir); }
69 
70  /// Set the memory type used for GeometricFactors and other large allocations
71  /// in PA extensions.
72  void SetPAMemoryType(MemoryType mt) { pa_mt = mt; }
73 
74  /// Get the integration rule of the integrator (possibly NULL).
75  const IntegrationRule *GetIntegrationRule() const { return IntRule; }
76 
77  /// Perform the local action of the NonlinearFormIntegrator
78  virtual void AssembleElementVector(const FiniteElement &el,
80  const Vector &elfun, Vector &elvect);
81 
82  /// @brief Perform the local action of the NonlinearFormIntegrator resulting
83  /// from a face integral term.
84  virtual void AssembleFaceVector(const FiniteElement &el1,
85  const FiniteElement &el2,
87  const Vector &elfun, Vector &elvect);
88 
89  /// Assemble the local gradient matrix
90  virtual void AssembleElementGrad(const FiniteElement &el,
92  const Vector &elfun, DenseMatrix &elmat);
93 
94  /// @brief Assemble the local action of the gradient of the
95  /// NonlinearFormIntegrator resulting from a face integral term.
96  virtual void AssembleFaceGrad(const FiniteElement &el1,
97  const FiniteElement &el2,
99  const Vector &elfun, DenseMatrix &elmat);
100 
101  /// Compute the local energy
102  virtual double GetElementEnergy(const FiniteElement &el,
104  const Vector &elfun);
105 
106  /// Method defining partial assembly.
107  /** The result of the partial assembly is stored internally so that it can be
108  used later in the methods AddMultPA(). */
109  virtual void AssemblePA(const FiniteElementSpace &fes);
110 
111  /** The result of the partial assembly is stored internally so that it can be
112  used later in the methods AddMultPA().
113  Used with BilinearFormIntegrators that have different spaces. */
114  virtual void AssemblePA(const FiniteElementSpace &trial_fes,
115  const FiniteElementSpace &test_fes);
116 
117  /** @brief Prepare the integrator for partial assembly (PA) gradient
118  evaluations on the given FE space @a fes at the state @a x. */
119  /** The result of the partial assembly is stored internally so that it can be
120  used later in the methods AddMultGradPA() and AssembleGradDiagonalPA().
121  The state Vector @a x is an E-vector. */
122  virtual void AssembleGradPA(const Vector &x, const FiniteElementSpace &fes);
123 
124  /// Compute the local (to the MPI rank) energy with partial assembly.
125  /** Here the state @a x is an E-vector. This method can be called only after
126  the method AssemblePA() has been called. */
127  virtual double GetLocalStateEnergyPA(const Vector &x) const;
128 
129  /// Method for partially assembled action.
130  /** Perform the action of integrator on the input @a x and add the result to
131  the output @a y. Both @a x and @a y are E-vectors, i.e. they represent
132  the element-wise discontinuous version of the FE space.
133 
134  This method can be called only after the method AssemblePA() has been
135  called. */
136  virtual void AddMultPA(const Vector &x, Vector &y) const;
137 
138  /// Method for partially assembled gradient action.
139  /** All arguments are E-vectors. This method can be called only after the
140  method AssembleGradPA() has been called.
141 
142  @param[in] x The gradient Operator is applied to the Vector @a x.
143  @param[in,out] y The result Vector: @f$ y += G x @f$. */
144  virtual void AddMultGradPA(const Vector &x, Vector &y) const;
145 
146  /// Method for computing the diagonal of the gradient with partial assembly.
147  /** The result Vector @a diag is an E-Vector. This method can be called only
148  after the method AssembleGradPA() has been called.
149 
150  @param[in,out] diag The result Vector: @f$ diag += diag(G) @f$. */
151  virtual void AssembleGradDiagonalPA(Vector &diag) const;
152 
153  /// Indicates whether this integrator can use a Ceed backend.
154  virtual bool SupportsCeed() const { return false; }
155 
156  /// Method defining fully unassembled operator.
157  virtual void AssembleMF(const FiniteElementSpace &fes);
158 
159  /** Perform the action of integrator on the input @a x and add the result to
160  the output @a y. Both @a x and @a y are E-vectors, i.e. they represent
161  the element-wise discontinuous version of the FE space.
162 
163  This method can be called only after the method AssembleMF() has been
164  called. */
165  virtual void AddMultMF(const Vector &x, Vector &y) const;
166 
168 
170  {
171  delete ceedOp;
172  }
173 };
174 
175 /** The abstract base class BlockNonlinearFormIntegrator is
176  a generalization of the NonlinearFormIntegrator class suitable
177  for block state vectors. */
179 {
180 public:
181  /// Compute the local energy
182  virtual double GetElementEnergy(const Array<const FiniteElement *>&el,
184  const Array<const Vector *>&elfun);
185 
186  /// Perform the local action of the BlockNonlinearFormIntegrator
189  const Array<const Vector *> &elfun,
190  const Array<Vector *> &elvec);
191 
192  virtual void AssembleFaceVector(const Array<const FiniteElement *> &el1,
193  const Array<const FiniteElement *> &el2,
195  const Array<const Vector *> &elfun,
196  const Array<Vector *> &elvect);
197 
198  /// Assemble the local gradient matrix
199  virtual void AssembleElementGrad(const Array<const FiniteElement*> &el,
201  const Array<const Vector *> &elfun,
202  const Array2D<DenseMatrix *> &elmats);
203 
204  virtual void AssembleFaceGrad(const Array<const FiniteElement *>&el1,
207  const Array<const Vector *> &elfun,
208  const Array2D<DenseMatrix *> &elmats);
209 
211 };
212 
213 
214 /// Abstract class for hyperelastic models
216 {
217 protected:
218  ElementTransformation *Ttr; /**< Reference-element to target-element
219  transformation. */
220 
221 public:
222  HyperelasticModel() : Ttr(NULL) { }
223  virtual ~HyperelasticModel() { }
224 
225  /// A reference-element to target-element transformation that can be used to
226  /// evaluate Coefficient%s.
227  /** @note It is assumed that Ttr_.SetIntPoint() is already called for the
228  point of interest. */
229  void SetTransformation(ElementTransformation &Ttr_) { Ttr = &Ttr_; }
230 
231  /** @brief Evaluate the strain energy density function, W = W(Jpt).
232  @param[in] Jpt Represents the target->physical transformation
233  Jacobian matrix. */
234  virtual double EvalW(const DenseMatrix &Jpt) const = 0;
235 
236  /** @brief Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
237  @param[in] Jpt Represents the target->physical transformation
238  Jacobian matrix.
239  @param[out] P The evaluated 1st Piola-Kirchhoff stress tensor. */
240  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const = 0;
241 
242  /** @brief Evaluate the derivative of the 1st Piola-Kirchhoff stress tensor
243  and assemble its contribution to the local gradient matrix 'A'.
244  @param[in] Jpt Represents the target->physical transformation
245  Jacobian matrix.
246  @param[in] DS Gradient of the basis matrix (dof x dim).
247  @param[in] weight Quadrature weight coefficient for the point.
248  @param[in,out] A Local gradient matrix where the contribution from this
249  point will be added.
250 
251  Computes weight * d(dW_dxi)_d(xj) at the current point, for all i and j,
252  where x1 ... xn are the FE dofs. This function is usually defined using
253  the matrix invariants and their derivatives.
254  */
255  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
256  const double weight, DenseMatrix &A) const = 0;
257 };
258 
259 
260 /** Inverse-harmonic hyperelastic model with a strain energy density function
261  given by the formula: W(J) = (1/2) det(J) Tr((J J^t)^{-1}) where J is the
262  deformation gradient. */
264 {
265 protected:
266  mutable DenseMatrix Z, S; // dim x dim
267  mutable DenseMatrix G, C; // dof x dim
268 
269 public:
270  virtual double EvalW(const DenseMatrix &J) const;
271 
272  virtual void EvalP(const DenseMatrix &J, DenseMatrix &P) const;
273 
274  virtual void AssembleH(const DenseMatrix &J, const DenseMatrix &DS,
275  const double weight, DenseMatrix &A) const;
276 };
277 
278 
279 /** Neo-Hookean hyperelastic model with a strain energy density function given
280  by the formula: \f$(\mu/2)(\bar{I}_1 - dim) + (K/2)(det(J)/g - 1)^2\f$ where
281  J is the deformation gradient and \f$\bar{I}_1 = (det(J))^{-2/dim} Tr(J
282  J^t)\f$. The parameters \f$\mu\f$ and K are the shear and bulk moduli,
283  respectively, and g is a reference volumetric scaling. */
285 {
286 protected:
287  mutable double mu, K, g;
290 
291  mutable DenseMatrix Z; // dim x dim
292  mutable DenseMatrix G, C; // dof x dim
293 
294  inline void EvalCoeffs() const;
295 
296 public:
297  NeoHookeanModel(double mu_, double K_, double g_ = 1.0)
298  : mu(mu_), K(K_), g(g_), have_coeffs(false) { c_mu = c_K = c_g = NULL; }
299 
301  : mu(0.0), K(0.0), g(1.0), c_mu(&mu_), c_K(&K_), c_g(g_),
302  have_coeffs(true) { }
303 
304  virtual double EvalW(const DenseMatrix &J) const;
305 
306  virtual void EvalP(const DenseMatrix &J, DenseMatrix &P) const;
307 
308  virtual void AssembleH(const DenseMatrix &J, const DenseMatrix &DS,
309  const double weight, DenseMatrix &A) const;
310 };
311 
312 
313 /** Hyperelastic integrator for any given HyperelasticModel.
314 
315  Represents @f$ \int W(Jpt) dx @f$ over a target zone, where W is the
316  @a model's strain energy density function, and Jpt is the Jacobian of the
317  target->physical coordinates transformation. The target configuration is
318  given by the current mesh at the time of the evaluation of the integrator.
319 */
321 {
322 private:
323  HyperelasticModel *model;
324 
325  // Jrt: the Jacobian of the target-to-reference-element transformation.
326  // Jpr: the Jacobian of the reference-to-physical-element transformation.
327  // Jpt: the Jacobian of the target-to-physical-element transformation.
328  // P: represents dW_d(Jtp) (dim x dim).
329  // DSh: gradients of reference shape functions (dof x dim).
330  // DS: gradients of the shape functions in the target (stress-free)
331  // configuration (dof x dim).
332  // PMatI: coordinates of the deformed configuration (dof x dim).
333  // PMatO: reshaped view into the local element contribution to the operator
334  // output - the result of AssembleElementVector() (dof x dim).
335  DenseMatrix DSh, DS, Jrt, Jpr, Jpt, P, PMatI, PMatO;
336 
337 public:
338  /** @param[in] m HyperelasticModel that will be integrated. */
340 
341  /** @brief Computes the integral of W(Jacobian(Trt)) over a target zone
342  @param[in] el Type of FiniteElement.
343  @param[in] Ttr Represents ref->target coordinates transformation.
344  @param[in] elfun Physical coordinates of the zone. */
345  virtual double GetElementEnergy(const FiniteElement &el,
347  const Vector &elfun);
348 
349  virtual void AssembleElementVector(const FiniteElement &el,
351  const Vector &elfun, Vector &elvect);
352 
353  virtual void AssembleElementGrad(const FiniteElement &el,
355  const Vector &elfun, DenseMatrix &elmat);
356 };
357 
358 /** Hyperelastic incompressible Neo-Hookean integrator with the PK1 stress
359  \f$P = \mu F - p F^{-T}\f$ where \f$\mu\f$ is the shear modulus,
360  \f$p\f$ is the pressure, and \f$F\f$ is the deformation gradient */
362 {
363 private:
364  Coefficient *c_mu;
365  DenseMatrix DSh_u, DS_u, J0i, J, J1, Finv, P, F, FinvT;
366  DenseMatrix PMatI_u, PMatO_u, PMatI_p, PMatO_p, Z, G, C;
367  Vector Sh_p;
368 
369 public:
371 
372  virtual double GetElementEnergy(const Array<const FiniteElement *>&el,
374  const Array<const Vector *> &elfun);
375 
376  /// Perform the local action of the NonlinearFormIntegrator
379  const Array<const Vector *> &elfun,
380  const Array<Vector *> &elvec);
381 
382  /// Assemble the local gradient matrix
383  virtual void AssembleElementGrad(const Array<const FiniteElement*> &el,
385  const Array<const Vector *> &elfun,
386  const Array2D<DenseMatrix *> &elmats);
387 };
388 
389 
391 {
392 private:
393  Coefficient *Q{};
394  DenseMatrix dshape, dshapex, EF, gradEF, ELV, elmat_comp;
395  Vector shape;
396  // PA extension
397  Vector pa_data;
398  const DofToQuad *maps; ///< Not owned
399  const GeometricFactors *geom; ///< Not owned
400  int dim, ne, nq;
401 
402 public:
404 
405  VectorConvectionNLFIntegrator() = default;
406 
407  static const IntegrationRule &GetRule(const FiniteElement &fe,
409 
410  virtual void AssembleElementVector(const FiniteElement &el,
412  const Vector &elfun,
413  Vector &elvect);
414 
415  virtual void AssembleElementGrad(const FiniteElement &el,
417  const Vector &elfun,
418  DenseMatrix &elmat);
419 
421 
422  virtual void AssemblePA(const FiniteElementSpace &fes);
423 
424  virtual void AssembleMF(const FiniteElementSpace &fes);
425 
426  virtual void AddMultPA(const Vector &x, Vector &y) const;
427 
428  virtual void AddMultMF(const Vector &x, Vector &y) const;
429 };
430 
431 
432 /** This class is used to assemble the convective form of the nonlinear term
433  arising in the Navier-Stokes equations \f$(u \cdot \nabla v, w )\f$ */
436 {
437 private:
438  Coefficient *Q{};
439  DenseMatrix dshape, dshapex, EF, gradEF, ELV, elmat_comp;
440  Vector shape;
441 
442 public:
444 
446 
447  virtual void AssembleElementGrad(const FiniteElement &el,
449  const Vector &elfun,
450  DenseMatrix &elmat);
451 };
452 
453 
454 /** This class is used to assemble the skew-symmetric form of the nonlinear term
455  arising in the Navier-Stokes equations
456  \f$.5*(u \cdot \nabla v, w ) - .5*(u \cdot \nabla w, v )\f$ */
459 {
460 private:
461  Coefficient *Q{};
462  DenseMatrix dshape, dshapex, EF, gradEF, ELV, elmat_comp;
463  Vector shape;
464 
465 public:
467 
469 
470  virtual void AssembleElementGrad(const FiniteElement &el,
472  const Vector &elfun,
473  DenseMatrix &elmat);
474 };
475 
476 }
477 
478 #endif
Abstract class for all finite elements.
Definition: fe_base.hpp:233
void trans(const Vector &u, Vector &x)
Definition: ex27.cpp:412
virtual void AddMultPA(const Vector &x, Vector &y) const
Method for partially assembled action.
Definition: nonlininteg.cpp:45
Class for an integration rule - an Array of IntegrationPoint.
Definition: intrules.hpp:96
virtual void AddMultMF(const Vector &x, Vector &y) const
ceed::Operator & GetCeedOp()
virtual double GetElementEnergy(const Array< const FiniteElement *> &el, ElementTransformation &Tr, const Array< const Vector *> &elfun)
Compute the local energy.
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:68
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:57
virtual double EvalW(const DenseMatrix &J) const
Evaluate the strain energy density function, W = W(Jpt).
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...
Structure for storing mesh geometric factors: coordinates, Jacobians, and determinants of the Jacobia...
Definition: mesh.hpp:2183
virtual void AssembleElementVector(const FiniteElement &el, ElementTransformation &trans, const Vector &elfun, Vector &elvect)
Perform the local action of the NonlinearFormIntegrator.
virtual void AddMultPA(const Vector &x, Vector &y) const
Method for partially assembled action.
virtual double GetElementEnergy(const Array< const FiniteElement *> &el, ElementTransformation &Tr, const Array< const Vector *> &elfun)
Compute the local energy.
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const =0
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
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 AssembleElementGrad(const FiniteElement &el, ElementTransformation &trans, const Vector &elfun, DenseMatrix &elmat)
Assemble the local gradient matrix.
virtual void AssemblePA(const FiniteElementSpace &fes)
Method defining partial assembly.
virtual void AssembleMF(const FiniteElementSpace &fes)
Method defining fully unassembled operator.
Definition: nonlininteg.cpp:63
virtual double EvalW(const DenseMatrix &J) const
Evaluate the strain energy density function, W = W(Jpt).
virtual bool SupportsCeed() const
Indicates whether this integrator can use a Ceed backend.
virtual void AssemblePA(const FiniteElementSpace &fes)
Method defining partial assembly.
Definition: nonlininteg.cpp:25
virtual void AddMultGradPA(const Vector &x, Vector &y) const
Method for partially assembled gradient action.
Definition: nonlininteg.cpp:51
virtual void AssembleGradDiagonalPA(Vector &diag) const
Method for computing the diagonal of the gradient with partial assembly.
Definition: nonlininteg.cpp:57
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_)
void SetPAMemoryType(MemoryType mt)
Definition: nonlininteg.hpp:72
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.
virtual void AddMultMF(const Vector &x, Vector &y) const
Definition: nonlininteg.cpp:69
Dynamic 2D array using row-major layout.
Definition: array.hpp:354
virtual void AssembleElementGrad(const FiniteElement &el, ElementTransformation &Tr, const Vector &elfun, DenseMatrix &elmat)
Assemble the local gradient matrix.
Definition: nonlininteg.cpp:91
virtual void AssembleFaceVector(const Array< const FiniteElement *> &el1, const Array< const FiniteElement *> &el2, FaceElementTransformations &Tr, const Array< const Vector *> &elfun, const Array< Vector *> &elvect)
virtual void AssembleFaceGrad(const Array< const FiniteElement *> &el1, const Array< const FiniteElement *> &el2, FaceElementTransformations &Tr, const Array< const Vector *> &elfun, const Array2D< DenseMatrix *> &elmats)
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Definition: fespace.hpp:219
void SetTransformation(ElementTransformation &Ttr_)
Base class Coefficients that optionally depend on space and time. These are used by the BilinearFormI...
Definition: coefficient.hpp:41
MemoryType
Memory types supported by MFEM.
Definition: mem_manager.hpp:31
const IntegrationRule * GetIntegrationRule() const
Get the integration rule of the integrator (possibly NULL).
Definition: nonlininteg.hpp:75
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:39
VectorConvectionNLFIntegrator(Coefficient &q)
virtual void AssembleElementGrad(const Array< const FiniteElement *> &el, ElementTransformation &Tr, const Array< const Vector *> &elfun, const Array2D< DenseMatrix *> &elmats)
Assemble the local gradient matrix.
HyperelasticNLFIntegrator(HyperelasticModel *m)
virtual void AssembleElementGrad(const Array< const FiniteElement *> &el, ElementTransformation &Tr, const Array< const Vector *> &elfun, const Array2D< DenseMatrix *> &elmats)
Assemble the local gradient matrix.
This class is used to express the local action of a general nonlinear finite element operator...
Definition: nonlininteg.hpp:27
void SetNURBSPatchIntRule(NURBSMeshRules *pr)
For patchwise integration, SetNURBSPatchIntRule must be called.
Definition: nonlininteg.hpp:62
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...
void EvalCoeffs() const
virtual double GetLocalStateEnergyPA(const Vector &x) const
Compute the local (to the MPI rank) energy with partial assembly.
Definition: nonlininteg.cpp:18
NeoHookeanModel(Coefficient &mu_, Coefficient &K_, Coefficient *g_=NULL)
NeoHookeanModel(double mu_, double K_, double g_=1.0)
virtual void EvalP(const DenseMatrix &J, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Abstract class for hyperelastic models.
virtual void EvalP(const DenseMatrix &J, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
virtual void AssembleElementVector(const FiniteElement &el, ElementTransformation &Ttr, const Vector &elfun, Vector &elvect)
Perform the local action of the NonlinearFormIntegrator.
Vector data type.
Definition: vector.hpp:58
NonlinearFormIntegrator(const IntegrationRule *ir=NULL)
Definition: nonlininteg.hpp:51
virtual double GetElementEnergy(const FiniteElement &el, ElementTransformation &Tr, const Vector &elfun)
Compute the local energy.
virtual void AssembleMF(const FiniteElementSpace &fes)
Method defining fully unassembled operator.
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...
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 NonlinearFormIntegrator.
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
Class for defining different integration rules on each NURBS patch.
Definition: intrules.hpp:275