residue

A fun application of Green’s functions and geometric algebra: Residue calculus

November 2, 2025 math and physics play No comments , , , , , , , , , , , , , , , , , ,

[Click here for a PDF version of this post]

Motivation.

A fun application of both Green’s functions and geometric algebra is to show how the Cauchy integral equation can be expressed in terms of the Green’s function for the 2D gradient. This is covered, almost as an aside, in [1]. I found that treatment a bit hard to understand, so I am going to work through it here at my own pace.

Complex numbers in geometric algebra.

Anybody who has studied geometric algebra is likely familiar with a variety of ways to construct complex numbers from geometric objects. For example, complex numbers can be constructed for any plane. If \( \Be_1, \Be_2 \) is a pair of orthonormal vectors for some plane in \(\mathbb{R}^N\), then any vector in that plane has the form
\begin{equation}\label{eqn:residueGreens:20}
\Bf = \Be_1 u + \Be_2 v,
\end{equation}
has an associated complex representation, by simply multiplying that vector one of those basis vectors. For example, if we pre-multiply \( \Bf \) by \( \Be_1 \), forming
\begin{equation}\label{eqn:residueGreens:40}
\begin{aligned}
z
&= \Be_1 \Bf \\
&= \Be_1 \lr{ \Be_1 u + \Be_2 v } \\
&= u + \Be_1 \Be_2 v.
\end{aligned}
\end{equation}

We may identify the unit bivector \( \Be_1 \Be_2 \) as an imaginary, designed by \( i \), since it has the expected behavior
\begin{equation}\label{eqn:residueGreens:60}
\begin{aligned}
i^2 &=
\lr{\Be_1 \Be_2}^2 \\
&=
\lr{\Be_1 \Be_2}
\lr{\Be_1 \Be_2} \\
&=
\Be_1 \lr{\Be_2
\Be_1} \Be_2 \\
&=
-\Be_1 \lr{\Be_1
\Be_2} \Be_2 \\
&=
-\lr{\Be_1 \Be_1}
\lr{\Be_2 \Be_2} \\
&=
-1.
\end{aligned}
\end{equation}

Complex numbers are seen to be isomorphic to even grade multivectors in a planar subspace. The imaginary is the grade-two pseudoscalar, and geometrically is an oriented unit area (bivector.)

Cauchy-equations in terms of the gradient.

It is natural to wonder about the geometric algebra equivalents of various complex-number relationships and identities. Of particular interest for this discussion is the geometric algebra equivalent of the Cauchy equations that specify required conditions for a function to be differentiable.

If a complex function \( f(z) = u(z) + i v(z) \) is differentiable, then we must be able to find the limit of
\begin{equation}\label{eqn:residueGreens:80}
\frac{\Delta f(z_0)}{\Delta z} = \frac{f(z_0 + h) – f(z_0)}{h},
\end{equation}
for any complex \( h \rightarrow 0 \), for any possible trajectory of \( z_0 + h \) toward \( z_0 \). In particular, for real \( h = \epsilon \),
\begin{equation}\label{eqn:residueGreens:100}
\lim_{\epsilon \rightarrow 0} \frac{u(x_0 + \epsilon, y_0) + i v(x_0 + \epsilon, y_0) – u(x_0, y_0) – i v(x_0, y_0)}{\epsilon}
=
\PD{x}{u(z_0)} + i \PD{x}{v(z_0)},
\end{equation}
and for imaginary \( h = i \epsilon \)
\begin{equation}\label{eqn:residueGreens:120}
\lim_{\epsilon \rightarrow 0} \frac{u(x_0, y_0 + \epsilon) + i v(x_0, y_0 + \epsilon) – u(x_0, y_0) – i v(x_0, y_0)}{i \epsilon}
=
-i\lr{ \PD{y}{u(z_0)} + i \PD{y}{v(z_0)} }.
\end{equation}
Equating real and imaginary parts, we see that existence of the derivative requires
\begin{equation}\label{eqn:residueGreens:140}
\begin{aligned}
\PD{x}{u} &= \PD{y}{v} \\
\PD{x}{v} &= -\PD{y}{u}.
\end{aligned}
\end{equation}
These are the Cauchy equations. When the derivative exists in a given neighbourhood, we say that the function is analytic in that region. If we use a bivector interpretation of the imaginary, with \( i = \Be_1 \Be_2 \), the Cauchy equations are also satisfied if the gradient of the complex function is zero, since
\begin{equation}\label{eqn:residueGreens:160}
\begin{aligned}
\spacegrad f
&=
\lr{ \Be_1 \partial_x + \Be_2 \partial_y } \lr{ u + \Be_1 \Be_2 v } \\
&=
\Be_1 \lr{ \partial_x u – \partial_y v } + \Be_2 \lr{ \partial_y u + \partial_x v }.
\end{aligned}
\end{equation}
We see that the geometric algebra equivalent of the Cauchy equations is simply
\begin{equation}\label{eqn:residueGreens:200}
\spacegrad f = 0.
\end{equation}
Roughly speaking, we may say that a function is analytic in a region, if the Cauchy equations are satisfied, or the gradient is zero, in a neighbourhood of all points in that region.

A special case of the fundamental theorem of geometric calculus.

Given an even grade multivector \( \psi \in \mathbb{R}^2 \) (i.e.: a complex number), we can show that
\begin{equation}\label{eqn:residueGreens:220}
\int_A \spacegrad \psi d^2\Bx = \oint_{\partial A} d\Bx \psi.
\end{equation}
Let’s get an idea why this works by expanding the area integral for a rectangular parameterization
\begin{equation}\label{eqn:residueGreens:240}
\begin{aligned}
\int_A \spacegrad \psi d^2\Bx
&=
\int_A \lr{ \Be_1 \partial_1 + \Be_2 \partial_2 } \psi I dx dy \\
&=
\int \Be_1 I dy \evalrange{\psi}{x_0}{x_1}
+
\int \Be_2 I dx \evalrange{\psi}{y_0}{y_1} \\
&=
\int \Be_2 dy \evalrange{\psi}{x_0}{x_1}
–
\int \Be_1 dx \evalrange{\psi}{y_0}{y_1} \\
&=
\int d\By \evalrange{\psi}{x_0}{x_1}
–
\int d\Bx \evalrange{\psi}{y_0}{y_1}.
\end{aligned}
\end{equation}
We took advantage of the fact that the \(\mathbb{R}^2\) pseudoscalar commutes with \( \psi \). The end result, is illustrated in fig. 1, shows pictorially that the remaining integral is an oriented line integral.

fig. 1. Oriented multivector line integral.

 

If we want to approximate a more general area, we may do so with additional tiles, as illustrated in fig. 2. We may evaluate the area integral using the line integral over just the exterior boundary using such a tiling, as any overlapping opposing boundary contributions cancel exactly.

fig. 2. A crude circular tiling approximation.

 

The reason that this is interesting is that it allows us to re-express a complex integral as a corresponding multivector area integral. With \( d\Bx = \Be_1 dz \), we have
\begin{equation}\label{eqn:residueGreens:260}
\oint dz\, \psi = \Be_1 \int \spacegrad \psi d^2\Bx.
\end{equation}

The Cauchy kernel as a Green’s function.

