On discrete ground states of rotating Bose-Einstein condensates

TL;DR

This paper analyzes the approximation of discrete ground states of rotating Bose-Einstein condensates using finite element methods, providing optimal error estimates and numerical validation despite the challenge of non-uniqueness.

Contribution

It develops an a priori error analysis for finite element approximations of ground states, addressing non-uniqueness via tangent space restrictions.

Findings

Optimal error estimates in L^2 and H^1 norms

Error bounds for ground state energy and chemical potential

Numerical experiments confirming theoretical results

Abstract

The ground states of Bose-Einstein condensates in a rotating frame can be described as constrained minimizers of the Gross-Pitaevskii energy functional with an angular momentum term. In this paper we consider the corresponding discrete minimization problem in Lagrange finite element spaces of arbitrary polynomial order and we investigate the approximation properties of discrete ground states. In particular, we prove a priori error estimates of optimal order in the $L^2$- and $H^1$-norm, as well as for the ground state energy and the corresponding chemical potential. A central issue in the analysis of the problem is the missing uniqueness of ground states, which is mainly caused by the invariance of the energy functional under complex phase shifts. Our error analysis is therefore based on an Euler-Lagrange functional that we restrict to certain tangent spaces in which we have local uniqueness of ground states. This gives rise to an error decomposition that is ultimately used to derive the desired a priori error estimates. We also present numerical experiments to illustrate various aspects of the problem structure.

Equations

$$E(w)=\frac{1}{2}{\int }_{\mathcal{D}}|\nabla w{|}^2+V|w{|}^2-\Omega \bar{w}{\mathcal{L}}_{3}w+\frac{\beta }{2}|w{|}^4\text{mbox}dx$$

$$(v,w)_{L^{2}(\mathcal{D})}:=\mbox{\rm Re}\Big{(}\int_{\mathcal{D}}v\,\overline{w}\hskip 2.0pt\mbox{d}x\Big{)},$$

$$\hskip 17.07164pt(v,w)_{H^{1}(\mathcal{D})}:=\mbox{\rm Re}\Big{(}\int_{\mathcal{D}}v\,\overline{w}\hskip 2.0pt\mbox{d}x+\int_{\mathcal{D}}\nabla v\cdot\overline{\nabla w}\hskip 2.0pt\mbox{d}x\Big{)}.$$

$$V(x)-\frac{1+\epsilon }{4}{\Omega }^2(x_1^2+x_2^2)\ge 0\text{for almost all }x\in \mathcal{D}.$$

$$E:H_0^1(\mathcal{D})\to \mathbb{R}$$

$$E(w):=\frac{1}{2}{\int }_{\mathcal{D}}|\nabla w{|}^2+V|w{|}^2-\Omega \bar{w}{\mathcal{L}}_{3}w+\frac{\beta }{2}|w{|}^4\text{mbox}dx.$$

$${\nabla }_{\text{mbox}\mathrm{R}}w:=\nabla w+\mathrm{i}\frac{\Omega }{2}{R}^{\top }w$$

$$\displaystyle(v,w)_{\mbox{\rm\tiny R}}:=\mbox{\rm Re}\Big{(}\int_{\mathcal{D}}\nabla_{\mbox{\rm\tiny R}}v\cdot\overline{\nabla_{\mbox{\rm\tiny R}}w}\hskip 2.0pt\mbox{d}x+\int_{\mathcal{D}}V_{\mbox{\rm\tiny R}}v\,\overline{w}\hskip 2.0pt\mbox{d}x\Big{)}\,$$

$$E(w)=\frac{1}{2}\Vert w{\Vert }_{\text{mbox}\mathrm{R}}^2+\frac{\beta }{4}{\int }_{\mathcal{D}}|w{|}^4\text{mbox}dx.$$

$$E(u)=\underset{v\in \mathbb{S}}{\inf }E(v)\text{mbox}where\mathbb{S}:=\{v\in H_0^1(\mathcal{D})\text{mbox}and\Vert v{\Vert }_{L^2(\mathcal{D})}=1\}.$$

$$E(v)\ge \frac{1}{2}\Vert v{\Vert }_{\text{mbox}\mathrm{R}}^2>0\text{mbox}forallv\in \mathbb{S}.$$

$$c(\epsilon )\Vert v{\Vert }_{H^1(\mathcal{D})}\le \Vert v{\Vert }_{\text{mbox}\mathrm{R}}\le C(V,\Omega )\Vert v{\Vert }_{H^1(\mathcal{D})}\text{mbox}forallv\in H_0^1(\mathcal{D}).$$

$$\mathcal{U}=\{u\in \mathbb{S}|E(u)=\underset{v\in \mathbb{S}}{\inf }E(v)\}$$

$$\hat{u}\in {\mathbb{B}}_{\delta }(u)\cap \mathcal{U}\Rightarrow \hat{u}=u{e}^{\mathrm{i}\theta }\text{mbox}forsome\theta \in [0,2\pi ).$$

$$T_u\mathbb{S}:=\{v\in H_0^1(\mathcal{D})|(u,v{)}_{L^2(\mathcal{D})}=0\}.$$

$$\frac{\text{mbox}{d}^2}{\text{mbox}d{t}^2}E(\gamma (t)){|}_{t=0}>0.$$

$${\mathcal{H}}_u^{\mathbb{C},\perp }:=T_u\mathbb{S}\cap T_{\mathrm{i}u}\mathbb{S}.$$

$${\mathbf{R}}_{\theta }=\left(\begin{array}{cc}\cos \theta & -\sin \theta \\ \sin \theta & \cos \theta \end{array}\right),$$

$$⟨E^{\prime }(u),v⟩=\lambda ⟨M^{\prime }(u),v⟩\text{mbox}forallv\in H_0^1(\mathcal{D}),$$

$$⟨M^{\prime }(u),v⟩=\text{mbox}\mathrm{Re}{\int }_{\mathcal{D}}u\overline{v}\text{mbox}dx=(u,v{)}_{L^2(\mathcal{D})}$$

$$⟨E^{\prime }(u),v⟩$$

$$\displaystyle\langle\mathcal{A}_{|u|}v,w\rangle:=(\nabla_{\mbox{\rm\tiny R}}v,\nabla_{\mbox{\rm\tiny R}}w)_{L^{2}(\mathcal{D})}+\big{(}(V_{\mbox{\rm\tiny R}}+\beta|u|^{2})v,w\big{)}_{L^{2}(\mathcal{D})}\qquad\mbox{for }v,w\in H^{1}_{0}(\mathcal{D}).$$

$$⟨{\mathcal{A}}_{|u|}v,v⟩\ge \alpha \Vert v{\Vert }_{H^1(\mathcal{D})}^2\text{mbox}forallv\in H_0^1(\mathcal{D}),$$

$$⟨{\mathcal{A}}_{|u|}u,v⟩=\lambda (u,v{)}_{L^2(\mathcal{D})}\text{mbox}forallv\in H_0^1(\mathcal{D}).$$

$$-\Delta u+Vu-\Omega {\mathcal{L}}_{3}u+\beta |u{|}^{2}u=\lambda u.$$

$$\underset{w\in {\mathcal{H}}_u^{\mathbb{C},\perp }}{\inf }\underset{v\in {\mathcal{H}}_u^{\mathbb{C},\perp }}{\sup }\frac{|⟨(E^{\prime \prime }(u)-\lambda \mathcal{I})v,w⟩|}{\Vert v{\Vert }_{H^1(\mathcal{D})}\Vert w{\Vert }_{H^1(\mathcal{D})}}⩾{\alpha }_1.$$

$$⟨E^{\prime \prime }(u)v,w⟩=⟨{\mathcal{A}}_{|u|}v,w⟩+2\beta (\text{mbox}\mathrm{Re}(u\overline{v})u,w{)}_{L^2(\mathcal{D})}.$$

$$-\Delta u=\lambda u-Vu+\Omega {\mathcal{L}}_{3}u-\beta |u{|}^{2}u=:g.$$

$$V_{h,k}=\{v\in H_0^1(\mathcal{D})\cap C^0(\overline{\mathcal{D}})|v{|}_K\in {\mathbb{P}}^k(K)\text{mbox}forallK\in {\mathcal{T}}_h\},$$

$$E(u_h)=\underset{v\in V_h\cap \mathbb{S}}{\inf }E(v).$$

Full text

1 Introduction

At extreme temperatures close to 0 Kelvin, dilute bosonic gases can condensate to a fascinating state of matter: a so-called Bose-Einstein condensate (BEC), cf. [11, 28, 39]. In such a condensate most of the particles occupy the same quantum state which makes it effectively behave like a macroscopically observable “super atom”. This super atom allows to study various quantum phenomena, where superfluidity, i.e., frictionless flow, is perhaps one of most intriguing ones [37]. Superfluidity can be verified by rotating the BEC with a stirring potential, triggering the formation of vortices with a quantized circulation if the rotation frequency is appropriately chosen.

The ground states of a rotating BEC can be described as the global minimizers $u\in H^{1}_{0}(\mathcal{D},\mathbb{C})$ of the Gross–Pitaevskii energy functional

[TABLE]

under the constraint that $\int_{\mathcal{D}}|w|^{2}\hskip 2.0pt\mbox{d}x=1$. Here, $\mathcal{D}\subset\mathbb{R}^{d}$ (for $d=2,3$) is the computational domain, $V\in L^{\infty}(\mathcal{D},\mathbb{R}_{\geq 0})$ models a trapping potential, $\beta\in\mathbb{R}_{\geq}0$ describes repulsive particle interactions, $\Omega\in\mathbb{R}$ is the angular velocity and $\mathcal{L}_{3}=-\mathrm{i}\left(x_{1}\partial_{x_{2}}-x_{2}\partial_{x_{1}}\right)$ the $x_{3}$-component of the angular momentum. The existence and non-existence of ground states in the above setting is e.g. studied in [10].

In practice, a numerical scheme for computing ground states involves two steps: 1. a space discretization which essentially requires a minimization of $E$ over a finite dimensional (discrete) space and 2. an appropriate iterative solver that allows to find the corresponding discrete minimizers. We will not discuss the second step here but just refer to [1, 2, 3, 4, 6, 7, 8, 9, 10, 13, 16, 25, 26, 33, 34] and the references therein for different approaches and complementary analytical results (as far as they exist). The focus of this paper is on the first step, i.e., the approximation of analytical ground states in discrete spaces and a corresponding error quantification. Typical space discretizations involve spectral and pseudo-spectral methods [6, 15], quasi-uniform finite element spaces [15] or adaptive finite elements [19, 23, 29]. Two-grid post-processing techniques in FE spaces are presented in [14, 31].

When it comes to the analysis of these space discretizations, and in particular the derivation of a priori error estimates, the case of rotating Bose–Einstein condensates has not been regarded yet. For non-rotating condensates (i.e. for $\Omega=0$) the first results were obtained by Zhou [41, 42] who considered arbitrary finite dimensional spaces and derived corresponding $H^{1}$-error estimates. The estimates were however not sharp yet due to a suboptimal contribution from the $L^{2}$-norm of the error. This gap could be closed by Cancès et al. [15] who revisited the problem and provided the first sharp error analysis, covering $L^{2}$- and $H^{1}$-error estimates together with estimates for the ground state energy and the so-called ground state eigenvalue. Shortly after, the results were generalized yet another time by Zhou and coworkers [21] by adding the effect of Hartree-type potentials to the analysis. Error estimates for finite element discretizations of the Kohn-Sham and Hartree-Fock equation were obtained in [17, 18, 36] and generalized finite element approximations for the standard Gross–Pitaevskii model were studied in [31, 32].

