MFEM  v4.0
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
tmop.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010, Lawrence Livermore National Security, LLC. Produced at
2 // the Lawrence Livermore National Laboratory. LLNL-CODE-443211. All Rights
3 // reserved. See file COPYRIGHT for details.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability see http://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the GNU Lesser General Public License (as published by the Free
10 // Software Foundation) version 2.1 dated February 1999.
11 
12 #ifndef MFEM_TMOP_HPP
13 #define MFEM_TMOP_HPP
14 
15 #include "../config/config.hpp"
16 #include "../linalg/invariants.hpp"
17 #include "nonlininteg.hpp"
18 
19 namespace mfem
20 {
21 
22 /** @brief Abstract class for local mesh quality metrics in the target-matrix
23  optimization paradigm (TMOP) by P. Knupp et al. */
25 {
26 protected:
27  const DenseMatrix *Jtr; /**< Jacobian of the reference-element to
28  target-element transformation. */
29 
30  /** @brief The method SetTransformation() is hidden for TMOP_QualityMetric%s,
31  because it is not used. */
33 
34 public:
35  TMOP_QualityMetric() : Jtr(NULL) { }
36  virtual ~TMOP_QualityMetric() { }
37 
38  /** @brief Specify the reference-element -> target-element Jacobian matrix
39  for the point of interest.
40 
41  The specified Jacobian matrix, #Jtr, can be used by metrics that cannot
42  be written just as a function of the target->physical Jacobian matrix,
43  Jpt. */
44  void SetTargetJacobian(const DenseMatrix &_Jtr) { Jtr = &_Jtr; }
45 
46  /** @brief Evaluate the strain energy density function, W = W(Jpt).
47  @param[in] Jpt Represents the target->physical transformation
48  Jacobian matrix. */
49  virtual double EvalW(const DenseMatrix &Jpt) const = 0;
50 
51  /** @brief Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
52  @param[in] Jpt Represents the target->physical transformation
53  Jacobian matrix.
54  @param[out] P The evaluated 1st Piola-Kirchhoff stress tensor. */
55  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const = 0;
56 
57  /** @brief Evaluate the derivative of the 1st Piola-Kirchhoff stress tensor
58  and assemble its contribution to the local gradient matrix 'A'.
59  @param[in] Jpt Represents the target->physical transformation
60  Jacobian matrix.
61  @param[in] DS Gradient of the basis matrix (dof x dim).
62  @param[in] weight Quadrature weight coefficient for the point.
63  @param[in,out] A Local gradient matrix where the contribution from this
64  point will be added.
65 
66  Computes weight * d(dW_dxi)_d(xj) at the current point, for all i and j,
67  where x1 ... xn are the FE dofs. This function is usually defined using
68  the matrix invariants and their derivatives.
69  */
70  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
71  const double weight, DenseMatrix &A) const = 0;
72 };
73 
74 
75 /// Metric without a type, 2D
77 {
78 protected:
80 
81 public:
82  // W = |J|^2.
83  virtual double EvalW(const DenseMatrix &Jpt) const;
84 
85  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
86 
87  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
88  const double weight, DenseMatrix &A) const;
89 };
90 
91 /// Skew metric, 2D.
93 {
94 public:
95  // W = 0.5 (1 - cos(angle_Jpr - angle_Jtr)).
96  virtual double EvalW(const DenseMatrix &Jpt) const;
97 
98  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
99  { MFEM_ABORT("Not implemented"); }
100 
101  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
102  const double weight, DenseMatrix &A) const
103  { MFEM_ABORT("Not implemented"); }
104 };
105 
106 /// Skew metric, 3D.
108 {
109 public:
110  // W = 1/6 (3 - sum_i cos(angle_Jpr_i - angle_Jtr_i)), i = 1..3.
111  virtual double EvalW(const DenseMatrix &Jpt) const;
112 
113  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
114  { MFEM_ABORT("Not implemented"); }
115 
116  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
117  const double weight, DenseMatrix &A) const
118  { MFEM_ABORT("Not implemented"); }
119 };
120 
121 /// Aspect ratio metric, 2D.
123 {
124 public:
125  // W = 0.5 (ar_Jpr/ar_Jtr + ar_Jtr/ar_Jpr) - 1.
126  virtual double EvalW(const DenseMatrix &Jpt) const;
127 
128  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
129  { MFEM_ABORT("Not implemented"); }
130 
131  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
132  const double weight, DenseMatrix &A) const
133  { MFEM_ABORT("Not implemented"); }
134 };
135 
136 /// Aspect ratio metric, 3D.
138 {
139 public:
140  // W = 1/3 sum [0.5 (ar_Jpr_i/ar_Jtr_i + ar_Jtr_i/ar_Jpr_i) - 1], i = 1..3.
141  virtual double EvalW(const DenseMatrix &Jpt) const;
142 
143  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
144  { MFEM_ABORT("Not implemented"); }
145 
146  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
147  const double weight, DenseMatrix &A) const
148  { MFEM_ABORT("Not implemented"); }
149 };
150 
151 /// Shape, ideal barrier metric, 2D
153 {
154 protected:
156 
157 public:
158  // W = 0.5|J|^2 / det(J) - 1.
159  virtual double EvalW(const DenseMatrix &Jpt) const;
160 
161  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
162 
163  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
164  const double weight, DenseMatrix &A) const;
165 };
166 
167 /// Shape & area, ideal barrier metric, 2D
169 {
170 protected:
172 
173 public:
174  // W = |J - J^-t|^2.
175  virtual double EvalW(const DenseMatrix &Jpt) const;
176 
177  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
178 
179  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
180  const double weight, DenseMatrix &A) const;
181 };
182 
183 /// Shape & area metric, 2D
185 {
186 protected:
188 
189 public:
190  // W = det(J) * |J - J^-t|^2.
191  virtual double EvalW(const DenseMatrix &Jpt) const;
192 
193  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
194 
195  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
196  const double weight, DenseMatrix &A) const;
197 };
198 
199 /// Shifted barrier form of metric 2 (shape, ideal barrier metric), 2D
201 {
202 protected:
203  double &tau0;
205 
206 public:
207  TMOP_Metric_022(double &t0): tau0(t0) {}
208 
209  // W = 0.5(|J|^2 - 2det(J)) / (det(J) - tau0).
210  virtual double EvalW(const DenseMatrix &Jpt) const;
211 
212  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
213 
214  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
215  const double weight, DenseMatrix &A) const;
216 };
217 
218 /// Shape, ideal barrier metric, 2D
220 {
221 protected:
223 
224 public:
225  // W = 0.5|J^t J|^2 / det(J)^2 - 1.
226  virtual double EvalW(const DenseMatrix &Jpt) const;
227 
228  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
229 
230  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
231  const double weight, DenseMatrix &A) const;
232 };
233 
234 /// Area metric, 2D
236 {
237 protected:
239 
240 public:
241  // W = (det(J) - 1)^2.
242  virtual double EvalW(const DenseMatrix &Jpt) const;
243 
244  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
245 
246  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
247  const double weight, DenseMatrix &A) const;
248 
249 };
250 
251 /// Area, ideal barrier metric, 2D
253 {
254 protected:
256 
257 public:
258  // W = 0.5( sqrt(det(J)) - 1 / sqrt(det(J)) )^2
259  // = 0.5( det(J) - 1 )^2 / det(J)
260  // = 0.5( det(J) + 1/det(J) ) - 1.
261  virtual double EvalW(const DenseMatrix &Jpt) const;
262 
263  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
264 
265  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
266  const double weight, DenseMatrix &A) const;
267 
268 };
269 
270 /// Shape, ideal barrier metric, 2D
272 {
273 protected:
275 
276 public:
277  // W = |J^t J|^2 / det(J)^2 - 2|J|^2 / det(J) + 2
278  // = I1b (I1b - 2).
279  virtual double EvalW(const DenseMatrix &Jpt) const;
280 
281  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
282 
283  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
284  const double weight, DenseMatrix &A) const;
285 
286 };
287 
288 /// Area, ideal barrier metric, 2D
290 {
291 protected:
293 
294 public:
295  // W = 0.5(det(J) - 1 / det(J))^2.
296  virtual double EvalW(const DenseMatrix &Jpt) const;
297 
298  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
299 
300  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
301  const double weight, DenseMatrix &A) const;
302 
303 };
304 
305 /// Untangling metric, 2D
307 {
308 protected:
309  const double eps;
311 
312 public:
313  TMOP_Metric_211(double epsilon = 1e-4) : eps(epsilon) { }
314 
315  // W = (det(J) - 1)^2 - det(J) + sqrt(det(J)^2 + eps).
316  virtual double EvalW(const DenseMatrix &Jpt) const;
317 
318  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
319 
320  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
321  const double weight, DenseMatrix &A) const;
322 };
323 
324 /// Shifted barrier form of metric 56 (area, ideal barrier metric), 2D
326 {
327 protected:
328  double &tau0;
330 
331 public:
332  /// Note that @a t0 is stored by reference
333  TMOP_Metric_252(double &t0): tau0(t0) {}
334 
335  // W = 0.5(det(J) - 1)^2 / (det(J) - tau0).
336  virtual double EvalW(const DenseMatrix &Jpt) const;
337 
338  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
339 
340  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
341  const double weight, DenseMatrix &A) const;
342 };
343 
344 /// Shape, ideal barrier metric, 3D
346 {
347 protected:
349 
350 public:
351  // W = |J| |J^-1| / 3 - 1.
352  virtual double EvalW(const DenseMatrix &Jpt) const;
353 
354  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
355 
356  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
357  const double weight, DenseMatrix &A) const;
358 };
359 
360 /// Shape, ideal barrier metric, 3D
362 {
363 protected:
365 
366 public:
367  // W = |J|^2 |J^-1|^2 / 9 - 1.
368  virtual double EvalW(const DenseMatrix &Jpt) const;
369 
370  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
371 
372  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
373  const double weight, DenseMatrix &A) const;
374 };
375 
376 /// Shape, ideal barrier metric, 3D
378 {
379 protected:
381 
382 public:
383  // W = |J|^2 / 3 * det(J)^(2/3) - 1.
384  virtual double EvalW(const DenseMatrix &Jpt) const;
385 
386  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
387 
388  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
389  const double weight, DenseMatrix &A) const;
390 };
391 
392 /// Volume metric, 3D
394 {
395 protected:
397 
398 public:
399  // W = (det(J) - 1)^2.
400  virtual double EvalW(const DenseMatrix &Jpt) const;
401 
402  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
403 
404  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
405  const double weight, DenseMatrix &A) const;
406 };
407 
408 /// Volume, ideal barrier metric, 3D
410 {
411 protected:
413 
414 public:
415  // W = 0.5( sqrt(det(J)) - 1 / sqrt(det(J)) )^2
416  // = 0.5( det(J) - 1 )^2 / det(J)
417  // = 0.5( det(J) + 1/det(J) ) - 1.
418  virtual double EvalW(const DenseMatrix &Jpt) const;
419 
420  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
421 
422  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
423  const double weight, DenseMatrix &A) const;
424 };
425 
426 /// Shape & volume, ideal barrier metric, 3D
428 {
429 protected:
431 
432 public:
433  // W = |J - J^-t|^2.
434  virtual double EvalW(const DenseMatrix &Jpt) const;
435 
436  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
437 
438  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
439  const double weight, DenseMatrix &A) const;
440 };
441 
442 /// Shifted barrier form of 3D metric 16 (volume, ideal barrier metric), 3D
444 {
445 protected:
446  double &tau0;
448 
449 public:
450  TMOP_Metric_352(double &t0): tau0(t0) {}
451 
452  // W = 0.5(det(J) - 1)^2 / (det(J) - tau0).
453  virtual double EvalW(const DenseMatrix &Jpt) const;
454 
455  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
456 
457  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
458  const double weight, DenseMatrix &A) const;
459 };
460 
461 
462 /// Base class for limiting functions to be used in class TMOP_Integrator.
463 /** This class represents a scalar function f(x, x0, d), where x and x0 are
464  positions in physical space, and d is a reference physical distance
465  associated with the point x0. */
467 {
468 public:
469  /// Returns the limiting function, f(x, x0, d).
470  virtual double Eval(const Vector &x, const Vector &x0, double d) const = 0;
471 
472  /** @brief Returns the gradient of the limiting function f(x, x0, d) with
473  respect to x. */
474  virtual void Eval_d1(const Vector &x, const Vector &x0, double dist,
475  Vector &d1) const = 0;
476 
477  /** @brief Returns the Hessian of the limiting function f(x, x0, d) with
478  respect to x. */
479  virtual void Eval_d2(const Vector &x, const Vector &x0, double dist,
480  DenseMatrix &d2) const = 0;
481 
482  /// Virtual destructor.
483  virtual ~TMOP_LimiterFunction() { }
484 };
485 
486 /// Default limiter function in TMOP_Integrator.
488 {
489 public:
490  virtual double Eval(const Vector &x, const Vector &x0, double dist) const
491  {
492  MFEM_ASSERT(x.Size() == x0.Size(), "Bad input.");
493 
494  return 0.5 * x.DistanceSquaredTo(x0) / (dist * dist);
495  }
496 
497  virtual void Eval_d1(const Vector &x, const Vector &x0, double dist,
498  Vector &d1) const
499  {
500  MFEM_ASSERT(x.Size() == x0.Size(), "Bad input.");
501 
502  d1.SetSize(x.Size());
503  subtract(1.0 / (dist * dist), x, x0, d1);
504  }
505 
506  virtual void Eval_d2(const Vector &x, const Vector &x0, double dist,
507  DenseMatrix &d2) const
508  {
509  MFEM_ASSERT(x.Size() == x0.Size(), "Bad input.");
510 
511  d2.Diag(1.0 / (dist * dist), x.Size());
512  }
513 
515 };
516 
517 
518 /** @brief Base class representing target-matrix construction algorithms for
519  mesh optimization via the target-matrix optimization paradigm (TMOP). */
520 /** This class is used by class TMOP_Integrator to construct the target Jacobian
521  matrices (reference-element to target-element) at quadrature points. It
522  supports a set of algorithms chosen by the #TargetType enumeration.
523 
524  New target-matrix construction algorithms can be defined by deriving new
525  classes and overriding the method ComputeElementTargets(). */
527 {
528 public:
529  /// Target-matrix construction algorithms supported by this class.
531  {
533  Ideal shape, unit size; the nodes are not used. */
535  Ideal shape, equal size/volume; the given nodes define the total target
536  volume; for each mesh element, the target volume is the average volume
537  multiplied by the volume scale, set with SetVolumeScale(). */
539  Ideal shape, given size/volume; the given nodes define the target
540  volume at all quadrature points. */
542  Given shape, given size/volume; the given nodes define the exact target
543  Jacobian matrix at all quadrature points. */
544  };
545 
546 protected:
547  // Nodes that are used in ComputeElementTargets(), depending on target_type.
548  const GridFunction *nodes; // not owned
549  mutable double avg_volume;
550  double volume_scale;
552 
553 #ifdef MFEM_USE_MPI
554  MPI_Comm comm;
555  bool Parallel() const { return (comm != MPI_COMM_NULL); }
556 #else
557  bool Parallel() const { return false; }
558 #endif
559 
560  // should be called only if avg_volume == 0.0, i.e. avg_volume is not
561  // computed yet
562  void ComputeAvgVolume() const;
563 
564 public:
565  /// Constructor for use in serial
567  : nodes(NULL), avg_volume(), volume_scale(1.0), target_type(ttype)
568  {
569 #ifdef MFEM_USE_MPI
570  comm = MPI_COMM_NULL;
571 #endif
572  }
573 #ifdef MFEM_USE_MPI
574  /// Constructor for use in parallel
575  TargetConstructor(TargetType ttype, MPI_Comm mpicomm)
576  : nodes(NULL), avg_volume(), volume_scale(1.0), target_type(ttype),
577  comm(mpicomm) { }
578 #endif
579  virtual ~TargetConstructor() { }
580 
581  /** @brief Set the nodes to be used in the target-matrix construction.
582 
583  This method should be called every time the target nodes are updated
584  externally and recomputation of the target average volume is needed. The
585  nodes are used by all target types except IDEAL_SHAPE_UNIT_SIZE. */
586  void SetNodes(const GridFunction &n) { nodes = &n; avg_volume = 0.0; }
587 
588  /// Used by target type IDEAL_SHAPE_EQUAL_SIZE. The default volume scale is 1.
589  void SetVolumeScale(double vol_scale) { volume_scale = vol_scale; }
590 
591  /** @brief Given an element and quadrature rule, computes ref->target
592  transformation Jacobians for each quadrature point in the element. */
593  virtual void ComputeElementTargets(int e_id, const FiniteElement &fe,
594  const IntegrationRule &ir,
595  DenseTensor &Jtr) const;
596 };
597 
598 class ParGridFunction;
599 
600 /** @brief A TMOP integrator class based on any given TMOP_QualityMetric and
601  TargetConstructor.
602 
603  Represents @f$ \int W(Jpt) dx @f$ over a target zone, where W is the
604  metric's strain energy density function, and Jpt is the Jacobian of the
605  target->physical coordinates transformation. The virtual target zone is
606  defined by the TargetConstructor. */
608 {
609 protected:
610  TMOP_QualityMetric *metric; // not owned
611  const TargetConstructor *targetC; // not owned
612 
613  // Weight Coefficient multiplying the quality metric term.
614  Coefficient *coeff1; // not owned, if NULL -> coeff1 is 1.
615  // Normalization factor for the metric term.
617 
618  // Nodes and weight Coefficient used for "limiting" the TMOP_Integrator.
619  // These are both NULL when there is no limiting.
620  // The class doesn't own nodes0 and coeff0.
623  // Limiting reference distance. Not owned.
625  // Limiting function. Owned.
627  // Normalization factor for the limiting term.
628  double lim_normal;
629 
630  // Jrt: the inverse of the ref->target Jacobian, Jrt = Jtr^{-1}.
631  // Jpr: the ref->physical transformation Jacobian, Jpr = PMatI^t DS.
632  // Jpt: the target->physical transformation Jacobian, Jpt = Jpr Jrt.
633  // P: represents dW_d(Jtp) (dim x dim).
634  // DSh: gradients of reference shape functions (dof x dim).
635  // DS: gradients of the shape functions in the target configuration,
636  // DS = DSh Jrt (dof x dim).
637  // PMatI: current coordinates of the nodes (dof x dim).
638  // PMat0: reshaped view into the local element contribution to the operator
639  // output - the result of AssembleElementVector() (dof x dim).
641 
643  double &metric_energy, double &lim_energy);
644 
645 public:
646  /** @param[in] m TMOP_QualityMetric that will be integrated (not owned).
647  @param[in] tc Target-matrix construction algorithm to use (not owned). */
649  : metric(m), targetC(tc),
650  coeff1(NULL), metric_normal(1.0),
651  nodes0(NULL), coeff0(NULL),
652  lim_dist(NULL), lim_func(NULL), lim_normal(1.0)
653  { }
654 
656 
657  /// Sets a scaling Coefficient for the quality metric term of the integrator.
658  /** With this addition, the integrator becomes
659  @f$ \int w1 W(Jpt) dx @f$.
660 
661  Note that the Coefficient is evaluated in the physical configuration and
662  not in the target configuration which may be undefined. */
663  void SetCoefficient(Coefficient &w1) { coeff1 = &w1; }
664 
665  /// Adds a limiting term to the integrator (general version).
666  /** With this addition, the integrator becomes
667  @f$ \int w1 W(Jpt) + w0 f(x, x_0, d) dx @f$,
668  where the second term measures the change with respect to the original
669  physical positions, @a n0.
670  @param[in] n0 Original mesh node coordinates.
671  @param[in] dist Limiting physical distances.
672  @param[in] w0 Coefficient scaling the limiting term.
673  @param[in] lfunc TMOP_LimiterFunction defining the limiting term f. If
674  NULL, a TMOP_QuadraticLimiter will be used. The
675  TMOP_Integrator assumes ownership of this pointer. */
676  void EnableLimiting(const GridFunction &n0, const GridFunction &dist,
677  Coefficient &w0, TMOP_LimiterFunction *lfunc = NULL);
678 
679  /** @brief Adds a limiting term to the integrator with limiting distance
680  function (@a dist in the general version of the method) equal to 1. */
681  void EnableLimiting(const GridFunction &n0,
682  Coefficient &w0, TMOP_LimiterFunction *lfunc = NULL);
683 
684  /// Update the original/reference nodes used for limiting.
685  void SetLimitingNodes(const GridFunction &n0) { nodes0 = &n0; }
686 
687  /** @brief Computes the integral of W(Jacobian(Trt)) over a target zone.
688  @param[in] el Type of FiniteElement.
689  @param[in] T Mesh element transformation.
690  @param[in] elfun Physical coordinates of the zone. */
691  virtual double GetElementEnergy(const FiniteElement &el,
693  const Vector &elfun);
694 
695  virtual void AssembleElementVector(const FiniteElement &el,
697  const Vector &elfun, Vector &elvect);
698 
699  virtual void AssembleElementGrad(const FiniteElement &el,
701  const Vector &elfun, DenseMatrix &elmat);
702 
703  /** @brief Computes the normalization factors of the metric and limiting
704  integrals using the mesh position given by @a x. */
705  void EnableNormalization(const GridFunction &x);
706 #ifdef MFEM_USE_MPI
708 #endif
709 };
710 
711 
712 /// Interpolates the @a metric's values at the nodes of @a metric_gf.
713 /** Assumes that @a metric_gf's FiniteElementSpace is initialized. */
714 void InterpolateTMOP_QualityMetric(TMOP_QualityMetric &metric,
715  const TargetConstructor &tc,
716  const Mesh &mesh, GridFunction &metric_gf);
717 
718 }
719 
720 #endif
Abstract class for Finite Elements.
Definition: fe.hpp:229
virtual void AssembleH(const DenseMatrix &Jpt, 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...
Definition: tmop.hpp:146
InvariantsEvaluator2D< double > ie
Definition: tmop.hpp:292
virtual void AssembleH(const DenseMatrix &Jpt, 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...
Definition: tmop.cpp:276
virtual void AssembleH(const DenseMatrix &Jpt, 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...
Definition: tmop.cpp:33
Shifted barrier form of metric 56 (area, ideal barrier metric), 2D.
Definition: tmop.hpp:325
Shifted barrier form of 3D metric 16 (volume, ideal barrier metric), 3D.
Definition: tmop.hpp:443
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:257
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 AssembleH(const DenseMatrix &Jpt, 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...
Definition: tmop.hpp:116
Class for an integration rule - an Array of IntegrationPoint.
Definition: intrules.hpp:85
const TargetType target_type
Definition: tmop.hpp:551
Class for grid function - Vector with associated FE space.
Definition: gridfunc.hpp:27
DenseMatrix PMatO
Definition: tmop.hpp:640
Shape &amp; volume, ideal barrier metric, 3D.
Definition: tmop.hpp:427
const DenseMatrix * Jtr
Definition: tmop.hpp:27
InvariantsEvaluator3D< double > ie
Definition: tmop.hpp:430
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:571
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.hpp:143
DenseMatrix Jrt
Definition: tmop.hpp:640
const GridFunction * lim_dist
Definition: tmop.hpp:624
Coefficient * coeff0
Definition: tmop.hpp:622
virtual void Eval_d2(const Vector &x, const Vector &x0, double dist, DenseMatrix &d2) const
Returns the Hessian of the limiting function f(x, x0, d) with respect to x.
Definition: tmop.hpp:506
InvariantsEvaluator2D< double > ie
Definition: tmop.hpp:238
void EnableNormalization(const GridFunction &x)
Computes the normalization factors of the metric and limiting integrals using the mesh position given...
Definition: tmop.cpp:1190
void SetSize(int s)
Resize the vector to size s.
Definition: vector.hpp:400
void InterpolateTMOP_QualityMetric(TMOP_QualityMetric &metric, const TargetConstructor &tc, const Mesh &mesh, GridFunction &metric_gf)
Interpolates the metric&#39;s values at the nodes of metric_gf.
Definition: tmop.cpp:1259
Shape, ideal barrier metric, 3D.
Definition: tmop.hpp:377
TMOP_QualityMetric * metric
Definition: tmop.hpp:610
InvariantsEvaluator2D< double > ie
Definition: tmop.hpp:255
virtual ~TMOP_LimiterFunction()
Virtual destructor.
Definition: tmop.hpp:483
Base class for limiting functions to be used in class TMOP_Integrator.
Definition: tmop.hpp:466
virtual double EvalW(const DenseMatrix &Jpt) const =0
Evaluate the strain energy density function, W = W(Jpt).
virtual void Eval_d2(const Vector &x, const Vector &x0, double dist, DenseMatrix &d2) const =0
Returns the Hessian of the limiting function f(x, x0, d) with respect to x.
Skew metric, 3D.
Definition: tmop.hpp:107
virtual void AssembleH(const DenseMatrix &Jpt, 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...
Definition: tmop.cpp:698
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:625
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:162
DenseMatrix Jpr
Definition: tmop.hpp:640
virtual void ComputeElementTargets(int e_id, const FiniteElement &fe, const IntegrationRule &ir, DenseTensor &Jtr) const
Given an element and quadrature rule, computes ref-&gt;target transformation Jacobians for each quadratu...
Definition: tmop.cpp:806
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:481
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const =0
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
InvariantsEvaluator2D< double > ie
Definition: tmop.hpp:79
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:561
Data type dense matrix using column-major storage.
Definition: densemat.hpp:23
Shape, ideal barrier metric, 3D.
Definition: tmop.hpp:361
int Size() const
Returns the size of the vector.
Definition: vector.hpp:150
void SetVolumeScale(double vol_scale)
Used by target type IDEAL_SHAPE_EQUAL_SIZE. The default volume scale is 1.
Definition: tmop.hpp:589
Shape, ideal barrier metric, 2D.
Definition: tmop.hpp:152
double metric_normal
Definition: tmop.hpp:616
Volume metric, 3D.
Definition: tmop.hpp:393
Area, ideal barrier metric, 2D.
Definition: tmop.hpp:289
virtual void AssembleH(const DenseMatrix &Jpt, 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...
Definition: tmop.cpp:493
Coefficient * coeff1
Definition: tmop.hpp:614
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.hpp:128
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.hpp:113
InvariantsEvaluator2D< double > ie
Definition: tmop.hpp:187
InvariantsEvaluator2D< double > ie
Definition: tmop.hpp:274
const GridFunction * nodes
Definition: tmop.hpp:548
Aspect ratio metric, 2D.
Definition: tmop.hpp:122
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:687
Shape &amp; area, ideal barrier metric, 2D.
Definition: tmop.hpp:168
Default limiter function in TMOP_Integrator.
Definition: tmop.hpp:487
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:364
Shifted barrier form of metric 2 (shape, ideal barrier metric), 2D.
Definition: tmop.hpp:200
TMOP_Metric_352(double &t0)
Definition: tmop.hpp:450
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:517
void SetTargetJacobian(const DenseMatrix &_Jtr)
Specify the reference-element -&gt; target-element Jacobian matrix for the point of interest.
Definition: tmop.hpp:44
Shape, ideal barrier metric, 2D.
Definition: tmop.hpp:271
virtual ~TMOP_QualityMetric()
Definition: tmop.hpp:36
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:67
virtual void AssembleH(const DenseMatrix &Jpt, 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...
Definition: tmop.cpp:607
virtual void AssembleH(const DenseMatrix &Jpt, 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...
Definition: tmop.cpp:380
Abstract class for local mesh quality metrics in the target-matrix optimization paradigm (TMOP) by P...
Definition: tmop.hpp:24
InvariantsEvaluator3D< double > ie
Definition: tmop.hpp:396
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:654
void EnableLimiting(const GridFunction &n0, const GridFunction &dist, Coefficient &w0, TMOP_LimiterFunction *lfunc=NULL)
Adds a limiting term to the integrator (general version).
Definition: tmop.cpp:867
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:451
void SetCoefficient(Coefficient &w1)
Sets a scaling Coefficient for the quality metric term of the integrator.
Definition: tmop.hpp:663
Area metric, 2D.
Definition: tmop.hpp:235
virtual void AssembleH(const DenseMatrix &Jpt, 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...
Definition: tmop.cpp:322
TargetConstructor(TargetType ttype)
Constructor for use in serial.
Definition: tmop.hpp:566
virtual void AssembleH(const DenseMatrix &Jpt, 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...
Definition: tmop.cpp:351
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:224
Volume, ideal barrier metric, 3D.
Definition: tmop.hpp:409
InvariantsEvaluator2D< double > ie
Definition: tmop.hpp:310
void SetNodes(const GridFunction &n)
Set the nodes to be used in the target-matrix construction.
Definition: tmop.hpp:586
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:595
const TargetConstructor * targetC
Definition: tmop.hpp:611
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:460
virtual void AssembleH(const DenseMatrix &Jpt, 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...
Definition: tmop.cpp:410
InvariantsEvaluator2D< double > ie
Definition: tmop.hpp:204
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:107
virtual ~TargetConstructor()
Definition: tmop.hpp:579
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:423
InvariantsEvaluator2D< double > ie
Definition: tmop.hpp:329
virtual void AssembleElementVector(const FiniteElement &el, ElementTransformation &T, const Vector &elfun, Vector &elvect)
Perform the local action of the NonlinearFormIntegrator.
Definition: tmop.cpp:996
virtual void AssembleH(const DenseMatrix &Jpt, 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...
Definition: tmop.cpp:199
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:27
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:191
void ComputeAvgVolume() const
Definition: tmop.cpp:768
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:676
virtual double Eval(const Vector &x, const Vector &x0, double dist) const
Returns the limiting function, f(x, x0, d).
Definition: tmop.hpp:490
virtual ~TMOP_QuadraticLimiter()
Definition: tmop.hpp:514
Shape, ideal barrier metric, 2D.
Definition: tmop.hpp:219
TMOP_Metric_211(double epsilon=1e-4)
Definition: tmop.hpp:313
virtual void AssembleElementGrad(const FiniteElement &el, ElementTransformation &T, const Vector &elfun, DenseMatrix &elmat)
Assemble the local gradient matrix.
Definition: tmop.cpp:1088
Untangling metric, 2D.
Definition: tmop.hpp:306
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.hpp:98
const double eps
Definition: tmop.hpp:309
Area, ideal barrier metric, 2D.
Definition: tmop.hpp:252
void SetLimitingNodes(const GridFunction &n0)
Update the original/reference nodes used for limiting.
Definition: tmop.hpp:685
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:430
Base class Coefficient that may optionally depend on time.
Definition: coefficient.hpp:31
InvariantsEvaluator3D< double > ie
Definition: tmop.hpp:447
void subtract(const Vector &x, const Vector &y, Vector &z)
Definition: vector.cpp:385
InvariantsEvaluator2D< double > ie
Definition: tmop.hpp:222
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:43
virtual void AssembleH(const DenseMatrix &Jpt, 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...
Definition: tmop.cpp:174
double DistanceSquaredTo(const double *p) const
Compute the square of the Euclidean distance to another vector.
Definition: vector.hpp:517
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:303
TMOP_Metric_252(double &t0)
Note that t0 is stored by reference.
Definition: tmop.hpp:333
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:335
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:473
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:394
Skew metric, 2D.
Definition: tmop.hpp:92
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:184
virtual void AssembleH(const DenseMatrix &Jpt, 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...
Definition: tmop.cpp:241
virtual void AssembleH(const DenseMatrix &Jpt, 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...
Definition: tmop.cpp:633
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:402
DenseMatrix PMatI
Definition: tmop.hpp:640
InvariantsEvaluator3D< double > ie
Definition: tmop.hpp:412
Shape &amp; area metric, 2D.
Definition: tmop.hpp:184
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:232
void SetTransformation(ElementTransformation &)
The method SetTransformation() is hidden for TMOP_QualityMetrics, because it is not used...
Definition: tmop.hpp:32
virtual double Eval(const Vector &x, const Vector &x0, double d) const =0
Returns the limiting function, f(x, x0, d).
InvariantsEvaluator3D< double > ie
Definition: tmop.hpp:348
virtual void AssembleH(const DenseMatrix &Jpt, 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...
Definition: tmop.cpp:579
void ComputeNormalizationEnergies(const GridFunction &x, double &metric_energy, double &lim_energy)
Definition: tmop.cpp:1208
void Diag(double c, int n)
Creates n x n diagonal matrix with diagonal elements c.
Definition: densemat.cpp:2455
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:372
virtual double GetElementEnergy(const FiniteElement &el, ElementTransformation &T, const Vector &elfun)
Computes the integral of W(Jacobian(Trt)) over a target zone.
Definition: tmop.cpp:903
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:21
virtual void Eval_d1(const Vector &x, const Vector &x0, double dist, Vector &d1) const =0
Returns the gradient of the limiting function f(x, x0, d) with respect to x.
Aspect ratio metric, 3D.
Definition: tmop.hpp:137
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:646
InvariantsEvaluator2D< double > ie
Definition: tmop.hpp:155
virtual void AssembleH(const DenseMatrix &Jpt, 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...
Definition: tmop.cpp:526
InvariantsEvaluator2D< double > ie
Definition: tmop.hpp:171
Abstract class for hyperelastic models.
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:127
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:266
TargetConstructor(TargetType ttype, MPI_Comm mpicomm)
Constructor for use in parallel.
Definition: tmop.hpp:575
virtual void AssembleH(const DenseMatrix &Jpt, 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...
Definition: tmop.cpp:662
double epsilon(const Vector &x)
Definition: maxwell.cpp:90
TMOP_Integrator(TMOP_QualityMetric *m, TargetConstructor *tc)
Definition: tmop.hpp:648
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:511
Vector data type.
Definition: vector.hpp:48
DenseMatrix DSh
Definition: tmop.hpp:640
DenseMatrix DS
Definition: tmop.hpp:640
TMOP_Metric_022(double &t0)
Definition: tmop.hpp:207
InvariantsEvaluator3D< double > ie
Definition: tmop.hpp:380
const GridFunction * nodes0
Definition: tmop.hpp:621
TargetType
Target-matrix construction algorithms supported by this class.
Definition: tmop.hpp:530
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:343
virtual void AssembleH(const DenseMatrix &Jpt, 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...
Definition: tmop.cpp:465
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:314
virtual void Eval_d1(const Vector &x, const Vector &x0, double dist, Vector &d1) const
Returns the gradient of the limiting function f(x, x0, d) with respect to x.
Definition: tmop.hpp:497
virtual void AssembleH(const DenseMatrix &Jpt, 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...
Definition: tmop.hpp:131
Class for parallel grid function.
Definition: pgridfunc.hpp:32
Base class representing target-matrix construction algorithms for mesh optimization via the target-ma...
Definition: tmop.hpp:526
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:168
Rank 3 tensor (array of matrices)
Definition: densemat.hpp:656
TMOP_LimiterFunction * lim_func
Definition: tmop.hpp:626
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:617
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:601
InvariantsEvaluator3D< double > ie
Definition: tmop.hpp:364
DenseMatrix Jpt
Definition: tmop.hpp:640
void ParEnableNormalization(const ParGridFunction &x)
Definition: tmop.cpp:1198
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:732
virtual void AssembleH(const DenseMatrix &Jpt, 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...
Definition: tmop.hpp:101
DenseMatrix P
Definition: tmop.hpp:640
virtual void AssembleH(const DenseMatrix &Jpt, 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...
Definition: tmop.cpp:744
bool Parallel() const
Definition: tmop.hpp:555
Shape, ideal barrier metric, 3D.
Definition: tmop.hpp:345
Metric without a type, 2D.
Definition: tmop.hpp:76
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:724
virtual void AssembleH(const DenseMatrix &Jpt, 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...
Definition: tmop.cpp:439
A TMOP integrator class based on any given TMOP_QualityMetric and TargetConstructor.
Definition: tmop.hpp:607