We’ve previously derived the Green’s function for the 2D Laplacian, and found
\begin{equation}\label{eqn:residueGreens:280}
\tilde{G}(\Bx, \Bx’) = \inv{2\pi} \ln \Abs{\lr{\Bx – \Bx’}},
\end{equation}
which satisfies
\begin{equation}\label{eqn:residueGreens:300}
\delta^2(\Bx – \Bx’) = \spacegrad^2 \tilde{G}(\Bx, \Bx’) = \spacegrad \lr{ \spacegrad \tilde{G}(\Bx, \Bx’) }.
\end{equation}
This means that \( G(\Bx, \Bx’) = \spacegrad \tilde{G}(\Bx, \Bx’) \) is the Green’s function for the gradient. That Green’s function is
\begin{equation}\label{eqn:residueGreens:320}
\begin{aligned}
G(\Bx, \Ba)
&= \inv{2 \pi} \frac{\spacegrad \Abs{\Bx – \Ba}}{\Abs{\Bx – \Ba}} \\
&= \inv{2 \pi} \frac{\Bx – \Ba}{\Abs{\Bx – \Ba}^2}.
\end{aligned}
\end{equation}
We may cast this Green’s function into complex form with \( z = \Be_1 \Bx, a = \Be_1 \Ba \). In particular
\begin{equation}\label{eqn:residueGreens:340}
\begin{aligned}
\inv{z – a}
&=
\frac{(z – a)^\conj}{\Abs{z – a}^2} \\
&=
\frac{(z – a)^\conj}{\Abs{z – a}^2} \\
&=
\frac{\Bx – \Ba}{\Abs{\Bx – \Ba}^2} \Be_1 \\
&=
2 \pi G(\Bx, \Ba) \Be_1.
\end{aligned}
\end{equation}

Cauchy’s integral.

With
\begin{equation}\label{eqn:residueGreens:360}
\psi = \frac{f(z)}{z – a},
\end{equation}
using \ref{eqn:residueGreens:260}, we can now evaluate
\begin{equation}\label{eqn:residueGreens:265}
\begin{aligned}
\oint dz\, \frac{f(z)}{z – a}
&= \Be_1 \int \spacegrad \frac{f(z)}{z – a} d^2\Bx \\
&= \Be_1 \int \lr{ \frac{\spacegrad f(z)}{z – a} + \lr{ \spacegrad \inv{z – a}} f(z) } I dA \\
&= \Be_1 \int f(z) \spacegrad 2 \pi G(\Bx – \Ba) \Be_1 I dA \\
&= 2 \pi \Be_1 \int \delta^2(\Bx – \Ba) \Be_1 f(\Bx) I dA \\
&= 2 \pi \Be_1^2 f(\Ba) I \\
&= 2 \pi I f(a),
\end{aligned}
\end{equation}
where we’ve made use of the analytic condition \( \spacegrad f = 0 \), and the fact that \( f \) and \( 1/(z-a) \), both even multivectors, commute.

The Cauchy integral equation
\begin{equation}\label{eqn:residueGreens:380}
f(a) = \inv{2 \pi I} \oint dz\, \frac{f(z)}{z – a},
\end{equation}
falls out naturally. This sort of residue calculation always seemed a bit miraculous. By introducing a geometric algebra encoding of complex numbers, we get a new and interesting interpretation. In particular,

  1. the imaginary factor in the geometric algebra formulation of this identity is an oriented unit area coming directly from the area element,
  2. the factor of \( 2 \pi \) comes directly from the Green’s function for the gradient,
  3. the fact that this particular form of integral picks up only the contribution at the point \( z = a \) is no longer mysterious seeming. This is directly due to delta-function filtering.

Also, if we are looking for an understanding of how to generalize the Cauchy equation to more general multivector functions, we now also have a good clue how that would be done.

References

[1] C. Doran and A.N. Lasenby. Geometric algebra for physicists. Cambridge University Press New York, Cambridge, UK, 1st edition, 2003.

Green’s function for the 3D wave equation.

October 22, 2025 math and physics play , , , , , , , ,

[Click here for a PDF version of this and related posts]

We’ve now evaluated the 1D Green’s function and 2D Green’s function for the wave equation.

For the sake of completeness, now let’s evaluate the Green’s function for the 3D wave equation operator. Again with \( \Br = \Bx – \Bx’, \tau = t – t’ \) we want the \( \epsilon \rightarrow 0 \) limit of
\begin{equation}\label{eqn:waveEquationGreens:1480}
G_\epsilon(\Br, \tau)
=
\inv{\lr{2 \pi}^4} \int d^3 \Bk d\omega \frac{e^{j \Bk \cdot \Br + j \omega \tau}}{(\omega/c – j \epsilon/c)^2 – \Bk^2}.
\end{equation}
For \(\epsilon > 0 \) this will presumably give us the retarded solution, with advanced for \( \epsilon < 0 \). We are using the nice pole displacement that leaves both poles on the same side of the upper or lower half plane, depending on the sign of \( \epsilon \). Let’s only do the \( \epsilon > 0 \) case by hand. Evaluating the \( \omega \) integral first with an upper half plane contour, we have
\begin{equation}\label{eqn:waveEquationGreens:1500}
\begin{aligned}
G_\epsilon(\Br, \tau)
&=
\frac{c^2}{\lr{2 \pi}^4} \int d^3 \Bk e^{j \Bk \cdot \Br}
\frac{e^{j \omega \tau}}{
\lr{\omega – \lr{ j \epsilon – \Abs{\Bk} c}}
\lr{\omega – \lr{ j \epsilon + \Abs{\Bk} c}}
} \\
&=
\frac{j c^2}{\lr{2 \pi}^3} \Theta(\tau) \int d^3 \Bk e^{j \Bk \cdot \Br}
\lr{
\evalbar{\frac{e^{j \omega \tau}}{\omega – \lr{ j \epsilon – \Abs{\Bk} c}}}{\omega = j \epsilon + \Abs{\Bk} c}
+
\evalbar{\frac{e^{j \omega \tau}}{\omega – \lr{ j \epsilon + \Abs{\Bk} c}}}{\omega = j \epsilon – \Abs{\Bk} c}
} \\
&=
\frac{j c^2}{\lr{2 \pi}^3} \Theta(\tau) e^{-\epsilon \tau} \int d^3 \Bk e^{j \Bk \cdot \Br}
\lr{
\frac{e^{j \Abs{\Bk} c \tau}}{2 \Abs{\Bk} c}
–
\frac{e^{-j \Abs{\Bk} c \tau}}{2 \Abs{\Bk} c}
} \\
&=
-\frac{c}{\lr{2 \pi}^3} \Theta(\tau) e^{-\epsilon \tau} \int \frac{d^3 \Bk}{\Abs{\Bk}} e^{j \Bk \cdot \Br}
\sin\lr{ \Abs{\Bk} c \tau }.
\end{aligned}
\end{equation}
We can evaluate the \( \epsilon \rightarrow 0 \) limit, and switch to spherical coordinates in k-space. Let \( \Br = r \Be_3 \)
\begin{equation}\label{eqn:waveEquationGreens:1520}
G(\Br, \tau)
=
-\frac{c}{\lr{2 \pi}^3} \Theta(\tau)
\int_{k = 0}^\infty \frac{k^2 dk}{k}
\int_{\phi = 0}^{2 \pi} d\phi
\int_{\theta = 0}^{\pi} \sin\theta d\theta
e^{j k r \cos\theta} \sin\lr{ k c \tau }.
\end{equation}
With \( u = \cos\theta \), this gives
\begin{equation}\label{eqn:waveEquationGreens:1540}
\begin{aligned}
G(\Br, \tau)
&=
\frac{c}{\lr{2 \pi}^2} \Theta(\tau)
\int_{k = 0}^\infty k dk \sin\lr{ k c \tau }
\int_{u = -1}^{1} du
e^{j k r u} \\
&=
\frac{c}{\lr{2 \pi}^2} \Theta(\tau)
\int_{k = 0}^\infty k dk \sin\lr{ k c \tau }
\lr{ \frac{e^{-j k r }}{j k r} – \frac{e^{j k r }}{j k r} } \\
&=
-\frac{c}{2 \pi^2 r} \Theta(\tau) \int_{k = 0}^\infty dk \sin\lr{ k c \tau } \sin\lr{ k r} \\
&=
-\frac{c}{4 \pi^2 r} \Theta(\tau) \int_{k = 0}^\infty dk
\lr{
\cos\lr{ k( c \tau – r ) }
–
\cos\lr{ k( c \tau + r ) }
} \\
&=
-\frac{c}{8 \pi^2 r} \Theta(\tau) \int_{k = -\infty}^\infty dk
\lr{
\cos\lr{ k( c \tau – r ) }
–
\cos\lr{ k( c \tau + r ) }
} \\
&=
-\frac{c}{8 \pi^2 r} \Theta(\tau) \int_{-\infty}^\infty dk
\lr{
e^{ j k( c \tau – r ) } – e^{ j k( c \tau + r ) }
} \\
&=
-\frac{c}{4 \pi r} \Theta(\tau)
\lr{
\delta( c \tau – r )
–
\delta( c \tau + r )
} \\
&=
-\frac{1}{4 \pi r} \Theta(\tau)
\lr{
\delta( \tau – r/c )
–
\delta( \tau + r/c )
}.
\end{aligned}
\end{equation}
Observe that the second delta function only has a value when \( \tau = -r/c \), but \( \Theta(-r/c) = 0 \). Similarly, the first delta function only has a value for \( \tau = r/c \ge 0 \), where the Heaviside step function is unity. That means we can simplify this to just
\begin{equation}\label{eqn:waveEquationGreens:1560}
\boxed{
G(\Br, \tau) = -\frac{1}{4 \pi \Abs{\Br}} \delta( \tau – \Abs{\Br}/c ),
}
\end{equation}
as expected.

