Simulation of bimolecular reactions: Numerical challenges with the graph Laplacian

An important framework for modelling and simulation of chemical reactions is a Markov process sometimes known as a master equation. Explicit solutions of master equations are rare; in general the explicit solution of the governing master equation for a bimolecular reaction remains an open question. We show that a solution is possible in special cases. One method of solution is diagonalization. The crucial class of matrices that describe this family of models are non-symmetric graph Laplacians. We illustrate how standard numerical algorithms for finding eigenvalues fail for the non-symmetric graph Laplacians that arise in master equations for models of chemical kinetics. We propose a novel way to explore the pseudospectra of the non-symmetric graph Laplacians that arise in this class of applications, and illustrate our proposal by Monte Carlo. Finally, we apply the Magnus expansion, which provides a method of simulation when rates change in time. Again the graph Laplacian structure presents some unique issues: standard numerical methods of more than second-order fail to preserve positivity. We therefore propose a method that achieves fourth-order accuracy, and maintain positivity. References A. Basak, E. Paquette, and O. Zeitouni. Regularization of non-normal matrices by gaussian noise–-the banded toeplitz and twisted toeplitz cases. In Forum of Mathematics, Sigma, volume 7. Cambridge University Press, 2019. doi:10.1017/fms.2018.29. S. Blanes, F. Casas, J. A. Oteo, and J. Ros. The magnus expansion and some of its applications. Phys. Rep., 470(5-6):151–238, 2009. doi:10.1016/j.physrep.2008.11.001. B. A. Earnshaw and J. P. Keener. Invariant manifolds of binomial-like nonautonomous master equations. SIAM J. Appl. Dyn. Sys., 9(2):568–588, 2010. doi10.1137/090759689. J. Gunawardena. A linear framework for time-scale separation in nonlinear biochemical systems. PloS One, 7(5):e36321, 2012. doi:10.1371/journal.pone.0036321. A. Iserles and S. MacNamara. Applications of magnus expansions and pseudospectra to markov processes. Euro. J. Appl. Math., 30(2):400–425, 2019. doi:10.1017/S0956792518000177. S. MacNamara. Cauchy integrals for computational solutions of master equations. ANZIAM Journal, 56:32–51, 2015. doi:10.21914/anziamj.v56i0.9345. S. MacNamara, A. M. Bersani, K. Burrage, and R. B. Sidje. Stochastic chemical kinetics and the total quasi-steady-state assumption: Application to the stochastic simulation algorithm and chemical master equation. J. Chem. Phys., 129:095105, 2008. doi:10.1063/1.2971036. S. MacNamara and K. Burrage. Stochastic modeling of naive T cell homeostasis for competing clonotypes via the master equation. SIAM Multiscale Model. Sim., 8(4):1325–1347, 2010. S. MacNamara, K. Burrage, and R. B. Sidje. Multiscale modeling of chemical kinetics via the master equation. SIAM Multiscale Model. and Sim., 6(4):1146–1168, 2008. doi:10.1137/060678154. S. MacNamara, Wi. McLean, and K. Burrage. Wider contours and adaptive contours, pages 79–98. Springer International Publishing, 2019. doi:10.1007/978-3-030-04161-8_7. M. J. Shon. Trapping and manipulating single molecules of DNA. PhD thesis, Harvard University, 2014. http://nrs.harvard.edu/urn-3:HUL.InstRepos:11744428. M. J. Shon and A. E. Cohen. Mass action at the single-molecule level. J. Am. Chem. Soc., 134(35):14618–14623, 2012. doi:10.1021/ja3062425. C. Timm. Random transition-rate matrices for the master equation. Phys. Rev. E, 80(2):021140, 2009. doi:10.1103/PhysRevE.80.021140. L. N. Trefethen and M. Embree. Spectra and pseudospectra: The behavior of nonnormal matrices and operators. Princeton University Press, 2005. https://press.princeton.edu/books/hardcover/9780691119465/spectra-and-pseudospectra.

Continuous spectra of the Batchelor vortex

The spectra of the Batchelor vortex are obtained by discretizing its linearized evolution operator using a modified Chebyshev polynomial approximation at a Reynolds number of 1000 and zero azimuthal wavenumber. Three types of eigenmodes are identified from the spectra: discrete modes, potential modes and free-stream modes. The discrete modes have been extensively documented but the last two modes have received little attention. A convergence study of the spectra and pseudospectra supports the classification that discrete modes correspond to discrete spectra while the other two modes correspond to continuous spectra. Free-stream modes have finite amplitude in the far field whilst potential modes decay to zero in the far field. The free-stream modes are therefore a limiting form of the potential modes when the radial decay rate of velocity components reduces to zero. The radial form of the free-stream modes with axial and radial wavenumbers is investigated and the penetration of the free-stream mode into the vortex core highlights the possibility of interaction between the potential region and the vortex core. A wavepacket pseudomode study confirms the existence of continuous spectra and predicts the locations and radial wavenumbers of the eigenmodes. The pseudomodes corresponding to the potential modes are observed to be in the form of one or two wavepackets while the free-stream modes are not observed to be in the form of wavepackets.

Infinite-dimensional numerical linear algebra: theory and applications

We present a new method for computing spectra and pseudospectra of bounded operators on separable Hilbert spaces. The core in this theory is a generalization of the pseudospectrum called the n -pseudospectrum.

Modal and non-modal linear stability of the plane Bingham–Poiseuille flow

The receptivity problem of plane Bingham–Poiseuille flow with respect to weak perturbations is addressed. The relevance of this study is highlighted by the linear stability analysis results (spectra and pseudospectra). The first part of the present paper thus deals with the classical normal-mode approach in which the resulting eigenvalue problem is solved using the Chebychev collocation method. Within the range of parameters considered, the Poiseuille flow of Bingham fluid is found to be linearly stable. The second part investigates the most amplified perturbations using the non-modal approach. At a very low Bingham number (B ≪ 1), the optimal disturbance consists of almost streamwise vortices, whereas at moderate or large B the optimal disturbance becomes oblique. The evolution of the obliqueness as function of B is determined. The linear analysis presented also indicates, as a first stage of a theoretical investigation, the principal challenges of a more complete nonlinear study.