# Constructive proofs for some semilinear PDEs on
$H^{2}(e^{|x|^{2}/4},\R^{d})$

Maxime Breden Hugo Chu CMAP, CNRS, École polytechnique, Institut Polytechnique
de Paris, 91120 Palaiseau, France<EMAIL_ADDRESS>of Mathematics, Imperial College London, London SW7 2AZ, United Kingdom.
<EMAIL_ADDRESS>
###### Abstract

We develop computer-assisted tools to study semilinear equations of the form

$-\Delta u-\frac{x}{2}\cdot\nabla{u}=f(x,u,\nabla u),\quad
x\in\mathbb{R}^{d}.$

Such equations appear naturally in several contexts, and in particular when
looking for self-similar solutions of parabolic PDEs. We develop a general
methodology, allowing us not only to prove the existence of solutions, but
also to describe them very precisely. We introduce a spectral approach based
on an eigenbasis of $\mathcal{L}:=-\Delta-\frac{x}{2}\cdot\nabla$ in spherical
coordinates, together with a quadrature rule allowing to deal with
nonlinearities, in order to get accurate approximate solutions. We then use a
Newton-Kantorovich argument, in an appropriate weighted Sobolev space, to
prove the existence of a nearby exact solution. We apply our approach to
nonlinear heat equations, to nonlinear Schrödinger equations and to a
generalised viscous Burgers equation, and obtain both radial and non-radial
self-similar profiles.

Keywords

Self-similar solutions $\cdot$ Elliptic PDEs $\cdot$ Computer-assisted proofs
$\cdot$ Weighted Sobolev spaces $\cdot$ Unbounded domains

Mathematics Subject Classification (2020)

35C06 $\cdot$ 35J61 $\cdot$ 65G20 $\cdot$ 65N15

##  Introduction

Partial differential equations (PDEs) are fundamental in the modelling of many
physical and biological systems, and the computation of their solutions is
often well understood on bounded domains. However, some PDEs can lose some key
features such as particular solutions or conservation laws when restricted to
a bounded domain. Furthermore, many problems, such as the ones studied in this
work, are more natural to consider on an unbounded domain, either because
their structure leads to simpler proofs or because they automatically give
rise to localised solutions whose decay is intrinsically prescribed by the
equation.

Among those, diffusion problems (both linear and nonlinear), many of which are
natural to consider on $\R^{d}$, have a simple treatment on bounded domains
based on standard Sobolev spaces but which fails to directly extend on
unbounded domains. Many alternative techniques have been developed or applied
to the study of diffusions on $\R^{d}$ such as variational formulations [19,
27, 53], Lyapunov (including energy/entropy) techniques and shooting methods
for radial solutions [11, 38]. In particular, one popular approach is to
recover some estimates which hold in the bounded case such as the Poincaré
inequality by working on weighted Sobolev spaces instead.

Slightly more recently, a new family of tools has emerged, namely computer-
assisted proofs (CAPs), which are complementary to the classical techniques.
One of the benefits of CAPs is that the existence results they provide always
come hand in hand with a very quantitative description of the obtained
solutions. Another advantage is that CAPs typically generalise in a
straightforward way to systems, which is not necessarily the case for some of
the above-mentioned techniques. Given a PDE, abstractly written as $F(u)=0$,
and a numerically obtained approximate solution $\bar{u}$, a common feature of
many CAPs is to consider an associated Newton-like operator of the form

$\displaystyle T(u)=u-AF(u),$ (1)

where $A$ is a well-chosen approximation of $DF(\bar{u})^{-1}$, and then to
prove that $T$ is a contraction on a small neighbourhood of $\bar{u}$, in a
suitable norm [47].

Following the seminal papers of Nakao [34] and Plum [39], many authors
contributed to the further development and popularisation of computer-assisted
proofs for elliptic PDEs (see e.g. [3, 18, 22, 36, 44, 48, 52, 55, 57]), but
the vast majority of these works can only be applied to problems on bounded
domains. Until recently, the main exception was the technique developed by
Plum, which takes $A$ equal to $DF(\bar{u})^{-1}$ and relies on eigenvalues
bounds obtained via the so-called homopty method in order to control
$\|DF(\bar{u})^{-1}\|$. This approach can handle a variety of problems on
unbounded domains, see [54] for a recent example and [35, Part II] for a more
general presentation. Regarding other CAP techniques on unbounded domains, but
for more specific problems, we also mention [31] which studies the spectrum of
coupled harmonic oscillators, and [32] which deals with eigenvalue exclusion
for localised perturbations of a one-dimensional periodic Schrödinger
operator. Lately, other general CAP methodologies were developed for (systems
of) PDEs on unbounded domains, first specifically for radial solutions [50,
51], and then for more general solutions in [13, 14]. In these recent
approaches, one goes back to a problem on a bounded domain, either by
compactifying the domain, or by dealing separately with the solutions at
infinity, or by introducing a periodised problem whose solutions approximate
those of the original problem. Such techniques are well-suited to study
localised solutions of equations having only bounded coefficients.

In this work, we study elliptic PDEs of the form

$\displaystyle\mathcal{L}u:=-\Delta u-\frac{x}{2}\cdot\nabla{u}=f(x,u,\nabla
u),\quad x\in\R^{d},$ (2)

and develop new computer-assisted tools for such equations. Equations of the
form (2) have already been studied extensively, as they naturally appear in
different contexts such as the ones we present below, or when studying self-
similarity in fluids problems (see [2], or [6] in the case of the Navier-
Stokes equations).

As a first example, consider the nonlinear heat equation

$\displaystyle\partial_{t}v=\Delta
v-\varepsilon|v|^{p-1}v,\qquad\varepsilon=\pm 1.$ (3)

If $u$ solves (2) with $f(x,u,\nabla u)=f(u)=u/(p-1)-\varepsilon|u|^{p-1}u$,
i.e.

$\mathcal{L}u=\frac{u}{p-1}-\varepsilon|u|^{p-1}u,$ (4)

then $v(t,x)=t^{-1/(p-1)}u(x/\sqrt{t})$ is a self-similar solution of (3). One
typically seeks “physical” solutions, having Gaussian decay or as critical
points of some energy functional. In the radial case, a thorough description
of solutions can then be achieved via shooting or ODE methods (see [11, 38,
56] and references therein). Alternatively, a variational approach [19, 53]
identifies such solutions by working directly on weighted Sobolev spaces.

As a second example, notice that, if $u(x)=e^{-|x|^{2}/8}\varphi(x)$ for an
arbitrary function $\varphi$, then

$\displaystyle\mathcal{L}u=e^{-|x|^{2}/8}\left(-\Delta\varphi+\left(\frac{|x|^{2}}{16}+\frac{d}{4}\right)\varphi\right),$
(5)

therefore equation (2) is also connected to Schrödinger equations. Indeed,
consider for instance the critical nonlinear Schrödinger equation

$i\partial_{t}v+\Delta v+\varepsilon|v|^{4/d}v=0,$ (6)

for $\varepsilon=\pm 1$. It is established (also via variational methods) in
[27] that (6) has self-similar solutions of the form

$v(t,x)=(1+t^{2})^{-d/4}\exp\left(\frac{it|x|^{2}}{4(1+t^{2})}\right)e^{i\omega\arctan{t}}\varphi\left(\frac{\sqrt{2}x}{(1+t^{2})^{1/2}}\right),$

where $\omega\in\R$ and $\varphi:\R^{d}\longrightarrow\mathbb{C}$ solves

$-\Delta\varphi+\frac{|x|^{2}}{16}\varphi+\frac{\omega}{2}\varphi=\frac{\varepsilon}{2}|\varphi|^{4/d}\varphi.$
(7)

Using the above-mentioned change of unknown $\varphi(x)=e^{|x|^{2}/8}u(x)$, we
get back to equation (2) for $u$, with

$\displaystyle f(x,u,\nabla
u)=f(x,u)=\left(\frac{d}{4}-\frac{\omega}{2}\right)u+\frac{\varepsilon}{2}e^{|x|^{2}/(2d)}|u|^{4/d}u.$

In this work, we develop a computer-assisted approach based on (1) to study
both radial and non-radial solutions of (2), which provides a very precise and
quantitative description of self-similar solutions of equations (3) and (6).
Defining

$H^{n}(\mu)=\left\\{f\in L^{2}(\mu)\mid\partial_{\beta}f\in L^{2}(\mu)\text{
for all }\beta\in\mathbb{N}^{d},|\beta|\leq
n\right\\},\qquad\mu(x)=\frac{\Gamma(d/2)}{(2\sqrt{\pi})^{d}}e^{|x|^{2}/4},$

we will rigorously enclose solutions in $H^{2}(\mu)$ with respect to the norm
$\|\cdot\|_{H^{2}(\mu)}:=\|\mathcal{L}\cdot\|_{L^{2}(\mu)}$. Here are examples
of results that can be obtained using this approach.

######  Theorem 1.

Let $d=2$, $p=5/3$, $\varepsilon=1$ and $\bar{u}:\R^{2}\to\R$ be the radial
function whose radial profile restricted to the interval $[0,4]$ is
represented in Figure 1, and whose precise description is available at [17].
There exists a positive radial solution $u^{\star}\in H^{2}(\mu)$ to Eq. (4)
such that

$\|u^{\star}-\bar{u}\|_{H^{2}(\mu)}\leq 1.6\times 10^{-23}.$

Figure 1: Plot of an approximate self-similar radial profile $\bar{u}$ to Eq.
(3) with $\varepsilon=+1$, $d=2$ and $p=5/3$.

######  Theorem 2.

Let $d=2$, $\varepsilon=-1$, $\omega=-5/2$, and $\bar{\varphi}:\R^{2}\to\R$ be
the function whose restriction to the disk $\\{(x,y)\mid x^{2}+y^{2}\leq
7^{2}\\}$ is represented in Figure 2, and whose precise description is
available at [17]. There exists a solution $\varphi^{\star}\in H^{2}(\R^{2})$
to Eq. (7) such that

$\|\varphi^{\star}-\bar{\varphi}\|_{H^{2}(\R^{2})}\leq 3.8\times 10^{-8}.$

Figure 2: Plot of a non-radial numerical solution
$\bar{\varphi}=e^{r^{2}/8}\bar{u}$ to Eq. (7) with $d=2$, $\varepsilon=-1$ and
$\omega=-5/2$.

More results with different self-similar profiles associated to the nonlinear
heat and Schrödinger equations can be found in Section 4 and Section 5.

######  Remark 3.

While the radial solutions we obtain were already known to exist thanks to
previous works [11, 38, 53], our results also provide a very precise
description of these solutions. Furthermore, our methods can also deal with
non-radial solutions, without relying on a perturbation argument. The only
other work we are aware of where non-radial solutions are found is [19, Remark
3.11], which relies on local bifurcation theory and therefore only holds for
$p$ close to specific values.

######  Remark 4.

Beyond the results themselves, we believe part of the interest of this work
lies in the computer-assisted tools that are developed in order to handle
equations of the form (2), and in the underlying ideas which may prove useful
for CAPs, and more largely for numerical analysis, in a broader context.

  * •

Our CAPs are based on a spectral expansion of the solution, which is a very
common approach on bounded domains, but has proven challenging on unbounded
domains. This is made possible by departing from the usual paradigm which
consists in choosing the basis only in terms of the leading order differential
operator (i.e. $\Delta$ in (2)), so that this operator becomes easy to invert
by hand. Indeed, the basis we use makes
$\mathcal{L}:=-\left(\Delta+\frac{x}{2}\cdot\nabla\right)$ diagonal (and
therefore trivial to invert), but the operator $\Delta$ alone does not have a
well-behaved inverse with respect to this basis.

  * •

Extracting some compactness is usually key for CAPs based on (1), and can be
especially challenging on unbounded domains. Here, compactness is retrieved by
working with appropriately weighted Sobolev spaces (see Section 2.2).

  * •

We use spherical coordinates basis functions of the Hilbert space $H^{2}(\mu)$
which are the counterparts of the spherical harmonics in the Hermite–Laguerre
family and allow us to handle our problems globally and with a unified
treatment. To the best of our knowledge, this is the first occurrence of the
use of these bases to solve partial differential equations, even from a purely
numerical perspective. Indeed, this differs from the previous literature (both
theoretical [19] and computational [21]) on this family of problems which have
exhibited (multivariate) Hermite functions as an eigenbasis for $\mathcal{L}$.
While at first glance, these bases might appear tricky to work with
(especially in nonlinear problems), we demonstrate that they are well-suited
to our problems. We show that they can be handled efficiently to both provide
a spectral method to find numerical solutions to problems of the form (2) and
then validate these solutions rigorously. On the one hand, these bases allow
us to handle non-radial solutions to radial equations. On the other hand, in
the case of radial solutions, these bases allow an approach based on the
analysis of the operator $\mathcal{L}$ while exploiting the dimensionality
reduction given by radial equations.

  * •

As already mentioned, spectral expansions are common in CAPs for nonlinear
equations, but only when they are associated to a Banach algebra structure
which facilitates the control of nonlinearities. The basis used in this work
does not enjoy this property, but we still manage to handle some nonlinear
terms, thanks to a combination of rigorous quadrature and Sobolev embeddings,
in a fashion which is more reminiscent of CAPs based on finite elements [35,
Part I].

  * •

Still regarding nonlinear terms, let us mention that we also handle non-
polynomial nonlinearities, namely a fractional power $p$ in (3). The existing
CAP techniques based on spectral methods are more convenient to use when only
polynomial nonlinearities occur, and therefore the usual approach in the
presence of non-polynomial terms is to first reformulate the equation in a
polynomial one, typically adding several new unknowns and associated equations
in the process [23, 28, 30]. For a different approach allowing to directly
treat some non-polynomial terms without adding extra equations and unknowns,
but relying heavily on a Banach algebra structure, we refer to [9]. In this
work, we directly deal with the fractional power in a somewhat ad-hoc manner
that also does not necessitate the introduction of extra equations and
unknowns.

Our approach can also handle some equations featuring first-order derivatives,
when $f$ in (2) really does depend on $\nabla u$. As an example, we also study
self-similar solutions of a generalised viscous Burgers equation on $\R_{+}$
(see, e.g. [2, 43])

$\partial_{t}v+v^{2}\partial_{x}v=\partial_{xx}v,$ (8)

of the form

$v(t,x)=t^{-1/4}u\left(\frac{x}{\sqrt{t}}\right).$

Then, one has to find a profile $u$ which solves an equation akin to (2) with
$f(x,u,\partial_{x}u)=u/4-u^{2}\partial_{x}u$:

$\displaystyle\partial_{xx}u+\frac{x}{2}\partial_{x}u+\frac{u}{4}-u^{2}\partial_{x}u=0,\quad
x\in\R_{+}.$ (9)

Finally, let us conclude the introduction by emphasising that there have been
previous occurrences of CAPs of somewhat different nature for self-similar
solutions, and that this topic has been attracting a lot of attention recently
for fluid flows [12, 15].

The remainder of the paper is organised as follows. We start by presenting the
main important properties of the operator $\mathcal{L}$, together with
appropriate weighted Sobolev spaces and bases of eigenfunctions in Section
2.1. We show how this setting allows us to recover a crucial Poincaré
inequality and projection error estimates in Section 2.2, and provide some
explicit constants for the Sobolev embeddings between these weighted spaces in
Section 2.3. The CAP methodology and the associated fixed point problem where
all these ingredients will be brought together is then introduced in Section
2.4. Proving the existence of a fixed point requires us to explicitly evaluate
(or at least upper-bound) several quantities associated to an approximate
solution, and we present in Section 3 a quadrature scheme allowing us to
compute these quantities rigorously. In Section 4, we then focus on self-
similar profiles for the nonlinear heat equation. Explicit bounds allowing to
use the CAP in the case of radial solutions and with an integer exponent $p$
are derived in Section 4.1, and used in Section 4.2 to obtain our first set of
results. In Section 4.3 we provide a sufficient condition which we use to
prove the positivity of the obtained solutions. We then show how to deal with
some fractional exponents $p$ in Section 4.4, and with non-radial solutions in
Section 4.5. Self-similar solutions of the nonlinear Schrödinger equations are
treated in an analogous manner in Section 5, and we end in Section 6 with a
generalised viscous Burgers equation. All the computations presented in this
paper have been implemented in Julia, using the IntervalArithmetic.jl library
[40] for interval arithmetic. The computer-assisted parts of the proofs can be
reproduced using the code available at [17].

##  Functional analytic preliminaries

###  Properties of the operator $\mathcal{L}$

One of the main ideas of this work is to make use of the operator

$\mathcal{L}=-\frac{x}{2}\cdot\nabla-\Delta,$

to play the role analogous to that of the Laplacian $-\Delta$ as it is usually
employed for computer-assisted proofs on a bounded domain. Indeed, consider
the weighted space $L^{2}(\mu)$, with scalar product

$\displaystyle\langle u,v\rangle=\int_{\R^{d}}u(x)v(x)\mu(x)\mathrm{d}x,$

where $\mu(x)=e^{|x|^{2}/4}/Z$, $Z=2^{d-1}\omega_{d-1}$ where
$\omega_{d-1}=2\pi^{d/2}/\Gamma(d/2)$ denotes the surface of the
$(d-1)$-dimensional sphere $\mathbb{S}^{d-1}$ (this choice of normalisation is
made in order to match the usual convention in numerical analysis). By
integration by parts, we obtain the Dirichlet form

$\displaystyle\langle\nabla u,\nabla v\rangle=\langle
u,\mathcal{L}v\rangle,\qquad\mbox{for all $u\in
H^{1}(\mu),v\in\mathcal{D}(\mathcal{L})=H^{2}(\mu)$,}$ (10)

which yields that $\mathcal{L}$ is a positive self-adjoint operator on
$L^{2}(\mu)$, and can be analysed in the same “Hilbertian” fashion as
$-\Delta$ is in a first course on elliptic PDEs, as will be seen below. We
refer to [26] for more properties of $\mathcal{L}$, and note that the above
Dirichlet form structure is also often exploited in the analysis of the
generators of reversible diffusions (see for instance [37]).

In the previous literature on the operator $\mathcal{L}$ [19, 21], its
eigenbasis is usually given in cartesian coordinates in terms of Hermite
functions

$\mathcal{L}\Psi_{\beta}=-\left(\frac{x}{2}\cdot\nabla+\Delta\right)\Psi_{\beta}=\frac{d+|\beta|}{2}\Psi_{\beta},\qquad\mbox{where
}\Psi_{\beta}(x)=e^{-|x|^{2}/4}\prod_{i=1}^{d-1}H_{\beta_{i}}(x_{i}/2),\qquad\beta\in\mathbb{N}^{d},$
(11)

where $H_{n}$ denotes the Hermite polynomial of degree $n$. This basis might
be suitable for some problems, but since the equations treated in this work
are radial, we will instead use bases in spherical coordinates (except for
$d=1$).

