> On Oct 30, 2014, at 3:36 AM, Florian Lindner <mailingli...@xgm.de> wrote:
>
> Am Mittwoch, 29. Oktober 2014, 13:00:17 schrieb Barry Smith:
>>
>> MatAssembly flushes out any preallocated space that was not used. Hence
>> when you try to put a diagonal entry in later it is as if you did not
>> preallocate. You should set zeros in the diagonal as part of your initial
>> setting values in the matrix before calling MatAssembly.
>
> Even with FLUSH_ASSEMBLY?
No, flush assembly does not clean out any unused preallocated space. final
assembly does
Barry
> I had this idea too, but was not able to reproduce the error I have in my
> application in a small example.
>
> Why does this program works? Shouldn't the MatDiagonalSet produce an
> error/warning that there were mallocs involved? Since the preallocation was
> lost?
>
> Thanks,
> Florian
>
> #include <iostream>
> #include "petscmat.h"
> #include "petscviewer.h"
> #include "petscksp.h"
>
> using namespace std;
>
> // Compiling with: mpic++ -g3 -Wall -I ~/software/petsc/include -I
> ~/software/petsc/arch-linux2-c-debug/include -L
> ~/software/petsc/arch-linux2-c-debug/lib -lpetsc test.cpp
>
> int main(int argc, char **args)
> {
> PetscInitialize(&argc, &args, "", NULL);
>
> PetscErrorCode ierr = 0;
> int N = 9;
> Mat matrix;
>
> ierr = MatCreate(PETSC_COMM_SELF, &matrix);
> ierr = MatSetType(matrix, MATSBAIJ); CHKERRQ(ierr);
> ierr = MatSetSizes(matrix, PETSC_DECIDE, PETSC_DECIDE, N, N); CHKERRQ(ierr);
> ierr = MatSetOption(matrix, MAT_SYMMETRY_ETERNAL, PETSC_TRUE); CHKERRQ(ierr);
>
> PetscInt nnz[N];
> nnz[0] = 2; nnz[1] = 1; nnz[2] = 1;
> nnz[3] = 1; nnz[4] = 1; nnz[5] = 1;
> nnz[6] = 1; nnz[7] = 1; nnz[8] = 1;
>
> for (int r = 0; r < N; r++) {
> cout << "Preallocated row " << r << " with " << nnz[r] << " elements." <<
> endl;
> }
> MatSeqSBAIJSetPreallocation(matrix, 1, PETSC_DEFAULT, nnz);
>
> MatSetValue(matrix, 0, 5, 10, INSERT_VALUES);
>
> ierr = MatAssemblyBegin(matrix, MAT_FLUSH_ASSEMBLY); CHKERRQ(ierr);
> ierr = MatAssemblyEnd(matrix, MAT_FLUSH_ASSEMBLY); CHKERRQ(ierr);
> Vec zeros;
> MatGetVecs(matrix, NULL, &zeros);
>
> MatDiagonalSet(matrix, zeros, ADD_VALUES);
> ierr = MatAssemblyBegin(matrix, MAT_FINAL_ASSEMBLY); CHKERRQ(ierr);
> ierr = MatAssemblyEnd(matrix, MAT_FINAL_ASSEMBLY); CHKERRQ(ierr);
>
> ierr = MatView(matrix, PETSC_VIEWER_STDOUT_WORLD); CHKERRQ(ierr);
>
> MatDestroy(&matrix);
> PetscFinalize();
> return 0;
> }
>
>
>
>
>>
>> Barry
>>
>>> On Oct 29, 2014, at 8:00 AM, Florian Lindner <mailingli...@xgm.de> wrote:
>>>
>>> Hello,
>>>
>>> I try to preallocate a sparse matrix like it was recommended in another
>>> posting, but get an error which kind of surprises me. Somehow I think it
>>> might be related to the order of assembly calls...
>>>
>>> My code creates the matrix:
>>>
>>> MatSetType(_matrixC.matrix, MATSBAIJ);
>>> MatSetSizes(_matrixC.matrix, PETSC_DECIDE, PETSC_DECIDE, n, n);
>>> MatSetOption(_matrixC.matrix, MAT_SYMMETRY_ETERNAL, PETSC_TRUE);
>>>
>>> Than I enter a loop to get the number of elements per row. The entire code
>>> is a bit too long, but right before the SetPreallocaiton I do:
>>>
>>> for (int r = 0; r < n; r++) {
>>> cout << "Preallocated row " << r << " with " << nnz[r] << " elements."
>>> << endl;
>>> }
>>>
>>> which results in the output:
>>>
>>> Matrix Size = 9
>>> Preallocated row 0 with 9 elements.
>>> Preallocated row 1 with 8 elements.
>>> Preallocated row 2 with 7 elements.
>>> Preallocated row 3 with 6 elements.
>>> Preallocated row 4 with 5 elements.
>>> Preallocated row 5 with 4 elements.
>>> Preallocated row 6 with 1 elements.
>>> Preallocated row 7 with 1 elements.
>>> Preallocated row 8 with 1 elements.
>>>
>>> at least in rows 0 to 5 it's a dense upper triangular matrix (incl. main
>>> diagonal).
>>>
>>> ierr = MatSeqSBAIJSetPreallocation(_matrixC.matrix, 1, PETSC_DEFAULT, nnz);
>>>
>>> Now more or less the same loops as above with MatSetValues call for each
>>> row.
>>>
>>> Since the KSPSolver requires that the diagonal is set, even if it's merely
>>> with zeros I do:
>>>
>>> _matrixC.assemble(MAT_FLUSH_ASSEMBLY) // my helper function, but I'm sure
>>> you get the point ;-)
>>> petsc::Vector zeros(_matrixC);
>>> MatDiagonalSet(_matrixC.matrix, zeros.vector, ADD_VALUES);
>>> ierr = MatAssemblyBegin(_matrixC.matrix, MAT_FINAL_ASSEMBLY);
>>> CHKERRV(ierr);
>>> [ some more code ]
>>> ierr = MatAssemblyEnd(_matrixC.matrix, MAT_FINAL_ASSEMBLY); CHKERRV(ierr);
>>>
>>> But when I run the code I get an error message:
>>>
>>> [0]PETSC ERROR: --------------------- Error Message
>>> --------------------------------------------------------------
>>> [0]PETSC ERROR: Argument out of range
>>> [0]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
>>> [0]PETSC ERROR: See http://www.mcs.anl.gov/petsc/documentation/faq.html for
>>> trouble shooting.
>>> [0]PETSC ERROR: Petsc Release Version 3.5.2, unknown
>>> [0]PETSC ERROR: ./petBench on a arch-linux2-c-debug named helium by
>>> lindnefn Wed Oct 29 13:48:49 2014
>>> [0]PETSC ERROR: Configure options --with-debugging=1
>>> [0]PETSC ERROR: #1 MatSetValues_SeqSBAIJ() line 977 in
>>> /data2/scratch/lindner/petsc/src/mat/impls/sbaij/seq/sbaij.c
>>> [0]PETSC ERROR: #2 MatSetValues() line 1135 in
>>> /data2/scratch/lindner/petsc/src/mat/interface/matrix.c
>>> [0]PETSC ERROR: #3 MatDiagonalSet_Default() line 203 in
>>> /data2/scratch/lindner/petsc/src/mat/utils/axpy.c
>>> [0]PETSC ERROR: #4 MatDiagonalSet() line 241 in
>>> /data2/scratch/lindner/petsc/src/mat/utils/axpy.c
>>>
>>> I'm surprised that I get this error message since I do have preallocated 9
>>> elements and I do not understand why inserting an element at (0, 0) could
>>> be a problem / needing a malloc. Has the preallocation got lost somewhere?
>>>
>>> Thanks once again....
>>>
>>> Florian
>>
