Numerical comparison of mass-conservative schemes for the Gross-Pitaevskii equation

Patrick Henning; Johan Wärnegård

doi:10.3934/krm.2019048

Article Contents

2019, Volume 12, Issue 6: 1247-1271. Doi: 10.3934/krm.2019048

This issue Previous Article Memory effects in measure transport equations Next Article A kinetic theory approach to model pedestrian dynamics in bounded domains with obstacles

Numerical comparison of mass-conservative schemes for the Gross-Pitaevskii equation

Patrick Henning and
Johan Wärnegård

Department of Mathematics, KTH Royal Institute of Technology, SE-100 44 Stockholm, Sweden

Received: October 2018

Revised: March 2019

Published: September 2019

The authors acknowledge the support by the Swedish Research Council (grant 2016-03339). P. Henning also acknowledges the support of the Göran Gustafsson Foundation

Abstract / Introduction Full Text(HTML) Figure(12) / Table(6) Related Papers Cited by

Abstract

In this paper we present a numerical comparison of various mass-conservative discretizations for the time-dependent Gross-Pitaevskii equation. We have three main objectives. First, we want to clarify how purely mass-conservative methods perform compared to methods that are additionally energy-conservative or symplectic. Second, we shall compare the accuracy of energy-conservative and symplectic methods among each other. Third, we will investigate if a linearized energy-conserving method suffers from a loss of accuracy compared to an approach which requires to solve a full nonlinear problem in each time-step. In order to obtain a representative comparison, our numerical experiments cover different physically relevant test cases, such as traveling solitons, stationary multi-solitons, Bose-Einstein condensates in an optical lattice and vortex pattern in a rapidly rotating superfluid. We shall also consider a computationally severe test case involving a pseudo Mott insulator. Our space discretization is based on finite elements throughout the paper. We will also give special attention to long time behavior and possible coupling conditions between time-step sizes and mesh sizes. The main observation of this paper is that mass conservation alone will not lead to a competitive method in complex settings. Furthermore, energy-conserving and symplectic methods are both reliable and accurate, yet, the energy-conservative schemes achieve a visibly higher accuracy in our test cases. Finally, the scheme that performs best throughout our experiments is an energy-conserving relaxation scheme with linear time-stepping proposed by C. Besse (SINUM, 42(3):934–952, 2004).

Keywords:

Mathematics Subject Classification: Primary: 35Q55, 65Y20, 65P10, 81Q05; Secondary: 65M60, 65Z99.

Citation:

Full Text(HTML)

Figure 1. $ L^1 $-errors of density and $ H^1 $-errors for all five methods for the single soliton test case

Download: Full-size image PowerPoint slide

Figure 2. We observe remarkably small $ L^2 $-errors for the LCN-FEM when the spatial and temporal discretizations are coupled as $ h/\tau = 2 $. This feature virtually disappears in the error plot of the density, implying that for this ratio the error in phase is minimal

Download: Full-size image PowerPoint slide

Figure 3. Figures (a)-(c) show the solution to the two stationary soliton test case

Download: Full-size image PowerPoint slide

Figure 4. Long time behaviour for Model Problem 3.2

Download: Full-size image PowerPoint slide

Figure 5. Energy evolution for the Two-Step FEM and the LCN-FEM for the time discretizations $ \tau = 2^{-8},2^{-9},\dots,2^{-20} $, note the different time scales. The Two-Step FEM is highly unstable for all ratios of $ h $ and $ \tau $. The LCN-FEM can be highly unstable depending on the ratio of $ h $ and $ \tau $

Download: Full-size image PowerPoint slide

Figure 6. Solution plots of density $ |u_{h,\tau}|^2 $ in Model Problem 4.1 for $ \tau = 2^{-8} $ and $ h = 0.06 $

Download: Full-size image PowerPoint slide

Figure 8. Solution plots of density $ |u_{h,\tau}|^2 $ in Model Problem 4.1 for $ \tau = 2^{-9} $ and $ h = 0.06 $

Download: Full-size image PowerPoint slide

Figure 9. For Model Problem 4.1, the $ L^1 $-errors in density and the $ H^1 $-errors in $ u $ for the four methods versus time-step size for the spatial discretization $ h = 0.06 $

Download: Full-size image PowerPoint slide