Extending the aforementioned error analysis to rotating BECs is however not straightforward. This is because the setting for $\Omega=0$ is fully real-valued, whereas the problem becomes complex-valued for $\Omega\not=0$. As a result, uniqueness of ground states is lost, cf. [10]. One particular issue are complex phase shifts of the form $e^{\mathrm{i}\theta}$ for angles $\theta\in[0,2\pi)$, i.e., if $u$ is an $L^{2}$-normalized minimizer of $E$ then $e^{\mathrm{i}\theta}u$ is also a normalized minimizer for any $\theta$. This means that ground states are not locally unique. Our analysis is therefore based on restricting the Euler-Lagrange functional for the constrained minimization problem to tangent spaces where the second Fréchet derivative of the Euler-Lagrange functional (evaluated at a ground state) has a bounded inverse. On these spaces, we can derive $H^{1}$-error estimates by applying abstract approximation results for discrete solutions to nonlinear equations. For estimating the $L^{2}$- and eigenvalue errors we change the approach and draw inspiration from the Kohn-Sham setting [17, 18, 36] (where also missing local uniqueness is faced) and consider an $L^{2}$-orthogonal decomposition of the error with contributions from the tangent spaces. This decomposition gives rise to strong error identities that we then use to derive different types of error estimates in $\mathbb{P}^{k}$-Lagrange finite element spaces $V_{h,k}$. In particular, we obtain convergence of order $\mathcal{O}(h^{k})$ for the $H^{1}$-error and convergence of order $\mathcal{O}(h^{k+1})$ for the $L^{2}$-error. We also derive optimal order estimates for the ground state energy and the ground state eigenvalue, both of order $\mathcal{O}(h^{2k})$. With this, our findings directly generalize the results of [15] to the Gross–Pitaevskii equation with rotation and to polynomial spaces of higher order.

Outline: The paper is structured as follows. In Section 2, we give an introduction to the analytical setting of the Gross-Pitaevskii equation with angular momentum, we discuss the issue of uniqueness and state our basic assumptions. The finite element discretization and our main results regarding the corresponding approximation properties are presented in Section 3 and afterwards illustrated by numerical experiments in Section 4. Section 5 is devoted to the proofs of the main results: Sections 5.1 and 5.2 provide an abstract asymptotic convergence analysis together with some preliminary estimates that allow to control energy- and eigenvalue-errors in terms of the $H^{1}$-error. The $H^{1}$-error itself is then estimated in Section 5.3, whereas the $L^{2}$-error estimate follows in Section 5.4. An error estimate of optimal order for the eigenvalue is then proved in Section 5.5. Complementary results on the inf-sup stability are given in the appendix.

2 Analytical setting

In the following we shall consider the Gross-Piteavskii equation on a computational domain $\mathcal{D}$ with homogeneous Dirichlet boundary conditions, where we assume that

1. (A1)

$\mathcal{D}\in\mathbb{R}^{d}$ is a bounded, convex domain for $d=2,3$ with polygonal boundary.

Since we are concerned with the minimization of a real-valued energy functional $E$, it is natural to consider all Hilbert spaces as real-linear spaces. To be precise, we interpret the Lebesgue space $L^{2}(\mathcal{D}):=L^{2}(\mathcal{D},\mathbb{C})$ as a real Hilbert space equipped with the inner product

[TABLE]

where $\overline{z}$ denotes the complex conjugate of a number $z\in\mathbb{C}$. Analogously, the Sobolev space $H_{0}^{1}(\mathcal{D}):=H_{0}^{1}(\mathcal{D},\mathbb{C})$ is also seen as a real Hilbert space equipped with the inner product

[TABLE]

