# Validation of EEG forward modeling approaches in the presence of anisotropy
in the source space

F. Drechsler1, J. Vorwerk2, J. Haueisen3, L. Grasedyck4, and C. H. Wolters2 1
Max Planck Institute for Mathematics in the Sciences, Inselstr. 22, D-04103
Leipzig, Germany 2 Institute for Biomagnetism and Biosignalanalysis,
University of Münster, 48149 Münster, Germany 3 Institute of Biomedical
Engineering and Informatics, TU Ilmenau, POB 100565, D-98684 Ilmenau, Germany
4 Institute for Geometry and Applied Mathematics, RWTH Aachen, Templergraben
55, D-52056 Aachen, Germany<EMAIL_ADDRESS>
###### Abstract

The quality of the inverse approach in electroencephalography (EEG) source
analysis is — among other things — depending on the accuracy of the forward
modeling approach, i.e., the simulation of the electric potential for a known
dipole source in the brain.
Here, we use multilayer sphere modeling scenarios to investigate the
performance of three different finite element method (FEM) based EEG forward
approaches – subtraction, Venant and partial integration – in the presence of
tissue conductivity anisotropy in the source space. In our studies, the effect
of anisotropy on the potential is related to model errors when ignoring
anisotropy and to numerical errors, convergence behavior and computational
speed of the different FEM approaches. Three different source space anisotropy
models that best represent adult, child and premature baby volume conduction
scenarios, are used.
Major findings of the study include (1) source space conductivity anisotropy
has a significant effect on electric potential computation: The effect
increases with increasing anisotropy ratio; (2) with numerical errors far
below anisotropy effects, all three FEM approaches are able to model source
space anisotropy accordingly, with the Venant approach offering the best
compromise between accuracy and computational speed; (3) FE meshes have to be
fine enough in the subdomain between the source and the sensors that capture
its main activity. We conclude that, especially for the analysis of cortical
development, but also for more general applications using EEG source analysis
techniques, source space conductivity anisotropy should be modeled and the FEM
Venant approach is an appropriate method.

###### pacs:

87.19.L, 87.19.R, 87.10.-e, 02.30.Dh

Keywords: EEG, source analysis, forward problem, finite element method,
anisotropy, multilayer sphere models, validation

##  1 Introduction

In EEG and magnetoencephalography (MEG) source analysis, the inverse solution,
i.e., the reconstruction of a current distribution in the brain from non-
invasive measurements of the electric potential at the head surface, depends —
among other things — on the accuracy of the forward problem [9, 22, 29]. For
an appropriate solution of the forward problem it is important to keep both
model and numerical errors as small as possible.

Model errors result, for example, from unrealistic or simplifying assumptions
on head volume conductor geometry and tissue conductivity structure. In this
paper, we will focus on the model error caused by ignoring tissue conductivity
anisotropy in the source space (the brain grey matter compartment). After the
volume conductor model has been determined, numerical errors of the chosen
numerical method result from an insufficient discretization and approximation
of the potential due to, for example, insufficient mesh resolution, bad
(deformed) elements or an inappropriate choice of the basis functions with
regard to the regularity of the solution function.

Different numerical methods have been proposed with the goal to solve the EEG
and MEG forward problem in the presence of brain tissue conductivity
anisotropy, under which are the finite difference method (FDM) [20, 35], the
finite volume method (FVM) [12] and the finite element method (FEM). This
paper focuses on the FEM which allows high accuracy for the numerical solution
of elliptic partial differential equations since it is specifically tailored
to the corresponding variational formulation [19, 8, 65, 16] and since it
allows high flexibility and accuracy in modeling the forward problem in
geometrically complicated inhomogeneous and anisotropic head volume conductors
[5, 2, 56, 10, 30, 38, 59, 47, 60, 55, 18, 43, 42].

Recent studies investigated the influence of brain anisotropy in realistic
head models [24, 62, 63, 21, 18]. It was found that, besides the non-
negligible effects of the so-called remote anisotropy (anisotropy in a certain
distance to the source), conductivity changes close or within the source space
have an especially significant effect on the forward problem [62, 23].
However, to the best of our knowledge, a validation of FEM-based forward
approaches in the presence of anisotropy in the source space has not yet been
presented, and the relationship between model errors due to anisotropies in
the source space and numerical errors remains unclear. The standard in source
analysis is still the use of a homogenized isotropically conducting source
space.

When using the effective medium approach, developed and validated in [3, 54,
36, 58], that relates conductivity tensors and water diffusion tensors in a
linear way, brain conductivity anisotropy can be determined non-invasively by
means of diffusion tensor magnetic resonance imaging (DT-MRI or DTI)
techniques.

Using modern 7T MRI machines, evidence for grey matter anisotropy in adults
was recently found in [25, 11]. In [25, Fig. 4], a fractional anisotropy (FA)
map with overlaid vector-orientation showed a clear normally-oriented (radial)
anisotropy in the cortex. In [26, Figure 39.7], this result is confirmed. It
has thus been measured that the largest diffusion tensor eigenvalue is the one
in radial direction. In [11, Fig. 4], FA values in the range of up to 0.15
were measured for the area of the central sulcus, which would relate to an
anisotropy ratio of about 1.3:1. FA-values were not constant in the region of
interest; higher FA was measured in Broadman area BA4 (primary motor cortex)
and lower in BA2 (sensory cortex). In [48], an anisotropy ratio of up to
1.41:1 (in the head of the caudate nucleus) was measured for adults. This
ratio would relate to an FA value of approximately 0.21. This motivates our
choice of a 1.41:1 (radial:tangential) conductivity anisotropy ratio for the
grey matter compartment in our adult model described in detail later in this
study. To the best of our knowledge, it is not likely that an even more
pronounced source space anisotropy is realistic for the adult case.

However, the human brain undergoes a dramatic maturation during the first
years of life. Both the grey and white matters change rapidly during the first
two years of life as indicated by the changes in T1 contrast of the grey
versus white matter around 6 to 9 months after birth [1]. In the period 26 to
29 weeks of gestation, a rapid dendritic differentiation in the cortical
plate, the appearance of a six-layered pattern and a maximal development of
the subplate zone has been reported [32]. As shown in [34, 33], the grey
matter anisotropy ratio is increasing with decreasing age. Since the
investigation of the cortical development using non-invasive EEG and MEG
techniques is of increasing interest [17, 41, 39, 37, 45, 40, 27], a correct
modeling of especially these higher levels of source space anisotropy might
significantly contribute to an understanding of the observed maturation
effects. In [34] and in [33, Fig. 2], anisotropy ratios have been measured in
the cortex of premature babies of 26 weeks of gestational age, which could
otherwise only been found in white matter tissue. Therefore, for our premature
baby model described in detail in the methods section of this paper, we assume
a source space anisotropy ratio of 5:1 as a representation of an upper limit
of grey matter anisotropy in humans and for our child model a ratio of 2.7:1
111unpublished maximal cortical anisotropy ratio in a DTI study of a child.

The three chosen cortical anisotropy ratios are upper limits of grey matter
anisotropy and surely only rough approximations of the realistic situations,
especially since they are known to vary strongly depending on the exact
location in the grey matter compartment [48, 11]. The chosen ratios have to be
seen in the light of the main interest of this study: We want to investigate
relationships of grey matter anisotropy modeling errors and corresponding
numerical errors. The goal is to gain deeper insight in the effects of
cortical conductivity anisotropy on the EEG forward problem and to present the
requirements with regard to a successful FEM modeling. Only if we can find FEM
parametrizations which lead to numerical errors that are below the model
errors, we can conclude that we model grey matter anisotropy accordingly.
Since numerical errors can only be studied in simplified volume conductors, we
restrict this study to multilayer-sphere models with constant radial and
tangential conductivity values for the different tissue compartments, where
quasi-analytical solutions exist [14].

In the next section, we will shortly introduce the EEG forward problem and
three different FEM approaches for the modeling of a current dipole in
inhomogeneous and anisotropic volume conductors, the subtraction approach
described in detail in [16, 65], the Venant approach introduced in [10] and
described in detail in [61] and the partial integration (PI) approach [59,
55]. In our previous work [28], we compared all three approaches in
tetrahedral volume conductors with isotropic source space, while, for the
first time, we validate in the presence of cortical conductivity anisotropy
here. We then describe our measures for model errors (or cortical anisotropy
effect errors) and numerical errors and detail the setup of the volume
conductor models. We optimize cortical conductivity values so that anisotropy
effects are purely due to anisotropy and not to a suboptimal choice of the
corresponding homogenized isotropic conductivity. Starting from a coarser
mesh, we then generate globally and locally refined tetrahedral models to
examine convergence behavior of the different FEM approaches. In the first two
studies, the refined mesh models are constraint to fulfill the so-called
homogeneity condition of the full subtraction approach, namely the condition
that conductivity is not allowed to change in a small subdomain around the
source (in our setup the tetrahedron of the coarsest mesh that contains the
source). The result section is divided in three studies: The first and the
second study, where globally refined and locally refined meshes are used,
resp., have the goal to present the effect of cortical conductivity anisotropy
in the three different volume conductor models, relate those effects to the
numerical errors and determine the FEM approach that performs best with regard
to numerical error, computational speed and constraints on mesh resolution and
practical applicability. In the third study, the model error is further
reduced by means of relaxing the homogeneity condition and the best-performing
FEM approach of studies 1 and 2 is used for computation. In the last section,
we discuss results and conclude.

##  2 Methods

###  2.1 Forward problem formulation

Assuming the quasi-static approximation of the Maxwell equations, the electric
potential $u$ in the head domain $\Omega$ with conductivity distribution
$\sigma$, evoked by a primary current $\mathbf{j}^{p}$, is given by a Poisson
equation with homogeneous Neumann boundary conditions on the head surface
$\Gamma=\partial\Omega$. This can be expressed by

$\nabla\cdot(\sigma\nabla u)=\nabla\cdot\mathbf{j}^{p}=J^{p}\mbox{ in
}\Omega,\qquad\langle\sigma\nabla u,\mathbf{n}\rangle=0\mbox{ on }\Gamma,$ (1)

with $\mathbf{n}$ being the unit surface normal. To achieve uniqueness of the
solution a reference electrode with given potential has to be fixed, i.e.,
$u(\mathbf{x}_{ref})=0$ [46, 22, 65]. A primary current at position
$\mathbf{y}\in\mathbb{R}^{3}$ and with moment $\mathbf{M}\in\mathbb{R}^{3}$ is
represented by a mathematical dipole

$J^{p}(\mathbf{x})=\nabla\cdot\mathbf{j}^{p}(\mathbf{x})=\nabla\cdot(\mathbf{M}\cdot\delta(\mathbf{x}-\mathbf{y})),$
(2)

where $\delta$ is the Dirac delta distribution [46, 15, 22, 65].

###  2.2 Finite element method

We discretize (1) using a variational formulation and a Galerkin FE approach.
We choose linear basis functions $\varphi_{i}$ at the vertex positions
$\xi_{i}$, $i\in\\{1,\ldots,N\\}$, of a discretization of $\Omega$. The finite
element approach leads to the following linear equation system

$Ku=b^{\mathbf{y}},\quad K\in\mathbb{R}^{N\times N},\quad
u,b^{\mathbf{y}}\in\mathbb{R}^{N}$ (3)

[65, 16] with the stiffness matrix

$K_{ij}:=\int_{\Omega}\langle\sigma(\mathbf{x})\nabla\varphi_{j}(\mathbf{x}),\nabla\varphi_{i}(\mathbf{x})\rangle\>d\mathbf{x},\quad
i,j\in\\{1,\ldots,N\\},$

and the right-hand side

$b^{\mathbf{y}}:=-\int_{\Omega}J^{p}(\mathbf{x})\varphi_{i}(\mathbf{x})\>d\mathbf{x}=-\int_{\Omega}\nabla\cdot\mathbf{j}^{p}(\mathbf{x})\varphi_{i}(\mathbf{x})\>d\mathbf{x}.$
(4)

###  2.3 Dipole models

As described in the introduction, we evaluate three different FE approaches
for the treatment of the singularity on the right hand side (4) with regard to
their performance in the presence of source space anisotropy. In the
following, we will give a short overview of their derivations.

####  2.3.1 Subtraction approach

Under the assumption that there exists a non-empty open neighborhood
$\Omega_{\infty}$ of the source position $\mathbf{y}$ with constant (possibly
anisotropic) conductivity
$\sigma^{\infty,\mathbf{y}}:=\sigma(\mathbf{y})\in\mathbb{R}^{3\times 3}$, we
split up the potential $u$ and the conductivity $\sigma$ into the parts

$\displaystyle u$ $\displaystyle=u^{corr,\mathbf{y}}+u^{\infty,\mathbf{y}},$
(5) $\displaystyle\sigma$
$\displaystyle=\sigma^{corr,\mathbf{y}}+\sigma^{\infty,\mathbf{y}}.$ (6)

Here, $u^{\infty,\mathbf{y}}$, the so-called singularity potential, is the
analytical solution for a dipole in an unbounded homogeneous conductor with
constant (possibly anisotropic) conductivity $\sigma^{\infty,\mathbf{y}}$,
which can be computed as

$u^{\infty,\mathbf{y}}(\mathbf{x}):=\frac{1}{4\pi\sqrt{{\rm
det}\,\sigma(\mathbf{y})}}\frac{\langle\mathbf{M}(\mathbf{y}),\sigma(\mathbf{y})^{-1}(\mathbf{x}-\mathbf{y})\rangle}{\langle\sigma(\mathbf{y})^{-1}(\mathbf{x}-\mathbf{y}),\mathbf{x}-\mathbf{y}\rangle^{3/2}}.$
(7)

$u^{corr,\mathbf{y}}$ is the so-called correction potential that solves the
equation

$\displaystyle\nabla\cdot(\sigma\nabla
u^{corr,\mathbf{y}})=f(\mathbf{x})\mbox{ in }\Omega,$ (8)
$\displaystyle\langle\sigma\nabla
u^{corr,\mathbf{y}},\mathbf{n}\rangle=g(\mathbf{x})\mbox{ on }\Gamma$ (9)

with

