Asymptotic solvers for second-order
differential equation systems with
multiple frequencies
Marissa CONDON
School of Electronic Engineering, Dublin City University
E-mail: marissa.condon@dcu.ie
Alfredo DEA
˜
NO
Depto.de Matem´aticas, Universidad Carlos III de Madrid
E-mail: alfredo.deanho@uc3m.es
Jing GAO
School of Mathematics and Statistics, Xi’an Jiaotong University
E-mail: jgao@mail.xjtu.edu.cn
Arieh ISERLES
DAMTP, Centre for Mathematical Sciences, University of Cambridge
E-mail: A.Iserles@damtp.cam.ac.uk
Abstract
In this paper, an asymptotic expansion is constructed to solve
second-order differential equation systems with highly oscillatory forc-
ing terms involving multiple frequencies. An asymptotic expansion
is derived in inverse of powers of the oscillatory parameter and its
truncation results in a very effective method of dicretizing the differ-
ential equation system in question. Numerical experiments illustrate
the effectiveness of the asymptotic method in contrast to the standard
Runge–Kutta method.
1 Introduction
The behaviour of signals comprising several non-commensurate frequen-
cies is very important in the design and analysis of electronic circuits. Non-
linearities in circuits can result in such signals giving rise to distortion and
resulting in degradation of performance. Hence, the accurate and efficient
1
simulation of circuit behaviour in the presence of such signals is essential.
It is to address this issue that the current paper is directed. As an exam-
ple of nonlinear circuits, the Van der Pol oscillator shall be considered. It
has numerous applications in science and engineering, for example, from de-
scribing the action potentials of biological neurons [5, 9] to the modelling
of resonant tunneling diode circuits [12]. A coupled Van der Pol–Duffing
system shall also be considered. Such coupled systems have applications in
secure communications, [6, 10]. In addition, these coupled systems can be
realised using analog circuitry [1].
Consider a second-order differential equation system of the form
y
00
(t) + f(y(t))y
0
(t) + g(y(t)) = F
ω
(t), t ≥ 0, (1.1)
where f : C
d
→ C
d
, g : C
d
→ C
d
are two analytical functions,
f(y) =
f
11
(y) f
12
(y) ··· f
1d
(y)
f
21
(y) f
22
(y) ··· f
2d
(y)
.
.
.
f
d1
(y) f
d2
(y) ··· f
dd
(y)
d×d
,
g =
g
1
(y)
g
2
(y)
.
.
.
g
d
(y)
d×1
,
where every entry f
jv
(y) : C
d
→ R and g
j
(y) : C
d
→ R is an analytic scalar
function for j, v = 1, 2, ··· , d. The initial conditions are y(0) = y
0
∈ C
d
and
y
0
(0) = y
0
0
∈ C
d
, and the forcing term F
ω
(t) is
F
ω
(t) =
M
X
m=1
a
m
(t)e
iω
m
t
,
in which a
1
, ··· , a
M
∈ R
+
→ C
d
are analytic functions. Note that we have
assumed that there is a finite set of frequencies ω
1
, ··· , ω
M
∈ R \{0} in the
forcing term. At least some of these frequencies are large which results in
a highly oscillatory solution and one which is very expensive to obtain with
classical discretization methods. Furthermore, assume that the functions
f(y(t)) and g(y(t)) are analytic to ensure the existence and uniqueness of
the solution y(t).
2
When d = 1, ω
2m−1
= mω, ω
2m
= −mω, m = 0, 1, ··· , b
M
2
c, ω 1, the
above multiple frequency case reduces to a single frequency case
y
00
(t) + f(y(t))y
0
(t) + g(y(t)) =
∞
X
k=−∞
b
k
(t)e
ikωt
, t ≥ 0,
which has been already analysed in [3].
2 Construction of the asymptotic expansion
We start with a set U
0
= {1, 2, ··· , M } and ω
j
= κ
j
ω, j = 1, 2, ··· , M,
where ω serves as the oscillatory parameter. Therefore, the original equation
(1.1) may be written in the form
y
00
(t) + f(y(t))y
0
(t) + g(y(t)) =
M
X
m=1
a
m
(t)e
iκ
m
ωt
(2.1)
=
X
m∈U
0
a
m
(t)e
iκ
m
ωt
, t ≥ 0.
Our basic ansatz is that the solution y(t) admits an expansion in inverse
powers of the oscillatory parameter ω,
y(t) ∼
∞
X
r=0
1
ω
r
X
m∈U
r
p
r,m
(t)e
iσ
m
ωt
,
The sets U
r
, the scalars σ
m
and the functions p
r,m
(t) will be described in
the sequel. The important point to note at this stage is that the functions
p
r,m
(t) are independent of ω and can be derived recursively: p
r,0
(t) by
solving a non-oscillatory ODE and p
r,m
(t) for m 6= 0 by recursion.
It is very important, however, to impose p
0,m
≡ 0 and p
1,m
≡ 0 for
m 6= 0 so that differentiation does not result in the presence of a positive
power of ω in the resultant equations for p
r,m
(t) involved in the asymptotic
method.
Therefore, the proposed solution y(t) in inverse powers of the oscillatory
parameter ω is
y(t) ∼ p
0,0
(t) +
1
ω
p
1,0
(t) +
∞
X
r=2
1
ω
r
X
m∈U
r
p
r,m
(t)e
iσ
m
ωt
(2.2)
3
As just stated it will be assumed in the ansatz that p
1,0
is the only non-zero
p
1,m
in U
1
. So U
1
= {0}.
Following the approach in [4] (cf. also [11]), the expression (2.2) for y(t) is
substituted into the second-order differential equation (2.1). The first-order
derivative of y(t) is
y
0
∼ p
0
0,0
+
1
ω
p
0
1,0
+
∞
X
r=2
1
ω
r
X
m∈U
r
p
0
r,m
e
iσ
m
ωt
+ iσ
m
ωp
r,m
e
iσ
m
ωt
= p
0
0,0
+
1
ω
p
0
1,0
+
∞
X
r=2
1
ω
r
X
m∈U
r
p
0
r,m
e
iσ
m
ωt
+
∞
X
r=2
1
ω
r−1
X
m∈U
r
iσ
m
p
r,m
e
iσ
m
ωt
= p
0
0,0
+
1
ω
p
0
1,0
+
∞
X
r=2
1
ω
r
X
m∈U
r
p
0
r,m
e
iσ
m
ωt
+
∞
X
r=1
1
ω
r
X
m∈U
r+1
iσ
m
p
r+1,m
e
iσ
m
ωt
,
= p
0
0,0
+
1
ω
p
0
1,0
+
X
m∈U
2
iσ
m
p
2,m
e
iσ
m
ωt
+
∞
X
r=2
1
ω
r
X
m∈U
r
p
0
r,m
e
iσ
m
ωt
+
X
m∈U
r+1
iσ
m
p
r+1,m
e
iσ
m
ωt
Likewise, the second-order derivative of y(t) is
y
00
(t) ∼ p
00
0,0
+
1
ω
p
00
1,0
+
X
m∈U
2
iσ
m
p
0
2,m
+ (iσ
m
ωp
2,m
)
e
iσ
m
ωt
+
∞
X
r=2
1
ω
r
"
X
m∈U
r
p
00
r,m
+ iσ
m
ωp
0
r,m
e
iσ
m
ωt
+
X
m∈U
r+1
iσ
m
p
0
r+1,m
+ iσ
m
ωp
r+1,m
e
iσ
m
ωt
,
= p
00
0,0
+
X
m∈U
2
(iσ
m
)
2
p
2,m
e
iσ
m
ωt
+
1
ω
p
00
1,0
+ 2
X
m∈U
2
(iσ
m
)p
0
2,m
e
iσ
m
ωt
+
X
m∈U
3
(iσ
m
)
2
p
3,m
e
iσ
m
ωt
+
∞
X
r=2
1
ω
r
"
X
m∈U
r
p
00
r,m
e
iσ
m
ωt
4
+ 2
X
m∈U
r+1
(iσ
m
)p
0
r+1,m
e
iσ
m
ωt
+
X
m∈U
r+2
(iσ
m
)
2
p
r+2,m
e
iσ
m
ωt
The function f(y(t))
d×d
is analytic and its Taylor expansion about
p
0,0
(t) may be determined. f
(n)
jv
(p
0,0
)[η
1
, ··· , η
n
] : C
d
×
n times
z }| {
C
d
× ··· × C
d
→ C
is the nth derivative operator which is linear in each of η
k
s such that
f
jv
(y
0
+ t) = f
jv
(y
0
) +
∞
X
n=1
t
n
n!
f
(n)
jv
(y
0
)[, ··· , ]
for sufficiently small |t| > 0. Hence,
f
jv
(y) = f
jv
p
0,0
(t) +
1
ω
p
1,0
(t) +
∞
X
r=2
1
ω
r
X
m∈U
r
p
r,m
(t)e
iσ
m
ωt
!
= f
jv
(p
0,0
) +
∞
X
n=1
1
n!
f
(n)
jv
(p
0,0
)
∞
X
`
1
=1
1
ω
`
1
X
k
1
∈U
`
1
p
`
1
,k
1
e
iσ
k
1
ωt
, ··· ,
∞
X
`
n
=1
1
ω
`
n
X
k
n
∈U
`
n
p
`
n
,k
n
e
iσ
k
n
ωt
= f
jv
(p
0,0
) +
∞
X
n=1
1
n!
∞
X
`
1
=1
···
∞
X
`
n
=1
1
ω
`
1
+···+`
n
X
k
1
∈U
`
1
···
X
k
n
∈U
`
n
f
(n)
jv
(p
0,0
)
p
`
1
,k
1
, ··· , p
`
n
,k
n
e
i(σ
k
1
+···+σ
k
n
)ωt
= f
jv
(p
0,0
) +
∞
X
n=1
1
n!
∞
X
r=n
1
ω
r
X
`∈I
o
n,r
X
k
1
∈U
`
1
···
X
k
n
∈U
`
n
f
(n)
jv
(p
0,0
)
p
`
1
,k
1
, ··· , p
`
n
,k
n
e
i(σ
k
1
+···+σ
k
n
)ωt
= f
jv
(p
0,0
) +
∞
X
r=1
1
ω
r
r
X
n=1
1
n!
X
`∈I
o
n,r
X
k
1
∈U
`
1
···
X
k
n
∈U
`
n
f
(n)
jv
(p
0,0
)
p
`
1
,k
1
, ··· , p
`
n
,k
n
e
i(σ
k
1
+···+σ
k
n
)ωt
,
where
I
o
n,r
= {` = (`
1
, ··· , `
n
)
T
∈ N
n
:
n
X
j=1
`
j
= r}, 1 ≤ n ≤ r.
5