Again, sort of sadly, we can skip all the fun and evaluate most of this in Mathematica. It needs only minor hand-holding to extract the delta function semantics. The retarded derivation is shown in fig. 1, and the advanced derivation in fig. 2.

fig. 1. Retarded 3D Green’s function for the wave equation.

fig. 2. Advanced 3D Green’s function for the wave equation.

Derivation of the 2D Green’s function for the wave equation operator.

October 22, 2025 math and physics play , , , , , , ,

[Click here for a PDF version of this post]

While it was difficult to attempt to verify the 2D Green’s function, it actually turns out to be fairly easy to derive it, provided we pick an alternate pole displacement from the 1D evaluation to make our lives easier.

With \( \Br = \Bx – \Bx’ \), and \( \tau = t – t’ \), and \( \epsilon > 0 \), we can form
\begin{equation}\label{eqn:waveEquationGreens:1340}
G_\epsilon(\Br, \tau) = \frac{c^2}{\lr{2 \pi}^3} \int d^2 \Bk d\omega \frac{ e^{j \Bk \cdot \Br + j \omega \tau}}{\lr{\omega -j \epsilon}^2 – \Bk^2 c^2 }
\end{equation}
This pole displacement has the nice property that both poles live in the upper half plane, so for \( \tau > 0 \), we have
\begin{equation}\label{eqn:waveEquationGreens:1360}
\begin{aligned}
G_\epsilon(\Br, \tau)
&= \frac{c^2}{\lr{2 \pi}^3} \int d^2 \Bk d\omega \frac{ e^{j \Bk \cdot \Br + j \omega \tau}}{
\lr{\omega -\lr{ \Abs{\Bk} c – j \epsilon}}
\lr{\omega -\lr{ -\Abs{\Bk} c – j \epsilon}}
} \\
&=
\Theta(\tau) \frac{c^2 j}{\lr{2 \pi}^2} \int d^2 \Bk e^{j \Bk \cdot \Br}
\lr{
\evalbar{
\frac{e^{ j \omega \tau}}{ \lr{\omega -\lr{ -\Abs{\Bk} c – j \epsilon}} }
}
{\omega = \Abs{\Bk} c – j \epsilon}
+
\evalbar{
\frac{e^{ j \omega \tau}}{ \lr{\omega -\lr{ \Abs{\Bk} c – j \epsilon}} }
}
{\omega = -\Abs{\Bk} c – j \epsilon}
}
\\
&=
\Theta(\tau) \frac{c^2 j}{\lr{2 \pi}^2} \int d^2 \Bk e^{j \Bk \cdot \Br} e^{ -j \epsilon \tau }
\lr{
\frac{e^{ j \Abs{\Bk} c \tau}}{ 2 \Abs{\Bk} c }
+
\frac{e^{ -j \Abs{\Bk} c \tau}}{ -2 \Abs{\Bk} c }
} \\
&=
\Theta(\tau) \frac{j^2 c}{\lr{2 \pi}^2} \int_{k=0}^\infty k dk \int_{\phi=0}^{2 \pi} d\phi e^{j k\Abs{\Br} \cos\phi } e^{ -j \epsilon \tau } \frac{\sin\lr{ k c \tau }}{k}.
\end{aligned}
\end{equation}
We’ve now successfully removed the singularity, and can evaluate the \(\epsilon \rightarrow 0 \) limit. We may also evaluate the \( \phi \) integral, remembering that
\begin{equation}\label{eqn:waveEquationGreens:1380}
\int_0^{2 \pi} e^{j \Abs{a} \cos\phi} d\phi = 2 \pi J_0(\Abs{a}),
\end{equation}
to find
\begin{equation}\label{eqn:waveEquationGreens:1400}
G(\Br, \tau) = -\Theta(\tau) \frac{c}{2 \pi} \int_{k=0}^\infty dk J_0(k\Abs{\Br}) \sin\lr{ k c \tau }.
\end{equation}
This integral yields easily to Mathematica, and we find
\begin{equation}\label{eqn:waveEquationGreens:1420}
G(\Br, \tau) = -\Theta(\tau) \frac{c}{2 \pi} \frac{\Theta(c \tau – \Abs{\Br})}{\sqrt{(c\tau)^2 – \Br^2}}.
\end{equation}
However, since \( \Theta(c \tau – \Abs{\Br}) = 1 \) only for \( \tau > \Abs{\Br}/c \), the \( \Theta(\tau) \) factor is redundant, and we find
\begin{equation}\label{eqn:waveEquationGreens:1440}
\boxed{
G(\Br, \tau) = – \frac{1}{2 \pi} \frac{\Theta(c \tau – \Abs{\Br})}{\sqrt{\tau^2 – \Br^2/c^2}},
}
\end{equation}
which matches the retarded Green’s function claimed by Grok.

Repeating this analysis for \( \tau < 0, \epsilon < 0 \), we find
\begin{equation}\label{eqn:waveEquationGreens:1460}
G(\Br, \tau) = -\Theta(-\tau) \frac{c}{2 \pi} \frac{\Theta(-c \tau – \Abs{\Br})}{\sqrt{(c\tau)^2 – \Br^2}},
\end{equation}
which we also see matches the Grok result for the advanced Green’s function. Both of these computations can be trivially performed in Mathematica following the same steps (taking all the fun from the story.) The advanced integral evaluation is shown in fig. 1 as an example.

fig. 1. Advanced 2D Green’s function for wave equation operator.

Deriving the Green’s functions for the 1D wave equation.

October 13, 2025 math and physics play , , , , , , , , ,

[Click here for a PDF version of this (and previous) post]

