A formula for the
n
-th Laplacian of a Gaussian function
t
Δ
-
2
r

2
L

=
-2t
L
-
2
r

2
L

t
∑
n=0
t-n
(-1)
n
4
(2t+1)!
(t-n)!(2n+1)!
2n
r
L
The current essay gives a closed form expression for the
n
-th Laplacian of a Gaussian function in three dimensions. To derive the results, some simple mathematical patterns are first observed, and induction is then used to guess the formula for the
n
-th Laplacian. To check the formula, induction is used, and as a more tangible test, its output is also compared by what we get by brute-force calculation. The need for the formula became apparent when trying to solve the wave equation
Δϕ-
2
∂
t
ϕ=F
with the external force
F
being a Gaussian function
F=
-
2
r
/
2
L

using a Taylor series
ϕ[t,r]=
∞
∑
n=0
n
t
a
n
[r]
.
Consider the function
f[r]=
-
2
r

2
L

(
1
)
Where
r≡
2
x
+
2
y
+
2
z
is the radial coordinate and
L
is some positive constant. The Laplacian operator
Δf≡(
2
∂
f/∂
2
x
)+
2
∂
f∂
2
y
+(
2
∂
f/∂
2
z
)
can be written in spherical coordinates as
Δf=(1/
2
r
)
∂
∂r
2
r
∂f
∂r
+O
∂
∂θ
+O
∂
∂ϕ

, where
O
∂
∂θ

indicates that there are some terms related to
∂f
∂θ
. Because the function
f
is spherically symmetric,
∂f
∂θ
=
∂f
∂ϕ
=0
, and
Δf=(1/
2
r
)
∂
∂r
2
r
∂f
∂r
. Differentiating the formula (1), we get
Δf=
-2
L
-
2
r
2
L

-6+4
2
r
2
L
We see that the result is a polynomial multiplying a Gaussian. By induction it can be proven that applying the Laplacian operator to the resulting expression would also yield a polynomial (a different one) multiplying a Gaussian. Moreover, the polynomial is only going to have even order terms
2n
r
, and the order of the highest term will be equal to
2t
. Thus, the
n
-th Laplacian of a Gaussian (1) can always be written in the form
t
Δ
-
2
r

2
L

=
-2t
L
t
∑
n=0
C
tn
n
(2
2
r
/
2
L
)
(
2
)
where
C
tn
is a constant. The formula (2) is written the way it is to make
C
tn
simpler. To find the formula for
C
tn
, we can use the fact that applying the Laplacian to a function zero times doesn’t change the function (because applying an operator zero times isn’t really doing anything), so we know that
t
Δ
-
2
r
/
2
L

=
-
2
r
/
2
L

. From here,
C
00
=1
(
3
)
We know that
t
Δ
-
2
r
/
2
L

=Δ
t-1
Δ
-
2
r
/
2
L


. Substituting the formula (2), we obtain a recursion relation for the coefficients
C
tn
:
C
tn
=(2n+1)(2n+2)
C
t-1,n+1
-(4n+3)
C
t-1,n
+
C
t-1,n-1
(
4
)
Combining the formulas (3) and (4), we can compute the coefficients
C
tn
for the first few
t
’s:
In[]:=
Column@NestList[Table[With[{c=Function[n,If[0≤n≤Length[#]-1,#[[n+1]],0]]},(2n+2)(2n+3)c[n+1]-(4n+3)c[n]+c[n-1]],{n,0,Length[#]}]&,{1},6]
Out[]=
{1}
{-3,1}
{15,-10,1}
{-105,105,-21,1}
{945,-1260,378,-36,1}
{-10395,17325,-6930,990,-55,1}
{135135,-270270,135135,-25740,2145,-78,1}
We can see that the rightmost coefficients
C
tt
are equal to 1 for all
t
’s. Also, the coefficients
C
t0
go like
1
,
3
,
15
,
105
,
945
; in the spirit of fourth-grade math problems “spot the sequence”, we can see that each next coefficient is the previous one times an odd number. In other words, the numbers go like
1=1
,
3=1*3
,
15=1*3*5
,
105=1*3*5*7
,
945=1*3*5*7*9
. This suggests that there is a simple pattern in the coefficients. Let us now look for a “horizontal pattern”. For that, let us divide each entry of the matrix by the number to the left of it. In other words, let’s divide
C
tn
by
C
t,n-1
:
In[]:=
ColumnTable
#[[nplus1]]
#[[nplus1-1]]
,{nplus1,2,Length[#]}&/@NestList[Table[With[{c=Function[n,If[0≤n≤Length[#]-1,#[[n+1]],0]]},(2n+2)(2n+3)c[n+1]-(4n+3)c[n]+c[n-1]],{n,0,Length[#]}]&,{1},6]
Out[]=
{}
-
1
3

-
2
3
,-
1
10

-1,-
1
5
,-
1
21

-
4
3
,-
3
10
,-
2
21
,-
1
36

-
5
3
,-
2
5
,-
1
7
,-
1
18
,-
1
55

-2,-
1
2
,-
4
21
,-
1
12
,-
2
55
,-
1
78

Amazingly, there is a simple pattern! The first column goes like
1/3,2/3,3/3,4/3,...
; the second column goes like
1/10,2/10,3/10,...
; the denominator stays constant for a given column. For the first few columns, the denominators are
3,10,21,36,55,78
. Taking the differences, we get
3,7,11,15,19,23
. Again we are faced with a miraculously simple primary-school-style problem. Isn’t it amazing how these things occur in calculus? It can be seen that the differences are of the form
4n-1
. Summing them, we deduce that the denominators are given by
n(2n+1)
. From here,
C
tn
C
t,n-1
=-
t+1-n
n(2n+1)
. Using the fact that
C
tt
=1
, we get that
C
tn
=
t-1
∏
m=n
-
(n+1)(2n+3)
t-m
. Computing this product, we obtain a (remarkably simple) formula for
C
tn
:
C
tn
=
t-n
-
1
2
(2t+1)!
(t-n)!(2n+1)!
(
5
)
And the
n
-th Laplacian of a Gaussian function
Exp[-
2
r
/
2
L
]
is thus given by
t
Δ
-
2
r

2
L

=
-2t
L
-
2
r

2
L

t
∑
n=0
t-n
(-1)
n
4
(2t+1)!
(t-n)!(2n+1)!
2n
r
L
(
6
)
To test this formula, let us compare the 5-th Laplacian first using brute force and then using the formula (6):
In[]:=
Withlapl=Functionexpr,
1
2
r
D[
2
r
D[expr,r],r],Expand@Nest[lapl,Exp[-
2
r
/
2
L
],5]
Out[]=
-
332640
-
2
r
2
L

10
L
+
1108800
-
2
r
2
L

2
r
12
L
-
887040
-
2
r
2
L

4
r
14
L
+
253440
-
2
r
2
L

6
r
16
L
-
28160
-
2
r
2
L

8
r
18
L
+
1024
-
2
r
2
L

10
r
20
L
In[]:=
With{t=5},Expand
-2t
L
-
2
r

2
L

Sum
t-n
(-1)
n
4
(2t+1)!
(t-n)!(2n+1)!
2n
r
L
,{n,0,5}
Out[]=
-
332640
-
2
r
2
L

10
L
+
1108800
-
2
r
2
L

2
r
12
L
-
887040
-
2
r
2
L

4
r
14
L
+
253440
-
2
r
2
L

6
r
16
L
-
28160
-
2
r
2
L

8
r
18
L
+
1024
-
2
r
2
L

10
r
20
L
As you can see, the results match exactly. A formal induction proof based on the recursion relation (4) also confirms the formula (6).