######  Definition 5 (Spherical coordinates basis).

For $d\in\\{2,3\\}$, define the following orthogonal basis of $L^{2}(\mu)$.

  * •

For $d=2$, $m\in\mathbb{N}$, $l\in\mathbb{Z}$, $x_{1}=r\cos\vartheta$,
$x_{2}=r\sin\vartheta$,

$\psi_{l,m}(r,\vartheta)=e^{-r^{2}/4}(r/2)^{|l|}L_{m}^{(|l|)}(r^{2}/4)\begin{cases}\cos
l\vartheta\qquad&\mbox{if $l\geq 0$}\\\ -\sin l\vartheta\qquad&\mbox{if
$l<0$,}\end{cases}$

where $L^{(\alpha)}_{m}$ denotes the generalised Laguerre polynomial of degree
$m$ and parameter $\alpha$ (see [1, p. 775]).

  * •

For $d=3$, $m\in\mathbb{N}$, $k\in\mathbb{Z}$, $l\in\mathbb{N}_{\geq|k|}$,
$x_{1}=r\sin\vartheta\cos\phi$, $x_{2}=r\sin\vartheta\sin\phi$,
$x_{3}=r\cos\vartheta$

$\psi_{k,l,m}(r,\vartheta,\phi)=e^{-r^{2}/4}(r/2)^{l}L_{m}^{(l+1/2)}(r^{2}/4)P^{k}_{l}(\cos\vartheta)\begin{cases}\cos
k\phi\qquad&\mbox{if $k\geq 0$}\\\ -\sin k\phi\qquad&\mbox{if
$k<0$,}\end{cases}$

where $P^{k}_{l}$ denotes the associated Legendre polynomial of order $k$ and
degree $l$ (see [1, Chapter 8]).

######  Remark 6.

For $d=2$, the $\psi_{l,m}$’s are weighted and reformulated versions of the
so-called Hermite polynomials of complex variables [25]. See [16] for the
derivations of some of these relations.

######  Lemma 7.

The families $\\{\psi_{l,m}\\}$ and $\\{\psi_{k,l,m}\\}$ from Definition 5 are
eigenfunctions of $\mathcal{L}$ for $d=2$ and $d=3$ respectively.

  * •

For $d=2$, $m\in\mathbb{N}$, $l\in\mathbb{Z}$,

$\mathcal{L}\psi_{l,m}=-\left(\frac{r}{2}\frac{\partial\psi_{l,m}}{\partial
r}+\frac{1}{r}\frac{\partial}{\partial
r}\left(r\frac{\partial\psi_{l,m}}{\partial
r}\right)+\frac{1}{r^{2}}\frac{\partial^{2}\psi_{l,m}}{\partial\vartheta^{2}}\right)=\left(\frac{2+|l|}{2}+m\right)\psi_{l,m}.$

  * •

For $d=3$, $m\in\mathbb{N}$, $k\in\mathbb{Z}$, $l\in\mathbb{N}_{\geq|k|}$,

$\displaystyle\mathcal{L}\psi_{k,l,m}$
$\displaystyle=-\left(\frac{r}{2}\frac{\partial\psi_{k,l,m}}{\partial
r}+\frac{1}{r^{2}}\frac{\partial}{\partial
r}\left(r^{2}\frac{\partial\psi_{k,l,m}}{\partial
r}\right)+\frac{1}{r^{2}\sin{\vartheta}}\frac{\partial}{\partial\vartheta}\sin{\vartheta}\left(\frac{\partial\psi_{k,l,m}}{\partial\vartheta}\right)+\frac{1}{r^{2}\sin^{2}{\vartheta}}\frac{\partial^{2}\psi_{k,l,m}}{\partial\phi^{2}}\right)$
$\displaystyle=\left(\frac{3+l}{2}+m\right)\psi_{k,l,m}.$

The construction of the bases given in Definition 5 hinges on the following
lemma.

######  Lemma 8.

Let $f$ be a hyperspherical function on $\mathbb{S}^{d-1}$ of order
$l\in\mathbb{N}$, i.e.

$\Delta_{\mathbb{S}^{d-1}}f=-l(l+d-2)f,$

where $\mathbb{S}^{d-1}$ denotes the Laplace-Beltrami operator on
$\mathbb{S}^{d-1}$. Then $e^{-r^{2}/4}(r/2)^{l}L_{m}^{(l+d/2-1)}(r^{2}/4)f$ is
an eigenfunction of $\mathcal{L}$ with eigenvalue $(d+l)/2+m$ on $\R^{d}$.

###### Proof.

This follows from separation of the radial variable and the spherical
variables and by identifying the Laguerre differential equation after the
change of variable $z=r^{2}/4$. ∎

In order to obtain the bases of Definition 5, it then suffices to have
orthogonal bases of the eigenspaces of $\Delta_{\mathbb{S}^{d-1}}$, and this
construction can also be used in higher dimensions (e.g. by using [24] for the
construction of an eigenbasis of $\Delta_{\mathbb{S}^{d-1}}$). Note that as a
consequence, this directly generalises the known corresponding bases for the
harmonic Schrödinger operator and for the generator of the radial
Ornstein–Uhlenbeck process. Furthermore, the linear relationship between
bivariate Hermite and Laguerre polynomials observed in [16] for $d=2$ (and the
one known for $d=1$) can be generalised by the above pattern.

######  Corollary 9 (Radial eigenbasis).

For $d\in\mathbb{N}$, let $\alpha=d/2-1$ and let
$\psi_{n}(r):=L^{(\alpha)}_{n}(r^{2}/4)e^{-r^{2}/4}$, then

$\mathcal{L}\psi_{n}=-\left(\frac{r}{2}+\frac{d-1}{r}\right)\frac{\partial\psi_{n}}{\partial
r}-\frac{\partial^{2}\psi_{n}}{\partial
r^{2}}=\left(\frac{d}{2}+n\right)\psi_{n}.$

Note that with these bases, we can handle radial equations globally without
treating either the behaviour at infinity or at the origin separately (see for
instance [51] for an example of CAP where the behaviours at the origin and at
infinity are handled separately).

###  Spaces, norms and Poincaré inequality

One of the main challenges of computer-assisted proofs on unbounded domains is
that the embedding $H^{1}(\Omega)\hookrightarrow L^{2}(\Omega)$ is not compact
anymore for $\Omega=\R^{d}$, as no Poincaré inequality even holds. Observe
that in our setting with weighted spaces, the eigenbasis of $\mathcal{L}$
together with (10) immediately implies a genuine Poincaré inequality

$\left\|u\right\|_{L^{2}(\mu)}\leq\sqrt{\frac{2}{d}\langle\nabla u,\nabla
u\rangle}=\sqrt{\frac{2}{d}}\left\|\nabla
u\right\|_{L^{2}(\mu)}\qquad\mbox{for all $u\in H^{1}(\mu)$},$ (12)

since $d/2$ is the smallest eigenvalue of $\mathcal{L}$. We can therefore
consider on $H^{1}(\mu)$ the scalar product
$\langle\cdot,\cdot\rangle_{H^{1}(\mu)}=\langle\nabla\cdot,\nabla\cdot\rangle$
and its associated norm, which is equivalent to the usual $H^{1}(\mu)$-norm.
Similarly, on the space $H^{2}(\mu)$ we consider the scalar product
$\langle\cdot,\cdot\rangle_{H^{2}(\mu)}=\langle\mathcal{L}\cdot,\mathcal{L}\cdot\rangle$
and associated norm, which is also equivalent to the usual $H^{2}(\mu)$-norm
(see Lemma 16).

Notice that these choices of norms yield, in a similar fashion to (12),

$\|u\|_{H^{1}(\mu)}\leq\sqrt{\frac{2}{d}}\|u\|_{H^{2}(\mu)}\qquad\mbox{for all
$u\in H^{2}(\mu)$}.$ (13)

######  Notation 10.

We denote by $\|u\|_{\infty}$ the essential supremum of $|u|$, i.e. for a
measure $\nu$ on $\Omega$,

$\|u\|_{\infty}:=\inf\left\\{M\mid|u(x)|\leq M\textrm{ for $\nu$-a.e.
$x\in\Omega$}\right\\}.$

Our notation does not depend on $\nu$, as all the measures we will be using
are equivalent (i.e. they have the same null sets). Finally, we adopt the
usual convention that $u$ is said to be continuous if it has a continuous
representative.

Another important feature of our setup is that that
$\mathcal{L}^{-1}:L^{2}(\mu)\to H^{1}(\mu)$ is well-defined and compact. More
quantitatively, if $P_{n}$ denotes the orthogonal projection in $L^{2}(\mu)$
onto $\mathrm{Span}\\{\psi\mid\mathcal{L}\psi=\lambda\psi\textrm{ with
}\lambda\leq d/2+n\\}$, (10) yields the following error estimates

$\displaystyle\|\mathcal{L}^{-1}(f-P_{n}f)\|_{H^{1}(\mu)}\leq\sqrt{\frac{2}{d+2n+1}}\|f-P_{n}f\|_{L^{2}(\mu)}\qquad\mbox{for
all $f\in L^{2}(\mu)$},$ (14)

and

$\displaystyle\left\|f-P_{n}f\right\|_{L^{2}(\mu)}\leq\frac{2}{d+2n+1}\left\|f-P_{n}f\right\|_{H^{2}(\mu)}\qquad\mbox{for
all $f\in H^{2}(\mu)$},$ (15)

which are critical for our CAPs, and in particular for controlling the quality
of the approximate inverse $A$ that will be defined in Section 2.4.

######  Remark 11.

Let us highlight two evident but useful properties of the orthogonal projector
$P_{n}$. First, $P_{n}$ commutes with $\mathcal{L}$ (and $\mathcal{L}^{-1}$).
Second $P_{n}$ is also an orthogonal projection with respect to the
$H^{1}(\mu)$ and $H^{2}(\mu)$ scalar products. In particular,
$P_{n}(L^{2}(\mu))=P_{n}(H^{1}(\mu))=P_{n}(H^{2}(\mu))$.

###  Sobolev embeddings

To control the nonlinearities in our setting, Sobolev embeddings play a
crucial role. We start by recalling some $H^{1}(\mu)$ embedding proved in [19,
29, 53].

######  Theorem 12.

We have the following Sobolev embeddings:

  * (i)

If $d=1$, $u\in H^{1}(\mu)\implies e^{x^{2}/8}u\in\mathcal{C}^{0,1/2}(\R)$.

  * (ii)

If $d=2$, $H^{1}(\mu)\subset L^{p}(\mu)$ for all $p\in[2,\infty)$.

  * (iii)

If $d\geq 3$, $H^{1}(\mu)\subset L^{p}(\mu)$ for all $p\in[2,2d/(d-2)]$.

While in principle, this shows the feasibility of a computer-assisted proof
using this functional framework, in practice we also require the corresponding
embedding constants that we derive below.

######  Notation 13.

In this section, we make use of classical embeddings between Sobolev spaces
without weights in order to derive the needed results on weighted spaces. When
we include $\R^{d}$ in the notation, e.g. $L^{2}(\R^{d})$, this means we refer
to the space without weight.

######  Lemma 14.

Let $\theta(x)=|x|^{2}/4$ and denoting $\|u\|_{H^{1}(e^{\theta})}:=\|\nabla
u\|_{L^{2}(e^{\theta})}$, then for all $u\in H^{1}(e^{\theta})$

$\|u\|_{L^{p}(e^{\theta})}\leq\|e^{\theta(1/2-1/p)}u\|_{L^{p}(e^{\theta})}\leq
C(d,p)\|u\|^{a}_{L^{2}(e^{\theta})}\|u\|^{1-a}_{H^{1}(e^{\theta})},$ (16)

where $a=1+d(1/p-1/2)$ and:

  * •

for $d=1$, $p\in[2,\infty]$ and $C(1,p)=2^{1/2-1/p}$;

  * •

for $d=2$, $p\in[2,\infty)$ and, writing $p=2(k+r)$ with $k\in\mathbb{N}^{*}$
and $r\in[0,1)$,

$C(2,2(k+r))=\sqrt{2}\left(\frac{k!}{\sqrt{2}}\right)^{1/(k+r)}(k+1)^{r/(k+r)};$

  * •

for $d\geq 3$, $p\in[2,2d/(d-2)]$ and

$C(d,p)=\left(\frac{1}{d(d-2)\pi}\right)^{(1-a)/2}\left\\{\frac{(d-1)!}{\Gamma(d/2)}\right\\}^{1/2-1/p}.$

###### Proof.

First, notice that $\theta(1/2-1/p)\geq 0$ since $p\geq 2$, which yields the
first inequality in (16). Note that by a density argument, as in [19, Lemma
4.17], for all $u\in H^{1}(e^{\theta})$ and $\varphi=e^{|x|^{2}/8}u$, we have
that

$-\frac{1}{2}\int_{\R^{d}}\left(x\cdot\nabla\varphi\right)\varphi\mathrm{d}x=\frac{d}{4}\int_{\R^{d}}\varphi^{2}\mathrm{d}x.$
(17)

And thus,

$\|\nabla
u\|_{L^{2}(e^{\theta})}^{2}=\int_{\R^{d}}\left|\nabla\varphi-\frac{x}{4}\varphi\right|^{2}\mathrm{d}x=\int_{\R^{d}}\left(\left|\nabla\varphi\right|^{2}-\frac{\varphi}{2}\left(x\cdot\nabla\varphi\right)+\frac{|x|^{2}}{16}\varphi^{2}\right)\mathrm{d}x=\int_{\R^{d}}\left[\left|\nabla\varphi\right|^{2}+\varphi^{2}\left(\frac{d}{4}+\frac{|x|^{2}}{16}\right)\right]\mathrm{d}x$

which yields

$\|\nabla\varphi\|_{L^{2}(\R^{d})}\leq\|\nabla u\|_{L^{2}(e^{\theta})}.$ (18)

  * •

For $d=1$, from the proof of [10, Theorem VIII.7], we have that

$\|\varphi\|_{L^{\infty}(\R)}\leq
2^{1/2}\|\varphi\|_{L^{2}(\R)}^{1/2}\|\varphi^{\prime}\|_{L^{2}(\R)}^{1/2}.$

Therefore,

$\|e^{\theta/2}u\|_{\infty}\leq
2^{1/2}\|u\|_{L^{2}(e^{\theta})}^{1/2}\|u\|_{H^{1}(e^{\theta})}^{1/2},$ (19)

which is the announced inequality for $p=\infty$. For $2\leq p<\infty$,
Hölder’s inequality yields

$\displaystyle\|e^{\theta(1/2-1/p)}u\|_{L^{p}(e^{\theta})}$
$\displaystyle\leq\left(\|e^{\theta/2}u\|_{\infty}^{p-2}\|u\|^{2}_{L^{2}(e^{\theta})}\right)^{1/p}$
$\displaystyle=\|e^{\theta/2}u\|_{\infty}^{1-2/p}\|u\|^{2/p}_{L^{2}(e^{\theta})}$
$\displaystyle\leq
2^{1/2-1/p}\|u\|^{1/2+1/p}_{L^{2}(e^{\theta})}\|u\|^{1/2-1/p}_{H^{1}(e^{\theta})}.$

  * •

Now, for $d=2$, we are in the “limiting case”, so from [10, Corollary IX.11],
for $k\in\mathbb{N}^{*}$

$\|\varphi\|_{L^{2k}(\mathbb{R}^{2})}\leq
k^{1/k}\|\varphi\|_{L^{2(k-1)}(\mathbb{R}^{2})}^{(k-1)/k}\|\nabla\varphi\|_{L^{2}(\mathbb{R}^{2})}^{1/k},$

and by induction on $k\in\mathbb{N}^{*}$

$\|\varphi\|_{L^{2k}(\R^{2})}\leq(k!)^{1/k}\|\varphi\|_{L^{2}(\R^{2})}^{1/k}\|\nabla\varphi\|_{L^{2}(\R^{2})}^{1-1/k}.$

Now, by interpolation, for $2k\leq 2(k+r)\leq 2(k+1)$, we get

$\|\varphi\|_{L^{2(k+r)}(\R^{2})}\leq(k!)^{1/(k+r)}(k+1)^{r/(k+r)}\|\varphi\|_{L^{2}(\R^{2})}^{1/(k+r)}\|\nabla\varphi\|^{1-1/(k+r)}_{L^{2}(\R^{2})},$

and thus

$\displaystyle\|e^{\theta(k+r-1)/(2(k+r))}u\|_{L^{2(k+r)}(e^{\theta})}=\|\varphi\|_{L^{2(k+r)}(\R^{2})}\leq\sqrt{2}\left(\frac{k!}{\sqrt{2}}\right)^{1/(k+r)}(k+1)^{r/(k+r)}\|u\|_{L^{2}(e^{\theta})}^{1/(k+r)}\|u\|^{1-1/(k+r)}_{H^{1}(e^{\theta})}.$

  * •

For $d\geq 3$, proceed similarly and use [45] (see also [19, Lemma 4.17]). ∎

######  Corollary 15.

For all $u\in H^{1}(\mu)$,

$\|u^{p}\|_{L^{2}(\mu)}\leq\|e^{(p-1)|x|^{2}/8}u^{p}\|_{L^{2}(\mu)}\leq
c(d,p)\|u\|^{p}_{H^{1}(\mu)},$

where

  * •

for $d=1$, $p\in[1,\infty)$, $c(1,p)=2^{(5p-3)/4}$;

  * •

for $d=2$, $p=k+r\in[1,\infty)$ with $k\in\mathbb{N}^{*}$ and $r\in[0,1)$,
$c(2,k+r)=(2\sqrt{2\pi})^{k+r-1}k!(k+1)^{r}$;

  * •

for $d\geq 3$, $p\in[1,d/(d-2)]$ and

$c(d,p)=\left(\frac{2}{d}\right)^{p/2}\left\\{\left(\frac{2}{d-2}\right)^{d/4}\frac{\Gamma(d)^{1/2}}{\Gamma(d/2)}\right\\}^{p-1}.$

###### Proof.

Noting that $\mu=Z^{-1}e^{\theta}$ and applying Lemma 14 followed by the
Poincaré inequality (12) yields

$\displaystyle\|u^{p}\|_{L^{2}(\mu)}=\frac{1}{\sqrt{Z}}\|u\|^{p}_{L^{2p}(e^{\theta})}$
$\displaystyle\leq\frac{\left(C(d,2p)\sqrt{Z}\right)^{p}}{\sqrt{Z}}\left(\frac{1}{\sqrt{\lambda_{0}}}\right)^{p+d(1/2-p/2)}\|u\|^{p}_{H^{1}(\mu)}$
$\displaystyle=\left(C(d,2p)\right)^{p}\left(\frac{2^{d}\pi^{d/2}}{\Gamma(d/2)}\right)^{\frac{p-1}{2}}\left(\frac{2}{d}\right)^{p/2+d(1/4-p/4)}\|u\|^{p}_{H^{1}(\mu)}.$

