Fundamental Solutions and Green's Functions for PDEs

A linear differential operator $\mathcal{L}$ turns an unknown field into a known source: $\mathcal{L}u=f$ . Solving it looks like an infinite problem — one equation for every point of the domain at once — but linearity offers a shortcut so powerful it organizes half of mathematical physics. Chop the source into an infinite dust of point sources, $f(\mathbf{x})=\int f(\mathbf{x}')\,\delta(\mathbf{x}-\mathbf{x}')\,\mathrm{d}\mathbf{x}'$ , solve the equation once for a single unit spike, and then add up the responses. That one universal response — the field a point source makes — is the Green's function $G$ , and the whole inhomogeneous problem collapses to a single integral $u=\int G f$ . This article is the hinge of the PDE chapter. Given $\mathcal{L}$ and a domain $\Omega$ , it asks the operational question: what is the inverse kernel, and how do the boundaries enter? Two ideas do all the work. First, one must separate the fundamental solution — the free-space response fixed only by decay or radiation at infinity — from the Green's function proper, which additionally obeys the homogeneous conditions on $\partial\Omega$ ; the two differ by a solution of $\mathcal{L}v=0$ , and steering that freedom is exactly the method of images and the eigenfunction expansion. Second, Green's identities fold the boundary data into a surface integral over the same kernel, and for self-adjoint $\mathcal{L}$ force the reciprocity $G(\mathbf{x},\mathbf{x}')=G(\mathbf{x}',\mathbf{x})$ . Each PDE type then stamps its physics onto the analytic shape of $G$ : the smooth $1/r$ potential of the ellipse, the spreading Gaussian of the parabola, the razor-sharp light-cone shell of the hyperbola.

Key takeaways

