# On the Application of Tschirnhaus Transformations to the Reduction of
Algebraic Equations

A Translation by Alexander J. Sutherland

##  1 Original Bibliographic Information

###  1.1 German

  * •

Author: Anders Wiman

  * •

Title: Über die Anwendung der Tschirnhausen-Transformation auf die Reduktion
algebraischer Gleichungen

  * •

Year: 1927

  * •

Language: German

  * •

Publisher: Nova Acta Regiae Societatis Scientiarum Upsaliensis (Nova Acta R.
Soc. scient. Uppsala)

  * •

Note: Der Königi. Societät der Wissenschaften Zu Uppsala Mitgeteilt

###  1.2 English

  * •

Author: Anders Wiman

  * •

Title: On the Application of Tschirnhaus Transformations to the Reduction of
Algebraic Equations

  * •

Year: 1927

  * •

Language: German

  * •

Publisher: New Proceedings of the Royal Society of Scientists of Uppsala

  * •

Note: Notice to the Royal Society of Scientists of Uppsala on 06 May 1927

This work was supported by the National Science Foundation under Grant No.
DMS-1944862.

##  2 The Translation

### Part 1

111Translator’s Note: This is a translation of the original mathematics. In
particular, errors in the text have not been fixed. The errors in question
come from considering intersections in affine spaces instead of in projective
spaces. Throughout this translation, there are additional footnotes with the
identifier “Translator’s Footnote:.” These footnotes refer to remarks in
Section 3 in which the translator provides additional mathematical commentary.

We consider a general equation of $n^{th}$ degree:

$x^{n}+c_{1}x^{n-1}+\cdots+c_{n}=0$ (1)

with roots $x_{1},\dotsc,x_{n}$. We then apply a Tschirnhaus transformation,
which has the general form

$y=a_{0}+a_{1}x+\cdots+a_{n-1}x^{n-1}.$ (2)

This will transform (1) into an equation

$y^{n}+C_{1}y^{n-1}+\cdots+C_{n}=0,$ (3)

where the coefficients $C_{i}$ are all homogeneous functions of degree $i$ in
the variables $a_{\nu}$ and include all such functions of weight up to $i$ in
the coefficents $c_{i}$. Tschirnhaus hoped to use this type of transformation
to convert equation (1) to the binomial form in such a way that the
determination of parameters $a_{\nu}$ should require the solution of equations
of degree at most $n-1$. As is well known, this is not the case, even though
there have many attempts for the degree 5; this will never be the case, just
like for [the problem of] trisecting an angle. However, the situation is
completely different if the problem is formulated in the following way: Is it
possible to satisfy the conditions

$C_{i}(a_{0},a_{1},\dotsc,a_{n-1})=0,\quad i=1,\dotsc,m$ (4)

in such a way that determining the parameters $a_{\nu}$ only requires
equations of degree up to $m$, when $n$ is sufficiently large? As a result of
the following treatment for the case $m=4$, it should not appear doubtful that
question should also be decided in the affirmative for larger $m$. However,
the general problem of determining the lower bound on $n$ associated to each
$m$ appears to be very complex.

### Section 2

Observe that

$\displaystyle C_{1}(a_{0},\dotsc,a_{n-1})=na_{0}+\cdots.$

If $C_{1}=0$, then $a_{0}$ is expressed linearly in the other parameters. The
coefficients $C_{i},(i=2,\dotsc,n-1)$ are then homogeneous functions of degree
$i$ in the parameters $a_{1},\dotsc,a_{n-1}$. In order to [find a point that
will] satisfy a single condition

$C_{x}=0,\quad(x>1)$

it is evident that it is only necessary to find an intersection of the
hypersurface $C_{x}=0$ with an arbitrary straight line to solve an equation of
degree $x$. If all the roots of (1) are real, then you cannot get a real
solution for $x=2$ because the hypersurface

$\sum\limits_{i=1}^{n}y_{i}^{2}=C_{2}=0$ (5)

has only the trivial 0. In contrast, there are always real points on the
surface $C_{3}=0$. Indeed, as you can easily see, the same is the case for all
surfaces $C_{x}=0\ (x>2)$ if $x$ is an even number.

If one has $n\geq 5$, then one obtains the solution of (4) for $m=3$ by the
well-known Bring-Jerrard transformation, which is illustrated geometrically by
F. Klein 222We refer to the in-depth treatment of F. Klein, Lectures on the
Icosahedron and the Solution of the Equation of Fifth Degree, Leipzig, 1884.
in the following way.333Translator’s footnote: See Remark 1 for the
translator’s summary of this argument. First, a point $P$ is obtained on the
surface $C_{2}=0$, which, as noted above, can be done using a square root. If
$n>5$, we then consider a three-dimensional space $R_{3}$ which is tangent [to
$C_{2}=0$] at $P$, which then has a hypersuface of degree 2 in common with
$C_{2}=0$. A second square root is now required in order to select one of the
two generators of this surface going through $P$. Determining an intersection
of one of these generators with $F_{3}=0$ requires the solution of a degree 3
equation. Although only one pair of imaginary roots need to occur in the real
equations $C_{1}=C_{2}=C_{3}=0$, at least two pairs of imaginary roots must
exist to actually execute this transformation. Otherwise, there is no real
line at $C_{2}=0$. This is due to the fact that if one converts

$\displaystyle C_{2}(a_{1},\dotsc,a_{n-1})=\sum\limits_{i=1}^{n}y_{i}^{2}$

to a sum of $n-1$ real squares, one gets the $\lambda$ with the sign $-$,
where $2\lambda$ denotes the number of imaginary roots.

### Section 3

Now, let $n>5$. We assume that a point $P$ that lies on both $C_{2}=0$ and
$C_{3}=0$ has been determined by the procedure given above. The associated
coordinates are $a_{1}^{(0)},\dotsc,a_{n-1}^{(0)}$. We write

$\displaystyle\alpha_{i}-\alpha_{i}^{(0)}=\beta_{i},\quad i=1,\dotsc,n-2$

and reconstruct $C_{2}$ and $C_{3}$ in terms of
$\beta_{1},\dotsc,\beta_{n-2}$. 444One can assume that $a_{n-1}^{(0)}\not=0$,
after possibly changing the indices, and then write $a_{n-1}=1$ In this
manner, by summing the terms with same total degree in the $\beta_{i}$, we
get:

$\displaystyle C_{2}$ $\displaystyle=\phi_{1}+\phi_{2}$ (6) $\displaystyle
C_{3}$ $\displaystyle=\psi_{1}+\psi_{2}+\psi_{3}$

We want to solve the present problem in such a way that we determine a
straight line going through $P$ that lies on both the surfaces $C_{2}=0$ and
$C_{3}=0$. If this is successful, then the further condition $C_{4}=0$ only
requires the solution of an equation of the fourth degree. We denote the space
with homogeneous coordinates $a_{1},\dotsc,a_{n-1}$ as a $R_{n-2}$, in
accordance with the number that is its dimension.

The first conditions to be introduced are

$\phi_{1}=\psi_{1}=0.$ (7)

The space $R_{n-2}$ is then reduced to a $R_{n-4}$.555Translator’s footnote:
See Remark 2 for the translator’s summary of this argument. If we consider the
straight lines through $P$ as elements of the space, the point space $R_{n-4}$
has only $n-5$ dimensions. From this perspective, we refer to it as
$\mathcal{L}_{n-5}$ and consider the subvarieties $\phi_{2}=0$, $\psi_{2}=0$,
and $\psi_{3}=0$ inside it. If one takes any plane in this line space
$\mathcal{L}_{n-5}$, then $\phi_{2}=0$ and $\psi_{2}=0$ have four elements in
common - that is, four straight line generators. If $n>7$, then a common
straight line of $\phi_{2}=0$ and $\psi_{2}=0$ can be determined by solving a
fourth degree equation. We denote such a straight line by $\ell_{1}$.

If we can complete the task in such a way that we get a plane tangent to
$\phi_{2}=0$ and $\psi_{2}=0$, then this result will be solved, as this plane
intersects $\psi_{3}=0$ in three straight lines, so that everything else comes
down to the solution of a degree three equation. The equations of the
hyperplanes, which meet the subvarieties $\phi_{2}=0$ and $\psi_{2}=0$ along
the straight line $\ell_{1}$ are given by

$\phi_{1}^{(1)}=0,\quad\psi_{1}^{(1)}=0.$ (8)

We assume the relations (8) are satisfied. The $R_{n-4}$ discussed above then
reduces to an $R_{n-6}$. If we restrict to this $\mathbb{P}^{n-6}$, then
elements of $\phi_{2}=0$ and $\psi_{2}=0$ appear as planes through $\ell_{1}$.
However, there are exceptions for when $n=7$, as $R_{n-6}$ coincides with
$\ell_{1}$, and for when $n=8$, as $\phi_{2}=0$ and $\psi_{2}=0$ are reduced
to the double-counted straight line $\ell_{1}$. We now replace the point space
$R_{n-6}$ with a space whose elements are the planes that contain the line
$\ell_{1}$. This plane space obviously has dimension $n-8$ and is denoted by
$\mathcal{E}_{n-8}$. We now have the intersection of the quadratic [algebraic]
manifolds $\phi_{2}=0$ and $\psi_{2}=0$ in $R_{n-6}$ in the plane space
$\mathcal{E}_{n-8}$. If $n\geq 10$, one concludes that only a fourth degree
equation has to be solved to determine a common plane of $\phi_{2}=0$ and
$\psi_{2}=0$ in $\mathcal{E}_{n-8}$. Therefore, we have the theorem:

If $n\geq 10$, the general equation (1) can be reduced to the form

$y^{n}+C_{5}y^{n-5}+\cdots+C_{n}=0$ (9)

by a transformation (2) in such a way that determining the parameters $a_{i}$
requires only the solution of a finitely many quartic, cubic, and quadratic
irrationalities.

Although an equation of the form (9) can be achieved with only two pairs of
imaginary roots, at least three pairs of imaginary roots of (1) must exist for
this transformation to be completed using only real numbers. It is only under
this condition that the hypersurface $C_{2}=0$ contains real planes. However,
if one wishes for not only necessary, but also sufficient conditions here,
this does seem to be possible without fairly in-depth discussions.

### Section 4

The case $n=9$ was examined by D. Hilbert in a recently published work.666Über
die Gleichungen neuten Grades, Math. Ann. 97, S. 243 (1926) It is first
demonstrated how one can determine a three-dimensional space $R_{3}$ which is
completely contained in the hypersurface $C_{2}=0$ and then one considers the
degree three surface $F_{3}$ in this $R_{3}$ which is cut out by the
hypersurface $C_{3}=0$. As is known, this surface $F_{3}$ only depends on four
fundamental parameters. This is also particularly evident when one transforms
only the left term of the equation to a sum of five cubes, for which it is
necessary to solve a fifth degree equation. Hilbert now puts the general
equation of ninth degree in the form

$y^{9}+C_{5}y^{4}+C_{6}y^{3}+C_{7}y^{2}+C_{8}y+C_{9}=0$ (10)

by first determining one of the 27 straight lines on the surface $F_{3}$ and
then intersecting one of these straight lines with the hypersurface $C_{4}=0$.
Both the equation (10), where one can easily set $C_{9}=1$, as well as the
equation of degree 27 are functions of only four parameters. The result of
this is that the solution of the general equation of ninth degree only
requires algebraic functions of four arguments in such a way that “one can get
by with functions of one argument, sums, and two special functions of four
arguments”.

It can be shown that, if the general equation of ninth degree is reduced to
the form (10), then there is no need to solve an equation of degree higher
than five, so that one of the special functions of four arguments above is
unnecessary. At the end, we generalize our task set at (4) by using auxiliary
equations larger than $m$, but with the restriction that each degree will
still always be $<n$.

