# [petsc-users] How to construct Jacobian for multicomponent problem ?

Matthew Knepley knepley at gmail.com
Wed Sep 25 06:46:31 CDT 2013

```On Wed, Sep 25, 2013 at 4:40 AM, Christophe Ortiz <
christophe.ortiz at ciemat.es> wrote:

> Hi Matt,
>
> On Wed, Sep 25, 2013 at 1:34 PM, Matthew Knepley <knepley at gmail.com>wrote:
>
>> On Wed, Sep 25, 2013 at 4:31 AM, Christophe Ortiz <
>> christophe.ortiz at ciemat.es> 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

> Christophe
>
>
>>
>>   Thanks,
>>
>>      Matt
>>
>>
>>> Is it a long vector like (u0, u1, ..., umx-1 , v0, v1, ..., vmx-1) ?
>>> Or is it (u0, v0, u1, v1, u2, v2, ....., umx-1, vmx-1) ?
>>>
>>> I don't understand the order of the components of the Jacobian and how
>>> to fill the matrix in that case.
>>>
>>> For instance, for the first row of the Jacobian, should I calculate
>>> dFu0/du0, dFu0/dFu1, ...., dFu0/dumx-1, dFu0/dv0, dFu0/dv1.....
>>>
>>> or
>>>
>>> dFu0/du0, dFu0/dv0, dFu0/du1, dFu0/dv1, ....
>>>
>>> And then, which subroutine should I use to fill in the Matrix ? Simply
>>> MatSetValues ?
>>>
>>> Many thanks in advance for your help.
>>>
>>> Christophe
>>>
>>
>>
>>
>> --
>> What most experimenters take for granted before they begin their
>> experiments is infinitely more interesting than any results to which their
>> -- Norbert Wiener
>>
>
>

--
What most experimenters take for granted before they begin their
experiments is infinitely more interesting than any results to which their