[petsc-users] MatSetValues() can't work right

Barry Smith bsmith at petsc.dev
Mon Mar 18 23:01:42 CDT 2024 


> On Mar 18, 2024, at 9:41 PM, Waltz Jan <jl2862237661 at gmail.com> wrote:
> 
> Thank you for your response. 
> However, even after reading the Notes in https://urldefense.us/v3/__https://petsc.org/release/manualpages/DM/DMCreateMatrix/__;!!G_uCfscf7eWS!fTO1ShsqXrxcXKmKrn7uXjX68PlSaKv4RBgRvwP9BUQpeowdAqyQyxq3cSp_3H231u74LG5cJRd24lnABMYgziE$ , I am still confused.
> According to the content in https://urldefense.us/v3/__https://petsc.org/release/manualpages/Mat/MatSetValues/__;!!G_uCfscf7eWS!fTO1ShsqXrxcXKmKrn7uXjX68PlSaKv4RBgRvwP9BUQpeowdAqyQyxq3cSp_3H231u74LG5cJRd24lnAAfForB0$ , idxm and idxn represent the global indices of the rows and columns, respectively. I created a matrix using DMDA with a size of 300x300, and I want to insert a value at the 101st row and 101st column. Shouldn't idxm and idxn be 100 in this case?

   It does go into the matrix at I,J of 101 and 101. But MatView() when the matrix comes from a DMDA prints the matrix out with a different ordering (the "natural" ordering for a 2 or 3d grid) thus printed out it is not at the 101, 101 location. 

  If you do PetscViewerPushFormat(viewer, PETSC_VIEWER_NATIVE); before calling MatView() then MatView will not print the matrix out in natural order, it will print it out as it is stored in parallel and you will see the location at 101, 101. 


>  But when NP=6, the inserted value appears at a different position. 
> 
  MatSetValuesStencil() allows each process to set entries for its points. You can use DMDAGetGhostCorners() to determine the allowed range of i,j,k that can be used in the MatStencil for each particular MPI process.  See for example src/snes/tests/ex20.c the function FormJacobian(). It is for finite difference matrices where points only speak to their neighbors in the matrix. It does not support setting arbitrary I, J locations into a matrix from arbitrary processes.