Figure 7. Energy evolution of the IM-FEM and RE-FEM for two time-step sizes and the spatial discretization $ h = 0.06 $. The corresponding solutions to Model Problem 4.1 are depicted in Fig.6a-8d

Download: Full-size image PowerPoint slide

Figure 11. Numerical approximations at $ T = 10 $ and mesh resolution $ h = 12\cdot2^{-6} $. The left picture depicts the reference solution. The middle figure shows the solution obtained with the Crank-Nicolson FE-method for $ \tau = T/256 $. We note that the approximations obtained with IM-FEM and RE-FEM for the same mesh size can visually not be distinguished from the CN-FEM result. However, the approximation obtained with LCN-FEM for $ \tau = T/256 $ as depicted in the right figure clearly suffers from an energy pollution and looks distorted compared to the other solutions

Download: Full-size image PowerPoint slide

Figure 10. Both figures depict a comparison between the isolines of the reference solution (black) and the isolines of a RE-FEM approximation (left, red lines) and a IM-FEM approximation (right, red), respectively. The snapshots are at maximum time $ T = 10 $ and the approximations are computed for $ \tau = T/256 $ and $ h = 12\cdot2^{-6} $

Download: Full-size image PowerPoint slide

Figure 12. Density plots, $ |u_{h,\tau}|^2 $, for $ \tau = 0.515\cdot2^{-11} $ and $ h = 0.01 $

Download: Full-size image PowerPoint slide

Table 1. Note that an optimal convergence rate for the error $ e_h = u_{h,\tau}^N - u^N $ in $ L^2 $ is of order $ \mathcal{O}(\tau^2+h^2) $ and an optimal convergence rate in $ H^1 $ is of order $ \mathcal{O}(\tau^2+h) $. As everywhere in the paper, this refers to the case of piecewise linear FEM (simplicial Langrange elements) for the space discretization

	IM-FEM	CN-FEM	RE-FEM	LCN-FEM	Two-Step FEM
Mass-Conservative	Yes	Yes	Yes	Yes	Yes
Energy-Conservative	No	Yes	Yes	No	No
Symplectic	Yes	No	No	No	No
Linear time-steps	No	No	Yes	Yes	Yes
$ L^2 $-convergence rates	Optimal	Optimal	-	Optimal	Optimal
	Conditional	Unconditional		Unconditional	Conditional
$ H^1 $-convergence rates	Optimal	-	-	Optimal	Optimal
	Conditional			Unconditional	Conditional

| Show Table

DownLoad: CSV

Table 4. Average CPU-times over five runs, the number of time-steps is 1024 and $ N $ denotes the degrees of freedom. The computational complexity is the same for the linear methods and, in 1D, roughly a tenth of that of the linear methods

	CPU [s] times 1D
$ N $	IM-FEM	CN-FEM	RE-FEM	LCN-FEM	TwoStep FEM
1600	47.92	46.58	5.78	5.40	5.50
3200	96.59	92.57	11.52	10.75	10.91
6400	197.89	192.53	23.38	21.83	22.18

| Show Table

DownLoad: CSV

Table 2. $ L^1 $-errors of the density for the two stationary soliton test case at time $ T = 2 $

$ \|\|\rho-\rho_{h,\tau}\|\|_{L^1} $
$ \tau $	IM-FEM	CN-FEM	RE-FEM	LCN-FEM	TwoStepFEM	$ \tau $	SP2
$ 2^{-3 } $	2.684e14	3.368e12	12.068	21.033	11.757	$ 2^{-7} $	11.582
$ 2^{-4 } $	23.911	22.546	13.774	18.639	11.0617	$ 2^{-8} $	6.2246
$ 2^{-5 } $	23.325	19.659	7.544	18.145	19.100	$ 2^{-9} $	2.1756
$ 2^{-6 } $	15.885	12.738	4.648	18.525	15.223	$ 2^{-10} $	0.5953
$ 2^{-7 } $	8.498	5.355	1.265	4.868	7.733	$ 2^{-11} $	0.1523
$ 2^{-8 } $	3.356	1.442	0.319	8.109	2.301	$ 2^{-12} $	0.0383
$ 2^{-9 } $	0.959	0.360	0.080	3.537	0.577	$ 2^{-13} $	0.0096
$ 2^{-10 } $	0.252	0.091	0.021	0.723	0.145	$ 2^{-14} $	0.0024

| Show Table

DownLoad: CSV

Table 3. $ H^1 $-errors for the two stationary soliton test case at time $ T = 2 $