Having had trouble verifying the 2D Green’s function, let’s try deriving them ourselves.

Setup.

Let’s try solving the forced wave equation
\begin{equation}\label{eqn:waveEquationGreens:860}
\lr{ \spacegrad^2 – \inv{c^2}\frac{\partial^2}{\partial t^2} } f(x,t) = g(x,t),
\end{equation}
using Fourier transform pairs
\begin{equation}\label{eqn:waveEquationGreens:880}
\begin{aligned}
F(\Bx, t) &= \inv{\lr{\sqrt{2 \pi}}^{N+1}} \int e^{j \Bk \cdot \Bx + j \omega t} \hat{F}(\Bk, \omega) d^N \Bk d\omega \\
\hat{F}(\Bk, \omega) &= \inv{\lr{\sqrt{2 \pi}}^{N+1}} \int e^{-j \Bk \cdot \Bx – j \omega t} F(\Bx, t) d^N \Bx dt.
\end{aligned}
\end{equation}
We can now transform \ref{eqn:waveEquationGreens:860}, expressing \(f, g\) in terms of their transforms
\begin{equation}\label{eqn:waveEquationGreens:900}
\lr{ \lr{ j \Bk}^2 – \lr{ j \omega }^2/c^2 } \hat{f} = \hat{g},
\end{equation}
or
\begin{equation}\label{eqn:waveEquationGreens:920}
\hat{f} = \frac{\hat{g}}{(\omega/c)^2 – \Bk^2},
\end{equation}
or
\begin{equation}\label{eqn:waveEquationGreens:940}
\begin{aligned}
f(\Bx, t)
&= \inv{\lr{\sqrt{2 \pi}}^{N+1}} \int e^{j \Bk \cdot \Bx + j \omega t} \frac{\hat{g}(\Bk, \omega)}{(\omega/c)^2 – \Bk^2} d^N \Bk d\omega \\
&= \inv{\lr{2 \pi}^{N+1}} \int e^{j \Bk \cdot \Bx + j \omega t} \frac{g(\Bx’, t’)}{(\omega/c)^2 – \Bk^2} d^N \Bk d\omega e^{-j \Bk \cdot \Bx’ – j \omega t’} d^N \Bx’ dt’ \\
&=
\int d^N \Bx’ dt’ g(\Bx’, t’) G(\Bx, \Bx’, t, t’),
\end{aligned}
\end{equation}
where
\begin{equation}\label{eqn:waveEquationGreens:960}
G(\Bx, \Bx’, t, t’)
=
\inv{\lr{2 \pi}^{N+1}} \int d^N \Bk d\omega \frac{e^{j \Bk \cdot (\Bx-\Bx’) + j \omega (t- t’)}}{(\omega/c)^2 – \Bk^2}.
\end{equation}

Evaluating the 1D Green’s function

For the 1D case we have
\begin{equation}\label{eqn:waveEquationGreens:980}
G(\Bx, \Bx’, t, t’)
=
\inv{\lr{2 \pi}^2} \int dk d\omega \frac{e^{j k (x-x’) + j \omega (t- t’)}}{(\omega/c)^2 – k^2}
\end{equation}
Let’s write \( u = x – x’ \), and \( \tau = t – t’ \), and displace the poles by an imaginary offset \( j \epsilon \)
\begin{equation}\label{eqn:waveEquationGreens:1000}
G_\epsilon(u, \tau)
=
-\inv{\lr{2 \pi}^2} \int dk d\omega \frac{e^{j k u + j \omega \tau }}{\lr{ k – \lr{ \omega/c + j \epsilon}}\lr{ k + \lr{ \omega/c + j \epsilon }}}.
\end{equation}

Let’s start by assuming that \( \epsilon > 0 \). When \( u > 0 \), we can use a upper half plane contour in the k-plane, enclosing \( \omega/c + j \epsilon \), to find
\begin{equation}\label{eqn:waveEquationGreens:1020}
\begin{aligned}
G_\epsilon(u, \tau)
&=
-\frac{2 \pi j}{\lr{2 \pi}^2} \int d\omega \evalbar{\frac{e^{j k u + j \omega \tau }}{k + \lr{ \omega/c + j \epsilon }}}{k = \omega/c + j \epsilon} \\
&=
\frac{1}{4 \pi j} \int d\omega \frac{e^{j \omega (\tau + u/c)}}{\omega/c + j \epsilon }.
\end{aligned}
\end{equation}
However, for \( u < 0 \) we need the lower half plane contour that encloses \( -\omega/c - j \epsilon \). Our residue calculation is \begin{equation}\label{eqn:waveEquationGreens:1040} \begin{aligned} G_\epsilon(u, \tau) &= -\frac{-2 \pi j}{\lr{2 \pi}^2} \int d\omega \evalbar{\frac{e^{j k u + j \omega \tau }}{k - \lr{ \omega/c + j \epsilon }}}{k = -\omega/c - j \epsilon} \\ &= \frac{1}{4 \pi j} \int d\omega \frac{e^{j \omega (\tau - u/c)}}{\omega/c + j \epsilon }. \end{aligned} \end{equation} Merging the two cases, we have \begin{equation}\label{eqn:waveEquationGreens:1060} \begin{aligned} G_\epsilon(u, \tau) &= \frac{1}{4 \pi j} \int d\omega \frac{e^{j \omega (\tau + \Abs{u}/c)}}{\omega/c + j \epsilon } \\ &= \frac{c}{4 \pi j} \int d\omega \frac{e^{j \omega (\tau + \Abs{u}/c)}}{\omega + j \epsilon c } \\ \end{aligned} \end{equation} This can be integrated in the \(\omega\)-plane, with the pole at \( -j \epsilon c \). For \( \tau + \Abs{u}/c > 0 \), we need an upper half plane infinite semicircular contour, but have no enclosed pole. For \( \tau + \Abs{u}/c < 0 \), we have \begin{equation}\label{eqn:waveEquationGreens:1080} \begin{aligned} G_\epsilon(u, \tau) &= \frac{c (-2 \pi j)}{4 \pi j} \evalbar{ e^{j \omega (\tau + \Abs{u}/c)}}{\omega = -j \epsilon c} \\ &= -\frac{c}{2}, \end{aligned} \end{equation} (in the limit.) Putting both pieces together, we have found the advanced Green's function for the 1D wave equation \begin{equation}\label{eqn:waveEquationGreens:1100} \boxed{ G(u, \tau) = -\frac{c}{2} \Theta(-\tau - \Abs{u}/c). } \end{equation} Having found the advanced solution with a positive pole displacement, it is reasonable to assume that we will get the retarded solution, with a negative pole displacement \( \epsilon < 0 \). This time, the upper half plane infinite semicircular contour encloses the \( -\omega/c -j \epsilon \) pole, and the lower half plane contour encloses the \( \omega/c + j \epsilon \) pole. This gives, us, for \( u > 0 \)
\begin{equation}\label{eqn:waveEquationGreens:1120}
\begin{aligned}
G_\epsilon(u, \tau)
&=
-\frac{2 \pi j}{\lr{2 \pi}^2} \int d\omega \evalbar{\frac{e^{j k u + j \omega \tau }}{k – \lr{ \omega/c + j \epsilon }}}{k = -\omega/c – j \epsilon} \\
&=
-\frac{1}{4 \pi j} \int d\omega \frac{e^{j \omega (\tau – u/c)}}{\omega/c + j \epsilon } \\
&=
-\frac{c}{4 \pi j} \int d\omega \frac{e^{j \omega (\tau – u/c)}}{\omega – (-j \epsilon c) } \\
&=
-\frac{2 \pi j c}{4 \pi j} \evalbar{e^{j \omega (\tau – u/c)}}{\omega = -j \epsilon c } \Theta(\tau – u/c) \\
&=
-\frac{c}{2} \Theta(\tau – u/c),
\end{aligned}
\end{equation}
and for \( u < 0 \) \begin{equation}\label{eqn:waveEquationGreens:1140} \begin{aligned} G_\epsilon(u, \tau) &= -\frac{-2 \pi j}{\lr{2 \pi}^2} \int d\omega \evalbar{\frac{e^{j k u + j \omega \tau }}{k + \lr{ \omega/c + j \epsilon }}}{k = \omega/c + j \epsilon} \\ &= \frac{1}{4 \pi j} \int d\omega \frac{e^{j \omega (\tau + u/c)}}{\omega/c + j \epsilon } \\ &= \frac{c}{4 \pi j} \int d\omega \frac{e^{j \omega (\tau + u/c)}}{\omega - (-j \epsilon c) } \\ &= -\frac{-2 \pi j c}{4 \pi j} \evalbar{e^{j \omega (\tau + u/c)}}{\omega = -j \epsilon c } \Theta(\tau + u/c) \\ &= -\frac{c}{2} \Theta(\tau + u/c). \end{aligned} \end{equation} Combining the two cases, we've found the retarded solution \begin{equation}\label{eqn:waveEquationGreens:1160} \boxed{ G(u, \tau) = -\frac{c}{2} \Theta(\tau - \Abs{u}/c). } \end{equation} This matches Grok's claim (which we also verified.)

