[ANN] NonlinearSchrodinger.jl: Solving nonlinear Schrodinger equations via higher-order integrators and Darboux Transformations

It is finally my pleasure to announce my very first Julia package, which was written with a ton of feedback from this community, so thank you everyone! I hope I can get some feedback on this package, my very first line of Julia code ever was in writing this, so it’s definitely not as well written as it could be. I still have much to learn from everyone here.



julia> ]add NonlinearSchrodinger


NonlinearSchrodinger.jl is a suite of tools for studying Nonlinear Schrodinger equations. The purpose of this package is to encourage the use of open source software in studying these equations, as most if not all works in this field rely on closed-source codes. This allows for easier reproduction of results, lowers the barrier to entry for new researchers, and alleviates the need to reinvent the wheel.


  1. Solving the cubic Nonlinear Schrodinger equation using 32 algorithms of order up to 8. Symplectic and RKN integrators are available.
  2. Solving the Hirota and Sasa-Satsuma equations using a combined split-step-finite-difference approach using a few different integrators. This is a preliminary feature and is not yet fully supported.
  3. Computing the integrals of motion (energy, momentum, and norm) and their errors.
  4. Computing the Darboux Transformation to study complicated analytical solutions. We currently support the breather and soliton seeds for extended nonlinear Schrodinger equations of order up to 5 (including cubic NLS, Hirota, LPD, Quintic, and arbitrary combinations thereof). We also support the cn and dn seeds for the cubic NLS.
  5. Easy visualization through Plots.jl recipes.
  6. A very simple interface that allows one to compute very complicated solutions via only a few lines of code.
  7. Many utilities for studying maximal intensity breather families on uniform and dnoidal backgrounds, pruning for Nonlinear talbot carpets, and breather to soliton conversion in extended NLSEs.

Supported Equations

The numerical integration algorithms currently only support the cubic Nonlinear Schrodinger equation

i \frac{\partial \psi}{\partial x} + \frac{1}{2} \frac{\partial^2 \psi}{\partial t^2} + |\psi|^2\psi = 0

Note that we use the confusing notation where x and t are flipped compared to the usual notation, we integrate in x and t is the transverse variable. This is in accordance with the literature in the field. We additionally have very preliminary support for the Hirota and Sasa-Satsuma equations, even though this is still under heavy development.

The Darboux Transformation scheme for studying analytical solutions supports the following extended nonlinear Schrodinger equation:

\begin{align} i{\psi _x} + S[\psi (x,t)] - i\alpha H[\psi (x,t)] + \gamma P[\psi (x,t)] - i\delta Q[\psi (x,t)] = 0, \end{align}


\begin{align} \begin{aligned} S[\psi (x,t)] &= \frac{1}{2}{\psi _{tt}} + {\left| \psi \right|^2}\psi, \\ H[\psi (x,t)] &= {\psi _{ttt}} + 6{\left| \psi \right|^2}{\psi _t}, \\ P[\psi (x,t)] &= {\psi _{tttt}} + 8{\left| \psi \right|^2}{\psi _{tt}} + 6{\left| \psi \right|^4}\psi + 4{\left| {{\psi _t}} \right|^2}\psi + 6{\psi _t}^2{\psi ^*} + 2{\psi ^2}\psi _{tt}^*, \\ Q[\psi (x,t)] &= {\psi _{ttttt}} + 10{\left| \psi \right|^2}{\psi _{ttt}} + 30{\left| \psi \right|^4}{\psi _t} + 10\psi {\psi _t}\psi _{tt}^* + 10\psi \psi _t^*{\psi _{tt}} + 20{\psi ^*}{\psi _t}{\psi _{tt}} + 10\psi _t^2\psi _t^*. \end{aligned} \end{align}

Special cases include the cubic nonlinear Schrodinger equation (\alpha = \gamma = \delta = 0), the Hirota equation (\alpha \neq 0, \gamma = \delta = 0), the Lakshmanan-Porsezian-Daniel (LPD) equation, (\gamma \neq 0, \alpha = \delta = 0) and the Quintic nonlinear Schrodinger equation (QNLSE) (\delta \neq 0, \alpha = \gamma = 0). The full equation with \alpha \neq 0, \gamma \neq 0, \delta \neq 0 is also sometimes called the QNLSE.

Supported Algorithms

See the documentation, I don’t want to make this announcement longer than it already is.

Example: 7 Soliton Collision

using NonlinearSchrodinger
using Plots, LaTeXStrings

xᵣ = -10=>10
T = 20
seed = "0"
box = Box(xᵣ, T, Nₓ=500, Nₜ = 500)
λ = [-0.45 + 0.775im, -0.35 + 0.8im, -0.25 + 0.825im, 0.85im, 0.25 + 0.875im, 0.35 + 0.9im, 0.45 + 0.925im]
xₛ = [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0]
tₛ = [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0]