As in the previous case, we determine a point $P$ on both $C_{2}=0$ and
$C_{3}=0$ and still suppose the conditions (8) hold. Since $n=9$ here, the
left terms of the subcone $\phi_{2}=0$ and $\psi_{2}=0$ can be written in five
homogeneous coordinates, such as $z_{1},z_{2},z_{3},z_{4},z_{5}$. We now look
for the self-conjugate pentahedron common to both $\phi_{2}$ and $\psi_{2}$,
which corresponds to solving a fifth degree equation.777Translator’s footnote:
See Remark 3 for more exposition. If this pentahedron is assumed to be a
coordinate pentahedron, we have expressions for $\phi_{2}$ and $\psi_{2}$ of
the form:

$\displaystyle\phi_{2}$
$\displaystyle=a_{1}z_{1}^{2}+a_{2}z_{2}^{2}+a_{3}z_{3}^{2}+a_{4}z_{4}^{2}+a_{5}z_{5}^{2}$
(11) $\displaystyle\psi_{2}$
$\displaystyle=b_{1}z_{1}^{2}+b_{2}z_{2}^{2}+b_{3}z_{3}^{2}+b_{4}z_{4}^{2}+b_{5}z_{5}^{2}.$

One can now eliminate any of the five variables from $\phi_{2}=0$ and
$\psi_{2}=0$ and thus obtain five relations, one of which we write in the form

$c_{1}z_{1}^{2}+c_{2}z_{2}^{2}+c_{3}z_{3}^{2}+c_{4}z_{4}^{2}=0.$ (12)

We take (12) as the equation of a degree two hypersurface and try to determine
the corresponding straight line generators, which only requires taking square
roots. The straight lines in this generating set are assigned to the values of
a a parameter $\lambda$ and, likewise, the points of a specified generator are
assigned to the values of another parameter $t_{1}$ and the their coordinates
$z_{1},z_{2},z_{3},z_{4}$ can be expressed linearly in both $\lambda$ and $t$.
According to (11), we get the relation for $z_{5}$:

$z_{5}^{2}=a_{2}(\lambda)t^{2}+b_{2}(\lambda)t+c_{2}(\lambda).$ (13)

Since the elements of (11) are straight lines through $P$, it can be seen that
a generator of (12) corresponds to a two-dimensional cone whose apex is $P$.
Since this cone must have six generators in common with the hypersurface
$\psi_{3}=0$, it follows that one need not use auxiliary equations of degree
more than six when reducing the general equation of ninth degree to the form
(10).

### Section 5

This matter can be simplified by looking for a value of $\lambda$ such that
the two-dimensional cone splits into two planes. We need only solve the degree
four equation

$\left[b_{2}(\lambda)\right]^{2}-4a_{2}(\lambda)c_{2}(\lambda)=0.$ (14)

The generating system of (12) has four other conjugate pairs of planes, so
that the total number of planes common to the subcones $\phi_{2}=0$ and
$\psi_{2}=0$ is 16. According to the five different relations in four
variables of the form (12), each of these planes can be paired with five
others. Two planes that can be paired together intersect each other in a
straight line. If not, they only intersect at $P$.

The theory of the intersection of the two cones $\phi_{2}=0$ and $\psi_{2}=0$
is indeed well known.888Translator’s footnote: See Remark 4 for the
translator’s summary of this argument. They system of equations (11) is often
used to study properties of a degree four surface with a double conic section.
The 16 common planes of the two subcones correspond to the 16 straight lines
lying on such a surface.

We have now demonstrated that the general equation of ninth degree can also be
converted to the form (9) without having to presuppose the solution of
auxiliary equations which each depend on more than one parameter. Among the
auxiliary equations, however, there is one of fifth-degree, so we cannot get
by with square roots and cube roots, as we can for $n>9$. However, it seems to
be difficult to prove strictly that the latter is not possible at all for
$n=9$.

##  3 Translator’s Notes

######  Remark 1.

(Main Argument of Section 2)
Let $\mathbb{A}^{n}$ be the affine space of Tschirnhaus transformations and
$\mathbb{P}^{n-1}$ its projectivization. $\mathbb{V}(C_{1})$ is a hyperplane
and thus $\mathbb{V}(C_{1})\cong\mathbb{P}^{n-2}$. A rational point $P$ of
$\mathbb{V}(C_{1})\cap\mathbb{V}(C_{2})$ can be determined over a quadratic
extension of the base field. Then, the tangent hyperplane $T$ at $P$ can be
computed rationally and will have dimension at least 4 in
$\mathbb{V}(C_{1})=\mathbb{P}^{n-2}$ if $n>6$. Hence
$\mathbb{V}(C_{1})\cap\mathbb{V}(C_{2})\cap T$ is a quadric that is singular
at $P$ in $\mathbb{V}(C_{1})\cap T\cong\mathbb{P}^{n-3}$. Consequently, a line
in this cone can be determined by solving a quadratic polynomial and it
suffices to intersect this line with $F_{3}=0$.

######  Remark 2.

(Main Argument of Section 3)
To re-state Wiman’s approach, this $\mathbb{P}^{n-4}$ is obtained by
projectivizing and then considering
$\mathbb{V}(C_{1})\cap\mathbb{V}(\phi_{1})\cap\mathbb{V}(\psi_{1})$ inside
$\mathbb{P}^{n-1}$. By shifting $P$ to the origin (e.g. $[0:\cdots:0:1]$),
Wiman then uses a classical correspondence to consider the pencil of lines
through $P$ and identifies it as $\mathcal{L}_{n-5}\cong\mathbb{P}^{n-5}$.
Note that $\phi_{2},\psi_{2},$ and $\psi_{3}$ induce hypersurfaces of the same
degree in $\mathcal{L}_{n-5}$. If $n\geq 7$, a point $Q$ of
$\mathbb{V}(\phi_{2})\cap\mathbb{V}(\psi_{2})\subseteq\mathcal{L}_{n-5}$ can
be determined by solving a quartic equation. Moreover, by construction, the
line $\ell_{1}$ determined by $P$ and $Q$ lies in
$\mathbb{V}(C_{1})\cap\mathbb{V}(C_{2})$ in the ambient space.

Now, consider the tangent hyperplanes of
$\mathbb{V}(\phi_{2}),\mathbb{V}(\psi_{2})\subseteq\mathbb{P}^{n-4}$ defined
by the polynomials $\phi_{1}^{(1)}$ and $\psi_{1}^{(1)}$. Consider the
$\mathbb{P}^{n-6}$ given by
$\mathbb{V}(C_{1})\cap\mathbb{V}(\phi_{1})\cap\mathbb{V}(\psi_{1})\cap\mathbb{V}(\phi_{1}^{(1)})\cap\mathbb{V}(\psi_{1}^{(1)})$.
Every point not on $\ell_{1}$ in
$\mathcal{E}_{n-8}=\mathbb{V}(\phi_{2})\cap\mathbb{V}(\psi_{2})\subseteq\mathbb{P}^{n-6}$
determines a plane on $\mathbb{V}(\phi_{2})\cap\mathbb{V}(\psi_{2})$ in the
ambient space; this is possible when $n\geq 9$. Determining such a point
$Q^{\prime}$ in $E_{n-8}$ determines a line
$\ell_{2}\subseteq\mathbb{V}(\phi_{2})\cap\mathbb{V}(\psi_{2})\subseteq\mathcal{L}_{n-5}$
and thus a point $Q^{\prime}$ of
$\ell_{2}\cap\mathbb{V}(\psi_{3})\subseteq\mathcal{L}_{n-5}$ can be determined
by solving a cubic equation. However, the line determined by $P$ and
$Q^{\prime}$ lies inside
$\mathbb{V}(C_{1})\cap\mathbb{V}(C_{2})\cap\mathbb{V}(C_{3})$ in the ambient
space and thus a point of
$\mathbb{V}(C_{1})\cap\mathbb{V}(C_{2})\cap\mathbb{V}(C_{3})\cap\mathbb{V}(C_{4})$
can be determined by solving a quartic equation.

######  Remark 3.

(Pencils of Quadrics)
The forms $\phi_{2}$ and $\psi_{2}$ define a pencil of quadratic forms in the
five variables $z_{1},\dotsc,z_{5}$. The singular fibers of the pencil are
given by the roots of the discriminant, which is a polynomial of degree 5.
Hence, determining a singular quadric in the pencil corresponds to solving a
degree 5 polynomial. Wiman then again uses the observation that a singular
quadric is a cone.

######  Remark 4.

(Degree 4 del Pezzo Surfaces)
Here Wiman is observing the fact that the intersection of two quadrics in
$\mathbb{P}^{4}$ is a degree 4 del Pezzo surface.


Safra Campus, Jerusalem 91904, Israel

E-mail address<EMAIL_ADDRESS>


# Supercurrent Noise in a Phase-Biased Superconductor-Normal Ring in Thermal
Equilibrium

Ziwei Dou1, Xavier Ballu1, Quan Dong2, Yong Jin3, Richard Deblock1, Sandrine
Autier-Laurent1, Sophie Guéron1,
Hélène Bouchiat1 , and Meydi Ferrier1 1\. Université Paris-Saclay, CNRS,
Laboratoire de Physique des Solides, 91405, Orsay, France
2\. CryoHEMT Electronics for Scientific Research, 91400, Orsay, France
3\. Université Paris-Saclay, CNRS, Centre de Nanosciences et de
Nanotechnologies, 91120, Palaiseau, France

###### Abstract

In superconductor-normal-superconductor (SNS) junctions, dissipationless
supercurrent is mediated via Andreev bound states (ABSs) controlled by the
phase difference between the two superconductors. Theory has long predicted
significant fluctuations, and thus a noise, of such supercurrent in
equilibrium, due to the thermal excitation between the ABSs. Via the
fluctuation-dissipation theorem (FDT), this leads paradoxically to a
dissipative conductance even in the zero frequency limit. In this article, we
directly measure the supercurrent noise in a phase-biased SNS ring inductively
coupled to a superconducting resonator. Using the same setup, we also measure
its admittance and quantitatively demonstrate the FDT for any phase. The
dissipative conductance shows an 1/T temperature dependence, in contrast to
the Drude conductance of unproximitized metal. This is true even at phase
$\pi$ where the Andreev spectrum is gapless. Supported by linear response
theory, we attribute this to the enhanced current correlation between the ABSs
symmetrical across the Fermi level. Our results provide insight into the
origin of noise in hydrid superconducting devices which may be useful for
probing topologically protected ABSs.

††preprint: APS/123-QED

To what extent is a dissipationless system noisy? This paradoxical question is
of fundamental importance for the ongoing research which combines
superconducting and topological materials with dissipationless features. A
paradigmatic example is a Josephson junction, where two superconducting
electrodes are coupled by a tunnel barrier. A well-known property of such a
junction is the dissipationless supercurrent that flows with no voltage bias.
As a consequence, no fluctuations can be detected in the supercurrent at
frequency or temperature much smaller than the superconducting gap $\Delta$
[1]. This is intuitively consistent with the fluctuation-dissipation theorem
[2, 3] (FDT) stating that dissipation and current noise are proportional in
equilibrium. In superconductor-normal-superconductor (SNS) junctions, the
tunnel barrier is replaced by a normal metal. There, a supercurrent is carried
by the phase-coherent Andreev bound states (ABSs), whose spectrum varies with
the phase difference $\varphi$ between the two superconducting electrodes. In
a long diffusive junction, the Andreev spectrum displays an induced “minigap”
$E_{g}\ll\Delta$ which is maximal at phase 0 and closes at phase $\pi$ [4, 5].
Greater noise has been predicted in such junctions since there exists a regime
where $E_{g}<k_{B}T\ll\Delta$ and transitions of the ABSs across the minigap
occur more frequently [6, 7, 8].

