On Wed, Sep 25, 2013 at 1:46 PM, Matthew Knepley <[email protected]> wrote:
> On Wed, Sep 25, 2013 at 4:40 AM, Christophe Ortiz < > [email protected]> wrote: > >> Hi Matt, >> Thanks for your prompt reply. >> >> On Wed, Sep 25, 2013 at 1:34 PM, Matthew Knepley <[email protected]>wrote: >> >>> On Wed, Sep 25, 2013 at 4:31 AM, Christophe Ortiz < >>> [email protected]> wrote: >>> >>>> Hi guys, >>>> >>>> Please apologize if this question has already been raised in previous >>>> posts. >>>> >>>> I am developing a code in fortran 90 to solve a 1D multicomponent >>>> problem. >>>> For instance, something like >>>> >>>> u_t = u_xx + R(u,v) >>>> v_t = v_xx + R(u,v) >>>> >>>> I have already implemented this for one component (u) and using DMDA. >>>> So far no problem. It works fine. >>>> >>>> Now I am trying to implement several components per node. Following >>>> example ex22f.F, for the residual function I would have for node i (1D) >>>> >>>> F(1,i) = .... >>>> F(2,i) = .... >>>> >>>> Now, I am not quite sure to understand how to construct the Jacobian >>>> when there are various components. In this case, how is the vector in PETSc >>>> ? >>>> >>> >>> The Jacobian is a matrix, and the rows are ordered in the same way as in >>> the residual function. The columns are the same as the rows. >>> >> >> Yes, I understand this for one component. My question is when you have >> dof > 1. For a multicomponent problem you can write the residual function >> as an array: >> >> F(1,i) =... (component 1) >> F(2,i) = ... (component 2) >> >> Since there is only one Jacobian for the whole system, how are the >> components of the Jacobian ordered in that case ? Is it first the >> components of u, then the components of v ? >> > > The ordering of the rows (and columns) of the Jacobian match the rows of > the residual EVEN in the multicomponent case. > > Matt > > Ok. But what I still don't understand is the ordering of the residual in the multicomponent case. Let's say dof=2, components u and v. Which is the correct ordering in the following ? Fu0, Fu1, Fu2,...Fumx-1, Fv0, Fv1, Fv2, ...Fvmx-1 or Fu0, Fv0, Fu1, Fv1, Fu2, Fv2, ... ? It is this point that is not clear to me. Also, which subroutine do you advice to use to fill in the Matrix ? MatSetValues ? MatSetValuesBlocked ? Christophe
