El Niño Temperature Anomalies Modeled by a Delay Differential Equation
Frank Wang, fwang@lagcc.cuny.edu
Delay differential equations are differential equations in which the derivative of the unknown function at a certain time depends on past values of the function and/or its derivatives. Max J. Suarez and Paul S. Schopf used such an equation to model the El Niño phenomenon. This worksheet demonstrate how Maple's dsolve command can be used to solve a delay differential equation numerically. This option was introduced in Maple 2015.
A Delay Differential Equation
Consider this equation
.
Let us assume the solution to be of the form .
restart:
y := t -> exp(w*t);
Eq1 := diff(y(t), t) = A*y(t-1);
Eq2 := expand(Eq1);
Soln := solve(Eq2, w);
We see that the solution is given as the Lambert W function (which was defined in Maple in its early years). Assume , we find the approximate value of to be
evalf(LambertW(1));
The graph of the solution curve is shown below.
plot(exp(LambertW(1)*t), t=0..5);
In applied mathematics and science, few differential equations admit closed form solutions, and often we need to employ numerical methods. The next example demonstrates the usage of the dsolve command for a delay differential equation.
A Delayed Action Oscillator for El Niño
Suarez and Schopf proposed a simple nonlinear model for the El Niño phenomenon. Their model relies on the existence of a strong positive feedback in the coupled ocean-atmosphere system and on some unspecified nonlinear mechanism invoked to limit the growth of unstable perturbations. Its key element is the inclusion of the effects of equatorially trapped oceanic waves propagating in a closed basin through a time delayed term. Let denote the temperature anomaly, that is, the deviation from a suitably defined long term average temperature, the time derivative is
,
where is the nondimensional delay (wave transit time) and measures the influence of the returning signal relative to that of the local feedback.
with(plots):
Dde1 := diff(T(t), t) = T(t) - T(t)^3 - a*T(t-d);
Let us assign the values for and to be 0.7 and 3, respectively.
Dde2 := eval(Dde1, {a = 0.7, d = 3});
Using an initial condition of , we solve the differential equation.
Soln2 := dsolve({Dde2, T(0) = 0.55}, T(t), numeric);
p1 := odeplot(Soln2, t = 0..100);
We try a different value for the delay.
Soln3 := dsolve({eval(Dde1, {a = 0.75, d = 6}), T(0) = 0.55}, T(t), numeric);
odeplot(Soln3, t = 0..100);
We observe that solutions of the delay differential equation are self-sustained oscillations. The reader can verify the conclusion of Suarez and Schopf that "the period of the oscillation is at least twice as long as the assumed delay." However, when the delay drops to a certain value, the oscillation ceases, as seen in the following.
Soln4 := dsolve({eval(Dde1, {a = 0.75, d = 1.5}), T(0) = 0.55}, T(t), numeric);
odeplot(Soln4, t = 0..100);
The delay is the time required for equatorially trapped oceanic waves to propagating across the Pacific Ocean. The existence of a threshold delay time offers a possible explanation for the lack of El Niño phenomenon in smaller water bodies such as the Atlantic and Indian Oceans.
If we consider the effect of global warming, with a rate of , the modified delay differential equation becomes:
Dde5 := diff(T(t), t) = T(t) - T(t)^3 - a*T(t-d) + b;
Let us try .
Soln6 := dsolve({eval(Dde5, {a = 0.7, d = 3, b = 0.05}), T(0) = 0.55}, T(t), numeric);
p2 := odeplot(Soln6, t = 0..100, color = blue);
display({p1, p2});
The above plot shows the comparison between the model with global warming (the blue curve, with ) and the original model. We see no noticeable difference in the steady state solutions. However, when the global warming rate gets too large, the oscillatory solution disappears, as seen below.
Soln7 := dsolve({eval(Dde5, {a = 0.7, d = 3, b = 0.5}), T(0) = 0.55}, T(t), numeric);
odeplot(Soln7, t = 0..100);
Drastic qualitative changes in the solution curves as one smoothly varies the parameter is called "bifurcation." In conclusion, Maple's newly added option to solve delay differential equations numerically allows users to investigate the delayed action oscillator for El Niño and investigate bifurcation due to the time delay and the global warming rate, and possibly other factors.
References
1. Ian Boutle, Richard H. S. Taylor, and Rudolf A. Römer, El Niño and the delay action oscillator, American Journal of Physics, 75, 15--24 (2007).
2. Nicholas J. Higham, The Princeton Companion to Applied Mathematics, Princeton University Press, 2015.
3. Max J. Suarez and Paul S. Schopf, A Delayed Action Oscillator for ENSO, Journal of the Atmospheric Sciences, 45, 3283--3287 (1988).