The supercurrent fluctuations are characterized by their noise power spectrum
$S_{I}(f)$, which is the Fourier transform of the time correlation of current.
From the FDT, a dissipation is expected: it takes the form of a conductance
$G$ and obeys the relation $S_{I}=4k_{B}TG$, where $T$ is the temperature. The
theoretical prediction of finite fluctuations down to zero frequency [6, 7, 8]
thus entails a non-zero $G$, in apparent contradiction with the
dissipationless character of the supercurrent. The resolution of this paradox
calls for a rigorous definition of the conductance. Naively, $G$ is the ratio
of the current to voltage in the linear regime and requires a voltage biasing
scheme. For a low impedance system such as an SNS junction, such a biasing
scheme is challenging, particularly because even an infinitely small dc
voltage bias drives the junction out of equilibrium due to the ac Josephson
effect [7]. A more suitable approach is to adopt a phase-biasing scheme using
an SNS ring. There the phase can be controlled by the dc magnetic flux via
$\varphi=2\pi\Phi/\Phi_{0}$ ($\Phi_{0}=h/2e$). Meanwhile, the conductance can
be accessed by adding a small ac flux $\delta\Phi\ll\Phi_{0}$ at frequency
$f$, which is produced for example by a superconducting resonator inductively
coupled to the ring. The magnetic susceptibility $\chi=\delta i/\delta\Phi$
can thus be measured. Generally, $\chi$ is complex due to the finite
relaxation time of the ABSs in response to microwave or thermal excitation. It
thus results in an admittance $Y=-i\chi/2\pi f$ whose conductance is the real
part $G=Re(Y)$ [9, 10, 11].

In this work, we have directly measured the long predicted supercurrent noise
of a phase-biased SNS ring. We have also measured its conductance, and
quantitatively verified the FDT over the full range of phase between 0 and
$2\pi$. Moreover, the experiment identifies an $1/T$ temperature dependence of
the conductance in contrast to the constant Drude conductance. Surprisingly,
it is true even at $\varphi=\pi$ when the minigap closes. Supported by linear
response theory, we show that such behavior is a manifestation of the enhanced
current correlation between the ABSs symmetrical across the Fermi level.

## Measurement method

Figure 1: Noise spectrum of the measurement system: (a) Schematic of the
measurement setup. Noise or transmission measurement is performed using the
same setup by connecting the input to a grounded $Z_{0}$ (port 2) at 300 K or
to an rf source (port 1), respectively. The superconducting resonator (blue)
is inductively coupled to the SNS ring (optical image, scale bar: 5 µm). Solid
rectangle: the normal metal wire. Dashed rectangle: the coupling inductance
$L_{c}$. The resonator (at $\sim$ 10 mK) is directly connected to a low-
temperature amplifier at 4 K with the gain $g$. (b) Effective circuit modeling
of (a). The total device (the resonator and the SNS ring, dashed box) is
modeled as a parallel-GLC circuit. See text. (c) Measured voltage noise
spectrum $S_{V}(f)$ at $T$ = 200 mK (red), fitted to Eq. 2 (black dashed). The
thermal contribution $S_{V,th}$ is plotted as the black continuous line. (d)
$S_{V}(f)$ at $\varphi$ = 0 around $f_{r}$, showing temperature dependence.
(e) $S_{V}(f_{r})$ extracted from (d) versus the mixing chamber temperature
$T_{MC}$. The linear relation (dashed line) is extrapolated using the high
temperature data.

Fig. 1(a) displays the schematic of the measurement setup, which can be
configured to measure both the noise spectrum and the transmission
coefficient. For transmission measurement, the input port is connected to an
rf source with sufficiently small power and the transmitted signal at the
output port is down-converted to obtain the transmission coefficient
$\Gamma(f)$. For noise measurement, the input port is simply grounded via
$Z_{0}$ = 50 and the voltage noise spectrum $g^{2}S_{V}(f)$ at the output is
recorded, where $g$ is the gain of the amplifier. Each spectrum is averaged
and filtered according to the procedure detailed in [12]. The optical image
shows the SNS ring, formed by the normal metal weak link on the left (solid
rectangle) enclosed by the coupling inductance $L_{c}$ on the right (dashed
rectangle). The dc magnetic flux $\Phi$ through the loop is used to phase-bias
the junction via $\varphi=2\pi\Phi/\Phi_{0}$. The normal wire is made by
evaporating a Ti/Au of 5nm/100nm thickness on an undoped silicon substrate,
with a length $L=$ 1.5 µm and a width $W=$ 100 nm. The rest of the ring is
formed by 80 nm thick superconductor molybdenum-rhenium (MoRe). The SNS ring
is probed via the superconducting resonator. The resonator is coupled to the
input transmission line by a coupling capacitance $C_{c}$, and is directly
connected to a home-made HEMT amplifier at 4 K [13]. The small $C_{c}$ ($\sim$
1 pF) and the large amplifier impedance ($>$ 1 G) are chosen to preserve the
resonator’s quality factor. To reduce the input capacitance caused by the
Miller effect, $g$ is limited around 1 [14]. Therefore, instead of
amplification, the amplifier is used to match the impedance between the
resonator and the external circuit.

From Fig. 1(b), the paralleled $G_{tot}$, $L_{tot}$ and $C_{tot}$ models the
total system of the SNS ring plus the resonator. $S_{I,th}$ is the total
current noise originated from the system. We note that, in the parallel GLC
model, larger $G$ corresponds to greater dissipation. The contributions from
the resonator and the ring can be separated [12]:

$\displaystyle\begin{split}&\frac{1}{L_{tot}}(\Phi)=\frac{1}{L_{reso}}+\kappa\frac{1}{L_{ring}}(\Phi)\\\
&G_{tot}(\Phi)=G_{reso}+\kappa G_{ring}(\Phi)\\\
&S_{I,th}(\Phi)=4k_{B}TG_{reso}+\kappa S_{I,ring}(\Phi)\\\ \end{split}$ (1)

where the resonator parameters $1/L_{reso}$ and $G_{reso}$ are constant while
$1/L_{ring}$ and $G_{ring}$ are phase-dependent. $\kappa=(L_{c}/L_{reso})^{2}$
is the inductive coupling coefficient and is a small factor in the order of
10-5. For the resonator, we can directly write its thermal noise as
$4k_{B}TG_{reso}$ [15]. It is important that we do not assume such relation
between $S_{I,ring}(\Phi)$ and $G_{ring}(\Phi)$. Indeed, the central goal of
this article is to demonstrate that these two quantities, both due to the ABS
dynamics, are also linked by the FDT.

As will be shown later, $1/L_{ring}(\Phi)$ and $G_{ring}(\Phi)$ are obtained
by the phase-dependent resonance frequency and quality factor of $\Gamma(f)$
[11]. $S_{I,th}$, on the other hand, is extracted from the voltage noise
$g^{2}S_{V}(f)$ measured at the amplifier output. From Fig. 1(b), the relation
between $S_{I,th}$ and $S_{V}(f)$ is [12]:

$\displaystyle\begin{split}&S_{V}(f)=S_{V,th}(f)+S_{V,amp}(f)\\\
&S_{V,th}=\frac{S_{I,th}}{|Y_{tot}|^{2}},S_{V,amp}=\frac{S_{I,amp}(Y_{tot})}{|Y_{tot}|^{2}}\\\
\end{split}$ (2)

Figure 2: Phase-dependent supercurrent noise and demonstration of FDT for the
SNS ring: (a) $1/L_{tot}$ and $G_{tot}$ measured by transmission for all
temperatures. (b) Measured $S_{V}(f)$ at $\varphi$ = 0 (black) and $\pi$
(red). (c) Extracted $S_{I,th}(f)$ at $\varphi$ = 0 (black) and $\pi$ (red),
with the mean value marked by the dashed lines. (d,e) Fluctuation-dissipation
relation of the total device, comparing the mean values of the extracted noise
$S_{I,th}$ (squares) and $4k_{B}TG_{tot}$ (solid lines) using the transmission
data in (a). The phase-independent background $4k_{B}TG_{reso}$ is the thermal
noise from the resonator. The phase-dependent part $\kappa S_{I,ring}$ is due
only to the supercurrent noise of the SNS ring.

where $Y_{tot}=G_{tot}+1/(j2\pi fL_{tot})+j2\pi fC_{tot}$ is the total
admittance of the device. The noise of the amplifier is modeled as
$S_{I,amp}$, which depends on both $Y_{tot}$ and the phase-independent
internal parameters [12]. From Eqs. 1 and 2, to accurately extract $S_{I,th}$,
both the resonator parameters $1/L_{reso}$ and $G_{reso}$ and the amplifier
noise $S_{I,amp}$ should be calibrated. Since we are predominantly interested
in the phase dependent part of the noise, we re-calibrate those parameters at
phase 0 for each temperature which is then kept constant for the other phases.
This is done by assuming the FDT for the phase-independent elements of the
setup $S_{I,th}(0)=4k_{B}TG_{th}(0)$ and fitting the measured $S_{V}(f)$ to
Eqs. 2 [12]. Such calibration scheme is required to compensate any drift of
the amplifier and resonator parameters with temperature or even time, which is
indeed seen later in Fig. 2(a). Nevertheless, we check that the maximum error
made on $G_{reso}$ in the order of $\kappa G_{ring}$ does not affect the
validity of our phase-dependent results [12].

The superb quality of the fit obtained in Fig. 1(c) confirms the relevance of
our model. By subtracting the amplifier contribution $S_{V,amp}$, the thermal
noise $S_{V,th}$ contributes to roughly 20% of $S_{V}$. We also see clearly
the temperature dependence of the noise on resonance in Figs 1(d) and (e).
Assuming a linear dependence at the highest $T$ [dashed line Fig. 1(e)], the
data below 50 mK digress slightly from the such relation, indicating a
different electronic temperature from that of the mixing chamber. Using the
extrapolated linear relation, the electronic temperatures are calibrated in
[12] whose values are noted in Fig. 1(d) and the later figures.

## Experimental Validation of FDT

The main result of this work is the experimental demonstration of the FDT for
an SNS ring. It consists of measuring independently $S_{I,ring}$ from the
voltage noise spectrum and $G_{ring}$ from the transmission coefficient and
check the relation $S_{I,ring}(\Phi)=4k_{B}TG_{ring}(\Phi)$ at all phases.
First we measure the admittance of the SNS ring via transmission coefficient
according to the following relations [12, 11]:

$\displaystyle\begin{split}&\frac{\kappa}{L_{ring}}(\Phi)=\frac{2}{L_{reso}}\frac{\delta
f_{r}(\Phi)}{f_{r}}\\\ &\kappa G_{ring}(\Phi)=\frac{1}{2\pi
f_{r}L_{reso}}\delta\left(\frac{1}{Q}\right)(\Phi)\\\ \end{split}$ (3)

where $f_{r}$ and $Q$ are the resonance frequency and the quality factor
respectively. Following the method described in [12, 11], as $\Phi$ is swept,
we simultaneously record $\delta f_{r}$ and $\delta(1/Q)$ using a phase-locked
feedback loop which maintains the resonator on resonance [11]. Since the above
method only determine the phase-dependence of the ring admittance, we can set
$\kappa/L_{ring}(0)$ and $\kappa G_{ring}(0)$ as 0. Adding to $1/L_{tot}(0)$
and $G_{tot}(0)$ from calibration, we thus obtain $1/L_{tot}(\Phi)$ and
$G_{tot}(\Phi)$ which are plotted in Fig. 2(a). At low temperatures,
$1/L_{tot}$ is reduced and $G_{tot}$ is enhanced at phase $\pi$, whereas such
effect is weakened at high temperatures, consistent with earlier results [16,
17, 18].