The convolution integrals.

Let’s write out the convolution integrals for fun. They are
\begin{equation}\label{eqn:waveEquationGreens:1180}
f(x,t) = -\frac{c}{2}
\int_{-\infty}^\infty dt’
\int_{-\infty}^\infty dx’
\Theta(\pm(t – t’) – \Abs{x – x’}/c) g(x’, t’).
\end{equation}

For the retarded case, we need only evaluate the step over the region
\begin{equation}\label{eqn:waveEquationGreens:1220}
t – t’ – \Abs{x – x’}/c > 0,
\end{equation}
or
\begin{equation}\label{eqn:waveEquationGreens:1240}
t – \Abs{x – x’}/c > t’.
\end{equation}
For the advanced case, we want the restriction
\begin{equation}\label{eqn:waveEquationGreens:1260}
-t + t’ – \Abs{x – x’}/c > 0,
\end{equation}
or
\begin{equation}\label{eqn:waveEquationGreens:1280}
t’ > t + \Abs{x – x’}/c,
\end{equation}
so the retarded convolution is
\begin{equation}\label{eqn:waveEquationGreens:1300}
f(x,t) = -\frac{c}{2}
\int_{-\infty}^\infty dx’
\int_{-\infty}^{t – \Abs{x – x’}/c} dt’
g(x’, t’),
\end{equation}
and the advanced convolution is
\begin{equation}\label{eqn:waveEquationGreens:1320}
f(x,t) = -\frac{c}{2}
\int_{-\infty}^\infty dx’
\int_{t + \Abs{x – x’}/c}^\infty dt’
g(x’, t’).
\end{equation}

TODO.

Next up will be an attempt to find the 2D Green’s function, and then for good measure, we should try to find the 3D case ourselves.

Correcting the errors: Green’s functions for the 1D Helmholtz and Laplacian operators.

September 28, 2025 math and physics play , , , , , , , , , , , , , , , ,

[Click here for a PDF version of this post, and others in this series]

The following recent posts explored 1D Green’s functions for the Helmholtz and Laplacian operators.  There was a sign error (wrong residue sign for a negatively oriented contour) that I made near the beginning that caused a lot of trouble.  Having found the error, I’ve now reworked all that exploratory content into a more coherent form.  That reworked content can be found in today’s blog post below (or in the PDF above, which includes all of this, plus the 2D and 3D derivations.)

Part 1. Green’s functions for the Helmholtz (wave equation) operator in various dimensions.

 

A trilogy in five+ parts: Confirming an error in the derived 1D Helmholtz Green’s function.

 

A trilogy in six+ parts: 1D Laplacian Green’s function

A trilogy in 7+ parts: A better check of the 1D Helmholtz Green’s function.

 

 

Helmholtz Green’s function in 1D.

Evaluating the integral.

For the one dimensional case, we want to evaluate
\begin{equation}\label{eqn:helmholtzGreensV2:200}
G(r) = -\inv{2 \pi} \int \inv{p^2 – k^2} e^{j p r} dp,
\end{equation}
an integral which is unfortunately non-convergent. Since we are dealing with delta functions, it is not surprising that we have convergence problems. The technique used in the book is to displace the pole slightly by a small imaginary amount, and then take the limit.

That is
\begin{equation}\label{eqn:helmholtzGreensV2:220}
G(r) = \lim_{\epsilon \rightarrow 0} G_\epsilon(r),
\end{equation}
where
\begin{equation}\label{eqn:helmholtzGreensV2:240}
\begin{aligned}
G_\epsilon(r)
&= -\inv{2 \pi} \int_{-\infty}^\infty \inv{p^2 – \lr{k + j \epsilon}^2} e^{j p r} dp \\
&= -\inv{2 \pi} \int_{-\infty}^\infty \frac{e^{j p r}}{\lr{ p – k -j \epsilon}\lr{p + k + j \epsilon}} dp.
\end{aligned}
\end{equation}
Our contours, for \( \epsilon > 0 \), are illustrated in fig 1.

fig 1. Contours for 3D Green’s function evaluation

For \( r > 0 \) we can use an upper half plane infinite semicircular contour integral, with \( R \rightarrow \infty \).

The residue calculation for this contour gives
\begin{equation}\label{eqn:helmholtzGreensV2:260}
\begin{aligned}
G_\epsilon(r)
&= -\frac{(+2 \pi j)}{2 \pi} \evalbar{\frac{e^{j p r}}{p + k + j \epsilon} }{p = k + j \epsilon} \\
&= -j \frac{e^{j \lr{k + j \epsilon} r}}{2\lr{k + j \epsilon}} \\
&= -j \frac{e^{j k r} e^{-\epsilon r}}{2\lr{k + j \epsilon}} \\
&\rightarrow -\frac{j}{2k} e^{j k r}.
\end{aligned}
\end{equation}

For \( r < 0 \) we use the lower half plane infinite semicircular contour For this contour, we find \begin{equation}\label{eqn:helmholtzGreensV2:280} \begin{aligned} G_\epsilon(r) &= -\frac{2 \pi j}{(-2 \pi)} \evalbar{\frac{e^{j p r}}{p – k – j \epsilon} }{p = -k – j \epsilon} \\ &= -j \frac{e^{-j \lr{k + j \epsilon} r}}{2\lr{k + j \epsilon}} \\ &= -j \frac{e^{-j k r} e^{\epsilon r}}{2\lr{k + j \epsilon}} \\ &\rightarrow -\frac{j}{2k} e^{-j k r}. \end{aligned} \end{equation} Combining both results, our Green’s function, after a positive pole displacement \( \epsilon > 0 \), is
\begin{equation}\label{eqn:helmholtzGreensV2:300}
G(r) = \frac{1}{2 j k} e^{j k \Abs{r}}.
\end{equation}