The corresponding $(\text{real})$ dual space is denoted by \big{(}H_{0}^{1}(\mathcal{D})\big{)}^{*}=H^{-1}(\mathcal{D}):=H^{-1}(\mathcal{D},\mathbb{C}), with $\langle\cdot,\cdot\rangle:=\langle\cdot,\cdot\rangle_{H^{-1}(\mathcal{D}),H^{1}_{0}(\mathcal{D}}$ being the canonical duality pairing in $H_{0}^{1}(\mathcal{D})$.

For the data, we assume that

1. (A2)

$V\in L^{\infty}(\mathcal{D},\mathbb{R}_{\geq 0})$ denotes a non-negative potential and 2. (A3)

the parameter $\beta\in\mathbb{R}_{\geq 0}$ characterizes the strength of (repulsive) particle interactions.

To model condensates that are rotating around the $x_{3}$-axis of the coordinate system, we require the third component of the angular momentum operator, which we denote by $\mathcal{L}_{3}:=-\mathrm{i}\left(x_{1}\partial_{x_{2}}-x_{2}\partial_{x_{1}}\right)$. The angular velocity of that component shall be given by $\Omega\in\mathbb{R}$. To guarantee existence of ground states, we formally need to add a balancing assumption that the trapping frequencies in $x_{1}$- and $x_{2}$-direction are larger than the angular frequency, i.e.

1. (A4)

there is $\varepsilon>0$ such that

[TABLE]

Physically speaking, assumption (A4) ensures that centrifugal forces do not become too strong in comparison to the strength of the trapping potential $V$.

With the above notation, we can define the Gross–Pitaevskii energy functional

[TABLE]

for $w\in H^{1}_{0}(\mathcal{D})$ by

[TABLE]

Since ground states will be defined as constrained minimizers of $E$, we briefly need to discuss the well-posedness. We start with the observation that $E$ is indeed real-valued and bounded from below: Using assumption (A4), we let $V_{\mbox{\rm\tiny R}}(x):=V(x)-\tfrac{1}{4}\,\Omega^{2}\,|x|^{2}\geq 0$ denote a non-negative modified potential and we let the covariant gradient $\nabla_{\mbox{\rm\tiny R}}$ be given by

[TABLE]

for the divergence-free vector field $R(x):=(x_{2},-x_{1},0)$ if $d=3$ and $R(x):=(x_{2},-x_{1})$ if $d=2$. With this

[TABLE]

defines an inner-product on $H^{1}_{0}(\mathcal{D})$ with induced norm $\|v\|_{\mbox{\rm\tiny R}}:=\sqrt{(v,v)_{\mbox{\rm\tiny R}}}$. It is now easy to see (cf. [24, 30]) that $E(w)$ can be equivalently rewritten as

[TABLE]

Hence, $E$ is a real functional that is bounded from below by zero.

A ground state of a Bose–Einstein condensate in the rotational frame is defined as a stable stationary state in the lowest possible energy level. Mathematically, it is characterized as a global minimizer $u$ of $E$ under the normalization condition that $\int_{\mathcal{D}}|u|^{2}\hskip 2.0pt\mbox{d}x=1$, i.e.

[TABLE]

The problem is a non-convex minimization problem on the $L^{2}$-sphere in $H^{1}_{0}(\mathcal{D})$. Exploiting the reformulation (4), the energy functional is positively bounded from below on the manifold $\mathbb{S}$ with

[TABLE]

Additionally, $E$ is weakly lower semi-continuous on $H^{1}_{0}(\mathcal{D})$ and it can be proved, cf. [27, Lemma 2.2], that $\|\cdot\|_{\mbox{\rm\tiny R}}$ is equivalent to the standard $H^{1}$-norm on $H^{1}_{0}(\mathcal{D})$. In particular, there exist constants $0<c(\varepsilon)<C(V,\Omega)<\infty$, depending on $V$, $\Omega$ and $\varepsilon$ respectively, such that

[TABLE]

It can be shown that $c(\varepsilon)=\mathcal{O}(\sqrt{\varepsilon})$, i.e. the lower constant degenerates for $\varepsilon\rightarrow 0$, which corresponds to a violation of the physical assumption (A4). Now, since $\|\cdot\|_{\mbox{\rm\tiny R}}$ is equivalent to the standard $H^{1}$-norm on $H^{1}_{0}(\mathcal{D})$, we easily see that $E$ fulfills the Palais–Smale condition on $\mathbb{S}$. As an even functional that fulfills the Palais–Smale condition we can apply classical Lusternik-Schnirelmann theory to $E$, cf. [35, Chapt. 4, Theorem 5.5], and conclude that $E$ has infinitely many critical values. The critical values are unbounded from above and bounded from below, where the smallest critical value is given by the ground state energy $E_{0}\in\mathbb{R}$. In particular, there exists a sequence of critical values $E_{k}\in\mathbb{R}$ with $-\infty<E_{0}\leq E_{1}\leq E_{2}\leq\cdots E_{k}\rightarrow\infty$ for $k\rightarrow\infty$ and corresponding critical points $u_{k}\in\mathbb{S}$ with $E_{k}=E(u_{k})$. A critical point to a smallest critical value $E_{k}=E_{0}$ is a ground state as defined in (5). Note that a critical point $u_{k}\in\mathbb{S}$ to $E_{k}$ is not unique due to the invariance of $E$ under rotations of $u_{k}$ on the complex unit circle. We summarize the previous discussion in the following corollary.

Corollary 2.1** (Existence of ground states).**

Under assumptions (A1)-(A4), there exists at least one (physically meaningful) ground state $u\in\mathbb{S}$ to problem (5) and it holds $E(u)>0$.

2.1 Uniqueness of ground states and quasi-isolation

Uniqueness of ground states can only be expected up to a constant phase factor $e^{\mathrm{i}\theta}$ with angle $\theta\in[0,2\pi)$. This is because for any $u\in\mathbb{S}$ we have $e^{\mathrm{i}\theta}u\in\mathbb{S}$ and $E(e^{\mathrm{i}\theta}u)=E(u)$. Uniqueness of the ground state density $|u|^{2}$ (which is independent of $e^{\mathrm{i}\theta}$) is possible, but can only be expected for frequencies $\Omega$ below a certain first critical frequency. For larger frequencies the ground state density $|u|^{2}$ can become non-unique, cf. [10, Section 6.2]. The phenomenon is not yet fully understood, but numerical experiments and theoretical considerations indicate that missing uniqueness of $|u|^{2}$ is a rather rare occasion that only happens for certain critical values, where a particular vortex state becomes energetically equally preferable to a ground state with less vortices (or without any vortices at all), cf. [10]. To account for the general lack of uniqueness, we define the set of ground states of (5) by

[TABLE]

and introduce the term of quasi-isolated ground states. Loosely speaking, we want that for every ground state $u\in\mathcal{U}$, there is a small $\delta$-neighborhood $\mathbb{B}_{\delta}(u):=\{v\in\mathbb{S}\,|\,\|v-u\|_{H^{1}(\mathcal{D})}<\delta\}$ of $u$ where all other ground states in $\mathbb{B}_{\delta}(u)$ are precisely phase shifts of $u$, i.e.,

[TABLE]

This condition can be expressed through curves $\gamma:t\mapsto\gamma(t)\in\mathbb{S}$ crossing a ground state with $\gamma(0)=u\in\mathcal{U}$. Since $\gamma$ remains on $\mathbb{S}$, we have $\gamma^{\prime}(0)\in T_{u}\mathbb{S}$ with the tangent space in $u$ given by

[TABLE]

Now, if condition (7) holds, the energy level $t\mapsto E(\gamma(t))$ has a strict local minimum in $t=0$ as long as $\gamma$ does not follow the direction $\mathrm{i}u$ of complex phase shifts (since $\tfrac{\mbox{\scriptsize d}}{\mbox{\scriptsize d}t}e^{\mathrm{i}t}u|_{t=0}=\mathrm{i}u$). This direction is blocked by imposing the condition $\gamma^{\prime}(0)\in T_{\mathrm{i}u}\mathbb{S}$. Hence, for all curves with $\gamma^{\prime}(0)\in T_{u}\mathbb{S}\cap T_{\mathrm{i}u}\mathbb{S}$, the necessary conditions for local minima imply $\tfrac{\mbox{\scriptsize d}^{\,2}}{\mbox{\scriptsize d}t^{2}}E(\gamma(t))|_{t=0}\geq 0$. The following definition formalizes these considerations, with the slight modification that the necessary condition for strict local minima is replaced by a sufficient condition.

Definition 2.2** (Quasi-isolated ground state).**

A ground state $u\in\mathcal{U}$ is called a quasi-isolated ground state if for any differentiable curve $\gamma:(-\varepsilon,\varepsilon)\rightarrow\mathbb{S}$ with $\gamma(0)=u$ and $\gamma^{\prime}(0)\in(T_{u}\mathbb{S}\cap T_{\mathrm{i}u}\mathbb{S})\setminus\{0\}$ it holds

[TABLE]

Since we will require the space $T_{u}\mathbb{S}\cap T_{\mathrm{i}u}\mathbb{S}$ often, we introduce for brevity

[TABLE]

The notation is motivated by the simple observation that $\mathscr{H}_{u}^{\mathbb{C},\perp}$ is the orthogonal complement of $u$ with respect to the complex $L^{2}$-inner product, whereas $T_{\mathrm{i}u}\mathbb{S}$ is the orthogonal complement of $\mathrm{i}u$ with respect to the real $L^{2}$-inner product. With Definition 2.2 we make the following assumption on the ground states:

1. (A5)

Every ground state $u\in\mathcal{U}$ is a quasi-isolated ground state.

With other words, two ground states are either well separated or complex phase shifts of each other.

In general, assumption (A5) is reasonable and its validity in a given setting can be always checked numerically. Even though the energy $E$ can be also invariant under rotations of the coordinate system in particular cases (such as for spherical domains $\mathcal{D}$ and radial symmetric potentials $V$), we can still expect (A5) to hold. The reason is that the corresponding critical direction is typically not an element of $H^{1}_{0}(\mathcal{D})$. For example, if $E(u)=E(u_{\theta})$ with $u_{\theta}(x):=u(\mathbf{R}_{\theta}x)$ for any rotation matrix

[TABLE]

then the critical direction of the curve (on which the energy level would remain the same) is given by $\tfrac{\mbox{\scriptsize d}}{\mbox{\scriptsize d}\theta}u(\mathbf{R}_{\theta}x)_{|\theta=0}=-(x\times\nabla)u$. Since we cannot expect $\nabla u$ to vanish on $\partial\mathcal{D}$, the function $-(x\times\nabla)u$ is not an element of $H^{1}_{0}(\mathcal{D})$ and, hence, $\gamma^{\prime}(0)=-(x\times\nabla)u$ is not admissible. In fact, this reasoning also implies that $-(x\times\nabla)u$ cannot be an eigenfunction of $E^{\prime\prime}(u)-\lambda\mathcal{I}$ to the eigenvalue [math] (see Proposition 2.4 below for a corresponding connection).

2.2 Characterization through the Gross–Pitaevskii equation

Since every ground state $u$ is a critical point of $E$ on $\mathbb{S}$, we can characterize it by the corresponding Euler–Lagrange equations seeking $u\in\mathbb{S}$ with Lagrange multiplier $\lambda\in\mathbb{R}_{>0}$ such that

[TABLE]

where $E^{\prime}:H^{1}_{0}(\mathcal{D})\rightarrow H^{-1}(\mathcal{D})$ denotes the Fréchet derivative of $E$ and $M:H^{1}_{0}(\mathcal{D})\rightarrow\mathbb{R}$ is the scaled constraint functional given by $M(v):=\tfrac{1}{2}\left(\int_{\mathcal{D}}|v|^{2}\hskip 2.0pt\mbox{d}x-1\right)$. Computing the Fréchet derivatives, we obtain

[TABLE]

and

[TABLE]

As a compact notation, we introduce for arbitrary $u\in H^{1}_{0}(\mathcal{D})$, the $u$-linearized approximation of $E^{\prime}$ by $\mathcal{A}_{|u|}:H^{1}_{0}(\mathcal{D})\rightarrow H^{-1}(\mathcal{D})$ with

[TABLE]

Note that $\mathcal{A}_{|u|}$ is a self-adjoint operator on $L^{2}(\mathcal{D})$ and that it is obviously elliptic with

[TABLE]

where $\alpha=c(\varepsilon)^{2}=\mathcal{O}(\varepsilon)$ is the square of the constant from (6). As $E^{\prime}(v)=\mathcal{A}_{|v|}(v)$, we can write the Euler–Lagrange equations for the critical points of $E$ equivalently as: find $u\in\mathbb{S}$ and $\lambda>0$ such that

[TABLE]

The problem can be interpreted as an eigenvalue problem with eigenvector nonlinearity and is known as the Gross–Pitaevskii eigenvalue problem (GPEVP). In the sense of distributions, problem (10) equivalently reads: find $\lambda>0$ and $u\in H_{0}^{1}(\mathcal{D})$ such that $\int_{\mathcal{D}}|u|^{2}dx=1$ and

[TABLE]

The equivalence of (10) and (11) is seen by noticing that the real part in front of the integrals can be dropped by using both $v$ and $\mathrm{i}v$ as test functions.

Any eigenfunction to (11) that is not a ground state is called a meta-stable excited state of the Bose–Einstein condensate.

Remark 2.3** (Ground state eigenvalue).**

An eigenvalue $\lambda$ that belongs to a ground state $u\in\mathbb{S}$ in the sense of (5) is called a ground state eigenvalue. Note however, that it is in general not possible to directly use the GPEVP (11) as a defining equation for ground states. The reason is that a ground state eigenvalue is not necessarily the smallest eigenvalue of the GPEVP (11) and an eigenfunction to the smallest eigenvalue can potentially be an excited state, cf. [1, Section 6.3]. Only for the case without rotation, i.e. $\Omega=0$, the smallest eigenvalue is guaranteed to coincide with the (unique) ground state eigenvalue [15].

With the above notation, we rephrase assumption (A5) (quasi-isolation of ground states) in terms of an inf-sup stability condition that we can exploit in the error analysis. The inf-sup condition is obtained from the property $\tfrac{\mbox{\scriptsize d}^{2}}{\mbox{\scriptsize d}t^{2}}E(\gamma(t))|_{t=0}>0$ by explicitly computing the derivative and exploiting the Euler–Lagrange equations. The argument is standard, but for completeness the proof of the following proposition is given in the appendix, Section A.

Proposition 2.4**.**

Assume (A1)-(A5), let $u\in\mathcal{U}$ be a ground state in the sense of (5) and $\lambda\in\mathbb{R}_{>0}$ the corresponding ground state eigenvalue given by (10). If $\mathcal{I}:H^{1}_{0}(\mathcal{D})\rightarrow H^{-1}(\mathcal{D})$ denotes the canonical identification $\mathcal{I}v:=(v,\cdot)_{L^{2}(\mathcal{D})}$, then the operator $E^{\prime\prime}(u)-\lambda\mathcal{I}$ is inf-sup stable on $\mathscr{H}_{u}^{\mathbb{C},\perp}$, namely there exists a constant $\alpha_{1}>0$ such that

[TABLE]

In particular, $E^{\prime\prime}(u)-\lambda\mathcal{I}$ has a bounded inverse on $\mathscr{H}_{u}^{\mathbb{C},\perp}$, which follows from standard theory for indefinite problems, cf. [5].

Similar assumptions such as the inf-sup condition (12) can be frequently found in the literature on the error analysis of finite element approximations for energy minimization problems in quantum physics, cf. [18, 19, 21, 36].

Using the notation in (8), $E^{\prime\prime}(u)$ in (12) can be computed as

[TABLE]

By selecting $v=\mathrm{i}u$ in (13) and exploiting that $\mathrm{i}u$ is an eigenfunction of $\mathcal{A}_{|u|}$ with eigenvalue $\lambda$, we have $\langle E^{\prime\prime}(u)\mathrm{i}u,w\rangle=\lambda\,(\mathrm{i}u,w)_{L^{2}(\mathcal{D})}$ for all $w\in H^{1}_{0}(\mathcal{D})$. Hence, $E^{\prime\prime}(u)-\lambda\mathcal{I}$ is singular on $\mbox{span}\{\mathrm{i}u\}$ (the directions of complex phase shifts) as previously predicted in Section 2.1.

We conclude this section with a natural regularity statement for ground states that we will exploit frequently in the convergence analysis.

Lemma 2.5** ($H^{2}$-regularity of ground states).**

Assume (A1)-(A3). Let $u$ be a ground state and $\lambda>0$ the corresponding ground state eigenvalue as in (11). Then, $u\in C^{0}(\overline{\mathcal{D}})\cap H^{2}(\mathcal{D}).$

Proof.

The proof is straightforward by writing (11) as

[TABLE]

Our assumptions on $V$ and the Sobolev embedding $H^{1}_{0}(\mathcal{D})\hookrightarrow L^{6}(\mathcal{D})$ for $d=2,3$, yield $g\in L^{2}(\mathcal{D})$. Hence, together with the convexity of $\mathcal{D}$, standard elliptic regularity theory (cf. [22, Theorem 3.2.2]) for the Poisson problem $-\Delta u=g$ guarantees $u\in C^{0}(\overline{\mathcal{D}})\cap H^{2}(\mathcal{D})$. ∎

3 Finite element discretization and main results

In the following we investigate the accuracy of discrete ground states to the minimization problem (5) in finite element spaces. For that, we consider a shape regular family $\{\mathcal{T}_{h}\}$ of nested conforming triangulations over the domain $\mathcal{D}$ with mesh size $h>0$. On each mesh $\mathcal{T}_{h}$, the corresponding $\mathbb{P}^{k}$-Lagrange finite element space is given by

[TABLE]

i.e., $V_{h,k}$ is the space of continuous functions in $H_{0}^{1}(\mathcal{D})$ such that their restriction to any element $K$ of $\mathcal{T}_{h}$ is a polynomial of degree not greater than $k$. For simplicity we shall mainly denote $V_{h}:=V_{h,k}\,$ as long as the specific value of $k$ is not crucial. If we want to stress the polynomial order, we will use the index “$h,k$” instead of “$h$”.

Let us now consider the finite element approximations for ground states of the energy minimization problem (5). A ground state in the finite element space $V_{h}$ is characterized by the minimization problem seeking $u_{h}\in V_{h}\cap\mathbb{S}$ such that

[TABLE]

In the following we shall call such minimizers discrete ground states. As in the continuous case, discrete ground states $u_{h}$ can be at most unique up to constant phase shifts. We denote the set of all discrete ground states by

[TABLE]

Analogously as in the continuous setting, any minimizer $u_{h}$ of (15) can be expressed through the corresponding Euler–Lagrange equations. Hence, there exists $\lambda_{h}\in\mathbb{R}_{>0}$ such that it holds

[TABLE]

Remark 3.1** (Same phase).**

Due to missing uniqueness caused by phase shifts, we can only make a reasonable comparison between an actual ground state $u$ and a discrete ground state $u_{h}$ if their phases are aligned. In the following we will call $u$ and $u_{h}$ to be in the same phase if their complex $L^{2}$-inner product is a value on the positive real axis, i.e. $\int_{\mathcal{D}}u_{h}\overline{u}\hskip 2.0pt\mbox{d}x\in\mathbb{R}_{\geq 0}$.

To formulate our main results and to keep the proofs compact, we write from now on “$a\lesssim b$” to abbreviate “$a\leq Cb$” for a constant $C>0$ that may depend on $\beta$, $\mathcal{D}$, $\Omega$, $V$ and $u$ but is independent of the mesh size $h$.

We are now prepared to present our main theorems. The first result guarantees that for each discrete ground state $u_{h}\in\mathcal{U}_{h}$ (on a sufficiently fine mesh) there is a $H^{1}$-neighbourhood of $u_{h}$ in which we will find a unique exact ground state $u\in\mathcal{U}$ that is in the same phase as $u_{h}$. Furthermore, such $u_{h}$ is a quasi-best approximation of $u$ in the $H^{1}$-norm.

Theorem 3.2**.**

Assume (A1)-(A4), then for every sufficiently small $\delta>0$, there is a $h_{\delta}>0$ such that for each $h\leq h_{\delta}$ and each $u_{h}\in\mathcal{U}_{h}$, there exists a ground state $u\in\mathcal{U}$ that is in the same phase, i.e., $\int_{\mathcal{D}}u_{h}\overline{u}\hskip 2.0pt\mbox{d}x\in\mathbb{R}_{\geq 0}$, and such that

[TABLE]

Furthermore, the error in energy is bounded by

[TABLE]

and distance between the ground state eigenvalue $\lambda$ to $u\in\mathcal{U}$ and the discrete ground state eigenvalue $\lambda_{h}$ to $u_{h}\in\mathcal{U}_{h}$ can be bounded by

[TABLE]

If additionally assumption (A5) is satisfied, then $u$ from the first part of the theorem is unique (in the $\delta$-neighborhood); $u_{h}$ is a $H^{1}$-quasi-best approximation to $u$

[TABLE]

and it holds the $L^{2}$-error estimate

[TABLE]

The proof of Theorem 3.2 is postponed to Section 5 where we will also state additional complementary results that address the general approximability of all ground states $u\in\mathcal{U}$ (which is important if the ground state density is not unqiue).

Theorem 3.2 directly yields optimal convergence rates (in $L^{2}(\mathcal{D})$, $H^{1}(\mathcal{D})$ and the energy error) for Lagrange spaces of arbitrary polynomial order $\mathbb{P}^{k}$, provided that the limiting ground state $u$ is sufficiently smooth. However, it does not contain a corresponding explicit estimate for the error in the ground state eigenvalue. For this we require an additional regularity assumption. The arising optimal convergence rates together with the necessary assumptions are summarized in the following result.

Theorem 3.3**.**

Assume (A1)-(A5) and let $(\lambda_{h,k},u_{h,k})$ denote a discrete ground state pair of (17) in $\mathbb{R}_{>0}\times\,\mathcal{U}_{h,k}$ for $k\in\mathbb{N}_{>0}$ and sufficiently small $h$. Let further $(\lambda,u)\in\mathbb{R}_{>0}\times\,\mathcal{U}$ be the corresponding (uniquely characterized) exact ground state pair with $u_{h,k}\in T_{\mathrm{i}u}\mathbb{S}$ from Theorem 3.2. If $u\in H^{k+1}(\mathcal{D})$, then it holds

[TABLE]

and for the error in energy

[TABLE]

Finally, let $w_{u|u|^{2}}\in\mathscr{H}_{u}^{\mathbb{C},\perp}$ denote the unique solution to the dual problem

[TABLE]

(which exists by Proposition 2.4). If $w_{u|u|^{2}}$ is sufficiently smooth in the sense that $w_{u|u|^{2}}\in H^{k+1}(\mathcal{D})$ with the corresponding regularity estimate

[TABLE]

then it also holds the following sharp error estimate for the eigenvalue

[TABLE]

The proof of Theorem 3.3 is presented in Section 5.5.

Remark 3.4** (Validity of the regularity assumptions in Theorem 3.3).**

Note that we always have $u\in H^{2}(\mathcal{D})$ and $w_{u|u|^{2}}\in H^{2}(\mathcal{D})$ with $\|w_{u|u|^{2}}\|_{H^{2}(\mathcal{D})}\lesssim\|u\|_{L^{6}(\mathcal{D})}^{3}$, which shows that the regularity assumptions of Theorem 3.3 are at least fulfilled with $k=1$ (cf. Lemma 5.11 below). For $k>1$, we would need additional regularity of the potential $V$, as well as a smooth domain $\mathcal{D}$ in order to formally guarantee the required regularity. However, smoothness of $\mathcal{D}$ would be in conflict with our assumption that $\mathcal{D}$ has a polygonal boundary. In practical situations, this is not an issue because the computational domain is selected sufficiently large (estimated by the Thomas-Fermi radius) so that the ground state decays quickly to zero towards the boundary. Hence, potential boundary effects on the regularity are negligible. The same holds for the solution $w_{u|u|^{2}}$ to the dual problem, which has the quickly decaying function $u|u|^{2}$ as a source term.

It is also worth to mention that in [15], additional regularity for rectangular domains could be rigorously proved with a prolongation by reflection argument. The argument does however no longer work in the case of rotation since the $\mathcal{L}_{3}u$-term in (14) prevents us from writing the equation as a Poisson problem with a source term in $H^{1}_{0}(\mathcal{D})$.

The rest of the paper is mainly dedicated to the proofs of the two above stated theorems. Before we turn our attention to these proofs, we present two numerical experiments to support the convergence rates (22), (23) and (25). Furthermore, we exemplarily computed the lower part of the spectrum of $E^{\prime\prime}(u)$ to verify that assumption (A5) is fulfilled in our test cases.

4 Numerical experiments

In this section we present numerical experiments for $\mathbb{P}^{1}$-FEM spaces to validate our theoretical convergence rates for ground states (in the $H^{1}$- and $L^{2}$-norm) and the respective eigenvalue and energy. Moreover, to illustrate the validity of assumption (A5) in our experiments, we compute the spectrum of $E^{\prime\prime}(u)$ for a ground state $u$ and verify that the ground state eigenvalue $\lambda$ is only appearing once and at the bottom of the spectrum.

In the following experiments we consider the energy minimization problem (5) on the square $\mathcal{D}=[-R,R]^{2}$, where the particular choice of $R$ will be specified in the two test cases. The problem is discretized with $\mathbb{P}^{1}$–FEM based on a uniform triangulation of $\mathcal{D}$ consisting of $2N^{2}_{h}$ triangles ($N_{h}\in\mathbb{N}$) with corresponding mesh size $h=12/N_{h}$. The values are specified in the test cases. The resulting discrete minimization problem (15) is solved with an energy-adaptive Riemannian gradient method suggested in [33], whose generalization to the GPE with rotation is formulated in [1, Section 6]. As a measure of accuracy we will stop the solver whenever the difference between two consecutive energies, i.e, $E(u^{n+1}_{h})-E(u^{n}_{h})$, falls below the tolerance of $\varepsilon=10^{-10}$. The starting value for the Riemannian gradient method is chosen, according to [10, Section 6.1], as the $L^{2}$-normalized interpolation of the function

[TABLE]

in $V_{h}$. In the experiments we chose the harmonic potential $V$ as

[TABLE]

where the precise trapping frequencies $\gamma_{x_{1}}$ and $\gamma_{x_{2}}$ are specified in the two test cases below.

4.1 Model problem 1

In the first model problem, we consider $\mathcal{D}=[-6,6]^{2}$, the angular velocity $\Omega=1.2$ and the repulsion parameter $\beta=100$. The trapping frequencies are selected as $\gamma_{x_{1}}=0.9$ and $\gamma_{x_{2}}=1.2$. The discrete ground states are computed on a mesh of decreasing mesh size with values $h=12\cdot 2^{-4}\,,12\cdot 2^{-5}\,,\,12\cdot 2^{-6}\,,\,12\cdot 2^{-7}\,,\,12\cdot 2^{-8}$. The reference ground state $u$ is computed on a fine mesh with size $h=12\cdot 2^{-9}$. The problem configuration is such that the ground state $u$ develops quantized vortices. The corresponding density $|u|^{2}$ of the rotating BEC is plotted (top-left) in Figure LABEL:gs-spectrum. Here we identified the ground state energy with $E(u)\approx 1.6440$ and the corresponding ground state eigenvalue with $\lambda\approx 4.4488$. The resulting convergence rates for the ground state (in the $H^{1}$- and in the $L^{2}$-norm), for the eigenvalue and the energy are plotted in Figure 2 and are in accordance with the predicted rates from Theorem 3.3.

The numerically computed spectrum of the second Frèchet derivative $E^{\prime\prime}(u)$, cf. (13), is presented in Figure LABEL:gs-spectrum. As expected, the ground state eigenvalue is a simple eigenvalue of $E^{\prime\prime}(u)$ that only appears at the very bottom of the spectrum. The corresponding eigenfunction is $\mathrm{i}u$. Note that this shows that assumption (A5) is fulfilled, because it implies (with the Courant–Fischer min-max principle)

[TABLE]

Hence, the smallest eigenvalue of $E^{\prime\prime}(u)-\lambda\mathcal{I}$ on $\mathscr{H}_{u}^{\mathbb{C},\perp}$ is positive. Together with the Gårding inequality, this implies that $E^{\prime\prime}(u)-\lambda\mathcal{I}$ has a bounded inverse on $\mathscr{H}_{u}^{\mathbb{C},\perp}$ and (A5) is fulfilled.

Remark 4.1** (Independence of the spectrum from the phase).**

Note that the spectrum of $E^{\prime\prime}(u)$ does not depend on the phase of $u$, i.e., for all $\theta\in[0,2\pi)$ the eigenvalues of $E^{\prime\prime}(u)$ and $E^{\prime\prime}(e^{\mathrm{i}\theta}u)$ are identical. This observation can be proved with the min-max principle for eigenvalues.

4.2 Model problem 2

In the second model problem we set $\mathcal{D}=[-8,8]^{2}$, the trapping frequencies of the potential to $\gamma_{x_{1}}=1.1$ and $\gamma_{x_{2}}=0.9$; the repulsion parameter to $\beta=400$ and the angular velocity to $\Omega=1.1$. Here the discrete ground states are computed on a mesh of decreasing mesh size with values $h=16\cdot 2^{-4}\,,16\cdot 2^{-5}\,,\,16\cdot 2^{-6}\,,\,16\cdot 2^{-7}\,,\,16\cdot 2^{-8}$. The reference ground state $u$ is computed on a fine mesh with size $h=16\cdot 2^{-9}$. The reference ground state has quantized vortices and we identified the ground state energy with $E(u)\approx 2.9107$ and the corresponding ground state eigenvalue with $\lambda\approx 8.2055$. The corresponding density $|u|^{2}$ (bottom-left ) is plotted in Figure LABEL:gs-spectrum. Again, the convergence rates for the $H^{1}$-error, the $L^{2}$-error, the energy error and the eigenvalue error, as plotted in Figure 3, are in accordance with Theorem 3.3.

The second Frèchet derivative $E^{\prime\prime}(u)$ is presented in Figure LABEL:gs-spectrum. Again, we find that the ground state eigenvalue is the smallest eigenvalue and appears only once in the spectrum, showing (with the same arguments as in Section 4.1) that assumption (A5) is fulfilled.

5 Proofs of Theorems 3.2 and 3.3

We now turn towards the proofs of our main results which are split into several subsections.

5.1 Asymptotic convergence results

The aim of this first proof section is to establish the existence of an actual unique quasi-isolated ground state in a sufficiently small $\delta$-neighborhood of each discrete ground state (for $h$ small enough).

First, note that with assumptions (A1)-(A4) and the norm equivalence (6) together with the fact that $u_{h}$ minimizes the energy of $E$ on the nested spaces $V_{h}$, we can conclude that a family of discrete ground states of (15) must be uniformly bounded in the $H^{1}$-norm, i.e., there exists a constant $C>0$ such that, for all sufficiently small $h$, it holds

[TABLE]

Now we will use the above uniform boundedness along with the lower semi-continuity of the energy functional $E$ in the weak topology of $\mathbb{S}$ to show the $H^{1}$-convergence of discrete ground states (up to subsequences) to an actual ground state $u\in\mathcal{U}$. The subsequent proposition also proves the first part of Theorem 3.2.

Proposition 5.1**.**

Assume (A1)-(A4) holds and let $\{(\lambda_{h},u_{h})\}$ be a family of discrete ground state pairs of (17) in $\mathbb{R}^{+}\times\mathcal{U}_{h}$. Then there exists a subsequence $\{(\lambda_{h_{j}},u_{h_{j}})\}_{j\in\mathbb{N}}\subset\{(\lambda_{h},u_{h})\}$ and an exact ground state $(\lambda_{0},u_{0})\in\mathbb{R}^{+}\times\mathcal{U}$ of $\eqref{eigen_value_problem_2}$ such that

[TABLE]

[TABLE]

Furthermore, for each sufficiently small $\delta>0$, there is a $h_{\delta}>0$ such that for each $h\leq h_{\delta}$ and each $u_{h}\in\mathcal{U}_{h}$, there exists a ground state $u\in\mathcal{U}$ with $\int_{\mathcal{D}}u_{h}\overline{u}\hskip 2.0pt\mbox{d}x\in\mathbb{R}_{\geq 0}$ and such that

[TABLE]

The above $u$ is unique if (A5) is fulfilled.

Proof.

The first part of the proof is based on standard compactness arguments, cf. [20, 33]. Let $\{u_{h}\}_{h>0}\subset\mathcal{U}_{h}$ be a sequence of discrete ground states of (15). Then the Rellich embedding theorem together with the uniform boundedness of $\{u_{h}\}$ in $H^{1}(\mathcal{D})$ ensures existence of a weakly convergent subsequence $\{u_{h_{j}}\}\subset\{u_{h}\}$ and a $u_{0}\in H_{0}^{1}(\mathcal{D})$ such that $u_{h_{j}}\rightharpoonup u_{0}$ in $H_{0}^{1}(\mathcal{D})$ and $u_{h_{j}}\rightarrow u_{0}$ strongly in $L^{p}(\mathcal{D})$ for $1\leq p\leq 4$. It remains to prove that $u_{0}$ is a ground state of (5) and that the weak-convergence in $H^{1}_{0}(\mathcal{D})$ is actually a strong convergence. To see this, first note that since the family of finite element spaces $\{V_{h_{j}}\}_{j\in\mathbb{N}}$ (based on regular mesh refinement) is dense in $H^{1}_{0}(\mathcal{D})$, there exists a minimising sequence $\{v_{h_{j}}\}$ for the energy functional $E$, i.e, $E(v_{h_{j}})\rightarrow\inf\limits_{v\in\mathbb{S}}E(v)$ where $v_{h_{j}}\in V_{h_{j}}\cap\mathbb{S}$. Hence, the lower semi-continuity of $E$ yields

[TABLE]

Since $\min\limits_{v\in\mathbb{S}}E(v)>E(u_{0})$ is not possible, we conclude $\lim\limits_{j\rightarrow 0}E(u_{h_{j}})=E(u_{0})=\min\limits_{v\in\mathbb{S}}E(v)$. The normalization of the limit, i.e., $\|u_{0}\|_{L^{2}(\mathcal{D})}=1$, follows from the strong $L^{2}$-convergence of the subsequence $\{u_{h_{j}}\}\subset\mathbb{S}$ to $u_{0}$. This implies $u_{0}$ is a ground state of (5). Consequently, the weak convergence of $\{u_{h_{j}}\}$ in $H_{0}^{1}(\mathcal{D})$ and the strong convergence in $L^{4}(\mathcal{D})$, together with $\lim\limits_{j\rightarrow 0}E(u_{h_{j}})=E(u_{0})$, yield $\lim\limits_{j\rightarrow 0}\|u_{h_{j}}\|_{\mbox{\rm\tiny R}}=\|u_{0}\|_{\mbox{\rm\tiny R}}$. The norm equivalence of $\|\cdot\|_{\mbox{\rm\tiny R}}$ and $\|\cdot\|_{H^{1}(\mathcal{D})}$ together with the Radon-Riesz property imply

[TABLE]

Also the convergence of the eigenvalue subsequence $\{\lambda_{h_{j}}\}$ to $\lambda_{0}$ follows from the strong convergence of $\{u_{h_{j}}\}$ to $u_{0}$ in $H^{1}(\mathcal{D})$ and $L^{p}(\mathcal{D})$ ($p=2,4$).

Next, we prove (26). The statement is easily seen with an argument by contradiction. Assume that (26) is wrong, then there is a $\delta>0$ and a subsequence $\{u_{\tilde{h}_{j}}\}\subset\mathcal{U}_{\tilde{h}_{j}}$ with $\tilde{h}_{j}\rightarrow 0$, but $\|u-u_{\tilde{h}_{j}}\|_{H^{1}(\mathcal{D})}>\delta$ for all ground states $u\in\mathcal{U}$. However, now we could again apply our previous arguments to the sequence $\{u_{\tilde{h}_{j}}\}$ to conclude $H^{1}$-convergence of a subsequence to a ground state. This is a contradiction to $\|u-u_{\tilde{h}_{j}}\|_{H^{1}(\mathcal{D})}>\delta$ for all $u\in\mathcal{U}$.

It remains to check that $u\in\mathcal{U}$ can be, without loss of generality, selected in the same phase as $u_{h}$. For that, let $u\in\mathcal{U}$ be such that (26) holds for a given sufficiently small $\delta$. Then

[TABLE]

Now consider the complex number $\boldsymbol{\alpha}:=\int_{\mathcal{D}}u_{h}\overline{u}\hskip 2.0pt\hskip 2.0pt\mbox{d}x\,,$ and the respective polar representation as $\boldsymbol{\alpha}=r_{h}e^{\mathrm{i}\theta_{h}},$ where $r_{h}>0$ and $\theta_{h}\in[0,2\pi)$. The inequalities in (28) imply that $\theta_{h}$ is in a small neighbourhood of [math], in particular $|\sin\theta_{h}|\,\lesssim\sqrt{\delta}$. With this, we consider the phase shift $e^{\mathrm{i}\theta_{h}}u$ of $u$ by the small angle $\theta_{h}$. The shifted function is still a ground state and obviously in the same phase as $u_{h}$, i.e. $\int_{\mathcal{D}}u_{h}\,\overline{e^{\mathrm{i}\theta_{h}}u}\hskip 2.0pt\hskip 2.0pt\mbox{d}x\in\mathbb{R}_{\geq 0}$. Furthermore, it holds

[TABLE]

where we used that $\theta_{h}$ is a small angle and where $\varepsilon(\delta)=\mathcal{O}(\sqrt{\delta})$. Note that we also exploited the fact that $\|u\|_{H^{1}(\mathcal{D})}$ (for any $u\in\mathcal{U}$) is uniformly bounded by a constant that only depends on the ground state energy and the data parameters. Therefore, for a given $u_{h}$ we have an actual (quasi-isolated) ground state $e^{\mathrm{i}\theta_{h}}u\in\mathcal{U}$ in a sufficiently small $\epsilon(\delta)$-neighbourhood of $u_{h}$ that is also in the same phase as $u_{h}$, i.e., we have

[TABLE]

This completes the proof as $\varepsilon(\delta)\rightarrow 0$ for $\delta\rightarrow 0$. ∎

Remark 5.2**.**

Note that if an actual ground state $u\in\mathcal{U}$ and a discrete ground state $u_{h}\in\mathcal{U}_{h}$ are in the same phase, then

[TABLE]

that is, $u_{h}\in T_{\mathrm{i}u}\mathbb{S}\,$. In this sense, $T_{\mathrm{i}u}\mathbb{S}\cap V_{h}$ becomes a natural space for seeking reasonable approximations of $u$.

With this we are prepared to turn towards the error estimates.

5.2 Abstract error bounds for eigenvalue and energy

In this subsection we derive a priori error estimates for the ground state eigenvalue and the ground state energy by bounding it by the $L^{2}$- and $H^{1}$-errors for the ground state.

We first present a lemma which gives a useful relation associated with a ground state eigenpair $(\lambda,u)$ of (5) and a discrete ground state pair $(\lambda_{h},u_{h})$ of (17). Moreover, we obtain an upper bound of the approximation error. The lemma gains interest if we realize that the estimate (30) helps us to derive the optimal rates for higher order finite element spaces.

Lemma 5.3**.**

Assume (A1)-(A4) holds. Let $(\lambda_{h},u_{h})\in\mathbb{R}_{>0}\times\mathcal{U}_{h}$ denote a discrete ground state pair of (15) and $(\lambda,u)\in\mathbb{R}_{>0}\times\mathcal{U}$ an arbitrary exact ground state pair. Then it holds

[TABLE]

and

[TABLE]

Note that (30) together with $\|u\|_{L^{6}(\mathcal{D})}\lesssim\|u\|_{H^{1}(\mathcal{D})}$ (for $d\leq 3$) directly implies

[TABLE]

which is however not sharp for higher order finite element spaces, whereas (30) allows to derive optimal rates.

Proof of Lemma 5.3.

Recall that the $L^{2}$-inner product only involves the real part. Hence, it is symmetric and we can write $\langle\mathcal{A}_{|u|}(v-u),v-u\rangle=\langle\mathcal{A}_{|u|}v,v\rangle+\langle\mathcal{A}_{|u|}u,u\rangle-2\langle\mathcal{A}_{|u|}u,v\rangle$. Using the definition of $\langle\mathcal{A}_{|u|}v,v\rangle$ and by rearranging the terms we get

[TABLE]

Adding $(\beta|v|^{2}v,v)_{L^{2}(\mathcal{D})}$ to both the sides of the above expression leads to the identity

[TABLE]

where the last equality follows from

[TABLE]

[TABLE]

Thus, by testing in (32) with $v=u_{h}$ and using the fact that $\|u_{h}\|_{L^{2}(\mathcal{D})}=1$, we see that (29) follows immediately. Now using assumptions (A1)-(A4) and the equivalence of the norms $\|\cdot\|_{\mbox{\rm\tiny R}}$ and $\|\cdot\|_{H^{1}(\mathcal{D})}$, we bound the first term on the right hand side of (29) as

[TABLE]

Since $u_{h}$ is uniformly bounded in $H^{1}(\mathcal{D})$ (by the fact that it is a discrete energy minimizer), the second term on the right hand side of (29) can be estimated as

[TABLE]

By combining the estimates (29), (33) and (34), we obtain the desired estimate (30). ∎

It is shown by Proposition 5.1 that for a sufficiently fine mesh $h$ and each discrete ground state $u_{h}$ to (15), there exists a (quasi-isolated) ground state of (5) that is in the same phase and in a small $H^{1}$-neighborhood of $u_{h}$. The following lemma (which also proves a part of Theorem 3.2) estimates the corresponding error in energy and shows that it converges to zero at least with the rate $\|u_{h}-u\|_{H^{1}(\mathcal{D})}^{2}$.

Lemma 5.4**.**

Assume (A1)-(A4), let $(\lambda_{h},u_{h})\in\mathbb{R}_{>0}\times\mathcal{U}_{h}$ denote a discrete ground state and $(\lambda,u)\in\mathbb{R}_{>0}\times\mathcal{U}$ an arbitrary exact ground state pair, then it holds

[TABLE]

Proof.

Using the definition of the energy $E$ and its first Fréchet derivative, we can express the difference of a discrete and the respective actual ground state energy by

[TABLE]

Noting that $\langle\mathcal{A}_{|u|}u_{h},u_{h}\rangle-\langle\mathcal{A}_{|u|}u,u\rangle=\langle\mathcal{A}_{|u|}(u_{h}-u),u_{h}\rangle+\langle\mathcal{A}_{|u|}(u_{h}-u),u\rangle$, the first term on the right hand side of equation (36) can be expressed as

[TABLE]

where the penultimate equality follows from the fact $\|u_{h}\|_{L^{2}(\mathcal{D})}=\|u\|_{L^{2}(\mathcal{D})}=1$ since $u,u_{h}\in\mathbb{S}$. With this, the continuity of $(\mathcal{A}_{|u|}-\lambda\mathcal{I})$ implies

[TABLE]

It remains to estimate the second term in (36). For this, we first note that

[TABLE]

which we can further estimate with the Cauchy-Schwarz inequality, the uniform boundedness of discrete ground states and the Sobolev embedding $H_{0}^{1}(\mathcal{D})\hookrightarrow L^{4}(\mathcal{D})$ for $(d\leq 3)$ by

[TABLE]

By putting (37) and (38) back in (36) we get the desired result. ∎

5.3 Proof of $H^{1}$-quasi optimality

In order to get explicit estimates for the error, we need to study the $H^{1}$-best-approximation properties of functions in $T_{\mathrm{i}u}\mathbb{S}$ by functions in $V_{h}\cap T_{\mathrm{i}u}\mathbb{S}$, as well as functions in $\mathscr{H}_{u}^{\mathbb{C},\perp}$ by functions in $V_{h}\cap\mathscr{H}_{u}^{\mathbb{C},\perp}$. The approximation properties in $T_{\mathrm{i}u}\mathbb{S}$ are directly required to estimate the $H^{1}$-error between $u\in T_{\mathrm{i}u}\mathbb{S}$ and $u_{h}\in V_{h}\cap T_{\mathrm{i}u}\mathbb{S}$, whereas the $L^{2}$-error estimates and the estimates for the eigenvalue error require to study dual problems posed in $\mathscr{H}_{u}^{\mathbb{C},\perp}$.

Lemma 5.5** ($H^{1}$-quasi best-approximations).**

Assume (A1)-(A5), let $h$ be sufficiently small. Then, for all $v\in T_{\mathrm{i}u}\mathbb{S}$, it holds

[TABLE]

and for all $w\in\mathscr{H}_{u}^{\mathbb{C},\perp}$ it holds

[TABLE]

Proof.

Exemplarily we only prove (40). The proof of (39) is analogous with obvious modifications. Now, consider the unconstrained $H^{1}$-projection $P_{H^{1},h}:H^{1}_{0}(\mathcal{D})\rightarrow V_{h}$ given by

[TABLE]

and the auxiliary projection $P_{L^{2},h}^{\perp}:V_{h}\rightarrow V_{h}\cap\mathscr{H}_{u}^{\mathbb{C},\perp}$ that maps a function of $V_{h}$ to the orthogonal complement of $u$ with respect to the complex $L^{2}$-inner product, i.e.,

[TABLE]

We obtain for every $w\in\mathscr{H}_{u}^{\mathbb{C},\perp}$:

[TABLE]

where the constant $C_{h}$ is uniformly bounded for all sufficiently small $h$ since

[TABLE]

The $H^{1}$-optimality of $P_{H^{1},h}$ on $V_{h}$ finishes the proof.

∎

In order to derive the desired $H^{1}$-error estimate (20) as stated in Theorem 3.2, we will make use of the following abstract result which can be e.g. found in [40, Theorem 4].

Theorem 5.6**.**

Let $X$ be a Hilbert space and $X_{h}\subset X$ a finite dimensional space. Further, let $\mathcal{J}:X\rightarrow X^{*}$ be a $C^{1}$-mapping and $u\in X$ such that

[TABLE]

We assume

1. (i)

the bilinear form $\langle\mathcal{J}^{\prime}(u)\,\cdot\,,\cdot\,\rangle$ defined on $X\times X$ is symmetric and $\mathcal{J}^{\prime}(u)$ is an isomorphism from $X$ to $X^{*}$, i.e., there exists a constant $\alpha>0$ such that

[TABLE] 2. (ii)

$\mathcal{J}^{\prime}(u)$* satisfies the discrete inf-sup condition, i.e., there is a constant $\tilde{\alpha}>0$ such that*

[TABLE] 3. (iii)

$\mathcal{J}^{\prime}$* is Lipschitz continuous at $u$, i.e., there exists a constant $\varepsilon_{0}>0$ and $L$ such that for all $v\in X$, $\|u-v\|_{X}\leq\varepsilon_{0}$, we have*

[TABLE]

Then for each sufficiently small $\delta>0$, there exists a $h_{\delta}>0$ such that for all $h\in(0,h_{\delta}]$ there is unique $u_{h}\in X_{h}$ satisfying

[TABLE]

Moreover we have the following error estimate:

[TABLE]

In order to write the GPEVP in a form that fits into the setting of Theorem 5.6, we can follow the idea applied in e.g. [18, 19, 21] to express the eigenvalue problem (10) in compact form based on the operator $\mathcal{J}:\mathbb{R}\times H^{1}_{0}(\mathcal{D})\rightarrow\mathbb{R}\times H^{-1}(\mathcal{D})$ given by

[TABLE]

for $(\sigma,v)$, $(\tau,w)\in\mathbb{R}\times H^{1}_{0}(\mathcal{D})$. The space $\mathbb{R}\times H^{1}_{0}(\mathcal{D})$ is equipped with the natural norm $\|(\sigma,v)\|_{\mathbb{R}\times H^{1}(\mathcal{D})}:=|\sigma|+\|v\|_{H^{1}(\mathcal{D})}$. Hence, the eigenvalue problem (10) can be equivalently expressed as: find $(\lambda,u)\in\mathbb{R}_{>0}\times H^{1}_{0}(\mathcal{D})$ such that

[TABLE]

The Frèchet derivative of $\mathcal{J}$ at $(\lambda,u)$, denoted by $\mathcal{J}^{\prime}(\lambda,u):\mathbb{R}\times H^{1}_{0}(\mathcal{D})\rightarrow\mathbb{R}\times H^{-1}(\mathcal{D})$, can be computed as

[TABLE]

for all directions $(\sigma,v)\in\mathbb{R}\times H^{1}_{0}(\mathcal{D})$ and all $(\tau,w)\in\mathbb{R}\times H^{1}_{0}(\mathcal{D})$. It is easy to see that $\mathcal{J}^{\prime}$ is Lipschitz continuous at $(\lambda,u).$

In order to apply the abstract result from Theorem 5.6, we select $X=\mathbb{R}\times T_{\mathrm{i}u}\mathbb{S}$ based on the tangent space in $\mathrm{i}u$ (noting that $(\lambda,u)\in\mathbb{R}\times T_{\mathrm{i}u}\mathbb{S}$). It remains to establish that $\mathcal{J}^{\prime}(\lambda,u)$ is an isomorphism on $\mathbb{R}\times T_{\mathrm{i}u}\mathbb{S}$ and the corresponding discrete inf-sup stability on $\mathbb{R}\times(V_{h}\cap T_{\mathrm{i}u}\mathbb{S})$. The first part is covered by the following lemma.

Lemma 5.7**.**

Assume (A1)-(A5) and let $(\lambda,u)\in\mathbb{R}_{>0}\times\mathcal{U}$ denote a ground state pair of (10), which hence solves (41). Then $\mathcal{J}^{\prime}(\lambda,u)$ is an isomorphism from $\mathbb{R}\times T_{\mathrm{i}u}\mathbb{S}$ to $\mathbb{R}\times(T_{\mathrm{i}u}\mathbb{S})^{*}$.

Proof.

The lemma is proved if we can show that the equation

[TABLE]

admits a unique solution for all $(\gamma,f)\in\mathbb{R}\times(T_{\mathrm{i}u}\mathbb{S})^{*}$. For that, we consider the saddle-point formulation of (43). With the representation of $\mathcal{J}^{\prime}(\lambda,u)$ in (42) we obtain:

[TABLE]

The proof is completed by showing that (44) is uniquely solvable. This holds (see [12, Theorem 1.1]) iff the variational problem

[TABLE]

has a unique solution for all $f\in(T_{\mathrm{i}u}\mathbb{S})^{*}$ and if there exists a constant $\alpha>0$ (independent of $v,\tau$) such that

[TABLE]

Since the kernel of $(\mathcal{I}u)_{|T_{\mathrm{i}u}\mathbb{S}}$ is $\mathscr{H}_{u}^{\mathbb{C},\perp}$, the unique solvability of problem (45) follows from Proposition 2.4. Furthermore, we obtain the inf-sup stability (46) for arbitrary $\tau\in\mathbb{R}$ with

[TABLE]

where we have used the fact that $u\in\mathbb{S}$ and that all ground states $u$ are uniformly bounded in the $H^{1}$-norm (exploiting energy minimization). Hence, $\mathcal{J}^{\prime}(\lambda,u)$ is an isomorphism from $\mathbb{R}\times T_{\mathrm{i}u}\mathbb{S}$ to $\mathbb{R}\times(T_{\mathrm{i}u}\mathbb{S})^{*}$. This completes the proof. ∎

On sufficiently fine meshes, discrete inf-sup stability can be concluded from the continuous inf-sup stability. It is easily seen that the corresponding result is analogous to the proof of Lemma 5.7 as soon as we have inf-sup stability of $E^{\prime\prime}(u)-\lambda\mathcal{I}$ on $V_{h}\cap\mathscr{H}_{u}^{\mathbb{C},\perp}$. The proof of the latter statement follows standard arguments and is elaborated in the appendix, see Proposition A.1. With this, we have the following lemma.

Lemma 5.8**.**

Assume $h$ to be sufficiently small and that (A1)-(A5) hold. Let $(\lambda,u)\in\mathbb{R}_{>0}\times\mathcal{U}$ be a ground state pair of (10), which hence solves (41). Then there exists a $h$-independent constant $\tilde{\alpha}>0$, such that

[TABLE]

With Lemma 5.7 and Lemma 5.8, we verified that $\mathcal{J}$ fulfills all the assumptions of the abstract Theorem 5.6 with $X=\mathbb{R}\times T_{\mathrm{i}u}\mathbb{S}$ and we can straightforwardly prove the following result.

Theorem 5.9**.**

Assume (A1)-(A5) holds and let $u\in\mathcal{U}$ be a (quasi-isolated) ground state. Then for each sufficiently small $\delta>0$ there is a mesh size $h_{\delta}$ such that for all $h\in(0,h_{\delta}]$ there exist a unique discrete stationary state $(\lambda_{h},u_{h})\in\mathbb{R}_{>0}\times T_{\mathrm{i}u}\mathbb{S}$ with $\langle E^{\prime}(u_{h}),v_{h}\rangle=\lambda_{h}\langle\mathcal{I}u_{h},v_{h}\rangle$ for all $v_{h}\in V_{h}$ and such that

[TABLE]

Moreover, $u_{h}$ is a $H^{1}$-quasi best approximation to u, i.e.,

[TABLE]

Proof.

Let $u\in\mathcal{U}$ be a given ground state with ground state eigenvalue $\lambda$, then it holds $\mathcal{J}(\lambda,u)=0$. Thanks to Lemma 5.7 and Lemma 5.8, we see that the assumptions of Theorem 5.9 are fulfilled for all sufficiently fines meshes. Hence, we get the existence of a unique discrete stationary state $(\lambda_{h},u_{h})\in\mathbb{R}_{>0}\times(V_{h}\cap T_{\mathrm{i}u}\mathbb{S})$ such that

[TABLE]

The proof is finished by recalling property (39) from Lemma 5.5 for $u\in T_{\mathrm{i}u}\mathbb{S}$. ∎

Note that Theorem 5.9 complements Theorem 3.2 in the sense that it states that every exact ground state $u\in\mathcal{U}$ can be (asymptotically) approximated with quasi-optimal accuracy. On the contrary, Theorem 3.2 does not address arbitrary ground states $u$, but it rather makes predictions about arbitrary discrete ground states $u_{h}$ and their individual relation to a specific exact ground state $u$. With this in mind, it remains to formally conclude estimate (20) in Theorem 3.2 from the above setting.

Proof of estimate (20) in Theorem 3.2.

Let the mesh size $h$ be small enough and let $u_{h}\in\mathcal{U}_{h}$ be a given discrete ground state. From Proposition 5.1 we know that in a sufficiently small $\delta$-neighborhood of $u_{h}$ there exists a unique exact ground state $u\in\mathcal{U}$ with $u_{h}\in T_{\mathrm{i}u}\mathbb{S}$. On the contrary, Theorem 5.9 guarantees the existence of unique discrete stationary state $\tilde{u}_{h}\in T_{\mathrm{i}u}\mathbb{S}$ in a sufficiently small neighborhood of $u$ for which (47) holds. Since every discrete ground state is in particular a stationary state, the uniqueness implies $\tilde{u}_{h}=u_{h}$ and the $H^{1}$-quasi best approximation property (47) is also valid for $u_{h}$. This proves (20). ∎

5.4 Proof of $L^{2}$-error estimates

Recalling that we already proved the estimates (18) (cf. Lemma 5.4) and (19) (cf. Lemma 5.3), as well as (20) (cf. Theorem 5.9), it only remains to check the $L^{2}$-error estimate (21) to finish the proof of Theorem 3.2. This will be done in this section.

For brevity, let

[TABLE]

denote the $H^{1}$-projection on $V_{h}\cap\mathscr{H}_{u}^{\mathbb{C},\perp}$ which is defined, for $w\in\mathscr{H}_{u}^{\mathbb{C},\perp}$, by

[TABLE]

Recalling Lemma 5.5, we know that for all sufficiently small mesh sizes $h$ and all $w\in\mathscr{H}_{u}^{\mathbb{C},\perp}$ it holds

[TABLE]

By the properties of the finite element space $V_{h}$, we have in particular

[TABLE]

In the following, we will also make use of an $L^{2}$-orthogonal decomposition of $H^{1}_{0}(\mathcal{D})$ into

[TABLE]

Here we note that the span refers to the real-valued span, i.e., $\mbox{span}\{u\}:=\{\alpha u\,|\,\alpha\in\mathbb{R}\}$ and $\mbox{span}\{\mathrm{i}u\}:=\{\alpha\mathrm{i}u\,|\,\alpha\in\mathbb{R}\}$. With this, we obtain a decomposition of the error $u_{h}-u$ as follows.

Lemma 5.10**.**

Let $u\in\mathcal{U}$ be a ground state and $u_{h}\in\mathcal{U}_{h}$ a corresponding discrete ground state such that $\int_{\mathcal{D}}u_{h}\overline{u}\hskip 2.0pt\mbox{d}x\in\mathbb{R}_{\geq 0}$. Then, there exists $v\in\mathscr{H}_{u}^{\mathbb{C},\perp}$ with $0\leq\|v\|^{2}_{L^{2}(\mathcal{D})}\leq 1$ such that

[TABLE]

In particular, with the constant $c(v):=\sqrt{1-\|v\|^{2}_{L^{2}(\mathcal{D})}}-1$, we can express the error as

[TABLE]

Proof.

From (51), there exist $c_{1},c_{2}\in\mathbb{R}$ and $v\in\mathscr{H}_{u}^{\mathbb{C},\perp}$ such that

[TABLE]

Since $u_{h},u,v\in T_{\mathrm{i}u}\mathbb{S}$, we have

[TABLE]

Since $v\in\mathscr{H}_{u}^{\mathbb{C},\perp}$, we conclude $\int_{\mathcal{D}}v\,\overline{\mathrm{i}u}\hskip 2.0pt\mbox{d}x=0$ and hence, $c_{2}=0$. Next, consider

[TABLE]

Since $c_{1}^{2}\in\mathbb{R}_{\geq 0}$, we conclude

[TABLE]

For the error we obtain $u_{h}-u=c_{1}u+v-u=c(v)u+v\,,$ where $c(v)=\sqrt{1-\|v\|^{2}_{L^{2}(\mathcal{D})}}-1$. Since $|\sqrt{1-x}-1|\leq x$ for all $0\leq x\leq 1$, we conclude $c(v)\leq\|v\|^{2}_{L^{2}(\mathcal{D})}$. Finally, we use $(u_{h}-u,v)_{L^{2}(\mathcal{D})}=(c(v)u,v)_{L^{2}(\mathcal{D})}+(v,v)_{L^{2}(\mathcal{D})}=\|v\|^{2}_{L^{2}(\mathcal{D})}$ to verify

[TABLE]

Combing the estimates to $c(v)\leq\|v\|^{2}_{L^{2}(\mathcal{D})}\leq\|u_{h}-u\|_{L^{2}(\mathcal{D})}^{2}$ yields (53) and completes the proof. ∎

In order to exploit $H^{2}$-regularity of solutions to the dual problem (54) below, we require a corresponding result for the operator $E^{\prime\prime}(u)-\lambda\mathcal{I}$ on $\mathscr{H}_{u}^{\mathbb{C},\perp}$. This is established by the following lemma.

Lemma 5.11** (Regularity of dual solutions).**

Let assumptions (A1)-(A5) hold. Then for any $\phi\in L^{2}(\mathcal{D})$, there exists a unique $w_{\phi}\in H^{2}(\mathcal{D})\cap\mathscr{H}_{u}^{\mathbb{C},\perp}$ such that

[TABLE]

and the following estimate holds true:

[TABLE]

Proof.

Since $E^{\prime\prime}(u)-\lambda\mathcal{I}$ fulfills a continuous inf-sup condition (12) on $\mathscr{H}_{u}^{\mathbb{C},\perp}$ guaranteed by (A5), we immediately conclude that (54) has a unique solution with the stability bound

[TABLE]

It remains to prove that $w_{\phi}\in H^{2}(\mathcal{D})$ with the regularity estimate (55). For this, we use the expression for the second Fréchet derivative given by (13) in combination with the definition of the $(\cdot,\cdot)_{\mbox{\rm\tiny R}}$-inner product in (3) to rewrite the equation (54) as

[TABLE]

where $f_{\phi}:=\phi+\lambda w_{\phi}-\beta|u|^{2}w_{\phi}-2\beta\,\mbox{\rm Re}(u\overline{w_{\phi}})u\in L^{2}(\mathcal{D})$. Hence, $\|w_{\phi}\|_{\mbox{\rm\tiny R}}^{2}\leq\|f_{\phi}\|_{L^{2}(\mathcal{D})}\|w_{\phi}\|_{L^{2}(\mathcal{D})}$. Since $\|\cdot\|_{\mbox{\rm\tiny R}}$ is equivalent to the standard $H^{1}$-norm, we conclude with the definition of $f_{\phi}$ that

[TABLE]

Now, let $v\in H^{1}_{0}(\mathcal{D})$. Using the decomposition (51), we can write $v$ uniquely in the form $v=c_{1}u+c_{2}\mathrm{i}u+w\,,$ where $c_{1}=(u,v)_{L^{2}(\mathcal{D})}\,,\,c_{2}=(\mathrm{i}u,v)_{L^{2}(\mathcal{D})}$ and $w\in\mathscr{H}_{u}^{\mathbb{C},\perp}$. This yields

[TABLE]

where $\tilde{f_{\phi}}:=(w_{\phi},u)_{\mbox{\rm\tiny R}}u+(w_{\phi},\mathrm{i}u)_{\mbox{\rm\tiny R}}\mathrm{i}u+f_{\phi}-(f_{\phi},u)_{L^{2}(\mathcal{D})}u-(f_{\phi},\mathrm{i}u)_{L^{2}(\mathcal{D})}\mathrm{i}u$. As $u,f_{\phi}\in L^{2}(\mathcal{D})$, this directly implies $\tilde{f_{\phi}}\in L^{2}(\mathcal{D})$. Furthermore, by the definition of $(\cdot\,\,\cdot)_{\mbox{\rm\tiny R}}$ (cf. (3)), we have

[TABLE]

As $\tilde{f_{\phi}},g_{\phi}\in L^{2}(\mathcal{D})$ and $v\in H^{1}_{0}(\mathcal{D})$ we conclude that $w_{\phi}\in\mathscr{H}_{u}^{\mathbb{C},\perp}\subset H^{1}_{0}(\mathcal{D})$ solves the Poisson equation $-\Delta w_{\phi}=\tilde{f_{\phi}}-g_{\phi}$ on the convex domain $\mathcal{D}$. Hence, standard elliptic regularity theory yields $w_{\phi}\in H^{2}(\mathcal{D})$ and $\|w_{\phi}\|_{H^{2}(\mathcal{D})}\lesssim\|\tilde{f_{\phi}}\|_{L^{2}(\mathcal{D})}+\|g_{\phi}\|_{L^{2}(\mathcal{D})}$. A direct estimation of these norms together with (58) gives $\|\tilde{f_{\phi}}\|_{L^{2}(\mathcal{D})}+\|g_{\phi}\|_{L^{2}(\mathcal{D})}\lesssim\|\phi\|_{L^{2}(\mathcal{D})}$. This establishes (55) and finishes the proof. ∎

Before we can turn towards the proof of the $L^{2}$-error estimate, we require a final lemma.

Lemma 5.12**.**

Assume (A1)(A5). Let $u\in\mathcal{U}$, $u_{h}\in\mathcal{U}_{h}$ and $\delta\in[0,1]$. Then, for any $v\in H_{0}^{1}(\mathcal{D})$ the following holds

[TABLE]

Proof.

By using the definition of the second Fréchet derivative, we can express

[TABLE]

We can now estimate the right-hand side of(60) by using Hölder’s inequality and the Sobolev embedding $H_{0}^{1}(\mathcal{D})$ $\hookrightarrow$ $L^{3}(\mathcal{D}),L^{6}(\mathcal{D})$ for $(d\leq 3)$ in the following way. To do so, we consider the first term of (60) and obtain

[TABLE]

We can now proceed analogously to bound the remaining terms in (60) as

[TABLE]

This completes the proof. ∎

Proof of estimate (21) in Theorem 3.2.

Noting that the rest of Theorem 3.2 is already proved, we assume that $h$ is sufficiently small and we let $u\in\mathcal{U}$ and $u_{h}\in\mathcal{U}_{h}\cap T_{\mathrm{i}u}\mathbb{S}$ denote a ground state pair such that $u_{h}$ is a $H^{1}$-quasi-best approximation of $u$ in the sense of (20). Using Lemma 5.10 there exists $v\in\mathscr{H}_{u}^{\mathbb{C},\perp}$ and $c(v)\in\mathbb{R}$ such that

[TABLE]

By Lemma 5.11 there exists a unique solution $w_{u_{h}-u}\in\mathscr{H}_{u}^{\mathbb{C},\perp}\cap H^{2}(\mathcal{D})$ such that

[TABLE]

and

[TABLE]

We obtain

[TABLE]

Next, we note that from the continuous and discrete eigenvalue problem, we have

[TABLE]

Here, the last term in (5.4) can be expressed, using $w_{u_{h}-u}\in\mathscr{H}_{u}^{\mathbb{C},\perp}$, as

[TABLE]

and the middle two terms in (5.4) can be written using Taylor expansion as

[TABLE]

for some $\delta\in[0,1]$. Inserting the last two identities into (5.4) and afterwards the result into (5.4) yields

[TABLE]

Next, we will estimate the terms appearing on the right hand side of the above equality (66) separately. The first term is bounded with estimate (53) as

[TABLE]

To estimate the second term in (66), we first use the continuity property of bilinear form $\langle(E^{\prime\prime}(u)-\lambda\mathcal{I})\,\cdot\,,\,\cdot\,\rangle$ and then the stability estimate $\|w_{u_{h}-u}\|_{H^{1}(\mathcal{D})}\lesssim\|u_{h}-u\|_{L^{2}(\mathcal{D})}$ to obtain

[TABLE]

For the third term in (66) we exploit $w_{u_{h}-u}\in H^{2}(\mathcal{D})$ and the corresponding regularity estimate (63) to conclude

[TABLE]

By using Lemma 5.12 with $v=P_{h}w_{u_{h}-u}$, we can estimate the fourth term in (66) as

[TABLE]

where the last step follows from the $H^{1}$-stability of $P_{h}$. Finally, the last term in (66) is again estimated using the approximation property of projection $P_{h}$ together with (63):

[TABLE]

Substituting the estimates (67)–(71) into (66) and using the (suboptimal) eigenvalue estimate (31), i.e., $|\lambda_{h}-\lambda|\lesssim\|u_{h}-u\|_{H^{1}(\mathcal{D})}^{2}+\|u_{h}-u\|_{L^{2}(\mathcal{D})}$, we get

[TABLE]

Finally, we can use the quasi-best-approximation property (20) together with the $H^{2}$-regularity from Lemma 2.5 to conclude $\|u_{h}-u\|_{H^{1}(\mathcal{D})}\lesssim h\,|u|_{H^{2}(\mathcal{D})}$. Plugging this into (5.4) and dividing by $\|u_{h}-u\|_{L^{2}(\mathcal{D})}$ finally yields

[TABLE]

∎

5.5 Proof of optimal convergence rates for the eigenvalue

Theorem 3.2 guarantees that for each discrete ground state pair $(\lambda_{h,k},u_{h,k})\in\mathbb{R}_{>0}\times\mathcal{U}_{h,k}$ to (17) (on a sufficiently fine mesh), there is a unique quasi-isolated ground state pair $(\lambda,u)\in\mathbb{R}_{>0}\times\mathcal{U}$ to (10) such that $u_{h,k}\in T_{\mathrm{i}u}\mathbb{S}$, as well as

[TABLE]

If the ground state $u$ has higher order regularity, i.e., $u\in H^{k+1}(\mathcal{D})$, then we directly get the following optimal order estimates for the $L^{2}$- and $H^{1}$-error (cf. [22, Theorem 3.2.2]):

[TABLE]

On top of that, we obtain from (35) an optimal order estimate for the energy with

[TABLE]

A sharp estimate for the eigenvalue is less straightforward. Here we recall from (30) in Lemma 5.3 that

[TABLE]

which we could crudely estimate (cf. (31)) by

[TABLE]

However, as the optimal rate claimed in Theorem 3.3 is of order $\mathcal{O}(h^{2k})$, the above estimate for the ground state eigenvalue is only optimal for $k=1$, i.e., for $\mathbb{P}^{1}$-FE approximations. In order to improve the result for arbitrary polynomial orders, we require additional regularity for the solution to a suitable dual problem with source term $u|u|^{2}$. In the following, we will elaborate the argument and prove (25) which then concludes Theorem 3.3.

Proof of Theorem 3.3.

As described above, it only remains to check the eigenvalue estimate (25). We start from (75) to observe that we only need to prove that $(u_{h,k}-u,u|u|^{2})_{L^{2}(\mathcal{D})}\lesssim h^{2k}$ to obtain the desired estimate of optimal order for $|\lambda-\lambda_{h,k}|$. Following a similar strategy as for the $L^{2}$-estimate in Theorem 3.2, we consider the unique solution $w_{u|u|^{2}}\in\mathscr{H}_{u}^{\mathbb{C},\perp}$ to the dual problem

[TABLE]

Recall here our regularity assumption that $w_{u|u|^{2}}\in H^{k+1}(\mathcal{D})\cap H^{1}_{0}(\mathcal{D})$ with $\|w_{u|u|^{2}}\|_{H^{k+1}(\mathcal{D})}\lesssim\|\,|u|^{3}\|_{H^{k-1}(\mathcal{D})}$. As in our previous estimates, it holds

[TABLE]

We can now proceed analogously as in the proof of Theorem 3.2 to estimate the terms in (78). For the first term we apply Lemma 5.10 to obtain

[TABLE]

In the same way, the second term in (78) is straightforwardly bounded by

[TABLE]

The approximation properties (49) of $P_{h}$ yield for the third term in (78) that

[TABLE]

For term four we get with Lemma 5.12

[TABLE]

For the last term in (78) we obtain

[TABLE]

Combining everything yields $|\lambda_{h,k}-\lambda|\lesssim h^{2k}$ and finishes the proof. ∎

Acknowledgements. The authors would like to thank the anonymous referees for their very insightful comments that greatly improved the contents of this paper.

Appendix A Appendix: inf-sup stability of $E^{\prime\prime}(u)-\lambda\mathcal{I}$

In the following we prove the inf-sup stability of the bilinear form $E^{\prime\prime}(u)-\lambda\mathcal{I}$ on $\mathscr{H}_{u}^{\mathbb{C},\perp}$ and $V_{h}\cap\mathscr{H}_{u}^{\mathbb{C},\perp}$. The proofs follow standard arguments in this context. We start with Proposition 2.4 which states that quasi-isolation in the sense of (A5) implies inf-sup stability on $\mathscr{H}_{u}^{\mathbb{C},\perp}$.

Proof of Proposition 2.4.

Fix a ground state pair $(\lambda,u)\in\mathbb{R}_{>0}\times H_{0}^{1}(\mathcal{D})$ of (11) and let $\gamma:(-\varepsilon,\varepsilon)\rightarrow\mathbb{S}$ be a smooth curve in the sense of Definition 2.2 on the $L^{2}$-unit sphere $\mathbb{S}$ with $\gamma(0)=u\in\mathcal{U}$ and $\gamma^{\prime}(0)=w\in(T_{u}\mathbb{S}\cap T_{\mathrm{i}u}\mathbb{S})\setminus\{0\}=\mathscr{H}_{u}^{\mathbb{C},\perp}\setminus\{0\}$. We define a Lagrange function for the constrained minimization problem (5) as $L(v):=E(v)-\lambda M(v)$ for $v\in H^{1}_{0}(\mathcal{D})$. A simple calculation gives

[TABLE]

Now since $u$ is a quasi-isolated ground state of the energy functional $E$, we have

[TABLE]

Hence [math] is not a eigenvalue of $E^{\prime\prime}(u)-\lambda\mathcal{I}$ on $\mathscr{H}_{u}^{\mathbb{C},\perp}$. Moreover, for any $w\in\mathscr{H}_{u}^{\mathbb{C},\perp}$, we have

[TABLE]

Thus the injectivity implied by (79) together with the Gårding inequality (80) yield the existence of a constant $\alpha_{1}>0$ such that

[TABLE]

∎

With the inf-sup stability on $\mathscr{H}_{u}^{\mathbb{C},\perp}$, we can conclude the inf-sup stability of $E^{\prime\prime}(u)-\lambda\mathcal{I}$ on $V_{h}\cap\mathscr{H}_{u}^{\mathbb{C},\perp}$ for sufficiently small $h$.

Proposition A.1** (Discrete inf-sup stability).**

Assume (A1)-(A5). Let $u\in\mathcal{U}$ denote a ground state in the sense of (5) and let $\lambda\in\mathbb{R}_{>0}$ denote the corresponding ground state eigenvalue given by (10). If $h$ is sufficiently small, then there is a constant $\alpha_{2}>0$ independent of $h$, such that the discrete inf-sup condition is valid, namely

[TABLE]

Proof.

The proof is standard and follows similarly as in [38] for the Helmholtz equation. Let $w\in V_{h}\cap\mathscr{H}_{u}^{\mathbb{C},\perp}$ be arbitrary but fixed. For $z\in\mathscr{H}_{u}^{\mathbb{C},\perp}$ we consider the expression

[TABLE]

By the continuous inf-sup stability (81) (implied by (A5)) there exists $z\in\mathscr{H}_{u}^{\mathbb{C},\perp}$ with

[TABLE]

and, as before, it holds the regularity estimate $\|z\|_{H^{2}(\mathcal{D})}\lesssim\|w\|_{L^{2}(\mathcal{D})}$. Plugging (83) in (82) we get $\langle(E^{\prime\prime}(u)-\lambda\mathcal{I})(w+z),w\rangle=\langle E^{\prime\prime}(u)w,w\rangle$ and hence with the coercivity of $E^{\prime\prime}(u)$:

[TABLE]

Now, we can choose an approximation $z_{h}\in V_{h}\cap\mathscr{H}_{u}^{\mathbb{C},\perp}$ of $z$ such that

[TABLE]

where we used Lemmas 5.5 and 5.11 here. Therefore, we have

[TABLE]

where $\tilde{\alpha}>0$, if $h$ is sufficiently small. We conclude

[TABLE]

Plugging (87) in (86), we get

[TABLE]

which implies

[TABLE]

Since $w\in\mathscr{H}_{u}^{\mathbb{C},\perp}\cap V_{h}$ was arbitrary, we can take infimum to obtain the desired result. ∎