I'm trying to program LU decomposition app in PETSc. My idea was, that the program will print the unfactorized matrix, then the factorized matrix and count time taken by a factorization itself.
I've written my code according to little information I've managed to find on the internet ( I used info from this post to initialize my matrix ), but, unfortunately, it isn't enough. My code compiles, but when I try to run it, it simply yells this error:
[0]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------
[0]PETSC ERROR: No support for this operation for this object type
[0]PETSC ERROR: Mat type mpiaij
[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, Sep, 08, 2014
[0]PETSC ERROR: ./petscLUFact on a arch-linux2-c-debug named martin-Aspire-E1-531 by martin Wed Oct 22 22:48:42 2014
[0]PETSC ERROR: Configure options
[0]PETSC ERROR: #1 MatLUFactor() line 2715 in /home/martin/petsc-3.5.2/src/mat/interface/matrix.c
[0]PETSC ERROR: #2 main() line 49 in petscLUFact.c
[0]PETSC ERROR: ----------------End of Error Message -------send entire error message to [email protected]
application called MPI_Abort(MPI_COMM_WORLD, 56) - process 0
===================================================================================
= BAD TERMINATION OF ONE OF YOUR APPLICATION PROCESSES
= EXIT CODE: 56
= CLEANING UP REMAINING PROCESSES
= YOU CAN IGNORE THE BELOW CLEANUP MESSAGES
===================================================================================
This error is caused by MatLUFactor function, which I want to use for in-place LU factorization. The problem is, that I don't know, what exactly is wrong with my code. I think, that the central problem may be in the matrix variable itself, maybe because of a bad allocation ( I think about MatMPIAIJSetPreallocation function ), but I'm not sure.
I tried to replace MatLUFactor functions with MatLUFactorNumeric and MatLUFactorSymbolic, but it works even worse than MatLUFactor, it's error is much "larger" :-)
Finally, I try to launch my program with this command:
mpiexec -n 4 ./petscLUFact
So, if you knew any solution, I'd greatly appreciate it ;-) Thank you!
P.s.: I was looking for the possible solution for a pretty long time, the closest article I've found so far was this maillist, but the questioner uses ParMETIS and SuperLU packages which, as far as I know, I don't use.
This is my source code:
static char help[] = "Reads a PETSc matrix and vector from a file and reorders it.\n\
f0 <input_file> : first file to load (small system)\n\
-f1 <input_file> : second file to load (larger system)\n\n";
#include <petscsys.h>
#include <petscmat.h>
int main( int argc, char **args ) {
Mat A; // 'main' matrix
PetscInt r = 2, c = 2; // matrix dimensions (row x col)
PetscInt i,j; // coordinates
PetscInt Istart, Iend;
PetscInt Ii; // counter
PetscScalar v; // 2-dimensional array ??? (I'm not sure)
PetscErrorCode ierr;
PetscInitialize( &argc, &args, (char*)0, help );
// Create matrix
ierr = MatCreate( PETSC_COMM_WORLD, &A );CHKERRQ( ierr );
ierr = MatSetSizes( A, PETSC_DECIDE, PETSC_DECIDE, r*c, r*c );CHKERRQ( ierr );
ierr = MatSetFromOptions(A);CHKERRQ( ierr );
ierr = MatMPIAIJSetPreallocation( A, 2, PETSC_NULL, 2, PETSC_NULL );CHKERRQ( ierr );
ierr = MatGetOwnershipRange(A,&Istart,&Iend);CHKERRQ( ierr );
// INIT matrix
ierr = MatSetValue( A, 0, 0, 1, INSERT_VALUES ); CHKERRQ( ierr );
ierr = MatSetValue( A, 0, 1, 2, INSERT_VALUES ); CHKERRQ( ierr );
ierr = MatSetValue( A, 1, 0, 3, INSERT_VALUES ); CHKERRQ( ierr );
ierr = MatSetValue( A, 1, 1, 1, INSERT_VALUES ); CHKERRQ( ierr );
ierr = MatAssemblyBegin( A, MAT_FINAL_ASSEMBLY ); CHKERRQ( ierr );
ierr = MatAssemblyEnd( A, MAT_FINAL_ASSEMBLY ); CHKERRQ( ierr );
// Print the matrix
ierr = MatView( A, PETSC_VIEWER_STDOUT_WORLD ); CHKERRQ( ierr );
// -----------------
// LU-decomposition
// -----------------
MatFactorInfo mfi;
// MatFactorInfo mfi INIT
ierr = MatFactorInfoInitialize( &mfi ); CHKERRQ( ierr );
mfi.fill = 2;
mfi.dtcol = 0;
IS rowPerm; // variable for row permutations
IS colPerm; // variable for column permutations
// Possible replace for MatLUFactor
/*
Mat Fact;
ierr = MatLUFactorSymbolic( Fact, A, rowPerm, colPerm, &mfi );
ierr = MatLUFactorNumeric( Fact, A, &mfi );
*/
// I've read somewhere, that zeros are enough for last three
// parameters, but it doesn't work too
//ierr = MatLUFactor( A, 0, 0, 0 ); CHKERRQ( ierr );
ierr = MatLUFactor( A, rowPerm, colPerm, &mfi );
ierr = MatView( A, PETSC_VIEWER_STDOUT_WORLD ); CHKERRQ( ierr );
MatDestroy(&A);
PetscFinalize();
return 0;
}
This is the error when I use MatLUFactorNumeric and MatLUFactorSymbolic functions:
[0]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------
[0]PETSC ERROR: No support for this operation for this object type
[0]PETSC ERROR: Matrix format mpiaij does not have a built-in PETSc LU
[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, Sep, 08, 2014
[0]PETSC ERROR: ./petscLUFact on a arch-linux2-c-debug named martin-Aspire-E1-531 by martin Wed Oct 22 23:40:55 2014
[0]PETSC ERROR: Configure options
[0]PETSC ERROR: #1 MatGetFactor() line 3961 in /home/martin/petsc-3.5.2/src/mat/interface/matrix.c
[0]PETSC ERROR: ------------------------------------------------------------------------
[0]PETSC ERROR: Caught signal number 11 SEGV: Segmentation Violation, probably memory access out of range
[0]PETSC ERROR: Try option -start_in_debugger or -on_error_attach_debugger
[0]PETSC ERROR: or see http://www.mcs.anl.gov/petsc/documentation/faq.html#valgrind[0]PETSC ERROR: or try http://valgrind.org on GNU/linux and Apple Mac OS X to find memory corruption errors
[0]PETSC ERROR: likely location of problem given in stack below
[0]PETSC ERROR: --------------------- Stack Frames ------------------------------------
[0]PETSC ERROR: Note: The EXACT line numbers in the stack are not available,
[0]PETSC ERROR: INSTEAD the line number of the start of the function
[0]PETSC ERROR: is given.
[0]PETSC ERROR: [0] MatLUFactorSymbolic line 2825 /home/martin/petsc-3.5.2/src/mat/interface/matrix.c
[0]PETSC ERROR: [0] MatGetFactor line 3944 /home/martin/petsc-3.5.2/src/mat/interface/matrix.c
[0]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------
[0]PETSC ERROR: Signal received
[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, Sep, 08, 2014
[0]PETSC ERROR: ./petscLUFact on a arch-linux2-c-debug named martin-Aspire-E1-531 by martin Wed Oct 22 23:40:55 2014
[0]PETSC ERROR: Configure options
[0]PETSC ERROR: #2 User provided function() line 0 in unknown file
application called MPI_Abort(MPI_COMM_WORLD, 59) - process 0
===================================================================================
= BAD TERMINATION OF ONE OF YOUR APPLICATION PROCESSES
= EXIT CODE: 59
= CLEANING UP REMAINING PROCESSES
= YOU CAN IGNORE THE BELOW CLEANUP MESSAGES
===================================================================================
Here is the error message:
[0]PETSC ERROR: No support for this operation for this object type
[0]PETSC ERROR: Mat type mpiaij
So there is no support for this operation, meaning LU factorization, for Mat type mpiaij. Now we can appeal to the online documentation, which does say that the built-in LU works only for sequential matrices (seqaij). It also shows the available packages, like SuperLU and MUMPS, that can be used with PETSc to do the LU factorization in parallel. These can both be automatically installed using
--download-superlu --download-mumps