A Green's function is the response of a linear operator to a unit point source, $\mathcal{L}_{\mathbf{x}}G(\mathbf{x},\mathbf{x}')=\delta(\mathbf{x}-\mathbf{x}')$ ; every inhomogeneous problem $\mathcal{L}u=f$ is then solved by superposition, $u(\mathbf{x})=\int G(\mathbf{x},\mathbf{x}')f(\mathbf{x}')\,\mathrm{d}\mathbf{x}'$ .
A fundamental solution is the free-space kernel (no boundaries); the Green's function proper adds a homogeneous solution so it obeys the boundary conditions — the content of the method of images.
Green's identities convert the volume PDE into a boundary integral: interior value = source term weighted by $G$ plus boundary data weighted by $G$ and $\partial_n G$ .
For self-adjoint $\mathcal{L}$ the kernel is symmetric, $G(\mathbf{x},\mathbf{x}')=G(\mathbf{x}',\mathbf{x})$ (reciprocity), and can be assembled from eigenfunctions as $G=\sum_n \phi_n(\mathbf{x})\phi_n^*(\mathbf{x}')/\lambda_n$ .
Each canonical PDE has its own kernel: $-1/(4\pi r)$ for Poisson in 3D, the Gaussian $(4\pi\kappa t)^{-n/2}e^{-r^2/4\kappa t}$ for diffusion, and the retarded shell $\delta(t-r/c)/(4\pi r)$ for the wave equation.
The contour or $i\epsilon$ prescription that inverts $\mathcal{L}$ in Fourier space selects retarded, advanced, or Feynman kernels — physically distinct solutions of one equation.
Key idea: Solve $\mathcal{L}G=\delta$ once; then $u=\int G f$ solves $\mathcal{L}u=f$ by superposition, and Green's identities fold the boundary data into a surface integral over the same kernel.

Where this lives today

Green's functions are the computational engine of much of modern physics. The boundary-element method (BEM), and its electromagnetic cousin the method of moments, discretize the representation (4) directly, solving for currents and potentials on surfaces alone rather than meshing all of space — the workhorse of antenna design and acoustic scattering. In nanoscale device physics the non-equilibrium Green's function (NEGF) formalism computes electron transport through a transistor channel, while the many-body Green's function underlies the GW approximation for the band gaps of real materials, and lattice QCD measures quark propagators — Green's functions of the Dirac operator — to extract hadron masses. Perhaps the most striking modern use is seismic interferometry: cross-correlating the ambient noise at two sensors recovers the Green's function between them, exactly as reciprocity predicts, so seismologists now image the crust and monitor volcanoes using only the background rumble of the oceans — turning noise into a source. And every gravitational-wave template LIGO matches is, at bottom, the retarded Green's function of the linearized Einstein equations for an inspiralling binary.

The physical picture

Think of $\mathcal{L}$ as an infinite-dimensional matrix acting on the vector of field values. The delta function is the "identity column" — a source concentrated at one coordinate $\mathbf{x}'$ — and the Green's function is literally the corresponding column of the inverse, $G=\mathcal{L}^{-1}$ , so $u=Gf$ is the continuum version of $u=A^{-1}b$ . A unit spike at $\mathbf{x}'$ radiates a disturbance, and $G(\mathbf{x},\mathbf{x}')$ records how much of it arrives at $\mathbf{x}$ . Swap emitter and detector and the reading is unchanged — reciprocity, the same statement behind Helmholtz reciprocity in acoustics and optics and the symmetry of a resistance matrix. The first figure makes "one source, one response" concrete, with the numbers of the 3D Coulomb kernel.

Figure

A unit δ-source at x′ makes the free-space response G(x,x′) at the observation point x. In 3D this is the Coulomb/Newton kernel G = −1/(4π r), falling as 1/r: −0.0796 at r = 1, −0.0398 at r = 2, −0.0265 at r = 3, and reciprocity means G(x,x′) = G(x′,x)

The analytic shape of $G$ reads off the physics, and this is where the three families of second-order PDE part ways. A $1/r$ tail with no time in it means an elliptic operator: an instantaneous, action-at-a-distance potential, smooth everywhere but the source. A spreading Gaussian means a parabolic operator: diffusive, time-irreversible smearing, nonzero everywhere the instant the source fires. A delta-shell riding $r=ct$ means a hyperbolic operator: a signal arriving sharply on the light cone, and not before. In Fourier space the story is algebraic, $\tilde{G}(\mathbf{k})=1/\tilde{\mathcal{L}}(\mathbf{k})$ , with the poles on the dispersion relation — the way one slips around them is exactly the choice between retarded, advanced, and Feynman. The second figure lines the three characters up.

Figure

The three characters of G by PDE type. Elliptic Poisson: a smooth 1/r potential. Parabolic heat: a Gaussian of width √(2κt) and peak (4πκt)^(−n/2) that spreads and never sharpens. Hyperbolic wave: a razor-thin shell δ(t − r/c) riding the light cone r = ct

The core equations

The whole method rests on a single defining equation: the Green's function is the field of a unit point source. Read $\mathcal{L}_{\mathbf{x}}$ as "the operator, differentiating in the unprimed variable $\mathbf{x}$ ", while $\mathbf{x}'$ marks where the source sits.

\mathcal{L}_{\mathbf{x}}\,G(\mathbf{x},\mathbf{x}')=\delta(\mathbf{x}-\mathbf{x}') \qquad (1)

where: $\mathcal{L}_{\mathbf{x}}$ — linear differential operator acting on the $\mathbf{x}$ variable; $G(\mathbf{x},\mathbf{x}')$ — Green's function, the response at $\mathbf{x}$ to a point source at $\mathbf{x}'$ ; $\delta(\mathbf{x}-\mathbf{x}')$ — Dirac delta, the unit point source. Equation (1) alone does not fix $G$ : one must also state the boundary or initial conditions. Once $G$ is known, every inhomogeneous problem with the same operator and conditions is solved by one integral. This is the payoff — a universal kernel, reused for any source.

u(\mathbf{x})=\int_\Omega G(\mathbf{x},\mathbf{x}')\,f(\mathbf{x}')\,\mathrm{d}^n x' \qquad (2)

where: $u$ — solution of $\mathcal{L}u=f$ ; $f$ — source or inhomogeneity; $\Omega$ — domain; $\mathrm{d}^n x'$ — volume element in the primed (source) variable. The formula says: add up the response to each point of the source, weighted by how strong the source is there. The bridge between the pointwise PDE and boundary data is Green's second identity — integration by parts done twice for the Laplacian. It trades a volume integral of $u\nabla^2 v-v\nabla^2 u$ for a pure surface integral.

\int_\Omega\!\big(u\,\nabla^2 v-v\,\nabla^2 u\big)\,\mathrm{d}V=\oint_{\partial\Omega}\!\big(u\,\partial_n v-v\,\partial_n u\big)\,\mathrm{d}S \qquad (3)

where: $u,v$ — any two twice-differentiable scalar fields; $\partial_n=\hat{\mathbf{n}}\cdot\nabla$ — outward normal derivative on the boundary $\partial\Omega$ ; $\mathrm{d}S$ — surface element. Setting $v=G$ is the master move that produces the representation formula. Doing exactly that yields the integral representation of the interior value: a volume integral of the source plus a surface integral carrying the boundary data. This is the equation that boundary-element methods discretize.

u(\mathbf{x})=\int_\Omega G\,f\,\mathrm{d}V'+\oint_{\partial\Omega}\!\big(u\,\partial_{n'}G-G\,\partial_{n'}u\big)\,\mathrm{d}S' \qquad (4)

where: $\nabla^2 u=f$ inside $\Omega$ ; $G$ — free-space or boundary Green's function, $\nabla^2 G=\delta$ . Choosing $G=0$ on $\partial\Omega$ (a Dirichlet $G$ ) kills the $G\,\partial_{n'}u$ term, so only the prescribed boundary values of $u$ are needed — not its unknown normal slope. For a formally self-adjoint operator with homogeneous conditions, the surface terms vanish and the kernel is forced symmetric. This reciprocity is a theorem, not a coincidence.