$ \|\|u-u_{h,\tau}\|\|_{H^1} $
$ \tau $	IM-FEM	CN-FEM	RE-FEM	LCN-FEM	TwoStepFEM	$ \tau $	SP2
$ 2^{-3 } $	3.795e10	4.270e9	11.445	24.614	10.629	$ 2^{-7} $	9.8755
$ 2^{-4 } $	10.368	11.108	15.194	22.414	11.157	$ 2^{-8} $	4.6047
$ 2^{-5 } $	9.544	11.474	13.620	29.591	10.877	$ 2^{-9} $	1.7293
$ 2^{-6 } $	10.400	9.935	3.768	36.725	10.542	$ 2^{-10} $	0.4877
$ 2^{-7 } $	6.175	4.662	1.022	9.990	6.601	$ 2^{-11} $	0.1258
$ 2^{-8 } $	2.704	1.261	0.260	12.261	2.043	$ 2^{-12} $	0.0317
$ 2^{-9 } $	0.806	0.315	0.066	3.277	0.514	$ 2^{-13} $	0.0079
$ 2^{-10 } $	0.214	0.080	0.018	0.632	0.129	$ 2^{-14} $	0.0020

| Show Table

DownLoad: CSV

Table 5. $ H^1 $-errors and CPU-times for different time-step sizes and mesh size $ h = 24/51200 $

$ \|\|u-u_{h,\tau}\|\|_{H^1} $, CPU [s]
$ \tau $	IM-FEM		CN-FEM		RE-FEM		LCN-FEM		Two-Step FEM
$ 2^{-9} $	2.401	1220	3.045	1130	3.045	130	95.175	130	9.338	140
$ 2^{-10} $	2.474	1960	0.980	1800	0.979	210	38.324	210	9.148	240
$ 2^{-11} $	0.710	3920	0.249	3540	0.249	370	1.948	370	8.942	450
$ 2^{-12} $	0.179	7800	0.062	7270	0.062	700	0.295	680	8.899	870
$ 2^{-13} $	0.044	10750	0.015	9800	0.015	1240	0.080	1190	8.683	1610

| Show Table

DownLoad: CSV

Table 6. $ L^2 $ and $ H^1 $-errors and $ L^1 $-errors of the density for the rotating condensate experiment (i.e. Model Problem 4.2) at time $ T = 10 $ and for the mesh size $ h = 12\cdot2^{-6} $

$ \|\|u-u_{h,\tau}\|\|_{L^2} $, $ \|\|\rho-\rho_{h,\tau}\|\|_{L^1} $, $ \|\|u-u_{h,\tau}\|\|_{H^1} $
$ T/\tau $	IM-FEM			CN-FEM			RE-FEM			LCN-FEM
64	0.6590	0.3560	1.9949	1.7527	0.3594	3.9750	1.7272	0.3553	3.9466	1.4095	1.4014	6.6062
128	1.1946	0.3496	2.1666	1.0690	0.3579	2.8822	1.0779	0.3579	2.9085	1.4500	1.2261	8.1410
256	1.1855	0.1893	3.1790	0.8746	0.1714	2.4329	0.8717	0.1707	2.4257	0.5402	0.1431	1.6532
512	0.3416	0.0515	0.9392	0.2362	0.0453	0.6699	0.2373	0.0456	0.6734	0.6026	0.4097	7.1779
1024	0.0847	0.0126	0.2333	0.0575	0.0110	0.1636	0.0599	0.0116	0.1706	0.0400	0.0118	0.2426
2048	0.0192	0.0027	0.0523	0.0110	0.0023	0.0348	0.0150	0.0029	0.0427	0.0111	0.0013	0.0239

| Show Table

DownLoad: CSV

Related Papers

Cited by

References