One then just has to check that

$c(d,p)=\left(C(d,2p)\right)^{p}\left(\frac{2^{d}\pi^{d/2}}{\Gamma(d/2)}\right)^{\frac{p-1}{2}}\left(\frac{2}{d}\right)^{p/2+d(1/4-p/4)}.\qed$

While the $H^{1}(\mu)\to L^{2}(\mu)$ embeddings will prove useful for some
intermediate estimates, our CAP will mainly be conducted in $H^{2}(\mu)$. We
derive the required embeddings with explicit constants below.

######  Lemma 16.

For all $f\in L^{2}(\mu)$, there exists a unique $u\in H^{2}(\mu)$ such that
$\mathcal{L}u=f$, and $\mathcal{L}^{-1}:L^{2}(\mu)\rightarrow H^{2}(\mu)$ is
bounded, so that $\|\cdot\|_{H^{2}(\mu)}:=\|\mathcal{L}\cdot\|_{L^{2}(\mu)}$
indeed induces the $H^{2}(\mu)$ topology. Furthermore, if $u\in H^{2}(\mu)$
and $\varphi=e^{|x|^{2}/8}u$, then $\varphi\in H^{2}(\R^{d})$ and

$\|D^{2}\varphi\|_{L^{2}(\R^{d})}=\|\Delta\varphi\|_{L^{2}(\R^{d})}\leq\sqrt{Z}\|u\|_{H^{2}(\mu)},$
(20)

where $\|D^{2}\varphi\|_{L^{2}(\R^{d})}$ denotes the Frobenius norm of the
Hessian $D^{2}\varphi$ of $\varphi$.

###### Proof.

The qualitative version of this lemma is the content of [26, Lemma 2.1.(vii)]
(The additional fact that $\varphi\in H^{2}(\R^{d})$ follows from [10, Theorem
IX.25]). Given the bounds obtained in the proof of [26, Lemma 2.1.(vii)], for
all $u\in H^{2}(\mu)$ the estimates below hold by density for
$\varphi=e^{|x|^{2}/8}u$. Using (5), we have that

$\displaystyle\int_{\R^{d}}\left(\mathcal{L}u\right)^{2}e^{|x|^{2}/4}\mathrm{d}x$
$\displaystyle=\int_{\R^{d}}\left[\left|\Delta\varphi\right|^{2}+\left(\frac{d}{4}+\frac{|x|^{2}}{16}\right)\varphi^{2}\right]\mathrm{d}x-\int_{\R^{d}}\varphi\Delta\varphi\left(\frac{d}{2}+\frac{|x|^{2}}{8}\right)\mathrm{d}x.$

Now, by integration by parts, the last term simplifies to

$\displaystyle-\int_{\R^{d}}\varphi\Delta\varphi\left(\frac{d}{2}+\frac{|x|^{2}}{8}\right)\mathrm{d}x$
$\displaystyle=\int_{\R^{d}}\left|\nabla\varphi\right|^{2}\left(\frac{d}{2}+\frac{|x|^{2}}{8}\right)\mathrm{d}x+\int_{\R^{d}}\varphi\left(\frac{x}{4}\cdot\nabla\varphi\right)\mathrm{d}x$
$\displaystyle=\int_{\R^{d}}\left|\nabla\varphi\right|^{2}\left(\frac{d}{2}+\frac{|x|^{2}}{8}\right)\mathrm{d}x-\frac{d}{8}\int_{\R^{d}}\varphi^{2}\mathrm{d}x\qquad\mbox{using\leavevmode\nobreak\
\eqref{eqn:useful-id}.}$

Adding back this term, we get

$\int_{\R^{d}}\left(\mathcal{L}u\right)^{2}e^{|x|^{2}/4}\mathrm{d}x\geq\int_{\R^{d}}\left|\Delta\varphi\right|^{2}\mathrm{d}x+\frac{d}{2}\int_{\R^{d}}\left|\nabla\varphi\right|^{2}\mathrm{d}x+\frac{d}{8}\int_{\R^{d}}\varphi^{2}\mathrm{d}x,$
(21)

and in particular

$\|D^{2}\varphi\|_{L^{2}(\R^{d})}=\|\Delta\varphi\|_{L^{2}(\R^{d})}\leq\|\mathcal{L}u\|_{L^{2}(e^{|x|^{2}/4})}\leq\sqrt{Z}\|u\|_{H^{2}(\mu)},$

where the first equality follows from [42, Lemma 2.5, Corollary 2.13], since
$\varphi\in H^{2}(\mathbb{R}^{d})$. ∎

Applying a similar argument as in Lemma 14, we then obtain the following
embeddings.

######  Corollary 17.

We have the following Sobolev embeddings:

  * (i)

If $d\in\\{1,2,3\\}$, $u\in H^{2}(\mu)\implies e^{|x|^{2}/8}u\in
L^{\infty}(\R^{d})$.

  * (ii)

If $d=4$, $H^{2}(\mu)\subset L^{p}(\mu)$ for all $p\in[2,\infty)$.

  * (iii)

If $d\geq 5$, $H^{2}(\mu)\subset L^{p}(\mu)$ for all $p\in[2,2d/(d-4)]$.

###### Proof.

The proof follows from the application of the Sobolev embeddings given by [10,
Corollary IX.13] on $\varphi=e^{|x|^{2}/8}u\in H^{2}(\R^{d})$. ∎

######  Remark 18.

Note that these embeddings are sufficient to cover the range of exponent $p$
for the semilinear heat equation for which rapidly decaying self-similar
solutions are known to exist [19].

We now provide an explicit embedding constant for the case that we will need
for our examples, namely (i) in Corollary 17, for $d\in\\{2,3\\}$ (for $d=1$,
one can use the $H^{1}(\mu)$-embedding from Corollary 15 together with (13)).

######  Lemma 19.

Let $d\in\\{2,3\\}$, then there exists $C(d)>0$ such that for all $u\in
H^{2}(\mu)$

$\|e^{|x|^{2}/8}u\|_{\infty}\leq C(d)\|u\|_{H^{2}(\mu)},$

where $C(d)$ can be chosen as

$C(d)=\sqrt{Z}\left(\frac{2}{d}C_{0}+\sqrt{\frac{2}{d}}C_{1}+C_{2}\right),$

with

$C_{0}=0.56419$, |  $C_{1}=0.79789$, |  $C_{2}=0.23033$, | if $d=2;$
---|---|---|---
$C_{0}=0.52319$, |  $C_{1}=1.0228$, |  $C_{2}=0.37467$, | if $d=3.$

###### Proof.

Letting once again $\varphi=e^{|x|^{2}/8}u$, we have that $\varphi\in
H^{2}(\R^{d})$ by Lemma 16. We can thus use the embedding from $H^{2}(\R^{d})$
to $L^{\infty}(\R^{d})$, with the explicit constants provided in [35, Example
6.12.(a)]:

$\displaystyle\|\varphi\|_{L^{\infty}(\R^{d})}\leq
C_{0}\|\varphi\|_{L^{2}(\R^{d})}+C_{1}\|\nabla\varphi\|_{L^{2}(\R^{d})}+C_{2}\|D^{2}\varphi\|_{L^{2}(\R^{d})}.$

Then, noticing that $\|u\|_{L^{2}(\mu)}\leq\frac{2}{d}\|u\|_{H^{2}(\mu)}$ (by
a direct computation or by combining (12) and (13)), we get

$\displaystyle\|\varphi\|_{L^{2}(\R^{d})}=\sqrt{Z}\|u\|_{L^{2}(\mu)}\leq\sqrt{Z}\frac{2}{d}\|u\|_{H^{2}(\mu)}.$

Similarly, starting from (18) and using once more (13), we get

$\displaystyle\|\nabla\varphi\|_{L^{2}(\R^{d})}\leq\sqrt{Z}\|u\|_{H^{1}(\mu)}\leq\sqrt{Z}\sqrt{\frac{2}{d}}\|u\|_{H^{2}(\mu)},$

and (20) provides the last missing inequality. ∎

###  The fixed-point theorem

As is common in the CAP literature, our approach will rely on a fixed-point
theorem allowing for the _a posteriori validation_ of an appropriate numerical
solution $\bar{u}$, that is, a proof of the existence of a solution
$u^{\star}$ to

$\mathcal{L}u=f(x,u,\nabla u),$ (22)

in an explicit neighbourhood of $\bar{u}$ in a Banach space $\mathcal{X}$.
This is achieved by constructing a well-chosen map
$T:\mathcal{X}\to\mathcal{X}$ with solutions of Eq. (22) as its fixed points
and which is contracting in a neighbourhood of $\bar{u}\in\mathcal{X}$. In
order to actually prove that $T$ is a contraction, we will use the following
statements (see [35, 49, 55] and the references therein for variants).

######  Theorem 20.

Let $(\mathcal{X},\|\cdot\|_{\mathcal{X}})$ be a Banach space, let
$T:\mathcal{X}\to\mathcal{X}$ be a $\mathcal{C}^{1}$-map and let
$\bar{u}\in\mathcal{X}$ be such that there exist a constant $Y\geq 0$ and an
increasing continuous function $\mathcal{Z}:\mathbb{R}_{+}\to\R_{+}$ such that

$\displaystyle\|T(\bar{u})-\bar{u}\|_{\mathcal{X}}$ $\displaystyle\leq Y,$
(23) $\displaystyle\|DT(u)\|_{\mathcal{X},\mathcal{X}}$
$\displaystyle\leq\mathcal{Z}(\|u-\bar{u}\|_{\mathcal{X}})\qquad\mbox{for all
$u\in\mathcal{X}$}$ (24)

If there exists $\delta>0$ such that

$\displaystyle Y+\int_{0}^{\delta}\mathcal{Z}(s)\mathrm{d}s$
$\displaystyle<\delta,$ (25) $\displaystyle\mathcal{Z}(\delta)<1,$ (26)

then $T$ has a unique fixed point $u^{\star}$ in $\bar{B}(\bar{u},\delta)$.

###### Proof.

See for instance [8, Theorem 2.1] ∎

If $T$ is polynomial in $u$ of degree $p$, we reformulate condition (24) into
bounding $D^{k}T(\bar{u})$ for $k\in\\{1,\ldots,p\\}$.

######  Corollary 21.

Let $T:\mathcal{X}\to\mathcal{X}$ be a polynomial of degree $p$ and let
$\bar{u}\in\mathcal{X}$ be such that there exist $Y,Z_{1},\ldots,Z_{p}\geq 0$
be such that

$\displaystyle\|T(\bar{u})-\bar{u}\|_{\mathcal{X}}$ $\displaystyle\leq Y,$
(27) $\displaystyle\|D^{k}T(\bar{u})\|_{\mathcal{X},\mathcal{X}}$
$\displaystyle\leq Z_{k}\qquad\mbox{for all $k\in\\{1,\ldots,p\\}$.}$ (28)

Introducing the radii polynomials

$\displaystyle P(\delta)$
$\displaystyle=Y-\delta+\sum_{k=1}^{p}\frac{Z_{k}}{k!}\delta^{k},$
$\displaystyle Q(\delta)$
$\displaystyle=-1+\sum_{k=1}^{p}\frac{Z_{k}}{(k-1)!}\delta^{k-1},$

if there exists $\delta>0$ such that $P(\delta)>0$, then denoting
$\underline{\delta}$ the smallest positive root of $P$ and $\bar{\delta}$ the
unique positive root of $Q$, for all
$\delta\in(\underline{\delta},\bar{\delta})$ $T$ has a unique fixed point
$u^{\star}\in\bar{B}(\bar{u},\delta)\subset\mathcal{X}$.

###### Proof.

See for instance [8, Corollary 4.4]. ∎

Provided $g:u\mapsto f(\cdot,u,\nabla u)$ maps $H^{2}(\mu)$ to $L^{2}(\mu)$,
our general problem

$\mathcal{L}u=f(x,u,\nabla u)$

can be reformulated as finding a zero of the map $F:H^{2}(\mu)\to H^{2}(\mu)$,
where

$F(u)=u-\mathcal{L}^{-1}f(x,u,\nabla u).$

Note that if $Dg(\bar{u})$ in fact maps $H^{2}(\mu)$ to $H^{1}(\mu)$, then
$DF(\bar{u})$ is a compact perturbation of the identity. This motivates the
introduction of the linear map $A:H^{2}(\mu)\to H^{2}(\mu)$ such that
$P_{n}AP_{n}=A_{n}$ is an approximate inverse of $P_{n}DF(\bar{u})P_{n}$
computed numerically, and where $P_{n}A(I-P_{n})=(I-P_{n})AP_{n}=O$ and
$(I-P_{n})A(I-P_{n})=I$. Recalling that $P_{n}$ denotes the orthogonal
projection in $L^{2}(\mu)$ onto
$\mathrm{Span}\\{\psi\mid\mathcal{L}\psi=\lambda\psi\textrm{ with }\lambda\leq
d/2+n\\}$, in view of (15) this operator $A$ should be an accurate approximate
inverse of $DF(\bar{u})$ if $n$ is large enough, and this statement will be
made quantitative in Section 4.1.2. We therefore consider the quasi-Newton
operator $T:H^{2}(\mu)\longrightarrow H^{2}(\mu)$ with

$T(u)=u-AF(u).$

Provided $A$ is invertible (this will be an automatic consequence of verifying
assumption (26)), zeros of $F$ are indeed in one-to-one correspondence with
fixed points of $T$. In Sections 4 to 6, we derive suitable estimates $Y$ and
$\mathcal{Z}$ (or $Z_{k}$) allowing us to apply Theorem 20 or Corollary 21 to
this fixed point $T$ for the different equations mentioned in the
introduction.

######  Remark 22.

One could choose to perform a computer-assisted proof in $H^{1}(\mu)$ instead
of $H^{2}(\mu)$, but the class of treatable problems would be more restricted
and this would not directly yield $L^{\infty}$-estimates on the solution in
dimension $d\in\\{2,3\\}$.

##  Handling nonlinearities via quadrature rules

For deriving the bounds $Y$ and $\mathcal{Z}$ needed for our CAP, we will need
to rigorously compute quantities like $\|\bar{u}^{p}\|_{L^{2}(\mu)}$. In the
more classical case of Fourier series, $\bar{u}$ is typically a finite Fourier
series (i.e. a trigonometric polynomial), and $\bar{u}^{p}$ is then still a
finite Fourier series (provided $p$ is an integer), whose integral or norm is
therefore straightforward to compute exactly. In contrast, in our setting,
whenever $u\neq 0$ and $v\neq 0$ belong to $P_{n}(H^{2}(\mu))$, the product
$uv$ does not belong to $P_{n^{\prime}}(H^{2}(\mu))$ for any
$n^{\prime}\in\mathbb{N}$. However, we can still compute exactly
$\|uv\|_{L^{2}(\mu)}$ (or $\|\bar{u}^{p}\|_{L^{2}(\mu)}$), and we explain how
this can be done in this section.

For clarity, we focus here on the case of the radial equation

$\mathcal{L}u=u/(p-1)-\varepsilon u^{p},\quad\varepsilon=\pm 1,\
p\in\mathbb{N},$

giving rise to self-similar solutions of the semilinear heat equation, but all
the constructions can be adapted to the other examples. We denote by
$\\{\hat{\psi}_{m}\\}_{m=0}^{\infty}$ the normalisation of the orthogonal
basis $\\{\psi_{m}\\}_{m=0}^{\infty}$ (see Corollary 9) with respect the
$L^{2}(\mu)$-inner product $\langle\cdot,\cdot\rangle$, i.e.

$\hat{\psi}_{m}(r)=\sqrt{\frac{m!}{\Gamma(m+\alpha+1)}}L^{(\alpha)}_{m}(r^{2}/4),\qquad\alpha=\frac{d}{2}-1.$

We will solve our problem with respect to this basis and look for our
numerical approximation $\bar{u}$ in
$\mathrm{Span}\\{\hat{\psi}_{m}\\}_{m=0}^{n}$ (which is the set of radial
functions in $P_{n}(H^{2}(\mu))$). We can then write $\bar{u}$ in terms of its
“Fourier” coefficients:

$\bar{u}=\sum_{m=0}^{n}\bar{u}_{m}\hat{\psi}_{m}\qquad\bar{\mathbf{u}}=(\bar{u}_{0},\ldots,\bar{u}_{n}).$
(29)

Since our computer-assisted proof will mostly rely on a Hilbertian approach,
one of the main challenges is to be able to compute many projections and
integrals (eventually rigorously) on unbounded domains. Such endeavour might
appear naive given the lack of simple error control on Gauss–Laguerre
quadrature and the notorious difficulties of working with this family of bases
[20]. Furthermore, we would ideally like to solve (via gradient descent or
Newton’s method) the Galerkin problem

$\left\langle\mathcal{L}\bar{u}-\bar{u}/(p-1)+\varepsilon\bar{u}^{p},\hat{\psi}_{j}\right\rangle=0\qquad\mbox{for
all $j\in\\{0,\ldots,n\\}$},$ (30)

without further approximating the problem. We show that this can in fact be
achieved relatively efficiently for a nonlinearity of moderate degree by
simply adapting the Gauss–Laguerre quadrature. Consider more generally the
problem of evaluating (eventually rigorously) integrals of the form

$\int_{\R^{d}}\prod_{k=1}^{m}f_{k}(x)\mu(x)\mathrm{d}x,$

for $f_{k}\in\mathrm{Span}\\{\hat{\psi}_{m}\\}_{m=0}^{n}$. For instance,
$\|\bar{u}^{p}\|_{L^{2}(\mu)}$ can be written in this form, but also the
coefficients of $P_{n}(\bar{u}^{p})$ given by
$\langle\bar{u}^{p},\hat{\psi}_{m}\rangle$ for $m=0,\ldots,n$. Since we are in
the radial setting, this boils down to evaluating integrals of the form

$\int_{0}^{\infty}\prod_{k=1}^{m}f_{k}(r)(r/2)^{d-1}e^{r^{2}/4}\mathrm{d}r,$
(31)

with $f_{k}(r)=p_{k}(r^{2}/4)e^{-r^{2}/4}$ for some polynomial $p_{k}$ of
degree at most $n$. Note that all calculations can be thought to be done with
respect to the variable $z=r^{2}/4$.

######  Remark 23.

When we have to deal with non-radial solutions, we also use spherical
coordinates, and the projections and integrals with respect to the spherical
variables only involve sines and cosines and thus do not require a quadrature.
This is another advantage of using spherical coordinates instead of cartesian
ones (i.e., of using the bases of Definition 5 instead of (11)), for which
nontrivial quadratures would be required in every direction.

We now show how to exactly compute integrals of the form (31).

######  Lemma 24.

Let $p_{1},\ldots,p_{m}$ be polynomials of degree at most $n$ and let
$N\in\mathbb{N}$ be such that $2N-1\geq mn$. Then

