Numerical solution master equation

To facilitate the discussion it is helpful to specify three of the numerous meanings of the word state . We shall call a site any value of the stochastic variable X or n. We shall call a macrostate any value of the macroscopic variable . A time-dependent macrostate is a solution of the macroscopic equation (X.3.1), a stationary macrostate is a solution of (X.3.3). We shall call a mesostate any probability distribution P. A time-dependent meso-state is a solution of the master equation, the stationary mesostate is the time-independent solution PS(X). 

In this Section we introduce a stochastic alternative model for surface reactions. As an application we will focus on the formation of NH3 which is described below, equations (9.1.72) to (9.1.76). It is expected that these stochastic systems are well-suited for the description via master equations using the Markovian behaviour of the systems under study. In such a representation an infinite set of master equations for the distribution functions describing the state of the surface and of pairs of surface sites (and so on) arises. As it was told earlier, this set cannot be solved analytically and must be truncated at a certain level. The resulting equations can be solved exactly in a small region and can be connected to a mean-field solution for large distances from a reference point. This procedure is well-suited for the description of surface reaction systems which includes such elementary steps as adsorption, diffusion, reaction and desorption.The numerical part needs only a very small amount of computer time compared to MC or CA simulations. 

We remark that the simulation scheme for master equation dynamics has a number of attractive features when compared to quantum-classical Liouville dynamics. The solution of the master equation consists of two numerically simple parts. The first is the computation of the memory function which involves adiabatic evolution along mean surfaces. Once the transition rates are known as a function of the subsystem coordinates, the sequential short-time propagation algorithm may be used to evolve the observable or density. Since the dynamics is restricted to single adiabatic surfaces, no phase factors... 

Stochastic analysis presents an alternative avenue for dealing with the inherently probabilistic and discontinuous microscopic events that underlie macroscopic phenomena. Many processes of chemical and physical interest can be described as random Markov processes.1,2 Unfortunately, solution of a stochastic master equation can present an extremely difficult mathematical challenge for systems of even modest complexity. In response to this difficulty, Gillespie3-5 developed an approach employing numerical Monte Carlo... 

For linear systems, the differential equation for the jth cumulant function is linear and it involves terms up to the jth cumulant. The same procedure will be followed subsequently with other models to obtain analogous differential equations, which will be solved numerically if analytical solutions are not tractable. Historically, numerical methods were used to construct solutions to the master equations, but these solutions have pitfalls that include the need to approximate higher-order moments as a product of lower moments, and convergence issues [383]. What was needed was a general method that would solve this sort of problem, and that came with the stochastic simulation algorithm. 

UNIMOL Calculation of Rate Coefficients for Unimolecular and Recombination Reactions. R. G. Gilbert, T. Jordan, and S. C. Smith, Department of Theoretical Chemistry, Sydney, NSW 2006, Australia, 1990. FORTRAN computer code for calculating the pressure and temperature dependence of unimolecular and recombination (association) rate coefficients. Theory based on RRKM and numerical solution of the master equation. See Theory of Unimolecular and Recombination Reactions, by R. G. Gilbert and S. C. Smith, Blackwell Scientific Publications, Oxford, 1990. 

Diezemann and co-workers have recently presented numerical solutions of a master equation for the orientational relaxation of rotational correlation functions in glass-formers such as salol [172,173]. They have considered both globally and locally connected models, where in the latter case only states similar in energy to the starting point can be reached in one hop. They report that the main features of their results are qualitatively the same for these two extremes of connectivity. 

The master equation (4.3) can be solved numerically in various ways [6,7]. A symmetric G matrix guarantees that the solutions are of the form... 

The above approach, based on the solution of Eq. (4.3), is numerically superior and more accurate than one based on conventional Monte Carlo simulations. For comparison, Pitsianis et al. [59b] performed Monte Carlo simulations using 100,000 random walks on a Sierpinski gasket of 29,526 sites and obtained a value of 1.354 for dj (the exact value is 1.365). In the approach elaborated above, the value 1.367 was obtained by solving the stochastic master equation on a gasket of only 366 sites. 

The mobility, fi, can be determined from a (numerical) solution of the master equation representing hopping of charge carriers between lattice sites ... 

The master equation, however, can only be solved analytically for very simple systems such as the gas-phase reaction A—>B. The analysis of these systems typically requires numerical simulation of a lattice-based kinetic Monte Carlo model. The lattice gas model can then be used to formulate the respective transition probabilities in order to solve the master equationThe groups of both Zhdanov[ ° ° ] and Kreuzerl ° l have been instrumental in demonstrating the application of lattice gas models to solve adsorption and desorption processed from surfaces. Once a lattice model has been formulated there are three types of solution ... 

In Fig. 4.3, [2], we plot a cut through the stationary solution of the master equation for selected parameters of the Selkov model vs. the variable Y the dotted line is a numerical solution of the probability distribution and the solid line is that distribution calculated from (4.15-4.18) with the approximation described in this paragraph and the same parameters for the Selkov model. The approximation gives a reasonable estimate. A different impression is gathered from the plot shown in Fig. 4.4 A most probable fluctuational trajectory obtained from numerical integration of the stationary solution of the master... 

For multi-variable systems this approach is more difficult the determination of the stochastic potential requires sufficient measurements to determine rate coefficients and then the numerical solution of the stationary form of the master equation. Details of this procedure are described in Appendix A of [1]. 

A direct test of the master equation for systems in non-equilibrium stationary states comes from the measurements of concentration fluctuations such measurements have not been made yet. Some other tests of the master equation are possible based on the earlier sections in this chapter, where we can compare measurements of the stochastic potential with numerical solutions of the master equation (which requires knowledge of rate coefficients and the reaction mechanism of the system). 

Other comparisons of the results of the linearized Fokker Planck equation and the numerical solutions of the master equation are shown in Fig. 19.3. 

Zeron ES, SantiUan M (2010) Distributions for negative-feedback-regulated stochastic gene expression Dimension reduction and numerical solution of the chemical master equation. J Theor Biol 264(2) 377-385... 

The numerical solution of the master equation (111) is straightforward, either in the matrix forms of equations (112), (113) or by means of direct numerical integration of the coupled equations, given that the column matrix (vector) p of the level populations is of modest order. We summarize here some of the most important considerations and steps, leading finally to the fluence, intensity, and time dependent rate coefficient. For constant intensity, one has the exponential solution given by equation (113). If the relevant part of the rate coefficient matrix can be written as proportional to radiation intensity /(f) with intensity independent K/, one finds equations (128) and (129) ... 

---