[petsc-users] How to get a matrix and vector for a coupled system of equations?

[Dave May]

On 28 November 2016 at 10:43, Rolf Kuiper <kuiper at mpia.de> wrote:

> Hi Dave,
>
> Thanks a lot for your prompt reply! This is even easier than I thought
> (and that is most likely the reason, why I could not think into this
> direction) :)
> And yes, I might/should upgrade our PETSc version used (3.1), but it will
> take me some days to check the full code package.
>
> Now, I will create my DA (now in 3.6 called DMDA) for the coupled system
> via
> DACreate3d(PETSC_COMM_WORLD, Periodicity, DA_STENCIL_BOX, Nx, Ny, Nz, Px,
> Py, Pz, 4, 1, lx, ly, lz, &MyDA);
> instead of
> DACreate3d(PETSC_COMM_WORLD, Periodicity, DA_STENCIL_BOX, Nx, Ny, Nz, Px,
> Py, Pz, 1, 1, lx, ly, lz, &MyDA);
> which I have used for the single equation.
>
> I would like to ask one follow-up question:
> Currently, I loop over my matrix (and vectors) via three for-loops over
> the 3 spatial directions, and set the columns indices in the following way:
> for(k){
>         row.k    = k;
>         col[0].k = k;
>         col[1].k = k;
>         col[2].k = k;
>         col[3].k = k;
>         col[4].k = k+1;
>         col[5].k = k-1;
>         for(j){
>             row.j    = j;
>             col[0].j = j;
>             col[1].j = j;
>             col[2].j = j+1;
>             col[3].j = j-1;
>             col[4].j = j;
>             col[5].j = j;
>             for(i){
>                 row.i    = i;
>                 col[0].i = i+1;
>                 col[1].i = i-1;
>                 col[2].i = i;
>                 col[3].i = i;
>                 col[4].i = i;
>                 col[5].i = i;
> Now for the coupled equations, I should overall additionally loop over the
> number of DOFs (just 2 in my previous email example).
>

So presumably you are using MatSetValuesStencil() to insert entries into
the matrix.
http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatSetValuesStencil.html

You can keep using MatSetValuesStencil(), for you multi-component example,
however you will additionally need to enter a value for the member "c"
within the MatStencil struct.
See
http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatStencil.html

"c" relates to the component or DOF index. So the code above will need to
be modified slightly to define values for row.c and col[0].c , col[1].c ...
etc


Alternatively you can use MatSetValuesStencilBlocked()
http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatSetValuesBlockedStencil.html#MatSetValuesBlockedStencil




> Could you give me an easy example or pseudo-code for the associated
> assignment of columns (I have the same stencil in each of the submatrices,
> so for the 2 DOFs in 3D, I would get 7+6=13 column entries per row)? Or can
> you link me to an existing example within the PETSc help?
>

Note that the manual pages I've sent links to contains links to example
codes (see bottom of the webpage) where you can see how to use these
functions.

Thanks,
  Dave


>
> Again, Thanks a lot!
> Rolf
>
>
>
>
> Am 24.11.2016 um 22:30 schrieb Dave May <dave.mayhem23 at gmail.com>:
>
> When you create the DMDA, set the number of DOFs (degrees of freedom) per
> point to be 2 instead of 1.
>
> You must be using and ancient version of petsc given the function names
> you quoted. Consider upgrading to 3.7
>
> Thanks,
> Dave
>
> On Thu, 24 Nov 2016 at 20:24, Rolf Kuiper <kuiper at mpia.de> wrote:
>
>> Dear PETSc users,
>>
>> maybe this is an easy question, but I can’t find the information right
>> away in the user’s guide nor online.
>>
>> What I am currently doing and which works fine:
>> To solve a partial differential equation for the quantity q on a parallel
>> distributed grid, which is represented by the Distributed Array MyDA, I am
>> currently creating the associated sparse matrix for the KSP solver via
>> DAGetMatrix(MyDA, MATMPIAIJ, &MyMatrix);
>>
>> The solution vector and right hand side vector are created via
>> DACreateGlobalVector(MyDA, &MyRightHandSideVector);
>> VecDuplicate(MyRightHandSideVector, &GlobalSolutionVector);
>>
>> The DA is constructed using DACreate3d() with the corresponding regular
>> structured grid information.
>>
>> And here is my problem:
>> Now, I would like to solve a coupled system of equations for the
>> quantities q1 and q2 on the same grid. I.e., the matrix should just get the
>> double number of rows and columns, the vectors contain twice the number of
>> entries (e.g. first all q1s and then all q2s). And I would like to be sure
>> that the entries of q1 and q2, which are associated with the same grid cell
>> are located on the same processor.
>> Is there already a pre-defined structures available (such as MATMPIAIJ)
>> within PETSc to enlarge such a single equation to store the entries of
>> coupled equations? Such as
>> DACreateTwiceTheGlobalVector()?
>>
>> The equation is (simplified) of the form
>> d/dt q1 + grad q2 = 0
>> d/dt q2 + f(q1) = 0
>> with an arbitrary function f depending on q1.
>>
>> Thanks a lot for your help in advance,
>> Rolf
>>
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20161128/3fd0aa8e/attachment.html>