$\int_{0}^{\infty}z^{\alpha}e^{-\beta
z}\prod_{k=1}^{m}p_{k}(z)\mathrm{d}z=\beta^{-\alpha-1}\sum_{i=1}^{N}W_{i}\prod_{k=1}^{m}p_{k}\left(\frac{z_{i}}{\beta}\right),\qquad\beta>0,$
(32)

where $z_{i}$ denotes the $i$th root of $L^{(\alpha)}_{N}$ and

$W_{i}=\frac{\Gamma(N+\alpha+1)z_{i}}{N!(N+1)^{2}[L^{(\alpha)}_{N+1}(z_{i})]^{2}}\qquad
i=1,\ldots,N.$

###### Proof.

This follows from using the chain rule for the change of variables $z\mapsto
z/\beta$ and the $N$th Gauss–Laguerre quadrature (see [20, Formula 4.1]) which
is exact for polynomials of degree $2N-1$ or less. ∎

In our implementation, $p_{k}$ will usually be written as

$p_{k}=\sum_{j=0}^{n}a^{k}_{j}L^{(\alpha)}_{j}/\|\psi_{j}\|_{L^{2}(\mu)}.$

Thus, the above quadrature rule would be written as

$\int_{0}^{\infty}z^{\alpha}e^{-\beta
z}\prod_{k=1}^{m}p_{k}(z)\mathrm{d}z=\beta^{-\alpha-1}\sum_{i=1}^{N}W_{i}\prod_{k=1}^{m}\sum_{j=0}^{n}a^{k}_{j}L^{(\alpha)}_{j}\left(\frac{z_{i}}{\beta}\right)=\beta^{-\alpha-1}\sum_{i=1}^{N}W_{i}\prod_{k=1}^{m}(V\mathbf{a}^{k})_{i},$

where $\mathbf{a}^{k}=(a^{k}_{j})_{0\leq j\leq n}$ and
$V=(V_{ij})=(L^{(\alpha)}_{j}(z_{i}/\beta)/\|\psi_{j}\|_{L^{2}(\mu)})$ is a
(pseudo)-Vandermonde matrix. However, this matrix has very large entries so
that (32) as written needs to be computed in high precision arithmetic. When
(32) needs to be evaluated fast and accurately, it proves better to use a
regularised Vandermonde matrix
$\bar{V}=\mathrm{Diag}(\beta^{-\alpha-1}W)^{1/m}V$ so that

$\int_{0}^{\infty}z^{\alpha}e^{-\beta
z}\prod_{k=1}^{m}p_{k}(z)\mathrm{d}z=\sum_{i=1}^{N}\prod_{k=1}^{m}(\bar{V}\mathbf{a}^{k})_{i}.$

In practice, $\bar{V}$ can be precomputed once (using high precision if needed
but then possibly stored in lower precision) and after this is done both the
rigorous and non-rigorous computations can be carried out in a manageable
amount of time.

######  Remark 25.

In practice, $\bar{V}$ proves much better behaved than $V$ and easier to work
with. It could be interesting to further investigate properties of $\bar{V}$.

######  Example 26.

In both the rigorous and non-rigorous parts, one needs to compute the matrix

$G_{ij}=\langle\hat{\psi}_{i},\bar{u}^{p-1}\hat{\psi}_{j}\rangle$

which can be written as
$G=\bar{V}^{T}\mathrm{Diag}(\bar{V}\bar{\mathbf{u}})^{p-1}\bar{V}$ for $m=p+1$
and $\beta=p$ (and $\bar{\mathbf{u}}$ as in (29)).

Note that here the nonlinearity in the Galerkin problem is not approximated by
its so-called _pseudospectral_ projection found in the previous related
literature [5, 21].

Like the validation procedure, the quadrature is implemented rigorously in
Julia using the package IntervalArithmetic.jl [40]. Our implementation also
relies on approximating the roots of $L^{(\alpha)}_{N}$ with high precision,
using initial guesses given by FastGaussQuadrature.jl [46] and refining them
using Newton’s method in BigFloat arithmetic. The roots are then enclosed
rigorously by combining the Fundamental Theorem of Algebra and the
Intermediate Value Theorem (recall that the roots of orthogonal polynomials
are always simple and real). $V$ is then computed by recursion. While this
quadrature might appear expensive, in our experience for problems with $N\leq
2000$, the enclosure of the roots together with the computation of $\bar{V}$
only took a few minutes on a standard laptop.

######  Remark 27.

Hermite and Laguerre bases are far less popular than Fourier and Chebyshev for
computer-assisted proofs, or in numerical analysis for this matter. One reason
is the lack of simple and global error control on interpolation and
quadrature. However, polynomials are well-controlled in the bulk, and the
tails of $\hat{\psi}_{m}$ can be controlled via the following estimate of
Szegö [1, formula 22.14.13]

$|\hat{\psi}_{m}(r)|\leq\hat{\psi}_{m}(0)e^{-r^{2}/8}=\frac{1}{\Gamma(\alpha+1)}\sqrt{\frac{\Gamma(m+\alpha+1)}{m!}}e^{-r^{2}/8},\qquad\alpha\geq
0.$ (33)

In particular, this estimate in theory allows for integration with rigorous
error control, say by bluntly integrating in the bulk via interval arithmetic.
Since our method will only rely on integration, it might handle more general
nonlinearities. Note that this estimate immediately provides
$L^{\infty}$-bounds on finite linear combinations of the $\hat{\psi}_{m}$’s.

##  The semilinear heat equation

In this section, we are interested in the semilinear heat equation

$\partial_{t}v=\Delta v-\varepsilon|v|^{p-1}v\qquad\mbox{on $\R^{d}$},$ (34)

with $p>1$. It is well-known that if $u$ solves the nonlinear elliptic problem

$-\Delta u-\frac{x}{2}\cdot\nabla{u}-\frac{u}{p-1}+\varepsilon|u|^{p-1}p=0,$
(35)

then $v(t,x)=t^{-1/(p-1)}u(x/\sqrt{t})$ is a self-similar solution of (34).
These solutions are relevant when studying the asymptotic behaviour of
solutions to Eq. (34). Indeed, provided that $v(0,\cdot)\in H^{1}(\mu)$, they
appear as (possibly unstable) “omega-limit points” [26] of $v(t,\cdot)$ in the
sense that there exists a solution $u$ to Eq. (35) and
$t_{n}\rightarrow\infty$ such that

$\lim_{n\to\infty}\left\|t_{n}^{1/(p-1)}v(t_{n},\cdot)-u\left(\frac{\cdot}{\sqrt{t_{n}}}\right)\right\|_{\infty}=0.$

Let us first focus on positive radial solutions of Eq. (35) (we will later
show in Section 4.3 that the solutions we obtain are indeed positive) which
have been the subject of [11, 38], so that Eq. (35) becomes

$\mathcal{L}u=-\partial_{r}^{2}u-\left(\frac{d-1}{r}+\frac{r}{2}\right)\partial_{r}u=\frac{u}{p-1}-\varepsilon
u^{p}.$ (36)

Furthermore, let us assume first for simplicity that $p\in\mathbb{N}$. An
example with a non-integer power $p$ will be presented in Section 4.4, and
non-radial solutions will be studied in Section 4.5.

###  Suitable bounds for the validation

In this section, we apply the general strategy outlined in Section 2.4 to Eq.
(36), and derive computable estimates $Y$ and $Z_{k}$ satisfying the
assumptions of Corollary 21. Since we are first looking for radial solutions,
we restrict ourselves to the Hilbert space

$\mathcal{H}=\overline{\mathrm{Span}\\{\psi_{m}\\}_{m=0}^{\infty}}^{H^{2}}\subset
H^{2}(\mu).$

Throughout this section, we then assume $\bar{u}$ is an element of
$\mathcal{H}_{n}:=\mathrm{Span}\\{\psi_{m}\\}_{m=0}^{n}\subset\mathcal{C}_{0}^{\infty}(\mathbb{R}_{+})$.
In particular, $\|\bar{u}\|_{\infty}$ is finite and can be computed rigorously
(or at least upper-bounded rigorously, say via the estimate (33)). This
quantity will appear in several of the estimates to come.

In order to study Eq. (36), we consider

$\displaystyle F:\mathcal{H}$ $\displaystyle\longrightarrow\mathcal{H},$
$\displaystyle u$ $\displaystyle\longmapsto
u-\mathcal{L}^{-1}\left(\frac{u}{p-1}-\varepsilon u^{p}\right),$ (37)

and $T$ defined as in Section 2.4. Since we work in $\mathcal{H}$, the
approximate inverse $A$ involved in $T$ is naturally taken to go from
$\mathcal{H}$ to $\mathcal{H}$, so that $T:\mathcal{H}\to\mathcal{H}$. We
obtain the following formulae for $DT^{k}(\bar{u})$

$\begin{cases}DT(\bar{u})h&=h-ADF(\bar{u})h=h-A(h-\mathcal{L}^{-1}(h/(p-1)-\varepsilon
p\bar{u}^{p-1}h))\\\
D^{k}T(\bar{u})h^{k}&=\dfrac{p!}{(p-k)!}A\mathcal{L}^{-1}(\bar{u}^{p-k}h^{k})\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad
1<k\leq p.\end{cases}$

######  Notation 28.

Recall Notation 10 and Remark 11. Our analysis will mostly be done in
$\mu$-weighted spaces, and we may shorten the notation by omitting the $\mu$.
Henceforth, $\|\cdot\|_{H^{2}}$ (resp. $\|\cdot\|_{H^{1}}$,
$\|\cdot\|_{L^{2}}$) will always denote $\|\cdot\|_{H^{2}(\mu)}$ (resp.
$\|\cdot\|_{H^{1}(\mu)}$, $\|\cdot\|_{L^{2}(\mu)}$). Furthermore, we denote by
$\\{\hat{\psi}_{m}\\}_{m=0}^{\infty}$ the normalisation of the radial
orthogonal basis $\\{\psi_{m}\\}_{m=0}^{\infty}$ (introduced in Corollary 9)
with respect the $L^{2}(\mu)$-inner product $\langle\cdot,\cdot\rangle$, and
let $\lambda_{m}=d/2+m$ so that
$\mathcal{L}\hat{\psi}_{m}=\lambda_{m}\hat{\psi}_{m}$. We denote
$P_{\scriptscriptstyle\infty}=I-P_{n}$ the complementary projection to $P_{n}$
(in $L^{2}(\mu)$ or $H^{2}(\mu)$). Note that $P_{n}$, and therefore also
$P_{\scriptscriptstyle\infty}$, leave $\mathcal{H}$ invariant, and we abuse
notation and still denote by $P_{n}$ and $P_{\scriptscriptstyle\infty}$ their
restriction to $\mathcal{H}$ so that $P_{n}:\mathcal{H}\to\mathcal{H}$ is the
orthogonal projection onto $\mathcal{H}_{n}$. We also denote by
$\mathcal{H}_{\scriptscriptstyle\infty}=P_{\scriptscriptstyle\infty}(\mathcal{H})$
the orthogonal complement of $\mathcal{H}_{n}$ in $\mathcal{H}$. Finally,
since our implementation will be with respect to the $L^{2}(\mu)$-basis, in
what follows, vectors in $\mathbb{R}^{n+1}$ are identified with the
corresponding decomposition with respect to $\\{\hat{\psi}_{m}\\}_{m=0}^{n}$
and we aim at stating final estimates with respect to the $L^{2}(\mu)$-norm.

####  4.1.1 The bound $Y$

It suffices to estimate

$\displaystyle\|T(\bar{u})-\bar{u}\|^{2}_{H^{2}}$
$\displaystyle=\|AF(\bar{u})\|_{H^{2}}^{2}$
$\displaystyle=\|P_{n}(AF(\bar{u}))\|^{2}_{H^{2}}+\|P_{\scriptscriptstyle\infty}(AF(\bar{u}))\|_{H^{2}}^{2}$
$\displaystyle=\|A_{n}(P_{n}F(\bar{u}))\|^{2}_{H^{2}}+\|P_{\scriptscriptstyle\infty}F(\bar{u})\|^{2}_{H^{2}}$
$\displaystyle=\|A_{n}(P_{n}F(\bar{u}))\|^{2}_{H^{2}}+\|P_{\scriptscriptstyle\infty}\mathcal{L}^{-1}(\bar{u}^{p})\|^{2}_{H^{2}}$
$\displaystyle=\|A_{n}(P_{n}F(\bar{u}))\|^{2}_{H^{2}}+\left(\|\mathcal{L}^{-1}\bar{u}^{p}\|^{2}_{H^{2}}-\|P_{n}\mathcal{L}^{-1}(\bar{u}^{p})\|^{2}_{H^{2}}\right).$

Thus, we may choose

$Y^{2}:=\|\mathcal{L}A_{n}(P_{n}F(\bar{u}))\|^{2}_{L^{2}}+\left(\|\bar{u}^{p}\|^{2}_{L^{2}}-\|P_{n}(\bar{u}^{p})\|^{2}_{L^{2}}\right).$

Note that all the terms in $Y$ can be evaluated rigorously by means of a
quadrature rule (see Lemma 24) and the output is expected to be small provided
the truncation level $n$ is large enough. In particular, the first term can
practically be made arbitrarily small in the non-rigorous computation of
$\bar{u}$.

######  Remark 29.

Here and in many subsequent places, we slightly abuse notation, in the sense
that the number $Y$ that we will actually compute is only an upper bound of
the quantity $Y$ defined above. Indeed, we use interval arithmetic to control
rounding errors, and in later estimates we also only upper-bounds quantities
like $\left\|\bar{u}\right\|_{\infty}$ rather than trying to compute them
_exactly_.

####  4.1.2 The bound $Z_{1}$

We choose to compute $Z_{1}$ by treating finite and infinite dimensional parts
separately, and therefore take

$Z_{1}:=\left\|\left(\begin{matrix}Z^{11}&Z^{12}\\\
Z^{21}&Z^{22}\end{matrix}\right)\right\|_{2}$

where the $Z^{ij}$ are real numbers satisfying

$\displaystyle\begin{cases}\|P_{n}DT(\bar{u})P_{n}\|_{H^{2},H^{2}}&\leq
Z^{11},\\\ \|P_{\scriptscriptstyle\infty}DT(\bar{u})P_{n}\|_{H^{2},H^{2}}&\leq
Z^{21},\\\ \|P_{n}DT(\bar{u})P_{\scriptscriptstyle\infty}\|_{H^{2},H^{2}}&\leq
Z^{12},\\\
\|P_{\scriptscriptstyle\infty}DT(\bar{u})P_{\scriptscriptstyle\infty}\|_{H^{2},H^{2}}&\leq
Z^{22}.\\\ \end{cases}$ (38)

Indeed, for $h\in\mathcal{H}$ with $\|h\|_{H^{2}}=1$, write
$h=h_{n}+h_{\scriptscriptstyle\infty}$ such that $h_{n}\in\mathcal{H}_{n}$,
$h_{\scriptscriptstyle\infty}\in\mathcal{H}_{\infty}$. Then,

$\displaystyle\|DT(\bar{u})h\|_{H^{2}}$
$\displaystyle=\sqrt{\|P_{n}DT(\bar{u})h\|_{H^{2}}^{2}+\|P_{n}DT(\bar{u})h\|_{H^{2}}^{2}}$
$\displaystyle\leq\sqrt{\left(\|P_{n}DT(\bar{u})h_{n}\|_{H^{2}}+\|P_{n}DT(\bar{u})h_{\scriptscriptstyle\infty}\|_{H^{2}}\right)^{2}+\left(\|P_{\scriptscriptstyle\infty}DT(\bar{u})h_{n}\|_{H^{2}}+\|P_{\scriptscriptstyle\infty}DT(\bar{u})h_{\scriptscriptstyle\infty}\|_{H^{2}}\right)^{2}}$
$\displaystyle\leq\sqrt{\left(Z^{11}\|h_{n}\|_{H^{2}}+Z^{12}\|h_{\scriptscriptstyle\infty}\|_{H^{2}}\right)^{2}+\left(Z^{21}\|h_{n}\|_{H^{2}}+Z^{22}\|h_{\scriptscriptstyle\infty}\|_{H^{2}}\right)^{2}}$
$\displaystyle=\left\|\begin{pmatrix}Z^{11}&Z^{12}\\\
Z^{21}&Z^{22}\end{pmatrix}\begin{pmatrix}\|h_{n}\|_{H^{2}}\\\
\|h_{\scriptscriptstyle\infty}\|_{H^{2}}\end{pmatrix}\right\|_{2}$
$\displaystyle\leq Z_{1}\qquad\mbox{since
$\|h\|_{H^{2}}^{2}+\|h_{\scriptscriptstyle\infty}\|_{H^{2}}^{2}=1$.}$

In the following subsections, we derive computable bounds $Z^{ij}$ satisfying
(38).

####  4.1.3 The bound $Z^{11}$

For $h_{n}\in\mathcal{H}_{n}$,

$P_{n}DT(\bar{u})h_{n}=P_{n}(I-ADF(\bar{u}))h_{n}=(I_{n}-A_{n}P_{n}DF(\bar{u})P_{n})h_{n}$

by construction of $A$. Then, we get the estimate

$\|P_{n}DT(\bar{u})P_{n}\|_{H^{2},H^{2}}=\|I_{n}-A_{n}P_{n}DF(\bar{u})P_{n}\|_{H^{2},H^{2}}=\|\mathcal{L}(I_{n}-A_{n}P_{n}DF(\bar{u})P_{n})\mathcal{L}^{-1}\|_{L^{2},L^{2}}=:Z^{11},$

which is expected to be small because
$A_{n}\approx(P_{n}DF(\bar{u})P_{n})^{-1}$. Using the isometric isomorphism
between $L^{2}(\mu)$ and $\ell^{2}$, this quantity $Z^{11}$ is in fact nothing
but the $2$-norm of an $(n+1)\times(n+1)$ matrix, whose entries can be
evaluated rigorously. For $M\in\R^{(n+1)\times(n+1)}$, we may bound above its
2-norm with the estimate $\|M\|_{2}\leq\sqrt{\|M\|_{1}\|M\|_{\infty}}$.

####  4.1.4 The bound $Z^{21}$

For $h_{n}\in\mathcal{H}_{n}$, we get that

$\displaystyle P_{\scriptscriptstyle\infty}DT(\bar{u})h_{n}$
$\displaystyle=-P_{\scriptscriptstyle\infty}DF(\bar{u})h_{n}$
$\displaystyle=P_{\scriptscriptstyle\infty}\mathcal{L}^{-1}(h_{n}/(p-1)-\varepsilon
p\bar{u}^{p-1}h_{n})$ $\displaystyle=-\varepsilon
pP_{\scriptscriptstyle\infty}\mathcal{L}^{-1}(\bar{u}^{p-1}h_{n}),$

and thus,

$\|P_{\scriptscriptstyle\infty}DT(\bar{u})h_{n}\|_{H^{2}}=p\|P_{\scriptscriptstyle\infty}(\bar{u}^{p-1}h_{n})\|_{L^{2}}.$

