A new Sobolev gradient method for direct minimization of the Gross-Pitaevskii energy with rotation

• Rotating condensates form quantized vortices, with increasing rotation producing arrays resembling Abrikosov lattices.Source: §1 · p. 1
• The algorithm-design problem has two connected parts: iteratively decreasing energy toward a local minimum and accounting for the mass constraint.Source: §1 · p. 2
• Prior approaches commonly apply unconstrained steepest descent and then normalize the result to satisfy the unit-norm constraint.The paper's projection method offers an alternative to this normalization approach.Source: §1 · pp. 2–3
• The new Sobolev gradient method outperforms classical L2 and H1 gradient methods, especially at high rotation rates.It is implemented with both finite differences and finite elements in two dimensions.Source: §1 · p. 3
• The paper proves existence and convergence in the asymptotic limit for trajectories of the Sobolev-gradient evolution equations.The introduction describes numerical tests computing physically relevant stationary states at high rotation and large interaction constants.Source: §1 · p. 3
• A projected Sobolev gradient preserves the initial state's norm during minimization, avoiding normalization and penalty-term approaches to the mass constraint.The authors derive an explicit projected gradient for Hilbert spaces beyond L2.Source: §1 · p. 3

2. Gradient descent methods using several gradients.

• The energy's reformulation motivates an inner product whose norm makes the functional coercive, a property suited to rotating cases.The descent direction is chosen to maximize the energy change at each step.Source: §2 · p. 4

2.1. Inner products and norms.

• The magnetic-gradient operator is ∇A = ∇ + iΩA^t; the resulting HA norm is equivalent to the H1 norm, and H1 and HA are equal as sets.The paper uses the equivalent HA norm to improve descent performance for the rotating GP energy.Source: §2.1 · pp. 4–5

2.2. Gradients.

• The energy derivative has a unique Hilbert-space representative for each inner product, defining the corresponding gradient.For finite differences the H1 and HA gradients use strong formulations; finite elements use weak formulations.Source: §2.2 · p. 6

3.1. Projection method for the mass constraint.

• Projecting the Sobolev gradient onto the constraint tangent space yields a descent trajectory that preserves the initial mass for all time.The paper derives an explicit projected-gradient expression to facilitate numerical implementation, including for Hilbert spaces other than L2.Source: §3.1 · pp. 6–7

3.2. Convergence results in an infinite dimensional Hilbert space.

• The convergence analysis extends earlier results for the non-rotating Gross-Pitaevskii energy to the rotating, constrained minimization problem.The authors note that general Sobolev-gradient theory does not provide the needed asymptotic-convergence criteria for constrained evolution.Source: §3.2 · p. 8
• Uniform and strict convexity in HA, and hence in H1 by norm equivalence, supports the minimizer and convergence result.The convexity argument uses Veff bounded away from zero; adding a multiple of the constraint to the energy can obtain this property without changing the minimizers.Source: §3.2 · p. 8
• Under the stated hypotheses, the evolution is uniquely defined for all time and converges to a constrained minimizer, with exponential convergence in the stated norm.The limit satisfies the constrained stationarity condition; the convexity argument assumes Veff is bounded away from zero, or that a multiple of the constraint is added to the energy.Source: §3.2 · pp. 8–9

4. Numerical implementation.

• The discrete algorithm advances a projected gradient trajectory by forward Euler in the finite-dimensional Hilbert space and can choose each time step by local energy minimization.The algorithm stops when the relative change in discrete energy falls below an imposed limit.Source: §4 · p. 9
• The projection method permits larger time steps than normalization alone, although the first-order time discretization preserves the norm only up to an error of order (δtn)^2∥∇zn,Y EG(zn)∥L2.Normalizing after iterations can improve constraint accuracy; convergence is based on the relative energy change falling below a threshold.Source: §4 · p. 9

4.1. Finite differences.

• The finite-difference implementation uses an equally spaced N by N grid and fourth-order centered approximations for first derivatives, which proved helpful for resolving vortex lattices on reasonably fine grids.A linesearch locally selects the time step, reducing the iterations needed for convergence; Sobolev gradients are computed by solving linear systems with preconditioned conjugate gradients.Source: §4.1 · pp. 10–11
• In the finite-difference setting, discrete L2, H, and HA inner products equip the grid-function space with corresponding finite-dimensional Hilbert-space structures.The discretized derivatives and their conjugate-transpose matrices provide the operators used to define the H and HA inner products.Source: §4.1 · p. 10

4.2. Finite elements.

• The finite-element implementation uses FreeFem++ and integral variational formulations to compute the ordinary, H1, and HA gradients.The HA gradient is solved directly from the weak formulation, avoiding an explicit L2-gradient solve; the implementation uses complex variables and P1 elements, with adaptive refinement for dense vortex lattices.Source: §4.2 · p. 12

5. Numerical experiments.

• The experiments progress from a manufactured-solution benchmark to metastable vortex states and then complex configurations, comparing the three gradients before using HA for the latter cases.The benchmark assesses convergence and computational cost in finite-difference and finite-element implementations.Source: §5 · p. 13, §5.1 · pp. 13–14

5.1. Test case with manufactured solutions.

• Approximately a factor-of-3 computational-time reduction was reported for HA over H1 at large Ω; L2 required many more time steps and much larger CPU times.The authors attribute L2's slow convergence to stability-limited small time steps, while Sobolev gradients allow larger steps.Source: §5.1 · pp. 15–16

5.2. Simulations of rotating Bose-Einstein condensates.

• The condensate simulations used Thomas-Fermi estimates to size the computational domain, with the rotated estimate based on the effective trapping potential.For strong interactions, the initial-state motivation is to neglect kinetic energy in the simplified energy functional.Source: §5.2 · p. 16

5.2.1. Off-center vortex case: harmonic trapping potential and small

• The off-center initial vortex converged to a single centered vortex; H1 and HA performed similarly at low Ω.A convergence criterion of ε = 10^-8 was needed to obtain a centered vortex and Lz = 1; a relaxed criterion left it slightly off-center.Source: §5.2.1 · p. 17, Figure 5.2 · p. 17

5.2.2. Vortex array case: harmonic-plus-quartic trapping potential and large Ω.

• Up to 36% CPU time gain was reported for HA in the high-rotation vortex-array case, and both discretizations produced the same converged state.HA also had better stability properties; the final state contains a central vortex and an array of seven vortices.Source: §5.2.2 · pp. 18–19, Figure 5.3 · p. 18

5.2.3. Giant vortex and Abrikosov vortex lattice.

• With HA and the projection method, increasing Ω from 2 to 4 produced a transition from a vortex lattice to a giant vortex in both discretizations.The giant vortex is a hole in the condensate with multiple phase defects; the simulations also formed an Abrikosov lattice with several hundred vortices, for which adaptive mesh refinement helped reduce time and capture positions.Source: §5.2.3 · pp. 19–20, Figure 5.4 · p. 20, Figure 5.5 · p. 20

6. Summary.

• The HA inner product yields a gradient that enhances descent-method performance over L2 or H1 gradients, with particularly important time gains at high rotation.The paper also introduces a projection method for the mass constraint and reports robust finite-difference and finite-element implementations for complex two-dimensional configurations.Source: §6 · pp. 20–21
• The projection method avoids penalty-term energy functionals and traditional normalization during descent.The authors suggest that the gradient and projection tools could also be combined with existing schemes, including sophisticated time stepping, for GP-type energy minimization.Source: §6 · p. 21