[1]	A. Aftalion, Vortices in Bose-Einstein Condensates, Progress in Nonlinear Differential Equations and their Applications, 67. Birkhäuser Boston, Inc., Boston, MA, 2006.
[2]	P. L. Christiansen, M. P. Sorensen and A. C. Scott, Nonlinear Science at the Dawn of the 21st Century, Lecture Notes in Physics, 542. Springer-Verlag, Berlin, 2000. doi: 10.1007/3-540-46629-0.
[3]	G. D. Akrivis, V. A. Dougalis and O. A. Karakashian, On fully discrete Galerkin methods of second-order temporal accuracy for the nonlinear Schrödinger equation, Numer. Math., 59 (1991), 31-53. doi: 10.1007/BF01385769.
[4]	T. Aktosun, T. Busse, F. Demontis and C. van der Mee, Exact solutions to the nonlinear Schrödinger equation, in Topics in Operator Theory, Systems and Mathematical Physics, Oper. Theory Adv. Appl., Birkhäuser Verlag, Basel, 2 (2010), 1–12. doi: 10.1007/978-3-0346-0161-0_1.
[5]	R. Altmann, P. Henning and D. Peterseim, Quantitative Anderson localization of Schrödinger eigenstates under disorder potentials, submitted, 2018, arXiv: 1803.09950.
[6]	P. W. Anderson, Absence of diffusion in certain random lattices, Phys. Rev., 109 (1958), 1492–1505, URL https://link.aps.org/doi/10.1103/PhysRev.109.1492. doi: 10.1142/9789812567154_0007.
[7]	X. Antoine, W. Z. Bao and C. Besse, Computational methods for the dynamics of the nonlinear Schrödinger/Gross-Pitaevskii equations, Comput. Phys. Commun., 184 (2013), 2621-2633. doi: 10.1016/j.cpc.2013.07.012.
[8]	D. N. Arnold, G. David, D. Jerison, S. Mayboroda and M. Filoche, Effective confining potential of quantum states in disordered media, Phys. Rev. Lett., 116 (2016), 056602. doi: 10.1103/PhysRevLett.116.056602.
[9]	W. Z. Bao and Y. Y. Cai, Mathematical theory and numerical methods for Bose-Einstein condensation, Kinet. Relat. Models, 6 (2013), 1-135. doi: 10.3934/krm.2013.6.1.
[10]	W. Z. Bao and Q. Du, Computing the ground state solution of Bose-Einstein condensates by a normalized gradient flow, SIAM J. Sci. Comput., 25 (2004), 1674-1697. doi: 10.1137/S1064827503422956.
[11]	W. Z. Bao, S. Jin and P. A. Markowich, Numerical study of time-splitting spectral discretizations of nonlinear Schrödinger equations in the semiclassical regimes, SIAM J. Sci. Comput., 25 (2003), 27-64. doi: 10.1137/S1064827501393253.
[12]	C. Besse, A relaxation scheme for the nonlinear Schrödinger equation, SIAM J. Numer. Anal., 42 (2004), 934-952. doi: 10.1137/S0036142901396521.
[13]	C. Besse, S. Descombes, G. Dujardin and I. Lacroix-Violet, Energy Preserving Methods for Nonlinear Schrödinger Equations, 2018, arXiv: 1812.04890.
[14]	T. Cazenave, Semilinear Schrödinger Equations, Courant Lecture Notes in Mathematics, 10. New York University, Courant Institute of Mathematical Sciences, New York, American Mathematical Society, Providence, RI, 2003. doi: 10.1090/cln/010.
[15]	Q. S. Chang, E. Jia and W. Sun, Difference schemes for solving the generalized nonlinear Schrödinger equation, J. Comput. Phys., 148 (1999), 397-415. doi: 10.1006/jcph.1998.6120.
[16]	D. Clément, A. F. Varón, J. A. Retter, L. Sanchez-Palencia, A Aspect and P. Bouyer, Experimental study of the transport of coherent interacting matter-waves in a 1d random potential induced by laser speckle, New Journal of Physics, 8 (2006), 165, URL http://stacks.iop.org/1367-2630/8/i=8/a=165.
[17]	L. Erdös, B. Schlein and H.-T. Yau, Derivation of the Gross-Pitaevskii equation for the dynamics of Bose-Einstein condensate, Ann. of Math., 172 (2010), 291-370. doi: 10.4007/annals.2010.172.291.
[18]	D. L. Feder, A. A. Svidzinsky, A. L. Fetter and C. W. Clark, Anomalous modes drive vortex dynamics in confined Bose-Einstein condensates, Phys. Rev. Lett., 86 (2001), 564-567. doi: 10.1103/PhysRevLett.86.564.
[19]	M. Greiner, O. Mandel, T. Esslinger, T. W. Hänsch and I. Bloch, Quantum phase transition from a superfluid to a Mott insulator in a gas of ultracold atoms, Nature, 415 (2002), 39–44, URL http://dx.doi.org/10.1038/415039a.
[20]	D. F. Griffiths, A. R. Mitchell and J. L. Morris, A numerical study of the nonlinear Schrödinger equation, Comput. Methods Appl. Mech. Engrg., 45 (1984), 177-215. doi: 10.1016/0045-7825(84)90156-7.
[21]	H. Hasimoto and H. Ono, Nonlinear modulation of gravity waves, Journal of the Physical Society of Japan, 33 (1972), 805-811. doi: 10.1143/JPSJ.33.805.
[22]	P. Henning and A. Målqvist, The finite element method for the time-dependent Gross-Pitaevskii equation with angular momentum rotation, SIAM J. Numer. Anal., 55 (2017), 923-952. doi: 10.1137/15M1009172.
[23]	P. Henning and D. Peterseim, Crank-Nicolson Galerkin approximations to nonlinear Schrödinger equations with rough potentials, Math. Models Methods Appl. Sci., 27 (2017), 2147-2184. doi: 10.1142/S0218202517500415.
[24]	E. Jarlebring, S. Kvaal and W. Michiels, An inverse iteration method for eigenvalue problems with eigenvector nonlinearities, SIAM J. Sci. Comput., 36 (2014), A1978-A2001. doi: 10.1137/130910014.
[25]	O. Karakashian and C. Makridakis, A space-time finite element method for the nonlinear Schrödinger equation: The continuous Galerkin method, SIAM J. Numer. Anal., 36 (1999), 1779-1807. doi: 10.1137/S0036142997330111.
[26]	C. A. Mülle and D. Delande, Disorder and interference: Localization phenomena, arXiv e-prints.
[27]	J. M. Sanz-Serna, Methods for the numerical solution of the nonlinear Schrödinger equation, Math. Comp., 43 (1984), 21–27, URL http://dx.doi.org/10.2307/2007397. doi: 10.1090/S0025-5718-1984-0744922-X.
[28]	J. M. Sanz-Serna, Runge-Kutta schemes for Hamiltonian systems, BIT, 28 (1988), 877-883. doi: 10.1007/BF01954907.
[29]	J. M. Sanz-Serna and J. G. Verwer, Conservative and nonconservative schemes for the solution of the nonlinear Schrödinger equation, IMA J. Numer. Anal., 6 (1986), 25-42. doi: 10.1093/imanum/6.1.25.
[30]	Y. Tourigny, Optimal H¹ estimates for two time-discrete Galerkin approximations of a nonlinear Schrödinger equation, IMA J. Numer. Anal., 11 (1991), 509-523. doi: 10.1093/imanum/11.4.509.
[31]	J. G. Verwer and J. M. Sanz-Serna, Convergence of method of lines approximations to partial differential equations, Computing, 33 (1984), 297-313. doi: 10.1007/BF02242274.
[32]	J. Wang, A new error analysis of Crank-Nicolson Galerkin FEMs for a generalized nonlinear Schrödinger equation, J. Sci. Comput., 60 (2014), 390-407. doi: 10.1007/s10915-013-9799-4.
[33]	H. C. Yuen and B. M. Lake, Instabilities of waves on deep water, Annual Review of Fluid Mechanics, 12 (1980), 303-334. doi: 10.1146/annurev.fl.12.010180.001511.
[34]	V. E. Zakharov, Stability of periodic waves of finite amplitude on a surface of deep fluid, Journal of Applied Mechanics and Technical Physics, 9 (1968), 190-194. doi: 10.1007/BF00913182.
[35]	V. E. Zakharov and A. B. Shabat, Exact theory of two-dimensional self-focusing and one-dimensional self-modulation of waves in nonlinear media, Ž. Éksper. Teoret. Fiz., 61 (1971), 118-134.
[36]	G. Zhong and J. E. Marsden, Lie-Poisson Hamilton-Jacobi theory and Lie-Poisson integrators, Physics Letters A, 133 (1988), 134-139. doi: 10.1016/0375-9601(88)90773-6.
[37]	G. E. Zouraris, On the convergence of a linear two-step finite element method for the nonlinear Schrödinger equation, M2AN Math. Model. Numer. Anal., 35 (2001), 389-405. doi: 10.1051/m2an:2001121.
[38]	G. E. Zouraris, Error Estimation of the Besse Relaxation Scheme for a Semilinear Heat Equation, 2018, arXiv: 1812.09273.