For the noise measurement, at each temperature, $S_{V}(f)$ is taken at a
series of phases from 0 to $2\pi$. Fig. 2(b) shows two examples at phase 0 and
$\pi$ for $T$ = 60 mK. The voltage noise on resonance decreases significantly
at phase $\pi$, consistent with the enhanced $G_{tot}$. Using the calibrated
$S_{I,amp}$ and the transmission measurement in Fig. 2(a) at each temperature,
$S_{I,th}(f)$ is extracted according to Eq. 2. As plotted in Fig. 2(c),
$S_{I,th}(f)$ for a given phase is frequency-independent and shows the
different mean value marked by the dashed line. To further reduce the
uncertainty of such mean value, around 100 data points shown in the figure are
averaged to calculate $S_{I,th}$ at the given phase. Such mean value (squares)
with the associated standard error of mean (error bar) for all phases and
temperatures are summarized in Figs. 2(d,e). Meanwhile, $4k_{B}TG_{tot}$ using
the transmission data is plotted as the solid lines. For the whole data set,
quantitative agreement between $S_{I,th}$ and $4k_{B}TG_{tot}$ is achieved.
Its phase-dependent part is thus a direct demonstration of the FDT for the SNS
ring alone. The very high precision of the measurement can be appreciated by
the fact that we are able to resolve a current noise variation of 20 fA2/Hz as
well as a temperature difference of a few milli-Kelvin. The simultaneous
satisfaction of the FDT for both the resonator and the SNS ring using the same
temperature also indicates that the ABSs of the SNS ring are well thermalized
with the electronic temperature of the resonator, which is otherwise hard to
be confirmed experimentally.

## Understanding temperature dependence of dissipation

Figure 3: Screening effect of SNS ring and temperature dependence of junction
conductance: (a) $1/L_{ring}$ (dashed lines) and $1/L_{j}$ (continuous lines)
after removing the screening effect. (b) Similar comparison for $G_{ring}$
(dashed) and $G_{j}$ (continuous). Only $T$ = 50 mK and 100 mK are shown for
clarity. The data for other $T$ are shown in [12]. (c) Critical current
$I_{c}$ versus $T$ (squares) is fitted to the long diffusive junction model
[19] (dashed line), giving $E_{TH}/k_{B}=$ 30 mK. (d) Normalized $G_{j}$
versus $T$ (squares). $G_{j}(T)$ is fitted to $G_{N}(1+T^{*}/T)$ (dashed line)
with $G_{N}=$ 54 mS and $T^{*}=$ 30 mK.

Having established the equivalence between the supercurrent noise and the
conductance of the SNS ring, we focus on the conductance and try to understand
its implied mechanism. The admittance of the SNS ring $Y_{ring}$ is reproduced
without the resonator contribution in Figs. 3(a, b) as the dashed lines.
However, due to the screening effect of the loop supercurrent, $Y_{ring}$ is
in general not equal to $Y_{j}$ [11] of an isolate junction. In order to
facilitate comparison with theory developed for the junction, we thus convert
from $Y_{ring}$ to $Y_{j}$ by
$1/L_{ring}=(1/L_{j})/(1+L_{l}\chi_{j}^{\prime})$ and
$G_{ring}=G_{j}/(1+L_{l}\chi_{j}^{\prime})^{2}$, where
$\chi^{\prime}_{j}=-1/L_{j}$ is the real component of the junction magnetic
susceptibility and $L_{l}=250$ pH is the loop inductance [11, 12]. Therefore,
$1/L_{ring}$ (or $G_{ring}$) and $1/L_{j}$ (or $G_{j}$) are identical only
when the screening coefficient $\beta(T)=L_{l}\chi_{j}^{\prime}\ll 1$. In our
case, $\beta\sim$ 0.1 at low $T$, and the screening effect is not negligible.
The junction $1/L_{j}$ and $G_{j}$ converted from the ring $1/L_{ring}$ and
$G_{ring}$ are plotted in Figs. 3(a, b) as the solid lines. For the inverse
inductance, the correction at low temperature is sizeable while at high
temperature it is negligible due to the reduced supercurrent. By integrating
$1/L_{j}=\partial I/\partial\Phi$, the critical current $I_{c}$ is obtained
and plotted in Fig. 3(c). From $I_{c}(T)$, we estimate the Thouless energy
$E_{TH}$ to be 30 mK [19, 12], and the temperature in the experiment is thus
always higher than $E_{TH}$. For conductance, on the other hand, $G_{j}$ is
almost phase-independent for all temperatures. We emphasize that such phase-
independent $G_{j}$ cannot be confused with $G_{reso}$ since $G_{reso}$ alone
does not produce the observed phase dependence in $\delta(1/Q)$[12, 11].
Instead, it is the screening effect that enables us to determine the flux-
independent $G_{j}$ which is otherwise difficult to be separated from
$G_{reso}$ obtained in Fig. 2(a) [12, 11]. In Fig. 3(d), $G_{j}(T)$ shows
strong enhancement at low temperature. This behavior is drastically different
from the Drude conductance of an incoherent normal metal, which is
temperature-independent in the milli-Kelvin range. $G_{j}(T)$ can also be
fitted to the relation $G_{N}(1+T^{*}/T)$ according to other works [11, 10,
7], giving $G_{N}$ = 54 mS close to the Drude conductance of a similar control
junction and $T^{*}$ = 30 mK close to the Thouless energy $E_{TH}$ estimated
from $I_{c}(T)$.

Figure 4: Theory to explain $G_{j}(\Phi,T)$: (a) Density of state versus flux.
$\Delta$ is the superconducting gap of the electrode. $E_{g}$ is the minigap.
(b) Current correlation function $\mathcal{C}(\epsilon^{+},\epsilon^{-})$ at
$\Phi/\Phi_{0}=$ 0.5. (c) $\mathcal{C}(\epsilon^{+},\epsilon^{-}=0)$ (blue
solid) and DOS (red dashed) at $\Phi/\Phi_{0}=0.5$. The peak of $\mathcal{C}$
at $\epsilon^{+}=0$ is essential to $G_{j}\sim 1/T$ (see text). (d) The total
conductance $G_{j}$ (solid) and its two components $G_{D}$ (dashed) and
$G_{ND}$ (dot-dashed), normalized by $G_{c}=G(0.5\Phi_{0})$. $G_{D}$ is offset
by $G_{c}$. $h\gamma=E_{g}$ and $k_{B}T=4E_{g}$. The aspect ratio $L/W$ = 8;
the superconducting coherence length and the mean-free path $\xi,l_{e}\approx
L/10$. $\Delta/E_{g}\approx 10$. See [12].

To explain the weak phase dependence as well as the strong temperature
dependence of $G_{j}$, we recall that by the linear response theory
$G_{j}=G_{D}+G_{ND}$ can be written in two terms [11, 20, 21]. In [12] we
justify that the resonance frequency is much less than the inelastic
scattering rate $\gamma$ and $k_{B}T/h$. Therefore, taking the zero-frequency
limit, we have:

$\displaystyle\begin{split}G_{D}&=-\frac{1}{h\gamma}\sum_{n}|J_{nn}|^{2}\frac{\partial
f_{n}}{\partial E_{n}}\\\ G_{ND}&=-\sum_{n\neq
m}|J_{nm}|^{2}\frac{h\gamma}{(E_{n}-E_{m})^{2}+(h\gamma)^{2}}\frac{f_{n}-f_{m}}{E_{n}-E_{m}}\\\
\end{split}$ (4)

where $E_{n}(\Phi)$ is the n-th Andreev level, $f_{n}=f(E_{n})$ is the Fermi-
Dirac function and $\gamma$ is the inelastic scattering rate.
$J_{nm}=(ie\hbar/m^{*})\langle n|\nabla|m\rangle$ is the matrix element of the
current operator. The term $G_{D}$ only involves the diagonal element of $J$
while the term $G_{ND}$ involves all the non-diagonal elements. $E_{n}$ and
$J_{nm}$ for long-diffusive junction can be numerically calculated by
diagonalizing the Bogoliubov-de Gennes Hamiltonian of the junction [22, 12].
The density of state (DOS) versus $\Phi$ calculated from $E_{n}(\Phi)$ is
plotted in Fig. 4(a), showing the phase-dependent minigap $E_{g}$ expected for
a long-diffusive junction (dashed line) [23]. We first explain the weak phase
dependence of $G_{j}$. In the experiment where $k_{B}T\gg E_{g}$, many ABSs
transitions across $E_{g}$ occur due to thermal excitation. Therefore, it is
justified to assume that all $|J_{nm}|^{2}$ with significant values are
included in $G_{ND}$ and
$G_{j}\propto\sum(|J_{nn}|^{2}+|J_{nm}|^{2})=Tr(|J|^{2})$ which does not
depend on the Aharanov-Bohm phase [11]. Indeed, the phase dependences of
$G_{D}$ and $G_{ND}$ cancel each other, resulting in a weakly phase-dependent
$G_{j}$ reproduced by calculation in Fig. 4(d). We also note that $h\gamma$
cannot be much smaller than $E_{g}$ to have a phase-independent $G_{j}$ [12].
With $E_{g}/k_{B}$ in the order of 10 mK, $\gamma$ is in the order of 100 MHz,
reasonable for similar set-ups [11, 24].

To explain the strong temperature dependence of $G_{j}$, we first approximate
$G_{ND}$ in the continuous spectrum limit, since the level spacing of the ABSs
in our junction is the smallest energy scale [12]:

$\begin{split}&G_{ND}\propto\iint\,d\epsilon^{+}\,d\epsilon^{-}\,\mathcal{C}(\epsilon^{+},\epsilon^{-})\mathcal{F}(\epsilon^{+},\epsilon^{-})\\\
&\mathcal{F}(\epsilon^{+},\epsilon^{-})=\frac{h\gamma}{(\epsilon^{-})^{2}+(h\gamma)^{2}}\frac{f(\epsilon^{+}+\frac{\epsilon^{-}}{2})-f(\epsilon^{+}-\frac{\epsilon^{-}}{2})}{\epsilon^{-}}\\\
\end{split}$ (5)

where $\epsilon^{+}$ and $\epsilon^{-}$ are the average and the difference
between two energy levels, respectively, and
$\mathcal{C}(\epsilon^{+},\epsilon^{-})$ is the current correlation function
convoluting the DOS and the matrix elements of the current operator
$|J_{nm}|^{2}$ (for the full definition see [12]). It is illustrative to
briefly examine the case of an unproximitized normal metal. Here $G_{D}\equiv
0$ since $E_{n}$ is phase-independent and $J_{nn}=\partial
E_{n}/\partial\varphi\equiv 0$. For $G_{ND}$, since
$\mathcal{C}(\epsilon^{+},\epsilon^{-})$ is constant in energy [12], the
integration involving the Fermi-Dirac distribution in Eq. 5 thus gives a
temperature-independent conductance. The result can be linked to the Drude
conductance if $\mathcal{C}$ includes the elastic scattering [12]. The
situation is drastically different for the proximitized normal metal. From Eq.
5, the key condition for a temperature-dependent $G_{ND}$ is an energy-
dependent $\mathcal{C}(\epsilon^{+},\epsilon^{-})$. As discussed in previous
works [7, 9], the presence of the minigap at phase 0 provides such energy-
dependence of $\mathcal{C}$, resulting $G_{j}\sim 1/T$. In our case, since
$G_{j}$ is phase-independent, it is thus surprising to see that such $1/T$
dependence is also true even at phase $\pi$ where the minigap closes. To
explain this, the maps of $\mathcal{C}(\epsilon^{+},\epsilon^{-})$ at phase
$\pi$ is shown in Fig. 4(b), with line cuts along $\epsilon^{-}=0$ shown in
Fig. 4(c) (blue solid). Although the DOS is indeed almost constant (red
dashed), a sharp peak in $\mathcal{C}$ appears at $\epsilon^{+}=0$ for a wide
range of $\epsilon^{-}$, revealing a strong current correlation between the
ABSs symmetrical across the Fermi level, even for levels far in energy [22].
In [12], we show that such feature does not exist for a normal metal without
superconducting contacts to impose the electron-hole symmetry for energy
smaller than $\Delta$. Since the peak width is much smaller than $h\gamma$ and
$k_{B}T$, it can be approximated as $\delta(\epsilon^{+})$. Inserting in Eq. 5
and using $G_{D}=0$ at phase $\pi$, we thus arrive at
$G_{j}=G_{ND}\propto\int\,d\epsilon^{-}\,\mathcal{F}(0,\epsilon^{-})$. In the
experiment where $k_{B}T\gg h\gamma$, the expression can be further
approximated by $G_{j}\sim\partial
f/\partial\epsilon^{\prime}(\epsilon^{\prime}=0)\sim 1/T$ thus giving the
observed temperature dependence. This shows that the electron-hole symmetry
imposed by the superconducting contacts makes a proximitized normal metal
fundamentally different from an unproximitized one even for phase $\pi$ where
the minigap closes.

