Pseudospectral discretization of delay differential equations in sun-star formulation: Results and conjectures

Odo Diekmann; Francesca Scarabel; Rossana Vermiglio

doi:10.3934/dcdss.2020196

Article Contents

2020, Volume 13, Issue 9: 2575-2602. Doi: 10.3934/dcdss.2020196

This issue Previous Article Controllability analysis of nonlinear fractional order differential systems with state delay and non-instantaneous impulsive effects Next Article Functionally-fitted block

$ \theta $

-methods for ordinary differential equations

Pseudospectral discretization of delay differential equations in sun-star formulation: Results and conjectures

1.
Department of Mathematics, Utrecht University, P.O. Box 80010, 3508 TA Utrecht, The Netherlands
2.
Department of Mathematics and Statistics, University of Helsinki, P.O. 68 (Pietari Kalmin katu 5), FI-00014 Helsinki, Finland
3.
Department of Mathematics and Statistics, York University, 4700 Keele St, Toronto, ON M3J 1P3, Canada
4.
Department of Mathematics, Computer Science and Physics, University of Udine, via delle Scienze 206, I-33100 Udine, Italy

^* Corresponding author: Rossana Vermiglio
^* Corresponding author: Rossana Vermiglio

Received: January 2019

Revised: May 2019

Early access: December 2019

Published: June 2020

The research of the second author was supported by Domast (Doctoral Programme in Mathematics and Statistics, University of Helsinki), and by the Centre of Excellence in Analysis and Dynamics Research, Academy of Finland. F.S. and R.V. are members of the INdAM Research group GNCS, and of CDLab (Computational Dynamics Laboratory), Department of Mathematics, Computer Science and Physics, University of Udine.

Abstract / Introduction Full Text(HTML) Figure(14) / Table(2) Related Papers Cited by

Abstract

In this paper we study the pseudospectral approximation of delay differential equations formulated as abstract differential equations in the $ \odot* $-space. This formalism also allows us to define rigorously the abstract variation-of-constants formula, where the $ \odot* $-shift operator plays a fundamental role. By applying the pseudospectral discretization technique we derive a system of ordinary differential equations, whose dynamics can be efficiently analyzed by existing bifurcation tools. To better understand to what extent the resulting finite-dimensional system "mimics" the dynamics of the original infinite-dimensional one, we study the pseudospectral approximations of the $ \odot* $-shift operator and of the $ \odot* $-generator in the supremum norm, which is the natural choice for delay differential equations, when the discretization parameter increases. In this context there are still open questions. We collect the most relevant results from the literature and we present some conjectures, supported by various numerical experiments, to illustrate the behavior w.r.t. the discretization parameter and to indicate the direction of ongoing and future research.

Keywords:

Mathematics Subject Classification: Primary: 65F35, 65L03, 65L07, 65M70, 65P99, 65Q20, 37M20X, 34K07, 34K20, 41A10; Secondary: 34K08, 34K28, 37K45.

Citation:

Full Text(HTML)

Figure 1. Left: plot of $ \sigma(D_M) $ for different $ M $, together with the line $ \mathrm{Re} \lambda = \log \epsilon \approx -36.0437 $ connected with the instability phenomenon studied in [63]. Right: plot of $ \frac{-\alpha(D_M)}{\log(M)} $, versus $ M $

Download: Full-size image PowerPoint slide

Figure 2. Error $ \| \psi_M - \psi\| $ for the resolvent operators applied to $ \beta = 0 $, $ \varphi = 1 $, versus $ M $. Left: $ \lambda = 1 $; right: $ \lambda = 10 $. Note the spectral accuracy and that the convergence is slower when $ \lambda $ has larger modulus

Download: Full-size image PowerPoint slide

Figure 3. Left: plot of $ \log \|{\rm e}^{D_M t}\| $ for $ t\in [0, 2] $ with different values of $ M $. Notice divergence in $ [0, 1] $ for increasing $ M $. Right: plot of $ \max_{t\in I} \frac{\log \|{\rm e}^{D_M t}\|}{\log(M)} $ versus $ M, $ for $ I $ specified in the legend. Notice that the convergence gains one order at every time interval

Download: Full-size image PowerPoint slide