G(\mathbf{x},\mathbf{x}')=G(\mathbf{x}',\mathbf{x}) \qquad (5)

where: $\mathcal{L}$ — formally self-adjoint (e.g. Sturm-Liouville or the Laplacian) with boundary conditions that make $\oint(G_1\partial_n G_2-G_2\partial_n G_1)\,\mathrm{d}S=0$ . Physically: source and detector are interchangeable. When the eigenfunctions of $\mathcal{L}$ are known, $G$ can be built directly by inverting the operator eigenvalue-by-eigenvalue — the spectral representation. It is the constructive face of " $G=\mathcal{L}^{-1}$ ".

G(\mathbf{x},\mathbf{x}')=\sum_n \frac{\phi_n(\mathbf{x})\,\phi_n^*(\mathbf{x}')}{\lambda_n} \qquad (6)

where: $\phi_n$ — orthonormal eigenfunctions of $\mathcal{L}$ obeying the boundary conditions; $\lambda_n$ — eigenvalue, $\mathcal{L}\phi_n=\lambda_n\phi_n$ . The sum diverges if any $\lambda_n=0$ (a zero mode), which signals the need for a modified Green's function and a solvability condition. Two further constructions recur below rather than as numbered cards. In Fourier space the inversion is algebraic, $\tilde G=1/[k^2-(\omega+i\epsilon)^2/c^2]$ : the ambiguity of $G$ becomes poles at $\omega=\pm ck$ , and the $i\epsilon$ shift chooses retarded, advanced, or Feynman (derivation 3). In simple geometries a boundary is met by an image source placed outside $\Omega$ ; for a grounded plane one mirror charge gives $G_{\mathrm{Dir}}=-\tfrac{1}{4\pi}\big(1/|\mathbf{x}-\mathbf{x}'|-1/|\mathbf{x}-\mathbf{x}'^{*}|\big)$ , which vanishes on the plane automatically (worked example 1).

Derivation AFundamental solution of the 3D Laplacian

The goal is to prove that the free-space kernel of Poisson's operator in three dimensions is $G=-1/(4\pi r)$ — the Coulomb and Newton potential — from only the defining equation $\nabla^2 G=\delta(\mathbf{x})$ and the symmetry of empty space. The strategy: use translational and rotational invariance to reduce the PDE to a radial ODE valid away from the origin, solve it, and fix the one undetermined constant not by a boundary condition but by the source strength, extracted with Gauss's theorem on a shrinking ball around the spike.

Step 1. Because free space has no preferred point or direction, the response to a source at the origin can depend only on $r=|\mathbf{x}|$ . For $r\neq0$ the source vanishes, so $G$ satisfies Laplace's equation, which in spherical coordinates is purely radial: $\dfrac{1}{r^2}\dfrac{\mathrm{d}}{\mathrm{d}r}\!\left(r^2\dfrac{\mathrm{d}G}{\mathrm{d}r}\right)=0$ for $r>0$ . Why: symmetry removes the angular derivatives, turning a PDE in three variables into one ODE in $r$ .

Step 2. Integrate once: $r^2\,\mathrm{d}G/\mathrm{d}r=\text{const}$ , hence $G=A/r+B$ . Demanding $G\to0$ at infinity kills $B$ , leaving $G=A/r$ . Why: one integration constant is fixed by behaviour at infinity, the other will be fixed by the source strength.

Step 3. To pin $A$ , integrate the defining equation over a small ball $V_\epsilon$ and convert to a flux: $\displaystyle\int_{V_\epsilon}\nabla^2 G\,\mathrm{d}V=\oint_{S_\epsilon}\nabla G\cdot\mathrm{d}\mathbf{S}=\int_{V_\epsilon}\delta\,\mathrm{d}V=1$ . Why: the delta is invisible to the pointwise ODE, but its unit strength is captured by integrating across it.

