R Adaptive Coordinate Transformations

• R-adaptive coordinate transformations are techniques that reposition computational mesh points to accurately capture sharp gradients and discontinuities without altering mesh connectivity.
• They utilize rigorous frameworks such as optimal transport, Monge–Ampère formulations, and variational principles to ensure controlled mesh alignment, density, and anisotropy.
• These methods are pivotal in adaptive PDE solvers, high-dimensional tensor PDE integrators, and neural operator learning frameworks like DeepONet for enhanced simulation accuracy.

R-adaptive coordinate transformations are a fundamental class of adaptive computational techniques in numerical analysis, scientific computing, and machine learning for capturing complex solution features—such as sharp gradients, discontinuities, or evolving fronts—by continuously redistributing computational points without altering mesh connectivity. The “r-adaptive” paradigm, differing from h-adaptivity (refinement/coarsening) and p-adaptivity (polynomial degree), is characterized by mesh point relocation, typically driven by monitor/density functions or metric tensors, and realized through coordinate transformations. They play a central role in mesh generation, adaptive PDE solvers, high-dimensional tensor PDE methods, neural operator learning (e.g., DeepONets), and sparse-grid regression.

1. Monge–Ampère and Optimal Transport Formulation

The optimal transport (OT) approach provides a rigorous geometric framework for r-adaptivity by reformulating mesh redistribution as a Monge–Ampère boundary value problem. Given computational (Ω_c) and physical (Ω_p) domains of equal volume, the goal is to construct an invertible mapping $x = x(\xi): \Omega_c \to \Omega_p$ that equidistributes a prescribed positive density $\rho(x)$ by enforcing

$$\rho(x) \, J(\xi) = \theta, \qquad J = \det D_\xi x,$$

where $$\theta = \frac{\int_{\Omega_p} \rho(x)\,dx}{\int_{\Omega_c} d\xi}$$ is the normalization constant, and the mapping $x(\xi)$ is Brenier-optimal, i.e., minimizes $\int_{\Omega_c} \|x(\xi)-\xi\|^2 d\xi$. The map is realized as the gradient of a convex potential $P(\xi)$, i.e., $x(\xi) = \nabla_\xi P(\xi)$, yielding a symmetric Jacobian $J = D_\xi x = D^2 P$. The equidistribution condition becomes a fully nonlinear Monge–Ampère PDE in $P$,

$\det(D^2P(\xi))\,\rho(\nabla P(\xi)) = \theta,$

subject to appropriate boundary conditions (Neumann, periodic, or mixed). The method generalizes naturally to three dimensions and remains scalar-valued at each time step, thus avoiding tensor-valued PDEs (Budd et al., 2014).

2. The Role of Metric Tensors and Implicit Mesh Metrics

Although the OT-based method does not explicitly reference a metric tensor, the induced mesh metric $M(x) = \theta J^{-2}$ is central for quantifying mesh alignment, size, and anisotropy. The symmetric positive-definite $J$ admits an orthonormal eigen-decomposition $J = \sum_{i=1}^2 \lambda_i v_i v_i^T$, leading to

$$M(x) = \sum_{i=1}^2 \mu_i v_i v_i^T,\quad \mu_i = \theta/\lambda_i^2.$$

The principal directions $v_i$ determine the axes of mesh stretching, and the eigenvalues $\mu_i$ quantify mesh tightness in those directions. For model densities, two canonical cases are analytically tractable:

• Linear features: With $\rho(x) = \rho_1(x \cdot e_1)\rho_2(x \cdot e_2)$, the metric eigenstructure aligns with $e_1, e_2$, and the mesh is strongly stretched tangential and compressed normal to features where $\rho_1 \gg 1$.
• Radially symmetric features: For $\rho = \rho(R)$, the eigenvectors are radial and tangential, with eigenvalues $dR/dr$ and $R/r$; the anisotropy ratio $(R/r)/(dR/dr)$ sharply increases near transition regions (e.g., rings).

The ratio $A = \lambda_2/\lambda_1$ quantifies local anisotropy, enabling predictive control of mesh stretching and element elongation.

3. Variational and Multisymplectic R-adaptivity

