On the usage of MatSetValuesBlocked
Hello Jed/Barry/Petsc friends I am trying to assemble a block matrix with 3x3 in 2D and 4x4 blocks in 3D coming from the fully coupled NS equation. I am not sure I am understanding the example provided here: https://petsc.org/main/docs/manualpages/Mat/MatSetValuesBlocked/ The description says that " *v -* a logically two-dimensional array of values", while in the example is passed as a 1D array. Should I pass a 2D array like v(1:nComp,1:nComp) or an array v(1:nComp*nComp) to MatSetValuesBlocked (I am using fortran that's why I start form 1 in my arrays) ? Are *idxm *and *idxn *the global index of each block right? I guess PETSc will correctly assign each block row and column as it knows the matrix has a given structure. Thank you as always!
On Thu, Nov 3, 2022 at 12:16 PM Edoardo alinovi <[email protected]> wrote:
Hello Jed/Barry/Petsc friends
I am trying to assemble a block matrix with 3x3 in 2D and 4x4 blocks in 3D coming from the fully coupled NS equation.
I am not sure I am understanding the example provided here: https://petsc.org/main/docs/manualpages/Mat/MatSetValuesBlocked/
The description says that " *v -* a logically two-dimensional array of values", while in the example is passed as a 1D array.
Should I pass a 2D array like v(1:nComp,1:nComp) or an array v(1:nComp*nComp) to MatSetValuesBlocked (I am using fortran that's why I start form 1 in my arrays) ?
No, you pass a contiguous chunk of memory, but it is _logically_ 2D in that the size is (idxm * bs) x (idxn * bs) Are *idxm *and *idxn *the global index of each block right? I guess PETSc
will correctly assign each block row and column as it knows the matrix has a given structure.
Yes. Thanks, Matt
Thank you as always!
-- What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead. -- Norbert Wiener https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>
Hello Matt, I see, so from an operation point of view, I should pass a 1D array... How should I unroll my nComp x nComp matrix? by row or by column? e.g. say that my block is: A=[ 1 2 3 4 5 6 7 8 9 ] than v = [1 2 3 4 5 6 7 8 9] or [1 4 7 2 5 8 3 6 9]?
You can find the current F90 interface definitions we support for MatSetValuesBlocked() in /src/mat/f90-mod/petscmat.h90
On Nov 3, 2022, at 12:16 PM, Edoardo alinovi <[email protected]> wrote:
Hello Jed/Barry/Petsc friends
I am trying to assemble a block matrix with 3x3 in 2D and 4x4 blocks in 3D coming from the fully coupled NS equation.
I am not sure I am understanding the example provided here: https://petsc.org/main/docs/manualpages/Mat/MatSetValuesBlocked/
The description says that " v - a logically two-dimensional array of values", while in the example is passed as a 1D array.
Should I pass a 2D array like v(1:nComp,1:nComp) or an array v(1:nComp*nComp) to MatSetValuesBlocked (I am using fortran that's why I start form 1 in my arrays) ?
We intend to support both approaches, whatever is most convenient for the rest of your code. (Perhaps more interface definitions are needed?) When we say "logically" two-dimensions this is intended to mean that you can pass a one dimensional array that contains all the nonzeros in the block in the order of the first column, followed by the second column etc. (How Fortran handles two dimensional arrays normally). But in most circumstances I would guess providing directly the two dimensal Fortran array is more natural.
Are idxm and idxn the global index of each block right? I guess PETSc will correctly assign each block row and column as it knows the matrix has a given structure.
Yes
Thank you as always!
Hi Barry, Thanks for popping in. This is my code: * blockValues = this%getDiagonalBlockValues(iElement=iElement) call MatSetValuesBlocked(this%A, 4-bdim, mesh%cellGlobalAddr(iElement)-1, 4-bdim, mesh%cellGlobalAddr(iElement)-1, blockValues, INSERT_VALUES, ierr)* *blockValues *is a 3x3 or 4x4 matrix and I am passing it straight away.
Also, just to be 100% sure, is m and equal to 3 in 2 and 4 in 3D if my blocks are 3x3 and 4x4 respectively?
Ah, I was forgetting the most important thing... Are the size of idxm and idxn equal to one if I insert 1 block or should I specify all the rows and columns in the block? I am getting some memory issues with unallocated non zero values so I must have made some mistake here... :( Sorry for the ton of questions!
On Nov 3, 2022, at 1:16 PM, Edoardo alinovi <[email protected]> wrote:
Ah, I was forgetting the most important thing... Are the size of idxm and idxn equal to one if I insert 1 block or should I specify all the rows and columns in the block?
Yes, for a single block they are one. The block size is set with MatSetBlockSize() or in the preallocation routines or automatically if you use DMCreateMatrix(), the block size is not passed in during MatSetValuesBlocked() the matrix uses its internal value. Barry
I am getting some memory issues with unallocated non zero values so I must have made some mistake here... :(
Sorry for the ton of questions!
Yes I am doing: call MatMPIBAIJSetPreallocation(this%A, 4-bdim, flubioSolvers%d_nz, mesh%d_nnz, flubioSolvers%o_nz, mesh%o_nnz, ierr) with d_nnz the number of diagonal blocks and o_nnz the number of off-diagonal blocks. However I am getting this: [0]PETSC ERROR: --------------------- Error Message -------------------------------------------------------------- [0]PETSC ERROR: Argument out of range [0]PETSC ERROR: New nonzero at (3,3) caused a malloc Use MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE) to turn off this check [0]PETSC ERROR: See https://petsc.org/release/faq/ for trouble shooting. [0]PETSC ERROR: Petsc Release Version 3.18.0, Sep 30, 2022 [0]PETSC ERROR: flubio_coupled on a gnu named alienware by edo Thu Nov 3 18:19:30 2022 [0]PETSC ERROR: Configure options PETSC_ARCH=gnu FOPTFLAGS=-O3 COPTFLAGS=-O3 CXXOPTFLAGS=-O3 -with-debugging=no -download-fblaslapack=1 -download-superlu_dist -download-mumps -download-hypre -download-metis -download-parmetis -download-scalapack -download-ml -download-slepc -download-hpddm -download-cmake -with-mpi-dir=/home/edo/software/openmpi-4.1.1/build/ [0]PETSC ERROR: #1 MatSetValuesBlocked_SeqBAIJ_Inlined() at /home/edo/software/petsc-3.18.0/src/mat/impls/baij/mpi/mpibaij.c:318 [0]PETSC ERROR: #2 MatSetValuesBlocked_MPIBAIJ() at /home/edo/software/petsc-3.18.0/src/mat/impls/baij/mpi/mpibaij.c:389 [1]PETSC ERROR: --------------------- Error Message -------------------------------------------------------------- [1]PETSC ERROR: Argument out of range [1]PETSC ERROR: New nonzero at (0,0) caused a malloc Use MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE) to turn off this check [1]PETSC ERROR: See https://petsc.org/release/faq/ for trouble shooting. [1]PETSC ERROR: Petsc Release Version 3.18.0, Sep 30, 2022 [1]PETSC ERROR: flubio_coupled on a gnu named alienware by edo Thu Nov 3 18:19:30 2022 [1]PETSC ERROR: Configure options PETSC_ARCH=gnu FOPTFLAGS=-O3 COPTFLAGS=-O3 CXXOPTFLAGS=-O3 -with-debugging=no -download-fblaslapack=1 -download-superlu_dist -download-mumps -download-hypre -download-metis -download-parmetis -download-scalapack -download-ml -download-slepc -download-hpddm -download-cmake -with-mpi-dir=/home/edo/software/openmpi-4.1.1/build/ [1]PETSC ERROR: #1 MatSetValuesBlocked_SeqBAIJ_Inlined() at /home/edo/software/petsc-3.18.0/src/mat/impls/baij/mpi/mpibaij.c:318 [1]PETSC ERROR: #2 MatSetValuesBlocked_MPIBAIJ() at /home/edo/software/petsc-3.18.0/src/mat/impls/baij/mpi/mpibaij.c:419 [1]PETSC ERROR: #3 MatSetValuesBlocked() at /home/edo/software/petsc-3.18.0/src/mat/interface/matrix.c:1978 [1]PETSC ERROR: #4 MatSetValuesBlocked_SeqBAIJ_Inlined() at /home/edo/software/petsc-3.18.0/src/mat/impls/baij/mpi/mpibaij.c:318 [1]PETSC ERROR: #5 MatSetValuesBlocked_MPIBAIJ() at /home/edo/software/petsc-3.18.0/src/mat/impls/baij/mpi/mpibaij.c:419 [1]PETSC ERROR: #6 MatSetValuesBlocked() at /home/edo/software/petsc-3.18.0/src/mat/interface/matrix.c:1978 [0]PETSC ERROR: #3 MatSetValuesBlocked() at /home/edo/software/petsc-3.18.0/src/mat/interface/matrix.c:1978 [0]PETSC ERROR: #4 MatSetValuesBlocked_SeqBAIJ_Inlined() at /home/edo/software/petsc-3.18.0/src/mat/impls/baij/mpi/mpibaij.c:318 [0]PETSC ERROR: #5 MatSetValuesBlocked_MPIBAIJ() at /home/edo/software/petsc-3.18.0/src/mat/impls/baij/mpi/mpibaij.c:419 [0]PETSC ERROR: #6 MatAssemblyEnd_MPIBAIJ() at /home/edo/software/petsc-3.18.0/src/mat/impls/baij/mpi/mpibaij.c:906 [0]PETSC ERROR: #7 MatAssemblyEnd() at /home/edo/software/petsc-3.18.0/src/mat/interface/matrix.c:5696
Just to give a bit more of context I am doing like this: call MatCreate(PETSC_COMM_WORLD, this%A, ierr) call MatSetSizes(this%A, lm, lm, M, M, ierr) ! lm is the local size of the matrix , M the global one call MatSetType(this%A, myType, ierr) ! myType is MATMPIBAIJ call MatSetBlockSize(this%A, 4-bdim, ierr) ! 4-bdim is equal to 3 in this case call MatMPIBAIJSetPreallocation(this%A, 4-bdim, flubioSolvers%d_nz, mesh%d_nnz, flubioSolvers%o_nz, mesh%o_nnz, ierr) ! d_nnz and o_nnz is the number of diagonal and off diagonal non zero blocks call MatSetUp(this%A, ierr) ... some non relevant code ..... call MatSetValuesBlocked(this%A, 4-bdim, mesh%cellGlobalAddr(iElement)-1, 4-bdim, mesh%cellGlobalAddr(iElement)-1, blockValues, INSERT_VALUES, ierr) mesh%cellGlobalAddr(iElement)-1 is an integer equal (not an array) to the block element number the block matrix *blockValues (3x3) *belongs to. Any evident errors?
Well, definitely using MatsetValuesBlocked in a bad way. Instead if seeing 11 in the first place 3 rows and 3 columns, I see all zeros and random numbers in row 9 10 and 11... [image: image.png]
The error indicates not enough nonzero blocks are preallocated for. Try something really simple, preallocate for one block and put in one block then call MatAssemblyBegin/End(), MatView(), if that works then work up to your full problem. Barry
On Nov 3, 2022, at 3:11 PM, Edoardo alinovi <[email protected]> wrote:
Just to give a bit more of context I am doing like this:
call MatCreate(PETSC_COMM_WORLD, this%A, ierr)
call MatSetSizes(this%A, lm, lm, M, M, ierr) ! lm is the local size of the matrix , M the global one
call MatSetType(this%A, myType, ierr) ! myType is MATMPIBAIJ
call MatSetBlockSize(this%A, 4-bdim, ierr) ! 4-bdim is equal to 3 in this case
call MatMPIBAIJSetPreallocation(this%A, 4-bdim, flubioSolvers%d_nz, mesh%d_nnz, flubioSolvers%o_nz, mesh%o_nnz, ierr) ! d_nnz and o_nnz is the number of diagonal and off diagonal non zero blocks
call MatSetUp(this%A, ierr)
... some non relevant code .....
call MatSetValuesBlocked(this%A, 4-bdim, mesh%cellGlobalAddr(iElement)-1, 4-bdim, mesh%cellGlobalAddr(iElement)-1, blockValues, INSERT_VALUES, ierr) mesh%cellGlobalAddr(iElement)-1 is an integer equal (not an array) to the block element number the block matrix blockValues (3x3) belongs to.
Any evident errors?
Barry, Can you please provide me with an example on how to use MatSetValuesBlocked? To play it easy, let's say that I want to insert a 3x3 block matrix b into the matrix A, rows 0-2, columns 0-2. Up to what I've understood (very few apparently XD ), I would do like this: b(3,3) = 11.0 call MatSetValuesBlocked(A, 3, 0, 3, 0, b, INSERT_VALUES, ierr). This does not work at all, I get this result that does not make any sense 😫 [image: image.png] It places 6 values instead of 9 and it they are in odd locations (0 1 2 9 10 11). Also I noted that I am getting different results if in place of the zero in red I use a fortran integer 😥 Super thanks for the help
You should pass 1 and 1 not 3, because you are setting one block. Regarding all the integer values passed in to PETSc routines To be completely portable you need to declare them as PetscInt and pass the variables. But if you do not use --with-64-bit-indices in ./configure and you do not use some Fortran compiler option to promote all integers to 64 bit integers just passing in directly the 0 etc is fine. The reason you get different results with a 0 or a Fortran integer is because of the 3 PETSc tries to three indices from 0 value you pass in so it is readying memory it is not suppose to be reading. Once you change the 3 to 1 it will likely be fine. Barry
On Nov 3, 2022, at 5:56 PM, Edoardo alinovi <[email protected]> wrote:
Barry,
Can you please provide me with an example on how to use MatSetValuesBlocked?
To play it easy, let's say that I want to insert a 3x3 block matrix b into the matrix A, rows 0-2, columns 0-2. Up to what I've understood (very few apparently XD ), I would do like this:
b(3,3) = 11.0 call MatSetValuesBlocked(A, 3, 0, 3, 0, b, INSERT_VALUES, ierr).
This does not work at all, I get this result that does not make any sense 😫
<image.png>
It places 6 values instead of 9 and it they are in odd locations (0 1 2 9 10 11).
Also I noted that I am getting different results if in place of the zero in red I use a fortran integer 😥
Super thanks for the help
Barry, Thanks, yes, I should pass 1 and not 3.... For some reason I have misunderstood the wording in the documentation and I have interpreted m and n as the number of row and columns of the block to insert, while I need to think of everything as divided by bs... Now I am with you!!! Many thanks and sorry to have been so dump! 😓
It is working like a charm now! Is it mandatory to use VecSetValuesBlocked to assemble the rhs? Does the Vec need to be of any other type than VECMPI? I am assembling it like this: brhs(1:3-bdim) = this%Ueqn%bC(iElement,1:3-bdim) brhs(4-bdim) = this%Peqn%bC(iElement,1) call VecSetValuesBlocked(this%rhs, 1, mesh%cellGlobalAddr(iElement)-1, brhs, INSERT_VALUES, ierr) But I am getting into troubles: [0]PETSC ERROR: --------------------- Error Message -------------------------------------------------------------- [0]PETSC ERROR: *PetscSegBufferAlloc_Private* [0]PETSC ERROR: See https://petsc.org/release/faq/ for trouble shooting. [0]PETSC ERROR: Petsc Release Version 3.18.0, Sep 30, 2022 [0]PETSC ERROR: flubio_coupled on a gnu named alienware by edo Fri Nov 4 09:31:03 2022 [0]PETSC ERROR: Configure options PETSC_ARCH=gnu FOPTFLAGS=-O3 COPTFLAGS=-O3 CXXOPTFLAGS=-O3 -with-debugging=no -download-fblaslapack=1 -download-superlu_dist -download-mumps -download-hypre -download-metis -download-parmetis -download-scalapack -download-ml -download-slepc -download-hpddm -download-cmake -with-mpi-dir=/home/edo/software/openmpi-4.1.1/build/ [0]PETSC ERROR: #1 PetscMallocAlign() at /home/edo/software/petsc-3.18.0/src/sys/memory/mal.c:55 [0]PETSC ERROR: #2 PetscSegBufferAlloc_Private() at /home/edo/software/petsc-3.18.0/src/sys/utils/segbuffer.c:31 [0]PETSC ERROR: #3 PetscSegBufferGet() at /home/edo/software/petsc-3.18.0/src/sys/utils/segbuffer.c:94 [1]PETSC ERROR: --------------------- Error Message -------------------------------------------------------------- [1]PETSC ERROR: General MPI error [1]PETSC ERROR: MPI error 1 MPI_ERR_BUFFER: invalid buffer pointer [1]PETSC ERROR: See https://petsc.org/release/faq/ for trouble shooting. [1]PETSC ERROR: Petsc Release Version 3.18.0, Sep 30, 2022 [1]PETSC ERROR: flubio_coupled on a gnu named alienware by edo Fri Nov 4 09:31:03 2022 [1]PETSC ERROR: Configure options PETSC_ARCH=gnu FOPTFLAGS=-O3 COPTFLAGS=-O3 CXXOPTFLAGS=-O3 -with-debugging=no -download-fblaslapack=1 -download-superlu_dist -download-mumps -download-hypre -download-metis -download-parmetis -download-scalapack -download-ml -download-slepc -download-hpddm -download-cmake -with-mpi-dir=/home/edo/software/openmpi-4.1.1/build/ [1]PETSC ERROR: #1 VecAssemblySend_MPI_Private() at /home/edo/software/petsc-3.18.0/src/vec/vec/impls/mpi/pbvec.c:133 [1]PETSC ERROR: #2 PetscCommBuildTwoSidedFReq_Reference() at /home/edo/software/petsc-3.18.0/src/sys/utils/mpits.c:314 [1]PETSC ERROR: #3 PetscCommBuildTwoSidedFReq() at /home/edo/software/petsc-3.18.0/src/sys/utils/mpits.c:526 [1]PETSC ERROR: [0]PETSC ERROR: #4 VecAssemblyRecv_MPI_Private() at /home/edo/software/petsc-3.18.0/src/vec/vec/impls/mpi/pbvec.c:164 [0]PETSC ERROR: #5 PetscCommBuildTwoSidedFReq_Reference() at /home/edo/software/petsc-3.18.0/src/sys/utils/mpits.c:320 [0]PETSC ERROR: #6 PetscCommBuildTwoSidedFReq() at /home/edo/software/petsc-3.18.0/src/sys/utils/mpits.c:526 [0]PETSC ERROR: #7 VecAssemblyBegin_MPI_BTS() at /home/edo/software/petsc-3.18.0/src/vec/vec/impls/mpi/pbvec.c:238 #4 VecAssemblyBegin_MPI_BTS() at /home/edo/software/petsc-3.18.0/src/vec/vec/impls/mpi/pbvec.c:238 [1]PETSC ERROR: #5 VecAssemblyBegin() at /home/edo/software/petsc-3.18.0/src/vec/vec/interface/vector.c:124 [1]PETSC ERROR: #6 VecAssemblyEnd_MPI_BTS() at /home/edo/software/petsc-3.18.0/src/vec/vec/impls/mpi/pbvec.c:337 [1]PETSC ERROR: #7 VecAssemblyEnd() at /home/edo/software/petsc-3.18.0/src/vec/vec/interface/vector.c:158 [1]PETSC ERROR: #8 VecView() at /home/edo/software/petsc-3.18.0/src/vec/vec/interface/vector.c:719 [0]PETSC ERROR: #8 VecAssemblyBegin() at /home/edo/software/petsc-3.18.0/src/vec/vec/interface/vector.c:124 Vec Object: 2 MPI processes
On Fri, Nov 4, 2022 at 4:32 AM Edoardo alinovi <[email protected]> wrote:
It is working like a charm now!
Is it mandatory to use VecSetValuesBlocked to assemble the rhs? Does the Vec need to be of any other type than VECMPI?
SetValuesBlocked() is never required. You can always use the normal versions, but you would have to supply all the indices, not just the block indices. Thanks, Matt
I am assembling it like this: brhs(1:3-bdim) = this%Ueqn%bC(iElement,1:3-bdim) brhs(4-bdim) = this%Peqn%bC(iElement,1) call VecSetValuesBlocked(this%rhs, 1, mesh%cellGlobalAddr(iElement)-1, brhs, INSERT_VALUES, ierr)
But I am getting into troubles:
[0]PETSC ERROR: --------------------- Error Message -------------------------------------------------------------- [0]PETSC ERROR: *PetscSegBufferAlloc_Private* [0]PETSC ERROR: See https://petsc.org/release/faq/ for trouble shooting. [0]PETSC ERROR: Petsc Release Version 3.18.0, Sep 30, 2022 [0]PETSC ERROR: flubio_coupled on a gnu named alienware by edo Fri Nov 4 09:31:03 2022 [0]PETSC ERROR: Configure options PETSC_ARCH=gnu FOPTFLAGS=-O3 COPTFLAGS=-O3 CXXOPTFLAGS=-O3 -with-debugging=no -download-fblaslapack=1 -download-superlu_dist -download-mumps -download-hypre -download-metis -download-parmetis -download-scalapack -download-ml -download-slepc -download-hpddm -download-cmake -with-mpi-dir=/home/edo/software/openmpi-4.1.1/build/ [0]PETSC ERROR: #1 PetscMallocAlign() at /home/edo/software/petsc-3.18.0/src/sys/memory/mal.c:55 [0]PETSC ERROR: #2 PetscSegBufferAlloc_Private() at /home/edo/software/petsc-3.18.0/src/sys/utils/segbuffer.c:31 [0]PETSC ERROR: #3 PetscSegBufferGet() at /home/edo/software/petsc-3.18.0/src/sys/utils/segbuffer.c:94 [1]PETSC ERROR: --------------------- Error Message -------------------------------------------------------------- [1]PETSC ERROR: General MPI error [1]PETSC ERROR: MPI error 1 MPI_ERR_BUFFER: invalid buffer pointer [1]PETSC ERROR: See https://petsc.org/release/faq/ for trouble shooting. [1]PETSC ERROR: Petsc Release Version 3.18.0, Sep 30, 2022 [1]PETSC ERROR: flubio_coupled on a gnu named alienware by edo Fri Nov 4 09:31:03 2022 [1]PETSC ERROR: Configure options PETSC_ARCH=gnu FOPTFLAGS=-O3 COPTFLAGS=-O3 CXXOPTFLAGS=-O3 -with-debugging=no -download-fblaslapack=1 -download-superlu_dist -download-mumps -download-hypre -download-metis -download-parmetis -download-scalapack -download-ml -download-slepc -download-hpddm -download-cmake -with-mpi-dir=/home/edo/software/openmpi-4.1.1/build/ [1]PETSC ERROR: #1 VecAssemblySend_MPI_Private() at /home/edo/software/petsc-3.18.0/src/vec/vec/impls/mpi/pbvec.c:133 [1]PETSC ERROR: #2 PetscCommBuildTwoSidedFReq_Reference() at /home/edo/software/petsc-3.18.0/src/sys/utils/mpits.c:314 [1]PETSC ERROR: #3 PetscCommBuildTwoSidedFReq() at /home/edo/software/petsc-3.18.0/src/sys/utils/mpits.c:526 [1]PETSC ERROR: [0]PETSC ERROR: #4 VecAssemblyRecv_MPI_Private() at /home/edo/software/petsc-3.18.0/src/vec/vec/impls/mpi/pbvec.c:164 [0]PETSC ERROR: #5 PetscCommBuildTwoSidedFReq_Reference() at /home/edo/software/petsc-3.18.0/src/sys/utils/mpits.c:320 [0]PETSC ERROR: #6 PetscCommBuildTwoSidedFReq() at /home/edo/software/petsc-3.18.0/src/sys/utils/mpits.c:526 [0]PETSC ERROR: #7 VecAssemblyBegin_MPI_BTS() at /home/edo/software/petsc-3.18.0/src/vec/vec/impls/mpi/pbvec.c:238 #4 VecAssemblyBegin_MPI_BTS() at /home/edo/software/petsc-3.18.0/src/vec/vec/impls/mpi/pbvec.c:238 [1]PETSC ERROR: #5 VecAssemblyBegin() at /home/edo/software/petsc-3.18.0/src/vec/vec/interface/vector.c:124 [1]PETSC ERROR: #6 VecAssemblyEnd_MPI_BTS() at /home/edo/software/petsc-3.18.0/src/vec/vec/impls/mpi/pbvec.c:337 [1]PETSC ERROR: #7 VecAssemblyEnd() at /home/edo/software/petsc-3.18.0/src/vec/vec/interface/vector.c:158 [1]PETSC ERROR: #8 VecView() at /home/edo/software/petsc-3.18.0/src/vec/vec/interface/vector.c:719 [0]PETSC ERROR: #8 VecAssemblyBegin() at /home/edo/software/petsc-3.18.0/src/vec/vec/interface/vector.c:124 Vec Object: 2 MPI processes
-- What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead. -- Norbert Wiener https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>
Thanks Matt, I have found out that setValuesblocked will work if I do: call MatCreateVecs(A, x, y, ierr) call setValuesBlocked(x, nblocks, varray, ierr) However, there is nogetValuesBlocked. Not the end of the world, it is handy to set and get stuff by block and not by single entry :) Cheers
On Nov 4, 2022, at 6:55 AM, Edoardo alinovi <[email protected]> wrote:
Thanks Matt,
I have found out that setValuesblocked will work if I do:
call MatCreateVecs(A, x, y, ierr) call setValuesBlocked(x, nblocks, varray, ierr)
Ah, likely the block size for the vector was not correct, leading to the memory corruption. MatCreateVecs() creates a vector compatible with the matrix, same block size and parallel layout so you don't need to worry about setting those values yourself. Barry
However, there is nogetValuesBlocked. Not the end of the world, it is handy to set and get stuff by block and not by single entry :)
Cheers
Yes, I did not set the block size for the vector... Missed it! I think I have nailed the way to handle block matrix/vectors, I am moving now on to solve the next facy error which is a ksp_diverged_its 😫😓 thanks gents for the support with this block madness :)
Matt, Barry, Should I do any particular trick to solve block matrices in ksp? I am doing a silly 3x3 cavity test case and I am struggling to converge using CG+bjacobi. It might be I have an error in the matrix, but just to be sure I am not missing something fundamental in the setup. This is my log that shows how norm residual is far to be happy: *10000 *KSP unpreconditioned resid norm 2.273088479279e+03 true resid norm 2.273088479279e+03 ||r(i)||/||b|| 6.561841227018e+02 Reason = -3 ERROR: KSP has not converged. Simulations stopped. [0]PETSC ERROR: --------------------- Error Message -------------------------------------------------------------- [0]PETSC ERROR: KSPSolve has not converged, reason DIVERGED_ITS Thank you
On 5 Nov 2022, at 10:31 AM, Edoardo alinovi <[email protected]> wrote:
Matt, Barry,
Should I do any particular trick to solve block matrices in ksp?
I am doing a silly 3x3 cavity test case and I am struggling to converge using CG+bjacobi.
This is far from the ideal preconditioner. First, you should check that your assembly is correct and stick to PCLU. Does the KSP report an error in this case? Does the solution « look » acceptable? Thanks, Pierre
It might be I have an error in the matrix, but just to be sure I am not missing something fundamental in the setup.
This is my log that shows how norm residual is far to be happy:
10000 KSP unpreconditioned resid norm 2.273088479279e+03 true resid norm 2.273088479279e+03 ||r(i)||/||b|| 6.561841227018e+02 Reason = -3 ERROR: KSP has not converged. Simulations stopped. [0]PETSC ERROR: --------------------- Error Message -------------------------------------------------------------- [0]PETSC ERROR: KSPSolve has not converged, reason DIVERGED_ITS
Thank you
Hello Pierre, Thank you for the suggestion. However, this is one of the cases where the error is done by the guy sitting in front of the screen! My matrix is no more symmetric and I was pretending to solve it with CG! At least the theory works! 😅
participants (4)
-
Barry Smith -
Edoardo alinovi -
Matthew Knepley -
Pierre Jolivet