Then, writing $h_{n}\in\mathcal{H}_{n}$ as

$h_{n}=\sum_{m=0}^{n}a_{m}\hat{\psi}_{m},$

we get

$P_{\scriptscriptstyle\infty}(\bar{u}^{p-1}h_{n})=\sum_{m=0}^{n}a_{m}P_{\scriptscriptstyle\infty}(\bar{u}^{p-1}\hat{\psi}_{m}).$

Thus,

$\displaystyle\|P_{\scriptscriptstyle\infty}(\bar{u}^{p-1}h_{n})\|_{L^{2}}$
$\displaystyle\leq\sum_{m=0}^{n}\|a_{m}P_{\scriptscriptstyle\infty}(\bar{u}^{p-1}\hat{\psi}_{m})\|_{L^{2}}\qquad\mbox{by
Minkowski inequality}$
$\displaystyle\leq\sum_{m=0}^{n}|a_{m}\lambda_{m}|\|\lambda^{-1}_{m}P_{\scriptscriptstyle\infty}(\bar{u}^{p-1}\hat{\psi}_{m})\|_{L^{2}}$
$\displaystyle\leq\left(\sum_{m=0}^{n}\lambda^{-2}_{m}\|P_{\scriptscriptstyle\infty}(\bar{u}^{p-1}\hat{\psi}_{m})\|^{2}_{L^{2}}\right)^{1/2}\left(\sum_{m=0}^{n}\lambda_{m}^{2}a_{m}^{2}\right)^{1/2}\qquad\mbox{by
Cauchy--Schwarz inequality}$
$\displaystyle=\left(\sum_{m=0}^{n}\lambda^{-2}_{m}\|P_{\scriptscriptstyle\infty}(\bar{u}^{p-1}\hat{\psi}_{m})\|^{2}_{L^{2}}\right)^{1/2}\|h_{n}\|_{H^{2}}.$

Introducing the vector $w\in\mathbb{R}^{n+1}$ with

$w_{m}:=\|P_{\scriptscriptstyle\infty}(\bar{u}^{p-1}\hat{\psi}_{m})\|_{L^{2}}=\left(\|\bar{u}^{p-1}\hat{\psi}_{m}\|^{2}_{L^{2}}-\|P_{n}(\bar{u}^{p-1}\hat{\psi}_{m})\|^{2}_{L^{2}}\right)^{1/2},$

and identifying $w$ with an element of $P_{n}(L^{2}(\mu))\subset L^{2}(\mu)$,
we get

$\|P_{\scriptscriptstyle\infty}DT(\bar{u})P_{n}\|_{H^{2},H^{2}}\leq
p\|\mathcal{L}^{-1}w\|_{L^{2}}=:Z^{21}.$

Note that $w_{m}$ is expected to be small, especially for small $m$, whereas
for larger $m$ the $\lambda_{m}^{-2}$ factors also contribute to making
$\|\mathcal{L}^{-1}w\|^{2}_{L^{2}}=\sum_{m=0}^{n}\lambda^{-2}_{m}w_{m}^{2}$
small.

####  4.1.5 The bound $Z^{12}$

Now, for
$h_{\scriptscriptstyle\infty}\in\mathcal{H}_{\scriptscriptstyle\infty}$,

$\displaystyle P_{n}DT(\bar{u})h_{\scriptscriptstyle\infty}$
$\displaystyle=-P_{n}ADF(\bar{u})h_{\scriptscriptstyle\infty}$
$\displaystyle=P_{n}A\mathcal{L}^{-1}h_{\scriptscriptstyle\infty}/(p-1)-\varepsilon
pP_{n}A\mathcal{L}^{-1}(\bar{u}^{p-1}h_{\scriptscriptstyle\infty})$
$\displaystyle=\mathcal{L}^{-1}P_{n}h_{\scriptscriptstyle\infty}/(p-1)-\varepsilon
pA_{n}P_{n}\mathcal{L}^{-1}(\bar{u}^{p-1}h_{\scriptscriptstyle\infty})$
$\displaystyle=-\varepsilon
pA_{n}\mathcal{L}^{-1}P_{n}(\bar{u}^{p-1}h_{\scriptscriptstyle\infty})\qquad\mbox{since
$h_{\scriptscriptstyle\infty}\in\mathcal{H}_{\scriptscriptstyle\infty}$.}$

Then, for $m\in\\{0,\ldots n\\}$, we estimate

$\displaystyle|\langle\bar{u}^{p-1}h_{\scriptscriptstyle\infty},\hat{\psi}_{m}\rangle|$
$\displaystyle=|\langle\bar{u}^{p-1}\hat{\psi}_{m},h_{\scriptscriptstyle\infty}\rangle|$
$\displaystyle=|\langle
P_{\scriptscriptstyle\infty}(\bar{u}^{p-1}\hat{\psi}_{m}),h_{\scriptscriptstyle\infty}\rangle|\qquad\mbox{since
$h_{\scriptscriptstyle\infty}\in\mathcal{H}_{\scriptscriptstyle\infty}$}$
$\displaystyle\leq\|P_{\scriptscriptstyle\infty}(\bar{u}^{p-1}\hat{\psi}_{m})\|_{L^{2}}\|h_{\scriptscriptstyle\infty}\|_{L^{2}}\qquad\mbox{by
Cauchy--Schwarz inequality}$
$\displaystyle\leq\frac{\|h_{\scriptscriptstyle\infty}\|_{H^{2}}}{\lambda_{n+1}}\|P_{\scriptscriptstyle\infty}(\bar{u}^{p-1}\hat{\psi}_{m})\|_{L^{2}}\qquad\mbox{by\leavevmode\nobreak\
\eqref{eq:projection_error_L2H2}}.$

Thus, reusing the vector $w$ introduced in the previous subsection,

$\|P_{n}DT(\bar{u})P_{\scriptscriptstyle\infty}\|_{H^{2},H^{2}}\leq{\frac{p}{\lambda_{n+1}}\left\|\,\left|A_{n}\mathcal{L}^{-1}\right|w\right\|_{H^{2}}=\frac{p}{\lambda_{n+1}}\left\|\,\left|\mathcal{L}A_{n}\mathcal{L}^{-1}\right|w\right\|_{L^{2}}}=:Z^{12},$

where for $M\in\R^{(n+1)\times(n+1)}$, $|M|=(|M|_{ij})=(|M_{ij}|)$.

####  4.1.6 The bound $Z^{22}$

For $h_{\scriptscriptstyle\infty}\in\mathcal{H}_{\scriptscriptstyle\infty}$,