Similarly, should we pick \( \epsilon < 0 \), the same sort of calculation yields an incoming wave solution \begin{equation}\label{eqn:helmholtzGreensV2:2240} G(r) = -\frac{1}{2 j k} e^{-j k \Abs{r}}. \end{equation} Allowing for either, we have Green’s functions for both the incoming and outgoing wave cases \begin{equation}\label{eqn:helmholtzGreensV2:2260} \boxed{ G_{\pm}(x – x’) = \pm \frac{1}{2 j k} e^{ \pm j k \Abs{x – x’}}. } \end{equation} With two Green’s functions, we can also make a linear combination. Specifically \begin{equation}\label{eqn:helmholtzGreensV2:2280} \begin{aligned} G(r) &= \inv{2}\lr{ G_{+}(r) + G_{-}(r) } \\ &= \inv{4 j k}\lr{ e^{ j k \Abs{r}} – e^{ – j k \Abs{r}} } \\ &= \inv{2 k} \sin\lr{ k \Abs{r} } \end{aligned} \end{equation} This real valued Green’s function is plotted in fig. 2.

fig. 2. Green’s function for 1D Helmholtz operator.

The convolution is now fully specified, providing a specific solution to the non-homogeneous equation \begin{equation}\label{eqn:helmholtzGreensV2:400} U(x) = \inv{2k} \int_{-\infty}^\infty \sin( k\Abs{r} ) V(x + r) dr. \end{equation}

The general solution may also include any solutions to the homogeneous Helmholtz equation \begin{equation}\label{eqn:helmholtzGreensV2:2300} U(x) = A e^{j k x} + B e^{-j k x} + \inv{2k} \int_{-\infty}^\infty \sin( k\Abs{r} ) V(x + r) dr. \end{equation}

A strictly causal solution.

We can split the convolution kernel a “causal” part, where only the spatially-“past” values of \( V \) contribute, and an “acausal” part \begin{equation}\label{eqn:helmholtzGreensV2:2320} U(x) = \inv{2k} \int_{-\infty}^0 \sin( k\Abs{r} ) V(x + r) dr + \inv{2k} \int_{0}^\infty \sin( k\Abs{r} ) V(x + r) dr. \end{equation} In a sense, we are averaging causal and acausal portions of the convolution. Suppose that we form a convolution with a built in cut-off, so that values of \( V(x’), x’ > x \) do not contribute to \( U(x) \). That is
\begin{equation}\label{eqn:helmholtzGreensV2:440}
f(x) = \int_{-\infty}^x \frac{\sin\lr{k \lr{x – x’}}}{k} V(x’) dx’.
\end{equation}
Here the one-half factor has been dropped, since we are no longer performing a QFT like average of causal and acausal terms.

Intuition suggests this should be a solution to the Helmholtz equation, but let’s test that guess. We start with the identity
\begin{equation}\label{eqn:helmholtzGreensV2:460}
\frac{d}{dx} \int_a^x g(x, x’) dx’
=
\evalbar{g(x, x’) }{x’ = x} + \int_a^x \frac{\partial g(x,x’)}{dx} dx’.
\end{equation}
Taking the first derivative of \( f(x) \), we find
\begin{equation}\label{eqn:helmholtzGreensV2:480}
\begin{aligned}
\frac{df}{dx}
&= \evalbar{ \frac{\sin\lr{k \lr{x – x’}}}{k} V(x’) }{x’ = x} + k \int_{-\infty}^x \frac{\cos\lr{k \lr{x – x’}}}{k} V(x’) dx’ \\
&= k \int_{-\infty}^x \frac{\cos\lr{k \lr{x – x’}}}{k} V(x’) dx’,
\end{aligned}
\end{equation}
where we have, somewhat lazily, treated the infinite limit as a constant. Effectively, this requires that the forcing function \( V(x) \) is zero at \( -\infty \). Taking the second derivative, we have
\begin{equation}\label{eqn:helmholtzGreensV2:500}
\begin{aligned}
\frac{d^2f}{dx^2}
&=
\evalbar{ k \frac{\cos\lr{k \lr{x – x’}}}{k} V(x’) }{x’ = x} – k^2 \int_{-\infty}^x \frac{\sin\lr{k \lr{x – x’}}}{k} V(x’) dx’ \\
&= V(x) – k^2 f(x),
\end{aligned}
\end{equation}
or
\begin{equation}\label{eqn:helmholtzGreensV2:520}
\frac{d^2}{dx^2} f(x) + k^2 f(x) = V(x).
\end{equation}
This verifies that \ref{eqn:helmholtzGreensV2:440} is also a specific solution to the wave equation, as expected and desired.

It appears that the general solution is likely of the following form
\begin{equation}\label{eqn:helmholtzGreensV2:380b}
U(x) =
A’ \cos\lr{ k x } +
B’ \sin\lr{ k x }
+\alpha \int_{-\infty}^x \frac{\sin\lr{k \lr{x – x’}}}{k} V(x’) dx’
+(1-\alpha)\int_x^\infty \frac{\sin\lr{k \lr{x – x’}}}{k} V(x’) dx’,
\end{equation}
with \( \alpha \in [0,1] \).

It’s pretty cool that we can completely solve the 1D forced wave equation, for any forcing function, from first principles. Yes, I took liberties that would make a mathematician cringe, but we are telling a story, and leaving the footnotes to somebody else.

Verification of the 1D Helmholtz Green’s function.

Let’s show that the outgoing Green’s function has the desired delta function semantics. That is
\begin{equation}\label{eqn:helmholtzGreensV2:2060}
\lr{ \spacegrad^2 + k^2 } G(x, x’) = \delta(x – x’),
\end{equation}
where
\begin{equation}\label{eqn:helmholtzGreensV2:2080}
G(x, x’) = \frac{e^{j k \Abs{x – x’}}}{2 j k}.
\end{equation}
Making a \( r = x – x’ \) change of variables gives
\begin{equation}\label{eqn:helmholtzGreensV2:2120}
\spacegrad^2 e^{j k \Abs{x – x’}} = \frac{d^2}{dr^2} e^{j k \Abs{r}}
\end{equation}

The function \( \Abs{r} \) formally has no derivative at the origin, but we may use the physics trick, rewriting the absolute in terms of the Heaviside theta function
\begin{equation}\label{eqn:helmholtzGreensV2:2340}
\Abs{r} = r \Theta(r) – r \Theta(-r).
\end{equation}
We then use the delta function identification for the derivative
\begin{equation}\label{eqn:helmholtzGreensV2:2360}
\Theta'(r) = \delta(r).
\end{equation}
In particular
\begin{equation}\label{eqn:helmholtzGreensV2:2380}
\begin{aligned}
\frac{d}{dr} \Abs{r}
&= \Theta(r) – \Theta(-r) + r \delta(r) – (-1)\delta(-r) \\
&= \Theta(r) – \Theta(-r) + 2 r \delta(r),
\end{aligned}
\end{equation}
using the symmetric property of the delta function \( \delta(-r) = \delta(r) \). The delta function contribution to this derivative is actually zero, as seen when we operate with \( r \delta(r) \) against a test function
\begin{equation}\label{eqn:helmholtzGreensV2:2400}
\begin{aligned}
\int r \delta(r) f(r) dr
&=
\evalbar{r f(r)}{r = 0} \\
&= 0.
\end{aligned}
\end{equation}
We’ve now found that
\begin{equation}\label{eqn:helmholtzGreensV2:2420}
\frac{d}{dr} \Abs{r} = \Theta(r) – \Theta(-r) = \mathrm{sgn}(r),
\end{equation}
the sign function. The derivative of the sign function is
\begin{equation}\label{eqn:helmholtzGreensV2:2440}
\begin{aligned}
\mathrm{sgn}'(r)
&= \lr{ \Theta(r) – \Theta(-r) }’ \\
&= \delta(r) -(-1)\delta(-r) \\
&= 2 \delta(r).
\end{aligned}
\end{equation}

