Before an oil company is drilling for oil, they try to estimate the magnitude of the excess pressure in underground layers. Knowledge of the excess pressure can help to predict the presence of oil and gas reservoirs. A mathematical model for the prediction of the excess pressure in a geological time scale is based on conservation of mass and Darcy's law. This leads to a time-dependent diffusion equation, where the region also changes in time. In order to solve this diffusion equation, the finite element method is applied. As a consequence in each time-step a linear system of equations has to be solved.

In practical applications we are faced with large regions in a three-dimensional space and as a consequence a large number of finite elements is necessary. The matrix itself is sparse, but due to fill-in a direct method requires to much memory to fit in core. So, only iterative methods are acceptable candidates for the solution of the linear systems of equations.

Since these equations are symmetric a preconditioned conjugate gradient method (ICCG) is a natural candidate. Unfortunately an extra complication of the physical problem we are dealing with, is that the underground consists of layers with very large differences in permeability. A large contrast in coefficients usually leads to a very ill-conditioned system of equations to be solved. Since the convergence rate of ICCG depends on the condition of the matrix one may expect a slow convergence rate. An even more alarming phenomenon is that numerical results suggest that ICCG has reached a certain accuracy but that the actual accuracy is in fact orders worse. This means that due to ill-conditioning the standard termination criterion is no longer reliable.

An analysis of the problem shows that without preconditioning there are many small eigenvalues in the matrix, but using an ILU preconditioned matrix this number is reduced to the number of sandstone layers that not reach the earth surface. This suggests a way of solving the problems mentioned. The convergence and reliability of the termination criterion is considerably improved by projecting the solution in each step onto the orthogonal complement of the space spanned by the eigenvectors corresponding to the very small eigenvalues of the preconditioned matrix. The idea is that, assuming that the vectors are written as linear combination of eigenvectors, the components corresponding to these specific eigenvectors do not play a role anymore. As a result one may expect a much faster convergence and a reliable termination criterion. A clear disadvantage of this method is of course that one has to compute the specific eigenvectors. It appears that eigenvectors can be easily approximated, based on physical arguments. Furthermore it is shown, that even approximate eigenvectors lead to a fast convergence.