calc = Calc(λ, tₛ, xₛ, seed, box)


Supporting and Citing

The paper will soon be made available on arXiv and submitted for publication.


Well done finishing a package!

I’m working on a package, Superfluids.jl, which also involves the nonlinear Schrodinger/Gross-Pitaevskii equation. It is nowhere near finished, but it works, and I should get round to registering it soon.

At first glance, there is surprisingly little overlap. You’ve implemented some sophisticated algorithms to solve the dynamics in free space. I’ve focused on fields in external potentials, relaxing the initial state, vortices, and the Bogoliubov-de Gennes eigenproblem for sound wave modes. Most of the algorithms are embarassingly inefficient, though the relaxation methods are a bit better, thanks to some help from the Optim maintainers.

We could think about compatibility: it would be nice if you could relax an initial state in one package, then solve it in the other.

If nothing else, a bit of competition is always good for motivation.


Thank you, but it is far from finished, there is still a lot to add! :smiley:

Unfortunately, I am not super familiar with the physics of nonlinear Schrodinger equations outside of nonlinear dynamics and optics, I have mostly studied them on the mathematical physics front. All I know is that the Gross-Pitaevskii equation reduces to a (1+1)D cubic NLSE in some special potentials, in which case all the methods in my package are pretty much applicable (see e.g. this paper: https://journals.jps.jp/doi/10.7566/JPSJ.83.104002, where they use the Lax system).

It seems that you are solving a 2+1 dimensional problem though, am I correct? I unfortunately have very little experience solving that problem numerically, but it’s more than doable using some similar algorithms I believe. But my package is not equipped to handle it.

On the Darboux Transformation front, unfortunately it is only applicable to 1+1 dimensional
Nonlinear Schrodinger problems as they must be completely integrable systems, a condition (2+1)D NLSEs do not satisfy.

I would love to somehow integrate our packages, assuming the physical problems are indeed compatible.

Competition is always a great way to drive more innovation!

Very cool! We also support nonlinear schrodinger equations (although of course with a different focus) over at GitHub - JuliaMolSim/DFTK.jl: Density-functional toolkit, so nice to see that there’s a community forming around these things :slight_smile:

1 Like

Huge fan of DFTK.jl! I actually do DFT and GW many-body perturbation theory for my dissertation work, this nonlinear dynamics stuff is just a hobby for me (I used to do research in this area as an undergrad).

I should have probably clarified in the announcement that the main purpose of this package is to study classically integrable Nonlinear Schrodinger systems in (1+1)D, for e.g. applications in nonlinear optics and dynamics, and not quantum mechanical systems at all (for the most part). So yes, a very different focus :smiley: .

I agree though, glad to see there’s a little community!

Oh, nice, that’s a pretty unusual background you have! MBPT always seemed pretty daunting to me, but we already have the basic building blocks in DFTK so I’m just waiting for a good excuse to get into it and implement it…

1 Like

My aim is for Superfluids to be a thin interface, which makes heavy use of numerical methods implemented by other packages, adds some methods that only physicists use, and wraps them in an API that is convenient for simulations. At the moment, the methods are imported from Optim and DifferentialEquations, but I’m very happy to have NonlinearSchrodinger available too. I’d also be very happy to have a Discretisations package, so I could stop rolling my own finite elements and spectral methods.

We both intend to solve the cubic NLSE in 1+1 dimensions, but we’re generalising from that point in quite different directions. So NLSEBase.jl would only have a few things to standardise at the moment.

  1. What precisely is the NLSE? We’re both using the same factor of 2 in the kinetic energy {1\over 2}{\partial^2\psi\over\partial t^2}. Physicists give the repulsion term a coefficient, C|\psi|^2\psi, because they often want to use length and time units associated with an external potential term V(t)\psi, such as the resonant frequency of a harmonic oscillator. But I’m happy to use C=1 internally and scale the axes as a plotting feature.

  2. How do you associate an array of sampled values with a domain and a discretisation, as in NonlinearSchrodinger.Box? This is an issue in Julia generally, not restricted to the NLSE. There is little prosect of a general answer, and we would have to choose something. In 1D, I guess it’s either equispaced grids or Chebyshev, but I have ambitions for polar coordinates and Fornberg spectral methods in 2D.

  3. How do you specify an initial configuration? In 1D, you presumably want to specify the positions and velocities of solitons, which I think is what the λ, xₛ and tₛ are doing in NonlinearSchrodinger. In 2D, you want to specify the positions and winding numbers of vortices.