The derivatives are
\begin{equation}\label{eqn:helmholtzGreensV2:2140}
\begin{aligned}
\frac{d}{dr} e^{j k \Abs{r} }
&=
j k e^{j k \Abs{r} } \frac{d\Abs{r}}{dr} \\
&=
j k e^{j k \Abs{r} } \mathrm{sgn}(r).
\end{aligned}
\end{equation}
and
\begin{equation}\label{eqn:helmholtzGreensV2:2160}
\frac{d^2}{dr^2} e{j k \Abs{r}}
=
\lr{ j k \mathrm{sgn}(r) }^2 e^{j k \Abs{r} } + 2 j k e^{j k \Abs{r} } \delta(r).
\end{equation}

We can identify \( e^{j k \Abs{r} } \delta(r) = \delta(r) \), just as we identified \( r \delta(r) = 0 \), by application to a test function. That is
\begin{equation}\label{eqn:helmholtzGreensV2:2180}
\begin{aligned}
\int e^{j k \Abs{r} } \delta(r) f(r) dr
&=
\evalbar{e^{j k \Abs{r} } f(r)}{r = 0} \\
&=
f(0) \\
&=
\int \delta(r) f(r) dr.
\end{aligned}
\end{equation}
With that identification
\begin{equation}\label{eqn:helmholtzGreensV2:2200}
\spacegrad^2 e^{j k \Abs{r} } = -k^2 e^{j k \Abs{r} } + 2 j k \delta(r),
\end{equation}
or
\begin{equation}\label{eqn:helmholtzGreensV2:2220}
\boxed{
\lr{ \spacegrad^2 + k^2 } \frac{e^{j k \Abs{x – x’} }}{2 j k} = \delta(x – x’).
}
\end{equation}

Verifying the Green’s function with convolution.

Avoiding the physics tricks, we may use a limiting argument to validate our Green’s function.

We first want to show that at points \( x’ \ne x \) the Helmholtz operator applied to the Green’s function is zero:
\begin{equation}\label{eqn:helmholtzGreensV2:1220}
\lr{ \spacegrad^2 + k^2} G(x,x’) = 0,
\end{equation}
Since we are avoiding the origin
\begin{equation}\label{eqn:helmholtzGreensV2:1240}
\lr{ k^2 + \frac{d^2}{dx^2} } e^{j k \Abs{x – x’}}.
\end{equation}
We expect that this will be zero. Making a change of variables \( r = x’ – x \), we want to evaluate
\begin{equation}\label{eqn:helmholtzGreensV2:1260}
\lr{ k^2 + \frac{d^2}{dr^2} } e^{j k \Abs{r}},
\end{equation}
assuming that we are omitting a neighbourhood around \( r = 0 \) where the absolute value causes trouble. For \( r > 0 \)
\begin{equation}\label{eqn:helmholtzGreensV2:1280}
\begin{aligned}
\frac{d}{dr} e^{j k \Abs{r}}
&= \frac{d}{dr} e^{j k r} \\
&= j k e^{j k r},
\end{aligned}
\end{equation}
and
\begin{equation}\label{eqn:helmholtzGreensV2:1300}
\begin{aligned}
\frac{d^2}{dr^2} e^{j k \Abs{r}}
&= \frac{d}{dr} j k e^{j k r} \\
&= (j k)^2 e^{j k r}.
\end{aligned}
\end{equation}
Similarly, for \( r < 0 \), we have
\begin{equation}\label{eqn:helmholtzGreensV2:1320}
\frac{d^2}{dr^2} e^{j k \Abs{r}} = (-jk)^2 e^{j k \Abs{r}}.
\end{equation}
In both cases, provided we are in a neighbourhood that omits \( r \ne 0 \), we have
\begin{equation}\label{eqn:helmholtzGreensV2:1340}
\lr{ k^2 + \frac{d^2}{dr^2} } e^{j k \Abs{r}} = 0,
\end{equation}
as desired.

The takeaway is that we have
\begin{equation}\label{eqn:helmholtzGreensV2:1360}
\begin{aligned}
\lr{ \spacegrad^2 + k^2 }\int_{-\infty}^\infty G(x, x’) V(x’) dx’
&=
\int_{\Abs{x’ – x} \le \epsilon} V(x’) \lr{ \frac{d^2}{d{x’}^2} + k^2 } G(x, x’) dx’ \\
&=
\int_{-\epsilon}^\epsilon V(x + r) \lr{ \frac{d^2}{dr^2} + k^2 } G(r) dr \\
\end{aligned}
\end{equation}
for some arbitrarily small value of \( \epsilon \). Observe that after bringing the operator into the integral, we also made a change of variables, first to \( x’ \) for the Laplacian, and then to \( r = x’ – x \).

We’d like the 1D equivalent of Green’s theorem to reduce this, so let’s work that out first.
\begin{equation}\label{eqn:helmholtzGreensV2:1380}
\begin{aligned}
\int dx\, v \frac{d^2 u}{dx^2} – \int dx\, u \frac{d^2 v}{dx^2}
&=
\int dx\,
\lr{
\frac{d}{dx} \lr{
v \frac{du}{dx}
}
– \frac{dv}{dx} \frac{du}{dx}
}
–
\int dx\,
\lr{
\frac{d}{dx} \lr{
u \frac{dv}{dx}
}
–
\frac{du}{dx} \frac{dv}{dx}
}
\\
&=
\int dx\,
\frac{d}{dx}
\lr{
v \frac{du}{dx}
}
–
\int dx\,
\frac{d}{dx} \lr{
u \frac{dv}{dx}
}
\\
&=
v \frac{du}{dx}
–
u \frac{dv}{dx},
\end{aligned}
\end{equation}
so
\begin{equation}\label{eqn:helmholtzGreensV2:1400}
\boxed{
\int_a^b dx\, v \frac{d^2 u}{dx^2}
=
\int_a^b dx\, u \frac{d^2 v}{dx^2}
+
\evalrange{v \frac{du}{dx}}{a}{b}
–
\evalrange{u \frac{dv}{dx}}{a}{b}.
}
\end{equation}

