Direct LU solver
Hong Zhang
hzhang at mcs.anl.gov
Thu Feb 28 14:57:24 CST 2008
or ~petsc/src/ksp/ksp/examples/tutorials/ex5.c to avoid using
matrix data file, e.g.
mpiexec -np 2 ./ex5 -ksp_type preonly -pc_type lu -mat_type superlu_dist
Hong
On Thu, 28 Feb 2008, Matthew Knepley wrote:
> I would recommend using KSP ex10 and then customizing with options.
> That way we know it should work. For SuperLU_dist,
>
> -ksp_type preonly -pc_type lu -mat_type superlu_dist
>
> Matt
>
> On Thu, Feb 28, 2008 at 1:07 PM, <Amit.Itagi at seagate.com> wrote:
>> Hi,
>>
>> I need to do direct LU solves (repeatedly, with the same matrix) in one of
>> my MPI applications. I am having trouble implementing the solver. To
>> identify the problem, I wrote a short toy code to run with 2 processes. I
>> can run it with either the spooles parallel matrix or the superlu_dist
>> matrix. I am using C++ and complex matrices.
>>
>> Here is the code listing:
>>
>> #include <iostream>
>> #include <complex>
>> #include <stdlib.h>
>> #include "petsc.h"
>> #include "petscmat.h"
>> #include "petscvec.h"
>> #include "petscksp.h"
>>
>> using namespace std;
>>
>> int main( int argc, char *argv[] ) {
>>
>> int rank, size;
>> Mat A;
>> PetscErrorCode ierr;
>> PetscInt loc;
>> PetscScalar val;
>> Vec x, y;
>> KSP solver;
>> PC prec;
>> MPI_Comm comm;
>>
>> // Number of non-zeros in each row
>> int d_nnz=1, o_nnz=1;
>>
>> ierr=PetscInitialize(&argc,&argv,PETSC_NULL,PETSC_NULL); CHKERRQ(ierr);
>> // Initialization (including MPI)
>>
>> comm=PETSC_COMM_WORLD;
>>
>> ierr=MPI_Comm_size(PETSC_COMM_WORLD,&size); CHKERRQ(ierr);
>> ierr=MPI_Comm_rank(PETSC_COMM_WORLD,&rank); CHKERRQ(ierr);
>>
>> // Assemble matrix A
>>
>> if(rank==0) {
>> ierr=MatCreateMPIAIJ(comm,1,1,2,2,0,&d_nnz,0,&o_nnz,&A); CHKERRQ(ierr);
>> val=complex<double>(1.0,0.0);
>> ierr=MatSetValue(A,0,0,val,INSERT_VALUES);CHKERRQ(ierr);
>> val=complex<double>(0.0,1.0);
>> ierr=MatSetValue(A,0,1,val,INSERT_VALUES);CHKERRQ(ierr);
>> }
>> else {
>> ierr=MatCreateMPIAIJ(PETSC_COMM_WORLD,1,1,2,2,0,&d_nnz,0,&o_nnz,&A);
>> CHKERRQ(ierr);
>> val=complex<double>(1.0,1.0);
>> ierr=MatSetValue(A,1,1,val,INSERT_VALUES);CHKERRQ(ierr);
>> val=complex<double>(0.0,-1.0);
>> ierr=MatSetValue(A,1,0,val,INSERT_VALUES);CHKERRQ(ierr);
>> }
>>
>> ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY); CHKERRQ(ierr);
>> ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY); CHKERRQ(ierr);
>>
>> cout << "============ Mat A ==================" << endl;
>> ierr=MatView(A,PETSC_VIEWER_STDOUT_WORLD); CHKERRQ(ierr);
>> cout << "======================================" << endl;
>>
>> // For spooles
>> //ierr=MatConvert(A,MATMPIAIJSPOOLES,MAT_REUSE_MATRIX,&A); CHKERRQ(ierr);
>>
>> // For superlu_dist
>> ierr=MatConvert(A,MATSUPERLU_DIST,MAT_REUSE_MATRIX,&A); CHKERRQ(ierr);
>>
>> // Direct LU solver
>> ierr=KSPCreate(comm,&solver); CHKERRQ(ierr);
>> ierr=KSPSetType(solver,KSPPREONLY); CHKERRQ(ierr);
>> ierr=KSPSetOperators(solver,A,A,SAME_NONZERO_PATTERN); CHKERRQ(ierr);
>> ierr=KSPGetPC(solver,&prec); CHKERRQ(ierr);
>> ierr=PCSetType(prec,PCLU); CHKERRQ(ierr);
>> ierr=KSPSetFromOptions(solver); CHKERRQ(ierr);
>>
>> //============ Vector assembly ========================
>>
>> if(rank==0) {
>> ierr=VecCreateMPI(PETSC_COMM_WORLD,1,2,&x); CHKERRQ(ierr);
>> val=complex<double>(1.0,0.0);
>> loc=0;
>> ierr=VecSetValues(x,1,&loc,&val,ADD_VALUES);CHKERRQ(ierr);
>> }
>> else {
>> ierr=VecCreateMPI(PETSC_COMM_WORLD,1,2,&x); CHKERRQ(ierr);
>> val=complex<double>(-1.0,0.0);
>> loc=1;
>> ierr=VecSetValues(x,1,&loc,&val,ADD_VALUES);CHKERRQ(ierr);
>> }
>>
>> ierr=VecAssemblyBegin(x); CHKERRQ(ierr);
>> ierr=VecAssemblyEnd(x); CHKERRQ(ierr);
>>
>> cout << "============== Vec x ==================" << endl;
>> ierr=VecView(x,PETSC_VIEWER_STDOUT_WORLD); CHKERRQ(ierr);
>> cout << "======================================" << endl;
>>
>> VecDuplicate(x,&y); // Duplicate the matrix storage
>>
>> // Solve the matrix equation
>> ierr=KSPSolve(solver,x,y); CHKERRQ(ierr);
>>
>> cout << "============== Vec y =================" << endl;
>> ierr=VecView(y,PETSC_VIEWER_STDOUT_WORLD); CHKERRQ(ierr);
>> cout << "======================================" << endl;
>>
>>
>> // Destructors
>> ierr=KSPDestroy(solver); CHKERRQ(ierr);
>> ierr=VecDestroy(x); CHKERRQ(ierr);
>> ierr=VecDestroy(y); CHKERRQ(ierr);
>> ierr=MatDestroy(A); CHKERRQ(ierr);
>>
>> // Finalize
>> ierr=PetscFinalize(); CHKERRQ(ierr);
>>
>>
>> return 0;
>>
>> }
>>
>>
>> When I run the program with spooles, I get the following output.
>>
>>
>> ============ Mat A ==================
>> ============ Mat A ==================
>> ======================================
>> row 0: (0, 1) (1, 0 + 1 i)
>> row 1: (0, 0 - 1 i) (1, 1 + 1 i)
>> ======================================
>> ============== Vec x ==================
>> ============== Vec x ==================
>> Process [0]
>> 1
>> ======================================
>> Process [1]
>> -1
>> ======================================
>>
>> fatal error in InpMtx_MPI_split()
>> firsttag = 0, tagbound = -1
>>
>> fatal error in InpMtx_MPI_split()
>> firsttag = 0, tagbound = -1
>> -----------------------------------------------------------------------------
>> One of the processes started by mpirun has exited with a nonzero exit
>> code. This typically indicates that the process finished in error.
>> If your process did not finish in error, be sure to include a "return
>> 0" or "exit(0)" in your C code before exiting the application.
>>
>> PID 22881 failed on node n0 (127.0.0.1) with exit status 255.
>> -----------------------------------------------------------------------------
>> mpirun failed with exit status 255
>>
>>
>> When run in the debugger, there is no stack trace. The error is
>>
>> Program exited with code 0377
>>
>>
>> With superlu_dist, the output is
>>
>> ============ Mat A ==================
>> ============ Mat A ==================
>> row 0: (0, 1) (1, 0 + 1 i)
>> row 1: (0, 0 - 1 i) (1, 1 + 1 i)
>> ======================================
>> ======================================
>> ============== Vec x ==================
>> ============== Vec x ==================
>> Process [0]
>> 1
>> Process [1]
>> -1
>> ======================================
>> ======================================
>> [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/petsc-as/documentation/troubleshooting.html#Signal[0]PETSC
>> ERROR: or try http://valgrind.org on linux or man libgmalloc on Apple to
>> find memory corruption errors
>> [0]PETSC ERROR: configure using --with-debugging=yes, recompile, link, and
>> run
>> [0]PETSC ERROR: to get more information on the crash.
>> [0]PETSC ERROR: --------------------- Error Message
>> ------------------------------------
>> [0]PETSC ERROR: Signal received!
>> [0]PETSC ERROR:
>> ------------------------------------------------------------------------
>> [0]PETSC ERROR: Petsc Release Version 2.3.3, Patch 8, Fri Nov 16 17:03:40
>> CST 2007 HG revision: 414581156e67e55c761739b0deb119f7590d0f4b
>> [0]PETSC ERROR: See docs/changes/index.html for recent updates.
>> [0]PETSC ERROR: See docs/faq.html for hints about trouble shooting.
>> [0]PETSC ERROR: See docs/index.html for manual pages.
>> [0]PETSC ERROR:
>> ------------------------------------------------------------------------
>> [0]PETSC ERROR: lu on a linux-gnu named tabla by amit Thu Feb 28 13:48:54
>> 2008
>> [0]PETSC ERROR: Libraries linked from
>> /home/amit/programs/ParEM/petsc-2.3.3-p8/lib/linux-gnu-c-debug
>> [0]PETSC ERROR: Configure run at Thu Feb 28 12:19:39 2008
>> [0]PETSC ERROR: Configure options --with-scalar-type=complex
>> --with-debugging=no --with-fortran-kernels=generic --with-clanguage=cxx
>> --with-metis=1 --download-metis=1 --with-parmetis=1 --download-parmetis=1
>> --with-superlu_dist=1 --download-superlu_dist=1 --with-spooles=1
>> --with-spooles-dir=/home/amit/programs/ParEM/spooles-2.2 COPTFLAGS="-O3
>> -march=p4 -mtune=p4 -ffast-math -malign-double -funroll-loops -pipe
>> -fomit-frame-pointer -finline-functions -msse2" CXXOPTFLAGS="-O3 -march=p4
>> -mtune=p4 -ffast-math -malign-double -funroll-loops -pipe
>> -fomit-frame-pointer -finline-functions -msse2" FOPTS="-O3 -qarch=p4
>> -qtune=p4" --with-shared=0
>> [0]PETSC ERROR:
>> ------------------------------------------------------------------------
>> [0]PETSC ERROR: User provided function() line 0 in unknown directory
>> unknown file
>> -----------------------------------------------------------------------------
>> One of the processes started by mpirun has exited with a nonzero exit
>> code. This typically indicates that the process finished in error.
>> If your process did not finish in error, be sure to include a "return
>> 0" or "exit(0)" in your C code before exiting the application.
>>
>> PID 22998 failed on node n0 (127.0.0.1) with exit status 1.
>> -----------------------------------------------------------------------------
>> mpirun failed with exit status 1
>>
>>
>> The debugger tracks the segmentation violation to
>>
>>
>> main -> KSPSolve -> KSPSetUp -> PCSetUp -> PCSetUp_LU -> MatLUFactorNumeric
>> -> MatLUFactorNumeric_SUperLU_DIST -> pzgssvx
>>
>>
>> Could someone kindly point out what I am missing ?
>>
>>
>> Thanks
>>
>> Rgds,
>> Amit
>>
>>
>
>
>
> --
> 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
>
>
More information about the petsc-users
mailing list