In conclusion, we have measured the supercurrent noise as well as the linear
admittance of a phase-biased SNS ring using a superconducting resonator. We
verify the FDT long-predicted for the SNS junction for all phases and
temperatures. In addition, the junction conductance $G_{j}$ is several times
larger than the Drude conductance and shows an $1/T$ dependence for all
phases, in particular at phase $\pi$. This highlights a key difference between
a proximitized and unproximitized metal, even when the minigap closes. Using
the linear response theory, we argue that $G_{j}\sim 1/T$ is a manifestation
of the enhanced current correlation between the ABSs symmetrical across the
Fermi level due to the superconducting proximity effect. Our results provides
insight to the nature of the noise in the hybrid superconducting systems, many
of which are promising candidates for ultrasensitive quantum devices. They may
also be applied to the more exotic systems such as the topological SNS
junction, whose supercurrent noise can reveal the evidence of a topologically
protected crossing between the Andreev levels [25].

## References

  * Rogovin and Scalapino [1974] D. Rogovin and D. Scalapino, Fluctuation phenomena in tunnel junctions, Annals of Physics 86, 1 (1974).
  * Callen and Welton [1951] H. B. Callen and T. A. Welton, Irreversibility and generalized noise, Phys. Rev. 83, 34 (1951).
  * Kubo [1966] R. Kubo, The fluctuation-dissipation theorem, Rep. Prog. Phys. 29, 255 (1966).
  * Ivanov  _et al._ [2002] D. A. Ivanov, R. von Roten, and G. Blatter, Minigap in a long disordered sns junction: Analytical results, Phys. Rev. B 66, 052507 (2002).
  * le Sueur  _et al._ [2008] H. le Sueur, P. Joyez, H. Pothier, C. Urbina, and D. Esteve, Phase controlled superconducting proximity effect probed by tunneling spectroscopy, Phys. Rev. Lett. 100, 197002 (2008).
  * Averin and Imam [1996] D. Averin and H. T. Imam, Supercurrent noise in quantum point contacts, Phys. Rev. Lett. 76, 3814 (1996).
  * Zhou and Spivak [1997] F. Zhou and B. Spivak, Resistance of superconductor-normal-metal-superconductor (sns) junctions, Journal of Experimental and Theoretical Physics Letters 65, 369 (1997).
  * Martín-Rodero  _et al._ [1996] A. Martín-Rodero, A. L. Yeyati, and F. J. García-Vidal, Thermal noise in superconducting quantum point contacts, Phys. Rev. B 53, R8891 (1996).
  * Virtanen  _et al._ [2011] P. Virtanen, F. S. Bergeret, J. C. Cuevas, and T. T. Heikkilä, Linear ac response of diffusive sns junctions, Phys. Rev. B 83, 144514 (2011).
  * Tikhonov and Feigel’man [2015] K. S. Tikhonov and M. V. Feigel’man, Admittance of a long diffusive sns junction, Phys. Rev. B 91, 054519 (2015).
  * Dassonneville  _et al._ [2018] B. Dassonneville, A. Murani, M. Ferrier, S. Guéron, and H. Bouchiat, Coherence-enhanced phase-dependent dissipation in long sns josephson junctions: Revealing andreev bound state dynamics, Phys. Rev. B 97, 184505 (2018).
  * [12] Supplementary Materials.
  * Jin  _et al._ [2016] Y. Jin, Q. Dong, U. Gennser, L. Couraud, A. Cavanna, and C. Ulysse, Ultra-low noise cryohemts for cryogenic high-impedance readout electronics: Results and applications, 2016 13th IEEE International Conference on Solid-State and Integrated Circuit Technology (ICSICT) , 342 (2016).
  * Pozar [2011] D. Pozar,  _Microwave Engineering, 4th Edition_ (Wiley, 2011).
  * Bonaldi  _et al._ [1999] M. Bonaldi, P. Falferi, M. Cerdonio, A. Vinante, R. Dolesi, and S. Vitale, Thermal noise in a high q cryogenic resonator, Review of Scientific Instruments 70, 1851 (1999).
  * Dou  _et al._ [2021] Z. Dou, T. Wakamura, P. Virtanen, N.-J. Wu, R. Deblock, S. Autier-Laurent, K. Watanabe, T. Taniguchi, S. Guéron, H. Bouchiat, and M. Ferrier, Microwave photoassisted dissipation and supercurrent of a phase-biased graphene-superconductor ring, Phys. Rev. Research 3, L032009 (2021).
  * Haller  _et al._ [2022] R. Haller, G. Fülöp, D. Indolese, J. Ridderbos, R. Kraft, L. Y. Cheung, J. H. Ungerer, K. Watanabe, T. Taniguchi, D. Beckmann, R. Danneau, P. Virtanen, and C. Schönenberger, Phase-dependent microwave response of a graphene josephson junction, Phys. Rev. Research 4, 013198 (2022).
  * Lake  _et al._ [2017] R. E. Lake, J. Govenius, R. Kokkoniemi, K. Y. Tan, M. Partanen, P. Virtanen, and M. Möttönen, Microwave admittance of gold-palladium nanowires with proximity-induced superconductivity, Advanced Electronic Materials 3, 1600227 (2017).
  * Dubos  _et al._ [2001] P. Dubos, H. Courtois, B. Pannetier, F. K. Wilhelm, A. D. Zaikin, and G. Schön, Josephson critical current in a long mesoscopic s-n-s junction, Phys. Rev. B 63, 064502 (2001).
  * Trivedi and Browne [1988] N. Trivedi and D. A. Browne, Mesoscopic ring in a magnetic field: Reactive and dissipative response, Phys. Rev. B 38, 9581 (1988).
  * Büttiker [1996] M. Büttiker, Dynamic conductance and quantum noise in mesoscopic conductors, Journal of Mathematical Physics 37, 4793 (1996).
  * Ferrier  _et al._ [2013] M. Ferrier, B. Dassonneville, S. Guéron, and H. Bouchiat, Phase-dependent andreev spectrum in a diffusive sns junction: Static and dynamic current response, Phys. Rev. B 88, 174505 (2013).
  * Zhou  _et al._ [1998] F. Zhou, P. Charlat, B. Spivak, and B. Pannetier, Density of states in superconductor-normal metal-superconductor junctions, Journal of Low Temperature Physics 110, 841 (1998).
  * Murani  _et al._ [2019] A. Murani, B. Dassonneville, A. Kasumov, J. Basset, M. Ferrier, R. Deblock, S. Guéron, and H. Bouchiat, Microwave signature of topological andreev level crossings in a bismuth-based josephson junction, Phys. Rev. Lett. 122, 076802 (2019).
  * Fu and Kane [2009] L. Fu and C. L. Kane, Josephson current and noise at a superconductor/quantum-spin-hall-insulator/superconductor junction, Phys. Rev. B 79, 161408 (2009).


# A unified approach to study the existence and numerical solution of
functional differential equation

Dang Quang A${}^{\text{a}}$, Dang Quang Long${}^{\text{b}}$
${}^{\text{a}}$ Center for Informatics and Computing, VAST
18 Hoang Quoc Viet, Cau Giay, Hanoi, Vietnam
Email<EMAIL_ADDRESS>
${}^{\text{b}}$ Institute of Information Technology, VAST,
18 Hoang Quoc Viet, Cau Giay, Hanoi, Vietnam
Email<EMAIL_ADDRESS>
###### Abstract

In this paper we consider a class of boundary value problems for third order
nonlinear functional differential equation. By the reduction of the problem to
operator equation we establish the existence and uniqueness of solution and
construct a numerical method for solving it. We prove that the method is of
second order accuracy and obtain an estimate for total error. Some examples
demonstrate the validity of the obtained theoretical results and the
efficiency of the numerical method. The approach used for the third order
nonlinear functional differential equation can be applied to functional
differential equations of any orders.

Keywords:  Third order boundary value problem; Functional differential
equation; Existence and uniqueness of solution; Iterative method; Total error.
AMS Subject Classification: 34B15, 65L10

##  1 Introduction

Functional differential equations have numerous applications in engineering
and sciences [6]. Therefore, for the last decades they have been studied by
many authors. There are many works concerning the numerical solution of both
initial and boundary value problems for them. The methods used are diverse
including collocation method [10], iterative methods [1, 8], neural networks
[7, 9], and so on. Below we mention some results of typical works.

First it is worthy to mention the work of Reutskiy in 2015 [10]. In this work
the author considered the linear pantograph functional differential equation
with proportional delay

$\displaystyle
u^{(n)}=\sum_{j=0}^{J}\sum_{k=0}^{n-1}p^{jk}(x)u^{(k)}(\alpha_{j}x)+f(x),\quad
x\in[0,T]$

associated with initial or boundary conditions. Here $\alpha_{j}$ are
constants ($0<\alpha_{j}<1$). The author proposed a method, where the initial
equation is replaced by an approximate equation which has an exact analytic
solution with a set of free parameters. These free parameters are determined
by the use of the collocation procedure. Many examples show the efficiency of
the method but no errors estimates are obtained.

In 2016 Bica et al. [1] considered the boundary value problem (BVP)

$\begin{split}x^{(2p)}(t)=f(t,x(t),x(\varphi(t))),\quad t\in[a,b],\\\
x^{(i)}(a)=a_{i},\;x^{(i)}(b)=b_{i},\quad i=\overline{0,p-1}\end{split}$ (1)

where $\varphi:[a,b]\rightarrow\mathbb{R},\;a\leq\varphi(t)\leq b,\forall
t\in[a,b]$. For solving the problem, the authors constructed successive
approximations for the equivalent integral equation with the use of cubic
spline interpolation at each iterative step. The error estimate was obtained
for the approximate solution under the very strong conditions including
$(\alpha+13\beta)(b-a)M_{G}<1$, where $\alpha$ and $\beta$ are the Lipshitz
coefficients of the function $f(s,u,v)$ in the variables $u$ and $v$,
respectively; $M_{G}$ is a number such that $|G(t,s)|\leq M_{G}\;\forall
t,s\in[a,b]$, $G(t,s)$ being the Green function for the above problem. Some
numerical experiments demonstrate the convergence of the proposed iterative
method. But is a regret that in the proof of the error estimate for fourth
order nonlinear BVP there is a vital mistake when the authors by default
considered that the partial derivatives $\frac{\partial^{3}G}{\partial
s^{3}},\frac{\partial^{4}G}{\partial s^{4}}$ are continuous in
$[a,b]\times[a,b]$. But it is invalid because $\frac{\partial^{3}G}{\partial
s^{3}}$ has discontinuity on the line $s=t$. Due to this mistake the authors
obtained that the error of the method for fourth order BVP is $O(h^{4})$.
Although in [1] the method is constructed for the general function
$\varphi(t)$ but in all numerical examples only the particular case
$\varphi(t)=\alpha t$ was considered and the conditions of convergence were
not verified.