$\displaystyle f(\mathbf{x})$
$\displaystyle:=\nabla(\sigma(\mathbf{y})-\sigma(\mathbf{x}))\nabla
u^{\infty,\mathbf{y}}(\mathbf{x})\quad\text{for}\;\mathbf{x}\in\Omega$ (10)
$\displaystyle g(\mathbf{x})$ $\displaystyle:=-\langle\sigma(\mathbf{x})\nabla
u^{\infty,\mathbf{y}}(\mathbf{x}),\mathbf{n}(\mathbf{x})\rangle\quad\text{for}\;\mathbf{x}\in\Gamma$
(11)

and is computed numerically. Applying the FE method, we obtain the right hand
side

$\displaystyle b_{i}^{\mathbf{y}}$
$\displaystyle:=\int_{\Omega}\langle(\sigma(\mathbf{y})-\sigma(\mathbf{x}))\nabla
u^{\infty,\mathbf{y}}(\mathbf{x}),\nabla\varphi_{i}(\mathbf{x})\rangle
d\mathbf{x}$
$\displaystyle-\int_{\Gamma}\varphi_{i}(\mathbf{x})\langle\mathbf{n}(\mathbf{x}),\sigma(\mathbf{y})\nabla
u^{\infty,\mathbf{y}}(\mathbf{x})\rangle d\Gamma,\quad i\in\\{1,\ldots,N\\}$
(12)

and, after solving (3) with right hand side (12) for the correction potential,
we obtain the approximated solution for the total potential

$u(\mathbf{x})=\sum_{i=1}^{N}(u_{i}+u^{\infty,\mathbf{y}}(\xi_{i}))\varphi_{i}(\mathbf{x}).$

We refer to [16, 65] for more detailed mathematics of the subtraction
approach, where a proof for existence and uniqueness and statements about
convergence properties can be found.

####  2.3.2 Partial integration approach

Applying partial integration to (4) gives us

$\displaystyle-\int_{\Omega}(\nabla\cdot\mathbf{j}^{p})\varphi_{i}\>d\mathbf{x}=\int_{\Omega}\langle\mathbf{j}^{p},\nabla\varphi_{i}\rangle\>d\mathbf{x}-\int_{\Gamma}\langle\mathbf{j}^{p},\mathbf{n}\rangle\>d\Gamma.$

We exploit the fact that the current density vanishes on the head surface and
achieve

$-\int_{\Omega}(\nabla\cdot\mathbf{j}^{p})\varphi_{i}\>d\mathbf{x}=\int_{\Omega}\langle\mathbf{j}^{p},\nabla\varphi_{i}\rangle\overset{(\ref{jpdelta})}{=}\langle\mathbf{M},\nabla\varphi_{i}(\mathbf{y})\rangle.$
(13)

This gives us the right hand side for the partial integration approach