Figure 4. $ \|{\rm e}^{D_M t}\| $ computed with different routines: built-in Matlab function $text{expm}$ and $\text{norm(., inf)}$ versus maximum over $ n = 100 $ random initial vectors of the solution of the ODE system. Left: $ \|{\rm e}^{D_M t}\| $ versus time. Right: $ \max_{t \in [0, 2]} \|{\rm e}^{D_M t}\| $ versus $ M $. The dotted line is the reference line $ \log(M) $. The fact that the effective norm of the exponential matrix diverges, while the norm computed by selecting a random set of vectors is uniformly bounded, suggests that the "bad" behavior of the norm is due to a small set of vectors

Download: Full-size image PowerPoint slide

Figure 5. Left: plot of the function $ \sum\limits_{i = 1}^M |\ell_i(\theta)| $ for $ M = 8 $. Right: plot of $ \|{\rm e}^{D_M t}\| $ (blue) and $ C_M(t) $ in (45) (red) versus time for $ M = 8 $

Download: Full-size image PowerPoint slide

Figure 6. Left: plot of $ \psi(\theta) = H(\theta) $ (black) and its interpolating polynomial for $ M = 8 $ (red). Right: $ \log \log $ plot of the error $ \|{\rm e}^{D_M t}\psi(\Theta_M)-R_M \mathcal{H}_{0}(t)\psi\| $ for $ t = 0.5, 1.5, 2.5 $ versus $ M $. Note the uniform bound for $ t = 0.5 $ and the convergence for $ t>1 $. The dotted lines are the reference lines $ M^{-k} $, $ k = 1, \dots, 4 $

Download: Full-size image PowerPoint slide

Figure 7. Same as Figure 6 for $ \psi(\theta) = 0.5-|\theta+0.5| $. Note that $ \psi \in Y $ and $ \psi^{(1)} $ has bounded variation. The convergence rate for $ t = 0.25, 0.5 $ is $ O(M^{-1}) $ (right)

Download: Full-size image PowerPoint slide

Figure 8. Same as Figure 6 for $ \psi(\theta) = {\rm e}^{\theta^2}-1 $. Note that the interpolating polynomial is indistinguishable from the function (left) and that $ \psi^{(2)} $ has bounded variation. The convergence rate for $ t = 0.5 $ is $ O(M^{-2}) $ (right)

Download: Full-size image PowerPoint slide

Figure 9. Same as Figure 6 for $ \psi(\theta) = -\theta^3. $ Note that the interpolating polynomial is indistinguishable from the function (left) and that $ \psi^{(3)} $ has bounded variation. The convergence rate for $ t = 0.5 $ is $ O(M^{-3}) $ (right)

Download: Full-size image PowerPoint slide

Figure 10. Same as Figure 6 for $ \psi(\theta) = {\rm e}^{-1/\theta^2}. $ Note that the interpolating polynomial is indistinguishable from the function (left) and the spectral convergence for $ t = 0.5, 1.5, 2.5 $ (right)

Download: Full-size image PowerPoint slide

Figure 11. Same as Figure 6 for $ \psi(\theta) = \sin\frac{1}{\theta} $, if $ \theta<0 $, $ \psi(0) = 0 $. Note there is no convergence for $ t = 0.5 $

Download: Full-size image PowerPoint slide

Figure 12. Same as Figure 6 for $ \psi(\theta) = \sin(6 \theta)+sign(\sin(\theta+{\rm e}^{2\theta})) $ [60,pag.10]. Note there is no convergence for $ t = 0.5 $

Download: Full-size image PowerPoint slide

Figure 13. Left: plot of the function $ z $ (saw function) Right: error $ \|\int_0^t e^{(t-s)D_M} z(s) {\rm d} s- w(t)(\Theta_M) \| $ for $ t = 0.5, 1, 1.5 $ versus $ M $. The dotted line is the reference line $ \log(M)/M $

Download: Full-size image PowerPoint slide

Figure 14. Same as Figures 13 for the function $ z $ built so that it is continuously differentiable with jumps in the second derivative at $ t = 0.5,1,1.5 $

Download: Full-size image PowerPoint slide

