MFEM  v4.2.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-2020, 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_TMOP_HPP
13 #define MFEM_TMOP_HPP
14 
15 #include "../linalg/invariants.hpp"
16 #include "nonlininteg.hpp"
17 
18 namespace mfem
19 {
20 
21 /** @brief Abstract class for local mesh quality metrics in the target-matrix
22  optimization paradigm (TMOP) by P. Knupp et al. */
24 {
25 protected:
26  const DenseMatrix *Jtr; /**< Jacobian of the reference-element to
27  target-element transformation. */
28 
29  /** @brief The method SetTransformation() is hidden for TMOP_QualityMetric%s,
30  because it is not used. */
32 
33 public:
34  TMOP_QualityMetric() : Jtr(NULL) { }
35  virtual ~TMOP_QualityMetric() { }
36 
37  /** @brief Specify the reference-element -> target-element Jacobian matrix
38  for the point of interest.
39 
40  The specified Jacobian matrix, #Jtr, can be used by metrics that cannot
41  be written just as a function of the target->physical Jacobian matrix,
42  Jpt. */
43  void SetTargetJacobian(const DenseMatrix &_Jtr) { Jtr = &_Jtr; }
44 
45  /** @brief Evaluate the strain energy density function, W = W(Jpt).
46  @param[in] Jpt Represents the target->physical transformation
47  Jacobian matrix. */
48  virtual double EvalW(const DenseMatrix &Jpt) const = 0;
49 
50  /** @brief Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
51  @param[in] Jpt Represents the target->physical transformation
52  Jacobian matrix.
53  @param[out] P The evaluated 1st Piola-Kirchhoff stress tensor. */
54  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const = 0;
55 
56  /** @brief Evaluate the derivative of the 1st Piola-Kirchhoff stress tensor
57  and assemble its contribution to the local gradient matrix 'A'.
58  @param[in] Jpt Represents the target->physical transformation
59  Jacobian matrix.
60  @param[in] DS Gradient of the basis matrix (dof x dim).
61  @param[in] weight Quadrature weight coefficient for the point.
62  @param[in,out] A Local gradient matrix where the contribution from this
63  point will be added.
64 
65  Computes weight * d(dW_dxi)_d(xj) at the current point, for all i and j,
66  where x1 ... xn are the FE dofs. This function is usually defined using
67  the matrix invariants and their derivatives.
68  */
69  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
70  const double weight, DenseMatrix &A) const = 0;
71 };
72 
73 
74 /// Metric without a type, 2D
76 {
77 protected:
79 
80 public:
81  // W = |J|^2.
82  virtual double EvalW(const DenseMatrix &Jpt) const;
83 
84  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
85 
86  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
87  const double weight, DenseMatrix &A) const;
88 };
89 
90 /// Skew metric, 2D.
92 {
93 public:
94  // W = 0.5 (1 - cos(angle_Jpr - angle_Jtr)).
95  virtual double EvalW(const DenseMatrix &Jpt) const;
96 
97  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
98  { MFEM_ABORT("Not implemented"); }
99 
100  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
101  const double weight, DenseMatrix &A) const
102  { MFEM_ABORT("Not implemented"); }
103 };
104 
105 /// Skew metric, 3D.
107 {
108 public:
109  // W = 1/6 (3 - sum_i cos(angle_Jpr_i - angle_Jtr_i)), i = 1..3.
110  virtual double EvalW(const DenseMatrix &Jpt) const;
111 
112  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
113  { MFEM_ABORT("Not implemented"); }
114 
115  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
116  const double weight, DenseMatrix &A) const
117  { MFEM_ABORT("Not implemented"); }
118 };
119 
120 /// Aspect ratio metric, 2D.
122 {
123 public:
124  // W = 0.5 (ar_Jpr/ar_Jtr + ar_Jtr/ar_Jpr) - 1.
125  virtual double EvalW(const DenseMatrix &Jpt) const;
126 
127  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
128  { MFEM_ABORT("Not implemented"); }
129 
130  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
131  const double weight, DenseMatrix &A) const
132  { MFEM_ABORT("Not implemented"); }
133 };
134 
135 /// Aspect ratio metric, 3D.
137 {
138 public:
139  // W = 1/3 sum [0.5 (ar_Jpr_i/ar_Jtr_i + ar_Jtr_i/ar_Jpr_i) - 1], i = 1..3.
140  virtual double EvalW(const DenseMatrix &Jpt) const;
141 
142  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
143  { MFEM_ABORT("Not implemented"); }
144 
145  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
146  const double weight, DenseMatrix &A) const
147  { MFEM_ABORT("Not implemented"); }
148 };
149 
150 /// Shape+Size+Orientation metric, 2D.
152 {
153 public:
154  // W = 0.5 (1 - cos(theta_Jpr - theta_Jtr)).
155  virtual double EvalW(const DenseMatrix &Jpt) const;
156 
157  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
158  { MFEM_ABORT("Not implemented"); }
159 
160  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
161  const double weight, DenseMatrix &A) const
162  { MFEM_ABORT("Not implemented"); }
163 };
164 
165 /// Shape, ideal barrier metric, 2D
167 {
168 protected:
170 
171 public:
172  // W = 0.5|J|^2 / det(J) - 1.
173  virtual double EvalW(const DenseMatrix &Jpt) const;
174 
175  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
176 
177  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
178  const double weight, DenseMatrix &A) const;
179 };
180 
181 /// Shape & area, ideal barrier metric, 2D
183 {
184 protected:
186 
187 public:
188  // W = |J - J^-t|^2.
189  virtual double EvalW(const DenseMatrix &Jpt) const;
190 
191  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
192 
193  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
194  const double weight, DenseMatrix &A) const;
195 };
196 
197 /// Shape & area metric, 2D
199 {
200 protected:
202 
203 public:
204  // W = det(J) * |J - J^-t|^2.
205  virtual double EvalW(const DenseMatrix &Jpt) const;
206 
207  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
208 
209  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
210  const double weight, DenseMatrix &A) const;
211 };
212 
213 /// Shifted barrier form of metric 2 (shape, ideal barrier metric), 2D
215 {
216 protected:
217  double &tau0;
219 
220 public:
221  TMOP_Metric_022(double &t0): tau0(t0) {}
222 
223  // W = 0.5(|J|^2 - 2det(J)) / (det(J) - tau0).
224  virtual double EvalW(const DenseMatrix &Jpt) const;
225 
226  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
227 
228  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
229  const double weight, DenseMatrix &A) const;
230 };
231 
232 /// Shape, ideal barrier metric, 2D
234 {
235 protected:
237 
238 public:
239  // W = 0.5|J^t J|^2 / det(J)^2 - 1.
240  virtual double EvalW(const DenseMatrix &Jpt) const;
241 
242  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
243 
244  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
245  const double weight, DenseMatrix &A) const;
246 };
247 
248 /// Area metric, 2D
250 {
251 protected:
253 
254 public:
255  // W = (det(J) - 1)^2.
256  virtual double EvalW(const DenseMatrix &Jpt) const;
257 
258  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
259 
260  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
261  const double weight, DenseMatrix &A) const;
262 
263 };
264 
265 /// Area, ideal barrier metric, 2D
267 {
268 protected:
270 
271 public:
272  // W = 0.5( sqrt(det(J)) - 1 / sqrt(det(J)) )^2
273  // = 0.5( det(J) - 1 )^2 / det(J)
274  // = 0.5( det(J) + 1/det(J) ) - 1.
275  virtual double EvalW(const DenseMatrix &Jpt) const;
276 
277  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
278 
279  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
280  const double weight, DenseMatrix &A) const;
281 
282 };
283 
284 /// Shape, ideal barrier metric, 2D
286 {
287 protected:
289 
290 public:
291  // W = |J^t J|^2 / det(J)^2 - 2|J|^2 / det(J) + 2
292  // = I1b (I1b - 2).
293  virtual double EvalW(const DenseMatrix &Jpt) const;
294 
295  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
296 
297  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
298  const double weight, DenseMatrix &A) const;
299 
300 };
301 
302 /// Area, ideal barrier metric, 2D
304 {
305 protected:
307 
308 public:
309  // W = 0.5(det(J) - 1 / det(J))^2.
310  virtual double EvalW(const DenseMatrix &Jpt) const;
311 
312  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
313 
314  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
315  const double weight, DenseMatrix &A) const;
316 
317 };
318 
319 /// Shape & orientation metric, 2D.
321 {
322 public:
323  // W = |T-T'|^2, where T'= |T|*I/sqrt(2).
324  virtual double EvalW(const DenseMatrix &Jpt) const;
325 
326  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
327  { MFEM_ABORT("Not implemented"); }
328 
329  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
330  const double weight, DenseMatrix &A) const
331  { MFEM_ABORT("Not implemented"); }
332 };
333 
334 /// Untangling metric, 2D
336 {
337 protected:
338  const double eps;
340 
341 public:
342  TMOP_Metric_211(double epsilon = 1e-4) : eps(epsilon) { }
343 
344  // W = (det(J) - 1)^2 - det(J) + sqrt(det(J)^2 + eps).
345  virtual double EvalW(const DenseMatrix &Jpt) const;
346 
347  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
348 
349  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
350  const double weight, DenseMatrix &A) const;
351 };
352 
353 /// Shifted barrier form of metric 56 (area, ideal barrier metric), 2D
355 {
356 protected:
357  double &tau0;
359 
360 public:
361  /// Note that @a t0 is stored by reference
362  TMOP_Metric_252(double &t0): tau0(t0) {}
363 
364  // W = 0.5(det(J) - 1)^2 / (det(J) - tau0).
365  virtual double EvalW(const DenseMatrix &Jpt) const;
366 
367  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
368 
369  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
370  const double weight, DenseMatrix &A) const;
371 };
372 
373 /// Shape, ideal barrier metric, 3D
375 {
376 protected:
378 
379 public:
380  // W = |J| |J^-1| / 3 - 1.
381  virtual double EvalW(const DenseMatrix &Jpt) const;
382 
383  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
384 
385  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
386  const double weight, DenseMatrix &A) const;
387 };
388 
389 /// Shape, ideal barrier metric, 3D
391 {
392 protected:
394 
395 public:
396  // W = |J|^2 |J^-1|^2 / 9 - 1.
397  virtual double EvalW(const DenseMatrix &Jpt) const;
398 
399  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
400 
401  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
402  const double weight, DenseMatrix &A) const;
403 };
404 
405 /// Shape, ideal barrier metric, 3D
407 {
408 protected:
410 
411 public:
412  // W = |J|^2 / 3 * det(J)^(2/3) - 1.
413  virtual double EvalW(const DenseMatrix &Jpt) const;
414 
415  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
416 
417  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
418  const double weight, DenseMatrix &A) const;
419 };
420 
421 /// Volume metric, 3D
423 {
424 protected:
426 
427 public:
428  // W = (det(J) - 1)^2.
429  virtual double EvalW(const DenseMatrix &Jpt) const;
430 
431  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
432 
433  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
434  const double weight, DenseMatrix &A) const;
435 };
436 
437 /// Volume, ideal barrier metric, 3D
439 {
440 protected:
442 
443 public:
444  // W = 0.5( sqrt(det(J)) - 1 / sqrt(det(J)) )^2
445  // = 0.5( det(J) - 1 )^2 / det(J)
446  // = 0.5( det(J) + 1/det(J) ) - 1.
447  virtual double EvalW(const DenseMatrix &Jpt) const;
448 
449  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
450 
451  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
452  const double weight, DenseMatrix &A) const;
453 };
454 
455 /// Shape & volume, ideal barrier metric, 3D
457 {
458 protected:
460 
461 public:
462  // W = |J - J^-t|^2.
463  virtual double EvalW(const DenseMatrix &Jpt) const;
464 
465  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
466 
467  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
468  const double weight, DenseMatrix &A) const;
469 };
470 
471 /// Shifted barrier form of 3D metric 16 (volume, ideal barrier metric), 3D
473 {
474 protected:
475  double &tau0;
477 
478 public:
479  TMOP_Metric_352(double &t0): tau0(t0) {}
480 
481  // W = 0.5(det(J) - 1)^2 / (det(J) - tau0).
482  virtual double EvalW(const DenseMatrix &Jpt) const;
483 
484  virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const;
485 
486  virtual void AssembleH(const DenseMatrix &Jpt, const DenseMatrix &DS,
487  const double weight, DenseMatrix &A) const;
488 };
489 
490 
491 /// Base class for limiting functions to be used in class TMOP_Integrator.
492 /** This class represents a scalar function f(x, x0, d), where x and x0 are
493  positions in physical space, and d is a reference physical distance
494  associated with the point x0. */
496 {
497 public:
498  /// Returns the limiting function, f(x, x0, d).
499  virtual double Eval(const Vector &x, const Vector &x0, double d) const = 0;
500 
501  /** @brief Returns the gradient of the limiting function f(x, x0, d) with
502  respect to x. */
503  virtual void Eval_d1(const Vector &x, const Vector &x0, double dist,
504  Vector &d1) const = 0;
505 
506  /** @brief Returns the Hessian of the limiting function f(x, x0, d) with
507  respect to x. */
508  virtual void Eval_d2(const Vector &x, const Vector &x0, double dist,
509  DenseMatrix &d2) const = 0;
510 
511  /// Virtual destructor.
512  virtual ~TMOP_LimiterFunction() { }
513 };
514 
515 /// Default limiter function in TMOP_Integrator.
517 {
518 public:
519  virtual double Eval(const Vector &x, const Vector &x0, double dist) const
520  {
521  MFEM_ASSERT(x.Size() == x0.Size(), "Bad input.");
522 
523  return 0.5 * x.DistanceSquaredTo(x0) / (dist * dist);
524  }
525 
526  virtual void Eval_d1(const Vector &x, const Vector &x0, double dist,
527  Vector &d1) const
528  {
529  MFEM_ASSERT(x.Size() == x0.Size(), "Bad input.");
530 
531  d1.SetSize(x.Size());
532  subtract(1.0 / (dist * dist), x, x0, d1);
533  }
534 
535  virtual void Eval_d2(const Vector &x, const Vector &x0, double dist,
536  DenseMatrix &d2) const
537  {
538  MFEM_ASSERT(x.Size() == x0.Size(), "Bad input.");
539 
540  d2.Diag(1.0 / (dist * dist), x.Size());
541  }
542 
544 };
545 
546 class FiniteElementCollection;
547 class FiniteElementSpace;
548 class ParFiniteElementSpace;
549 
551 {
552 protected:
553  // Owned.
556 
557 #ifdef MFEM_USE_MPI
558  // Owned.
561 #endif
562 
563  int dim, ncomp;
564 
565 public:
566  AdaptivityEvaluator() : mesh(NULL), fes(NULL)
567  {
568 #ifdef MFEM_USE_MPI
569  pmesh = NULL;
570  pfes = NULL;
571 #endif
572  }
573  virtual ~AdaptivityEvaluator();
574 
575  /** Specifies the Mesh and FiniteElementCollection of the solution that will
576  be evaluated. The given mesh will be copied into the internal object. */
577  void SetSerialMetaInfo(const Mesh &m,
578  const FiniteElementCollection &fec, int num_comp);
579 
580 #ifdef MFEM_USE_MPI
581  /// Parallel version of SetSerialMetaInfo.
582  void SetParMetaInfo(const ParMesh &m,
583  const FiniteElementCollection &fec, int num_comp);
584 #endif
585 
586  // TODO use GridFunctions to make clear it's on the ldofs?
587  virtual void SetInitialField(const Vector &init_nodes,
588  const Vector &init_field) = 0;
589 
590  virtual void ComputeAtNewPosition(const Vector &new_nodes,
591  Vector &new_field) = 0;
592 };
593 
594 /** @brief Base class representing target-matrix construction algorithms for
595  mesh optimization via the target-matrix optimization paradigm (TMOP). */
596 /** This class is used by class TMOP_Integrator to construct the target Jacobian
597  matrices (reference-element to target-element) at quadrature points. It
598  supports a set of algorithms chosen by the #TargetType enumeration.
599 
600  New target-matrix construction algorithms can be defined by deriving new
601  classes and overriding the methods ComputeElementTargets() and
602  ContainsVolumeInfo(). */
604 {
605 public:
606  /// Target-matrix construction algorithms supported by this class.
608  {
610  Ideal shape, unit size; the nodes are not used. */
612  Ideal shape, equal size/volume; the given nodes define the total target
613  volume; for each mesh element, the target volume is the average volume
614  multiplied by the volume scale, set with SetVolumeScale(). */
616  Ideal shape, given size/volume; the given nodes define the target
617  volume at all quadrature points. */
619  Given shape, given size/volume; the given nodes define the exact target
620  Jacobian matrix at all quadrature points. */
622  Full target tensor is specified at every quadrature point. */
623  };
624 
625 protected:
626  // Nodes that are used in ComputeElementTargets(), depending on target_type.
627  const GridFunction *nodes; // not owned
628  mutable double avg_volume;
629  double volume_scale;
631 
632 #ifdef MFEM_USE_MPI
633  MPI_Comm comm;
634  bool Parallel() const { return (comm != MPI_COMM_NULL); }
635 #else
636  bool Parallel() const { return false; }
637 #endif
638 
639  // should be called only if avg_volume == 0.0, i.e. avg_volume is not
640  // computed yet
641  void ComputeAvgVolume() const;
642 
643 public:
644  /// Constructor for use in serial
646  : nodes(NULL), avg_volume(), volume_scale(1.0), target_type(ttype)
647  {
648 #ifdef MFEM_USE_MPI
649  comm = MPI_COMM_NULL;
650 #endif
651  }
652 #ifdef MFEM_USE_MPI
653  /// Constructor for use in parallel
654  TargetConstructor(TargetType ttype, MPI_Comm mpicomm)
655  : nodes(NULL), avg_volume(), volume_scale(1.0), target_type(ttype),
656  comm(mpicomm) { }
657 #endif
658  virtual ~TargetConstructor() { }
659 
660  /** @brief Set the nodes to be used in the target-matrix construction.
661 
662  This method should be called every time the target nodes are updated
663  externally and recomputation of the target average volume is needed. The
664  nodes are used by all target types except IDEAL_SHAPE_UNIT_SIZE. */
665  void SetNodes(const GridFunction &n) { nodes = &n; avg_volume = 0.0; }
666 
667  /// Used by target type IDEAL_SHAPE_EQUAL_SIZE. The default volume scale is 1.
668  void SetVolumeScale(double vol_scale) { volume_scale = vol_scale; }
669 
670  /// Checks if the target matrices contain non-trivial size specification.
671  virtual bool ContainsVolumeInfo() const;
672 
673  /** @brief Given an element and quadrature rule, computes ref->target
674  transformation Jacobians for each quadrature point in the element.
675  The physical positions of the element's nodes are given by @a elfun. */
676  virtual void ComputeElementTargets(int e_id, const FiniteElement &fe,
677  const IntegrationRule &ir,
678  const Vector &elfun,
679  DenseTensor &Jtr) const;
680 
681  virtual void ComputeElementTargetsGradient(const IntegrationRule &ir,
682  const Vector &elfun,
684  DenseTensor &dJtr) const;
685 };
686 
688 {
689 public:
690  explicit TMOPMatrixCoefficient(int dim) : MatrixCoefficient(dim, dim) { }
691 
692  /** @brief Evaluate the derivative of the matrix coefficient with respect to
693  @a comp in the element described by @a T at the point @a ip, storing the
694  result in @a K. */
695  virtual void EvalGrad(DenseMatrix &K, ElementTransformation &T,
696  const IntegrationPoint &ip, int comp) = 0;
697 
699 };
700 
702 {
703 protected:
704  // Analytic target specification.
708 
709 public:
711  : TargetConstructor(ttype),
712  scalar_tspec(NULL), vector_tspec(NULL), matrix_tspec(NULL) { }
713 
714  virtual void SetAnalyticTargetSpec(Coefficient *sspec,
715  VectorCoefficient *vspec,
716  TMOPMatrixCoefficient *mspec);
717 
718  /** @brief Given an element and quadrature rule, computes ref->target
719  transformation Jacobians for each quadrature point in the element.
720  The physical positions of the element's nodes are given by @a elfun. */
721  virtual void ComputeElementTargets(int e_id, const FiniteElement &fe,
722  const IntegrationRule &ir,
723  const Vector &elfun,
724  DenseTensor &Jtr) const;
725 
726  virtual void ComputeElementTargetsGradient(const IntegrationRule &ir,
727  const Vector &elfun,
729  DenseTensor &dJtr) const;
730 };
731 
732 #ifdef MFEM_USE_MPI
733 class ParGridFunction;
734 #endif
735 
737 {
738 protected:
739  // Discrete target specification.
740  // Data is owned, updated by UpdateTargetSpecification.
742  Vector tspec; //eta(x)
744  Vector tspec_pert1h; //eta(x+h)
745  Vector tspec_pert2h; //eta(x+2*h)
746  Vector tspec_pertmix; //eta(x+h,y+h)
747  // The order inside these perturbation vectors (e.g. in 2D) is
748  // eta1(x+h,y), eta2(x+h,y) ... etan(x+h,y), eta1(x,y+h), eta2(x,y+h) ...
749  // same for tspec_pert2h and tspec_pertmix.
750 
751  // Components of Target Jacobian at each quadrature point of an element. This
752  // is required for computation of the derivative using chain rule.
754 
755  // Note: do not use the Nodes of this space as they may not be on the
756  // positions corresponding to the values of tspec.
759 
760  // These flags can be used by outside functions to avoid recomputing the
761  // tspec and tspec_perth fields again on the same mesh.
763 
764  // Evaluation of the discrete target specification on different meshes.
765  // Owned.
767 
768  void SetDiscreteTargetBase(const GridFunction &tspec_);
769  void SetTspecAtIndex(int idx, const GridFunction &tspec_);
771 #ifdef MFEM_USE_MPI
772  void SetTspecAtIndex(int idx, const ParGridFunction &tspec_);
774 #endif
775 
776 public:
778  : TargetConstructor(ttype),
779  ncomp(0),
780  sizeidx(-1), skewidx(-1), aspectratioidx(-1), orientationidx(-1),
782  tspec_fes(NULL), tspec_fesv(NULL),
783  good_tspec(false), good_tspec_grad(false), good_tspec_hess(false),
784  adapt_eval(NULL) { }
785 
787  {
788  delete adapt_eval;
789  delete tspec_fes;
790  delete tspec_fesv;
791  }
792 
793  /** @name Target specification methods.
794  The following methods are used to specify geometric parameters of the
795  targets when these parameters are given by discrete FE functions.
796  Note that every GridFunction given to the Set methods must use a
797  H1_FECollection of the same order. The number of components must
798  correspond to the type of geometric parameter and dimension.
799 
800  @param[in] tspec_ Input values of a geometric parameter. Note that
801  the methods in this class support only functions that
802  use H1_FECollection collection of the same order. */
803  ///@{
804  virtual void SetSerialDiscreteTargetSpec(const GridFunction &tspec_);
805  virtual void SetSerialDiscreteTargetSize(const GridFunction &tspec_);
806  virtual void SetSerialDiscreteTargetSkew(const GridFunction &tspec_);
807  virtual void SetSerialDiscreteTargetAspectRatio(const GridFunction &tspec_);
808  virtual void SetSerialDiscreteTargetOrientation(const GridFunction &tspec_);
809 #ifdef MFEM_USE_MPI
810  virtual void SetParDiscreteTargetSpec(const ParGridFunction &tspec_);
811  virtual void SetParDiscreteTargetSize(const ParGridFunction &tspec_);
812  virtual void SetParDiscreteTargetSkew(const ParGridFunction &tspec_);
813  virtual void SetParDiscreteTargetAspectRatio(const ParGridFunction &tspec_);
814  virtual void SetParDiscreteTargetOrientation(const ParGridFunction &tspec_);
815 #endif
816  ///@}
817 
818  /// Used in combination with the Update methods to avoid extra computations.
821 
822  /** Used to update the target specification after the mesh has changed. The
823  new mesh positions are given by new_x. If @a use_flags is true, repeated
824  calls won't do anything until ResetUpdateFlags() is called. */
825  void UpdateTargetSpecification(const Vector &new_x, bool use_flag = false);
826 
827  void UpdateTargetSpecification(Vector &new_x, Vector &IntData);
828 
831  int nodenum, int idir,
832  const Vector &IntData);
833 
835 
836  /** Used for finite-difference based computations. Computes the target
837  specifications after a mesh perturbation in x or y direction.
838  If @a use_flags is true, repeated calls won't do anything until
839  ResetUpdateFlags() is called. */
840  void UpdateGradientTargetSpecification(const Vector &x, double dx,
841  bool use_flag = false);
842  /** Used for finite-difference based computations. Computes the target
843  specifications after two mesh perturbations in x and/or y direction.
844  If @a use_flags is true, repeated calls won't do anything until
845  ResetUpdateFlags() is called. */
846  void UpdateHessianTargetSpecification(const Vector &x, double dx,
847  bool use_flag = false);
848 
850  {
851  if (adapt_eval) { delete adapt_eval; }
852  adapt_eval = ae;
853  }
854 
855  const Vector &GetTspecPert1H() { return tspec_pert1h; }
856  const Vector &GetTspecPert2H() { return tspec_pert2h; }
857  const Vector &GetTspecPertMixH() { return tspec_pertmix; }
858 
859  /** @brief Given an element and quadrature rule, computes ref->target
860  transformation Jacobians for each quadrature point in the element.
861  The physical positions of the element's nodes are given by @a elfun.
862  Note that this function assumes that UpdateTargetSpecification() has
863  been called with the position vector corresponding to @a elfun. */
864  virtual void ComputeElementTargets(int e_id, const FiniteElement &fe,
865  const IntegrationRule &ir,
866  const Vector &elfun,
867  DenseTensor &Jtr) const;
868 
869  virtual void ComputeElementTargetsGradient(const IntegrationRule &ir,
870  const Vector &elfun,
872  DenseTensor &dJtr) const;
873 };
874 
875 class TMOPNewtonSolver;
876 
877 /** @brief A TMOP integrator class based on any given TMOP_QualityMetric and
878  TargetConstructor.
879 
880  Represents @f$ \int W(Jpt) dx @f$ over a target zone, where W is the
881  metric's strain energy density function, and Jpt is the Jacobian of the
882  target->physical coordinates transformation. The virtual target zone is
883  defined by the TargetConstructor. */
885 {
886 protected:
887  friend class TMOPNewtonSolver;
888  friend class TMOPComboIntegrator;
889 
890  TMOP_QualityMetric *metric; // not owned
891  const TargetConstructor *targetC; // not owned
892 
893  // Custom integration rules.
896 
897  // Weight Coefficient multiplying the quality metric term.
898  Coefficient *coeff1; // not owned, if NULL -> coeff1 is 1.
899  // Normalization factor for the metric term.
901 
902  // Nodes and weight Coefficient used for "limiting" the TMOP_Integrator.
903  // These are both NULL when there is no limiting.
904  // The class doesn't own nodes0 and coeff0.
907  // Limiting reference distance. Not owned.
909  // Limiting function. Owned.
911  // Normalization factor for the limiting term.
912  double lim_normal;
913 
914  // Adaptive limiting.
915  const GridFunction *zeta_0; // Not owned.
916  GridFunction *zeta; // Owned. Updated by adapt_eval.
917  Coefficient *coeff_zeta; // Not owned.
919 
921 
922  // Parameters for FD-based Gradient & Hessian calculation.
923  bool fdflag;
924  double dx;
925  double dxscale;
926  // Specifies that ComputeElementTargets is being called by a FD function.
927  // It's used to skip terms that have exact derivative calculations.
929  // Compute the exact action of the Integrator (includes derivative of the
930  // target with respect to spatial position)
932 
935 
936  // Jrt: the inverse of the ref->target Jacobian, Jrt = Jtr^{-1}.
937  // Jpr: the ref->physical transformation Jacobian, Jpr = PMatI^t DS.
938  // Jpt: the target->physical transformation Jacobian, Jpt = Jpr Jrt.
939  // P: represents dW_d(Jtp) (dim x dim).
940  // DSh: gradients of reference shape functions (dof x dim).
941  // DS: gradients of the shape functions in the target configuration,
942  // DS = DSh Jrt (dof x dim).
943  // PMatI: current coordinates of the nodes (dof x dim).
944  // PMat0: reshaped view into the local element contribution to the operator
945  // output - the result of AssembleElementVector() (dof x dim).
947 
949  double &metric_energy, double &lim_energy);
950 
951 
954  const Vector &elfun, Vector &elvect);
955 
958  const Vector &elfun, DenseMatrix &elmat);
959 
960  void AssembleElementVectorFD(const FiniteElement &el,
962  const Vector &elfun, Vector &elvect);
963 
964  // Assumes that AssembleElementVectorFD has been called.
965  void AssembleElementGradFD(const FiniteElement &el,
967  const Vector &elfun, DenseMatrix &elmat);
968 
969  void AssembleElemVecAdaptLim(const FiniteElement &el, const Vector &weights,
971  const IntegrationRule &ir, DenseMatrix &m);
972  void AssembleElemGradAdaptLim(const FiniteElement &el, const Vector &weights,
974  const IntegrationRule &ir, DenseMatrix &m);
975 
976  double GetFDDerivative(const FiniteElement &el,
978  Vector &elfun, const int nodenum,const int idir,
979  const double baseenergy, bool update_stored);
980 
981  /** @brief Determines the perturbation, h, for FD-based approximation. */
982  void ComputeFDh(const Vector &x, const FiniteElementSpace &fes);
983 #ifdef MFEM_USE_MPI
984  void ComputeFDh(const Vector &x, const ParFiniteElementSpace &pfes);
985 #endif
986  void ComputeMinJac(const Vector &x, const FiniteElementSpace &fes);
987 
988  void UpdateAfterMeshChange(const Vector &new_x);
989 
991  {
992  nodes0 = NULL; coeff0 = NULL; lim_dist = NULL; lim_func = NULL;
993  }
994 
996  {
997  if (IntegRules)
998  {
999  return IntegRules->Get(el.GetGeomType(), integ_order);
1000  }
1001  return (IntRule) ? *IntRule
1002  /* */ : IntRules.Get(el.GetGeomType(), 2*el.GetOrder() + 3);
1003  }
1005  {
1006  // TODO the energy most likely needs less integration points.
1007  return EnergyIntegrationRule(el);
1008  }
1010  {
1011  // TODO the action and energy most likely need less integration points.
1012  return EnergyIntegrationRule(el);
1013  }
1014 
1015 public:
1016  /** @param[in] m TMOP_QualityMetric that will be integrated (not owned).
1017  @param[in] tc Target-matrix construction algorithm to use (not owned). */
1019  : metric(m), targetC(tc), IntegRules(NULL), integ_order(-1),
1020  coeff1(NULL), metric_normal(1.0),
1021  nodes0(NULL), coeff0(NULL),
1022  lim_dist(NULL), lim_func(NULL), lim_normal(1.0),
1023  zeta_0(NULL), zeta(NULL), coeff_zeta(NULL), adapt_eval(NULL),
1024  discr_tc(dynamic_cast<DiscreteAdaptTC *>(tc)),
1025  fdflag(false), dxscale(1.0e3), fd_call_flag(false), exact_action(false)
1026  { }
1027 
1028  ~TMOP_Integrator();
1029 
1030  /// Prescribe a set of integration rules; relevant for mixed meshes.
1031  /** This function has priority over SetIntRule(), if both are called. */
1032  void SetIntegrationRules(IntegrationRules &irules, int order)
1033  {
1034  IntegRules = &irules;
1035  integ_order = order;
1036  }
1037 
1038  /// Sets a scaling Coefficient for the quality metric term of the integrator.
1039  /** With this addition, the integrator becomes
1040  @f$ \int w1 W(Jpt) dx @f$.
1041 
1042  Note that the Coefficient is evaluated in the physical configuration and
1043  not in the target configuration which may be undefined. */
1044  void SetCoefficient(Coefficient &w1) { coeff1 = &w1; }
1045 
1046  /** @brief Limiting of the mesh displacements (general version).
1047 
1048  Adds the term @f$ \int w_0 f(x, x_0, d) dx @f$, where f is a measure of
1049  the displacement between x and x_0, given the max allowed displacement d.
1050 
1051  @param[in] n0 Original mesh node coordinates (x0 above).
1052  @param[in] dist Allowed displacement in physical space (d above).
1053  @param[in] w0 Coefficient scaling the limiting integral.
1054  @param[in] lfunc TMOP_LimiterFunction defining the function f. If
1055  NULL, a TMOP_QuadraticLimiter will be used. The
1056  TMOP_Integrator assumes ownership of this pointer. */
1057  void EnableLimiting(const GridFunction &n0, const GridFunction &dist,
1058  Coefficient &w0, TMOP_LimiterFunction *lfunc = NULL);
1059 
1060  /** @brief Adds a limiting term to the integrator with limiting distance
1061  function (@a dist in the general version of the method) equal to 1. */
1062  void EnableLimiting(const GridFunction &n0, Coefficient &w0,
1063  TMOP_LimiterFunction *lfunc = NULL);
1064 
1065  /** @brief Restriction of the node positions to certain regions.
1066 
1067  Adds the term @f$ \int c (z(x) - z_0(x_0))^2 @f$, where z0(x0) is a given
1068  function on the starting mesh, and z(x) is its image on the new mesh.
1069  Minimizing this, means that a node at x0 is allowed to move to a
1070  position x(x0) only if z(x) ~ z0(x0).
1071  Such term can be used for tangential mesh relaxation.
1072 
1073  @param[in] z0 Function z0 that controls the adaptive limiting.
1074  @param[in] coeff Coefficient c for the above integral.
1075  @param[in] ae AdaptivityEvaluator to compute z(x) from z0(x0). */
1076  void EnableAdaptiveLimiting(const GridFunction &z0, Coefficient &coeff,
1077  AdaptivityEvaluator &ae);
1078 #ifdef MFEM_USE_MPI
1079  /// Parallel support for adaptive limiting.
1080  void EnableAdaptiveLimiting(const ParGridFunction &z0, Coefficient &coeff,
1081  AdaptivityEvaluator &ae);
1082 #endif
1083 
1084  /// Update the original/reference nodes used for limiting.
1085  void SetLimitingNodes(const GridFunction &n0) { nodes0 = &n0; }
1086 
1087  /** @brief Computes the integral of W(Jacobian(Trt)) over a target zone.
1088  @param[in] el Type of FiniteElement.
1089  @param[in] T Mesh element transformation.
1090  @param[in] elfun Physical coordinates of the zone. */
1091  virtual double GetElementEnergy(const FiniteElement &el,
1093  const Vector &elfun);
1094 
1095  virtual void AssembleElementVector(const FiniteElement &el,
1097  const Vector &elfun, Vector &elvect);
1098 
1099  virtual void AssembleElementGrad(const FiniteElement &el,
1101  const Vector &elfun, DenseMatrix &elmat);
1102 
1104 
1105  /** @brief Computes the normalization factors of the metric and limiting
1106  integrals using the mesh position given by @a x. */
1107  void EnableNormalization(const GridFunction &x);
1108 #ifdef MFEM_USE_MPI
1109  void ParEnableNormalization(const ParGridFunction &x);
1110 #endif
1111 
1112  /** @brief Enables FD-based approximation and computes dx. */
1113  void EnableFiniteDifferences(const GridFunction &x);
1114 #ifdef MFEM_USE_MPI
1116 #endif
1117 
1118  void SetFDhScale(double _dxscale) { dxscale = _dxscale; }
1119  bool GetFDFlag() const { return fdflag; }
1120  double GetFDh() const { return dx; }
1121 
1122  /** @brief Flag to control if exact action of Integration is effected. */
1123  void SetExactActionFlag(bool flag_) { exact_action = flag_; }
1124 };
1125 
1127 {
1128 protected:
1129  // Integrators in the combination. Owned.
1131 
1132 public:
1134 
1136  {
1137  for (int i = 0; i < tmopi.Size(); i++) { delete tmopi[i]; }
1138  }
1139 
1140  /// Adds a new TMOP_Integrator to the combination.
1141  void AddTMOPIntegrator(TMOP_Integrator *ti) { tmopi.Append(ti); }
1142 
1144 
1145  /// Adds the limiting term to the first integrator. Disables it for the rest.
1146  void EnableLimiting(const GridFunction &n0, const GridFunction &dist,
1147  Coefficient &w0, TMOP_LimiterFunction *lfunc = NULL);
1148 
1149  /** @brief Adds the limiting term to the first integrator. Disables it for
1150  the rest (@a dist in the general version of the method) equal to 1. */
1151  void EnableLimiting(const GridFunction &n0, Coefficient &w0,
1152  TMOP_LimiterFunction *lfunc = NULL);
1153 
1154  /// Update the original/reference nodes used for limiting.
1155  void SetLimitingNodes(const GridFunction &n0);
1156 
1157  virtual double GetElementEnergy(const FiniteElement &el,
1159  const Vector &elfun);
1160  virtual void AssembleElementVector(const FiniteElement &el,
1162  const Vector &elfun, Vector &elvect);
1163  virtual void AssembleElementGrad(const FiniteElement &el,
1165  const Vector &elfun, DenseMatrix &elmat);
1166 
1167  /// Normalization factor that considers all integrators in the combination.
1168  void EnableNormalization(const GridFunction &x);
1169 #ifdef MFEM_USE_MPI
1170  void ParEnableNormalization(const ParGridFunction &x);
1171 #endif
1172 };
1173 
1174 /// Interpolates the @a metric's values at the nodes of @a metric_gf.
1175 /** Assumes that @a metric_gf's FiniteElementSpace is initialized. */
1176 void InterpolateTMOP_QualityMetric(TMOP_QualityMetric &metric,
1177  const TargetConstructor &tc,
1178  const Mesh &mesh, GridFunction &metric_gf);
1179 }
1180 
1181 #endif
Abstract class for all finite elements.
Definition: fe.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.hpp:145
InvariantsEvaluator2D< double > ie
Definition: tmop.hpp:306
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:294
void SetFDhScale(double _dxscale)
Definition: tmop.hpp:1118
virtual void SetParDiscreteTargetSize(const ParGridFunction &tspec_)
Definition: tmop.cpp:1049
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:34
Shifted barrier form of metric 56 (area, ideal barrier metric), 2D.
Definition: tmop.hpp:354
Shifted barrier form of 3D metric 16 (volume, ideal barrier metric), 3D.
Definition: tmop.hpp:472
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:275
virtual void SetInitialField(const Vector &init_nodes, const Vector &init_field)=0
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 EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.hpp:157
void SetSerialMetaInfo(const Mesh &m, const FiniteElementCollection &fec, int num_comp)
Definition: tmop.cpp:1849
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:115
Class for an integration rule - an Array of IntegrationPoint.
Definition: intrules.hpp:90
VectorCoefficient * vector_tspec
Definition: tmop.hpp:706
void EnableFiniteDifferences(const GridFunction &x)
Enables FD-based approximation and computes dx.
Definition: tmop.cpp:2747
const TargetType target_type
Definition: tmop.hpp:630
Class for grid function - Vector with associated FE space.
Definition: gridfunc.hpp:30
void UpdateAfterMeshChange(const Vector &new_x)
Definition: tmop.cpp:2723
DenseMatrix PMatO
Definition: tmop.hpp:946
Shape &amp; volume, ideal barrier metric, 3D.
Definition: tmop.hpp:456
const DenseMatrix * Jtr
Definition: tmop.hpp:26
InvariantsEvaluator3D< double > ie
Definition: tmop.hpp:459
const IntegrationRule & Get(int GeomType, int Order)
Returns an integration rule for given GeomType and Order.
Definition: intrules.cpp:915
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:607
Base class for vector Coefficients that optionally depend on time and space.
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.hpp:142
DenseMatrix Jrt
Definition: tmop.hpp:946
const GridFunction * lim_dist
Definition: tmop.hpp:908
virtual void SetSerialDiscreteTargetOrientation(const GridFunction &tspec_)
Definition: tmop.cpp:1156
virtual double GetElementEnergy(const FiniteElement &el, ElementTransformation &T, const Vector &elfun)
Compute the local energy.
Definition: tmop.cpp:2804
double epsilon
Definition: ex25.cpp:136
Coefficient * coeff0
Definition: tmop.hpp:906
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:535
virtual void EvalGrad(DenseMatrix &K, ElementTransformation &T, const IntegrationPoint &ip, int comp)=0
Evaluate the derivative of the matrix coefficient with respect to comp in the element described by T ...
InvariantsEvaluator2D< double > ie
Definition: tmop.hpp:252
void EnableNormalization(const GridFunction &x)
Computes the normalization factors of the metric and limiting integrals using the mesh position given...
Definition: tmop.cpp:2619
void SetSize(int s)
Resize the vector to size s.
Definition: vector.hpp:459
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:2881
Shape, ideal barrier metric, 3D.
Definition: tmop.hpp:406
TMOP_QualityMetric * metric
Definition: tmop.hpp:890
void SetAdaptivityEvaluator(AdaptivityEvaluator *ae)
Definition: tmop.hpp:849
InvariantsEvaluator2D< double > ie
Definition: tmop.hpp:269
virtual ~TMOP_LimiterFunction()
Virtual destructor.
Definition: tmop.hpp:512
Base class for limiting functions to be used in class TMOP_Integrator.
Definition: tmop.hpp:495
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.
virtual void ComputeElementTargetsGradient(const IntegrationRule &ir, const Vector &elfun, IsoparametricTransformation &Tpr, DenseTensor &dJtr) const
Definition: tmop.cpp:1422
Skew metric, 3D.
Definition: tmop.hpp:106
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:734
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:661
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:180
DenseMatrix Jpr
Definition: tmop.hpp:946
IntegrationRules * IntegRules
Definition: tmop.hpp:894
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:517
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:78
int GetOrder() const
Returns the order of the finite element. In the case of anisotropic orders, returns the maximum order...
Definition: fe.hpp:319
virtual void ComputeElementTargetsGradient(const IntegrationRule &ir, const Vector &elfun, IsoparametricTransformation &Tpr, DenseTensor &dJtr) const
Definition: tmop.cpp:984
Container class for integration rules.
Definition: intrules.hpp:309
ParFiniteElementSpace * pfes
Definition: tmop.hpp:560
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:597
Data type dense matrix using column-major storage.
Definition: densemat.hpp:23
Shape, ideal barrier metric, 3D.
Definition: tmop.hpp:390
int Size() const
Returns the size of the vector.
Definition: vector.hpp:160
void SetVolumeScale(double vol_scale)
Used by target type IDEAL_SHAPE_EQUAL_SIZE. The default volume scale is 1.
Definition: tmop.hpp:668
void AssembleElemGradAdaptLim(const FiniteElement &el, const Vector &weights, IsoparametricTransformation &Tpr, const IntegrationRule &ir, DenseMatrix &m)
Definition: tmop.cpp:2374
Shape, ideal barrier metric, 2D.
Definition: tmop.hpp:166
void FinalizeSerialDiscreteTargetSpec()
Definition: tmop.cpp:1165
double metric_normal
Definition: tmop.hpp:900
void EnableAdaptiveLimiting(const GridFunction &z0, Coefficient &coeff, AdaptivityEvaluator &ae)
Restriction of the node positions to certain regions.
Definition: tmop.cpp:1921
Volume metric, 3D.
Definition: tmop.hpp:422
Area, ideal barrier metric, 2D.
Definition: tmop.hpp:303
void SetIntegrationRules(IntegrationRules &irules, int order)
Prescribe a set of integration rules; relevant for mixed meshes.
Definition: tmop.hpp:1032
virtual ~DiscreteAdaptTC()
Definition: tmop.hpp:786
void EnableNormalization(const GridFunction &x)
Normalization factor that considers all integrators in the combination.
Definition: tmop.cpp:2848
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:529
Coefficient * coeff1
Definition: tmop.hpp:898
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.hpp:127
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:470
Abstract parallel finite element space.
Definition: pfespace.hpp:28
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.hpp:112
InvariantsEvaluator2D< double > ie
Definition: tmop.hpp:201
void SetExactActionFlag(bool flag_)
Flag to control if exact action of Integration is effected.
Definition: tmop.hpp:1123
InvariantsEvaluator2D< double > ie
Definition: tmop.hpp:288
void AssembleElementVectorFD(const FiniteElement &el, ElementTransformation &T, const Vector &elfun, Vector &elvect)
Definition: tmop.cpp:2460
const GridFunction * nodes
Definition: tmop.hpp:627
virtual void SetSerialDiscreteTargetSkew(const GridFunction &tspec_)
Definition: tmop.cpp:1139
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:329
void UpdateGradientTargetSpecification(const Vector &x, double dx, bool use_flag=false)
Definition: tmop.cpp:1769
Aspect ratio metric, 2D.
Definition: tmop.hpp:121
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:723
Shape &amp; area, ideal barrier metric, 2D.
Definition: tmop.hpp:182
Default limiter function in TMOP_Integrator.
Definition: tmop.hpp:516
virtual void ComputeElementTargets(int e_id, const FiniteElement &fe, const IntegrationRule &ir, const Vector &elfun, DenseTensor &Jtr) const
Given an element and quadrature rule, computes ref-&gt;target transformation Jacobians for each quadratur...
Definition: tmop.cpp:1237
void AssembleElementGradExact(const FiniteElement &el, ElementTransformation &T, const Vector &elfun, DenseMatrix &elmat)
Definition: tmop.cpp:2229
Array< Vector * > ElemDer
Definition: tmop.hpp:933
const Vector & GetTspecPertMixH()
Definition: tmop.hpp:857
void EnableLimiting(const GridFunction &n0, const GridFunction &dist, Coefficient &w0, TMOP_LimiterFunction *lfunc=NULL)
Adds the limiting term to the first integrator. Disables it for the rest.
Definition: tmop.cpp:2775
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:382
double GetFDh() const
Definition: tmop.hpp:1120
Shifted barrier form of metric 2 (shape, ideal barrier metric), 2D.
Definition: tmop.hpp:214
TMOP_Metric_352(double &t0)
Definition: tmop.hpp:479
virtual void ComputeElementTargets(int e_id, const FiniteElement &fe, const IntegrationRule &ir, const Vector &elfun, DenseTensor &Jtr) const
Given an element and quadrature rule, computes ref-&gt;target transformation Jacobians for each quadratur...
Definition: tmop.cpp:950
GridFunction * zeta
Definition: tmop.hpp:916
AdaptivityEvaluator * adapt_eval
Definition: tmop.hpp:918
Coefficient * coeff_zeta
Definition: tmop.hpp:917
Coefficient * scalar_tspec
Definition: tmop.hpp:705
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:553
virtual void SetSerialDiscreteTargetAspectRatio(const GridFunction &tspec_)
Definition: tmop.cpp:1147
void UpdateTargetSpecification(const Vector &new_x, bool use_flag=false)
Definition: tmop.cpp:1188
void SetTargetJacobian(const DenseMatrix &_Jtr)
Specify the reference-element -&gt; target-element Jacobian matrix for the point of interest.
Definition: tmop.hpp:43
TMOPMatrixCoefficient(int dim)
Definition: tmop.hpp:690
const FiniteElementSpace * tspec_fesv
Definition: tmop.hpp:758
Shape, ideal barrier metric, 2D.
Definition: tmop.hpp:285
virtual ~TMOP_QualityMetric()
Definition: tmop.hpp:35
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:68
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:643
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:398
virtual void ComputeElementTargetsGradient(const IntegrationRule &ir, const Vector &elfun, IsoparametricTransformation &Tpr, DenseTensor &dJtr) const
Definition: tmop.cpp:930
Geometry::Type GetGeomType() const
Returns the Geometry::Type of the reference element.
Definition: fe.hpp:312
Abstract class for local mesh quality metrics in the target-matrix optimization paradigm (TMOP) by P...
Definition: tmop.hpp:23
virtual ~TMOPMatrixCoefficient()
Definition: tmop.hpp:698
virtual ~AdaptivityEvaluator()
Definition: tmop.cpp:1875
InvariantsEvaluator3D< double > ie
Definition: tmop.hpp:425
void SetLimitingNodes(const GridFunction &n0)
Update the original/reference nodes used for limiting.
Definition: tmop.cpp:2796
void SetDiscreteTargetBase(const GridFunction &tspec_)
Definition: tmop.cpp:1090
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:690
DiscreteAdaptTC * discr_tc
Definition: tmop.hpp:920
void EnableLimiting(const GridFunction &n0, const GridFunction &dist, Coefficient &w0, TMOP_LimiterFunction *lfunc=NULL)
Limiting of the mesh displacements (general version).
Definition: tmop.cpp:1896
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:487
void SetCoefficient(Coefficient &w1)
Sets a scaling Coefficient for the quality metric term of the integrator.
Definition: tmop.hpp:1044
Area metric, 2D.
Definition: tmop.hpp:249
void UpdateTargetSpecificationAtNode(const FiniteElement &el, ElementTransformation &T, int nodenum, int idir, const Vector &IntData)
Definition: tmop.cpp:1206
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:340
TargetConstructor(TargetType ttype)
Constructor for use in serial.
Definition: tmop.hpp:645
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:369
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:242
Volume, ideal barrier metric, 3D.
Definition: tmop.hpp:438
DiscreteAdaptTC(TargetType ttype)
Definition: tmop.hpp:777
InvariantsEvaluator2D< double > ie
Definition: tmop.hpp:339
void SetNodes(const GridFunction &n)
Set the nodes to be used in the target-matrix construction.
Definition: tmop.hpp:665
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:631
void AssembleElementVectorExact(const FiniteElement &el, ElementTransformation &T, const Vector &elfun, Vector &elvect)
Definition: tmop.cpp:2091
const TargetConstructor * targetC
Definition: tmop.hpp:891
Array< Vector * > ElemPertEnergy
Definition: tmop.hpp:934
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:496
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:428
InvariantsEvaluator2D< double > ie
Definition: tmop.hpp:218
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:108
virtual void SetParDiscreteTargetAspectRatio(const ParGridFunction &tspec_)
Definition: tmop.cpp:1065
virtual ~TargetConstructor()
Definition: tmop.hpp:658
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:441
void ParEnableNormalization(const ParGridFunction &x)
Definition: tmop.cpp:2864
InvariantsEvaluator2D< double > ie
Definition: tmop.hpp:358
virtual void AssembleElementVector(const FiniteElement &el, ElementTransformation &T, const Vector &elfun, Vector &elvect)
Perform the local action of the NonlinearFormIntegrator.
Definition: tmop.cpp:2062
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:217
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:28
virtual bool ContainsVolumeInfo() const
Checks if the target matrices contain non-trivial size specification.
Definition: tmop.cpp:846
void RestoreTargetSpecificationAtNode(ElementTransformation &T, int nodenum)
Definition: tmop.cpp:1223
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:209
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:164
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.hpp:326
void SetTspecAtIndex(int idx, const GridFunction &tspec_)
Definition: tmop.cpp:1118
void ComputeAvgVolume() const
Definition: tmop.cpp:804
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:712
virtual double Eval(const Vector &x, const Vector &x0, double dist) const
Returns the limiting function, f(x, x0, d).
Definition: tmop.hpp:519
virtual ~TMOP_QuadraticLimiter()
Definition: tmop.hpp:543
Shape, ideal barrier metric, 2D.
Definition: tmop.hpp:233
TMOP_Metric_211(double epsilon=1e-4)
Definition: tmop.hpp:342
void UpdateHessianTargetSpecification(const Vector &x, double dx, bool use_flag=false)
Definition: tmop.cpp:1795
virtual void AssembleElementGrad(const FiniteElement &el, ElementTransformation &T, const Vector &elfun, DenseMatrix &elmat)
Assemble the local gradient matrix.
Definition: tmop.cpp:2076
Untangling metric, 2D.
Definition: tmop.hpp:335
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.hpp:97
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Definition: fespace.hpp:87
const double eps
Definition: tmop.hpp:338
Area, ideal barrier metric, 2D.
Definition: tmop.hpp:266
void SetLimitingNodes(const GridFunction &n0)
Update the original/reference nodes used for limiting.
Definition: tmop.hpp:1085
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:448
Base class Coefficients that optionally depend on space and time. These are used by the BilinearFormI...
Definition: coefficient.hpp:39
InvariantsEvaluator3D< double > ie
Definition: tmop.hpp:476
Collection of finite elements from the same family in multiple dimensions. This class is used to matc...
Definition: fe_coll.hpp:26
void subtract(const Vector &x, const Vector &y, Vector &z)
Definition: vector.cpp:408
TMOPMatrixCoefficient * matrix_tspec
Definition: tmop.hpp:707
virtual void ComputeElementTargets(int e_id, const FiniteElement &fe, const IntegrationRule &ir, const Vector &elfun, DenseTensor &Jtr) const
Given an element and quadrature rule, computes ref-&gt;target transformation Jacobians for each quadratur...
Definition: tmop.cpp:860
const IntegrationRule & ActionIntegrationRule(const FiniteElement &el) const
Definition: tmop.hpp:1004
Base class for Matrix Coefficients that optionally depend on time and space.
InvariantsEvaluator2D< double > ie
Definition: tmop.hpp:236
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:44
void AssembleElemVecAdaptLim(const FiniteElement &el, const Vector &weights, IsoparametricTransformation &Tpr, const IntegrationRule &ir, DenseMatrix &m)
Definition: tmop.cpp:2335
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:192
void DisableLimiting()
Definition: tmop.hpp:990
virtual void AssembleElementGrad(const FiniteElement &el, ElementTransformation &T, const Vector &elfun, DenseMatrix &elmat)
Assemble the local gradient matrix.
Definition: tmop.cpp:2832
double DistanceSquaredTo(const double *p) const
Compute the square of the Euclidean distance to another vector.
Definition: vector.hpp:589
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:321
DiscreteAdaptTC * GetDiscreteAdaptTC() const
Definition: tmop.hpp:1103
TMOP_Metric_252(double &t0)
Note that t0 is stored by reference.
Definition: tmop.hpp:362
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:353
DenseTensor Jtrcomp
Definition: tmop.hpp:753
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:509
const IntegrationRule * IntRule
Definition: nonlininteg.hpp:29
void ResetUpdateFlags()
Used in combination with the Update methods to avoid extra computations.
Definition: tmop.hpp:819
AnalyticAdaptTC(TargetType ttype)
Definition: tmop.hpp:710
void AddTMOPIntegrator(TMOP_Integrator *ti)
Adds a new TMOP_Integrator to the combination.
Definition: tmop.hpp:1141
Vector tspec_pertmix
Definition: tmop.hpp:746
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:412
Skew metric, 2D.
Definition: tmop.hpp:91
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:202
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:259
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:669
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:420
DenseMatrix PMatI
Definition: tmop.hpp:946
InvariantsEvaluator3D< double > ie
Definition: tmop.hpp:441
Shape &amp; area metric, 2D.
Definition: tmop.hpp:198
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:250
void SetTransformation(ElementTransformation &)
The method SetTransformation() is hidden for TMOP_QualityMetrics, because it is not used...
Definition: tmop.hpp:31
void ComputeFDh(const Vector &x, const FiniteElementSpace &fes)
Determines the perturbation, h, for FD-based approximation.
Definition: tmop.cpp:2729
virtual double Eval(const Vector &x, const Vector &x0, double d) const =0
Returns the limiting function, f(x, x0, d).
Class for integration point with weight.
Definition: intrules.hpp:25
AdaptivityEvaluator * adapt_eval
Definition: tmop.hpp:766
InvariantsEvaluator3D< double > ie
Definition: tmop.hpp:377
Array< TMOP_Integrator * > tmopi
Definition: tmop.hpp:1130
A standard isoparametric element transformation.
Definition: eltrans.hpp:348
void ComputeMinJac(const Vector &x, const FiniteElementSpace &fes)
Definition: tmop.cpp:2690
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:615
void ComputeNormalizationEnergies(const GridFunction &x, double &metric_energy, double &lim_energy)
Definition: tmop.cpp:2638
void Diag(double c, int n)
Creates n x n diagonal matrix with diagonal elements c.
Definition: densemat.cpp:1348
This class is used to express the local action of a general nonlinear finite element operator...
Definition: nonlininteg.hpp:26
const IntegrationRule & EnergyIntegrationRule(const FiniteElement &el) const
Definition: tmop.hpp:995
const IntegrationRule & GradientIntegrationRule(const FiniteElement &el) const
Definition: tmop.hpp:1009
double GetFDDerivative(const FiniteElement &el, ElementTransformation &T, Vector &elfun, const int nodenum, const int idir, const double baseenergy, bool update_stored)
Definition: tmop.cpp:2438
void SetParMetaInfo(const ParMesh &m, const FiniteElementCollection &fec, int num_comp)
Parallel version of SetSerialMetaInfo.
Definition: tmop.cpp:1862
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:390
int dim
Definition: ex24.cpp:53
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:1955
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:22
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:160
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.
const GridFunction * zeta_0
Definition: tmop.hpp:915
Aspect ratio metric, 3D.
Definition: tmop.hpp:136
const Vector & GetTspecPert2H()
Definition: tmop.hpp:856
const Vector & GetTspecPert1H()
Definition: tmop.hpp:855
void AssembleElementGradFD(const FiniteElement &el, ElementTransformation &T, const Vector &elfun, DenseMatrix &elmat)
Definition: tmop.cpp:2524
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:682
InvariantsEvaluator2D< double > ie
Definition: tmop.hpp:169
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:562
virtual void SetParDiscreteTargetSpec(const ParGridFunction &tspec_)
Definition: tmop.cpp:1083
InvariantsEvaluator2D< double > ie
Definition: tmop.hpp:185
Abstract class for hyperelastic models.
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:128
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:284
TargetConstructor(TargetType ttype, MPI_Comm mpicomm)
Constructor for use in parallel.
Definition: tmop.hpp:654
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
TMOP_Integrator(TMOP_QualityMetric *m, TargetConstructor *tc)
Definition: tmop.hpp:1018
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:547
void FinalizeParDiscreteTargetSpec(const ParGridFunction &tspec_)
Definition: tmop.cpp:1018
bool GetFDFlag() const
Definition: tmop.hpp:1119
Vector data type.
Definition: vector.hpp:51
DenseMatrix DSh
Definition: tmop.hpp:946
DenseMatrix DS
Definition: tmop.hpp:946
TMOP_Metric_022(double &t0)
Definition: tmop.hpp:221
Shape+Size+Orientation metric, 2D.
Definition: tmop.hpp:151
InvariantsEvaluator3D< double > ie
Definition: tmop.hpp:409
const GridFunction * nodes0
Definition: tmop.hpp:905
const FiniteElementSpace * tspec_fes
Definition: tmop.hpp:757
TargetType
Target-matrix construction algorithms supported by this class.
Definition: tmop.hpp:607
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:361
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:501
virtual void SetSerialDiscreteTargetSpec(const GridFunction &tspec_)
Definition: tmop.cpp:1181
virtual void SetParDiscreteTargetSkew(const ParGridFunction &tspec_)
Definition: tmop.cpp:1057
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:332
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:526
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:130
virtual void SetParDiscreteTargetOrientation(const ParGridFunction &tspec_)
Definition: tmop.cpp:1074
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:603
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:186
virtual void SetAnalyticTargetSpec(Coefficient *sspec, VectorCoefficient *vspec, TMOPMatrixCoefficient *mspec)
Definition: tmop.cpp:941
Rank 3 tensor (array of matrices)
Definition: densemat.hpp:728
Class for parallel meshes.
Definition: pmesh.hpp:32
virtual void SetSerialDiscreteTargetSize(const GridFunction &tspec_)
Definition: tmop.cpp:1130
IntegrationRules IntRules(0, Quadrature1D::GaussLegendre)
A global object with all integration rules (defined in intrules.cpp)
Definition: intrules.hpp:378
TMOP_LimiterFunction * lim_func
Definition: tmop.hpp:910
const Array< TMOP_Integrator * > & GetTMOPIntegrators() const
Definition: tmop.hpp:1143
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:653
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:637
InvariantsEvaluator3D< double > ie
Definition: tmop.hpp:393
virtual void AssembleElementVector(const FiniteElement &el, ElementTransformation &T, const Vector &elfun, Vector &elvect)
Perform the local action of the NonlinearFormIntegrator.
Definition: tmop.cpp:2816
DenseMatrix Jpt
Definition: tmop.hpp:946
Shape &amp; orientation metric, 2D.
Definition: tmop.hpp:320
void ParEnableNormalization(const ParGridFunction &x)
Definition: tmop.cpp:2627
virtual void EvalP(const DenseMatrix &Jpt, DenseMatrix &P) const
Evaluate the 1st Piola-Kirchhoff stress tensor, P = P(Jpt).
Definition: tmop.cpp:768
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:100
DenseMatrix P
Definition: tmop.hpp:946
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:780
virtual void ComputeAtNewPosition(const Vector &new_nodes, Vector &new_field)=0
bool Parallel() const
Definition: tmop.hpp:634
Shape, ideal barrier metric, 3D.
Definition: tmop.hpp:374
Metric without a type, 2D.
Definition: tmop.hpp:75
virtual double EvalW(const DenseMatrix &Jpt) const
Evaluate the strain energy density function, W = W(Jpt).
Definition: tmop.cpp:760
FiniteElementSpace * fes
Definition: tmop.hpp:555
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:457
A TMOP integrator class based on any given TMOP_QualityMetric and TargetConstructor.
Definition: tmop.hpp:884