Step 4. Evaluate the flux with $\nabla G=-\dfrac{A}{r^2}\hat{\mathbf{r}}$ on the sphere $S_\epsilon$ (area element $r^2\,\mathrm{d}\Omega$ ): $\displaystyle\oint\!\left(-\frac{A}{r^2}\right)r^2\,\mathrm{d}\Omega=-4\pi A=1$ , giving $A=-1/4\pi$ . Why: the area $4\pi r^2$ cancels the $1/r^2$ exactly, making the answer independent of $\epsilon$ — the hallmark of a genuine point source.

Step 5. Therefore $G(\mathbf{x})=-\dfrac{1}{4\pi r}$ . The same argument in $n$ dimensions gives $G=-\Gamma(n/2-1)/(4\pi^{n/2}r^{n-2})$ for $n\ge3$ and the logarithm $G=\tfrac{1}{2\pi}\ln r$ in two dimensions, where no power-law decays. Why: the dimension enters through the surface area of the unit sphere, so the potential's falloff is a statement about the geometry of space itself. *As a check: the surface flux $\oint\nabla(-1/4\pi r)\cdot\mathrm{d}\mathbf{S}=1$ evaluates to unity, confirming the normalization. Numerically $-1/(4\pi)=-0.0796$ at $r=1$ .*

Go deeper · Derivation B · 1 minThe heat kernel by Fourier transform

Derivation BThe heat kernel by Fourier transform

The aim is to construct the Green's function of the diffusion operator, $(\partial_t-\kappa\nabla^2)G=\delta(\mathbf{x})\delta(t)$ in $n$ space dimensions, and show it is a spreading Gaussian normalized to unit mass. The method: Fourier transform in space so $\nabla^2\to-k^2$ , reducing the PDE to a first-order ODE in time for each mode; solve that ODE with the causal condition that nothing exists before the source fires; then invert, which factorizes into identical Gaussian integrals, one per dimension.

Step 1. Transform in space, $\nabla^2\to-k^2$ and $\delta(\mathbf{x})\to1$ , giving one ODE per wavevector: $\partial_t\tilde{G}+\kappa k^2\tilde{G}=\delta(t)$ . Why: diffusion does not couple different Fourier modes, so the transform diagonalizes the operator and each mode evolves independently.

Step 2. For $t>0$ the right side vanishes and the mode decays, $\tilde{G}(\mathbf{k},t)=\Theta(t)\,e^{-\kappa k^2 t}$ , where the $\delta(t)$ imposes the jump $\tilde{G}(0^+)=1$ from the causal past $\tilde{G}(t<0)=0$ . Why: high- $k$ (short-wavelength) modes decay fastest — diffusion erases fine structure first, which is why it smooths and cannot be run backward.

Step 3. Invert: $\displaystyle G(\mathbf{x},t)=\Theta(t)\!\int\!\frac{\mathrm{d}^n k}{(2\pi)^n}\,e^{i\mathbf{k}\cdot\mathbf{x}-\kappa k^2 t}$ , a product of $n$ identical one-dimensional Gaussian integrals. Why: the exponent separates over Cartesian components, so the integral is the $n$ th power of a single line integral.

Step 4. Each factor, by completing the square ( $\int e^{-ak^2+bk}\,\mathrm{d}k=\sqrt{\pi/a}\,e^{b^2/4a}$ ), gives $\displaystyle\int\frac{\mathrm{d}k}{2\pi}\,e^{ikx-\kappa k^2 t}=\frac{1}{\sqrt{4\pi\kappa t}}\,e^{-x^2/4\kappa t}$ . Why: a Gaussian in $k$ transforms to a Gaussian in $x$ whose width grows as $\sqrt{t}$ — the diffusive scaling.

Step 5. Assemble the $n$ factors: $\displaystyle G=\frac{\Theta(t)}{(4\pi\kappa t)^{n/2}}\,e^{-|\mathbf{x}|^2/4\kappa t}$ , which integrates to $\int G\,\mathrm{d}^n x=1$ for every $t>0$ and tends to $\delta(\mathbf{x})$ as $t\to0^+$ . Why: unit mass at all times is conservation of heat (or of probability — this is also the Brownian transition density); and $G>0$ everywhere the instant $t>0$ , so diffusion has infinite propagation speed. *As a check: $\int_{\mathbb{R}^n}G\,\mathrm{d}^n x=1$ for all $t>0$ and $G\to\delta$ as $t\to0^+$ , confirming unit-mass conservation and the correct initial condition. In 1D with $\kappa=10^{-4}\ \mathrm{m^2/s}$ the RMS spread $\sqrt{2\kappa t}=1.41\ \mathrm{cm}$ at $t=1\ \mathrm{s}$ .*

Go deeper · Derivation C · 1 minThe retarded Green's function of the wave equation

Derivation CThe retarded Green's function of the wave equation