Recently, in 2018 Khuri and Sayfy [8] proposed a Green function based
iterative method for functional differential equations of arbitrary orders.
But the scope of application of the method is very limited due to the
difficulty in calculation of integrals at each iteration.

For solving functional differential equations, besides analytical and
numerical methods recently computational intelligence algorithms also are used
(see, e.g.,[7, 9]), where feed-forward artificial neural networks of different
architecture are applied. These algorithms are heuristic, so no errors
estimates are obtained and they require large computational efforts.

In this paper we propose a new approach to functional differential equations
(FDE). Although this approach can be applied to functional differential
equations of any orders with nonlinear terms containing derivatives but for
simplicity we consider the FDE of the form

$u^{\prime\prime\prime}=f(t,u(t),u(\varphi(t))),\quad t\in[0,a]$ (2)

associated with the general boundary conditions

$\begin{split}B_{1}[u]=\alpha_{1}u(0)+\beta_{1}u^{\prime}(0)+\gamma_{1}u^{\prime\prime}(0)=b_{1},\\\
B_{2}[u]=\alpha_{2}u(0)+\beta_{2}u^{\prime}(0)+\gamma_{2}u^{\prime\prime}(0)=b_{2},\\\
B_{3}[u]=\alpha_{3}u(1)+\beta_{3}u^{\prime}(1)+\gamma_{3}u^{\prime\prime}(1)=b_{3},\\\
\end{split}$ (3)

or

$\begin{split}B_{1}[u]=\alpha_{1}u(0)+\beta_{1}u^{\prime}(0)+\gamma_{1}u^{\prime\prime}(0)=b_{1},\\\
B_{2}[u]=\alpha_{2}u(1)+\beta_{2}u^{\prime}(1)+\gamma_{2}u^{\prime\prime}(1)=b_{2},\\\
B_{3}[u]=\alpha_{3}u(1)+\beta_{3}u^{\prime}(1)+\gamma_{3}u^{\prime\prime}(1)=b_{3},\\\
\end{split}$ (4)

such that

$Rank\begin{pmatrix}\alpha_{1}&\beta_{1}&\gamma_{1}&0&0&0\\\
\alpha_{2}&\beta_{2}&\gamma_{2}&0&0&0\\\
0&0&0&\alpha_{3}&\beta_{3}&\gamma_{3}\\\ \end{pmatrix}=3.$

As in (1), the function $\varphi(t)$ is assumed to be continuous and maps
$[0,a]$ into itself.

Developing the unified approach for fully third order nonlinear differential
equation

