/ tutorials / laplace-dirichlet
Sign in
tutorials / laplace-dirichlet
laplace-dirichlet / src / chapters / ch10 / text.md
60 lines · 6.6 KBPreviewCodeBlameHistoryRaw
1Every method up to Chapter 9 discretised the *area*: a grid of unknowns filling the domain, an operator relating each to its neighbours, and an error controlled by how finely the area is covered. This chapter does the opposite. The domain is a smooth closed curve
3$$r(t) = 1 + a\cos(n_{\text{lobe}}\, t),$$
5which no tensor-product grid fits, and the unknown lives on the boundary alone.
7### Reformulating on the boundary
9Look for the solution as the field of a dipole layer spread over the curve with density $\sigma$:
11$$u(x) \;=\; \frac{1}{2\pi}\int_\Gamma \frac{n_y \cdot (x-y)}{|x-y|^2}\,\sigma(y)\,ds(y).$$
13This is the double-layer potential. It is harmonic at every $x$ off the curve, whatever $\sigma$ is, because the integrand is harmonic in $x$; the equation is therefore satisfied *exactly*, and only the boundary condition remains. Letting $x$ approach the curve from inside, the potential jumps by $-\sigma/2$, and matching the data $g$ gives a Fredholm integral equation of the second kind,
15$$-\tfrac{1}{2}\sigma(x) + \int_\Gamma K(x,y)\,\sigma(y)\,ds(y) \;=\; g(x), \qquad x \in \Gamma.$$
17Two properties of this equation are what the chapter is about.
19**It is well conditioned.** The operator is the identity plus a compact operator, so its eigenvalues cluster at $-1/2$ and stay away from zero. The condition number of the discrete system is a modest constant, independent of the discretisation. Contrast the five-point matrix, whose condition number grows like $N^2$ and whose accuracy floor followed from that.
21**Its kernel is not singular.** $K(x,y)$ looks like it blows up as $y \to x$, and for a general kernel it would, but for the double layer on a smooth curve the numerator vanishes at the same rate as the denominator and the limit is finite:
23$$\lim_{y \to x} K(x,y) \;=\; -\frac{\kappa(x)}{4\pi},$$
25with $\kappa$ the curvature. So the integrand is a smooth periodic function of the parameter $t$, and for smooth periodic integrands the plain trapezoid rule is spectrally accurate. Discretising with $n$ equally spaced points in $t$ and the trapezoid rule gives an $n \times n$ dense system whose solution converges geometrically. This is the Nyström method, and on this geometry it is very hard to beat.
27For comparison the lab also runs the **method of fundamental solutions**: place $n$ point sources on a circle of radius $R$ outside the curve and fit their combined field to the boundary data by least squares. It is three lines of code, every basis function is harmonic, and it also converges geometrically.
29<!--lab:field-->
31<!--lab:sweep-->
33### What the numbers say
35On the default geometry, measured at interior points out to three quarters of the radius:
37| method | error | boundary unknowns | time |
38|---|---|---|---|
39| double-layer Nyström | $5 \times 10^{-16}$ | 512 | 13 ms |
40| fundamental solutions | $2 \times 10^{-14}$ | 128 | 12 ms |
41| five-point on the square (Ch. 6) | $1 \times 10^{-8}$ | $10^{6}$ | 440 ms |
43The first row is the point of the whole tutorial in one line: full double precision, on a curved domain, from a few hundred unknowns and a few milliseconds, because the formulation put the unknown where the information is and the quadrature matched the smoothness of the integrand.
45### Three things this does not do
47**Close evaluation.** The dashed curve in the convergence plot is the error at points at 98 per cent of the radius, and it is catastrophically worse: $10^{-1}$ where the interior error is $10^{-16}$. The reason is that the integrand, smooth for a target well inside, develops a near-singularity as the target approaches the curve, and the trapezoid rule stops being spectrally accurate. Nothing is wrong with the density $\sigma$: it is the *evaluation* that fails. Fixing it needs a different quadrature near the boundary, and the literature on that is substantial. The method of fundamental solutions does not have this problem, because its basis functions are smooth all the way to the curve, and the lab shows its near-boundary error tracking its interior error exactly.
49**Conditioning of the fitting approach.** The method of fundamental solutions wins on unknowns here, but its least-squares matrix is exponentially ill conditioned in the number of sources, and its accuracy depends on the radius $R$: too small and the sources cannot represent the field, too large and the matrix is numerically rank deficient. Move the slider and watch the error change by orders of magnitude for the same amount of work. There is no comparable knob on the Nyström method. That difference, a method with a free parameter that must be tuned against one without, matters more as problems get harder.
51**Density and dimension.** The Nyström matrix is dense, so the solve is $O(n^3)$ and the matrix is $O(n^2)$ in memory. In two dimensions with $n$ in the hundreds this is irrelevant. In three dimensions, with $n$ in the hundreds of thousands, it is fatal, and boundary integral methods only became practical there with fast summation: the fast multipole method and its relatives, which apply the operator in $O(n)$ and turn the dense solve into a preconditioned iteration. The well-conditioning of the second-kind equation is what makes that iteration converge in a handful of steps.
53### Where this leaves the tutorial
55The methods in these ten chapters differ by fourteen orders of magnitude in the accuracy they reach for a comparable amount of work, on problems whose exact solutions were known throughout. The pattern behind that spread is worth stating plainly.
57- **Accuracy is bounded by the smoothness of the solution, and only then by the method.** Chapters 4 and 9 are the same lesson twice: a spectral method on a solution with a $|x-\tfrac12|^{5/2}$ kink is no better than a three-point difference, and second-order differences on a re-entrant corner are no better than order $4/3$.
58- **Where the solution is smooth, high order is nearly free and enormously profitable.** The compact nine-point stencil costs the same per grid point as the five-point stencil and is ten orders of magnitude more accurate at $N = 32$.
59- **Structure, when it exists, is worth more than any amount of brute force.** The sine transform, the Chebyshev matrix decomposition and the second-kind integral equation are all the same move: identify what makes the problem special and let the solver use it.
60- **Every one of these statements has conditions, and the conditions are the useful part.** A separable operator on a rectangle, an analytic solution, a smooth boundary, a target not too close to it. Knowing which condition a method needs is what tells you when to reach for it, and the labs above are set up so that each condition can be switched off and the consequences watched.
↑↓ move↵openescclose