using Plots, ComplexPhasePortrait, ApproxFun, SingularIntegralEquations, OscillatoryIntegrals, QuadGK
gr();
Dr Sheehan Olver
We want to solve ∫1−1log|x−t|u(t)dt=1x2+1
This is an inverse Hilbert problem so we know u(x)=−1√1−x2H[√1−t2f(t)]+C√1−x2
We determine the Cauchy transform of √1−x2f(x) using the usual methods: start with the ansatz ϕ(z)=√z−1√z+12i2zπ(1+z2)2
In other words, for ˜u(x)=2i√i−1√i+1π√1−x2(−14(x−i)2+2i8(x−i))−2i√−i−1√1−iπ√1−x2(14(x+i)2+i8(x+i))
To find C we impose the condition that ∫1−1log|t|u(t)dt=1
Check derivation
Let's check the derivation. First we can calculate u numerically:
L = SingularIntegral(0) : JacobiWeight(-0.5,-0.5,Chebyshev())
x = Fun()
g = (1/(x^2+1))
u = (π*L) \ g
π*logkernel(u, 0.1) - g(0.1)
-2.220446049250313e-16
Differentiating it satisfies Hu=g′/(−π):
hilbert(u,0.1) ≈ g'(0.1)/(-π) ≈ (2*0.1)/(π*(1+0.1^2)^2)
true
The Cauchy transform also works:
f = g'/(-π)
z = 1+im
ψ = z -> sqrt(z-1)sqrt(z+1)/(2im) * 2z/(π*(1+z^2)^2) -
sqrt(im-1)sqrt(im+1)/π * (-1/(4*(z-im)^2) + im/(8(z-im))) -
sqrt(-im-1)sqrt(1-im)/π * (1/(4*(z+im)^2) + im/(8(z+im)))
cauchy(sqrt(1-x^2)*f, z) ≈ ψ(z)
true
Therefore the Hilbert transform is given by:
H = x-> - 2im*sqrt(im-1)sqrt(im+1)/π * (-1/(4*(x-im)^2) + im/(8(x-im))) -
2im*sqrt(-im-1)sqrt(1-im)/π * (1/(4*(x+im)^2) + im/(8(x+im)))
hilbert(sqrt(1-x^2)*f)(0.1) ≈ H(0.1)
true
We thus have an expression for u, we are just missing the constant:
ũ = x -> -H(x)/sqrt(1-x^2)
C = u(0.0) + H(0.0)
u(0.1) ≈ ũ(0.1) + C/sqrt(1-0.1^2)
true
We can find C interms of the relevant integral, which we call D:
D = quadgk(x -> x == 0 ? 0.0 : ũ(x)*log(abs(x)),-1,1)[1]
C, (1-D)/(π*log(1/2))
(-0.3247204711377926 + 0.0im, -0.32472047136213555 + 0.0im)
Since ∫1−1xdx=0, we have no logarithmic growth at infinity. Thus by the formula in Lecture 19 we find for U(x)=∫1xtdt=1−x22
So we just have to work out the Cauchy transform. Try as an ansatz (1−z22)log(z−1)−log(z+1)2πi
We can confirm the formula:
z = 2+im; π*logkernel(x,z) ≈ real((1-z^2)/2 * (log(z-1)-log(z+1)) - z)
true
From Lecture 19 we know that we want to solve
We write v(x,y)=∫1−1u(t)log|z−t|dt+log|z−i|
We can actually solve this numerically:
ũ = SingularIntegral(0) \ (-log(abs(x-im))/π)
u = ũ - sum(ũ)/(π*sqrt(1-x^2)) # ensure integrates to zero
v = z -> π*logkernel(u, z) + log(abs(z-im))
D = v(0)
0.18822640645958963
xx = yy = -4:0.011:4
V = v.(xx' .+ im*yy)
contour(xx, yy, V; nlevels=50)
plot!(domain(x); color=:black, legend=false)
As usual for logarithmic singular integral equations we want to solve ∫−1−1u(t)t−xdt=f′(x)
Verification Let's check our work. First we see that the Hilbert transform of u does indeed satisfy the specified equation (using the numerical u
as calculated above):
fp = x -> x/(x^2+1)
π*hilbert(u, 0.1) ≈ fp(0.1)
true
And the Cauchy transform of √1−x2f′(x) satisfies the derived formula:
C = z -> z * sqrt(z-1)sqrt(z+1)/(2im*(z^2+1)) - 1/(2im) + im*sqrt(im-1)sqrt(im+1)/(4(z-im)) +
im*sqrt(-im-1)sqrt(-im+1)/(4(z+im))
cauchy(sqrt(1-x^2)fp(x), z) ≈ C(z)
true
Therefore we can invert to the Hilbert transform for f′:
w = -hilbert(sqrt(1-x^2)fp(x))/sqrt(1-x^2)
w̃ = sqrt(im-1)sqrt(im+1)/(2(x-im) * sqrt(1-x^2)) +
sqrt(-im-1)sqrt(-im+1)/(2(x+im) * sqrt(1-x^2))
hilbert(w,0.1) ≈ hilbert(w̃,0.1) ≈ fp(0.1)
true
And we have recovered u up to C/√1−x2:
C = 1/π
u(0.1) ≈ w̃(0.1)/π + C/sqrt(1-0.1^2)
true
It has exponential decay in the right-half plane, therefore eγxf(x)=eγx1+ex
We can verify this by exact computation using Residue calculus: for 0<ℑs<1, we can integrate over a rectangle to get: $$ \left(\int_{-R}^R + \int_R^{2\I \pi + R} + \int_{2 \I \pi + R}^{2\I \pi - R} + \int_{2 \I \pi - R}^{-R} \right) {\E^{-\I s x} \over 1 + \E^x} \dx = 2 \pi \I \Res_{z = \I \pi } {\E^{-\I s z} \over 1 + \E^z} =
{\E^{-\I s (R + \I t)} \over 1 + \E^{R + \I t}} = {\E^{-\I R \Re s + R \Im s + t} \over 1 + \E^{R + \I t}} \rightarrow 0 and
\int_{2 \I \pi + \infty}^{2\I \pi - \infty} {\E^{-\I s t} \over 1 + \E^t} \dt = \int_{\infty}^{-\infty} {\E^{-\I s (x+2 \I \pi)} \over 1 + \E^x} \dx = -\E^{2 \pi s} \int_{-\infty}^\infty {\E^{-\I s x} \over 1 + \E^x} \dx Therefore,wehave
phaseplot(-3..3, -10..10, z -> 1/(1+exp(z))) #integrand
phaseplot(-3..3, -3..3, z -> im*π*csch(π*z)) # transform
Here eγxf(x)=e(γ+2)x has decay at +∞ proved γ<−2, hence we have the strip ℑs<−2.
Indeed, its Fourier transform is −i2i+s
Here it's ℑs>0: unlike 1.1, we now have decay at x→∞ since fL(x) is identically zero.
It's Fourier transform is determinable by integration-by-parts: ˆf(s)=∫0−∞xe−isxdx=1is∫0−∞e−isxdx=1s2
The Fourier transforms are given above.
It's actually an entire function, but non-decaying. This is hinting at the relationship between smoothness of a function and decay of its Fourier transform, and vice-versa: since δ(x) "decays" to all orders, we expect its Fourier transform to be entire, but since its n ot smooth at all, we expect no decay, so on a formal level we can predict the analyticity properties.
α = 0.3
x = Fun(0..100)
f = 1 + α*x
h = s -> (im/s + α/s^2)
γ = -0.5 # we take the Fourier transform on R + im*γ
s = -0.5 + im*γ
-sum(f*exp(-im*s*x)) - h(s)
5.773159728050814e-15 + 7.549516567451064e-15im
Transforming the equation, we have Φ+(s)−(1+ˆK(s))Φ−(s)=is+αs2
g = s -> (4+s^2)/(1+s^2)
κ = z -> imag(z) > γ ? (z+im*2)/(z+im) :
(z-im)/(z-im*2)
phaseplot(-3..3, -3..3, κ)
s = 0.1 + γ*im
κ₊ = κ(s + eps()*im)
κ₋ = κ(s - eps()*im)
κ₊ - κ₋*g(s)
-1.3322676295501878e-15 - 2.7755575615628914e-16im
We thus get the RH problem Y+(s)−Y−(s)=h(s)/κ+(s)=(is+αs2)s+is+2i
s = 0.1 + γ*im
Y = z -> imag(z) > γ ? im*(2+α)/(4*(z+2im)) :
- α/(2z^2) + im*(α-2)/(4z)
Y₊ = Y(s + eps()*im)
Y₋ = Y(s - eps()*im)
Y₊ - Y₋ - h(s)/κ₊
5.551115123125783e-16 - 4.440892098500626e-16im
We therefore have
Φ(z)=κ(z)Y(z)={i(2+α)4(z+i)ℑz>γ(−α2z2+i(α−2)4z)z−iz−2iℑz<γΦ = z -> imag(z) > γ ? im*(2+α)/(4*(z+im)) :
(-α/(2z^2) + im*(α-2)/(4z))*(z-im)/(z-2im)
Φ₊ = Φ(s+eps()im)
Φ₋ = Φ(s-eps()im)
Φ₊ - Φ₋*g(s) - h(s)
8.881784197001252e-16 - 1.1102230246251565e-15im
Finally, we recover the solution by inverting Φ−, using Residue calculus in the upper half plane: for x>0 we have u(x)=12π∫∞+iγ−∞+iγ(−α2z2+i(α−2)4z)z−iz−2ieizxdz=i(Resz=0+Resz=2i)(−α2z2+i(α−2)4z)z−iz−2ieizx=1+xα4−α+14e−2x
Did it work? yes:
t = Fun(0 .. 50)
u = (1+t*α)/4 - (α-1)/4*exp(-2t)
x = 0.1
u(x) + 3/2*sum(exp(-abs(t-x))*u) , f(x)
(1.0300000000000011, 1.03)
Setting up the problem as above, we arrive at a degenerate RH problem: Φ+(s)−g(s)Φ−(s)=h(s)
Suppose we allow κ−(s)∼s to have growth, then we can write κ(z)={1z+iαℑz>γz−iα2αℑz<γ
α = 0.3
g = s -> (2α)/(α^2+s^2)
h = s -> (im/s + α/s^2)
κ = z -> imag(z) > γ ? 1/(z + im*α) :
(z-im*α)/(2α)
phaseplot(-3..3, -3..3, κ)
s = 0.1 + γ*im
κ₊ = κ(s + eps()*im)
κ₋ = κ(s - eps()*im)
κ₊ - κ₋*g(s)
2.6645352591003757e-15 + 1.7763568394002505e-15im
Then we have h(s)/κ+(s)=is2+α2s2=i+iα2s2
s = 0.1 + γ*im
Y = z -> imag(z) > γ ? im :
-im*α^2/s^2
Y₊ = Y(s + eps()*im)
Y₋ = Y(s - eps()*im)
Y₊ - Y₋ - h(s)/κ₊
1.3877787807814457e-16 + 7.771561172376096e-16im
Putting things together, we get Φ(z)=κ(z)Y(z)={iz+iαℑz>γ−iα2z2z−iα2αℑz<γ
Φ = z -> imag(z) > γ ? im/(z + im*α) :
-im*α^2/z^2* (z-im*α)/(2α)
Φ₊ = Φ(s+eps()im)
Φ₋ = Φ(s-eps()im)
Φ₊ - Φ₋*g(s) - h(s)
-3.9968028886505635e-15 + 4.218847493575595e-15im
We now invert the Fourier transform of Φ−(s) using Jordan's lemma: u(x)=12π∫∞+iγ−∞+iγΦ−(s)eisxds=α2Resz=0z−iαz2eizx=α2(1+xα)
t = Fun(0 .. 200)
u = α*(1+t*α)/2
x = 0.1
sum(exp(-α*abs(t-x))*u) - (1 + α*x)
1.199040866595169e-14
where g(s)=1−2λs2+1=s2+1−2λs2+1=(s−iγ)(s+iγ)(s+i)(s−i)
t = Fun(0 .. 200)
λ = 0.1
γ = sqrt(1-2λ)
u = t/γ^2 - exp(-t*γ)*(γ-1)/γ^2
x = 0.1
u(x) - λ*sum(exp(-abs(t-x))*u) - x
-3.239075674343894e-14
Oddly, this is definitely a solution, but not in the form the question asked for. To get the other solution, consider now the bad winding number case of −1<δ<−γ. Motivated by 2.2, what if we allow κ to have different behaviour? Consider κ(z)={1z+iℑz>δ(z−i)(z−iγ)(z+iγ)ℑz<δ
Thus we want to solve Y+(s)−Y−(s)=h(s)κ+(s)−1=s+is2=1s+is2
What's the moral of the story?
where p(x)=725∫∞−∞e−5|x−t|uR(t)dt=725∫∞0e−5|x−t|uR(t)dt.
κ = z -> imag(z) > γ ?
(z+3im)*(z+4im)/(z+5im) :
(z-5im)/((z-3im)*(z-4im))
γ = -1.0
s = 0.1+γ*im
g = s -> (s^2+9)*(s^2+16)/(s^2+25)
κ(s+eps()im) - g(s)κ(s-eps()im)
-5.551115123125783e-17 - 2.220446049250313e-16im
Writing Φ(z)=κ(z)Y(z) we get the subtractive RH problem Y+(s)−Y−(s)=h(s)κ+(s)=(α+1is)s+5i(s+3i)(s+4i)
We can now invert the Fourier transform of Φ−(s)=s−5i(s−3i)(s−4i)512s
Here's we check the solution:
t = Fun(0 .. 200)
u = -25/144 - 5exp(-4t)/48 + 5exp(-3t)/18
x = 1.1
u''(x) - 72/5*sum(exp(-5abs(x-t))*u)
1.0000000000000018
Here we check the jump of Y:
α = u'(0)
h = s -> α + 1/(im*s)
Y = z -> imag(z) > γ ?
(-(α+1/4)/(z+4im) + (2/3 + 2α)/(z+3im)) :
5/(12z)
Y(s+eps()im) - Y(s-eps()im) - h(s)/κ(s + eps()im)
2.7755575615628914e-17 + 5.551115123125783e-17im
Here we check the jump of Φ:
γ = -1.0
Φ = z -> imag(z) > γ ?
(z+3im)*(z+4im)/(z+5im) * (-(α+1/4)/(z+4im) + (2/3 + 2α)/(z+3im)) :
(z-5im)/((z-3im)*(z-4im)) * 5/(12z)
Φ(s + eps()*im) - g(s)*Φ(s - eps()*im) - h(s)
-1.1102230246251565e-16 + 1.3877787807814457e-17im