A third-order linear ordinary differential equation in general has the form , where is the unknown function and are known functions of . These equations are typically very difficult to solve explicitly, so instead, people try to at least understand their behavior through simplifications and tricks. One such trick is presented in this essay: how to convert a general 3rd order linear ODE to the (relatively) simple form via a change of variables - in other words, how to set the coefficients before and to zero. The computation is fairly nontrivial: it involves solving another ODE and taking a couple of difficult integrals; however, more often than not it is doable and thus makes a valuable addition to the mathematical toolbox.
u'''+au''+bu'+cu=f
u=u[t]
a,b,c,f
t
v'''+Zv=F
v''
v'
w''[t]=[t]w[t]
2
p
Deriving the formula
Deriving the formula
A general third-order linear ordinary differential equation has the form
u'''[t]+a[t]u''[t]+b[t]u'[t]+c[t]u[t]=f[t]
(
1
) The equation is currently written in terms of , however, we can make a change of variables and re-write it in terms of some other function . There are many ways in which a change of variables could be done. Consider the following transformation:
u[t]
v[t]
u[t]=A[t]v[ϕ[t]]
(
2
) Effectively, this change of variables replaces with and multiplies the expression by some scaling function . Substituting this into (1) yields
t
ϕ
A[t]
Av+(aA+3A'+3Aϕ'ϕ'')v+(Abϕ'+2aA'ϕ'+3ϕ'A''+aAϕ''+3A'ϕ''+Aϕ''')+(Ac+bA'+aA''+A''')v=f[t]
3
(ϕ')
3
d
3
dϕ
2
(ϕ')
2
(ϕ')
2
d
2
dϕ
dv
dϕ
(
3
) That is how the ODE is written in terms of the variables and . Effectively, it is in the same form as the initial ODE (1), the only difference being that the coefficients have changed. Here, denotes a derivative and is reserved for the time derivative. At this point, the functions and could be anything, but for the change of variables to be useful, let us try to find such and that the coefficients before v and are zero:
v
ϕ
d
dϕ
ϕ
'
A
ϕ
A
ϕ
2
d
2
dϕ
dv
dϕ
aA+3A'+3Akk'=0Abk+2aA'k+3kA''+aAk'+3A'k'+Ak''=0
2
k
2
k
(
4
) Here, is a short-hand notation for . If we solve these equations, the ODE (3) will have the desired simple form. Dividing the first equation by yields . Integrating this with respect to time, we get a[]+3Log[A]+3Log[k]=. Solving for yields
k≡ϕ'[t]
dϕ/dt
A
2
k
a+3+3=0
A'
A
k'
k
t
∫
t
1
dt
1
∼
A[t]
A=Exp-a[]
1
1
3
t
∫
t
1
dt
1
1
k
(
5
) Plugging this in to the second equation in (4), we get . The equation is nonlinear which makes it difficult to solve; however, it can be made linear by introducing
3-2=+a'-b
2
k'
k
k''
k
1
3
2
a
w[t]≡
1
k[t]
(
6
) The equation is written in terms of as
3-2=+a'-b
2
k'
k
k''
k
1
3
2
a
w
w''[t]=[t]w[t]
2
p
(
7
) Where [t]≡+. Unlike the previous equation, (7) is linear, which makes it tremendously easier to solve, both exactly and approximately. From (6), can then be computed as . The definition requires that . Coming back to the definition and requiring to be smooth, we see that the condition simply demands a one-to-one relationship between and . This is very natural because in the change of variables (2) is the “transformed time”, and it only makes sense if there is a one-to-one relationship between time and “transformed time”. Substituting this into (5) and setting =1, we get the final formula for :
2
p
2
a[t]
12
a'[t]-b
4
w[t]≡1
k[t]
k
k=
1
2
w
w[t]≡1
k[t]
k!=0
k[t]≡ϕ'[t]
ϕ[t]
k!=0
ϕ
t
ϕ
1
A[t]
A[t]=Exp-a[][t]
1
3
t
∫
t
1
dt
1
2
w
(
8
) Making equal to any other nonzero number (for instance, 2) wouldn’t really change much because it would just correspond to rescaling by some constant. Since , is the integral of :
1
A
k[t]≡ϕ'[t]
ϕ[t]
k[t]=1/
2
w
ϕ[t]=[]
t
∫
dt
1
2
w
t
1
(
9
) And this is it! The problem is solved. If we substitute , where and are given by (8) and (9), into the ODE (3), we would see that the coefficients before and are indeed zero. The ODE can be written in the form
u[t]=A[t]v[ϕ[t]]
A
ϕ
v''
v'
3
d
3
dϕ
(
10
) Where and are given by
Z
F
Z[ϕ]≡[t](4-18a[t](b[t]-[t])+9(6c[t]-3[t]+[t]))F[ϕ]=[t]f[t]
1
54
6
w
3
a[t]
′
a
′
b
′′
a
6
w
A[t]
(
11
)Code used to do the calculations
Code used to do the calculations
Computer implementation
Computer implementation
The function below executes the algorithm described above. It first computes and tries to find . If no solutions were found, it just returns the general formulas for . If one solution was found, it substitutes the formula for into (8) and (9) and outputs concrete formulas for . If multiple solutions were found, it outputs the general formulas for and gives the user the list of possible .
2
p
w[t]
ϕ,A,Z,F
w
A,ϕ,Z
ϕ,A,Z,F
w
In[]:=
ClearAll[findϕAZF];findϕAZF[a_Function, b_Function, c_Function, f_Function] := Block t, p2, w, possiblew, wIndicator, Aval , p2 = + ; possiblew = DSolve[w''[t] p2 w[t], w, t]; Which[ SameQ[Head[possiblew], DSolve], wIndicator = Evaluate[w''[t] p2 w[t]], (* if couldn't solve *) Length[possiblew] 1, w = possiblew[[1, 1, 2]]; wIndicator = Nothing, (* if found exactly one solution *) Length[possiblew] > 1, wIndicator = (w ∈ possiblew[[All, 1, 2]]) (* if found several solutions *) ]; Aval = Exp-Integrate[a[t], t] ; Simplify@ϕ Integrate, t, wIndicator, A Aval, Z (4 -18 a[t] (b[t]-[t])+9 (6 c[t]-3 [t]+[t])), F -> f[t];
2
a[t]
12
a'[t] - b[t]
4
1
3
2
w[t]
1
2
w[t]
1
54
6
w[t]
3
a[t]
′
a
′
b
′′
a
6
w[t]
Aval
In[]:=
ClearAll[getZFform];getZFform[{ϕ -> ϕval_, A -> Aval_, Z -> Zval_, F -> Fval_}] := Block[{texpr}, texpr = Solve[ϕval == ϕ, t]; Which[ SameQ[Head[texpr], Solve], Return[$Failed], Length[texpr] >= 1, texpr = texpr[[1, 1, 2]]; u[t] == Aval v[ϕ], t == texpr, Z -> Zval /. t -> texpr, F -> Fval /. t -> texpr ]];
Here is how to use this method:
1
.Get the ODE to the form ; identify the functions .
u'''[t]+a[t]u''[t]+b[t]u'[t]+c[t]u[t]=f[t]
a[t],b[t],c[t],f[t]
2
.Run findϕAZF to find .
ϕ[t],A[t],Z[t],F[t]
3
.The expressions for are going to contain two arbitrary constants. Choose values for the constants, preferably to make the expressions simpler.
ϕ[t],A[t],Z[t],F[t]