Variational formulations provide robust structure-preserving r-adaptive methods, especially for time-dependent Lagrangian field theories. Tyranowski and Desbrun present two canonical strategies in 1+1D: (i) post-discretization control-theoretic coupling, where mesh equations are treated as algebraic constraints after discretizing the physical field, and (ii) a pseudo-particle (vakonomic) approach, where the mesh points $X_i(t)$ are full dynamical variables and enforced to follow adaptation constraints $g_i(y, X)=0$ via Lagrange multipliers in the variational principle. The resulting discretizations yield high-order index-1 or index-3 DAEs, which, when solved via constrained symplectic Runge-Kutta methods (e.g., 4th-order Lobatto IIIA–IIIB), rigorously conserve discrete multisymplectic structure and ensure long-time energy stability (Tyranowski et al., 2013).

A typical monitor function for mesh density is

$\rho = \sqrt{1 + \alpha^2 \phi_X^2},$

favoring refinement near steep gradients. For continuous r-adaptation, one solves the weighted equidistribution PDE

$$\partial_x(\rho(X)\, X_x) = 0 \implies \rho(X)\, X_x = \tfrac{1}{\sigma},$$

with $\sigma$ the integral normalization. In unsteady problems, relaxation via a moving mesh PDE,

$X_t = \frac{1}{\tau} \partial_x(\rho(X) X_x),$

smoothly tracks evolving features.

4. Variational R-adaptation for High-Order Meshes

A variational mesh optimization approach for high-order r-adaptation seeks to minimize a global deformation energy

$$E[\varphi] = \sum_e \int_{\Omega_e^c} W(J_e(x)) d\Omega,$$

where $J_e$ is the local Jacobian of the mapping, and $W$ is an element-wise strain energy (e.g., hyperelastic). Elementwise target metric tensors $M_e$ specify desired size, shape, and anisotropy. Each element's reference-to-target affine map $T_e$ satisfies $T_e^T T_e = M_e$; the full mapping $\varphi_M = \varphi \circ T_e$ ensures element alignment with prescribed metric shape. Optimization proceeds via Newton–Raphson with line search, using per-element residuals and tangents defined by the first Piola–Kirchhoff stress $P = \partial W/\partial J$ and its derivatives. This framework naturally accommodates high-order, curved elements through polynomial Lagrange or modal shape functions, allowing exact CAD geometry projection and r-adaptation without remeshing (Marcon et al., 2019).

5. Coordinate Transformations in Machine Learning and High-dimensional PDEs

R-adaptive transformations are integrated into high-dimensional approximation and machine learning frameworks:

a) Coordinate-adaptive tensor PDE integrators

To control tensor rank in high-dimensional PDEs (e.g., for functional tensor train or hierarchical Tucker formats), one introduces a time-dependent diffeomorphic map $\Phi_t: x \mapsto y$ to deform computational coordinates. The evolution of $\Phi_t$ is chosen to minimize the normal component of the PDE operator relative to the tensor manifold (i.e., ensure the projected dynamics are as tangential as possible), leading to a strictly convex variational principle for $\dot{\Phi}_t$. For linear flows ($\Phi_t(x) = \Gamma(t)x$), the optimal flow generator $\dot{\Gamma}(t)$ is obtained by solving a global linear system, ensuring globally unique solutions. The method markedly reduces rank growth, as demonstrated on Liouville and Fokker–Planck equations, achieving lower tensor ranks for a fixed error compared to non-convex Riemannian optimization (Dektor et al., 2023).

b) R-adaptive DeepONet

For neural operator learning with discontinuous solution operators, DeepONet's linear reconstruction structure imposes a fundamental approximation barrier. The R-adaptive DeepONet overcomes this by learning coordinate transformations ($a \mapsto x(\xi)$) via the equidistribution principle, then composing with a second DeepONet in computational coordinates to yield smooth target representations. The overall solution operator is

$$u(x) = [\,\tilde{\mathcal{G}}(a) \circ \mathcal{T}(a)^{-1}\,](x),$$

where both $\mathcal{T}$ and $\tilde{\mathcal{G}}$ are learned via DeepONet parameterizations. Weighted loss functions reflect Jacobian scaling and solution gradients, penalizing mesh tangling. Theoretical results show superlinear convergence (e.g., $O(p^{-3/2})$ for 1D advection and Burgers) surpassing the $O(p^{-1})$ rate of vanilla DeepONet, as confirmed for linear advection, low-viscosity Burgers, and gas-dynamics equations (Zhu et al., 8 Aug 2024).

c) ANOVA-effective dimension reduction via optimally rotated coordinates

