> Let me rewrite the expression you wrote as > > J*v = (F(u+epsilon*v)-F(u)) / epsilon > > Where epsilon is a small perturbation and F(u) is the nonlinear residual > function and J is the jacobian matrix of the nonlinear system. The above > formula compute the action of a Jacobian on a given vector, or more > specifically the Krylov vector if you are using say GMRES or CG to solve > your system.
That's absolutely true. To expound a little more, by definition F = F(U) J = dF/dU J*v = [dF/dU]*v (1) which is simply the directional derivative of F (in the direction of v. As specified above, you can approximate this matrix-vector product as J*v ~ (F(u+epsilon*v)-F(u)) / epsilon (2) In Krylov subspace methods the operation (1) occurs repeatedly at each iteration. So, if you are willing to approximate it, you can use (2) instead. There are a couple of reasons why you might want to do this: - obviously, this finite difference of vectors eliminates the need to store J, hence "matrix free" - not to whine, but sometimes computing an accurate J is *hard* and/or computationally intensive. It can be error-prone. Using (2) you get Now of course any hard problem will need preconditioning in the Krylov solver. You can accomplish this in several ways... One would be to build an approximate Jacobian and use ILU or the like. Now you have one matrix instead of two. Also, this could be a much simpler matrix (block-diagonal for certain cases), thus limiting the requirements. But there are other options too, allowing for a matrix free preconditioner as well. These include Gauss-Siedel, Geometric multigrid, another Krylov technique, etc... -Ben ------------------------------------------------------------------------- This SF.net email is sponsored by: Microsoft Defy all challenges. Microsoft(R) Visual Studio 2008. http://clk.atdmt.com/MRT/go/vse0120000070mrt/direct/01/ _______________________________________________ Libmesh-users mailing list [email protected] https://lists.sourceforge.net/lists/listinfo/libmesh-users