$b_{i}^{\mathbf{y}}:=\left\\{\begin{array}[]{ll}\langle\mathbf{M},\nabla\varphi_{i}(\mathbf{y})\rangle&\mbox{if
}i\in\textit{NodesOfEle}(\mathbf{y}),\\\
0&\textit{otherwise}.\end{array}\right.$ (14)

where $\textit{NodesOfEle}(\mathbf{y})$ gives the nodes of the finite element
that contains the dipole.

####  2.3.3 Venant approach

Exploiting Saint Venant’s principle, we approximate a point dipole by a
distribution of electrical monopoles [10, 61]. We place monopoles at the
finite element node closest to the source position and those sharing an edge
with this node. We adjust the charges of the monopoles so that the dipole
moment of this configuration matches the dipole moment of the source as good
as possible, while the sum of charges and higher moments remain zero. To avoid
numerical instabilities as a consequence of huge source loads, an additional
regularization is needed [10, 61]. If we denote the set of neighboring FE
nodes for a dipole at position $\mathbf{y}$ by
$\textit{Neighbours}(\mathbf{y})$ and the strength of the monopole at FE node
$\xi_{i}\in\textit{Neighbours}(\mathbf{y})$ by $q_{i}$, we obtain the right-
hand side

$b_{i}^{\mathbf{y}}:=\left\\{\begin{array}[]{ll}-q_{i}&\mbox{if
}\xi_{i}\in\textit{Neighbours}(\mathbf{y}),\\\
0&\textit{otherwise}.\end{array}\right.$ (15)

###  2.4 Computational efficiency

The two direct approaches Venant and PI, which approximate the dipolar source
locally by a distribution of electrical monopoles, have a sparse right hand
side and allow the calculation of a single solution in a few milliseconds when
applying the transfer matrix approach [64]. In contrast, the subtraction
approach has a much higher computational demand since both the right hand side
computation (12) and also its multiplication to the transfer matrix are
computationally demanding since (12) is a fully populated vector (see [16] for
more details).

###  2.5 Quasi-analytical solution, difference criterion and electrode
configuration

In [14], a quasi-analytical solution was derived for a mathematical dipole in
an anisotropic multilayer sphere model that will be used here to compute
modeling errors (isotropic versus anisotropic) and numerical errors (numerical
versus quasi-analytical). We use the relative error (RE)

$RE(\mathbf{u}^{1},\mathbf{u}^{2}):=\frac{\|\mathbf{u}^{1}-\mathbf{u}^{2}\|_{2}}{\|\mathbf{u}^{2}\|_{2}}$

as a difference measure throughout our study.
$\mathbf{u}^{1},\mathbf{u}^{2}\in\mathbb{R}^{m}$ with $m$ the number of
surface electrodes might, for example, denote the numerical and the quasi-
analytical solutions or the forward solution in anisotropic and isotropic
brain compartment models, resp..

To achieve error measures which are independent of the choice of the sensor
configuration, we distribute the electrodes in a most-regular way over the
outer sphere surface. In this way we generated a configuration with $m=748$
electrodes on the surface of the outer sphere of our four layer sphere model.

###  2.6 Volume conductor models

Three four-layer sphere volume conductor models are used in our study, an
adult model, a child model and a premature baby model.

Table 1: Four layer sphere models used for the isotropic reference scenarios Compartment | Outer radius | Isotropic conductivity (S/m)
---|---|---
Skin | 92mm | 0.43
Skull | 88mm | 0.01
CSF | 80mm | 1.79
Brain (adult) | 76mm | 0.33
Brain (child) | 76mm | 0.51
Brain (premature baby) | 76mm | 0.59
Table 2: Conductivities in radial and tangential directions used for the anisotropic brain scenarios Model | radial (S/m) | tangential (S/m)
---|---|---
Adult | 0.41 | 0.29
Child | 0.783 | 0.29
Premature baby | 1.2 | 0.24

Table 1 presents the radii and the conductivity parameters for the isotropic
reference scenarios of the four spherical layers that approximate the
compartments skin, skull, cerebrospinal fluid (CSF) and brain. While the radii
of all compartments and the isotropic conductivities of the compartments skin
[44], skull [13] and CSF [4] were fixed throughout our studies, the
conductivity of the brain compartment differs between the three models in the
isotropic (Table 1) and anisotropic (Table 2) modeling situations.

For the adult model in Table 1, we used the standard isotropic brain
conductivity of $0.33S/m$ as generally used in source analysis [20, 35, 44].
Furthermore, we fixed the ratios between radial and tangential conductivities
in Table 2 for the anisotropic brain models of 1.41:1 (adult), 2.7:1 (child)
and 5:1 (premature baby) following the motivation given in the introduction.

The values for radial and tangential conductivity for the adult model in Table
2 were computed by means of an optimization procedure. We computed forward
solutions using the quasi-analytical forward method for 16,338 dipole
locations and three dipole moments in the carthesian directions, overall thus
$49,014$ forward solutions, first for the isotropic adult reference scenario
and, in a second step, for a series of corresponding anisotropic scenarios
with fixed ratio of 1.41:1 but differing magnitudes for radial and tangential
conductivities. For each of those anisotropic scenarios, we summed up the RE’s
between isotropic and anisotropic electrode potentials over all $49,014$
sources. The radial and tangential conductivity parameters with the minimal RE
were then chosen as the corresponding anisotropic brain adult model as
indicated in Table 2. We can therefore state that the presented effects in the
results section are anisotropy effects and not due to suboptimally chosen
conductivity values.

We now motivate the chosen optimized conductivity parameters for the child and
premature baby models. The grey matter cytoarchitecture in premature babies is
dominated by the apical dendrites of pyramidal cells that have few branches
and are arrayed radially [34, 33]. The decrease of relative anisotropy in the
grey matter compartment (in our model a decreasing anisotropy ratio from 5:1
to 2.7:1 and down to 1.41:1) was accompanied by a decrease in the average
apparent diffusion coefficient values (in our model an assumed decreasing
radial conductivity in Table 2 and an assumed decreasing homogenized isotropic
conductivity in Table 1) when development progresses [34, 33]. The authors
furthermore report that during maturation pyramidal cells elaborate highly-
branched basal dendrites which tend to be arranged parallel to the cortical
surface and, in addition, thalamocortical afferent fibers are added that are
also oriented parallel to the cortical surface [34, 33]. We therefore assumed
a slightly increasing tangential conductivity between premature baby and child
model (Table 2). Under these constraints, we performed again an optimization
as explained for the adult model to ensure that the chosen homogenized
isotropic conductivities for child and premature baby model in Table 1 have
minimal RE’s to their corresponding anisotropic models shown in Table 2 so
that the presented effects in the results section are again anisotropy effects
and not due to suboptimally chosen corresponding conductivity values.

###  2.7 Tetrahedral mesh generation

Figure 1: Visualization of the basic FE mesh mesh142k (left) and the globally-refined meshes mesh1.1m (middle) and mesh9.5m (right). See also Table 3. Table 3: Globally refined tetrahedral FE meshes: Generation procedure and number of nodes and elements. Mesh | Generated with | Nodes | Elements
---|---|---|---
mesh142k | Tetgen | 142,563 | 900,372
mesh1.1m | Refinement | 1,188,329 | 7,202,976
mesh9.5m | Refinement | 9,590,961 | 57,623,808

For the generation of our coarsest tetrahedral FE mesh, namely mesh142k (Table
3 and Figure 1), we used TetGen [49], which implements a constrained Delaunay
tetrahedralization (CDT) approach [51, 50]. The meshing procedure started with
the preparation of a suitable boundary discretization of the model. For each
of the four layers and for a given triangle edge length, nodes were
distributed in a most-regular way and connected through triangles. This
yielded a valid triangular surface mesh for each of the four layers. Meshes of
different layers were not intersecting each other. The CDT approach was then
used to construct a tetrahedralization conforming to the surface meshes. It
first built a Delaunay tetrahedralization of the vertices of the surface
meshes. It then used a local degeneracy removal algorithm combining vertex
perturbation and vertex insertion to construct a new set of vertices which
included the input set of vertices. In the last step, a fast facet recovery
algorithm was used to construct the CDT [52, 51]. This approach is combined
with two further constraints to the size and shape of the tetrahedra. The
first constraint can be used to restrict the volume of the generated
tetrahedra in a certain compartment, the so-called volume constraint. The
second constraint is important for the generation of quality tetrahedra. If
$R$ denotes the radius of the unique circumsphere of a tetrahedron and $L$ its
shortest edge length, the so-called radius-edge ratio of the tetrahedron can
be defined as $Q=\frac{R}{L}$. The radius-edge ratio can detect almost all
badly-shaped tetrahedra except one type of tetrahedra, so-called slivers. A
sliver is a very flat tetrahedron which has no small edges, but can have
arbitrarily large dihedral angles (close to $\pi$). For this reason, an
additional mesh smoothing and optimization step was used to remove the slivers
and improve the overall mesh quality.

Table 4: Locally refined tetrahedral FE meshes: Generation procedure, refinement area around the fixed source location (in mm) and number of nodes and elements. Mesh | Generated | Ref. area | Nodes | Elements
---|---|---|---|---
| with | around |  |
|  | source |  |
mesh142k | Tetgen | 0mm | 142,563 | 900,372
mesh142k_10 | Refinement | 10mm | 144,137 | 909,747
mesh142k_20 | Refinement | 20mm | 153,325 | 964,809
mesh142k_30 | Refinement | 30mm | 174,468 | 1,093,566
mesh142k_40 | Refinement | 40mm | 211,729 | 1,315,947
mesh142k_50 | Refinement | 50mm | 266,499 | 1,646,015
mesh142k_60 | Refinement | 60mm | 339,627 | 2,087,099
mesh142k_70 | Refinement | 70mm | 430,432 | 2,634,985
mesh142k_100 | Refinement | 100mm | 774,189 | 4,710,009
mesh1.1m | Refinement | globally refined | 1,188,329 | 7,202,976

Based upon mesh142k, we created further globally and locally refined meshes.
Refinement is carried out using regular refinement which splits every
tetrahedron of the mesh in eight new tetrahedra [7, 6]. The specialty of this
approach is that the set of nodes of the coarser meshes are a subset of the
nodes of each finer mesh, needed later for maintaining the homogeneity
condition of the subtraction approach. In this way, two regularly and globally
refined meshes were created (Table 3 and Figure 1). Furthermore, we are
interested in the effect of local in comparison to global mesh refinement.
Therefore, we fixed the source at

$(183.044;169.393;149.969)^{T}$

having a distance of 73.93mm to the midpoint and of 2.07mm to the CSF. We then
refined mesh142k inside a given radius around this position. First, a regular
refinement, also called red refinement, of each tetrahedron inside the radius
is applied. Secondly, in the green refinement, a closure is created so that
the locally refined mesh has no hanging nodes [7, 6] (Table 4). We refine only
_mesh142k_ and denote the locally refined meshes with _mesh142k_x_ whereby x
is the radius of the refinement area around the dipole. For example, if we
refine all tetrahedra that lie in at most 10mm distance to the fixed position,
we denote the corresponding model _mesh142k_10_.

For mesh142k, the conductivity tensors indicated in Tables 1 and 2 were then
computed in the barycenters of the elements. First, with regard to the studies
in the results-sections 3.1 and 3.2, where we want to maintain the homogeneity
condition for the subtraction approach, we assigned the conductivity of the
original tetrahedron to each of the eight new tetrahedra in the refinement,
that is, we accept that the refinement will only improve the FE approximation,
while the model error is not reduced. Secondly, for the study in results-
section 3.3, we calculate adapted conductivity tensors for the refined meshes,
that is, we obtain meshes where also the conductivity tensors of the refined
brain elements change from element to element.

###  2.8 Dipole positions

Dipoles were positioned in the brain compartment between the radii 68mm and
75mm, that is, in a distance of minimally 1mm and maximally 8mm to the CSF
(the next bigger change in conductivity). We chose this intervall, because it
is well-known that numerical errors increase with decreasing distance of the
source to the next bigger conductivity jump (i.e., the CSF) and a mathematical
reasoning for this phenomenon was given in [65]. We thus used on purpose a
source space intervall where we have to be aware of bigger numerical forward
modeling errors [65, 16, 28]. 10 different dipole positions were chosen
randomly in each of the seven (68-69mm,…,74-75mm) 1mm intervals, that is,
altogether 70 different source positions. Sources were placed in the
barycenters of elements of mesh142k. Numerical errors were shown to be minimal
for sources in element barycenters for both subtraction and PI approaches,
while for the Venant approach, best source positions are the FE nodes [57].
Since these errors are on an overall low scale [57], our choice is slightly in
favor of the subtraction and PI approaches. Since normally-oriented dipoles in
the grey matter compartment are the anatomically and physiologically realistic
sources in combination with the used radial conductivity anisotropy, we mainly
only present results for radially-oriented sources. However, for academic
interest, we also present one result figure for tangentially-oriented sources.

###  2.9 Computational platform and used software

All simulations ran on an IBM server with 2 Xeon X5650 processors (2.67 GHz)
using the SimBio software environment [53].

##  3 Results

For studies 3.1 and 3.2, the refined mesh models are constrained to fulfill
the homogeneity condition of the full subtraction approach. The studies have
the goal to present the effect of cortical conductivity anisotropy in the
three different volume conductor models, relate those effects to the numerical
errors and determine the FEM approach that overall performs best. In the last
study 3.3, the model error is further reduced by means of relaxing the
homogeneity condition, and the best-performing FEM approach of studies 3.1 and
3.2 is used for computation.

For the presentation of our results we use boxplots.

Volume conductor model (adult, child or premature baby, see 2.6) and the FE
dipole modeling approach (Venant, PI or subtraction, see 2.3) are denoted in
the title of each figure. The x-axis shows source eccentricity and the y-axis
the used error measure.

###  3.1 Globally refined meshes that fulfill the homogeneity condition

Figure 2: RE for normally-oriented sources and the Venant FE approach for
adult (left column), child (middle column) and premature baby (right column)
model: Boxplots are shown for the effect of cortical anisotropy (cyan), for
the numerical anisotropy effect (green) (see description in 3.1) and for the
numerical errors in the coarsest mesh142k (blue) and the once (mesh1.1m, red)
and twice (mesh9.5m, grey) refined meshes.

Figure 3: RE for normally-oriented (upper row) and tangentially-oriented
(lower row) sources for the _premature baby_ model for the coarsest mesh142k
(left column) and the once (mesh1.1m, middle column) and twice (mesh9.5m,
right column) refined meshes. Boxplots are shown for the effect of cortical
anisotropy (cyan) and for the numerical errors for subtraction (blue), Venant
(red) and partial integration (grey) approach.

Figures 2 and 3 show the RE results for normally-oriented sources, for
different volume conductor models and different dipole modeling approaches.
The figures contain

  * •

the effect of cortical anisotropy (_aniso effect_ , cyan), that is, the
comparison between quasi-analytical forward results in the anisotropic and in
the corresponding isotropic volume conductor model (see 2.6, Tables 1 and 2),

  * •

the numerical anisotropy effect (_num. aniso effect_ , green), defined here as
the results of the Venant approach in mesh9.5m for isotropic conductivities
versus the quasi-analytical results in the corresponding anisotropic model and

  * •

the numerical errors of the three different FE approaches in the coarsest
mesh142k and in the once (mesh1.1m) and twice (mesh9.5m) refined meshes, that
is, the numerical forward modeling results versus the quasi-analytical
results, both for anisotropic conductivity in the cortex.

Figure 2 presents those results for the Venant approach and for the adult
model (left), for the child model (middle) and for the premature baby model
(right). The _aniso effect_ (cyan) is increasing proportionally to the degree
of cortical anisotropy in the model, that is, it increases from the adult (in
intervall 68-69mm, an RE median of 0.2067) to the premature baby model (an RE
median of 1.0821). As expected, the _num. aniso effect_ (green) closely
follows the _aniso effect_ , only on a slightly higher error scale, which
shows that the error is mainly due to ignoring cortical anisotropy, while the
additional numerical error is small. Most importantly, for mesh142k (blue),
mesh1.1m (red) and mesh9.5m (grey) we obtain a clear convergence behaviour,
that is, the numerical errors are still significant for the coarsest mesh142k,
while for the finest mesh mesh9.5m they are in the area of 0.01$\leq$ RE
$\leq$ 0.04 and thus far below the _aniso effect_.

Figure 3 shows a direct comparison of the three FE approaches subtraction
(blue), Venant (red) and partial integration (grey) in mesh142k (left column),
mesh1.1m (middle column) and mesh9.5m (right column) for both normally-
oriented (upper row) as well as tangentially-oriented (lower row) sources in
the premature baby model. All three approaches show a clear convergence
behavior with increasing mesh size. However, it is clearly visible that
although the subtraction is by far the most computationally expensive FE
approach and although we chose source positions in favor of the subtraction
(and the PI) approach, it shows highest numerical errors for all three meshes.
Both Venant and PI approach are at about the same error level with a sligh
advantage of the Venant approach, even if source positions were chosen in
favor of the PI.

###  3.2 Locally refined meshes that fulfill the homogeneity condition

Table 5: Relative error (RE) for all three dipole models and locally refined meshes. | normally-oriented sources | tangentially-oriented sources
---|---|---
Model | sub. | venant | partial | sub. | venant | partial
mesh142k | 0.16899 | 0.07985 | 0.10579 | 0.09028 | 0.04279 | 0.04030
mesh142k_10 | 0.16844 | 0.07979 | 0.09093 | 0.09034 | 0.04286 | 0.04049
mesh142k_20 | 0.17036 | 0.07861 | 0.08995 | 0.09112 | 0.04446 | 0.04191
mesh142k_30 | 0.18633 | 0.07101 | 0.07591 | 0.13560 | 0.03878 | 0.03075
mesh142k_40 | 0.05510 | 0.04554 | 0.04650 | 0.05632 | 0.02629 | 0.02507
mesh142k_50 | 0.06458 | 0.04153 | 0.04306 | 0.03438 | 0.02540 | 0.02354
mesh142k_60 | 0.06310 | 0.04132 | 0.04284 | 0.03377 | 0.02530 | 0.02318
mesh142k_70 | 0.06205 | 0.04154 | 0.04303 | 0.03578 | 0.02523 | 0.02305
mesh142k_100 | 0.06124 | 0.04198 | 0.04348 | 0.03354 | 0.02492 | 0.02278
mesh1.1m | 0.06048 | 0.04202 | 0.04344 | 0.03278 | 0.02492 | 0.02275

In the second study, we focus our interest on the numerical RE errors in the
locally-refined meshes around the fixed source position as described in 2.7
(see, especially Table 4). Table 5 shows the results for all three dipole
models and for normally- and tangentially-oriented sources. While the dipole
has a distance of 18.07mm to the surface of the multisphere model, we observe
a clear error decrease for a radius of 40mm and nearly no more improvement for
further refinement. Nearly the same numerical accuracy could have thus been
achieved with less than one fifth of the FE nodes. Table 5 confirms that the
subtraction approach can not compete against the direct approaches Venant and
PI and, when taking into account that the source position was chosen in the
element barycenter, that is, in favor of the PI approach, we can conclude
again that the Venant approach seems to be the best choice out of the three
tested FE approaches.

###  3.3 Globally refined meshes with adapted conductivities

Figure 4: RE for normally-oriented sources and the Venant FE approach for the
_premature baby_ model and adapted conductivities: Boxplots are shown for the
effect of cortical anisotropy (cyan) and for the numerical errors in the
coarsest mesh142k (blue) and the once (mesh1.1m, red) and twice (mesh9.5m,
grey) refined meshes without adapted conductivities from study 3.1 and for the
once (mesh1.1m reduced, yellow) and twice (mesh9.5m reduced, pink) refined
meshes with adapted conductivities.

For the last study, we will only use the best performing FE approach from
studies 3.1 and 3.2, namely the Venant FE approach, and the model with highest
cortical anisotropy, the premature baby model. The results of the last
simulation study are shown in Figure 4. The presented boxplots for the _aniso
effect_ (cyan) and for _mesh142k_ (blue) as well as for _mesh1.1m_ (red) and
_mesh9.5m_ (grey) are the known numerical accuracy results for the globally-
refined meshes from study 3.1, that is, _mesh1.1m_ and _mesh9.5m_ fulfill the
homogeneity condition. Additionally, we present the results in _reduced
mesh1.1m_ (yellow) and _reduced mesh9.5m_ (pink), where the homogeneity
condition was relaxed, that is, where the conductivities in the refined models
were adapted to the real model. A comparison of the results _mesh1.1m_ versus
_reduced mesh1.1m_ and _mesh9.5m_ versus _reduced mesh9.5m_ shows only a
slight overall RE numerical error decrease for the two adapted models. The
slight improvements in the modeling of the cortical conductivities thus does
not show a significant improvement with regard to the numerical errors.

##  4 Discussion and conclusion

In this paper, we studied the effect of cortical conductivity anisotropy
versus modeling and numerical errors of three different FE approaches for the
EEG forward problem, namely the subtraction, the partial integration (PI) and
the Venant approach. To do so, we used multi-layer sphere volume conductors
and examined maximal cortical conductivity anisotropy for three different
scenarios, reflecting the increasing anisotropy ratios from an adult over a
child to a premature baby modeling situation. We assumed maximal anisotropy
ratios of 1.41:1, 2.7:1 and 5:1, resp., as motivated by the measurements of
[48, 34, 33].

We used a relative error (RE) based optimization procedure to find the best-
approximating conductivity parameters for the isotropy-anisotropy pairs that
were compared. This procedure allows us to state that the presented RE-effects
are due to cortical anisotropy and not due to suboptimally chosen conductivity
values. For other choices of the conductivities, the RE-errors would thus even
be larger.

We showed that cortical conductivity anisotropy has a significant effect on
the forward EEG problem and should thus be modeled. The effect increased with
increasing anisotropy ratio so that the significance of modeling it increases
with decreasing age of the studied subjects. Since the investigation of the
cortical development is of high interest [17, 45], a correct modeling of
source space anisotropy might thus significantly contribute to an
understanding of cortical maturation effects. Realistic measurement of grey
matter conductivity anisotropy can be performed using a combination of newest
MRI technology [25, 11] and the linear effective medium approach model for
relating water diffusion and conductivity [54, 36, 58], so that the
application of the proposed and evaluated methodology to realistic situations
seems straightforward.

In our first two studies, we used globally- and locally-refined meshes that
fulfilled the so-called homogeneity condition. This condition is important for
the subtraction approach, where the mathematical theory suggests that a
sufficiently (sufficient relative to the mesh resolution) large neighborhood
around the source position with constant conductivity is needed for an
appropriate treatment of the source singularity [65, 16]. In isotropically
conducting brain compartments, this can be easily and well fulfilled [65, 16],
but the situation is more complicated in volume conductors with anisotropic
grey matter conductivity, as evaluated for the first time in this paper. We
found that the smaller the neighborhood (note that a smaller neighborhood
offers the possibility of reducing the model error by means of adapting the
conductivities), the larger the numerical errors of the subtraction approach,
resulting in larger overall errors. As a consequence, for the subtraction
approach, we found that it is better not to adapt conductivity tensors in
refined areas around the source. While the subtraction approach is very
sensitive to this condition, it is of less importance for both direct
approaches PI and Venant, so that the performance of those approaches could
even slightly (even though, in our study with the chosen parameters, not
significantly) be improved (see study 3.3). This effect was only small, since
the elements in the grey matter compartment were already relatively small for
_mesh142k_ , so that already the coarsest mesh captured the real conductivity
structure in the volume conductor quite well. More importantly, we showed FE
convergence for all three FE approaches, so that all of them are appropriate
methods for the modeling of source space anisotropy. Relative errors (RE) were
shown to be far below anisotropy effect errors, especially for the highest FE
resolutions. We also investigated the relative difference measure (RDM) and
the magnification factor (MAG) [31] and found that the RE errors decomposed
into both of them with, however, a more distinct MAG component. In our study
with the chosen meshes, we found the Venant approach to be the overall best-
performing method with regard to accuracy, computational complexity and
practicability. We can thus state that cortical anisotropy should be modeled
and, from the results of this study, the Venant FE approach should be used for
it.

In study 3.2, we showed that for accurate FE modeling, mesh-refinement is not
only needed around the source, but between the source and those EEG sensors,
that capture the main source activity. Since we did not observe an error
reduction for larger refinement radii, we can also conclude that in areas far
away from the source (for superficial sources, e.g., deeper areas in the
brain), coarser mesh resolutions seem sufficient to reduce computational
complexity. However, since our FEM approach is completely linear in the number
of nodes [64], even high resolutions are practicable, as shown in this study,
where the finest resolution had nearly 10 Million FE nodes.

Finally, we would like to mention that from our former experience when using
realistic head models [23, 62, 63] and from the results of this study, we
expect that the proposed modeling is not only relevant for EEG, but also for
MEG, and, because of the higher anisotropy ratios, especially for modeling of
newborn MEG [39, 37, 40, 27] and fetal MEG [17, 41].

## Acknowledgment

This work was supported by the Deutsche Forschungsgemeinschaft (WO1425/3-1,
GR3179/3-1, HA2899/14-1).

## References

## References

  * [1] C.R. Almli, M.J. Rivkin, and R.C. McKinstry.  The NIH MRI study of normal brain development (objective-2): newborns, infants, toddlers, and preschoolers.  Neuroimage, 35:308–325, 2007.
  * [2] K.A. Awada, D.R. Jackson, J.T. Williams, D.R. Wilton, S.B. Baumann, and A.C. Papanicolaou.  Computational aspects of finite element modeling in EEG source localization.  IEEE Trans Biomed. Eng., 44(8):736–751, 1997.
  * [3] P.J. Basser, D. Le Bihan, and J. Mattiello.  Estimation of the effective self-diffusion-tensor from the NMR spin echo.  J. Magn. Res. B, 103:247–254, 1994.
  * [4] S.B. Baumann, D.R. Wozny, S.K. Kelly, and F.M. Meno.  The electrical conductivity of human cerebrospinal fluid at body temperature.  IEEE Trans Biomed. Eng., 44(3):220–223, 1997.
  * [5] O. Bertrand, M. Thévenet, and F. Perrin.  3D finite element method in brain electrical activity studies.  In J. Nenonen, H.M. Rajala, and T. Katila, editors, Biomagnetic Localization and 3D Modelling, pages 154–171. Report of the Dep. of Tech.Physics, Helsinki University of Technology, 1991.
  * [6] J. Bey.  Tetrahedral grid refinement.  Computing, 55:355–378, 1995.  10.1007/BF02238487.
  * [7] J. Bey.  Finite-Volumen- und Mehrgitter-Verfahren für elliptische Randwertprobleme.  Teubner Verlag, 1998.
  * [8] D. Braess.  Finite Elements: Theory, Fast Solvers and Applications in Solid Mechanics. Cambridge University Press, 2007.
  * [9] R. Brette and A. Destexhe.  Handbook of Neural Activity Measurement.  Cambridge University Press, 2012.
  * [10] H. Buchner, G. Knoll, M. Fuchs, A. Rienäcker, R. Beckmann, M. Wagner, J. Silny, and J. Pesch.  Inverse localization of electric dipole current sources in finite element models of the human head.  Electroenc. Clin. Neurophysiol., 102:267–278, 1997.
  * [11] J. Cohen-Adad, J.R. Polimeni, K.G. Helmer, T. Benner, J.A. McNab, L.L. Wald, B.R. Rosen, and C. Mainero.  T${}_{2}^{*}$ mapping and B0 orientation-dependence at 7T reveal cyto- and myeloarchitecture organization of the human cortex.  NeuroImage, 60:1006–1014, 2012.
  * [12] M.J.D. Cook and Z.J. Koles.  The effect of tissue anisotropy on the EEG inverse problem.  In Proc. of the 30th Annual Int. Conf. of the IEEE Engineering in Medicine and Biology Society, Vancouver, British Columbia, Canada, August 20-24, pages 4563–4566, 2008.
  * [13] M. Dannhauer, B. Lanfer, C.H. Wolters, and T.R. Knösche.  Modeling of the human skull in EEG source analysis.  Human Brain Mapping, 32(9):1383–1399, 2011.  DOI: 10.1002/hbm.21114, PMID: 20690140.
  * [14] J.C. de Munck and M. Peters.  A fast method to compute the potential in the multi sphere model.  IEEE Trans Biomed. Eng., 40(11):1166–1174, 1993.
  * [15] J.C. de Munck, B.W. van Dijk, and H. Spekreijse.  Mathematical dipoles are adequate to describe realistic generators of human brain activity.  IEEE Trans Biomed. Eng., 35(11):960–966, 1988.
  * [16] F. Drechsler, C.H. Wolters, T. Dierkes, H. Si, and L. Grasedyck.  A full subtraction approach for finite element method based source analysis using constrained delaunay tetrahedralisation.  NeuroImage, 46(4):1055–1065, 2009.  http://dx.doi.org/10.1016/j.neuroimage.2009.02.024, PMID: 19264145.
  * [17] H. Eswaran, J. Wilson, H. Preissl, S. Robinson, J. Vrba, P. Murphy, D. Rose, and C. Lowery.  Magnetoencephalographic recordings of visual evoked brain activity in the human fetus.  Lancet, 360:779–780, 2002.
  * [18] D. Güllmar, J. Haueisen, and J.R. Reichenbach.  Influence of anisotropic electrical conductivity in white matter tissue on the EEG/MEG forward and inverse solution. a high-resolution whole head simulation study.  NeuroImage, 2010.  doi:10.1016/j.neuroimage.2010.02.014.
  * [19] W. Hackbusch.  Elliptic Differential Equations. Theory and Numerical Treatment.  Springer-Verlag, Berlin, 1992.
  * [20] H. Hallez, B. Vanrumste, P. Van Hese, Y. D’Asseler, I. Lemahieu, and R. Van de Walle.  A finite difference method with reciprocity used to incorporate anisotropy in electroencephalogram dipole source localization.  Phys.Med.Biol., 50:3787–3806, 2005.
  * [21] H. Hallez, B. Vanrumste, P. Van Hese, S. Delputte, and I. Lemahieu.  Dipole estimation errors due to differences in modeling anisotropic conductivities in realistic head models for EEG source analysis.  Phys.Med.Biol., 53:1877–1894, 2008.
  * [22] M.S. Hämäläinen, R. Hari, R.Ilmoniemi, J. Knuutila, and O.Lounasmaa.  Magnetoencephalography: theory, instrumentation, and applications to noninvasive studies of the working human brain.  Rev.Mod.Phys., 65:413–497, 1993.
  * [23] J. Haueisen, C. Ramon, H. Brauer, and H. Nowak.  The influence of local conductivity changes on MEG and EEG.  Biomedizinische Technik, 45(7-8):211–214, 2000.
  * [24] J. Haueisen, D.S. Tuch, C. Ramon, P.H. Schimpf, V.J. Wedeen, J.S. George, and J.W. Belliveau.  The influence of brain tissue anisotropy on human EEG and MEG.  NeuroImage, 15:159–166, 2002.
  * [25] R.M. Heidemann, D.A. Porter, A. Anwander, T. Feiweier, K. Heberlein, T. Knösche, and R. Turner.  Diffusion imaging in humans at 7T using readout-segmented EPI and GRAPPA.  Magnetic Resonance in Medicine, 64(1):9–14, 2010.
  * [26] D.K. Jones.  Diffusion MRI: Theory, Methods, and Applications.  Oxford University Press, USA, ISBN-13: 9780195369779, December 2010.
  * [27] S. Lew, M. Hämäläinen, C.H. Wolters, D.D. Silva, M.-S. Choe, E. Grant, and Y. Okada.  Effect of the fontanel on MEG and EEG source analysis using a finite element model of an infant head.  In Human Brain Mapping, Quebec City, Canada, June 26-30, http://www.humanbrainmapping.org/hbm2011/, 2011.
  * [28] S. Lew, C.H. Wolters, T. Dierkes, C. Röer, and R.S. MacLeod.  Accuracy and run-time comparison for different potential approaches and iterative solvers in finite element method based EEG source analysis.  Applied Numerical Mathematics, 59(8):1970–1988, 2009.  NIHMSID 120338, PMCID: PMC2791331.
  * [29] F. Lucka, S. Pursiainen, M. Burger, and C.H. Wolters.  Hierarchical bayesian inference for the EEG inverse problem using realistic FE head models: Depth localization and source separation for focal primary currents.  NeuroImage, 61(4):1364 – 1382, 2012.
  * [30] G. Marin, C. Guerin, S. Baillet, L. Garnero, and Meunier G.  Influence of skull anisotropy for the forward and inverse problem in EEG: simulation studies using the FEM on realistic head models.  Human Brain Mapping, 6:250–269, 1998.
  * [31] J.W.H Meijs, O.W. Weier, M.J. Peters, and A. van Oosterom.  On the numerical accuracy of the boundary element method.  IEEE Trans Biomed. Eng., 36:1038–1049, 1989.
  * [32] L. Mrzljak, H.B. Uylings, I. Kostovic, and C.G. Van Eden.  Prenatal development of neurons in the human prefrontal cortex: I. a qualitative golgi study.  J Comp Neurol, 271(3):355–386, 1988.
  * [33] J. Neil, J. Miller, P. Mukherjee, and P. S. Hüppi.  Diffusion tensor imaging of normal and injured developing human brain \- a technical review.  NMR in Biomedicine, 15:543–552, 2002.
  * [34] J.J. Neil, R.C. McKinstry, B.L. Schlaggar, G. Schefft, S.I. Shiran, J.S. Shimony, A.Z. Snyder, C.R. Almli, E. Akbudak, and T.E. Conturo.  Evaluation of diffusion anisotropy during human cortical grey matter development.  In Proceedings of ISMRM, 8th Annual Meeting and Exhibition, page 591, 2000.
  * [35] L.A. Neilson, M. Kovalyovc, and Z.J. Koles.  A computationally efficient method for accurately solving the EEG forward problem in a finely discretized head model.  Clin.Neurophysiol., 116:2302–2314, 2005.
  * [36] S.H. Oh, S.Y. Lee, M.H. Cho, T.-S. Kim, and I.H. Kim.  Electrical conductivity estimation from diffusion tensor and t2: a silk yarn phantom study.  In Proc.Intl.Soc.Mag.Reson.Med., volume 14, page 3034, 2006.
  * [37] Y. Okada, K. Pratt, C. Atwood, A. Mascarenas, R. Reineman, J. Nurminen, and D. Paulson.  Babysquid: a mobile, high-resolution multichannel MEG system for neonatal brain assessment.  Rev. Sci. Instrum., 77:24301–24309, 2006.
  * [38] J. Ollikainen, M. Vauhkonen, P.A. Karjalainen, and J.P. Kaipio.  Effects of local skull inhomogeneities on EEG source estimation.  Med.Eng.Phys., 21:143–154, 1999.
  * [39] E. Pihko, L. Lauronen, H. Wikstrom, S. Taulu, J. Nurminen, S. Kivitie-Kallio, and Y. Okada.  Somatosensory evoked potentials and magnetic fields elicited by tactile stimulation of the hand during active and quiet sleep in newborns.  Clin Neurophysiol, 115:448–455, 2004.
  * [40] E. Pihko, P. Nevalainen, J. Stephen, Y. Okada, and L. Lauronen.  Maturation of somatosensory cortical processing from birth to adulthood revealed by magnetoencephalography.  Clin Neurophysiol, 120:1552–1561, 2009.
  * [41] H. Preissl, C.L. Lowery, and H. Eswaran.  Fetal magnetoencephalography: current progress and trends.  Exp Neurol., 190(Suppl 1):S28–36, 2004.
  * [42] S. Pursiainen.  Raviart-thomas-type sources adapted to applied EEG and MEG: implementation and results.  Inverse Problems, 28:16pp., 2012.
  * [43] S. Pursiainen, A. Sorrentino, C. Campi, and M. Piana.  Forward simulation and inverse dipole localization with the lowest order raviart-thomas elements for electroencephalography.  Inverse Problems, 27(4), 2011.
  * [44] C. Ramon, P. Schimpf, and J. Haueisen.  Influence of head models on EEG simulations and inverse source localizations.  BioMedical Engineering OnLine, 5(10), 2006.
  * [45] N. Roche-Labarbe, A.Aarabi, G. Kongolo, C. Gondry-Jouet, M Dümpelmann, R. Grebe, and F. Wallois.  High-resolution electroencephalography and source localization in neonates.  Human Brain Mapping, 29:167–176, 2008.
  * [46] J. Sarvas.  Basic mathematical and electromagnetic concepts of the biomagnetic inverse problem.  Phys.Med.Biol., 32(1):11–22, 1987.
  * [47] P.H. Schimpf, C.R. Ramon, and J. Haueisen.  Dipole models for the EEG and MEG.  IEEE Trans Biomed. Eng., 49(5):409–418, 2002.
  * [48] J. S. Shimony, R.C. McKinstry, E. Akbudak, J.A. Aronovitz, A.Z. Snyder, N.F. Lori, T.S. Cull, and T.E. Conturo.  Quantitative diffusion-tensor anisotropy brain MR imaging: Normative human data and anatomic analysis.  Radiology, 212(3):770–784, 1999.
  * [49] H. Si.  TetGen, a quality tetrahedral mesh generator and three-dimensional delaunay triangulator, v1.3, user’s manual.  Technical Report 9, Weierstrass Institute for Applied Analysis and Stochastics, 2004.  http://tetgen.berlios.de.
  * [50] H. Si.  TetView.  http://tetgen.berlios.de/tetview.html, 2007.
  * [51] H. Si.  Adaptive tetrahedral mesh generation by constrained Delaunay refinement.  International Journal for Numerical Methods in Engineering, 75(7):856–880, 2008.
  * [52] H. Si and K. Gärtner.  Meshing piecewise linear complexes by constrained Delaunay tetrahedralizations.  In Proc. 14th International Meshing Roundtable, pages 147–163. Sandia National Laboratories, 2005.
  * [53] SimBio Development Group.  SimBio: A generic environment for bio-numerical simulations.  online, https://www.mrt.uni-jena.de/simbio, accessed June 15, 2012.
  * [54] D.S. Tuch, V.J. Wedeen, A.M. Dale, J.S. George, and J.W. Belliveau.  Conductivity tensor mapping of the human brain using diffusion tensor MRI.  Proc.Natl.Acad.Sci.USA, 98(20):11697–11701, 2001.
  * [55] S. Vallaghé and T. Papadopoulo.  A trilinear immersed finite element method for solving the electroencephalography forward problem.  SIAM J. Sci. Comput., 32(4):2379–2394, 2010.
  * [56] S.P. van den Broek.  Volume Conduction Effects in EEG and MEG.  PhD thesis, Proefschrift Universiteit Twente Enschede, The Netherlands, 1997.
  * [57] J. Vorwerk.  Comparison of numerical approaches to the EEG forward problem.  Diploma thesis in Mathematics, Institut für Biomagnetismus und Biosignalanalyse, Westfälische Wilhelms-Universität Münster, april 2011.
  * [58] K. Wang, J. Li, S. Zhu, B. Mueller, K. Lim, Z. Liu, and B. He.  A new method to derive white matter conductivity from diffusion tensor MRI.  IEEE Trans. Biomed. Eng., 55(10):2481–2486, 2008.
  * [59] D. Weinstein, L. Zhukov, and C. Johnson.  Lead-field bases for electroencephalography source imaging.  Annals of Biomed.Eng., 28(9):1059–1066, 2000.
  * [60] K. Wendel, N.G. Narra, M. Hannula, P. Kauppinen, and J. Malmivuo.  The influence of CSF on EEG sensitivity distributions of multilayered head models.  IEEE Trans. Biomed. Eng., 55(4):1454–1456, 2008.
  * [61] C.H. Wolters, A. Anwander, G. Berti, and U. Hartmann.  Geometry-adapted hexahedral meshes improve accuracy of finite element method based EEG source analysis.  IEEE Trans.Biomed.Eng., 54(8):1446–1453, 2007.
  * [62] C.H. Wolters, A. Anwander, X. Tricoche, S. Lew, and C.R. Johnson.  Influence of local and remote white matter conductivity anisotropy for a thalamic source on EEG/MEG field and return current computation.  Int.Journal of Bioelectromagnetism, 7(1):203–206, 2005.
  * [63] C.H. Wolters, A. Anwander, D. Weinstein, M. Koch, X. Tricoche, and R.S. MacLeod.  Influence of tissue conductivity anisotropy on EEG/MEG field and return current computation in a realistic head model: A simulation and visualization study using high-resolution finite element modeling.  NeuroImage, 30(3):813–826, 2006.
  * [64] C.H. Wolters, L. Grasedyck, and W. Hackbusch.  Efficient computation of lead field bases and influence matrix for the FEM-based EEG and MEG inverse problem.  Inverse Problems, 20(4):1099–1116, 2004.
  * [65] C.H. Wolters, H. Köstler, C. Möller, J. Härtlein, L. Grasedyck, and W. Hackbusch.  Numerical mathematics of the subtraction method for the modeling of a current dipole in EEG source reconstruction using finite element head models.  SIAM J. on Scientific Computing, 30(1):24–45, 2007.


# Turbulent Threshold for Continuum Calogero-Moser Models

James Hogan  Department of Mathematics, University of California, Los Angeles,
CA 90095, USA<EMAIL_ADDRESS>and  Matthew Kowalski  Department of
Mathematics, University of California, Los Angeles, CA 90095, USA
<EMAIL_ADDRESS>
###### Abstract.

We determine the sharp mass threshold for Sobolev norm growth for the focusing
continuum Calogero–Moser model. It is known that below the mass of $2\pi$,
solutions to this completely integrable model enjoy uniform-in-time $H^{s}$
bounds for all $s\geq 0$. In contrast, we show that for arbitrarily small
$\varepsilon>0$ there exists initial data $u_{0}\in H^{\infty}_{+}$ of mass
$2\pi+\varepsilon$ such that the corresponding maximal lifespan solution
$u:(T_{-},T_{+})\times\mathbb{R}\to\mathbb{C}$ satisfies $\lim_{t\to
T_{\pm}}\|u(t)\|_{H^{s}}=\infty$ for all $s>0$. As part of our proof, we
demonstrate an orbital stability statement for the soliton and a dispersive
decay bound for solutions with suitable initial data.

##  1\. Introduction

This paper is devoted to the focusing continuum Calogero–Moser model

(CCM) $iu_{t}+u_{xx}-2iu\Pi^{+}\partial_{x}(|u|^{2})=0,$

known in some sources as the Calogero–Moser derivative nonlinear Schrödinger
equation. Here $\Pi^{+}$ is the Cauchy–Szegő projection; see (1.7). Solutions
are complex-valued functions of spacetime that belong to a Hardy–Sobolev space

(1.1) $H^{s}_{+}(\mathbb{R}):=\\{f\in
H^{s}(\mathbb{R}):\widehat{f}(\xi)=0\text{ for }\xi<0\\},\quad s\geq 0,$

for each time $t$.

A defocusing version of (CCM) appeared in [16] as a model for internal waves
in a stratified fluid, with the nonlinear Schrödinger equation emerging as a
small-depth limit. The focusing version studied here is derived in [1] as a
formal hydrodynamic limit of the Calogero-Moser system — a Hamiltonian system
of particles on a line interacting with an inverse square potential:

(CM) $V(x_{j})=\sum_{k\neq j}\frac{1}{(x_{j}-x_{k})^{2}}.$

This particle system, and its generalizations, are well-studied; it is solved
in the quantum mechanical case [4, 5] and is integrable in the classical case
[15]. The discrete periodic case was studied by Sutherland [18, 19] and the
periodic continuum limit has been studied under the name of
Calogero–Sutherland derivative NLS in [2, 3].

As a continuum limit of an integrable system, it is unsurprising that (CCM) is
completely integrable. It has a Lax pair given by

$\mathcal{L}_{u}=-i\partial_{x}-u\Pi^{+}\overline{u}\quad\text{and}\quad\mathcal{P}_{u}=i\partial_{x}^{2}+2u\Pi^{+}\partial_{x}\overline{u}$

where $u,\overline{u}$ are understood as multiplication by $u,\overline{u}$,
respectively; see [10]. This Lax pair generates an infinite hierarchy of
conserved quantities [10]

(1.2) $I_{k}(u):=\langle u,\mathcal{L}_{u}^{k}u\rangle.$

The first and most relevant of these are the mass

$M(u):=I_{0}(u)=\int|u|^{2},$

the momentum

$P(u):=I_{1}(u)=\int-i\overline{u}u_{x}-\tfrac{1}{2}|u|^{4},$

and the Hamiltonian

$E(u):=\tfrac{1}{2}I_{2}(u)=\tfrac{1}{2}\int\big{|}u_{x}-iu\Pi^{+}(|u|^{2})\big{|}^{2}.$

This Hamiltonian generates the flow of (CCM), albeit with a nonstandard
Poisson structure arising from a gauge transformation in [1].

Given the infinite family of conserved quantities (1.2), one might expect that
solutions to (CCM) satisfy a priori bounds in $H^{s}$ for all $s\geq 0$.
Indeed, for initial data $u_{0}\in H^{\infty}_{+}$ with $M(u_{0})<2\pi$,
Gérard and Lenzmann showed in [10] that the corresponding maximal lifespan
solution $u(t)$ satisfies

(1.3) $\|u\|_{L^{\infty}_{t}H_{x}^{n}}\leq C(I_{0}(u),I_{1}(u),...,I_{2n}(u))$

for all integers $n\geq 0$. Using (1.3) they showed that for initial data with
mass less than $2\pi$, (CCM) is globally well-posed in $H^{n}_{+}$ for all
integers $n\geq 1$.

On the other hand, Gérard and Lenzmann also showed that (CCM) has multisoliton
solutions for which these a priori bounds fail. Specifically, for $N\geq 2$,
the $N$-soliton solutions, which have mass $2\pi N$, exhibit a growth of
Sobolev norms of the form $\|u(t)\|_{H^{s}}\sim|t|^{2s}$ for large $|t|$. This
raised the question of what is the true threshold for the growth of Sobolev
norms. We answer this question by proving the following:

######  Theorem 1.1 (Turbulent behaviour).

For any sufficiently small $\varepsilon>0$, there exist initial data $u_{0}\in
H^{\infty}_{+}(\mathbb{R})$ with $M(u_{0})=2\pi+\varepsilon$, a time
$T\in(0,\infty]$, and a maximal lifespan solution $u(t)$ to (CCM) such that
for all $s>0$, $u\in C_{t}H^{s}_{x}([0,T)\times\mathbb{R})$ and

$\lim_{t\nearrow T}\|u(t)\|_{H^{s}}=+\infty.$

In particular, if $T=\infty$, then we have the bounds

$\|u(t)\|_{H^{s}}\gtrsim t^{s}.$

An analogous statement holds backwards in time.

This theorem shows that the global-in-time $H^{s}$ bounds (1.3) fail above the
mass of $2\pi$. This presents a significant obstacle to existing techniques
for proving global well-posedness of (CCM) above the mass of $2\pi$.

The $2\pi$ mass threshold appears in this analysis through the sharp
inequality

(1.4) $\|\Pi^{+}\overline{u}f\|_{L^{2}}^{2}\leq\tfrac{1}{2\pi}M(u)\langle
f,-i\partial_{x}f\rangle$

valid for $u\in L^{2}_{+}$ and $f\in H^{\frac{1}{2}}_{+}$; see [10]. This
inequality is saturated by the function

$Q(x)=\frac{\sqrt{2}}{x+i}$

which is the unique nonzero optimizer up to phase, translation, and
$L^{2}$-scaling symmetries. It is known as the soliton solution to (CCM) and
does not evolve in time. The soliton satisfies $M(Q)=2\pi$ and $I_{k}(Q)=0$
for all $k\geq 1$.

Within our proof, this mass threshold also occurs prominently in the existence
of eigenvalues for the Lax operator. In [10], it was shown that

(1.5) $2\pi N\leq M(u)$

where $N$ is the number of eigenvalues of $\mathcal{L}_{u}$. Our proof
requires initial data $u_{0}$ for which the Lax operator $\mathcal{L}_{u_{0}}$
has a negative eigenvalue, which is consistent with our solutions having mass
at least $2\pi$. The soliton $Q$ and the multisolitons are introduced in [10]
as optimizers of (1.5).

Let us briefly review the available well-posedness results on this model.
Local well-posedness in $H^{s}_{+}(\mathbb{R})$ for $s>\tfrac{1}{2}$ was
established by Gérard and Lenzmann in [10]. This was shown for $s>\frac{3}{2}$
via a Kato-type iterative scheme and standard energy methods and extended to
$s>\frac{1}{2}$ following an argument in [7]. For $M(u)<2\pi$, Gérard and
Lenzmann used (1.4) to control $\|u\|_{H^{\frac{1}{2}}}$ in terms of $P(u)$
and $M(u)$. They obtained the bounds (1.3) on $H^{s}$ for all $s\geq 0$
iteratively using the conserved quantities $I_{k}(u)$ defined in (1.2). This
extended local well-posedness in $H^{s}_{+}$ to global well-posedness in
$H^{n}_{+}$ for integers $n\geq 1$, so long as $M(u)<2\pi$. They also showed
that $H^{s}_{+}$ solutions for $s\geq 1$ to (CCM) with mass exactly $2\pi$
cannot blow up in finite time, which extended global well-posedness to mass at
most $2\pi$ for such regular solutions.

Recent work by Laurens, Killip, and Vişan [12] extended global well-posedness
below the mass of $2\pi$ to $H^{s}_{+}(\mathbb{R})$ for $s\geq 0$. Central to
their analysis was the concept of $L^{2}$-equicontinuity:

######  Definition 1.2 ($L^{2}$ Equicontinuity).

A bounded set $U\subset L^{2}(\mathbb{R})$ is said to be equicontinuous in the
$L^{2}(\mathbb{R})$ topology if

$\limsup_{\delta\to 0}\sup_{\begin{subarray}{c}u\in
U,\;|y|<\delta\end{subarray}}\|u(\cdot+y)-u(\cdot)\|_{L^{2}}=0,$

or equivalently,

$\limsup_{\kappa\to\infty}\,\sup_{u\in
U}\int_{|\xi|\geq\kappa}|\widehat{u}(\xi)|^{2}d\xi=0.$

The argument in [12] relied on showing that any $L^{2}$-equicontinuous family
$U\subset H^{\infty}_{+}$ of initial data leads to a family of partial orbits

$U_{T}^{*}=\\{u(t):u\text{ solves \eqref{CM-DNLS}},u(0)\in U,t\in[0,T)\\}$

which is $L^{2}$-equicontinuous for any choice of $T>0$. They defined the
sharp mass threshold $M^{*}$ such that if $U\subset H^{\infty}_{+}$ is
$L^{2}$-equicontinuous and $\sup_{u_{0}\in U}M(u_{0})<M^{*}$, then $U_{T}^{*}$
is also $L^{2}$-equicontinuous for all $T>0$. They showed that $M^{*}\geq
2\pi$ and that (CCM) is globally well-posed in $H^{s}_{+}$ for all $s\geq 0$
below the mass $M^{*}$. The exact value of $M^{*}$ is currently unknown, but
Theorem 1.1 shows that global-in-time $L^{2}$-equicontinuity fails in general
for mass greater than $2\pi$. Specifically, for any $\varepsilon>0$
sufficiently small, we find initial data $u_{0}\in H^{\infty}_{+}$ with
$M(u_{0})=2\pi+\varepsilon$ for which the maximal orbit fails to be
$L^{2}$-equicontinuous.

The local theory in [10] allows us to consider $H^{\infty}_{+}$ solutions of
(CCM), which are sufficient for our analysis. For $H^{\infty}_{+}$ solutions,
Laurens, Killip, and Vişan showed in [12] that finite-time blowup must be
accompanied by a loss of $L^{2}$-equicontinuity and hence simultaneous blowup
in $H^{s}$ for all $s>0$. The well-posedness necessary for our main theorem is
then stated as follows:

######  Proposition 1.3.

For initial data $u_{0}\in H^{\infty}_{+}$ there exists a maximal time
$T\in(0,\infty]$ and a unique solution $u(t)$ such that $u\in
C_{t}H^{s}_{+}([0,T)\times\mathbb{R})$ for all $s\geq 0$. Moreover, if
$T<\infty$ then

$\lim_{t\nearrow T}\|u(t)\|_{H^{s}}=+\infty$

for all $s>0$. An analogous statement holds backwards in time.

###  1.1. Overview of proof

We start by proving an orbital stability result for the soliton for a
restricted class of initial data in Theorem 2.1. Specifically, we show that
for any $\varepsilon>0$, if $u_{0}$ is sufficiently close to $Q$ in $H^{1}$
and $u(t)$ is uniformly bounded below in $\dot{H}^{1}$, then

$\inf_{\lambda>0;\;\theta,y\in\mathbb{R}}\|e^{i\theta}\lambda^{\frac{1}{2}}u(t,\lambda
x+y)-Q(x)\|_{H^{1}_{x}}<\varepsilon,$

for all times of existence $t$. This allows us to choose parameters
$\lambda(t)>0$ and $\theta(t),y(t)\in\mathbb{R}$ such that

(1.6)
$\|e^{i\theta(t)}\lambda(t)^{\frac{1}{2}}u(t,\lambda(t)x+y(t))-Q(x)\|_{H^{1}_{x}}<\varepsilon.$

The scaling parameter $\lambda(t)$ then represents the characteristic width of
$u(t)$. We verify that this set of solutions is non-empty by showing that
initial data $u_{0}\in H^{\infty}_{+}$ for which $\mathcal{L}_{u_{0}}$ has a
negative eigenvalue lead to solutions that satisfy a uniform lower bound in
$\dot{H}^{1}$; see Lemma 2.3. An example of such initial data is
$u_{0}=(1+\varepsilon)Q$ for $\varepsilon>0$ as shown in Lemma 2.4.

We then prove a dispersive decay bound in Theorem 3.1 that shows that global-
in-time $H^{\infty}_{+}$ solutions to (CCM) satisfy $|u(t,z)|\to 0$ for
$\operatorname{Im}z>0$ as $t\to\infty$, where $u(t,z)$ is defined by
holomorphic extension to the upper half-plane. For solutions satisfying (1.6),
this dispersive decay is inconsistent with $\lambda(t)\sim 1$ uniformly in
$t$; we must either have that $\lambda(t)\to 0$ or $\lambda(t)\to\infty$ as
$t\to\infty$. The lower bound on $\|u(t)\|_{\dot{H}^{1}}$ rules out the case
of $\lambda(t)\to\infty$ and so we must have $\lambda(t)\to 0$ as
$t\to\infty$. This forces $\|u(t)\|_{H^{s}}$ to grow to $\infty$ as
$t\to\infty$, as desired. The lower bound $\|u(t)\|_{H^{s}}\gtrsim t^{s}$ then
follows from the rate of dispersive decay that we obtain.

###  1.2. Acknowledgements

J.H. and M.K. were supported in part by NSF grant DMS-2154022 and NSF grant
DMS-2054194. We would like to thank our advisors, Rowan Killip and Monica
Vişan, for their guidance and support.

###  1.3. Notation

We use the standard notation $A\lesssim B$ to indicate that $A\leq CB$ for
some universal constant $C>0$. If both $A\lesssim B$ and $B\lesssim A$ then we
use the notation $A\sim B$. When the implied constant fails to be universal,
the relevant dependencies will be indicated within the text.

Our conventions for the Fourier transform are

$\widehat{f}(\xi)=\tfrac{1}{\sqrt{2\pi}}\int e^{-i\xi
x}f(x)dx\quad\text{so}\quad f(x)=\tfrac{1}{\sqrt{2\pi}}\int e^{i\xi
x}\widehat{f}(\xi)d\xi.$

These Fourier transforms are unitary on $L^{2}$ and hence yield the standard
Plancherel identities. With these conventions, we define the Cauchy–Szegö
projector $\Pi^{+}:L^{2}\to L^{2}_{+}$ as

(1.7)
$[\Pi^{+}f](x)=\tfrac{1}{\sqrt{2\pi}}\int_{0}^{\infty}e^{ix\xi}\widehat{f}(\xi)d\xi.$

For $s\geq 0$, we define the Sobolev spaces $H^{s}$ as the completion of
$\mathcal{S}(\mathbb{R})$ with respect to the norm

$\|f\|^{2}_{H^{s}}=\int\big{(}1+|\xi|^{2}\big{)}^{s}|\widehat{f}(\xi)|^{2}d\xi.$

The homogeneous Sobolev spaces $\dot{H}^{s}$ are defined similarly by the norm

$\|f\|^{2}_{\dot{H}^{s}}=\int|\xi|^{2s}|\widehat{f}(\xi)|^{2}d\xi.$

We work in the Hardy-Sobolev spaces, denoted
$H_{+}^{s}=H_{+}^{s}(\mathbb{R})\subset H^{s}(\mathbb{R})$, comprised of those
functions whose Fourier transforms are supported on $[0,\infty)$ as in (1.1).
Equivalently, $H^{s}_{+}$ is comprised of $f\in H^{s}$ whose Poisson integral

$f(z)=\int\frac{\operatorname{Im}z}{\pi|x-z|^{2}}f(x)dx\quad\text{defined
for}\quad\operatorname{Im}z>0$

is holomorphic in the upper half-plane. This formulation of $H^{s}_{+}$ allows
us to work interchangeably with $f(x)$ for $x\in\mathbb{R}$ and $f(z)$ for
$\operatorname{Im}z>0$. With the Poisson formulation, Hölder’s inequality
implies

(1.8) $|f(z)|\lesssim(\operatorname{Im}z)^{-1/p}\|f\|_{L^{p}(\mathbb{R})}$

for all $1\leq p\leq\infty$. In this paper, we mainly work in the spaces
$L^{2}_{+}=H^{0}_{+}$ and $H^{\infty}_{+}=\cap_{s}H^{s}_{+}$.

We observe that (CCM) is invariant under phase rotation, translation, and
$L^{2}$-critical scaling

$u(t,x)\mapsto e^{i\theta}\lambda^{\frac{1}{2}}u(\lambda^{2}t,\lambda x+y)$

for $\lambda>0$ and $y,\theta\in\mathbb{R}$. We adopt the notation

(1.9) $u_{\lambda,y,\theta}(t,x)=e^{i\theta}\lambda^{\frac{1}{2}}u(t,\lambda
x+y)$

for $\lambda>0$ and $y,\theta\in\mathbb{R}$. Note that only the spatial
component is altered here.

##  2\. Orbital stability of the soliton

In this section we seek to prove

######  Theorem 2.1 (Restricted orbital stability of $Q$).

Fix $c>0$. For all $\varepsilon>0$ there exists $\delta>0$ such that if
$u_{0}\in H^{1}$ satisfies

$\|u_{0}-Q\|_{H^{1}}<\delta$

and the corresponding maximal lifespan solution $u$ satisfies
$\|u(t)\|_{\dot{H}^{1}}\geq c$, then

$\inf_{\lambda>0;\;\theta,y\in\mathbb{R}}\|u_{\lambda,\theta,y}(t)-Q\|_{H^{1}}<\varepsilon$

for all times of existence $t$, using the notation from (1.9).

Note that in our orbital distance, the symmetries are applied to $u$ rather
than $Q$. This is reminiscent of modulation theory for gKdV and NLS, [13, 14,
20]. In the case of $s=0$, scaling-criticality implies

$\inf_{\lambda>0;\;\theta,y\in\mathbb{R}}\|u_{\lambda,\theta,y}-Q\|_{L^{2}}=\inf_{\lambda>0;\;\theta,y\in\mathbb{R}}\|u-Q_{\lambda,\theta,y}\|_{L^{2}},$

which is the standard distance from $u$ to the manifold of solitons. In
contrast, for $s\neq 0$, the scaling term prevents the symmetries from being
moved to $Q$ without accruing additional factors of $\lambda$. Because of
this, our orbital distance for $s\neq 0$ differs from the standard notion.

We argue via a variational approach akin to Cazenave and Lions [6]. We
therefore recall the characterization of $Q$ as a minimizer of the energy, as
given by Gérard and Lenzmann in [10].

######  Lemma 2.2 (Minimal mass bubble, [10]).

Suppose that $(v^{n})_{n\in\mathbb{N}}\subset H^{1}$ is a $L^{2}$-bounded
sequence which satisfies

$\displaystyle\lim_{n\to\infty}E(v^{n})=0,\;\;\;\liminf_{n\to\infty}\|v^{n}\|_{L^{2}}^{2}=\|Q\|_{L^{2}}^{2},\;\;\;\text{and}\;\;\;\|v^{n}\|_{\dot{H}^{1}}=\mu>0\quad\text{for
all }n\in\mathbb{N}.$

Then after possibly passing to a subsequence,

$\|v^{n}_{\lambda,\theta,x_{n}}-Q\|_{L^{2}}\to 0$

for some constants $\lambda>0$, $\theta\in\mathbb{R}$, and some sequence
$(x_{n})\subset\mathbb{R}$.

As $Q$ is the unique energy minimizer up to symmetries (see [10]), the
preceding lemma informally asserts that convergence of mass and energy to
those of $Q$ implies convergence in $L^{2}$ to $Q$ up to symmetries. Since
mass and energy are conserved and $Q$ is stationary, this suggests the flow of
(CCM) should preserve closeness to $Q$ in the $L^{2}$ topology. To upgrade
this to the $H^{1}$ topology, we use the conservation of energy. This is
formalized in what follows.

###### Proof of Theorem 2.1.

We argue by contradiction. Suppose that there exist some $\varepsilon>0$ and a
sequence of initial data $u^{n}_{0}\to Q$ in $H^{1}$ with maximal lifespan
solutions $u^{n}$ such that $\|u^{n}(t)\|_{\dot{H}^{1}}\geq c$ and

$\inf_{\lambda>0;\;\theta,y\in\mathbb{R}}\|u^{n}_{\lambda,\theta,y}(t_{n})-Q\|_{H^{1}}\geq\varepsilon$

for some time of existence $t_{n}$.

To reach a contradiction, we aim to apply Lemma 2.2 to $u^{n}(t_{n})$. We
first normalize $u^{n}(t_{n})$ in $\dot{H}^{1}$ and define

$\lambda_{n}=\|u^{n}(t_{n})\|_{\dot{H}^{1}}^{-1}\leq
c^{-1}\quad\text{and}\quad v^{n}(x)=u^{n}_{\lambda_{n},0,0}(t_{n}).$

We note that

$\sup_{n}\|v^{n}\|_{L^{2}}=\sup_{n}\|u^{n}_{0}\|_{L^{2}}<\infty\quad\text{and}\quad\|v^{n}\|_{\dot{H}^{1}}=\lambda_{n}\|u^{n}(t_{n})\|_{\dot{H}^{1}}=1.$

Furthermore,

$0\leq E(v^{n})=\lambda_{n}^{2}E(u_{0}^{n})\leq c^{-2}E(u_{0}^{n})\to 0$

because $u_{0}^{n}\to Q$ in $H^{1}$. Passing to a subsequence, Lemma 2.2 then
guarantees that

$\|v^{n}_{\lambda,\theta,x_{n}}-Q\|_{L^{2}}\to 0$

for some $\lambda>0$, $\theta\in\mathbb{R}$, and a sequence
$(x_{n})\subset\mathbb{R}$.

Next, we upgrade this $L^{2}$-convergence to $H^{1}$-convergence using an
energy argument. Define

$w^{n}=v^{n}_{\lambda,\theta,x_{n}}-Q.$

To show that $v^{n}_{\lambda,\theta,x_{n}}\to Q$ in $H^{1}$ it then suffices
to show that $\|w^{n}\|_{\dot{H}^{1}}\to 0$. Because

$\mathcal{L}_{Q}Q=-iQ_{x}-\Pi_{+}(|Q|^{2})Q=0$

we compute that

$\displaystyle E(v^{n}_{\lambda,\theta,x_{n}})$
$\displaystyle=\tfrac{1}{2}\int|\partial_{x}w^{n}-i\Pi_{+}(|Q|^{2})w^{n}-i\Pi_{+}(|w^{n}|^{2}+2\mathrm{Re}\,\overline{Q}w^{n})(Q+w^{n})|^{2}$
$\displaystyle=\tfrac{1}{2}\int|\partial_{x}w^{n}|^{2}+O(\|w^{n}\|_{L^{2}}+\|w^{n}\|_{L^{2}}^{2}).$

We note that the preceding calculations repeatedly use the fact that
$\|w^{n}\|_{L^{\infty}}\lesssim 1$ uniformly in $n$, which follows from
Sobolev embedding.

Because $w^{n}\to 0$ in $L^{2}$ and
$E(v^{n}_{\lambda,\theta,x_{n}})=\lambda^{2}E(v^{n})\to 0$, it then follows
that $\|w^{n}\|_{\dot{H}^{1}}\to 0$ and so $v^{n}_{\lambda,\theta,x_{n}}\to Q$
in $H^{1}$. This implies that

$\lim_{n\to\infty}\inf_{\lambda>0;\;\theta,y\in\mathbb{R}}\|u^{n}_{\lambda,\theta,y}(t_{n})-Q\|_{H^{1}}=0$

which is a contradiction. ∎

To verify that this stability result is not vacuous, we show that a solution
$u(t)$ for which $\mathcal{L}_{u}$ has a negative eigenvalue satisfies
$\|u(t)\|_{\dot{H}^{1}}\gtrsim 1$ uniformly over the times of existence. As
$\mathcal{L}_{u},\mathcal{P}_{u}$ form a Lax pair, we expect that the map
$t\mapsto\mathcal{L}_{u(t)}$ is an isospectral deformation and thus that the
negative eigenvalue of $\mathcal{L}_{u}$ is preserved in time. This was shown
in [15] in the case that $\mathcal{L}_{u}$ and $\mathcal{P}_{u}$ are matrices
and applied to (CM). To extend this result to the case where $\mathcal{L}_{u}$
and $\mathcal{P}_{u}$ are unbounded operators, more careful analysis is
required. For an example of this analysis for (CCM), see [12], where it is
shown that $t\mapsto\mathcal{L}_{u(t)}$ is an isospectral deformation provided
that $u(t)\in H^{\infty}_{+}$.

######  Lemma 2.3 (Uniform lower bound in $\dot{H}^{1}$).

Suppose there exist $\varepsilon>0$ and $u_{0}\in H^{\infty}_{+}$ such that
$\|u_{0}\|_{L^{2}}^{2}\leq 2\pi+\varepsilon$ and such that
$\mathcal{L}_{u_{0}}$ has an eigenvalue in $(-\infty,-c\varepsilon]$ for some
$c>0$. Let $u$ denote the corresponding maximal lifespan solution to (CCM).
Then

$\|u(t)\|_{\dot{H}^{1}}\gtrsim c$

uniformly for all times of existence $t$.

###### Proof.

Since $\mathcal{L}_{u(t)}$ is isospectral, there exists an eigenvector
$\psi(t)\in H^{\frac{1}{2}}_{+}$ for all times of existence $t$ with
$\|\psi(t)\|_{L^{2}}=1$ such that

$\langle\psi(t),\mathcal{L}_{u(t)}\psi(t)\rangle\leq-c\varepsilon.$

Expanding $\mathcal{L}_{u(t)}$ yields

$\displaystyle-c\varepsilon$
$\displaystyle\geq\langle\psi(t),-i\partial_{x}\psi(t)\rangle-\Bigg{[}\frac{2\pi}{\|u(t)\|_{L^{2}}^{2}}+\frac{\|u(t)\|_{L^{2}}^{2}-2\pi}{\|u(t)\|_{L^{2}}^{2}}\Bigg{]}\|\Pi^{+}\overline{u(t)}\psi(t)\|_{L^{2}}^{2}.$

Inequality (1.4) then implies

$\displaystyle c\varepsilon$
$\displaystyle\leq\frac{\|u(t)\|_{L^{2}}^{2}-2\pi}{\|u(t)\|_{L^{2}}^{2}}\|\Pi^{+}\overline{u(t)}\psi(t)\|_{L^{2}}^{2}$
$\displaystyle\leq\frac{\varepsilon}{2\pi+\varepsilon}\|u(t)\|_{\infty}^{2}$
$\displaystyle\leq\varepsilon\|u(t)\|_{\infty}^{2}.$

Therefore $c\leq\|u\|_{\infty}^{2}$ and so
$c\lesssim\|u\|_{L^{2}}\|u\|_{\dot{H}^{1}}\lesssim\|u\|_{\dot{H}^{1}}$ as
desired. ∎

We now show that initial data of the form $u_{0}=(1+\varepsilon)Q$ for
$\varepsilon>0$ has a uniformly negative eigenvalue. This data will be the
starting point for our proof of the turbulent threshold.

######  Lemma 2.4 (Existence of negative eigenvalue).

For any $\varepsilon>0$, the Lax operator $\mathcal{L}_{(1+\varepsilon)Q}$ has
a negative eigenvalue in $(-\infty,-\varepsilon]$.

###### Proof.

Because the Lax operator is self-adjoint,

$\displaystyle 2\pi\inf\sigma(\mathcal{L}_{(1+\varepsilon)Q})$
$\displaystyle=\min_{\psi\neq
0}\frac{2\pi}{\|\psi\|_{L^{2}}^{2}}\left\langle\psi,\mathcal{L}_{(1+\varepsilon)Q}\psi\right\rangle$
$\displaystyle\leq\left\langle Q,\mathcal{L}_{(1+\varepsilon)Q}Q\right\rangle$
$\displaystyle=\left\langle
Q,-i\partial_{x}Q\right\rangle-(1+\varepsilon)^{2}\|\Pi_{+}|Q|^{2}\|_{L^{2}}^{2}$
$\displaystyle=\left\langle
Q,\mathcal{L}_{Q}Q\right\rangle-(2\varepsilon+\varepsilon^{2})\pi$
$\displaystyle\leq-2\pi\varepsilon,$

where the final equality follows from $\mathcal{L}_{Q}Q=0$. Therefore
$\inf\sigma(\mathcal{L}_{(1+\varepsilon)Q})\leq-\varepsilon$. Because
$\sigma_{\text{ess}}(\mathcal{L}_{(1+\varepsilon)Q})=[0,\infty)$ (see [10,
12]), this concludes the proof. ∎

##  3\. Dispersive decay

Our goal in this section is to establish the following dispersive decay
estimate for solutions to (CCM) with suitable initial data:

######  Theorem 3.1 (Dispersive decay).

Given a set of initial data $U\subset H^{\infty}_{+}$ which is bounded and
equicontinuous in $L^{2}_{+}$ and satisfies $\langle x\rangle u_{0}\in L^{2}$
for all $u_{0}\in U$,

(3.1)
$|u(t,z)|\lesssim|t|^{-\frac{1}{2}}\|u_{0}\|_{L^{1}}\big{[}1+\|u_{0}\|_{L^{2}}^{2}(1+(\operatorname{Im}z)^{-1})\big{]}$

uniformly for $\operatorname{Im}z>0$, $u_{0}\in U$, and all times of existence
$t$.

Using Theorem 2.1, we will consider solutions $u(t)$ for which

$\|u_{\lambda(t),\theta(t),y(t)}(t)-Q\|_{H^{1}}<\varepsilon,$

for some $\lambda(t)>0$ and $\theta(t),y(t)\in\mathbb{R}$. Theorem 3.1 then
implies that the characteristic width $\lambda(t)$ of $u(t)$ cannot satisfy
$\lambda(t)\sim 1$ uniformly in $t$. This deviates from many classical
stability results for which the scaling remained constant; e.g. [21].

To prove Theorem 3.1, we use the explicit formula established in [12] for
$H^{\infty}_{+}$ solutions to (CCM). This explicit formula is of the type
developed by Gérard and collaborators; see [8, 9, 11]. A similar explicit
formula was established for an intertwined matrix version of (CCM) in [17].
The explicit formula for (CCM) was also shown to hold for initial data in
$L^{2}_{+}$ below the mass of $2\pi$ in [12], which can be used to expand the
application of Theorem 3.1. As we are largely interested in $H^{\infty}_{+}$
solutions, we recall the explicit formula in that case.

######  Theorem 3.2 (Explicit formula).

For any $H^{\infty}_{+}(\mathbb{R})$ solution $u(t)$ to (CCM) with initial
data satisfying $\langle x\rangle u_{0}\in L^{2}$,

$u(t,z)=\tfrac{1}{2\pi
i}I_{+}\Big{\\{}(X+2t\mathcal{L}_{u_{0}}-z)^{-1}u_{0}\Big{\\}}$

for all $\operatorname{Im}z>0$ and all times of existence $t$.

Here $X$ and $I_{+}$ are defined as follows:

######  Definition 3.3.

The unbounded linear functional $I_{+}$ is defined by

$I_{+}(f):=\lim_{y\to\infty}\sqrt{2\pi}\int_{0}^{\infty}ye^{-y\xi}\widehat{f}(\xi)d\xi$

with domain $D(I_{+})$ given by the set of those $f\in L^{2}_{+}(\mathbb{R})$
for which the limit exists.

######  Definition 3.4.

The unbounded operator $X:L^{2}_{+}\to L^{2}_{+}$ is given by

$\widehat{Xf}(\xi):=i\tfrac{d\widehat{f}}{d\xi}(\xi)\quad\text{for}\quad f\in
D(X)$

where $D(X)$ is the domain $D(X)=\\{f\in L^{2}_{+}(\mathbb{R}):\widehat{f}\in
H^{1}([0,\infty])$.

To simplify notation, let

$A(t,z;u):=(X+2t\mathcal{L}_{u}-z)^{-1},$

so that the explicit formula in Theorem 3.2 becomes

$u(t,z)=\tfrac{1}{2\pi i}I_{+}A(t,z;u_{0})u_{0}.$

Using the resolvent identity as was done in [12], we may write

(3.2)
$A(t,z;u_{0})=A_{0}(t,z)+2tA_{0}(t,z)u_{0}\Pi^{+}\overline{u_{0}}A(t,z;u_{0})$

where $A_{0}(t,z)=A(t,z;0)$.

We will employ the following bounds on $A_{0}$ and $A$:

######  Proposition 3.5 (Inequalities (4.21) and (4.41) from [12]).

The operator $I_{+}A_{0}$ is bounded on $L^{1}$, namely

$|I_{+}A_{0}(t,z)f|\lesssim|t|^{-\frac{1}{2}}\|f\|_{L^{1}},$

uniformly for $t\in\mathbb{R}$ and $\operatorname{Im}z>0$.

Given a bounded and equicontinuous set $U\subset L^{2}_{+}(\mathbb{R})$, the
bound

$\|A(t,z;u)f\|_{L^{\infty}}\lesssim|t|^{-1}\big{[}1+(\operatorname{Im}z)^{-1}\big{]}\|f\|_{L^{1}}$

holds uniformly for $u\in U$, $t\in\mathbb{R}$, $f\in L^{1}$, and
$\operatorname{Im}z>0$.

Combining these, we establish a dispersive decay for $H^{\infty}_{+}$
solutions.

###### Proof of Theorem 3.1.

Fix $u_{0}\in U$. A simple Grönwall argument ensures that $\langle x\rangle
u(t)\in L^{2}$ for all times of existence and Hölder’s inequality then
guarantees that $U\subset L^{1}$. Combining the explicit formula in Theorem
3.2, the resolvent identity (3.2), and the bounds from Proposition 3.5, we may
estimate

$\displaystyle|u(t,z)|$
$\displaystyle\lesssim\left|I_{+}\\{A_{0}(t,z)u_{0}\\}\right|+|t|\left|I_{+}\\{A_{0}(t,z)u_{0}\Pi^{+}\overline{u_{0}}A(t,z;u_{0})u_{0}\\}\right|$
$\displaystyle\lesssim|t|^{-\frac{1}{2}}\|u_{0}\|_{L^{1}}+|t|^{\frac{1}{2}}\|u_{0}\Pi^{+}\overline{u_{0}}A(t,z;u_{0})u_{0}\|_{L^{1}}$
$\displaystyle\lesssim|t|^{-\frac{1}{2}}\|u_{0}\|_{L^{1}}\big{[}1+\|u_{0}\|_{L^{2}}^{2}(1+(\operatorname{Im}z)^{-1})\big{]}$

uniformly for $\operatorname{Im}z>0$, $u_{0}\in U$, and all times of existence
$t$. ∎

##  4\. Growth of Sobolev norms

In this section, we find $H^{\infty}_{+}$ solutions arbitrarily close to, but
above, the mass threshold of $2\pi$ that exhibit global Sobolev growth:

######  Theorem 4.1 (Turbulent behaviour).

For any sufficiently small $\varepsilon>0$, there exist initial data $u_{0}\in
H^{\infty}_{+}(\mathbb{R})$ with $M(u_{0})=2\pi+\varepsilon$, a time
$T\in(0,\infty]$, and a solution $u(t)$ to (CCM) such that for all $s>0$,
$u\in C_{t}H^{s}_{x}([0,T)\times\mathbb{R})$ and

$\lim_{t\nearrow T}\|u(t)\|_{H^{s}}=+\infty.$

In particular, if $T=\infty$, then we have the bounds

$\|u(t)\|_{H^{s}}\gtrsim t^{s}.$

An analogous statement holds backwards in time.

Our initial data will be a perturbation of the soliton $Q$ belonging to
$\mathcal{S}(\mathbb{R})\cap L^{2}_{+}$. This is more than sufficient to
satisfy the regularity constraints arising in Theorem 3.1. Our proof will
require that the implicit constants which come from Theorem 3.1 are uniform;
see the progression from (4.2) to (4.3). We therefore construct an
$L^{2}$-equicontinuous family of initial data as follows:

######  Lemma 4.2.

There exists an $L^{2}$-equicontinuous family
$U=\\{u_{0}^{\varepsilon}:0<\varepsilon\leq 1\\}\subset
L^{2}_{+}\cap\mathcal{S}(\mathbb{R})$ such that
$\|u_{0}^{\varepsilon}-Q\|_{L^{2}}\lesssim\varepsilon$ and the corresponding
solutions $u^{\varepsilon}(t)$ satisfy

$\|u^{\varepsilon}(t)\|_{\dot{H}^{1}}\gtrsim
1\quad\text{and}\quad\inf_{\lambda>0;\;\theta,y\in\mathbb{R}}\|u_{\lambda,\theta,y}^{\varepsilon}(t)-Q\|_{H^{1}}<\varepsilon$

uniformly for all times of existence $t$ and initial data
$u^{\varepsilon}_{0}\in U$.

###### Proof.

To define $u_{0}^{\varepsilon}$, we begin with $(1+\varepsilon)Q$ which has a
negative eigenvalue in $(-\infty,-\varepsilon]$ as shown in Lemma 2.4. We
recall that $Q$ has the Fourier transform

(4.1) $\widehat{Q}(\xi)=2\sqrt{\pi}ie^{-\xi}\chi_{\\{\xi\geq 0\\}}.$

Therefore the only obstruction to $(1+\varepsilon)Q$ belonging to
$\mathcal{S}(\mathbb{R})$ is that $(1+\varepsilon)\widehat{Q}$ has a
discontinuity at $\xi=0$.

We define $u^{\varepsilon}_{0}\in\mathcal{S}\cap L^{2}_{+}$ to be a smoothing
of $(1+\varepsilon)Q$ on the Fourier side which removes the discontinuity at
$\xi=0$. This can be done with sufficiently small perturbation in $H^{1}$ to
ensure that $\mathcal{L}_{u^{\varepsilon}_{0}}$ has a negative eigenvalue in
$(-\infty,\frac{\varepsilon}{2}]$.

Let $u^{\varepsilon}(t)$ denote the maximal solution to (CCM) with initial
data $u_{0}^{\varepsilon}$. The bounds
$\|u^{\varepsilon}(t)\|_{\dot{H}^{1}}\gtrsim 1$ and
$\inf_{\lambda>0;\theta,y\in\mathbb{R}}\|u_{\lambda,\theta,\varepsilon}^{\varepsilon}(t)-Q\|_{H^{1}}<\varepsilon$
then follow from Theorem 2.1 and Lemma 2.3. Since $L^{2}$-equicontinuity is
equivalent to tightness on the Fourier side, the equicontinuity of $U$ follows
from the exponential decay of $\widehat{Q}$. ∎

With the initial data established, we now prove the main result.

###### Proof of Theorem 4.1.

Fix initial data $u^{\varepsilon}_{0}\in U$ from Lemma 4.2 for
$0<\varepsilon\leq 1$ to be chosen later. Let $u^{\varepsilon}(t)$ denote the
maximal solution to (CCM) with initial data $u^{\varepsilon}_{0}$.

Choose $\lambda(t)>0$ and $y(t),\theta(t)\in\mathbb{R}$ for all times of
existence $t$ so that

$\|u^{\varepsilon}_{\lambda(t),\theta(t),y(t)}(t)-Q\|_{H^{1}}<\varepsilon.$

Define

$q(t,x)=e^{-i\theta(t)}\lambda(t)^{-\frac{1}{2}}Q\Big{(}\frac{x-y(t)}{\lambda(t)}\Big{)}.$

By scaling criticality, we note that
$\|u^{\varepsilon}(t)-q(t)\|_{L^{2}}<\varepsilon$. We will start by showing
that $\lambda(t)\to 0$ and thus $\|q(t)\|_{\dot{H}^{s}}\to\infty$ for all
$s>0$.

We initially show a uniform upper bound for $\lambda(t)$. Since
$\|u^{\varepsilon}(t)\|_{\dot{H}^{1}}\gtrsim 1$, we find

$\displaystyle 1$
$\displaystyle\lesssim\|u^{\varepsilon}(t)-q(t)\|_{\dot{H}^{1}}+\|q(t)\|_{\dot{H}^{1}}$
$\displaystyle=\lambda(t)^{-1}\|u^{\varepsilon}_{\lambda(t),\theta(t),y(t)}(t)-Q\|_{\dot{H}^{1}}+\lambda(t)^{-1}\|Q\|_{\dot{H}^{1}}$
$\displaystyle\lesssim\lambda(t)^{-1}\varepsilon+\lambda(t)^{-1}\|Q\|_{\dot{H}^{1}}.$

Therefore $\lambda(t)\lesssim\varepsilon+\|Q\|_{\dot{H}^{1}}\lesssim 1$ for
$0<\varepsilon\leq 1$.

To conclude $\lambda(t)\to 0$, we use the dispersive decay bounds. By Theorem
3.1 and inequality (1.8), for all $\operatorname{Im}z>0$,

$\displaystyle|q(t,z)|$
$\displaystyle\leq|q(t,z)-u^{\varepsilon}(t,z)|+|u^{\varepsilon}(t,z)|$
$\displaystyle\lesssim\|q(t)-u^{\varepsilon}(t)\|_{L^{2}}(\operatorname{Im}z)^{-\frac{1}{2}}+|t|^{-\frac{1}{2}}\|u^{\varepsilon}_{0}\|_{L^{1}}\big{[}1+\|u^{\varepsilon}_{0}\|_{L^{2}}^{2}(1+(\operatorname{Im}z)^{-1})\big{]}$
$\displaystyle\lesssim\varepsilon(\operatorname{Im}z)^{-\frac{1}{2}}+|t|^{-\frac{1}{2}}\|u^{\varepsilon}_{0}\|_{L^{1}}\big{[}1+\|u^{\varepsilon}_{0}\|_{L^{2}}^{2}(1+(\operatorname{Im}z)^{-1})\big{]}.$

Note that the implicit constants are uniform in $\varepsilon$ by the
$L^{2}$-equicontinuity of $U$.

We now evaluate this inequality at the point $z(t)=y(t)+i\lambda(t)$. Direct
calculation shows

(4.2)
$(2\lambda(t))^{-\frac{1}{2}}=|q(z(t))|\lesssim\varepsilon\lambda(t)^{-\frac{1}{2}}+|t|^{-\frac{1}{2}}\|u^{\varepsilon}_{0}\|_{L^{1}}\big{[}1+\|u^{\varepsilon}_{0}\|_{L^{2}}^{2}(1+\lambda(t)^{-1})\big{]}.$

Because the implicit constants are uniform in $0<\varepsilon\leq 1$, we may
take $\varepsilon>0$ sufficiently small to obtain

(4.3)
$\lambda(t)^{\frac{1}{2}}\lesssim|t|^{-\frac{1}{2}}\|u^{\varepsilon}_{0}\|_{L^{1}}\big{[}1+\|u^{\varepsilon}_{0}\|_{L^{2}}^{2}(1+\lambda(t)^{-1})\big{]},$

and so

$|t|^{\frac{1}{2}}\lesssim\|u^{\varepsilon}_{0}\|_{L^{1}}\Big{[}\big{(}1+\|u^{\varepsilon}_{0}\|_{L^{2}}^{2}\big{)}\lambda(t)^{\frac{1}{2}}+\|u^{\varepsilon}_{0}\|^{2}_{L^{2}}\lambda(t)^{-\frac{1}{2}}\Big{]}.$

The upper bounds on $\lambda(t)$ and $\|u^{\varepsilon}_{0}\|_{L^{2}}$ then
imply

$|t|^{\frac{1}{2}}\lesssim\|u^{\varepsilon}_{0}\|_{L^{1}}\big{[}1+\lambda(t)^{-\frac{1}{2}}\big{]}.$

Since $\|u_{0}^{\varepsilon}\|_{L^{1}}<\infty$, it follows that
$\lambda(t)\lesssim|t|^{-1}$ for sufficiently large $|t|$.

Using the formula (4.1) for $\widehat{Q}$, direct calculation yields

$\displaystyle\|u^{\varepsilon}(t)\|_{\dot{H}^{s}}^{2}$
$\displaystyle\gtrsim|\lambda(t)|^{-2s}\int_{1/\lambda(t)}^{\infty}|\widehat{u^{\varepsilon}}(t,\xi)|^{2}d\xi$
$\displaystyle\gtrsim|\lambda(t)|^{-2s}\bigg{(}\int_{1/\lambda(t)}^{\infty}|\widehat{q}(t,\xi)|^{2}d\xi-\varepsilon\bigg{)}$
$\displaystyle=|\lambda(t)|^{-2s}\bigg{(}\int_{1/\lambda(t)}^{\infty}4\pi\lambda(t)e^{-2\lambda(t)\xi}d\xi-\varepsilon\bigg{)}$
$\displaystyle=|\lambda(t)|^{-2s}\big{(}2\pi e^{-2}-\varepsilon\big{)}.$

Choosing $\varepsilon>0$ sufficiently small, this implies that
$\|u^{\varepsilon}(t)\|_{\dot{H}^{s}}\gtrsim|t|^{s}$ for sufficiently large
$|t|$. Because $u^{\varepsilon}$ is either global in time or exhibits finite-
time blowup, this concludes the proof of the theorem. ∎

A natural further question is whether there is interesting behaviour below the
$2\pi$ mass threshold. From the a priori $H^{s}$ bounds, we cannot expect
Sobolev norm growth, but there is potential for some decay in the homogeneous
Sobolev spaces $\dot{H}^{s}$ for $s>0$. This is formalized in the following
corollary.

######  Corollary 4.3.

For all $c,\varepsilon,s>0$ there exists initial data $u_{0}\in
H^{\infty}_{+}$ such that $M(u_{0})<2\pi$ and
$\|u_{0}-Q\|_{L^{2}}<\varepsilon$ for which the global solution $u(t)$
satisfies

$\inf_{t\in\mathbb{R}}\|u(t)\|_{\dot{H}^{s}}\leq c.$

###### Proof.

We begin with the case $s=1$. Suppose for the sake of contradiction that there
exist $c_{0},\varepsilon_{0}>0$ such that for all initial data $u_{0}\in
H^{\infty}_{+}$ with $M(u_{0})<2\pi$ and
$\|u_{0}-Q\|_{L^{2}}<\varepsilon_{0}$, the corresponding global solution
$u(t)$ satisfies

$\inf_{t}\|u(t)\|_{\dot{H}^{1}}>c_{0}.$

Starting with initial data of the form $(1-\varepsilon)Q$ for
$0<\varepsilon<\varepsilon_{0}$, we follow the proof of Lemma 4.2 to construct
an $L^{2}$-equicontinuous family
$U=\\{u_{0}^{\varepsilon}:0<\varepsilon<\varepsilon_{0}\\}$ satisfying all the
conditions in Lemma 4.2 along with $M(u_{0}^{\varepsilon})<2\pi$ for all
$0<\varepsilon<\varepsilon_{0}$. By only smoothing
$(1-\varepsilon)\widehat{Q}$ near $\xi=0$, we further ensure for all $r>0$,

(4.4) $\lim_{\varepsilon\to 0}\|u_{0}^{\varepsilon}-Q\|_{H^{r}}=0.$

This will be necessary to establish the case when $0\leq s<1$. Note that in
the proof of Lemma 4.2, we used initial data with a negative eigenvalue in
order to gain a lower bound in $\dot{H}^{1}$. In contrast, here we have
supposed for the sake of contradiction that such a lower bound exists. This is
necessary because (1.5) implies that eigenvalues cannot exist below the mass
of $2\pi$.

Fix $u_{0}^{\varepsilon}\in U$ for $0<\varepsilon<\varepsilon_{0}$. Because
$M(u_{0}^{\varepsilon})<2\pi$, there exists a global solution
$u^{\varepsilon}(t)$ to (CCM) with initial data $u_{0}^{\varepsilon}$.
Following the same computations as in the proof of Theorem 4.1, we choose
sufficiently small $\varepsilon>0$ and find that
$\|u^{\varepsilon}(t)\|_{\dot{H}^{1}}\gtrsim|t|$ for sufficiently large $|t|$.
This contradicts the a priori $H^{1}$ bounds from [10], concluding the proof
of the corollary in the case $s=1$.

For $0<s<1$, we use the bound

$\|u^{\varepsilon}\|_{\dot{H}^{s}}\lesssim\|u^{\varepsilon}\|_{\dot{H}^{1}}^{s}\|u^{\varepsilon}\|_{L^{2}}^{1-s}$

along with conservation of mass to prove the result.

For $s>1$, we choose an integer $n>2s-1$. Since $u_{0}^{\varepsilon}\to Q$ in
$H^{n}$ by (4.4), it follows that the conserved quantities
$I_{0}(u^{\varepsilon}),...,I_{2n}(u^{\varepsilon})$ are uniformly bounded and
thus $\|u^{\varepsilon}(t)\|_{H^{n}}$ satisfy global-in-time $H^{s}$ bounds
uniformly in both $t$ and $\varepsilon$; see [10]. We then use the inequality

$\|u^{\varepsilon}\|_{\dot{H}^{s}}\lesssim\|u^{\varepsilon}\|_{\dot{H}^{1}}^{\frac{1}{2}}\|u^{\varepsilon}\|_{\dot{H}^{2s-1}}^{\frac{1}{2}},$

to conclude the proof. ∎

Despite this corollary, it is impossible for individual $H^{\infty}_{+}$
solutions to decay to zero in $\dot{H}^{1}$. A quick calculation shows that

$E(u)\lesssim\big{[}1+M(u)^{2}\big{]}\|u(t)\|^{2}_{\dot{H}^{1}}$

and so energy provides a non-trivial lower bound for $\|u(t)\|_{\dot{H}^{1}}$
as soon as $E(u)>0$. In [10], it was shown that $\mathcal{L}_{u}u=0$, and
hence $E(u)=0$, if and only if either $u=0$ or $u=Q_{\lambda,\theta,y}$ for
some $\lambda>0$ and $\theta,y\in\mathbb{R}$. Since these solutions do not
evolve in time, this rules out the possibility of an individual nonzero
solution satisfying $\inf_{t}\|u(t)\|_{\dot{H}^{1}}=0$.

For $s>1$, the bound

$\|f\|_{\dot{H}^{1}}\lesssim\|f\|^{1-\frac{1}{s}}_{L^{2}}\|f\|^{\frac{1}{s}}_{\dot{H}^{s}}$

shows that decay in $\dot{H}^{s}$ implies decay in $\dot{H}^{1}$. Thus there
cannot exist an individual solution $u$ in $H^{\infty}_{+}$ which satisfies
$\inf_{t}\|u(t)\|_{\dot{H}^{s}}=0$. For $0\leq s<1$, the global-in-time bounds
on $H^{s}$ and the bound

$\|f\|_{\dot{H}^{1}}\lesssim\|f\|^{\frac{1}{2}}_{\dot{H}^{s}}\|f\|^{\frac{1}{2}}_{\dot{H}^{2-s}}$

imply that the same statement holds for $H^{\infty}_{+}$ solutions with mass
below $2\pi$.

## References

  * [1] A. G. Abanov, E. Bettelheim, and P. Wiegmann.  Integrable hydrodynamics of Calogero–Sutherland model: bidirectional Benjamin–Ono equation.  J. Phys. A, 42(13):135201, 24, 2009.
  * [2] R. Badreddine.  On the global well-posedness of the Calogero–Sutherland derivative nonlinear Schrödinger equation.  Preprint arXiv:2303.01087.
  * [3] R. Badreddine.  Traveling waves & finite gap potentials for the Calogero–Sutherland derivative nonlinear Schrödinger equation.  Preprint arXiv:2307.01592.
  * [4] F. Calogero.  Solution of the one-dimensional $N$-body problems with quadratic and/or inversely quadratic pair potentials.  J. Mathematical Phys., 12:419–436, 1971.
  * [5] F. Calogero and C. Marchioro.  Exact solution of a one-dimensional three-body scattering problem with two-body and/or three-body inverse-square potentials.  J. Mathematical Phys., 15:1425–1430, 1974.
  * [6] T. Cazenave and P.-L. Lions.  Orbital stability of standing waves for some nonlinear Schrödinger equations.  Comm. Math. Phys., 85(4):549–561, 1982.
  * [7] R. P. de Moura and D. Pilod.  Local well posedness for the nonlocal nonlinear Schrödinger equation below the energy space.  Adv. Differential Equations, 15(9-10):925–952, 2010.
  * [8] P. Gérard.  An explicit formula for the Benjamin–Ono equation.  Tunis. J. Math., 5(3):593–603, 2023.
  * [9] P. Gérard and S. Grellier.  An explicit formula for the cubic Szegő equation.  Trans. Amer. Math. Soc., 367(4):2979–2995, 2015.
  * [10] P. Gérard and E. Lenzmann.  The Calogero–Moser derivative nonlinear Schrödinger equation.  Preprint arXiv:2208.04105.
  * [11] P. Gérard and A. Pushnitski.  The cubic Szegő equation on the real line: explicit formula and well-posedness on the Hardy class.  Preprint arXiv:2307.06734.
  * [12] R. Killip, T. Laurens, and M. Vişan.  Scaling-critical well-posedness for continuum Calogero–Moser models.  Preprint arXiv:2311.12334.
  * [13] Y. Martel and F. Merle.  Instability of solitons for the critical generalized Korteweg-de Vries equation.  Geom. Funct. Anal., 11(1):74–123, 2001.
  * [14] Y. Martel and F. Merle.  Blow up in finite time and dynamics of blow up solutions for the $L^{2}$-critical generalized KdV equation.  J. Amer. Math. Soc., 15(3):617–664, 2002.
  * [15] J. Moser.  Three integrable Hamiltonian systems connected with isospectral deformations.  In Surveys in applied mathematics (Proc. First Los Alamos Sympos. Math. in Natural Sci., Los Alamos, N.M., 1974), pages 235–258. Academic Press, New York-London, 1976.
  * [16] D. Pelinovsky.  Intermediate nonlinear Schrödinger equation for internal waves in a fluid of finite depth.  Physics Letters A, 197(5-6):401–406, 1995.
  * [17] Ruoci Sun.  The intertwined derivative Schrödinger system of Calogero–Moser–Sutherland type.  Preprint HAL:04227081.
  * [18] B. Sutherland.  Exact results for a quantum many-body problem in one dimension.  Physical Review A, 4(5):2019, 1971.
  * [19] B. Sutherland.  Exact results for a quantum many-body problem in one dimension. II.  Physical Review A, 5(3):1372, 1972.
  * [20] M. I. Weinstein.  Modulational stability of ground states of nonlinear Schrödinger equations.  SIAM J. Math. Anal., 16(3):472–491, 1985.
  * [21] M. I. Weinstein.  Lyapunov stability of ground states of nonlinear dispersive evolution equations.  Comm. Pure Appl. Math., 39(1):51–67, 1986.