The target is the causal kernel of the three-dimensional wave operator, $(\partial_t^2/c^2-\nabla^2)G=\delta(\mathbf{x})\delta(t)$ , and the claim is that it is a delta-function shell supported strictly on the forward light cone $r=ct$ . The route: Fourier transform in space and time, invert the operator algebraically, and confront the two poles it places on the real $\omega$ axis; causality is injected by pushing those poles off the axis with the $+i\epsilon$ prescription, after which the time integral is a residue sum and the angular $\mathbf{k}$ integral collapses to a single delta.

Step 1. Transform in space and time. The operator becomes $-\omega^2/c^2+k^2$ , so $\tilde{G}(\mathbf{k},\omega)=\dfrac{-c^2}{\omega^2-c^2k^2}$ , with simple poles at $\omega=\pm ck$ on the contour. Why: the poles are the free waves of the homogeneous equation (the dispersion relation), which is exactly why $G$ is ambiguous until we say how to handle them.

Step 2. Enforce causality by shifting $\omega\to\omega+i\epsilon$ , moving both poles into the lower half-plane. For $t<0$ the factor $e^{-i\omega t}$ forces the contour closed in the upper half-plane, which now encloses no pole, so $G=0$ before the source fires. Why: "no effect before its cause" is not automatic — it is a choice of contour, and this is the choice that makes it true.

Step 3. For $t>0$ close below and sum the two residues: $\displaystyle\int\frac{\mathrm{d}\omega}{2\pi}\frac{-c^2e^{-i\omega t}}{(\omega+i\epsilon)^2-c^2k^2}=\Theta(t)\,\frac{c}{k}\sin(ckt)$ . Why: the two poles combine into an oscillation, the free-field response of a single mode kicked at $t=0$ .

Step 4. Do the angular integral with $\int e^{i\mathbf{k}\cdot\mathbf{x}}\,\mathrm{d}\Omega=4\pi\sin(kr)/(kr)$ , leaving $\displaystyle G_{\mathrm{ret}}=\frac{\Theta(t)}{2\pi^2 r}\int_0^\infty c\,\sin(kr)\sin(ckt)\,\mathrm{d}k$ . Why: rotational symmetry reduces the 3D $\mathbf{k}$ -integral to a radial one, at the cost of one power of $r$ — the origin of the $1/r$ amplitude.

Step 5. Convert the product of sines to a difference of cosines and use $\int_0^\infty\cos(k\xi)\,\mathrm{d}k=\pi\delta(\xi)$ to obtain $\displaystyle G_{\mathrm{ret}}=\frac{1}{4\pi r}\,\delta\!\left(t-\frac{r}{c}\right)$ . Why: the kernel is nonzero only when $t=r/c$ — a sharp wavefront travelling at exactly $c$ , Huygens' principle in odd space dimensions. *As a check: the support lies strictly on $r=ct$ , confirming the light-cone structure. Flipping the sign of $i\epsilon$ gives the advanced kernel $\delta(t+r/c)/(4\pi r)$ .*

Limiting cases

The single defining equation $\mathcal{L}G=\delta$ produces qualitatively different kernels as the operator changes type and the boundary or radiation condition changes. The table collects the four canonical cases.