Table 1. $ \|D_M\| $ and estimation of the order $ \log_2\frac{\|D_{2M}\|}{\|D_M\|} $ varying $ M. $

$ M $	$ \\|D_M\\| $	order	$ \\|D_M\\|/M^2 $
4	31		1.9375
8	127	2.0345	1.9844
16	511	2.0085	1.9961
32	2047	2.0021	1.9990
64	8191	2.0005	1.9998
128	32767	2.0001	1.9999
256	131071	2.0000	2.0000

| Show Table

DownLoad: CSV

Table 2. $ \mu(D_M) $ and estimation of the order $ \log_2\frac{\mu(D_{2M})}{\mu(D_M)} $ varying $ M. $

$ M $	$ \mu(D_M) $	order	$ \mu(D_M)/M^2 $
4	6.6569e+00		4.1605e-01
8	3.8921e+01	2.5476	6.0814e-01
16	1.6698e+02	2.1011	6.5227e-01
32	6.7900e+02	2.0237	6.6308e-01
64	2.7270e+03	2.0058	6.6577e-01
128	1.0919e+04	2.0015	6.6644e-01
256	4.3687e+04	2.0004	6.6661e-01

| Show Table

DownLoad: CSV

Related Papers

Cited by

References

[1]	A. H. Al-Mohy and N. J. Higham, A new scaling and squaring algorithm for the matrix exponential, SIAM J. Matrix Anal. Appl., 31 (2009), 970-989. doi: 10.1137/09074721X.
[2]	A. Andò, D. Breda, L. Davide, S. Maset, F. Scarabel and R. Vermiglio, 15 years or so of pseudospectral collocation methods for stability and bifurcation of delay equations, in Advances on Delays and Dynamics, Springer, New York, 2019.
[3]	C. T. H. Baker, The Numerical Treatment of Integral Equations, Monographs on Numerical Analysis, Clarendon Press, Oxford, 1977.
[4]	C. Baker, Numerical analysis of Volterra functional and integral equations, in The State of the Art in Numerical Analysis, Math. Appl. Conf. Ser. New Ser., 63, Oxford Univ. Press, New York, 1997.
[5]	H. T. Banks, J. A. Burns and E. M. Cliff, Spline-based approximation methods for control and identification of hereditary systems, in International Symposium on Systems Optimization and Analysis, Lecture Notes in Control and Information Sci., 14, Springer, Berlin-New York, 1979, 314–320. doi: 10.1007/BFb0002662.
[6]	A. Batkai and S. Piazzera, Semigroup for Delay Equations, Research Notes in Mathematics, 10, A K Peters, Ltd., Wellesley, MA, 2005.
[7]	A. Bellen and M. Zennaro, Numerical Methods for Delay Differential Equations, Numerical Mathematics and Scientific Computation, Oxford University Press, 2003.
[8]	D. Breda, O. Diekmann, M. Gyllenberg, F. Scarabel and R. Vermiglio, Pseudospectral discretization of nonlinear delay equations: New prospects for numerical bifurcation analysis, SIAM J. Appl. Dyn. Syst., 15 (2016), 1-23. doi: 10.1137/15M1040931.
[9]	D. Breda, O. Diekmann, D. Liessi and F. Scarabel, Numerical bifurcation analysis of a class of nonlinear renewal equations, Electron. J. Qual. Theory Differ. Equ., 2016, 1–24. doi: 10.14232/ejqtde.2016.1.65.
[10]	D. Breda, P. Getto, J. Sánchez Sanz and R. Vermiglio, Computing the eigenvalues of realistic Daphnia models by pseudospectral methods, SIAM J. Sci. Comput., 37 (2015), A2607–A2629. doi: 10.1137/15M1016710.
[11]	D. Breda, S. Maset and R. Vermiglio, Pseudospectral differencing methods for characteristic roots of delay differential equations, SIAM J. Sci. Comput., 27 (2005), 482-495. doi: 10.1137/030601600.
[12]	D. Breda, S. Maset and R. Vermiglio, Pseudospectral approximation of eigenvalues of derivative operators with non-local boundary conditions, Appl. Numer. Math., 56 (2006), 318-331. doi: 10.1016/j.apnum.2005.04.011.
[13]	D. Breda, S. Maset and R. Vermiglio, Approximation of eigenvalues of evolution operators for linear retarded functional differential equations, SIAM J. Numer. Anal., 50 (2012), 1456-1483. doi: 10.1137/100815505.
[14]	D. Breda, S. Maset and R. Vermiglio, Stability of Linear Delay Differential Equations. A Numerical Approach with MATLAB, SpringerBriefs in Electrical and Computer Engineering, Springer, New York, 2015. doi: 10.1007/978-1-4939-2107-2.
[15]	H. Brunner, Collocation Methods for Volterra Integral and Related Functional Differential Equations, Cambridge Monographs on Applied and Computational Mathematics, 15, Cambridge University Press, Cambridge, 2004. doi: 10.1017/CBO9780511543234.
[16]	H. Brunner, Volterra Integral Equations. An Introduction to Theory and Applications, Cambridge Monographs on Applied and Computational Mathematics, 30, Cambridge University Press, Cambridge, 2017. doi: 10.1017/9781316162491.
[17]	H. Brunner and P. van der Houwen, The Numerical Solution of Volterra Equations, CWI Monographs, 3, North-Holland Publishing Co., Amsterdam, 1986.
[18]	C. Canuto, Y. M. Hussaini, A. Quarteroni and T. A. J. Zang, Spectral Methods in Fluid Dynamics, Springer Series in Computational Physics, Springer-Verlag, New York, 1988. doi: 10.1007/978-3-642-84108-8.
[19]	M. D. Chekroun, M. Ghil, H. Liu and S. Wang, Low-dimensional Galerkin approximations of nonlinear delay differential equations, Discrete Contin. Dyn. Syst., 36 (2016), 4133-4177. doi: 10.3934/dcds.2016.36.4133.
[20]	P. Clément, O. Diekmann, M. Gyllenberg, H. Heijmans and H. Thieme, Perturbation theory for dual semigroups. Ⅲ: Nonlinear Lipschitz continuous perturbations in the sun-reflexive case, in Volterra Integrodifferential Equations in Banach Spaces and Applications, Pitman Res. Notes Math. Ser., 190, Longman Sci. Tech., Harlow, 1989, 67–89.
[21]	P. J. Davis, Interpolation and Approximation, Dover Publications, Inc. New York, 1975.
[22]	A. Dhooge, W. Govaerts, Y. A. Kuznetsov, H. G. E. Meijer and B. Sautois, New features of the software MatCont for bifurcation analysis of dynamical systems, Math. Comput. Model. Dyn. Syst., 14 (2008), 147-175. doi: 10.1080/13873950701742754.
[23]	O. Diekmann, P. Getto and M. Gyllenberg, Stability and bifurcation analysis of Volterra functional equations in the light of suns and stars, SIAM J. Math. Anal., 39 (2007/08), 1023-1069. doi: 10.1137/060659211.
[24]	O. Diekmann and M. Gyllenberg, Equations with infinite delay: Blending the abstract and the concrete, J. Differential Equations, 252 (2012), 819-851. doi: 10.1016/j.jde.2011.09.038.
[25]	O. Diekmann, S. A. van Gils, S. M. Verduyn Lunel and H.-O. Walther, Delay Equations. Functional, Complex, and Nonlinear Analysis, Applied Mathematical Sciences, 110, Springer-Verlag, New York, 1995. doi: 10.1007/978-1-4612-4206-2.
[26]	M. Dubiner, Asymptotic analysis of spectral methods, J. Sci. Comput., 2 (1987), 3-31. doi: 10.1007/BF01061510.
[27]	K. Engelborghs, T. Luzyanina, G. Samaey, D. Roose and K. Verheyden, DDE-BIFTOOL: A MATLAB package for bifurcation analysis of delay differential equations., Available from: http://ddebiftool.sourceforge.net.
[28]	B. Fornberg, A Practical Guide to Pseudospectral Methods, Cambridge Monographs on Applied and Computational Mathematics, 1, Cambridge University Press, Cambridge, 1996. doi: 10.1017/CBO9780511626357.
[29]	D. Funaro, A preconditioning matrix for the Chebyshev differencing operator, SIAM J. Numer. Anal., 24 (1987), 1024-1031. doi: 10.1137/0724067.
[30]	D. Funaro, Some results about the spectrum of the Chebyshev differencing operator, in Numerical Approximation of Partial Differential Equations, North-Holland Math. Stud., 133, North-Holland, Amsterdam, 1987, 271–284. doi: 10.1016/S0304-0208(08)71738-9.
[31]	D. Funaro, Polynomial Approximation of Differential Equations, Lecture Notes in Physics. New Series M: Monographs, 8, Springer-Verlag, Berlin, 1992. doi: 10.1007/978-3-540-46783-0.
[32]	P. Getto, M. Gyllenberg, Y. Nakata and F. Scarabel, Stability analysis of a state-dependent delay differential equation for cell maturation: Analytical and numerical methods, J. Math. Biol., 79 (2019), 281-328. doi: 10.1007/s00285-019-01357-0.
[33]	D. Gottlieb, The stability of pseudospectral-Chebyshev methods, Math. Comp., 36 (1981), 107-118. doi: 10.1090/S0025-5718-1981-0595045-1.
[34]	D. Gottlieb and L. Lustman, The spectrum of the Chebyshev collocation operator for the heat equation, SIAM J. Numer. Anal., 20 (1983), 909-921. doi: 10.1137/0720063.
[35]	D. Gottlieb and S. A. Orszag, Numerical Analysis of Spectral Methods: Theory and Applications, CBMS-NSF Regional Conference Series in Applied Mathematics, SIAM, Philadelphia, PA, 1977. doi: 10.1137/1.9781611970425.
[36]	D. Gottlieb and E. Turkel, Topics in spectral methods, in Numerical Methods in Fluid Dynamics, Lecture Notes in Math., 1127, Springer, Berlin, 1985, 115–155. doi: 10.1007/BFb0074530.
[37]	M. Gyllenberg, F. Scarabel and R. Vermiglio, Equations with infinite delay: Numerical bifurcation analysis via pseudospectral discretization, Appl. Math. Comput., 333 (2018), 490-505. doi: 10.1016/j.amc.2018.03.104.
[38]	E. Hairer, S. P. Nørsett and G. Wanner, Solving Ordinary Differential Equations. I: Nonstiff Problems, Springer Series in Computational Mathematics, 8, Springer-Verlag, Berlin, 1993. doi: 10.1007/978-3-540-78862-1.
[39]	J. K. Hale and S. M. Verduyn Lunel, Introduction to Functional Differential Equations, Applied Mathematical Sciences, 99, Springer-Verlag, New York, 1993. doi: 10.1007/978-1-4612-4342-7.
[40]	J. S. Hesthaven, S. Gottlieb and D. Gottlieb, Spectral Methods for Time-Dependent Problems, Cambridge Monographs on Applied and Computational Mathematics, 21, Cambridge University Press, Cambridge, 2007. doi: 10.1017/CBO9780511618352.
[41]	N. J. Higham, The scaling and squaring method for the matrix exponential revisited, SIAM J. Matrix Anal. Appl., 26 (2005), 1179-1193. doi: 10.1137/04061101X.
[42]	N. J. Higham and A. H. Al-Mohy, Computing matrix functions, Acta Numer., 19 (2010), 159-208. doi: 10.1017/S0962492910000036.
[43]	K. Ito and F. Kappel, A uniformly differentiable approximation scheme for delay systems using splines, Appl. Math. Optim., 23 (1991), 217-262. doi: 10.1007/BF01442400.
[44]	K. Ito and F. Kappel, Two families of approximation schemes for delay systems, Results Math., 21 (1992), 93-137. doi: 10.1007/BF03323074.
[45]	Z. Jackiewicz and B. D. Welfert, Stability of Gauss-Radau pseudospectral approximations of the one-dimensional wave equation, J. Sci. Comput., 18 (2003), 287-313. doi: 10.1023/A:1021121008091.
[46]	R. Kress, Linear Integral Equations, Applied Mathematical Sciences, 82, Springer-Verlag, Berlin, 1989. doi: 10.1007/978-3-642-97146-4.
[47]	D. Liessi, Pseudospectral Methods for Stability of Periodic Solutions of Delay Models, Ph.D thesis, University of Udine, 2018.
[48]	S. Maset, The collocation method in the numerical solution of boundary value problems for neutral functional differential equations. Part Ⅰ: Convergence results, SIAM J. Numer. Anal., 53 (2015), 2771-2793. doi: 10.1137/130935550.
[49]	S. Maset, The collocation method in the numerical solution of boundary value problems for neutral functional differential equations. Part Ⅱ: Differential equations with deviating arguments, SIAM J. Numer. Anal., 53 (2015), 2794-2821. doi: 10.1137/140979022.
[50]	S. Maset, An abstract framework in the numerical solution of boundary value problems for neutral functional differential equations, Numer. Math., 133 (2016), 525-555. doi: 10.1007/s00211-015-0754-1.
[51]	G. Mastroianni and G. V. Milovanović, Interpolation Processes. Basic Theory and Applications, Springer Monographs in Mathematics, Springer-Verlag, Berlin, 2008. doi: 10.1007/978-3-540-68349-0.
[52]	W. Michiels, I. Boussaada and S.-I. Niculescu, An explicit formula for the splitting of multiple eigenvalues for nonlinear eigenvalue problems and connections with the linearization for the delay eigenvalue problem, SIAM J. Matrix Anal. Appl., 38 (2017), 599-620. doi: 10.1137/16M107774X.
[53]	W. Michiels and S.-I. Niculescu, Stability, Control, and Computation for Time-Delay Systems. An Eigenvalue-Based Approach, Advances in Design and Control, 27, SIAM, Philadelphia, PA, 2014. doi: 10.1137/1.9781611973631.
[54]	S. C. Reddy and L. N. Trefethen, Lax-stability of fully discrete spectral methods via stability regions and pseudo-eigenvalues, Comput. Methods Appl. Mech. Engrg., 80 (1990), 147-164. doi: 10.1016/0045-7825(90)90019-I.
[55]	T.-J. Rivlin, Chebyshev Polynomials. From Approximation Theory to Algebra and Number Theory, Pure and Applied Mathematics, 10, John Wiley & Sons, Inc., New York, 1990.
[56]	F. Scarabel, Capture the Past to Portray the Future: Numerical Bifurcation Analysis of Delay Equations, with a Focus on Population Dynamics, Ph.D thesis, University of Helsinki, 2018.
[57]	D. M. Sloan, On the norms of inverses of pseudospectral differentiation matrices, SIAM J. Numer. Anal., 42 (2004), 30-48. doi: 10.1137/S0036142902414542.
[58]	A. Solomonoff and E. Turkel, Global properties of pseudospectral methods, J. Comput. Phys., 81 (1989), 239-276. doi: 10.1016/0021-9991(89)90208-8.
[59]	R. Szalai, Knut., Available from: http://rs1909.github.io/knut/.
[60]	L. N. Trefethen, Approximation Theory and Approximation Practice, SIAM, Philadelphia, PA, 2013.
[61]	L. N. Trefethen, Inverse yogiisms, Notices Amer. Math. Soc., 63 (2016), 1281-1285. doi: 10.1090/noti1446.
[62]	L. N. Trefethen and M. Embree, Spectra and Pseudospectra: The Behavior of Nonnormal Matrices And Operators, Princeton University Press, Princeton, NJ, 2005.
[63]	L. N. Trefethen and M. R. Trummer, An instability phenomenon in spectral methods, SIAM J. Numer. Anal., 24 (1987), 1008-1023. doi: 10.1137/0724066.
[64]	J. Wang and F. Waleffe, The asymptotic eigenvalues of first-order spectral differentiation matrices, J. Appl. Math. Phys., 2 (2014), 176-188. doi: 10.4236/jamp.2014.25022.
[65]	G. Webb, Functional differential equations and nonlinear semigroup in $L_p$-spaces, J. Differential Equations, 20 (1976), 71-89. doi: 10.1016/0022-0396(76)90097-8.
[66]	Z. Wu and W. Michiels, Reliably computing all characteristic roots of delay differential equations in a given right half plane using a spectral method, J. Comput. Appl. Math., 236 (2012), 2499-2514. doi: 10.1016/j.cam.2011.12.009.