$\displaystyle
u^{\prime\prime\prime}=f(t,u(t),u^{\prime}(t),u^{\prime\prime}(t)$

in the previous works [2, 3], in this paper we establish the existence and
uniqueness of solution of the problem (2)-(3) and propose an iterative method
for finding the solution at both continuous and discrete levels. Some examples
demonstrate the validity of obtained theoretical results and the efficiency of
the proposed numerical method.

##  2 Existence and uniqueness of solution

Following the approach in [2, 3] (see also [4, 5]) to investigate the problem
(2)-(3) we introduce the nonlinear operator $A$ defined in the space of
continuous functions $C[0,a]$ by the formula:

$(A\psi)(t)=f(t,u(t),u(\varphi(t))),$ (5)

where $u(t)$ is the solution of the problem

$\begin{split}u^{\prime\prime\prime}(t)&=\psi(t),\quad 0<t<1\\\
B_{1}[u]&=b_{1},B_{2}[u]=b_{2},B_{3}[u]=b_{3},\end{split}$ (6)

where $B_{1}[u],B_{2}[u],B_{3}[u]$ are defined by (3). It is easy to verify
the following

###### Proposition 2.1

If the function $\psi$ is a fixed point of the operator $A$, i.e., $\psi$ is
the solution of the operator equation

$A\psi=\psi,$ (7)

where $A$ is defined by (5)-(6) then the function $u(t)$ determined from the
BVP (6) is a solution of the BVP (2)-(3). Conversely, if the function $u(x)$
is the solution of the BVP (2)-(3) then the function

$\psi(t)=f(t,u(t),u(\varphi(t)))$

satisfies the operator equation (7).

Now, let $G(t,s)$ be the Green function of the problem (6). Then the solution
of the problem can be represented in the form

$u(t)=g(t)+\int_{0}^{a}G(t,s)\psi(s)ds,$ (8)

where $g(t)$ is the polynomial of second degree satisfying the boundary
conditions

$B_{1}[g]=b_{1},B_{2}[g]=b_{2},B_{3}[g]=b_{3},$ (9)

Denote

$M_{0}=\max_{0\leq t\leq a}\int_{0}^{1}|G(t,s)|ds.$ (10)

For any positive number $M$ define the domain

$\mathcal{D}_{M}=\Big{\\{}(t,u,v)\mid 0\leq t\leq
a;|u|\leq\|g\|+M_{0}M;|v|\leq\|g\|+M_{0}M\Big{\\}},$ (11)

where $\|g\|=\max_{0\leq t\leq a}|g(t)|$.
As usual, we denote by $B[0,M]$ the closed ball of the radius $M$ centered at
$0$ in the space of continuous functions $C[0,a]$.

###### Theorem 2.2

Assume that:

(i) The function $\varphi(t)$ is a continuous map from $[0,a]$ to $[0,a]$.

(ii) The function $f(t,u,v)$ is continuous and bounded by $M$ in the domain
$\mathcal{D_{M}}$, i.e.,

$|f(t,u,v)|\leq M\quad\forall(t,u,v)\in\mathcal{D_{M}}.$ (12)

(iii) The function $f(t,u,v)$ satisfies the Lipshitz conditions in the
variables $u,v$ with the coefficients $L_{1},L_{2}\geq 0$ in
$\mathcal{D}_{M}$, i.e.,

$\begin{split}|f(t,u_{2},v_{2})-f(t,u_{1},v_{1})|\leq
L_{1}|u_{2}-u_{1}|+L_{2}|v_{2}-v_{1}|\quad\\\
\forall(t,u_{i},v_{i})\in\mathcal{D}_{M}\;(i=1,2)\end{split}$ (13)

(iv)

$q:=(L_{1}+L_{2})M_{0}<1.$ (14)

The the problem (2)-(3) has a unique solution $u(t)\in C^{3}[0,a]$, satisfying
the estimate

$|u(t)|\leq\|g\|+M_{0}M\quad\forall t\in[0,a].$ (15)

Proof. The proof of the theorem will be done in the following steps:
First we show that the operator $A$ is a mapping $B[0,M]\rightarrow B[0,M]$.
Indeed, for any $\psi\in B[0,M]$ we have $\|\psi\|\leq M$. Let $u(t)$ be the
solution of the problem (6). From (8) it follows

$|u(t)|\leq\|g\|+M_{0}M\quad\forall t\in[0,a].$ (16)

Since $0\leq\varphi(t)\leq a$ we also have

$|u(\varphi(t)|\leq\|g\|+M_{0}M\quad\forall t\in[0,a].$

Therefore, if $t\in[0,a]$ then $(t,u(t),u(\varphi(t)))\in\mathcal{D}_{M}$. By
the assumption (12) we have $|f(t,u(t),u(\varphi(t)))|\leq M\quad\forall
t\in[0,a]$. In view of (5) we have $|(A\psi)(t)|\leq M\quad\forall t\in[0,a]$.
It means $\|A\psi\|\leq M$ or $A\psi\in B[0,M]$.
Next, we prove that $A$ is a contraction in $B[0,M]$. Let
$\psi_{1},\psi_{2}\in B[0,M]$ and $u_{1}(t),u_{2}(t)$ be the solutions of the
problem (6), respectively. Then from the assumption (13) we obtain

$|A\psi_{2}-A\psi_{1}|\leq
L_{1}|u_{2}(t)-u_{1}(t)|+L_{2}|u_{2}(\varphi(t))-u_{1}(\varphi(t))|.$ (17)

From the representations

$u_{i}(t)=g(t)+\int_{0}^{a}G(t,s)\psi_{i}(s)ds,\quad(i=1,2)$

and (10) it is easy to obtain

$\displaystyle|u_{2}(t)-u_{1}(t)|$ $\displaystyle\leq
M_{0}\|\psi_{2}-\psi_{1}\|,$
$\displaystyle|u_{2}(\varphi(t))-u_{1}(\varphi(t))|$ $\displaystyle\leq
M_{0}\|\psi_{2}-\psi_{1}\|$

Combining the above estimates and (17), in view of the assumption (14) we
obtain

$\displaystyle\|A\psi_{2}-A\psi_{1}\|\leq q\|\psi_{2}-\psi_{1}\|,\quad q<1.$

Thus, $A$ is a contraction mapping in $B[0,M]$.

Therefore, the operator equation (7) has a unique solution $\psi\in B[0,M]$.
By Proposition 2.1 the solution of the problem (6) for this right-hand side
$\psi(t)$ is the solution of the original problem (2)-(3).

##  3 Solution method and its convergence

Consider the following iterative method:

  1. 1.

Given $\psi_{0}\in B[0,M]$, for example,

$\psi_{0}(t)=f(t,0,0).$ (18)

  2. 2.

Knowing $\psi_{k}(t)$ $(k=0,1,...)$ compute

$\begin{split}u_{k}(t)&=g(t)+\int_{0}^{a}G(t,s)\psi_{k}(s)ds,\\\
v_{k}(t)&=g(\varphi(t))+\int_{0}^{1}G(\varphi(t),s)\psi_{k}(s)ds,\end{split}$
(19)

  3. 3.

Update

$\psi_{k+1}(t)=f(t,u_{k}(t),v_{k}(t)).$ (20)

Set

$p_{k}=\dfrac{q^{k}}{1-q},\;d=\|\psi_{1}-\psi_{0}\|.$ (21)

######  Theorem 3.1 (Convergence)

Under the assumptions of Theorem 2.2 the above iterative method converges and
there holds the estimate

$\|u_{k}-u\|\leq M_{0}p_{k}d,$

where $u$ is the exact solution of the problem (2)-(3) and $M_{0}$ is given by
(10).

This theorem follows straightforward from the convergence of the successive
approximation method for finding fixed point of the operator $A$ and the
representations (8) and the first equation in (19).

To numerically realize the above iterative method we construct the
corresponding discrete iterative method. For this purpose cover the interval
$[0,a]$ by the uniform grid
$\bar{\omega}_{h}=\\{t_{i}=ih,\;h=a/N,i=0,1,...,N\\}$ and denote by
$\Phi_{k}(t),U_{k}(t),V_{k}(t)$ the grid functions, which are defined on the
grid $\bar{\omega}_{h}$ and approximate the functions
$\psi_{k}(t),u_{k}(t),v_{k}(t)$ on this grid, respectively.

Below we describe the discrete iterative method:

  1. 1.

Given

$\Psi_{0}(t_{i})=f(t_{i},0,0),\ i=0,...,N.$ (22)

  2. 2.

Knowing $\Psi_{k}(t_{i}),\;k=0,1,...;\;i=0,...,N,$ compute approximately the
definite integrals (19) by the trapezoidal rule

$\begin{split}U_{k}(t_{i})&=g(t_{i})+\sum_{j=0}^{N}h\rho_{j}G(t_{i},t_{j})\Psi_{k}(t_{j}),\\\
V_{k}(t_{i})&=g(\xi_{i})+\sum_{j=0}^{N}h\rho_{j}G(\xi_{i},t_{j})\Psi_{k}(t_{j}),i=0,...,N,\end{split}$
(23)

where $\rho_{j}$ are the weights

$\rho_{j}=\begin{cases}1/2,\;j=0,N\\\ 1,\;j=1,2,...,N-1\end{cases}$

and $\xi_{i}=\varphi(t_{i})$.

  3. 3.

Update

$\Psi_{k+1}(t_{i})=f(t_{i},U_{k}(t_{i}),V_{k}(t_{i})).$ (24)

Now study the convergence of the above discrete iterative method. For this
purpose we need some auxiliary results.

###### Proposition 3.2

If the function $f(t,u,v)$ has all partial derivatives continuous up to second
order and the function $\varphi(t)$ also has continuous derivatives up to
second order then the functions $\psi_{k}(t),u_{k}(t),v_{k}(t)$ constructed by
the iterative method (18)-(20) also have continuous derivatives up to second
order.

This proposition is obvious.

###### Proposition 3.3

For any function $\psi(t)\in C^{2}[0,a]$ there hold the estimates

$\displaystyle\int_{0}^{a}G(t_{i},s)\psi(s)ds=\sum_{j=0}^{N}h\rho_{j}G(t_{i},s_{j})\psi(s_{j})+O(h^{2}),$
(25)
$\displaystyle\int_{0}^{a}G(\xi_{i},s)\psi(s)ds=\sum_{j=0}^{N}h\rho_{j}G(\xi_{i},s_{j})\psi(s_{j})+O(h^{2}),$
(26)

where in order to avoid possible confusion we denote $s_{j}=t_{j}$.

Proof. The validity of (25) is guaranteed by [3, Proposition 3]. Here we
notice that (25) is not automatically deduced from the estimate for the
composite trapezoidal rule because the function
$\frac{\partial^{2}G(t_{i},s)}{\partial s^{2}}$ has discontinuity at
$s=t_{i}$.

Now we prove the estimate (26). Since $0\leq\xi_{i}=\varphi(t_{i})\leq a$,
there are two cases.
Case 1: $\xi_{i}$ coincides with one node $s_{j}$ of the grid
$\bar{\omega_{h}}$, i.e., there exists $s_{j}\in\bar{\omega_{h}}$ such that
$\xi_{i}=s_{j}$. Because the Green function $G(t,s)$ as a function of $s$, it
is continuous function at $s=\xi_{i}$ and is a polynomial of $s$ in the
intervals $[0,\xi_{i}$ and $[\xi_{i},a]$, we have

$\displaystyle\int_{0}^{a}G(\xi_{i},s)\psi(s)ds=\int_{0}^{\xi_{i}}G(\xi_{i},s)\psi(s)ds+\int_{\xi_{i}}^{a}G(\xi_{i},s)\psi(s)ds$
$\displaystyle=h\Big{(}\frac{1}{2}G(\xi_{i},s_{0})\psi(s_{0})+\sum_{m=1}^{j-1}G(\xi_{i},s_{m})\psi(s_{m})+\frac{1}{2}G(\xi_{i},s_{j})\psi(s_{j})\Big{)}+O(h^{2})$
$\displaystyle+h\Big{(}\frac{1}{2}G(\xi_{i},s_{j})\psi(s_{j})+\sum_{m=j+1}^{N-1}G(\xi_{i},s_{m})\psi(s_{m})+\frac{1}{2}G(\xi_{i},s_{N})\psi(s_{N})\Big{)}+O(h^{2})$
$\displaystyle=\sum_{j=0}^{N}h\rho_{j}G(t_{i},s_{j})\psi(s_{j})+O(h^{2}).$

Thus, (26) is proved for Case 1.
Case 2: $\xi_{i}$ lies between $s_{l}$ and $s_{l+1}$, i.e.,
$s_{l}<\xi_{i}<s_{l+1}$ for some $l=\overline{0,N-1}$. In this case, we
represent

$\int_{0}^{a}G(\xi_{i},s)\psi(s)ds=\int_{0}^{s_{l}}F(s)ds+\int_{s_{l}}^{\xi_{i}}F(s)ds+\int_{\xi_{i}}^{s_{l+1}}F(s)ds+\int_{s_{l+1}}^{a}F(s)ds.$
(27)

Here, for short we denote $F(s)=G(\xi_{i},s)\psi(s)$. Note that $F(s)\in
C^{2}$ in $[s_{l},\xi_{i}]$ and $[\xi_{i},s_{l+1}]$. Applying the composite
trapezoidal rule to the first and the last integrals in the right-hand side of
(27) we obtain

$\displaystyle T_{1}:$
$\displaystyle=\int_{0}^{s_{l}}F(s)ds+\int_{s_{l+1}}^{a}F(s)ds$ (28)
$\displaystyle=\sum_{j=0}^{l}\rho_{j}^{(l-)}F(s_{j})+\sum_{j=l+1}^{N}\rho_{j}^{(l+)}F(s_{j})+O(h^{2}),$

where

$\displaystyle\rho_{j}^{(l-)}=\left\\{\begin{array}[]{ll}\frac{1}{2},\quad
j=0,l\\\
1,\quad,1<j<l\end{array}\right.,\quad\rho_{j}^{(l+)}=\left\\{\begin{array}[]{ll}\frac{1}{2},\quad
j=l+1,N\\\ 1,\quad l+1<j<N\end{array}\right.$

For calculating the second and the third integrals in the right-hand side of
(27) we use the trapezoidal rule

$\displaystyle T_{2}:$
$\displaystyle=\int_{s_{l}}^{\xi_{i}}F(s)ds+\int_{\xi_{i}}^{s_{l+1}}F(s)ds$
(29)
$\displaystyle=\frac{1}{2}\big{[}(F(s_{l})+F(\xi_{i}))(\xi_{i}-s_{l})+(F(\xi_{i})+F(s_{l+1}))(s_{l+1}-\xi_{i})\big{]}+O(h^{2}).$

Using the points $s_{l}$ and $s_{l+1}$ for linearly interpolating $F(s)$ in
the point $\xi_{i}$ we have

$\displaystyle
F(\xi_{i})=F(s_{l})\frac{\xi_{i}-s_{l+1}}{s_{l}-s_{l+1}}+F(s_{l+1})\frac{\xi_{i}-s_{l}}{s_{l+1}-s_{l}}+O(h^{2}).$

From here we obtain

$F(\xi_{i})(s_{l+1}-s_{l})=F(s_{l})(s_{l+1}-\xi_{i})+F(s_{l+1})(\xi_{i}-s_{l})+O(h^{3}).$
(30)

Now, transforming $T_{2}$ we have

$\displaystyle
T_{2}=\frac{1}{2}\big{[}F(s_{l})(\xi_{i}-s_{l})+F(s_{l+1})(s_{l+1}-\xi_{i})\big{]}+F(\xi_{i})(s_{l+1}-s_{l})+O(h^{2})$

Further, in view of (30) it is easy to obtain

$T_{2}=\frac{1}{2}h(F(s_{l})+F(s_{l+1}))+O(h^{3}).$

Taking into account the above estimate, (28) and (27) we have

$\int_{0}^{a}G(\xi_{i},s)\psi(s)ds=\sum_{j=0}^{N}h\rho_{j}G(\xi_{i},s_{j})\psi(s_{j})+O(h^{2}).$

Thus, (26) is proved for Case 2 and the proof of Proposition 3.3 is complete.

###### Remark 3.4

If in Proposition 3.3 replace $G(t_{i},s)$ and $G(\xi_{i},s)$ by
$|G(t_{i},s)|$ and $|G(\xi_{i},s)|$, respectively then we obtain the analogous
estimates

$\displaystyle\int_{0}^{a}|G(t_{i},s)|\psi(s)ds=\sum_{j=0}^{N}h\rho_{j}|G(t_{i},s_{j})|\psi(s_{j})+O(h^{2}),$
(31)
$\displaystyle\int_{0}^{a}|G(\xi_{i},s)|\psi(s)ds=\sum_{j=0}^{N}h\rho_{j}|G(\xi_{i},s_{j})|\psi(s_{j})+O(h^{2}),$
(32)

###### Proposition 3.5

Under the assumptions of Theorem 2.2 we have the estimates

$\displaystyle\|\Psi_{k}-\psi_{k}\|_{\bar{\omega_{h}}}=O(h^{2}),$ (33)
$\displaystyle\|U_{k}-u_{k}\|_{\bar{\omega_{h}}}=O(h^{2}),,$ (34)

where $\|.\|_{\bar{\omega_{h}}}$ is the max-norm of grid function defined on
the grid $\bar{\omega_{h}}$.

Proof. We prove the proposition by induction. For $k=0$ we have at once
$\|\Psi_{0}-\psi_{0}\|_{\bar{\omega_{h}}}$ because
$\Psi_{0}(t_{i})=f(t_{i},0,0)$ and
$\psi_{0}(t_{i})=f(t_{i},0,0),\;i=\overline{0,N}$, too. Next, by (19) and
Proposition 3.3 we have

$\displaystyle u_{0}(t_{i})$
$\displaystyle=g(t_{i})+\int_{0}^{a}G(t_{i},s)\psi_{0}(s)ds$
$\displaystyle=g(t_{i})+\sum_{j=0}^{N}h\rho_{j}G(t_{i},s_{j})\psi_{0}(s_{j})+O(h^{2}).$

On the other hand, by (23) we have

$\displaystyle U_{0}(t_{i})$
$\displaystyle=g(t_{i})+\sum_{j=0}^{N}h\rho_{j}G(t_{i},s_{j})\Psi_{0}(s_{j}).$

Therefore,

$\displaystyle|U_{0}(t_{i})-u_{0}(t_{i})|=O(h^{2}).$

It implies $\|U_{0}-u_{0}\|_{\bar{\omega_{h}}}=O(h^{2})$. Thus, the estimates
(33) and (34) are valid for $k=0$.

Now, suppose that these estimates are valid for $k\geq 0$. We shall show that
they are valid for $k+1$. Indeed, from (20), (24) and the Lipshitz conditions
for the function $f(t,u,v)$ we have

$\begin{split}|\Psi_{k+1}(t_{i})-\psi_{k+1}(t_{i})|&=|f(t_{i},U_{k}(t_{i}),V_{k}(t_{i}))-f(t_{i},u_{k}(t_{i}),v_{k}(t_{i}))|\\\
&\leq
L_{1}|U_{k}(t_{i})-u_{k}(t_{i})|+L_{2}|V_{k}(t_{i})-v_{k}(t_{i})|.\end{split}$
(35)

Now estimate $|V_{k}(t_{i})-v_{k}(t_{i})|$. We have by Proposition 3.3

$\displaystyle v_{k}(t_{i})$
$\displaystyle=g(\varphi(t_{i}))+\int_{0}^{a}G(\varphi(t_{i}),s)\psi_{k}(s)ds$
$\displaystyle=g(\xi_{i})+\sum_{j=0}^{N}h\rho_{j}G(\xi_{i},s_{j})\psi_{k}(s_{j})+O(h^{2}).$

In view of (23) we have

$\begin{split}|V_{k}(t_{i})-v_{k}(t_{i})|&=|\sum_{j=0}^{N}h\rho_{j}G(\xi_{i},s_{j})(\Psi_{k}(s_{j})-\psi_{k}(s_{j}))|+O(h^{2})\\\
&\leq\sum_{j=0}^{N}h\rho_{j}|G(\xi_{i},s_{j})|\|\Psi_{k}-\psi_{k}\|_{\bar{\omega_{h}}}+O(h^{2}).\end{split}$
(36)

Notice that (32) for $\psi(s)=1$ gives

$\displaystyle\int_{0}^{a}|G(\xi_{i},s)|ds=\sum_{j=0}^{N}h\rho_{j}|G(\xi_{i},s_{j})|+O(h^{2}).$

From here it follows that

$\begin{split}&\sum_{j=0}^{N}h\rho_{j}|G(\xi_{i},s_{j})|=\int_{0}^{a}|G(\xi_{i},s)|ds+O(h^{2})\\\
&\leq\max_{0\leq t\leq
a}\int_{0}^{1}|G(t,s)|ds+O(h^{2})=M_{0}+O(h^{2}).\end{split}$

Thanks to this estimate, from (36) we obtain

$|V_{k}(t_{i})-v_{k}(t_{i})|\leq\|\Psi_{k}-\psi_{k}\|_{\bar{\omega_{h}}}+O(h^{2}).$

So, due to the induction hypothesis it implies

$\|V_{k}-v_{k}\|_{\bar{\omega_{h}}}=O(h^{2}).$ (37)

Combining the induction hypothesis
$\|U_{k}-u_{k}\|_{\bar{\omega_{h}}}=O(h^{2})$ and (37), from (35) we obtain

$\|\Psi_{k+1}-\psi_{k+1}\|_{\bar{\omega_{h}}}=O(h^{2}).$ (38)

In order to prove

$\|U_{k+1}-u_{k+1}\|_{\bar{\omega_{h}}}=O(h^{2}).$ (39)

we take into account that

$\displaystyle|U_{k+1}(t_{i})-u_{k+1}(t_{i})|\leq\sum_{j=0}^{N}h\rho_{j}|G(t_{i},s_{j})||\Psi_{k+1}(s_{j})-\psi_{k+1}(s_{j})|+O(h^{2}).$

Doing the similar argument as above and using the proved estimate (38) it is
easy to obtain

$\displaystyle|U_{k+1}(t_{i})-u_{k+1}(t_{i})|=O(h^{2}),$

or (39).
Thus, Proposition 3.5 is proved.

Now combining Proposition 3.5 with Theorem 3.1 we obtain the following result.

###### Theorem 3.6

Under the assumptions of Theorem 2.2 for the approximate solution of the
problem (2)-(3) obtained by the discrete iterative method (22)-(24) we have
the estimate

$\|U_{k}-u\|_{\bar{\omega_{h}}}\leq M_{0}p_{k}d+O(h^{2}),$

where $p_{k}$ and $d$ are defined by (21).

###### Remark 3.7

The results in Section 2 and 3 are obtained for the nonlinear third order FDE
with nonlinear term $f=f(t,u(t),u(\varphi(t)))$. Analogously, it is possible
to obtain similar results of existence and convergence of the iterative method
at continuous level for the general case

$f=f(t,u(t),u(\varphi(t)),u^{\prime}(\varphi_{1}(t)),u^{\prime\prime}(\varphi_{2}(t))).$

But for numerical realization of the iterative method it is needed to take
attention that the second derivative $\frac{\partial^{2}G(t,s)}{\partial
t^{2}}$ of the Green function has discontinuity at the line $s=t$. In this
case for computing integrals containing $\frac{\partial G(t,s)}{\partial t}$
and $\frac{\partial^{2}G(t,s)}{\partial t^{2}}$ it is needed to use the
formulas constructed in our previous work [3].

###### Remark 3.8

The technique developed in this paper for the third order nonlinear FDE can be
applied to nonlinear FDE of any order.

##  4 Examples

In all numerical examples below we perform the iterative method (22)-(24)
until $\|\Psi_{k}-\Psi_{k-1}\|_{\bar{\omega_{h}}}\leq 10^{-10}$. In the tables
of results for the convergence of the iterative method
$Error=\|U_{K}-u\|_{\bar{\omega_{h}}}$, $K$ is the number of iterations
performed.

Example 1. Consider the following problem

$\displaystyle u^{\prime\prime\prime}(t)$
$\displaystyle=e^{t}-\frac{1}{4}u(t)+\frac{1}{4}u^{2}(\frac{t}{2}),\quad
0<t<1,$ (40) $\displaystyle u(0)$
$\displaystyle=1,\;u^{\prime}(0)=1,\;u^{\prime}(1)=e$

with the exact solution $u(t)=e^{t}.$ The Green function for the above problem
is

$\displaystyle
G(t,s)=\left\\{\begin{array}[]{ll}\dfrac{s}{2}(t^{2}-2t+s),\quad 0\leq s\leq
t\leq 1,\\\ \,\,\dfrac{t^{2}}{2}(s-1),\quad 0\leq t\leq s\leq 1.\\\
\end{array}\right.$

So, we have

$M_{0}=\max_{0\leq t\leq a}\int_{0}^{1}|G(t,s)|ds=\frac{1}{2}.$

The second degree polynomial satisfying the boundary conditions of the problem
is

$g(t)=1+t+\frac{e-1}{2}t^{2}.$

Therefore, $\|g\|=2+\dfrac{e-1}{2}=2.7183$. In this example
$f(t,u,v)=e^{t}-\frac{1}{4}u+\frac{1}{4}v^{2}$. It is easy to verify that for
$M=6.5$ we have $|f(t,u,v)|\leq M$ in the domain $\mathcal{D}_{M}$ defined by
(11). Moreover, in this domain the function $f(t,u,v)$ satisfies the Lipshitz
conditions in $u$ and $v$ with the coefficients $L_{1}=\frac{1}{4}$ and
$L_{2}=1.7004$. Therefore, $q:=(L_{1}+L_{2})M_{0}=0.16$. Thus, all the
assumptions of Theorem 2.2 are satisfied. By the theorem, the problem (40) has
a unique solution. This is the above exact solution.

The results of convergence of the iterative method (22)-(24) are given in
Table 1.

Table 1: The convergence in Example 1.  $N$ | $h^{2}$ | $K$ | $Error$
---|---|---|---
50 | 4.0000e-04 | 3 | 6.1899e-05
100 | 1.0000e-04 | 3 | 1.5475e-05
150 | 4.4444e-05 | 3 | 6.877 -06
200 | 2.5000e-05 | 3 | 3.8688e-06
300 | 1.1111e-05 | 3 | 1.7195e-06
400 | 6.2500e-06 | 3 | 9.6721e-07
500 | 4.0000e-06 | 3 | 6.1901e-07
800 | 1.5625e-06 | 3 | 6.1901e-07
1000 | 1.0000e-06 | 3 | 1.5475e-07

From this table we see that the results of computation support the conclusion
that the accuracy of the iterative method is $O(h^{2})$.

###### Remark 4.1

Theorem 3.6 gives sufficient conditions for convergence of the iterative
method (22)-(24). In the cases when these conditions are not satisfied the
iterative also can converge to some solution. For example, for the case
$f(t,u,v)=e^{t}+u^{2}+v^{2}+1$ with the same boundary conditions as in (40)
the iterative method converges after 15 iterations. And for the case
$f(t,u,v)=e^{2t}-u^{3}+v^{2}+5$ after $16$ iterations the iterative process
reaches the $TOL=10^{-10}$. Notice that the number of iterations do not depend
on the grid size as in Example 1.

Example 2. Consider the following problem

$\displaystyle u^{\prime\prime\prime}(t)$
$\displaystyle=\sin(u^{2}(t))+\cos(u^{2}(t^{2})),\quad 0<t<1,$ (41)
$\displaystyle u(0)$
$\displaystyle=0,\;u^{\prime}(0)=\pi,\;u^{\prime}(1)=-\pi.$

For this problem $f(t,u,v)=\sin(u^{2})+\cos(v^{2})$ and $\varphi(t)=t^{2}$. It
is easy to verify that all the conditions of Theorem 3.1 are satisfied,
therefore the problem has a unique solution. Also, by Theorem 3.6 the
iterative method (22)-(24) converges. With $TOL=10^{-10}$ the iterative method
for any number of grid points stops after $8$ iterations. The graph of the
approximate solution is depicted in Figure 1.

Figure 1: The graph of the approximate solution in Example $2$.

Example 3. (Example 5 in [8]) Consider the problem

$\displaystyle u^{\prime\prime\prime}(t)$ $\displaystyle=-1+2u^{2}(t/2),\quad
0<t<\pi,$ (42) $\displaystyle u(0)$
$\displaystyle=0,\;u^{\prime}(0)=1,\;u(\pi)=0.$

which has exact solution $u(t)=\sin(t)$. For this problem the Green function
has the form

$\displaystyle
G(t,s)=\left\\{\begin{array}[]{ll}-\dfrac{t^{2}(\pi-s)^{2}}{2\pi^{2}}+\dfrac{(t-s)^{2}}{2},\quad
0\leq s\leq t\leq\pi,\\\ \,\,-\dfrac{t^{2}(\pi-s)^{2}}{2\pi^{2}},\quad 0\leq
t\leq s\leq\pi\\\ \end{array}\right.$

and the function $f(t,u,v)=-1+2v^{2}$.

The results of convergence of the iterative method (22)-(24) for this example
are given in Table 2

Table 2: The convergence in Example 3.  $N$ | $h^{2}$ | $K$ | $Error$
---|---|---|---
50 | 4.0000e-04 | 25 | 1.4455e-04
100 | 1.0000e-04 | 25 | 3.6142e-05
150 | 4.4444e-05 | 25 | 1.6063e-05
200 | 2.5000e-05 | 25 | 9.0345e-06
300 | 1.1111e-05 | 25 | 4.0155e-06
400 | 6.2500e-06 | 25 | 2.2587e-06
500 | 4.0000e-06 | 25 | 1.4456e-06
800 | 1.5625e-06 | 25 | 5.6467e-07
1000 | 1.0000e-06 | 25 | 3.6139e-07

From Table 2 it is seen also that the numerical method has the accuracy
$O(h^{2})$.

##  5 Conclusion

In this paper we have proposed a unified approach to nonlinear functional
differential equations via boundary value problems for nonlinear third order
functional differential equations as a particular case. We have established
the existence and uniqueness of solution and proved the convergence of the
discrete iterative method for finding the solution. Some examples demonstrate
the validity of the theoretical results and the efficiency of the numerical
method.

The proposed approach can be applied to boundary value problems for nonlinear
functional differential equations of any order. It also can be applied to
integro-differential equations.

## References

  * [1] Bica, A. M., Curila, M., Curila, S., Two-point boundary value problems associated to functional differential equations of even order solved by iterated splines, Applied Numerical Mathematics 110 (2016) 128-147.
  * [2] Dang, Q.A, Dang, Q.L., A unified approach to fully third order nonlinear boundary value problems, J. Nonlinear Funct. Anal. 2020 (2020), Article ID 9, https://doi.org/10.23952/jnfa.2020.9
  * [3] Dang, Q.A., Dang, Q.L., Simple numerical methods of second and third-order convergence for solving a fully third-order nonlinear boundary value problem, Numer Algor (2020). https://doi.org/10.1007/s11075-020-01016-2
  * [4] Dang, Q.A, Ngo, T.K.Q.: Existence results and iterative method for solving the cantilever beam equation with fully nonlinear term, Nonlinear Anal. Real World Appl. 36, 56-68 (2017)
  * [5] Dang, Q.A, Dang, Q.L., Ngo, T.K.Q.: A novel efficient method for nonlinear boundary value problems, Numer. Algor. 76, 427-439 (2017)
  * [6] Hale, Jack K., Theory of functional differential equations, Springer-Verlag, New York, 1977,
  * [7] Hou C-C, Simos T.E., Famelis I.T., Neural network solution of pantograph type differential equations. Math Meth Appl Sci. 2020;1-6. https://doi.org/10.1002/mma.6126
  * [8] Khuri, S.A. , Sayfy, A. , Numerical solution of functional differential equations: a Green’s function-based iterative approach, International Journal of Computer Mathematics Volume 95, 2018 - Issue 10 , 1937-1949.
  * [9] Raja, M. A. Z. , Numerical treatment for boundary value problems of Pantograph functional differential equation using computational intelligence algorithms, Applied Soft Computing, Volume 24, 2014, Pages 806-821.
  * [10] Reutskiy, S.Yu., A new collocation method for approximate solution of the pantograph functional differential equations with proportional delay, Applied Mathematics and Computation 266 (2015) 642-655.
  * [11] Yang, C., Hou, J. Jacobi spectral approximation for boundary value problems of nonlinear fractional pantograph differential equations. Numer Algor (2020). https://doi.org/10.1007/s11075-020-00924-7