On Apr 17, 2013, at 1:53 PM, "Jin, Shuangshuang" <Shuangshuang.Jin at pnnl.gov>
wrote:
> Hello, Barry, I got a new question here.
>
> I have integrated the AX=B solver into my code using the example from
> ex125.c. The only major difference between my code and the example is the
> input matrix.
>
> In the example ex125.c, the A matrix is created by loading the data from a
> binary file. In my code, I passed in the computed A matrix and the B matrix
> as arguments to a user defined solvingAXB() function:
>
> static PetscErrorCode solvingAXB(Mat A, Mat B, PetscInt n, PetscInt nrhs, Mat
> X, const int me);
>
> I noticed that A has to be stored as an MATAIJ matrix and the B as a MATDENSE
> matrix. So I converted their storage types inside the solvingAXB as below
> (They were MPIAIJ type before):
Note that MATAIJ is MATMPIAIJ when rank > 1 and MATSEQAIJ when rank == 1.
>
> ierr = MatConvert(A, MATAIJ, MAT_REUSE_MATRIX, &A); // A has to be MATAIJ
> for built-in PETSc LU!
> MatView(A, PETSC_VIEWER_STDOUT_WORLD);
>
> ierr = MatConvert(B, MATDENSE, MAT_REUSE_MATRIX, &B); // B has to be a
> SeqDense matrix!
> MatView(B, PETSC_VIEWER_STDOUT_WORLD);
>
> With this implementation, I can run the code with expected results when 1
> processor is used.
>
> However, when I use 2 processors, I get a run time error which I guess is
> still somehow related to the matrix format but cannot figure out how to fix
> it. Could you please take a look at the following error message?
PETSc doesn't have a parallel LU solver, you need to install PETSc with
SuperLU_DIST or MUMPS (preferably both) and then pass in one of those names
when you call MatGetFactor().
Barry
>
> [d3m956 at olympus ss08]$ mpiexec -n 2 dynSim -i 3g9b.txt
> Matrix Object: 1 MPI processes
> type: mpiaij
> row 0: (0, 0 - 16.4474 i) (3, 0 + 17.3611 i)
> row 1: (1, 0 - 8.34725 i) (7, 0 + 16 i)
> row 2: (2, 0 - 5.51572 i) (5, 0 + 17.0648 i)
> row 3: (0, 0 + 17.3611 i) (3, 0 - 39.9954 i) (4, 0 + 10.8696 i) (8, 0 +
> 11.7647 i)
> row 4: (3, 0 + 10.8696 i) (4, 0.934 - 17.0633 i) (5, 0 + 5.88235 i)
> row 5: (2, 0 + 17.0648 i) (4, 0 + 5.88235 i) (5, 0 - 32.8678 i) (6, 0 +
> 9.92063 i)
> row 6: (5, 0 + 9.92063 i) (6, 1.03854 - 24.173 i) (7, 0 + 13.8889 i)
> row 7: (1, 0 + 16 i) (6, 0 + 13.8889 i) (7, 0 - 36.1001 i) (8, 0 + 6.21118 i)
> row 8: (3, 0 + 11.7647 i) (7, 0 + 6.21118 i) (8, 1.33901 - 18.5115 i)
> Matrix Object: 1 MPI processes
> type: mpidense
> 0.0000000000000000e+00 + -1.6447368421052630e+01i 0.0000000000000000e+00 +
> 0.0000000000000000e+00i 0.0000000000000000e+00 + 0.0000000000000000e+00i
> 0.0000000000000000e+00 + 0.0000000000000000e+00i 0.0000000000000000e+00 +
> -8.3472454090150254e+00i 0.0000000000000000e+00 + 0.0000000000000000e+00i
> 0.0000000000000000e+00 + 0.0000000000000000e+00i 0.0000000000000000e+00 +
> 0.0000000000000000e+00i 0.0000000000000000e+00 + -5.5157198014340878e+00i
> 0.0000000000000000e+00 + 0.0000000000000000e+00i 0.0000000000000000e+00 +
> 0.0000000000000000e+00i 0.0000000000000000e+00 + 0.0000000000000000e+00i
> 0.0000000000000000e+00 + 0.0000000000000000e+00i 0.0000000000000000e+00 +
> 0.0000000000000000e+00i 0.0000000000000000e+00 + 0.0000000000000000e+00i
> 0.0000000000000000e+00 + 0.0000000000000000e+00i 0.0000000000000000e+00 +
> 0.0000000000000000e+00i 0.0000000000000000e+00 + 0.0000000000000000e+00i
> 0.0000000000000000e+00 + 0.0000000000000000e+00i 0.0000000000000000e+00 +
> 0.0000000000000000e+00i 0.0000000000000000e+00 + 0.0000000000000000e+00i
> 0.0000000000000000e+00 + 0.0000000000000000e+00i 0.0000000000000000e+00 +
> 0.0000000000000000e+00i 0.0000000000000000e+00 + 0.0000000000000000e+00i
> 0.0000000000000000e+00 + 0.0000000000000000e+00i 0.0000000000000000e+00 +
> 0.0000000000000000e+00i 0.0000000000000000e+00 + 0.0000000000000000e+00i
> PETSC LU:
> [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:
> ------------------------------------------------------------------------
> [0]PETSC ERROR: Petsc Development HG revision:
> 6e0ddc6e9b6d8a9d8eb4c0ede0105827a6b58dfb HG Date: Mon Mar 11 22:54:30 2013
> -0500
> [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: dynSim on a arch-complex named olympus.local by d3m956 Wed
> Apr 17 11:38:12 2013
> [0]PETSC ERROR: Libraries linked from
> /pic/projects/ds/petsc-dev/arch-complex/lib
> [0]PETSC ERROR: Configure run at Tue Mar 12 14:32:37 2013
> [0]PETSC ERROR: Configure options --with-scalar-type=complex
> --with-clanguage=C++ PETSC_ARCH=arch-complex --with-fortran-kernels=generic
> [0]PETSC ERROR:
> ------------------------------------------------------------------------
> [0]PETSC ERROR: MatGetFactor() line 3949 in src/mat/interface/matrix.c
> [0]PETSC ERROR: solvingAXB() line 880 in "unknowndirectory/"reducedYBus.C
> End!
> [1]PETSC ERROR:
> ------------------------------------------------------------------------
> [1]PETSC ERROR: Caught signal number 15 Terminate: Somet process (or the
> batch system) has told this process to end
> [1]PETSC ERROR: Try option -start_in_debugger or -on_error_attach_debugger
> [1]PETSC ERROR: or see
> http://www.mcs.anl.gov/petsc/documentation/faq.html#valgrind[1]PETSC ERROR:
> or try http://valgrind.org on GNU/linux and Apple Mac OS X to find memory
> corruption errors
> [1]PETSC ERROR: likely location of problem given in stack below
> [1]PETSC ERROR: --------------------- Stack Frames
> ------------------------------------
> [1]PETSC ERROR: Note: The EXACT line numbers in the stack are not available,
> [1]PETSC ERROR: INSTEAD the line number of the start of the function
> [1]PETSC ERROR: is given.
> [1]PETSC ERROR: [1] PetscSleep line 35 src/sys/utils/psleep.c
> [1]PETSC ERROR: [1] PetscTraceBackErrorHandler line 172
> src/sys/error/errtrace.c
> [1]PETSC ERROR: [1] PetscError line 361 src/sys/error/err.c
> [1]PETSC ERROR: [1] MatGetFactor line 3932 src/mat/interface/matrix.c
> [1]PETSC ERROR: --------------------- Error Message
> ------------------------------------
> [1]PETSC ERROR: Signal received!
> [1]PETSC ERROR:
> ------------------------------------------------------------------------
> [1]PETSC ERROR: Petsc Development HG revision:
> 6e0ddc6e9b6d8a9d8eb4c0ede0105827a6b58dfb HG Date: Mon Mar 11 22:54:30 2013
> -0500
> [1]PETSC ERROR: See docs/changes/index.html for recent updates.
> [1]PETSC ERROR: See docs/faq.html for hints about trouble shooting.
> [1]PETSC ERROR: See docs/index.html for manual pages.
> [1]PETSC ERROR:
> ------------------------------------------------------------------------
> [1]PETSC ERROR: dynSim on a arch-complex named olympus.local by d3m956 Wed
> Apr 17 11:38:12 2013
> [1]PETSC ERROR: Libraries linked from
> /pic/projects/ds/petsc-dev/arch-complex/lib
> [1]PETSC ERROR: Configure run at Tue Mar 12 14:32:37 2013
> [1]PETSC ERROR: Configure options --with-scalar-type=complex
> --with-clanguage=C++ PETSC_ARCH=arch-complex --with-fortran-kernels=generic
> [1]PETSC ERROR:
> ------------------------------------------------------------------------
> [1]PETSC ERROR: User provided function() line 0 in unknown directory unknown
> file
> --------------------------------------------------------------------------
>
> Thanks,
> Shuangshuang
>
>
> -----Original Message-----
> From: petsc-users-bounces at mcs.anl.gov [mailto:petsc-users-bounces at
> mcs.anl.gov] On Behalf Of Barry Smith
> Sent: Tuesday, April 16, 2013 5:50 PM
> To: PETSc users list
> Subject: Re: [petsc-users] ksp for AX=B system
>
>
> Shuangshuang
>
> This is what I was expecting, thanks for the confirmation. For these size
> problems you definitely want to use a direct solver (often parallel but not
> for smaller matrices) and solve multiple right hand sides. This means you
> actually will not use the KSP solver that is standard for most PETSc work,
> instead you will work directly with the MatGetFactor(), MatGetOrdering(),
> MatLUFactorSymbolic(), MatLUFactorNumeric(), MatMatSolve() paradigm where the
> A matrix is stored as an MATAIJ matrix and the B (multiple right hand side)
> as a MATDENSE matrix.
>
> An example that displays this paradigm is src/mat/examples/tests/ex125.c
>
> Once you have something running of interest to you we would like to work
> with you to improve the performance, we have some "tricks" we haven't yet
> implemented to make these solvers much faster than they will be by default.
>
> Barry
>
>
>
> On Apr 16, 2013, at 7:38 PM, "Jin, Shuangshuang" <Shuangshuang.Jin at
> pnnl.gov> wrote:
>
>> Hi, Barry, thanks for your prompt reply.
>>
>> We have various size problems, from (n=9, m=3), (n=1081, m = 288), to
>> (n=16072, m=2361) or even larger ultimately.
>>
>> Usually the dimension of the square matrix A, n is much larger than the
>> column dimension of B, m.
>>
>> As you said, I'm using the loop to deal with the small (n=9, m=3) case.
>> However, for bigger problem, I do hope there's a better approach.
>>
>> This is a power flow problem. When we parallelized it in OpenMP previously,
>> we just parallelized the outside loop, and used a direct solver to solve it.
>>
>> We are now switching to MPI and like to use PETSc ksp solver to solve it in
>> parallel. However, I don't know what's the best ksp solver I should use
>> here? Direct solver or try a preconditioner?
>>
>> I appreciate very much for your recommendations.
>>
>> Thanks,
>> Shuangshuang
>>
>>