$P_{\scriptscriptstyle\infty}DT(\bar{u})h_{\scriptscriptstyle\infty}=-P_{\scriptscriptstyle\infty}DF(\bar{u})h_{\scriptscriptstyle\infty}=P_{\scriptscriptstyle\infty}\mathcal{L}^{-1}((1/(p-1)-\varepsilon
p\bar{u}^{p-1})h_{\scriptscriptstyle\infty}.$

Thus,

$\displaystyle\|P_{\scriptscriptstyle\infty}DT(\bar{u})h_{\scriptscriptstyle\infty}\|_{H^{2}}$
$\displaystyle=\|(1/(p-1)-\varepsilon
p\bar{u}^{p-1})h_{\scriptscriptstyle\infty}\|_{L^{2}}$
$\displaystyle\leq\|1/(p-1)-\varepsilon
p\bar{u}^{p-1}\|_{\infty}\|h_{\scriptscriptstyle\infty}\|_{L^{2}}\qquad\mbox{by
H\"{o}lder's inequality }$
$\displaystyle\leq\frac{(1/(p-1)+p\|\bar{u}\|^{p-1}_{\infty})}{\lambda_{n+1}}\|h_{\scriptscriptstyle\infty}\|_{H^{2}}\qquad\mbox{by\leavevmode\nobreak\
\eqref{eq:projection_error_L2H2}}.$

Thus, we may choose

$Z^{22}:=\frac{1/(p-1)+p\|\bar{u}\|^{p-1}_{\infty}}{\lambda_{n+1}},$

where $\|\bar{u}\|_{\infty}$ can be computed (or at least upper-bounded) using
(33). Note that while $T:\mathcal{H}\rightarrow\mathcal{H}$ may not be
compact, $DT(\bar{u}):\mathcal{H}\rightarrow\mathcal{H}$ is, provided
$\bar{u}\in L^{\infty}(\R^{d})$.

####  4.1.7 The bounds $Z_{k}$ for $1<k\leq p$

If $d\leq 2$, or $p\leq d/(d-2)$ for $d\geq 3$, we can use the $H^{1}$-Sobolev
embeddings given by Lemma 14, which yield, for $h\in\mathcal{H}$,

$\displaystyle\left\|D^{k}T(\bar{u})h^{k}\right\|_{H^{2}}$
$\displaystyle=\frac{p!}{(p-k)!}\left\|A\mathcal{L}^{-1}(\bar{u}^{p-k}h^{k})\right\|_{H^{2}}$
$\displaystyle\leq\frac{p!}{(p-k)!}\|A\mathcal{\|}_{H^{2},H^{2}}\|\mathcal{L}^{-1}(\bar{u}^{p-k}h^{k})\|_{H^{2}}$
$\displaystyle\leq\frac{p!}{(p-k)!}\|\mathcal{L}A\mathcal{L}^{-1}\|_{L^{2},L^{2}}\|\bar{u}^{p-k}h^{k}\|_{L^{2}}$
$\displaystyle\leq\frac{p!\|\bar{u}\|_{\infty}^{p-k}}{(p-k)!}\|\mathcal{L}A\mathcal{L}^{-1}\|_{L^{2},L^{2}}\|h^{k}\|_{L^{2}}\qquad\mbox{by
H\"{o}lder's inequality}$
$\displaystyle\leq\frac{p!\|\bar{u}\|_{\infty}^{p-k}}{(p-k)!}\|\mathcal{L}A\mathcal{L}^{-1}\|_{L^{2},L^{2}}c(d,k)\|h\|_{H^{1}}^{k}$
$\displaystyle\leq\frac{p!\|\bar{u}\|_{\infty}^{p-k}}{(p-k)!}\|\mathcal{L}A\mathcal{L}^{-1}\|_{L^{2},L^{2}}c(d,k)\lambda_{0}^{-k/2}\|h\|_{H^{2}}^{k},$

where $c(d,k)$ is given by Corollary 15 and we have used (13). Thus, we may
choose

$Z_{k}:=\frac{p!\|\bar{u}\|_{\infty}^{p-k}}{(p-k)!}\lambda_{0}^{-k/2}c(d,k)\|\mathcal{L}A\mathcal{L}^{-1}\|_{L^{2},L^{2}},$
(39)

where

$\|\mathcal{L}A\mathcal{L}^{-1}\|_{L^{2},L^{2}}\leq\max\left\\{\|\mathcal{L}A_{n}\mathcal{L}^{-1}\|_{L^{2},L^{2}},1\right\\}.$

Note that instead of using $\|\bar{u}\|_{\infty}$, for $k<p$, we may derive a
bound only using the embedding $H^{1}(\mu)\hookrightarrow L^{2p}(\mu)$, but
this is not expected to be as sharp.

For $d\in\\{2,3\\}$, we may instead use Lemma 19 relying on $H^{2}$-Sobolev
embeddings

$\displaystyle\left\|D^{k}T(\bar{u})h^{k}\right\|_{H^{2}}$
$\displaystyle\leq\frac{p!\|\bar{u}\|_{\infty}^{p-k}}{(p-k)!}\|\mathcal{L}A\mathcal{L}^{-1}\|_{L^{2},L^{2}}\|h^{k}\|_{L^{2}}$
$\displaystyle\leq\frac{p!\|\bar{u}\|_{\infty}^{p-k}}{(p-k)!}\|\mathcal{L}A\mathcal{L}^{-1}\|_{L^{2},L^{2}}\|h\|_{\infty}^{k-1}\|h\|_{L^{2}}$
$\displaystyle\leq\frac{p!\|\bar{u}\|_{\infty}^{p-k}}{\lambda_{0}(p-k)!}\|\mathcal{L}A\mathcal{L}^{-1}\|_{L^{2},L^{2}}C(d)^{k-1}\|h\|_{H^{2}}^{k},$

where $C(d)$ is given by Lemma 19. Thus, we may also choose

$Z_{k}:=\frac{p!\|\bar{u}\|_{\infty}^{p-k}}{\lambda_{0}(p-k)!}C(d)^{k-1}\|\mathcal{L}A\mathcal{L}^{-1}\|_{L^{2},L^{2}}.$
(40)

In practice, when $d\in\\{2,3\\}$, if both estimates are valid, we take for
$Z_{k}$ the minimum between (39) and (40).

###  Results

We first focus on the case $\varepsilon=-1$ in (36), such that the equation
reads

$\mathcal{L}u=-\partial_{r}^{2}u-\left(\frac{d-1}{r}+\frac{r}{2}\right)\partial_{r}u=\frac{u}{p-1}+u^{p}.$
(41)

In the cases $(d,p)=(2,3)$ and $(d,p)=(3,2)$, we find the numerical solutions
plotted below with $n=500$.

(a) $d=2$, $p=3$

(b) $d=3$, $p=2$

Figure 3: Representation of the numerical approximate solutions $\bar{u}$ to Eq. (41) which are validated in Theorem 30. Note that the spatial domain of course has to be truncated for the plot, but these approximate solutions are actually defined and computed on $[0,\infty)$ thanks to our choice of basis (see Corollary 9). |  $\varepsilon=-1$, $d=2$, $p=3$ |  $\varepsilon=-1$, $d=3$, $p=2$
---|---|---
$Y$ | $6.830397272197097\times 10^{-13}$ | $1.4114627225937143\times 10^{-21}$
$Z_{1}$ | $0.009916643027658332$ | $0.008621565176356588$
$Z_{2}$ | $170.320283777825$ | $17.04079355792355$
$Z_{3}$ | $719.2983917297151$ | 0
$\underline{\delta}$ | $6.898810312173614\times 10^{-13}$ | $1.4237375688374744\times 10^{-21}$
$\bar{\delta}$ | $0.0057434127215418$ | $0.05817677630175133$
Table 1: Results for the a posteriori validation with the bounds given in
Section 4.1
and $\bar{u}$ represented in Figure 3.

######  Theorem 30.

For $(d,p)=(2,3)$ (resp. $(d,p)=(3,2)$), let $\bar{u}:\R_{+}\to\R$ be the
radial function whose restriction to the interval is represented in Figure
3(a) (resp. Figure 3(b)) and whose precise description in terms of “Fourier”
coefficients in the basis $\\{\hat{\psi}_{m}\\}_{m=0}^{n}$ is available at
[17], then there exists a solution $u^{\star}\in\mathcal{H}$ to Eq. (41) such
that

$\|u^{\star}-\bar{u}\|_{H^{2}(\mu)}\leq\underline{\delta},$

where $\underline{\delta}$ is as in Table 1. It is the only solution to Eq.
(41) in $\mathcal{H}$ such that
$\|u^{\star}-\bar{u}\|_{H^{2}(\mu)}\leq\bar{\delta}$.

###### Proof.

We consider $F:\mathcal{H}\to\mathcal{H}$ defined in (37). In the file
Heat/d2_p3/proof.ipynb (resp.
Heat/d3_p2/proof.ipynb) at [17], we construct $A_{n}$ (and thus indirectly $A$
and $T$) as outlined in Section 2.4. Following the estimates of Section 4.1,
we obtain the upper bounds in Table 1 for $Y$ and the $Z_{k}$’s to conclude
with Corollary 21. ∎

###  Proving positivity a posteriori

Similarly to the previous literature on the semilinear heat equation, we also
wish to verify the (strict) positivity of our radial solutions in this work.
Note that even if $\bar{u}$ could be proved to be strictly positive, for any
$\delta>0$ $B(\bar{u},\delta)$ necessarily contains functions that are
eventually negative at $+\infty$, e.g. $\bar{u}-(\delta/2)\hat{\psi}_{m}$ for
large enough $m$. Thus, the positivity of $u^{\star}$ cannot be deduced
directly. Since $\bar{u}$ decays faster than our $L^{\infty}$-control, a
direct comparison of $\bar{u}$ and $u^{\star}$ is only sufficient in a
neighbourhood of zero. We thus resort to a maximum principle argument to show
the (strict) positivity of $u^{\star}$ away from zero with the lemma below.

######  Lemma 31.

Let $b:\R^{d}\to\R^{d}$, $c:\R^{d}\to\R$ two $\mathcal{C}^{0}$ functions and
$\varphi:\R^{d}\to\R$ a $\mathcal{C}^{2}$ function such that

$\displaystyle-\Delta\varphi+b\cdot\nabla\varphi+c\varphi=0.$ (42)

Assume there exists $r_{0}>0$ such that $c(x)>0$ for all $|x|\geq r_{0}$.
Assume further that $\varphi(x)>0$ for all $|x|\leq r_{0}$, and that
$\varphi(x)\underset{|x|\to\infty}{\longrightarrow}0$. Then, $\varphi(x)>0$
for all $x\in\R^{d}$.

This statement is essentially due to the maximum principle, which can be used
even though we are on an unbounded domain because we already control the
behaviour of the solution at infinity.

###### Proof.

We argue by contradiction and assume that $\varphi$ is nonpositive somewhere.
Then, since $\varphi(x)\underset{|x|\to\infty}{\longrightarrow}0$, there
exists $x_{0}\in\R^{d}$ such that $\varphi(x_{0})=\min_{\R^{d}}\varphi\leq 0$,
and we must have $|x_{0}|>r_{0}$. Let $U$ be a connected open subset of
$\R^{d}$ containing $x_{0}$, such that $|x|\geq r_{0}$ for all $x\in U$, and
such that there exists $x$ in $\partial U$ with $|x|=r_{0}$. Since $\varphi$
reaches its (nonpositive) minimum over $\bar{U}$ at the interior point
$x_{0}$, while $c$ is positive in $U$, the strong maximum principle yields
that $\varphi$ must be constant on $\bar{U}$. In particular, $\varphi$ is
nonpositive on $\partial U$, which contradicts the fact the $\varphi(x)$ is
positive for all $|x|\leq r_{0}$. ∎

Note that with $\varepsilon=-1$ in (36), we are in the “bad-sign” case of the
semilinear heat equation, so that the strong maximum principle and the above
lemma cannot be applied directly to equation (35). However,
$\varphi=e^{|x|^{2}/8}u^{\star}$ satisfies an equation of the form (42),
therefore its (strict) positivity can be verified via Lemma 31, which in turn
is equivalent to the (strict) positivity of $u^{\star}$. For clarity, let us
denote below a radial solution $u^{\star}$ as a function of $r=|x|$.

######  Corollary 32.

Let $u^{\star}\in H^{1}(\mu)$ be a solution of Eq. (36) and $r_{0}>0$ be such
that

$\frac{d}{4}+\frac{r^{2}}{16}-\frac{1}{p-1}+\varepsilon
u^{\star}(r)^{p-1}>0\qquad\mbox{for all $r\geq r_{0}$}.$ (43)

If in addition, $u^{\star}(r)>0$ for all $r\leq r_{0}$, then $u^{\star}(r)>0$
for all $r\in\mathbb{R}_{+}$.

###### Proof.

Define

$\displaystyle\varphi:\R_{+}$ $\displaystyle\longrightarrow\R$ $\displaystyle
r$ $\displaystyle\longmapsto e^{r^{2}/8}u^{\star}(r).$

If $u^{\star}\in H^{1}(\mu)$ is a solution to Eq. (36), then by [19, Theorem
3.12, Proposition 3.15] it is also a classical solution in the sense that
$u^{\star}\in\mathcal{C}^{2}(\R^{d})$. Since
$u^{\star}(|\cdot|)\in\mathcal{C}^{2}(\R^{d})$, we observe that $u^{\star}\in
H^{1}(\R_{+},e^{r^{2}/4})$ and thus from the proof of Lemma 14 in the case
$d=1$, we have that $\varphi\in H^{1}(\R_{+})$. Therefore from [10, Corollary
VIII.8], $\lim_{r\to+\infty}\varphi(r)=0$. Now, letting $b=0$ and

$c(x)=\frac{d}{4}+\frac{|x|^{2}}{16}-\frac{1}{p-1}+\varepsilon
u^{\star}(|x|)^{p-1},$

we have that $\varphi(|\cdot|)$ solves the equation

$-\Delta\varphi+c\varphi=0.$

Since $c,\varphi(|\cdot|)\in\mathcal{C}^{2}(\R^{d})$ fulfill the assumptions
of Lemma 31, under condition (43), we obtain that $\varphi(|x|)>0$ for all
$x\in\R^{d}$ and thus $u^{\star}(r)>0$ for all $r\in\mathbb{R}_{+}$. ∎

In the case $d\in\\{1,2,3\\}$, it is straightforward to find some $r_{0}$ such
that condition (43) holds by using the $L^{\infty}$-bound obtained from Lemma
19 and using the estimate (33). Similarly, checking $u^{\star}$ is positive on
$[0,r_{0}]$ can be achieved efficiently by using Lemma 19 and by using
interval arithmetic to evaluate $P_{m}\bar{u}$ on $[0,r_{0}]$ for some small
$m$ and using the estimate (33) to bound $|(I-P_{m})\bar{u}(r)|$. Let us
remark that alternatively (in particular for $d\geq 4$),
$L^{\infty}$-estimates on $u^{\star}$ could be obtained by bootstrap arguments
along the line of the proof of [19, Theorem 3.12]. By iteratively applying
such arguments, one may also be able to obtain $\mathcal{C}^{k}$-bounds on
$u^{\star}-\bar{u}$.

######  Corollary 33.

Let $u^{\star}\in H^{2}(\mu)$ be as given by Theorem 30, $u^{\star}$ is
strictly positive on $\R^{d}$.

######  Remark 34.

In the event that the positivity of $u^{\star}$ is verified, in the case of
$\varepsilon=-1$, since $u^{\star}\in C^{2}(\R_{+})$, we immediately obtain
that it is a classical solution (in the sense of [38]) of

$\partial_{r}^{2}u+\left(\frac{d-1}{r}+\frac{r}{2}\right)\partial_{r}u+\frac{u}{p-1}+|u|^{p-1}u=0,$
(44)

and since $u^{\star}$ is rapidly decaying, the results of [38] can be applied
and yield $u^{\star}(r)=\mathcal{O}(r^{2/(p-1)-d}e^{-r^{2}/4})$. Furthermore,
[56, Theorem 1] and [33] ensure that $u^{\star}$ is the unique rapidly
decaying and strictly positive solution of Eq. (36) and of Eq. (35)
respectively. The analogous result for the case $\varepsilon=+1$ is given by
[11].

###  An example with a fractional exponent

In the case $\varepsilon=+1$, rapidly decaying self-similar solutions only
exist for $p<1+2/d$ [11, 19], and therefore integer exponents are excluded as
soon as $d>1$. We now explain how our method can be generalised to handle
cases when the exponent $p$ is a rational number. We focus on enclosing the
positive self-similar profile of the nonlinear heat equation with
$\varepsilon=+1$, $d=2$ and $p=5/3$. We thus aim to solve

$\mathcal{L}u-\frac{3}{2}u+u^{5/3}=0,$ (45)

where we choose the real branch $x\in\R\mapsto x^{5/3}:=(\sqrt[3]{x})^{2}x$.
We first look for an approximate solution $\bar{u}=\bar{U}^{3}$ such that

$\bar{U}(r)=q(r)e^{-r^{2}/12},$

where $q$ is a polynomial. The fractional exponent can then be _rigorously_
handled via a quadrature as in Lemma 24. For instance, here, we first start by
finding an approximate solution $\tilde{u}$ via a pseudospectral method [21],
and compute for $n\in 3\mathbb{N}$

$q(r)=\sum_{m=0}^{n/3}a_{m}L_{m}(r^{2}/12)\qquad\text{where
}a_{m}\approx\frac{1}{6}\int_{0}^{\infty}L_{m}(r^{2}/12)\sqrt[3]{\tilde{u}(r)}r\mathrm{d}r,$

so that $\bar{U}\approx\sqrt[3]{\tilde{u}}$. Finally, we define
$\bar{u}:=\bar{U}^{3}$ to apply Theorem 20. From here, $\bar{u}$ can be
identified as a vector of $n+1$ intervals as in (29) and computations can be
carried out in the same manner as in the other examples.

######  Remark 35.

It is clear that this approach can work with any rational exponent with an odd
denominator. In the case of an exponent with an even denominator, this
approach can only work if $\bar{U}$ can be verified to be nonnegative. Beyond
the fact that this verification most likely requires a lengthy computation,
finding such nonnegative $\bar{U}$ in the first place probably needs some
fine-tuning. A more sophisticated approach would be to bypass this issue by
reformulating the problem as a differential-algebraic equation as in [7], but
this would probably require additional bounds with an adapted functional
analysis.

Now, analysing the function $x\mapsto x^{2/3}:=(\sqrt[3]{x})^{2}$, one finds
that for all $x,y\in\R$,

$\left|\frac{\mathrm{d}y^{5/3}}{\mathrm{d}y}-\frac{\mathrm{d}x^{5/3}}{\mathrm{d}x}\right|=\frac{5}{3}\left|y^{2/3}-x^{2/3}\right|\leq\frac{5}{3}|x-y|^{2/3}.$

Thus, writing
$\|DT(u)\|_{H^{2},H^{2}}\leq\|DT(\bar{u})\|_{H^{2},H^{2}}+\|DT(u)-DT(\bar{u})\|_{H^{2},H^{2}}$
and estimating that for $h\in H^{2}(\mu)$

$\displaystyle\|(DT(u)-DT(\bar{u}))h\|_{H^{2}}$
$\displaystyle\leq\frac{5}{3}\|A\mathcal{L}^{-1}(u^{2/3}-\bar{u}^{2/3})h\|_{H^{2}}$
$\displaystyle\leq\frac{5}{3}\|\mathcal{L}A\mathcal{L}^{-1}\|_{L^{2},L^{2}}\|(u^{2/3}-\bar{u}^{2/3})h\|_{L^{2}}$
$\displaystyle\leq\frac{5}{3}\|\mathcal{L}A\mathcal{L}^{-1}\|_{L^{2},L^{2}}\|(u^{2/3}-\bar{u}^{2/3})\|_{\infty}\|h\|_{L^{2}}$
$\displaystyle\leq\frac{5}{3}\|\mathcal{L}A\mathcal{L}^{-1}\|_{L^{2},L^{2}}\|(u-\bar{u})\|^{2/3}_{\infty}\|h\|_{L^{2}}\qquad\mbox{by\leavevmode\nobreak\
\eqref{eqn:5-3-bound}}$
$\displaystyle\leq\frac{5}{3}C(2)^{2/3}\|\mathcal{L}A\mathcal{L}^{-1}\|_{L^{2},L^{2}}\|(u-\bar{u})\|^{2/3}_{H^{2}}\|h\|_{H^{2}}\qquad\mbox{by
Lemma\leavevmode\nobreak\ \ref{lem:infty-bound},}$

we may choose

$\mathcal{Z}(\delta):=Z_{1}+Z_{5/3}\delta^{2/3}$

satisfying (24), where $Z_{1}$ is as in Section 4.1 and

$Z_{5/3}=\frac{5}{3}C(2)^{2/3}\|\mathcal{L}A\mathcal{L}^{-1}\|_{L^{2},L^{2}}.$

For this example, we choose $n=1200$ and find $\bar{u}$ depicted in Figure 1.

###### Proof of Theorem 1.

We consider $F:\mathcal{H}\to\mathcal{H}$ defined in (37). In the file
Heat/fractional/proof.ipynb at [17], we load a precomputed $\bar{U}$ and
compute $\bar{u}=\bar{U}^{3}$, then we construct $A_{n}$ (and thus indirectly
$A$ and $T$) as outlined in Section 2.4, and we find the following sufficient
bounds:

$\displaystyle Y$ $\displaystyle=1.5917976375189734\times 10^{-23},$
$\displaystyle Z_{1}$ $\displaystyle=0.002929002447260958,$ $\displaystyle
Z_{5/3}$ $\displaystyle=19.560776011689786,$ $\displaystyle\underline{\delta}$
$\displaystyle=1.5964737129311592\times 10^{-23}$ $\displaystyle\bar{\delta}$
$\displaystyle=0.011508269024627918.$

We conclude with Theorem 20. Finally, positivity is deduced from Corollary 32.
∎

###  An example of a non-radial solution

Let us now turn to the problem of finding (real-valued) non-radial solutions
to problem (35). We limit ourselves to the case $\varepsilon=-1,d=2,p=3$ here
(though, such construction can be generalised), so that the problem remains
polynomial while allowing for non-necessarily positive solutions. In polar
coordinates, the equation for
$u(r,\vartheta)=e^{-r^{2}/8}\varphi(r\cos\vartheta,r\sin\vartheta)$ reads

$\mathcal{L}u=\frac{1}{2}u+u^{3},$ (46)

where

$\mathcal{L}u=-\frac{r}{2}\frac{\partial u}{\partial
r}-\Delta=-\frac{r}{2}\frac{\partial u}{\partial
r}-\frac{1}{r}\frac{\partial}{\partial r}\left(r\frac{\partial u}{\partial
r}\right)-\frac{1}{r^{2}}\frac{\partial^{2}u}{\partial\vartheta^{2}}.$

Observe that since the equation itself is radial, if
$u:\R_{+}\times\mathbb{T}\to\R$ is a solution to Eq. (46), then so is
$u(\,\cdot\,,\,\cdot+\vartheta_{0})$ for any phase
$\vartheta_{0}\in\mathbb{T}$. Thus, our approach which relies on a contraction
principle, and hence requires local uniqueness, cannot possibly work directly
on $H^{2}(\mu)$. By treating the problem in polar coordinates, instead of
using cartesian coordinates with the Hermite basis (11), we get a natural way
of dealing with the translation invariance in the angular variable.

We thus restrict the problem to an appropriate subspace, namely

$\mathcal{H}=\overline{\mathrm{Span}\left\\{\psi_{k,m}\,\mid\,k\in\mathbb{N},m\in\mathbb{N}\right\\}}^{H^{2}}\subset
H^{2}(\mu),$

where $\psi_{k,m}$ are as in Definition 5, and, as before, we denote by
$\hat{\psi}_{k,m}$ the normalisation of $\psi_{k,m}$ with respect to the
$L^{2}(\mu)$-inner product $\langle\cdot,\cdot\rangle$. Since we restricted
ourselves to nonnegative $k$’s, $\mathcal{H}$ is nothing but the space of
$H^{2}(\mu)$-functions which are even with respect to $\vartheta$. Solving for
Eq. (46) in this space addresses the above symmetry “issue” and it is
straightforward to check that

$\displaystyle F:\mathcal{H}$ $\displaystyle\longrightarrow\mathcal{H}$
$\displaystyle u$ $\displaystyle\longmapsto
u-\mathcal{L}^{-1}\left(\frac{u}{2}+u^{3}\right)$ (47)

is well-defined and compact. In this respect, our strategy goes along the same
line as the one employed in [4] to find non-radial solutions of analogous
semilinear equations on the disc. In the most general setting, the above-
defined space $\mathcal{H}$ shall be used, but one can simplify the proof
further if the solution has a zero projection onto certain frequencies. In our
example, we will perform our computer-assisted proof in the subspace

$\tilde{\mathcal{H}}=\overline{\mathrm{Span}\left\\{\psi_{2k+1,m}\,\mid\,k\in\mathbb{N},m\in\mathbb{N}\right\\}}^{H^{2}}\subset\mathcal{H}\subset
H^{2}(\mu).$

Working on such subspaces can also automatically discard other non-zero
solutions which are not of interest such as radial ones. Here, we have that

$\mathcal{L}\psi_{2k+1,m}=\lambda_{k+m}\psi_{2k+1,m},\qquad\lambda_{j}:=\frac{3}{2}+j,$

and the $\psi_{2k+1,m}$ for $k,m\in\mathbb{N}$ form a Hilbert basis of
$\tilde{\mathcal{H}}$. Note that the definition of $\lambda_{j}$ in that
context is slightly different from the one introduced in Notation 28. In
particular, this change reflects the fact that the first eigenvalue of
$\mathcal{L}$ restricted to $\tilde{\mathcal{H}}$ is strictly larger than
$d/2=1$, which slightly improves the Poincaré inequalities (12) and (13) on
$\tilde{\mathcal{H}}$.

The estimates for the radii polynomials coefficients are essentially the same
as in Section 4.1, with $P_{n}$ now denoting the projection onto
$\mathrm{Span}\left\\{\psi_{2k+1,m}\,\mid\,k\in\mathbb{N},m\in\mathbb{N}\text{
with }k+m\leq n\right\\}$, and

$w_{k,m}:=\|P_{\scriptscriptstyle\infty}(\bar{u}^{2}\hat{\psi}_{2k+1,m})\|_{L^{2}(\mu)}=\left(\|\bar{u}^{2}\hat{\psi}_{2k+1,m}\|^{2}_{L^{2}(\mu)}-\|P_{n}(\bar{u}^{2}\hat{\psi}_{2k+1,m})\|^{2}_{L^{2}(\mu)}\right)^{1/2}.$

The computation of integrals described in Section 3 is adapted as follows:

  * •

the integral is first reformulated as a sum of integrals of product of
functions, separable in $r$ and $\vartheta$ components so as to apply Fubini’s
Theorem,

  * •

the integrals with respect to $\vartheta$ are treated via the multiplication
formula for cosines and orthogonality,

  * •

the integrals with respect to $r$ are of the form (32) after the change of
variable $z=r^{2}/4$.

As for the $L^{\infty}$-bounds on $\hat{\psi}_{2k+1,m}$, we have that
$|\hat{\psi}_{2k+1,m}(r,\vartheta)|\leq|\hat{\psi}_{2k+1,m}(r,0)|$, and we
compute the local extrema of $r\mapsto\partial_{r}\hat{\psi}_{2k+1,m}(r,0)$ by
enclosing its zeros (which must all be real), this time using
PolynomialRoots.jl, the Julia implementation of [41].

######  Remark 36.

In practice, while we might do the computer-assisted proof with many basis
functions, we set some of the entries of $\bar{u}$ to zero to eliminate
unnecessary computational costs so that the validation can be performed within
a manageable amount of time on a personal computer.

We choose $n=150$ (i.e. 11,476 bases functions) and obtain the approximate
solution $\bar{u}$ depicted in Figure 4.

Figure 4: Plot of a non-radial numerical solution $\bar{u}$ to Eq. (46).

######  Theorem 37.

Let $\bar{u}:\mathbb{R}^{2}\to\mathbb{R}$ be the function whose restriction to
the disk $\\{(x,y)\mid x^{2}+y^{2}\leq 4^{2}\\}$ is represented on Fig 4 and
whose description in terms of “Fourier” coefficients in
$\\{\hat{\psi}_{2k+1,m}\\}_{k\geq 0,k+m\leq n}$ is available at [17]. There
exists a solution $u^{\star}\in\tilde{\mathcal{H}}$ to Eq. (46) such that

$\|u^{\star}-\bar{u}\|_{H^{2}(\mu)}\leq 0.0015107162514175158.$

###### Proof.

We consider $F:\tilde{\mathcal{H}}\to\tilde{\mathcal{H}}$ defined in (4.5). In
the file Heat/asymmetric/proof.ipynb at [17], we construct $A_{n}$ (and thus
indirectly $A$ and $T$) as outlined in Section 2.4, we find the following
sufficient bounds:

$\displaystyle Y$ $\displaystyle=0.0011482317939412424,$ $\displaystyle Z_{1}$
$\displaystyle=0.09887537542580953,$ $\displaystyle Z_{2}$
$\displaystyle=236.2041502678645,$ $\displaystyle Z_{3}$
$\displaystyle=456.22972624236917,$ $\displaystyle\bar{\delta}$
$\displaystyle=0.003801071181766125,$

and we conclude with Corollary 21. ∎

##  Self-similar solutions of the nonlinear Schrödinger equation

Let us now turn to the critical nonlinear Schrödinger equation

$i\partial_{t}v+\Delta v+\varepsilon|v|^{4/d}v=0,$ (48)

for $\varepsilon=\pm 1$. It is known [27] that this equation admits self-
similar solutions of the form

$v(t,x)=(1+t^{2})^{-d/4}\exp\left(\frac{it|x|^{2}}{4(1+t^{2})}\right)e^{i\omega\arctan{t}}\varphi\left(\frac{\sqrt{2}x}{(1+t^{2})^{1/2}}\right),$
(49)

for any $\omega\in\R$, as soon as $\varphi:\R^{d}\longrightarrow\mathbb{C}$
solves

$-\Delta\varphi+\frac{|x|^{2}}{16}\varphi+\frac{\omega}{2}\varphi=\frac{\varepsilon}{2}|\varphi|^{4/d}\varphi.$
(50)

Assuming $\varphi(x)=e^{|x|^{2}/8}u(x)$, we get the equation

$\mathcal{L}u=\left(\frac{d}{4}-\frac{\omega}{2}\right)u+\frac{\varepsilon}{2}e^{(p-1)|x|^{2}/8}|u|^{4/d}u,$
(51)

where $p=4/d+1$.

######  Remark 38.

In [27], Eq. (50) is solved via variational methods for $\varphi$ in the
Hilbert space

$\mathcal{X}=\left\\{\varphi\in H^{1}(\R^{d}):|x|\varphi\in
L^{2}(\R^{d})\right\\},$

but this is equivalent to solving Eq. (51) for $u\in H^{1}(e^{|x|^{2}/4})$.
Furthermore, an enclosure of a solution $u^{\star}$ to (51) around $\bar{u}$
with respect to the $H^{2}(\mu)$-norm will automatically yield an enclosure
for a solution $\varphi^{\star}:=e^{|x|^{2}/8}u^{\star}$ to (50) around
$\bar{\varphi}:=e^{|x|^{2}/8}\bar{u}$ with respect to the $H^{2}(\R^{d})$-norm
by inequality (20).

In the radial case, we thus aim to find a non-trivial zero of

$F(u)=u-\mathcal{L}^{-1}\left(\left(\frac{d}{4}-\frac{\omega}{2}\right)u+\frac{\varepsilon}{2}e^{(p-1)r^{2}/8}u^{p}\right).$
(52)

We get similar estimates for the radii polynomials coefficients as in Section
4.1 and take

$\displaystyle Y$
$\displaystyle=\left\\{\|\mathcal{L}A_{n}(P_{n}F(\bar{u}))\|^{2}_{L^{2}(\mu)}+\frac{1}{4}\left(\|e^{(p-1)r^{2}/8}\bar{u}^{p}\|^{2}_{L^{2}(\mu)}-\|P_{n}(e^{(p-1)r^{2}/8}\bar{u}^{p})\|^{2}_{L^{2}(\mu)}\right)\right\\}^{1/2}$
$\displaystyle Z^{11}$
$\displaystyle=\|\mathcal{L}(I_{n}-A_{n}P_{n}DF(\bar{u})P_{n})\mathcal{L}^{-1}\|_{L^{2}(\mu),L^{2}(\mu)}$
$\displaystyle Z^{21}$
$\displaystyle=\frac{p}{2}\|\mathcal{L}^{-1}w\|_{L^{2}(\mu)}$ $\displaystyle
Z^{12}$
$\displaystyle=\frac{p}{2\lambda_{n+1}}\left\|\,\left|\mathcal{L}A_{n}\mathcal{L}^{-1}\right|w\right\|_{L^{2}(\mu)}$
$\displaystyle Z^{22}$ $\displaystyle=\frac{\|d/4-\omega/2+\varepsilon
p(e^{r^{2}/8}\bar{u})^{p-1}/2\|_{\infty}}{\lambda_{n+1}}$ $\displaystyle
Z_{k}$
$\displaystyle=\frac{p!}{d(p-k)!}C(d)^{k-1}\|\mathcal{L}A\mathcal{L}^{-1}\|_{L^{2}(\mu),L^{2}(\mu)}\|e^{r^{2}/8}\bar{u}\|_{\infty}^{p-k}$
$\displaystyle\text{or }Z_{k}$
$\displaystyle=\frac{p!}{2(p-k)!}\left(\frac{2}{d}\right)^{k/2}c(d,k)\|\mathcal{L}A\mathcal{L}^{-1}\|_{L^{2}(\mu),L^{2}(\mu)}\|e^{r^{2}/8}\bar{u}\|_{\infty}^{p-k},$

where

$w_{m}:=\|P_{\scriptscriptstyle\infty}((e^{r^{2}/8}\bar{u})^{p-1}\hat{\psi}_{m})\|_{L^{2}(\mu)}=\left(\|(e^{r^{2}/8}\bar{u})^{p-1}\hat{\psi}_{m}\|^{2}_{L^{2}(\mu)}-\|P_{n}((e^{r^{2}/8}\bar{u})^{p-1}\hat{\psi}_{m})\|^{2}_{L^{2}(\mu)}\right)^{1/2}.$

Similarly, the quadrature rule of Section 3 can be adapted. For instance, the
matrix

$G_{ij}=\langle\hat{\psi}_{i},e^{(p-1)r^{2}/8}\bar{u}^{p-1}\hat{\psi}_{j}\rangle$

can be written as
$G=\bar{V}^{T}\mathrm{Diag}(\bar{V}\bar{\mathbf{u}})^{p-1}\bar{V}$ for
$\bar{V}$ as in Section 3, with $m=p+1$ and $\beta=(p+1)/2$. One could also
study the positivity of our radial solutions by applying Lemma 31 to Eq. (50).

Furthermore, one can find non-radial solutions by adapting the setup and the
above bounds analogously to Section 4.5.

(a) $\varepsilon=+1$, $d=2$, $\omega=2$

(b) $\varepsilon=-1$, $d=2$, $\omega=-2$

Figure 5: Plots of the numerical approximations $\bar{\varphi}=e^{|x|^{2}/8}\bar{u}$ to Eq. (50). |  $\varepsilon=+1$, $d=2$, $\omega=2$ |  $\varepsilon=-1$, $d=2$, $\omega=-2$ |  $\varepsilon=-1$, $d=2$, $\omega=-5/2$
---|---|---|---
$Y$ | $1.1097624368731403\times 10^{-7}$ | $1.510794398394117\times 10^{-33}$ | $7.709839958677216\times 10^{-9}$
$Z_{1}$ | $0.015947448213879025$ | $0.009444194874814784$ | $0.04065143211428217$
$Z_{2}$ | $77.79632835336696$ | $468.0575302222208$ | $109.71390801884884$
$Z_{3}$ | $134.26940590853712$ | $1264.711434525654$ | $486.14503594089473$
$\underline{\delta}$ | $1.1277557297872954\times 10^{-7}$ | $1.5251987173285673\times 10^{-33}$ | $8.036540516245275\times 10^{-9}$
$\bar{\delta}$ | $0.007787422383319831$ | $0.0021102955124448573$ | $0.008580958099091766$
$\eta$ | $1.8917745266839649\times 10^{-6}$ | $7.152376831616066\times 10^{-33}$ | $3.768713253012379\times 10^{-8}$
$\mathrm{dim}(P_{n}\mathcal{H})$ | $501$ | $501$ | $5151$
Figure | 5(a) | 5(b) | 2
Table 2: Results for the a posteriori validation of the solutions to Eq. (50)
in Theorem 39.

######  Theorem 39.

For $d=2$ and $(\varepsilon,\omega)=(+1,2)$, (resp.
$(\varepsilon,\omega)=(-1,-2)$, resp. $(\varepsilon,\omega)=(-1,-5/2)$), let
$\bar{\varphi}$ be the function depicted on a compact set in Figure 5(a)
(resp. 5(b), resp. 2) and whose description in terms of “Fourier” coefficients
is available at [17]. There exists a solution $\varphi^{\star}\in
H^{2}(\R^{d})$ to Eq. (50) such that

$\|\varphi^{\star}-\bar{\varphi}\|_{H^{2}(\R^{2})}\leq\eta,$

where $\eta$ is as in Table 2.

###### Proof.

We consider $F$ defined in (52). In the file Schrodinger/plus/proof.ipynb
(resp.
Schrodinger/plus/proof.ipynb, resp. Schrodinger/asymmetric/proof.ipynb) at
[17], we construct $A_{n}$ (and thus indirectly $A$ and $T$) as outlined in
Section 2.4 and following the estimates above, we obtain the upper bounds in
Table 2 for $Y$ and the $Z_{k}$’s to conclude with Corollary 21. Now, observe
that from combining (21) with the Poincaré inequalities (12) and (13), for
$d\geq 2$, for all $u\in H^{2}(\mu)$ and $\varphi=e^{|x|^{2}/8}u$, we have
that

$\|\varphi\|_{H^{2}(\R^{d})}^{2}:=\|D^{2}\varphi\|_{L^{2}(\R^{d})}^{2}+\|\nabla\varphi\|_{L^{2}(\R^{d})}^{2}+\|\varphi\|_{L^{2}(\R^{d})}^{2}\leq
Z\left[1+\max\left(0,1-\frac{d}{8}\right)\frac{4}{d^{2}}\right]\|u\|_{H^{2}(\mu)}^{2}.$

And thus, for $d=2$, we may choose $\eta=(\sqrt{7Z}/2)\underline{\delta}$. ∎

Note that the third case of Theorem 39 corresponds to Theorem 2.

##  An example with first order terms

In this final section, we consider a toy problem whose interest lies in the
fact that it leads to an equation of the form (2) with a first-order term,
i.e. a right-hand side of the form $f(x,u,\nabla u)$ with an actual dependence
on $\nabla u$. We focus on a generalised viscous Burgers equation on $\R_{+}$
with Neumann boundary condition

$\begin{cases}\partial_{t}v+v^{2}\partial_{x}v=\partial_{xx}v,\\\
\partial_{x}v(t,0)=0\qquad\qquad\quad\mbox{for all $t>0.$}\end{cases}$ (53)

Then, under the similarity transformation

$v(t,x)=t^{-1/4}u\left(\frac{x}{\sqrt{t}}\right),$

$v$ solves Eq. (53) for $t>0$ if $u$ solves

$\mathcal{L}u-\frac{u}{4}+u^{2}\frac{\partial u}{\partial x}=0,\qquad x\geq
0.$ (54)

Similar problems have been studied in [2, 43] and the references therein. This
equation can be approached with our method, as $u\in H^{2}(\mu)$ implies that
$u\in\mathcal{C}^{1}_{0}(\R)$ and thus $u^{2}\partial_{x}u\in L^{2}(\mu)$.

######  Remark 40.

We deal with a one dimensional example here for simplicity, but terms like
$(u\cdot\nabla)u$ could in principle also be handled in dimension
$d\in\\{2,3\\}$, as $u\in H^{2}(\mu)$ is then still enough to guarantee that
$(u\cdot\nabla)u\in L^{2}(\mu)$ since $u\in L^{\infty}(\R^{d})$.

In order to enforce the Neumann boundary condition, we symmetrise the problem
and work on the Hilbert space
$\mathcal{H}=\overline{\mathrm{Span}\\{\hat{\psi}_{m}\\}_{m=0}^{\infty}}^{H^{2}}$,
i.e. we make use of “even half-Hermite series”. Note that working on
$\mathcal{H}$ is equivalent to working on $H^{2}(\R_{+},e^{x^{2}/4})$ with the
additional Neumann boundary condition (see also [53, Section 3]). Furthermore,
for clarity, in what follows, all elements of $L^{2}(\mu)$ that we will
consider should be treated as functions on $\R_{+}$ having an even extension
on $\R$. In particular, for such functions $u$ and $v$, we have that

$\langle
u,v\rangle=\int_{-\infty}^{\infty}uv\mu(x)\mathrm{d}x=\frac{1}{2}\int_{-\infty}^{\infty}uve^{x^{2}/4}\mathrm{d}x=\int_{0}^{\infty}uve^{x^{2}/4}\mathrm{d}x.$

We consider

$\displaystyle F:\mathcal{H}$ $\displaystyle\longrightarrow\mathcal{H}$
$\displaystyle u$ $\displaystyle\longmapsto
u+\mathcal{L}^{-1}\left(-\frac{u}{4}+u^{2}\partial_{x}u\right),$ (55)

and construct our fixed point operator $T:\mathcal{H}\to\mathcal{H}$ as
outlined in 2.4. Note that $T$ is still polynomial with respect to $u$ in this
case and most bounds for the validation can be derived similarly as in Section
4.1:

$\displaystyle Y$
$\displaystyle=\left\\{\|\mathcal{L}A_{n}(P_{n}F(\bar{u}))\|^{2}_{L^{2}(\mu)}+\left(\|\bar{u}^{2}\partial_{x}\bar{u}\|^{2}_{L^{2}(\mu)}-\|P_{n}(\bar{u}^{2}\partial_{x}\bar{u})\|^{2}_{L^{2}}\right)\right\\}^{1/2},$
$\displaystyle Z^{11}$
$\displaystyle=\|\mathcal{L}(I_{n}-A_{n}P_{n}DF(\bar{u})P_{n})\mathcal{L}^{-1}\|_{L^{2}(\mu),L^{2}(\mu)},$
$\displaystyle Z^{21}$ $\displaystyle=\|\mathcal{L}^{-1}w\|_{L^{2}(\mu)},$
$\displaystyle Z^{22}$
$\displaystyle=\frac{\|1/4+2\bar{u}\partial_{x}\bar{u}\|_{\infty}}{\lambda_{n+1}}+\frac{\|\bar{u}\|_{\infty}^{2}}{\sqrt{\lambda_{n+1}}},$
$\displaystyle Z_{2}$
$\displaystyle=2^{15/4}\|\mathcal{L}A\mathcal{L}^{-1}\|_{L^{2}(\mu),L^{2}(\mu)}(\sqrt{2}\|\bar{u}\|_{\infty}+\|\partial_{x}\bar{u}\|_{\infty}),$
$\displaystyle Z_{3}$
$\displaystyle=96\|\mathcal{L}A\mathcal{L}^{-1}\|_{L^{2}(\mu),L^{2}(\mu)},$

with

$\displaystyle w_{m}:$
$\displaystyle=2\|P_{\infty}(\bar{u}\hat{\psi}_{m}\partial_{x}\bar{u})\|_{L^{2}}+\|P_{\infty}(\bar{u}^{2}\partial_{x}\hat{\psi}_{m})\|_{L^{2}}$
$\displaystyle=2\left(\|\bar{u}\hat{\psi}_{m}\partial_{x}\bar{u}\|^{2}_{L^{2}}-\|P_{n}(\bar{u}\hat{\psi}_{m}\partial_{x}\bar{u})\|^{2}_{L^{2}}\right)^{1/2}+\left(\|\bar{u}^{2}\partial_{x}\hat{\psi}_{m}\|^{2}_{L^{2}}-\|P_{n}(\bar{u}^{2}\partial_{x}\hat{\psi}_{m})\|^{2}_{L^{2}}\right)^{1/2}.$

The only estimate which requires substantial adaptations is $Z^{12}$. For
$h_{\scriptscriptstyle\infty}\in\mathcal{H}_{\scriptscriptstyle\infty}$, we
get

$\displaystyle\langle
2\bar{u}\partial_{x}\bar{u}h_{\scriptscriptstyle\infty}+\bar{u}^{2}\partial_{x}h_{\scriptscriptstyle\infty},\hat{\psi}_{m}\rangle$
$\displaystyle=\langle\partial_{x}(\bar{u}^{2}h_{\scriptscriptstyle\infty}),\hat{\psi}_{m}\rangle$
$\displaystyle=\int_{0}^{\infty}\partial_{x}(\bar{u}^{2}h_{\scriptscriptstyle\infty})\hat{\psi}_{m}e^{x^{2}/4}\mathrm{d}x$
$\displaystyle=\bar{u}(0)^{2}h_{\scriptscriptstyle\infty}(0)\hat{\psi}_{m}(0)-\int_{0}^{\infty}\bar{u}^{2}h_{\scriptscriptstyle\infty}\partial_{x}(\hat{\psi}_{m}e^{x^{2}/4})\mathrm{d}x\qquad\mbox{by
parts.}$ (56)

We deal with both terms appearing in (6) separately. First, writing
$h_{\scriptscriptstyle\infty}$ as

$h_{\scriptscriptstyle\infty}=\sum_{m=n+1}^{\infty}a_{m}\hat{\psi}_{m},$

we have that

$\displaystyle h_{\scriptscriptstyle\infty}(0)^{2}$
$\displaystyle=\left(\sum_{m=n+1}^{\infty}a_{m}\hat{\psi}_{m}(0)\right)^{2}$
$\displaystyle\leq\left(\sum_{m=n+1}^{\infty}\lambda_{m}^{2}a_{m}^{2}\right)\left(\sum_{m=n+1}^{\infty}\frac{\hat{\psi}_{m}(0)^{2}}{\lambda_{m}^{2}}\right)\qquad\mbox{by
Cauchy--Schwarz inequality}$
$\displaystyle=\frac{\|h_{\scriptscriptstyle\infty}\|_{H^{2}(\mu)}^{2}}{\Gamma(1/2)^{2}}\sum_{m=n+1}^{\infty}\frac{\Gamma(m+1/2)}{(m+1/2)^{2}m!}\qquad\mbox{by\leavevmode\nobreak\
\eqref{eqn:tail-estimate}}$
$\displaystyle=\frac{\|h_{\scriptscriptstyle\infty}\|_{H^{2}(\mu)}^{2}}{\pi}\left(\sum_{m=0}^{\infty}\frac{\Gamma(m+1/2)}{(m+1/2)^{2}m!}-\sum_{m=0}^{n}\frac{\Gamma(m+1/2)}{(m+1/2)^{2}m!}\right)$
$\displaystyle\leq\frac{\|h_{\scriptscriptstyle\infty}\|_{H^{2}(\mu)}^{2}}{\pi}\left(\pi^{3/2}\log{4}-\sum_{m=0}^{n}\frac{\Gamma(m+1/2)}{(m+1/2)^{2}m!}\right)$
$\displaystyle\leq\|h_{\scriptscriptstyle\infty}\|_{H^{2}(\mu)}^{2}S_{n},$

where

$S_{n}=\frac{1}{\pi}\left(\pi^{3/2}\log{4}-\sum_{m=0}^{n}\frac{\Gamma(m+1/2)}{(m+1/2)^{2}m!}\right).$

Now, for the integral term in (6),

$\displaystyle\left|\int_{0}^{\infty}\bar{u}^{2}h_{\scriptscriptstyle\infty}\partial_{x}(\hat{\psi}_{m}e^{x^{2}/4})\mathrm{d}x\right|$
$\displaystyle=\left|\int_{0}^{\infty}\bar{u}^{2}h_{\scriptscriptstyle\infty}\left(\partial_{x}\hat{\psi}_{m}+\frac{x}{2}\hat{\psi}_{m}\right)e^{x^{2}/4}\mathrm{d}x\right|$
$\displaystyle\leq\left|\left\langle
h_{\scriptscriptstyle\infty},\bar{u}^{2}\left(\partial_{x}\hat{\psi}_{m}+\frac{x}{2}\hat{\psi}_{m}\right)\right\rangle\right|$
$\displaystyle\leq\left|\left\langle
h_{\scriptscriptstyle\infty},P_{\scriptscriptstyle\infty}\left(\bar{u}^{2}\left(\partial_{x}\hat{\psi}_{m}+\frac{x}{2}\hat{\psi}_{m}\right)\right)\right\rangle\right|\qquad\mbox{as
$h_{\scriptscriptstyle\infty}\in\mathcal{H}_{\scriptscriptstyle\infty}$}$
$\displaystyle\leq\|h_{\scriptscriptstyle\infty}\|_{L^{2}(\mu)}\left\|P_{\scriptscriptstyle\infty}\left(\bar{u}^{2}\left(\partial_{x}\hat{\psi}_{m}+\frac{x}{2}\hat{\psi}_{m}\right)\right)\right\|_{L^{2}(\mu)}$
$\displaystyle\leq\frac{\|h_{\scriptscriptstyle\infty}\|_{H^{2}(\mu)}}{\lambda_{n+1}}\left(\left\|\bar{u}^{2}\left(\partial_{x}\hat{\psi}_{m}+\frac{x}{2}\hat{\psi}_{m}\right)\right\|_{L^{2}(\mu)}^{2}-\left\|P_{n}\left(\bar{u}^{2}\left(\partial_{x}\hat{\psi}_{m}+\frac{x}{2}\hat{\psi}_{m}\right)\right)\right\|_{L^{2}(\mu)}^{2}\right)^{1/2}.$

Thus, with

$\tilde{w}_{m}=\bar{u}(0)^{2}\hat{\psi}_{m}(0)\sqrt{S_{n}}+\frac{1}{\lambda_{n+1}}\left(\left\|\bar{u}^{2}\left(\partial_{x}\hat{\psi}_{m}+\frac{x}{2}\hat{\psi}_{m}\right)\right\|_{L^{2}(\mu)}^{2}-\left\|P_{n}\left(\bar{u}^{2}\left(\partial_{x}\hat{\psi}_{m}+\frac{x}{2}\hat{\psi}_{m}\right)\right)\right\|_{L^{2}(\mu)}^{2}\right)^{1/2},$

we can choose

$Z^{12}=\left\|\,\left|\mathcal{L}A_{n}\mathcal{L}^{-1}\right|\tilde{w}\right\|_{L^{2}(\mu)}.$

An estimate based on Parseval’s identity on $H^{1}(\mu)$ is also possible, but
gives very similar results.

######  Remark 41.

Note that if $u(x)=q(x)e^{-x^{2}/4}$ for a polynomial $q$, then
$\partial_{x}u(x)=(\partial_{x}q(x)+x/2q(x))e^{-x^{2}/4}$ and thus the
integrals necessary to compute our bounds can also be evaluated with
quadrature rules as in Section 3. Concerning $L^{\infty}$-bounds, from [1,
Formulae 22.14.15–16], we have that for all $m\in\mathbb{N}$

$\displaystyle\|\hat{\psi}_{m}e^{r^{2}/8}\|_{\infty}$
$\displaystyle\leq\frac{m!}{\pi^{1/4}2^{m}\sqrt{(2m)!}}\left[2^{2m+1}-\left(\begin{matrix}2m\\\
m\\\ \end{matrix}\right)\right],$
$\displaystyle\|\partial_{x}\hat{\psi}_{m}\|_{\infty}$
$\displaystyle\leq\frac{\sqrt{2m+1}\sqrt{(2m+1)!}}{\sqrt{e}\pi^{1/4}m!}.$

With $n=1500$, we find the approximate solution $\bar{u}$ represented in
Figure 6.

Figure 6: Plot of the numerical solution $\bar{u}$ to Eq. (54).

######  Theorem 42.

Let $\bar{u}:\R_{+}\to\R$ be the function represented on the interval $[0,4]$
in Fig 6 and whose description in terms of “Fourier” coefficients in
$\\{\hat{\psi}_{m}\\}_{m=0}^{n}$ is available at [17], then there exists a
positive solution $u^{\star}\in\mathcal{H}$ to Eq. (54) such that

$\|u^{\star}-\bar{u}\|_{H^{2}(\mu)}\leq 0.0009686954373860599.$

###### Proof.

We consider $F:\mathcal{H}\to\mathcal{H}$ defined in (55). In the file
Burger/proof.ipynb at [17], we construct $A_{n}$ (and thus indirectly $A$ and
$T$) as outlined in Section 2.4 and following the estimates above, we find the
following sufficient bounds:

$\displaystyle Y$ $\displaystyle=0.0007563639234949774,$ $\displaystyle Z_{1}$
$\displaystyle=0.08518340287046657,$ $\displaystyle Z_{2}$
$\displaystyle=276.49885111394093,$ $\displaystyle Z_{3}$
$\displaystyle=556.4779967366688,$ $\displaystyle\bar{\delta}$
$\displaystyle=0.002894449581015444,$

and we conclude with Corollary 21 to the existence of a solution
$u^{\star}\in\mathcal{H}$ to Eq. (54). As Eq. (54) can be seen as an ordinary
differential equation with smooth vector field, we immediately get that
$u^{\star}\in\mathcal{C}^{\infty}(\R_{+})$ and thus we can apply Lemma 31
(again to $\varphi=e^{|x|^{2}/8}u^{\star}$) with $b(x)=u^{2}(x)$ and
$c(x)=x^{2}/16-xu^{2}(x)/4$, to get that $u^{\star}$ is positive. ∎

## Acknowledgments

HC thanks Jeroen Lamb and Martin Rasmussen for their support and guidance. The
authors thank Miguel Escobedo, Olivier Hénot, Marc Nualart Batalla and Sheehan
Olver for useful discussions. HC was supported by a scholarship from the
Imperial College London EPSRC DTP in Mathematical Sciences (EP/W523872/1) and
by the EPSRC Centre for Doctoral Training in Mathematics of Random Systems:
Analysis, Modelling and Simulation (EP/S023925/1). MB was supported by the ANR
project CAPPS: ANR-23-CE40-0004-01, and also acknowledges the hospitality of
the Department of Mathematics of Imperial College London, and the support of
the CNRS-Imperial _Abraham de Moivre_ IRL.

## References

  * [1] Milton Abramowitz and I A Stegun  “Handbook of Mathematical Functions: With Formulas, Graphs, and Mathematical Tables Applied mathematics series”  In _National Bureau of Standards, Washington, DC_ , 1970
  * [2] J. Aguirre, M. Escobedo and E. Zuazua  “Self-Similar Solutions of a Convection Diffusion Equation And Related Semilinear Elliptic Problems”  In _Communications in Partial Differential Equations_ 15.2, 1990, pp. 139–157  DOI: 10.1080/03605309908820681
  * [3] G. Arioli and H. Koch  “Computer-assisted methods for the study of stationary solutions in dissipative systems, applied to the Kuramoto-Sivashinski equation”  In _Arch. Ration. Mech. Anal._ 197.3, 2010, pp. 1033–1051
  * [4] Gianni Arioli and Hans Koch  “Non-radial solutions for some semilinear elliptic equations on the disk”  In _Nonlinear Analysis, Theory, Methods and Applications_ 179 Elsevier Ltd, 2019, pp. 294–308  DOI: 10.1016/j.na.2018.09.001
  * [5] Ben-yu Guo and Cheng-long Xu  “Hermite pseudospectral method for nonlinear partial differential equations”  In _ESAIM: Mathematical Modelling and Numerical Analysis_ 34.4, 2000, pp. 859–872  DOI: 10.1051/m2an:2000100
  * [6] Zachary Bradshaw and Tai-Peng Tsai  “Self-similar solutions to the Navier-Stokes equations: a survey of recent results”, 2018, pp. 1–25  URL: http://arxiv.org/abs/1802.00038
  * [7] Maxime Breden  “Computer-assisted proofs for some nonlinear diffusion problems”  In _Communications in Nonlinear Science and Numerical Simulation_ 109, 2022, pp. 106292  DOI: 10.1016/j.cnsns.2022.106292
  * [8] Maxime Breden and Christian Kuehn  “Rigorous validation of stochastic transition paths”  In _Journal de Mathématiques Pures et Appliquées_ 131 Elsevier Masson SAS, 2019, pp. 88–129  DOI: 10.1016/j.matpur.2019.04.012
  * [9] Maxime Breden and Maxime Payan  “Computer-assisted proofs for the many steady states of a chemotaxis model with local sensing”  In _arXiv preprint arXiv:2311.13896_ , 2023
  * [10] Haim Brezis  “Functional Analysis, Sobolev Spaces and Partial Differential Equations”  New York, NY: Springer New York, 2011  DOI: 10.1007/978-0-387-70914-7
  * [11] Haim Brezis, Lambertus Adrianus Peletier and David Terman  “A very singular solution of the heat equation with absorption”  In _Archive for Rational Mechanics and Analysis_ 95.3, 1986, pp. 185–209  DOI: 10.1007/BF00251357
  * [12] Tristan Buckmaster, Gonzalo Cao-Labora and Javier Gómez-Serrano  “Smooth imploding solutions for 3D compressible fluids”  In _arXiv preprint arXiv:2208.09445_ , 2022
  * [13] Matthieu Cadiot, Jean-Philippe Lessard and Jean-Christophe Nave  “Rigorous computation of solutions of semi-linear PDEs on unbounded domains via spectral methods”  In _arXiv preprint arXiv:2302.12877_ , 2023
  * [14] Matthieu Cadiot, Jean-Philippe Lessard and Jean-Christophe Nave  “Stationary non-radial localized patterns in the planar Swift-Hohenberg PDE: constructive proofs of existence”  In _arXiv preprint arXiv:2403.10450_ , 2024
  * [15] Jiajie Chen and Thomas Y Hou  “Stable nearly self-similar blowup of the 2D Boussinesq and 3D Euler equations with smooth data II: Rigorous Numerics”  In _arXiv preprint arXiv:2305.05660_ , 2023
  * [16] Yong Chen and Yong Liu  “On the eigenfunctions of the complex Ornstein–Uhlenbeck operators”  In _Kyoto Journal of Mathematics_ 54.3, 2014, pp. 577–596  DOI: 10.1215/21562261-2693451
  * [17] Hugo Chu  “Codes associated to “Constructive proofs for some semilinear PDEs on $H^{2}(e^{|x|^{2}/4},\mathbb{R}^{d})$”” Available at https://github.com/Huggzz/Hermite-Laguerre_proofs/tree/main, 2024
  * [18] S. Day, J.-P. Lessard and K. Mischaikow  “Validated continuation for equilibria of PDEs”  In _SIAM J. Numer. Anal._ 45.4 SIAM, 2007, pp. 1398–1424
  * [19] Miguel Escobedo and Otared Kavian  “Variational problems related to self-similar solutions of the heat equation”  In _Nonlinear Analysis: Theory, Methods & Applications_ 11.10, 1987, pp. 1103–1133  DOI: 10.1016/0362-546X(87)90001-0
  * [20] Daniele Funaro  “Computational aspects of pseudospectral Laguerre approximations”  In _Applied Numerical Mathematics_ 6.6, 1990, pp. 447–457  DOI: 10.1016/0168-9274(90)90003-X
  * [21] Daniele Funaro and Otared Kavian  “Approximation of Some Diffusion Evolution Equations in Unbounded Domains by Hermite Functions”  In _Mathematics of Computation_ 57.196, 1991, pp. 597  DOI: 10.2307/2938707
  * [22] Javier Gómez-Serrano  “Computer-assisted proofs in PDE: a survey”  In _SeMA Journal_ 76.3 Springer, 2019, pp. 459–484
  * [23] Olivier Hénot  “On polynomial forms of nonlinear functional differential equations”  In _Journal of Computational Dynamics_ 8.3 Journal of Computational Dynamics, 2021, pp. 307–323
  * [24] Atsushi Higuchi  “Symmetric tensor spherical harmonics on the N -sphere and their application to the de Sitter group $SO(N,1)$”  In _Journal of Mathematical Physics_ 28.7, 1987, pp. 1553–1566  DOI: 10.1063/1.527513
  * [25] Kiyosi Itô  “Complex Multiple Wiener Integral”  In _Japanese journal of mathematics :transactions and abstracts_ 22, 1952, pp. 63–86  DOI: 10.4099/jjm1924.22.0–“˙˝63
  * [26] Otared Kavian  “Remarks on the large time behaviour of a nonlinear diffusion equation”  In _Annales de l’Institut Henri Poincaré C, Analyse non linéaire_ 4.5, 1987, pp. 423–452  DOI: 10.1016/s0294-1449(16)30358-4
  * [27] Otared Kavian and Fred B. Weissler  “Self-similar solutions of the pseudo-conformally invariant nonlinear Schrödinger equation”  In _Michigan Mathematical Journal_ 41.1, 1994  DOI: 10.1307/mmj/1029004922
  * [28] Shane Kepley and J.. Mireles James  “Chaotic motions in the restricted four body problem via Devaney’s saddle-focus homoclinic tangle theorem”  In _Journal of Differential Equations_ Elsevier, 2018
  * [29] J.. Kurtz  “Weighted Sobolev spaces with applications to singular nonlinear boundary value problems”  In _Journal of Differential Equations_ 49.1, 1983, pp. 105–123  DOI: 10.1016/0022-0396(83)90021-9
  * [30] J.-P. Lessard, J.. Mireles James and J. Ransford  “Automatic differentiation for Fourier series and the radii polynomial approach”  In _Physica D: Nonlinear Phenomena_ 334 Elsevier, 2016, pp. 174–186
  * [31] K Nagatou, MT Nakao and M Wakayama  “Verified numerical computations for eigenvalues of non-commutative harmonic oscillators”  In _Numerical Functional Analysis and Optimization_ 23.5-6 Taylor & Francis, 2002, pp. 633–650
  * [32] Kaori Nagatou, Michael Plum and Mitsuhiro T Nakao  “Eigenvalue excluding for perturbed-periodic one-dimensional Schrödinger operators”  In _Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences_ 468.2138 The Royal Society Publishing, 2012, pp. 545–562
  * [33] Yūki Naito and Takashi Suzuki  “Radial Symmetry of Self-Similar Solutions for Semilinear Heat Equations”  In _Journal of Differential Equations_ 163.2, 2000, pp. 407–428  DOI: 10.1006/jdeq.1999.3742
  * [34] Mitsuhiro T. Nakao  “A numerical approach to the proof of existence of solutions for elliptic problems”  In _Japan Journal of Applied Mathematics_ 5.2, 1988, pp. 313–332  DOI: 10.1007/BF03167877
  * [35] Mitsuhiro T. Nakao, Michael Plum and Yoshitaka Watanabe  “Numerical Verification Methods and Computer-Assisted Proofs for Partial Differential Equations” 53, Springer Series in Computational Mathematics  Singapore: Springer Singapore, 2019  DOI: 10.1007/978-981-13-7669-6
  * [36] Shin’ichi Oishi  “Numerical verification of existence and inclusion of solutions for nonlinear operator equations”  In _Journal of Computational and Applied Mathematics_ 60.1-2 Elsevier, 1995, pp. 171–185
  * [37] Grigorios A. Pavliotis  “Stochastic Processes and Applications” 60, Texts in Applied Mathematics  New York, NY: Springer New York, 2014  DOI: 10.1007/978-1-4939-1323-7
  * [38] L A Peletier, D Terman and F B Weissler  “On the equation $\Delta u+\tfrac{x}{2}\cdot\nabla u+f(u)=0$”  In _Archive for Rational Mechanics and Analysis_ 94.1, 1986, pp. 83–99  DOI: 10.1007/BF00278244
  * [39] Michael Plum  “Explicit $H^{2}$-estimates and pointwise bounds for solutions of second-order elliptic boundary value problems”  In _Journal of Mathematical Analysis and Applications_ 165.1 Elsevier, 1992, pp. 36–61
  * [40] David P. Sanders et al.  “JuliaIntervals/IntervalArithmetic.jl: v0.16.7”, 2020  DOI: 10.5281/ZENODO.3727070
  * [41] Jan Skowron and Andrew Gould  “General Complex Polynomial Root Solver and Its Further Optimization for Binary Microlenses”, 2012  URL: http://www.astrouw.edu.pl/
  * [42] C Spina, V Manco and Giorgio Metafune  “Equazioni Ellittiche del Secondo Ordine, Parte Seconda: Teoria $L^{p}$”  Università di Lecce - Coordinamento SIBA, 2005
  * [43] Ch. Srinivasa Rao, P.L. Sachdev and Mythily Ramaswamy  “Self-similar solutions of a generalized Burgers equation with nonlinear damping”  In _Nonlinear Analysis: Real World Applications_ 4.5, 2003, pp. 723–741  DOI: 10.1016/S1468-1218(02)00083-4
  * [44] Akitoshi Takayasu, Xuefeng Liu and Shin’ichi Oishi  “Verified computations to semilinear elliptic boundary value problems on arbitrary polygonal domains”  In _Nonlinear Theory and Its Applications, IEICE_ 4.1 The Institute of Electronics, InformationCommunication Engineers, 2013, pp. 34–61
  * [45] Giorgio Talenti  “Best constant in Sobolev inequality”  In _Annali di Matematica Pura ed Applicata_ 110.1, 1976, pp. 353–372  DOI: 10.1007/BF02418013
  * [46] Alex Townsend and Sheehan Olver  “FastGaussQuadrature.jl”  GitHub, 2014  URL: https://github.com/JuliaApproximation/FastGaussQuadrature.jl
  * [47] J.. Berg and J.-P. Lessard  “Rigorous numerics in dynamics”  In _Notices Amer. Math. Soc._ 62.9, 2015
  * [48] J.. Berg and J.. Williams  “Validation of the bifurcation diagram in the 2D Ohta–Kawasaki problem”  In _Nonlinearity_ 30.4 IOP Publishing, 2017, pp. 1584
  * [49] Jan Bouwe Berg, Maxime Breden, Jean-Philippe Lessard and Lennaert Veen  “Spontaneous periodic orbits in the Navier–Stokes flow”  In _Journal of Nonlinear Science_ 31.2 Springer, 2021, pp. 1–64
  * [50] Jan Bouwe Berg, Chris M Groothedde and JF Williams  “Rigorous computation of a radially symmetric localized solution in a Ginzburg–Landau problem”  In _SIAM Journal on Applied Dynamical Systems_ 14.1 SIAM, 2015, pp. 423–447
  * [51] Jan Bouwe Berg, Olivier Hénot and Jean-Philippe Lessard  “Constructive proofs for localised radial solutions of semilinear elliptic systems on $\mathbb{R}^{d}$”  In _Nonlinearity_ 36.12 IOP Publishing, 2023, pp. 6476
  * [52] Thomas Wanner  “Computer-assisted equilibrium validation for the diblock copolymer model”  In _Discrete & Continuous Dynamical Systems-A_ 37.2 American Institute of Mathematical Sciences, 2017, pp. 1075
  * [53] Fred B. Weissler  “Rapidly decaying solutions of an ordinary differential equation with applications to semilinear elliptic and parabolic partial differential equations”  In _Archive for Rational Mechanics and Analysis_ 91.3, 1986, pp. 247–266  DOI: 10.1007/BF00250744
  * [54] Jonathan Wunderlich  “Computer-assisted Existence Proofs for Navier-Stokes Equations on an Unbounded Strip with Obstacle”, 2022
  * [55] N. Yamamoto  “A numerical verification method for solutions of boundary value problems with local uniqueness by Banach’s fixed-point theorem”  In _SIAM J. Numer. Anal._ 35, 1998, pp. 2004–2013
  * [56] Eiji Yanagida  “Uniqueness of Rapidly Decaying Solutions to the Haraux–Weissler Equation”  In _Journal of Differential Equations_ 127.2, 1996, pp. 561–570  DOI: 10.1006/jdeq.1996.0083
  * [57] Piotr Zgliczynski and Konstantin Mischaikow  “Rigorous numerics for partial differential equations: The Kuramoto—Sivashinsky equation”  In _Foundations of Computational Mathematics_ 1.3 Springer, 2001, pp. 255–288