> I tried using MatSetValuesStencil, but the code threw an error with the following message:
> CODES:
> #include <petscdmda.h>
> #include <petscsys.h>
> #include <petscdm.h>
> #include <petscdmda.h>
> #include <petscsys.h>
> #include <petscdm.h>
> #include <petsc/private/snesimpl.h>
> 
> int main()
> {
>     PetscInitialize(NULL, NULL, NULL, NULL);
>     DM da;
>     DMDACreate3d(PETSC_COMM_WORLD, DM_BOUNDARY_GHOSTED, DM_BOUNDARY_GHOSTED, DM_BOUNDARY_GHOSTED, DMDA_STENCIL_STAR,
>                  10, 1, 10, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, 3, 1, NULL, NULL, NULL, &da);
>     DMSetFromOptions(da);
>     DMSetUp(da);
>     Mat Jac;
>     DMCreateMatrix(da, &Jac);
> 
>     double val = 1.;
> 
>     MatStencil st{5,1,5,0};
>     MatSetValuesStencil(Jac, 1, &st, 1, &st, &val, INSERT_VALUES);
>     MatAssemblyBegin(Jac, MAT_FINAL_ASSEMBLY);
>     MatAssemblyEnd(Jac, MAT_FINAL_ASSEMBLY);
> 
>     PetscViewer viewer;
>     PetscViewerASCIIOpen(PETSC_COMM_WORLD, "./jacobianmatrix.m", &viewer);
>     PetscViewerPushFormat(viewer, PETSC_VIEWER_ASCII_MATLAB);
>     MatView(Jac, viewer);
>     PetscViewerDestroy(&viewer);
> 
>     PetscFinalize();
> }
> 
> [0]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------
> [0]PETSC ERROR: Argument out of range
> [0]PETSC ERROR: Local index 438 too large 377 (max) at 0
> [0]PETSC ERROR: See https://urldefense.us/v3/__https://petsc.org/release/faq/__;!!G_uCfscf7eWS!fTO1ShsqXrxcXKmKrn7uXjX68PlSaKv4RBgRvwP9BUQpeowdAqyQyxq3cSp_3H231u74LG5cJRd24lnAY72xwtI$  for trouble shooting.
> [0]PETSC ERROR: Petsc Release Version 3.20.4, Jan 29, 2024 
> [0]PETSC ERROR: Unknown Name on a arch-linux-cxx-debug named TUF-Gaming by lei Tue Mar 19 09:31:01 2024
> [0]PETSC ERROR: Configure options --with-cc=gcc --with-cxx=g++ --with-fc=gfortran --download-fblaslapack --download-hypre --with-debugging=yes --download-mpich --with-clanguage=cxx
> [0]PETSC ERROR: #1 ISLocalToGlobalMappingApply() at /home/lei/Software/PETSc/petsc-3.20.4/src/vec/is/utils/isltog.c:789
> [0]PETSC ERROR: #2 MatSetValuesLocal() at /home/lei/Software/PETSc/petsc-3.20.4/src/mat/interface/matrix.c:2408
> [0]PETSC ERROR: #3 MatSetValuesStencil() at /home/lei/Software/PETSc/petsc-3.20.4/src/mat/interface/matrix.c:1762
> [1]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------
> [1]PETSC ERROR: Argument out of range
> [1]PETSC ERROR: Local index 423 too large 377 (max) at 0
> [1]PETSC ERROR: See https://urldefense.us/v3/__https://petsc.org/release/faq/__;!!G_uCfscf7eWS!fTO1ShsqXrxcXKmKrn7uXjX68PlSaKv4RBgRvwP9BUQpeowdAqyQyxq3cSp_3H231u74LG5cJRd24lnAY72xwtI$  for trouble shooting.
> [1]PETSC ERROR: Petsc Release Version 3.20.4, Jan 29, 2024 
> [1]PETSC ERROR: Unknown Name on a arch-linux-cxx-debug named TUF-Gaming by lei Tue Mar 19 09:31:01 2024
> [1]PETSC ERROR: Configure options --with-cc=gcc --with-cxx=g++ --with-fc=gfortran --download-fblaslapack --download-hypre --with-debugging=yes --download-mpich --with-clanguage=cxx
> [1]PETSC ERROR: #1 ISLocalToGlobalMappingApply() at /home/lei/Software/PETSc/petsc-3.20.4/src/vec/is/utils/isltog.c:789
> [1]PETSC ERROR: #2 MatSetValuesLocal() at /home/lei/Software/PETSc/petsc-3.20.4/src/mat/interface/matrix.c:2408
> [1]PETSC ERROR: #3 MatSetValuesStencil() at /home/lei/Software/PETSc/petsc-3.20.4/src/mat/interface/matrix.c:1762
> 
> Is it not possible to set values across processors using MatSetValuesStencil? If I want to set values of the matrix across processors, what should I do? 
> I am really confused, and I would greatly appreciate your help.
> 
> On Mon, Mar 18, 2024 at 9:28 PM Barry Smith <bsmith at petsc.dev <mailto:bsmith at petsc.dev>> wrote:
>> 
>>    The output is correct (only confusing). For PETSc DMDA by default viewing a parallel matrix converts it to the "natural" ordering instead of the PETSc parallel ordering.
>> 
>>    See the Notes in https://urldefense.us/v3/__https://petsc.org/release/manualpages/DM/DMCreateMatrix/__;!!G_uCfscf7eWS!fTO1ShsqXrxcXKmKrn7uXjX68PlSaKv4RBgRvwP9BUQpeowdAqyQyxq3cSp_3H231u74LG5cJRd24lnABMYgziE$ 
>> 
>>   Barry
>> 
>> 
>>> On Mar 18, 2024, at 8:06 AM, Waltz Jan <jl2862237661 at gmail.com <mailto:jl2862237661 at gmail.com>> wrote:
>>> 
>>> This Message Is From an External Sender
>>> This message came from outside your organization.
>>> PETSc version: 3.20.4
>>> Program:
>>> #include <petscdmda.h>
>>> #include <petscsys.h>
>>> #include <petscdm.h>
>>> #include <petscdmda.h>
>>> #include <petscsys.h>
>>> #include <petscdm.h>
>>> #include <petsc/private/snesimpl.h>
>>> 
>>> int main()
>>> {
>>>     PetscInitialize(NULL, NULL, NULL, NULL);
>>>     DM da;
>>>     DMDACreate3d(PETSC_COMM_WORLD, DM_BOUNDARY_GHOSTED, DM_BOUNDARY_GHOSTED, DM_BOUNDARY_GHOSTED, DMDA_STENCIL_STAR,
>>>                  10, 1, 10, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, 3, 1, NULL, NULL, NULL, &da);
>>>     DMSetFromOptions(da);
>>>     DMSetUp(da);
>>>     Mat Jac;
>>>     DMCreateMatrix(da, &Jac);
>>>     int row = 100, col = 100;
>>>     double val = 1.;
>>>     MatSetValues(Jac, 1, &row, 1, &col, &val, INSERT_VALUES);
>>>     MatAssemblyBegin(Jac, MAT_FINAL_ASSEMBLY);
>>>     MatAssemblyEnd(Jac, MAT_FINAL_ASSEMBLY);
>>> 
>>>     PetscViewer viewer;
>>>     PetscViewerASCIIOpen(PETSC_COMM_WORLD, "./jacobianmatrix.m", &viewer);
>>>     PetscViewerPushFormat(viewer, PETSC_VIEWER_ASCII_MATLAB);
>>>     MatView(Jac, viewer);
>>>     PetscViewerDestroy(&viewer);
>>> 
>>>     PetscFinalize();
>>> }
>>> 
>>> When I ran the program with np = 6, I got the result as the below
>>> <image.png>
>>> It's obviously wrong.
>>> When I ran the program with np = 1 or 8, I got the right result as
>>> <image.png>
>> 

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20240319/c8307f7a/attachment-0001.html>

More information about the petsc-users mailing list