I am solving a coupled system of non-linear PDEs in 1D. Something like,

$$ u_t = F_1(u,v,w) \\ v_t = F_2(u,v,w) \\ w_t = F_3(u,v,w) $$

where each variable is a function of $x$ (the spatial dimension).

So for each point in space I wish to calculate the Jacobian of the vector,

$$ \boldsymbol{F} = \pmatrix{F_1(u,v,w) \\ F_2(u,v,w) \\ F_3(u,v,w)} $$


Should I use the PDEs themselves as a starting point to calculate the Jacobian, or is it better to use the discretised form of the equations?

I ask because the equations are very unstable without using exponential fitting, so the discretised form is quite non-trivial. My intuition is that Jacobian should follow from the discretised form as this seems more consistent. Is there a standard approach?

For example, if I were to start with the discretised form, I would simply partially differentiate each term with respect to the appropriate variable.


If you need Jacobian matrix information for a numerical method, you should calculate the Jacobian matrix of the discretized form of the equations, since that will be consistent with the discretized equations you are solving.

  • 1
    $\begingroup$ As a comment, you can see this as a larger issue for schemes with stabilization (shock capturing terms in CFD as an example). I'm not sure if your specific monotone/exponentially fitted scheme falls under this umbrella, but if it does, it's often advantageous to include these additional terms in the Jacobian as well due to their nonlinear effects. Jacobian-free Newton-Krylov methods (approximating a matrix-vector product by a difference scheme) often cite this as an advantage, that they automatically account for all such nonlinearities in a discrete scheme. $\endgroup$ – Jesse Chan Aug 12 '13 at 17:37

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.