  In this paper, we develop a new extrapolation cascadic multigrid
(ECMG$_{jcg}$) method, which makes it possible to solve 3D elliptic boundary
value problems on rectangular domains of over 100 million unknowns on a desktop
computer in minutes. First, by combining Richardson extrapolation and
tri-quadratic Serendipity interpolation techniques, we introduce a new
extrapolation formula to provide a good initial guess for the iterative
solution on the next finer grid, which is a third order approximation to the
finite element (FE) solution. And the resulting large sparse linear system from
the FE discretization is then solved by the Jacobi-preconditioned Conjugate
Gradient (JCG) method. Additionally, instead of performing a fixed number of
iterations as cascadic multigrid (CMG) methods, a relative residual stopping
criterion is used in iterative solvers, which enables us to obtain conveniently
the numerical solution with the desired accuracy. Moreover, a simple Richardson
extrapolation is used to cheaply get a fourth order approximate solution on the
entire fine grid. Test results are reported to show that ECMG$_{jcg}$ has much
better efficiency compared to the classical MG methods. Since the initial guess
for the iterative solution is a quite good approximation to the FE solution,
numerical results show that only few number of iterations are required on the
finest grid for ECMG$_{jcg}$ with an appropriate tolerance of the relative
residual to achieve full second order accuracy, which is particularly important
when solving large systems of equations and can greatly reduce the
computational cost. It should be pointed out that when the tolerance becomes
smaller, ECMG$_{jcg}$ still needs only few iterations to obtain fourth order
extrapolated solution on each grid, except on the finest grid. Finally, we
present the reason why our ECMG algorithms are so highly efficient for solving
such problems.
