We've recently introduced some new functionality in Chebfun for computing the nullspace of differential operators. Let's explore this with a couple of simple examples.

### 1. Simple example #1

Let's start as simply as we can, and take

$$ (Lu)(x) := u''(x), \quad x\in [-1, 1]. $$

L = chebop(@(u) diff(u, 2));

Clearly the nullspace of this operator -- that is, the space of functions $v$ for which $L(v)=0$ -- is spanned the two functions

v = [1, chebfun('x')]; norm(L(v))

ans = 0

Supposing we didn't know this, we could compute a basis for the nullspace with the `null`

method:

LW = 'LineWidth'; lw = 1.6; V = null(L) plot(V, LW, lw) V'*V norm(L(V))

V = chebfun column1 (1 smooth piece) interval length endpoint values [ -1, 1] 4 -1.4 0.44 Epslevel = 4.009500e-16. Vscale = 1.384490e+00. chebfun column2 (1 smooth piece) interval length endpoint values [ -1, 1] 4 -0.29 1.3 Epslevel = 4.132711e-16. Vscale = 1.343214e+00. ans = 1.000000000000000 -0.000000000000000 -0.000000000000000 1.000000000000000 ans = 2.644251170154956e-13

where we find that $V^T V = I$ and $LV \approx 0$ as required.

Clearly `V`

doesn't correspond directly to $1$ and $x$, since there is some freedom in how we orthogonalise the basis. However, we can check that `V`

and ${1, x}$ correspond to the same spaces by computing the angle between the spaces with the `subspace`

command.

subspace(v, V)

ans = 1.385544576431027e-14

### 2. Incomplete boundary conditions

Now let's consider the more complicated 2nd-order operator

$$ Lu = u'' + 0.1x(1-x^2)u' + \sin(x)u, \quad x\in [-\pi, \pi]. \qquad (*) $$

dom = [-pi, pi]; L = chebop(@(x, u) diff(u, 2) + .1*x.*(1-x.^2).*diff(u) + sin(x).*u, dom);

As before, it has a nullspace of rank 2.

V = null(L) plot(V, LW, lw) V'*V norm(L(V))

V = chebfun column1 (1 smooth piece) interval length endpoint values [ -3.1, 3.1] 64 1.8 0.043 Epslevel = 1.753455e-16. Vscale = 1.786120e+00. chebfun column2 (1 smooth piece) interval length endpoint values [ -3.1, 3.1] 64 -0.011 1.3 Epslevel = 2.468482e-16. Vscale = 1.268748e+00. ans = 0.999999999999998 0.000000000000000 0.000000000000000 1.000000000000000 ans = 3.096261776816828e-11

However, now suppose we impose one boundary condition, say, a Dirichlet condition at the left. This removes one degree of freedom, and we are left with a rank 1 nullspace.

L.lbc = 0; L.rbc = []; v = null(L) plot(v, LW, lw), shg v'*v norm(L(v))

v = chebfun column (1 smooth piece) interval length endpoint values [ -3.1, 3.1] 64 -1.7e-16 1.3 Epslevel = 1.708133e-15. Vscale = 1.268983e+00. ans = 1.000000000000000 ans = 1.427155766422003e-11

Clearly this null vector must satisfy the given condition $v(-\pi) = 0.$

v(-pi)

ans = -1.665334536937735e-16

### 3. An application

Where might these ideas be useful? Well, suppose we were interested in equation $(*)$ with a homogeneous Dirichlet condition at the left, and wanted to know what inhomogeneous Dirichlet condition gave the minimal 2-norm of the solution to $Lu = 1$. Rather than solving the linear system for a number of different boundary conditions (which would be computationally expensive) we could simply solve for one, say again a homogeneous Dirichlet condition,

L.rbc = 0; u = L\1; hold on, plot(u, '--r', LW, lw), hold off

and compute the rest by adding a scalar multiple of the null-function $v$.

E = chebfun(@(c) norm(u + c*v, 2), [-10, 10], 'vectorize', 'splitting', 'on'); plot(E,LW,lw)

We compute the 2-norm as a chebfun in the unknown variable $c$, which we can then minimise to obtain the minimal energy solution

[minE, c_star] = min(E) u_star = u + c_star*v plot(u_star,LW,lw)

minE = 4.121950420615883 c_star = 3.143771420957320 u_star = chebfun column (1 smooth piece) interval length endpoint values [ -3.1, 3.1] 64 -4e-15 4 Epslevel = 1.776357e-15. Vscale = 3.989391e+00.

So the condition we require is that $u(\pi)$ = `bc_star`

, where

bc_star = u_star(pi)

bc_star = 3.989391428267542

### 4. Exotic constraints

The Chebfun `null`

function can also handle the more exotic types of boundary conditions that can be imposed in Chebfun (see [1]). For example, suppose we wish to compute the nullspace of the 3rd-order piecewise-smoooth ODE

$$ Lu := 0.1u''' + \sin(x)u'' + u, \quad x\in[-1,1] $$

with the 'boundary' condition

$$ \int(u) = u(0). $$

dom = [-1, 1]; L = chebop(@(x, u) .1*diff(u, 3) + sin(x).*diff(u, 2) + u, dom); L.bc = @(x, u) sum(u) - u(0);

Here `null`

has no problems!

V = null(L) plot(V, LW, lw), shg V'*V

V = chebfun column1 (1 smooth piece) interval length endpoint values [ -1, 1] 64 1.1 -1.2 Epslevel = 2.269192e-16. Vscale = 1.223148e+00. chebfun column2 (1 smooth piece) interval length endpoint values [ -1, 1] 64 0.77 -0.15 Epslevel = 2.567754e-16. Vscale = 1.080928e+00. ans = 1.000000000000000 -0.000000000000000 -0.000000000000000 1.000000000000001

sum(V) - V(0,:) norm(L(V), 1)

ans = 1.0e-15 * -0.111022302462516 0.666133814775094 ans = 9.372923955371351e-09

### 5. References

- Chebfun Example ode-linear/NonstandardBCs