In adaptive sparse-grid regression, the curse of dimensionality is mitigated by finding an orthogonal coordinate transformation $Q$ that minimizes the mean ANOVA dimension of the problem; the optimal $Q$ is found by optimizing the variance captured by low-order ANOVA terms in the rotated coordinates. A polynomial surrogate is fit to the data, closed-form partial variances are computed, and a conjugate gradient on the Stiefel manifold solves the optimization. This rotation can reduce the number of sparse-grid points by factors of 2–10 for fixed error targets (Bohn et al., 2018).

6. Practical Implementation and Algorithmic Aspects

In the OT/Monge–Ampère paradigm, the scalar PDE is solved by finite difference, Newton–Raphson, or parabolic relaxation (PMA) methods; Neumann or periodic boundary conditions maintain the invertibility of the mapping. The r-adaptive mapping is computed as the solution to (MA), then used to move mesh nodes, with mesh metric tensors, eigenstructures, and anisotropy ratios guiding further adaptation or error control. For variational mesh optimization, elementwise optimization is performed using assembled sparse global systems, with high-order node movement constrained (for CAD models, nodes slide on the surface). In variational–multisymplectic integrators, post-discretization control or Lagrange-multiplier coupling is enforced via (constrained) symplectic integrators, preserving geometric structure.

In high-dimensional settings, coordinate adaptation is performed stepwise: at each time, the optimal flow is computed as the unique minimizer of quadratic cost functional, advanced with the tensor PDE solver, and the coordinate map updated (e.g., via a matrix ODE for linear flows). Non-convex Riemannian optimization or static linear transformations can provide initial low-dimensional subspaces, but convex continuous adaptation yields globally optimal, feedback-driven flows with better rank control and stability.

For machine learning applications, preprocessing (e.g., in R) involves surrogate fitting (least squares polynomial), computation of ANOVA variances, objective/gradient evaluation, manifold optimization (e.g., with manifoldOptim), and transformation of inputs prior to adaptive regression.

7. Numerical Behavior, Anisotropy, and Applications

Mesh alignment, stretching, and skewness in OT-generated r-adaptive meshes are analytically characterized: strong anisotropy (large $A$) occurs next to sharp features (lines, rings), mesh elements align tangentially, and transition smoothly away from features. Numerical experiments verify close agreement between analytical predictions and computed anisotropy ratios, smooth transitions between linear and curved arcs, and modulation of mesh quality near boundaries. In variational high-order methods, mesh adaptation achieves strong local anisotropy (e.g., radial shrinking), and deformed high-order nodes maintain CAD geometrical validity. For multisymplectic field theories, structure-preserving r-adaptivity yields stable energy evolution on challenging PDEs such as Sine–Gordon, with conservation up to soliton collisions.

Coordinate-adaptive tensor PDE solvers halve tensor rank requirements compared to fixed-coordinates methods for the same error threshold. R-adaptive DeepONet achieves superlinear error decay for discontinuous solutions, exceeding the rates attainable by meshless or uniform network architectures. In high-dimensional sparse-grid regression, optimal coordinate rotations substantially accelerate accuracy improvements.

A summary table of principal r-adaptive coordinate transformation strategies appears below:

Method/Class	Primary Equation/Principle	Key Application Domain
Monge–Ampère/OT mesh adaptivity	$\det D^2P\cdot\rho(\nabla P) = \theta$	Unstructured mesh PDE solvers (Budd et al., 2014)
Variational mesh optimization	Minimize $\sum_e \int W(J_e(x))d\Omega$	High-order/curved meshes (Marcon et al., 2019)
Multisymplectic variational methods	Variational action with mesh constraints	Structure-preserving time-dependent PDEs (Tyranowski et al., 2013)
Tensor PDE coordinate adaptivity	Minimize normal component of PDE operator	High-dimensional PDEs (FTT/HT) (Dektor et al., 2023)
Machine Learning/DeepONet R-adapt.	Learn a$\mapsto$x($\xi$) via equidistribution	Neural operator learning (Zhu et al., 8 Aug 2024)
ANOVA-optimal rotation	Minimize mean effective ANOVA dimension	Sparse-grid regression (Bohn et al., 2018)

The diversity and rigor of r-adaptive coordinate transformations underline their indispensability across computational mathematics, geometric integration, tensor-based algorithms, and scientific machine learning. Their ability to encode mesh alignment, maintain global smoothness, and provide predictive control of anisotropy renders them essential for accurate, efficient, and structure-preserving simulation of PDEs and data-driven models.