This gives us
\begin{equation}\label{eqn:helmholtzGreensV2:1420}
\begin{aligned}
\lr{ \spacegrad^2 + k^2 }\int_{\Abs{x’- x} \le \epsilon} e^{j k \Abs{x – x’}} V(x’) dx’
&=
\int_{-\epsilon}^\epsilon e^{j k \Abs{r}} \lr{ k^2 + \frac{d^2}{dr^2} } V(x + r) dr \\
&\quad +
\evalrange{ V(x + r) \frac{d}{dr} e^{j k \Abs{r}} }{-\epsilon}{\epsilon}
–
\evalrange{ e^{j k \Abs{r}} \frac{d}{dr} V(x+r) }{-\epsilon}{\epsilon}.
\end{aligned}
\end{equation}
If we can assume that \( V \) and it’s first and second derivatives are all continuous over this small interval, then the first integral is approximately
\begin{equation}\label{eqn:helmholtzGreensV2:1440}
\begin{aligned}
\int_{-\epsilon}^\epsilon e^{j k \Abs{r}} \lr{ k^2 + \frac{d^2}{dr^2} } V(x + r) dr
&\sim
\lr{ \lr{ k^2 + \frac{d^2}{dr^2} } V(x + r) }
\int_{-\epsilon}^\epsilon e^{j k \Abs{r}} dr \\
&=
\frac{2 j}{k} \lr{ 1 – e^{ j k \epsilon} }
\lr{ \lr{ k^2 + \frac{d^2}{dr^2} } V(x + r) } \\
&\rightarrow 0.
\end{aligned}
\end{equation}
Similarly, with \( dV/dr \) continuity condition, that last term is also zero. We are left, for \( \epsilon \) sufficiently small, we are left with
\begin{equation}\label{eqn:helmholtzGreensV2:1460}
\lr{ \spacegrad^2 + k^2 }\int_{\Abs{x’- x} \le \epsilon} e^{j k \Abs{x – x’}} V(x’) dx’
=
V(x)
\evalrange{ \frac{d}{dr} e^{j k \Abs{r}} }{-\epsilon}{\epsilon}.
\end{equation}
Because we are evaluating this derivative only at points \( r = \pm \epsilon \ne 0 \), that derivative is
\begin{equation}\label{eqn:helmholtzGreensV2:2460}
\frac{d}{dr} e^{j k \Abs{r}}
=
j k e^{j k \Abs{r}} \mathrm{sgn}(r),
\end{equation}
leaving us with
\begin{equation}\label{eqn:helmholtzGreensV2:2480}
\begin{aligned}
\evalrange{ \frac{d}{dr} e^{j k \Abs{r}} }{-\epsilon}{\epsilon}
&=
j k e^{j k \Abs{\epsilon}} \mathrm{sgn}(\epsilon) – j k e^{j k \Abs{-\epsilon}} \mathrm{sgn}(-\epsilon) \\
&=
j k e^{j k \Abs{\epsilon}} \lr{ 1 – (-1) } \\
2 j k e^{j k \Abs{\epsilon}}.
\end{aligned}
\end{equation}
or
\begin{equation}\label{eqn:helmholtzGreensV2:2500}
\lr{ \spacegrad^2 + k^2 }\int_{-\infty}^\infty e^{j k \Abs{x – x’} } V(x’) dx’
=
2 j k e^{j k \Abs{\epsilon}} V(x).
\end{equation}
Taking limits and dividing through by \( 2 j k \) proves the result.

1D Laplacian Green’s function.

Having blundered our way to what appears to be the correct Green’s function for the 1D Helmholtz operator, let’s further validate that by deriving the Green’s function for the 1D Laplacian. We should also be able to verify that it has the correct delta function semantics.

Expansion in series and taking the limit.

Expanding the Helmholtz Green’s function in series around \( k \Abs{r} \) we have
\begin{equation}\label{eqn:helmholtzGreensV2:1780}
\begin{aligned}
G(r)
&= -\frac{j}{2k} \lr{ 1 + j k \Abs{r} + O((k \Abs{r})^2) } \\
&= -\frac{j}{2} \lr{ \inv{k} + j \Abs{r} + \inv{k} O((k \Abs{r})^2) } \\
\end{aligned}
\end{equation}
This means that to first order in \( k \), we have
\begin{equation}\label{eqn:helmholtzGreensV2:1800}
G(r) + \frac{j}{2k} = \frac{\Abs{r}}{2}.
\end{equation}
As before, we are free to add constant terms to the Green’s function for the Laplacian, and we conclude that the 1D Green’s function for the Laplacian is
\begin{equation}\label{eqn:helmholtzGreensV2:1820}
\boxed{
G(r) = \frac{\Abs{r}}{2}.
}
\end{equation}

Alternatively, we may use the real sine form of the Green’s function, which has a nice expansion around \( k = 0 \), and arrive at the same result.

Observe that
\begin{equation}\label{eqn:helmholtzGreensV2:2520}
\frac{d^2}{dr^2} \frac{\Abs{r}}{2} =
\frac{d}{dr} \frac{\mathrm{sgn}(r)}{2} =
\delta(r),
\end{equation}
which verifies that this is a valid Green’s function for the 1D Laplacian.

Verifying the Green’s function with convolution.

We can also operate on the convolution with the Laplacian, to verify correctness. We are interested in evaluating
\begin{equation}\label{eqn:helmholtzGreensV2:1840}
\spacegrad^2 \int \frac{\Abs{x – x’}}{2} V(x’) dx’ = \int V(x + r) \frac{d^2}{dr^2} \frac{\Abs{r}}{2} dr.
\end{equation}
If all goes well, this should evaluate to \( V(x) \), indicating that \( \spacegrad^2 \Abs{x – x’}/2 = \delta(x – x’) \). As a first step, we expect \( \spacegrad^2 G = 0 \), for \( x \ne x’ \). Consider first \( r > 0 \), where
\begin{equation}\label{eqn:helmholtzGreensV2:1860}
\frac{d}{dr} \Abs{r}
=
\frac{d}{dr} r
= 1,
\end{equation}
and for \( r < 0 \) where
\begin{equation}\label{eqn:helmholtzGreensV2:1880}
\frac{d}{dr} \Abs{r}
=
\frac{d}{dr} (-r)
= -1.
\end{equation}
This means that, away from the origin \( d\Abs{r}/dr = \mathrm{sgn}(r) \), and \( d^2 \Abs{r}/dr^2 = 0\). We can conclude that, for some non-zero positive epsilon that we will eventually let approach zero, we have
\begin{equation}\label{eqn:helmholtzGreensV2:1900}
\begin{aligned}
\spacegrad^2 \int \frac{\Abs{x – x’}}{2} V(x’) dx’
&= \int_{-\epsilon}^\epsilon V(x + r) \frac{d^2}{dr^2} \frac{\Abs{r}}{2} dr \\
&= \int_{-\epsilon}^\epsilon \lr{
\frac{d}{dr} \lr{ V(x + r) \frac{d}{dr} \frac{\Abs{r}}{2} }
– \frac{dV(x + r)}{dr} \frac{d}{dr} \frac{\Abs{r}}{2}
}
dr \\
&=
\evalrange{ V(x + r) \frac{d}{dr} \frac{\Abs{r}}{2} }{-\epsilon}{\epsilon}
– \int_{-\epsilon}^\epsilon
\lr{
\frac{d}{dr} \lr{ \frac{dV(x + r)}{dr} \frac{\Abs{r}}{2} }
–
\frac{d^2V(x + r)}{dr^2} \frac{\Abs{r}}{2}
} dr \\
&=
\evalrange{ V(x + r) \frac{d}{dr} \frac{\Abs{r}}{2} }{-\epsilon}{\epsilon}
-\evalrange{ \frac{dV(x + r)}{dr} \frac{\Abs{r}}{2} }{-\epsilon}{\epsilon}
+
\int_{-\epsilon}^\epsilon \frac{d^2V(x + r)}{dr^2} \frac{\Abs{r}}{2} dr \\
&=
\inv{2} \lr{ V(x + \epsilon) + V(x – \epsilon) } \\
&-\quad \frac{\epsilon}{2}\lr{
\frac{dV(x + \epsilon)}{dr}
–
\frac{dV(x – \epsilon)}{dr}
} \\
&+
\frac{\epsilon^2}{2} \lr{
\frac{d^2V(x + \epsilon)}{dr^2}
+
\frac{d^2V(x + \epsilon)}{dr^2}
}.
\end{aligned}
\end{equation}
In the limit we have
\begin{equation}\label{eqn:helmholtzGreensV2:1920}
\boxed{
\spacegrad^2 \int G(x, x’) V(x’) dx’ = \inv{2} \lr{ V(x^+) + V(x^-) }.
}
\end{equation}

If the test (or driving) function is continuous at \( x’ = x \), then this is exactly the delta-function semantics that we expect of a Green’s function. It’s interesting that this check provides us with precise semantics for the Green’s function for discontinuous functions too.