On Finding the Smallest Generalized Eigenpair Using Markov Chain Monte Carlo Algorithm ()
1. Introduction
Monte Carlo and quasi Monte Carlo methods comprise that branch of experimental mathematics which is concerned with experiments on random numbers. In the last decade they have found extensive use in the fields of operational research and nuclear physics, where there are a variety of problems beyond the available resources of theoretical mathematics [1,2]. The problem of calculating the largest or smallest generalized eigenvalue problem is one of the most important problems in science and engineering. This problem arises naturally in many applications. Mathematically it is a generalization of the symmetric eigenvalue problem and it can be reduced to an equivalent symmetric eigenvalue problem. Recently, evaluating the smallest and largest eigenpair of large scale matrices using Monte Carlo and quasi Monte Carlo methods has been studied [3-5].
Let are real symmetric matrices and the matrix B is a positive definite matrix. Consider the problem of evaluating the smallest eigenpair of the pencil (A, B) i.e. the values for which
(1)
A generalized eigenvalue problem (1) is said to be symmetric positive definite (S/PD) if A is symmetric and B is positive definite.
2. Markov Chain Monte Carlo Algorithm
Suppose that the matrix and two vectors are given. Consider the following Markov chain with length
(2)
where for,. The statistical nature of constructing the chain (2) follow as
(3)
where and show the probability of starting chain at and transition probability from state to, respectively. In other words, we have
Now, define the following random variable
where, and
Theorem 1. Consider the linear system. Let us the nonsingular matrix, such that , then the system can be presented in the following form
(4)
where f = Mb. Then under condition, the random variable is unbiased estimator, i.e.
Proof [3].
Suppose that is the iterative solution of the following recursion relation with Now, we consider the random variable
then
By simulating N random paths with length i:
we can write
The Monte Carlo estimation can be evaluated by, which is an approximation of
Now, consider the following choice for the initial density vector and the transition probability matrix leads to an Almost Optimal Monte Carlo (AOMC) algorithm.
Theorem 2. Using the above choice and the variance of the unbiased estimator for obtaining the inverse matrix is minimized.
3. Inverse Monte Carlo Iterative (IMCI)
Inverse Monte Carlo iterative algorithm can be applied when A is a nonsingular matrix. In this method, we calculate the following functional in each.
It is more efficient that we first evaluate the inverse matrix using the Monte Carlo algorithm [3,4]. The algorithm can be realized as follows
The Proposed Algorithm
1) Inputs,.
2) Starting from initial vector.
3) For
4) Using global algorithm [4], Calculate the sequence of Monte Carlo iterations by solving the following System
Set
5) Output: Smallest eigenvector and the corresponding eigenvector.
Table 1. Computing the smallest eigenpair for different matrices using MCMC algorithm.
Figure 1. Relative error based on number of Markov chains.
6) End of Algorithm.
4. Numerical Results
In this section, the results listed in Table 1 are relative errors obtained by MCMC. In all of computations we used the MATLAB function eig (A, B) results in contrast to MCMC algorithm results.
5. Conclusion
A new numerical algorithm based on the inverse Markov chain Monte Carlo algorithm to calculate the smallest generalized eigenpair proposed in this paper. Also, we saw that by increasing the number of Markov chains the relative error is decreased (Figure 1).
6. Acknowledgements
The author is very grateful to two anonymous referees for giving very useful and valuable comments.