Case	Result	Meaning
Poisson (elliptic), 3D free space	$G=-\dfrac{1}{4\pi\vert\mathbf{x}-\mathbf{x}'\vert}$	Smooth, instantaneous, $1/r$ decay; the Coulomb/Newton potential, $\nabla^2 G=\delta$ .
Heat (parabolic), $n$ dimensions	$G=\dfrac{\Theta(t-t')}{[4\pi\kappa(t-t')]^{n/2}}\,e^{-\vert\mathbf{x}-\mathbf{x}'\vert^2/4\kappa(t-t')}$	Causal ( $t>t'$ ) but instantly nonzero everywhere — infinite propagation speed, irreversible.
Wave (hyperbolic), 3D retarded	$G_{\mathrm{ret}}=\dfrac{\delta\!\big(t-t'-\vert\mathbf{x}-\mathbf{x}'\vert/c\big)}{4\pi\vert\mathbf{x}-\mathbf{x}'\vert}$	Supported strictly on the forward light cone; a sharp signal at speed $c$ (Huygens).
Helmholtz (oscillatory), 3D outgoing	$G_+=-\dfrac{e^{ik\vert\mathbf{x}-\mathbf{x}'\vert}}{4\pi\vert\mathbf{x}-\mathbf{x}'\vert}$	Solves $(\nabla^2+k^2)G=\delta$ ; the Sommerfeld radiation condition picks outgoing waves.

Common trap"The Green's function of an operator is a single, well-defined object."

It is not — not until you have named the boundary or initial conditions. Equation (1) fixes $G$ only up to a solution of the homogeneous equation $\mathcal{L}v=0$ , and that freedom is where all the physics lives. The free-space fundamental solution and the Dirichlet Green's function for the same Laplacian differ by a harmonic function chosen to vanish on the boundary. For the wave operator the retarded, advanced, and Feynman kernels all satisfy $\mathcal{L}G=\delta$ , yet one is a signal going forward in time, one backward, and one the time-ordered quantum propagator — utterly different physics from one equation. So "the Green's function" is meaningless shorthand: always state the conditions, or equivalently the contour, before you write $G$ .

Worked examples

1 · Point charge above a grounded plane, by images (fully worked)

A charge $q=1.0\ \mathrm{nC}$ sits a distance $d=1.0\ \mathrm{cm}$ above an infinite grounded conducting plane. Use the Dirichlet Green's function to find the force on the charge, the peak induced surface-charge density, and the total induced charge.

Step 1. The grounded plane demands $\phi=0$ on the surface, achieved by placing an image $-q$ at the mirror point $\mathbf{x}'^{*}$ , a distance $d$ below the plane. The potential above is then charge plus image, $\phi(\mathbf{x})=\dfrac{q}{4\pi\varepsilon_0}\!\left(\dfrac{1}{r_+}-\dfrac{1}{r_-}\right)$ , with $r_\pm$ the distances to each; on the plane $r_+=r_-$ , so $\phi=0$ as required.

Step 2. The force on the real charge is the Coulomb attraction to its image, separated by $2d$ : $F=\dfrac{1}{4\pi\varepsilon_0}\dfrac{q^2}{(2d)^2}=\dfrac{q^2}{16\pi\varepsilon_0 d^2}$ . Numerically, $F=\dfrac{(10^{-9})^2}{16\pi(8.854\times10^{-12})(0.01)^2}=\mathbf{2.25\times10^{-5}\ N}$ , toward the plane — a real pull with no second charge present, supplied by the conductor.

Step 3. The induced surface-charge density follows from $\sigma=-\varepsilon_0\,\partial_n\phi$ evaluated on the conductor, giving the standard result $\sigma(\rho)=-\dfrac{q\,d}{2\pi(\rho^2+d^2)^{3/2}}$ , where $\rho$ is the in-plane distance from the foot of the charge. It peaks directly beneath the charge, $\sigma(0)=-\dfrac{q}{2\pi d^2}=-\dfrac{10^{-9}}{2\pi(0.01)^2}=\mathbf{-1.59\ \mu C/m^2}$ .

Step 4. Integrate over the plane: $\displaystyle Q_{\mathrm{ind}}=\int_0^\infty\!\!\sigma(\rho)\,2\pi\rho\,\mathrm{d}\rho=-q\!\int_0^\infty\!\frac{d\,\rho\,\mathrm{d}\rho}{(\rho^2+d^2)^{3/2}}=-q=\mathbf{-1.0\ nC}$ . The conductor draws up exactly $-q$ — every field line from $q$ terminates on the plane, which is why the image trick works. F = 2.25×10⁻⁵ N toward the plane, σ(0) = −1.59 μC/m², total induced charge = −q exactly.

2 · Diffusion length from the heat kernel

A thin rod of thermal diffusivity $\kappa=1.0\times10^{-4}\ \mathrm{m^2/s}$ receives a sharp pulse of heat at one point at $t=0$ . How far has the heat spread after 1 s and after 100 s, and how long until it spreads 10 cm?

Step 1. The 1D heat kernel $G=(4\pi\kappa t)^{-1/2}e^{-x^2/4\kappa t}$ is a Gaussian of variance $\langle x^2\rangle=2\kappa t$ , so the characteristic spread is the RMS width $\sqrt{2\kappa t}$ — growing as $\sqrt{t}$ , not $t$ , the signature of diffusion.

Step 2. At $t=1\ \mathrm{s}$ : $\sqrt{2\kappa t}=\sqrt{2\times10^{-4}}=\mathbf{1.41\ cm}$ . At $t=100\ \mathrm{s}$ : $\sqrt{2\times10^{-2}}=\mathbf{14.1\ cm}$ . A hundredfold longer wait spreads the heat only tenfold farther — the $\sqrt{t}$ law, which is why diffusion moves heat so slowly over long distances.

Step 3. Invert for the time to reach $L=0.10\ \mathrm{m}$ : $t=\dfrac{L^2}{2\kappa}=\dfrac{(0.10)^2}{2\times10^{-4}}=\mathbf{50\ s}$ . The quadratic $t\propto L^2$ is why cooking scales badly with size while micron-scale chips reach steady state almost instantly. Spread 1.41 cm at 1 s, 14.1 cm at 100 s; 50 s to reach 10 cm.

3 · Green's function of a string, built by matching

Solve $-\dfrac{\mathrm{d}^2 G}{\mathrm{d}x^2}=\delta(x-x')$ on $[0,L]$ with $G(0)=G(L)=0$ — the static deflection of a stretched string pulled at one point.

Step 1. Away from $x'$ the source vanishes, so $G$ is linear in each region. Impose the boundary conditions: $G=a\,x$ for $x<x'$ (vanishes at $0$ ) and $G=b\,(L-x)$ for $x>x'$ (vanishes at $L$ ) — a straight segment on either side of the pull.

Step 2. Continuity at $x=x'$ requires $a\,x'=b\,(L-x')$ ; integrating the equation across the spike gives the slope jump $G'(x'^+)-G'(x'^-)=-1$ , i.e. $-b-a=-1$ . Solving the pair: $a=(L-x')/L$ , $b=x'/L$ .

Step 3. Hence $G(x,x')=\dfrac{x_<\,(L-x_>)}{L}$ , with $x_<=\min(x,x')$ , $x_>=\max(x,x')$ — manifestly symmetric, the reciprocity of this self-adjoint problem. For $L=1$ pulled at the centre $x'=0.5$ , the peak deflection is $G(0.5,0.5)=0.25$ ; expanding this same $G$ in the sine eigenfunctions reproduces equation (6) term by term. G = x_<(L − x_>)/L, symmetric, peak 0.25 at the centre of a unit string.

Physics of the portrait

The Green's function carries the name of a man who had almost no formal schooling and whose defining work went unread for two decades — a reminder that the deepest tools of physics sometimes come from the most unlikely places.

Origins · 1 minGeorge Green and how this was discovered

George Green

1793–1841, English mathematician and physicist · potential theory and boundary-value problems

George Green was a miller's son from Nottingham with barely one year of schooling, who taught himself mathematics while working in the family windmill. In 1828 he published, by private subscription in an edition of about fifty copies, An Essay on the Application of Mathematical Analysis to the Theories of Electricity and Magnetism. In it he introduced the "potential function", the identities now called Green's theorem, and the auxiliary function — the Green's function — that solves a boundary-value problem by an integral over its source and boundary data. The essay was so far ahead of its time, and so obscurely circulated, that it sank almost without trace. Green entered Gonville and Caius College, Cambridge, as an undergraduate in 1833 at the age of forty, graduated in 1837, and died in 1841, still largely unknown. Four years later a young William Thomson (later Lord Kelvin) found a copy and had it republished, from where it spread through the whole of nineteenth-century physics — the mathematical common denominator of "how does a source make a field", worked out beside a grinding millstone.

History

1828

— George Green's Essay introduces the potential function and the function that bears his name

1845

— William Thomson rediscovers Green's Essay and, with the method of images, revives potential theory

1867

— Ludvig Lorenz writes the retarded potentials, the wave equation's Green's function at work

1882

— Kirchhoff casts diffraction as a boundary integral over the free-space wave kernel

1948

— Feynman's iε propagator makes the Green's function the backbone of quantum field theory

Going further

Texts. Stone & Goldbart, Mathematics for Physics, Ch. 11 (Green Functions) · Arfken, Weber & Harris, Mathematical Methods for Physicists, 7th ed., Ch. 10 · Hassani, Mathematical Physics, Ch. 21–22 · Byron & Fuller, Mathematics of Classical and Quantum Physics, Ch. 7 · Jackson, Classical Electrodynamics, 3rd ed., Ch. 1 & 6 · Morse & Feshbach, Methods of Theoretical Physics, Ch. 7. The canonical kernels and numerical results collected here are: the Gauss-flux normalization giving $-1/(4\pi r)$ in 3D, equal to $-0.0796$ at $r=1$ ; the heat kernel normalized to $1$ with $G\to\delta$ as $t\to0^+$ and characteristic spread $\sqrt{2\kappa t}$ ; the retarded wave kernel with light-cone support obtained by contour integration; the image-charge force $F=q^2/16\pi\varepsilon_0 d^2=2.25\times10^{-5}\ \mathrm{N}$ and total induced charge $-q$ ; and the string Green's function $G=x_<(L-x_>)/L$ consistent with both the matching conditions and the sine-series expansion. Constants are exact.

Prerequisites. Classification of Second-Order Linear PDEs · Boundary Conditions and Well-Posedness · Sturm-Liouville Theory and Eigenfunction Expansions · Laplace's and Poisson's Equations. Leads to. The Helmholtz Equation · The Wave Equation · The Heat (Diffusion) Equation · Method of Characteristics.

Ask the tutor: Point me at any operator and boundary — a drum, a grounded sphere, a diffusing gas — and I will build its Green's function step by step and show you the field it makes.

What is the difference between a fundamental solution and a Green's function?

Both solve $\mathcal{L}G=\delta$ , but they answer to different conditions. A fundamental solution is the response in unbounded free space, fixed only by how it must decay or radiate at infinity — for example $-1/(4\pi r)$ for the 3D Laplacian. A Green's function proper lives on a bounded domain and additionally satisfies the homogeneous boundary conditions there, say $G=0$ on $\partial\Omega$ . The two differ by a solution of $\mathcal{L}v=0$ chosen to enforce the boundary data — exactly what the method of images and the eigenfunction expansion construct. So a fundamental solution is the raw building block; the Green's function is that block tailored to fit a particular room.

What does a "point source" actually mean physically?

The Dirac delta $\delta(\mathbf{x}-\mathbf{x}')$ idealizes a source concentrated at one point with unit total strength: a unit point charge, a unit impulse of heat, a hammer-tap on a membrane. No real source is truly pointlike, but any smooth source is a continuous sum of them, $f(\mathbf{x})=\int f(\mathbf{x}')\delta(\mathbf{x}-\mathbf{x}')\,\mathrm{d}\mathbf{x}'$ , so the point source is the natural basis. Because $\mathcal{L}$ is linear, the response to the whole source is the same sum of point responses — the integral $u=\int Gf$ . The delta is thus a bookkeeping device: it isolates the response to one location so everything else follows by superposition.

How can the heat kernel have infinite propagation speed without breaking relativity?

The heat kernel $G=(4\pi\kappa t)^{-n/2}e^{-r^2/4\kappa t}$ is strictly positive for every $r$ the instant $t>0$ , so a pulse at the origin registers, however faintly, arbitrarily far away immediately. This does not violate relativity because the diffusion equation is not a relativistic law; it is a coarse-grained description of many random microscopic collisions, valid only for speeds far below $c$ . The "signal" at large $r$ is exponentially tiny and carries no usable information, and a relativistic treatment (the telegrapher's equation) restores a finite front. So the infinite speed is an artifact of the approximation — a reminder that a Green's function inherits the assumptions of its operator.

Why are there several Green's functions for one wave operator?

Because $\mathcal{L}G=\delta$ fixes $G$ only up to a homogeneous solution, and the wave operator has a rich space of free waves to add. In Fourier space this is poles of $\tilde G$ sitting right on the contour at $\omega=\pm ck$ ; the integral is undefined until you say how to go around them. Pushing both poles below the real axis gives the retarded kernel, nonzero only for $t>0$ — cause precedes effect. Pushing both above gives the advanced kernel, its time-reverse. Splitting them — positive-frequency pole below, negative above — gives the Feynman propagator, forward for particles and backward for antiparticles, the object every QFT diagram is built from. Same equation, three contours, three physics.

Why does self-adjointness force reciprocity?

Apply Green's second identity with $u=G(\mathbf{x},\mathbf{x}_1)$ , $v=G(\mathbf{x},\mathbf{x}_2)$ . The volume side becomes $\int(G_1\,\mathcal{L}G_2-G_2\,\mathcal{L}G_1)\,\mathrm{d}V=G(\mathbf{x}_2,\mathbf{x}_1)-G(\mathbf{x}_1,\mathbf{x}_2)$ by sifting. If $\mathcal{L}$ is formally self-adjoint and both kernels obey the same homogeneous conditions, the surface term vanishes, so the difference is zero: $G(\mathbf{x}_1,\mathbf{x}_2)=G(\mathbf{x}_2,\mathbf{x}_1)$ . Physically, source and detector are interchangeable — the principle exploited in seismic interferometry.

When does the eigenfunction expansion fail, and what is the fix?

The spectral formula $G=\sum_n\phi_n(\mathbf{x})\phi_n^*(\mathbf{x}')/\lambda_n$ diverges the moment any eigenvalue $\lambda_n=0$ — a zero mode. The classic case is the Neumann Laplacian on a closed domain, where the constant is a zero mode: $\nabla^2 u=f$ has no unique solution (add any constant) and none at all unless the source integrates to zero. The remedy is a modified (pseudo) Green's function: solve $\mathcal{L}G=\delta-\phi_0(\mathbf{x})\phi_0^*(\mathbf{x}')$ , which requires the source to be orthogonal to the zero mode. That orthogonality is the Fredholm solvability condition — physically a conservation law, such as "total flux through a closed insulated boundary must vanish".