|
restrict_anisotropy=True, |
|
): |
|
""" |
|
Apply :math:`L^p` normalisation to the metric. |
|
|
|
See |
|
https://mesh-adaptation.github.io/docs/animate/1-metric-based.html#metric-normalisation |
|
for more details. |
|
|
|
:kwarg global_factor: pre-computed global normalisation factor |
|
:type global_factor: :class:`float` |
|
:kwarg boundary: is the normalisation to be done over the boundary? |
|
:type boundary: :class:`bool` |
|
:kwarg restrict_sizes: should minimum and maximum metric magnitudes be enforced? |
|
:type restrict_sizes: :class:`bool` |
|
:kwarg restrict_anisotropy: should maximum anisotropy be enforced? |
|
:type restrict_anisotropy: :class:`bool` |
|
:return: the normalised metric, modified in-place |
|
:rtype: :class:`~.RiemannianMetric` |
|
""" |
|
d = self._tdim - 1 if boundary else self._tdim |
|
p = self.metric_parameters.get("dm_plex_metric_p", 1.0) |
|
target = self.metric_parameters.get("dm_plex_metric_target_complexity") |
|
if target is None: |
|
raise ValueError("dm_plex_metric_target_complexity must be set.") |
|
|
|
# Enforce that the metric is SPD |
|
self.enforce_spd( |
|
restrict_sizes=False, restrict_anisotropy=self._restrict_anisotropy_first |
|
) |
|
|
|
# Compute global normalisation factor |
|
detM = ufl.det(self) |
|
if global_factor is None: |
|
dX = (ufl.ds if boundary else ufl.dx)(domain=self._mesh) |
|
exponent = 0.5 if np.isinf(p) else (p / (2 * p + d)) |
|
integral = firedrake.assemble(pow(detM, exponent) * dX) |
|
global_factor = firedrake.Constant(pow(target / integral, 2 / d)) |
|
|
|
# Normalise the metric |
|
if boundary: |
|
raise NotImplementedError( |
|
"Normalisation on the boundary not yet implemented." |
|
) |
|
determinant = 1 if np.isinf(p) else pow(detM, -1 / (2 * p + d)) |
|
self.interpolate(global_factor * determinant * self) |
|
|
|
# Enforce element constraints |
|
return self.enforce_spd( |
|
restrict_sizes=restrict_sizes, |
|
restrict_anisotropy=not self._restrict_anisotropy_first, |
|
) |
animate/animate/metric.py
Lines 536 to 587 in cc40cd8
The
restrict_anisotropyargument in the :function:enforce_spdmethod of the :class:RiemannianMetricconflicts with theself._restrict_anisotropy_firstvariable in the :function:normalisemethod. It is unclear whetherrestrict_anisotropyis still necessary or has been replaced byself._restrict_anisotropy_first.