[petsc-users] parallel computing error
Matthew Knepley
knepley at gmail.com
Fri May 5 04:21:37 CDT 2023
On Fri, May 5, 2023 at 3:49 AM 권승리 / 학생 / 항공우주공학과 <ksl7912 at snu.ac.kr>
wrote:
> Dear Barry Smith
>
> Thanks to you, I knew the difference between MATAIJ and MATDENSE.
>
> However, I still have some problems.
>
> There is no problem when I run with a single core. But, MatGetFactor error
> occurs when using multi-core.
>
> Could you give me some advice?
>
> The error message is
>
> [0]PETSC ERROR: --------------------- Error Message
> --------------------------------------------------------------
> [0]PETSC ERROR: See https://petsc.org/release/overview/linear_solve_table/
> for possible LU and Cholesky solvers
> [0]PETSC ERROR: MatSolverType petsc does not support matrix type mpidense
>
PETSc uses 3rd party packages for parallel dense factorization. You would
need to reconfigure with either ScaLAPACK
or Elemental.
Thanks,
Matt
> [0]PETSC ERROR: See https://petsc.org/release/faq/ for trouble shooting.
> [0]PETSC ERROR: Petsc Release Version 3.18.5, unknown
> [0]PETSC ERROR: ./app on a arch-linux-c-opt named ubuntu by ksl Fri May 5
> 00:35:23 2023
> [0]PETSC ERROR: Configure options --download-mpich --with-debugging=0
> COPTFLAGS="-O3 -march=native -mtune=native" CXXOPTFLAGS="-O3 -march=native
> -mtune=native" FOPTFLAGS="-O3 -march=native -mtune=native" --download-mumps
> --download-scalapack --download-parmetis --download-metis
> --download-parmetis --download-hpddm --download-slepc
> [0]PETSC ERROR: #1 MatGetFactor() at
> /home/ksl/petsc/src/mat/interface/matrix.c:4757
> [0]PETSC ERROR: #2 main() at /home/ksl/Downloads/coding_test/coding/a1.c:66
> [0]PETSC ERROR: No PETSc Option Table entries
> [0]PETSC ERROR: ----------------End of Error Message -------send entire
> error message to petsc-maint at mcs.anl.gov----------
> application called MPI_Abort(MPI_COMM_SELF, 92) - process 0
>
> My code is below:
>
> int main(int argc, char** args)
> {
> Mat A, E, A_temp, A_fac;
> int n = 15;
> PetscInitialize(&argc, &args, NULL, NULL);
> PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
>
> PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
> PetscCall(MatSetType(A,MATDENSE));
> PetscCall(MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, n, n));
> PetscCall(MatSetFromOptions(A));
> PetscCall(MatSetUp(A));
> // Insert values
> double val;
> for (int i = 0; i < n; i++) {
> for (int j = 0; j < n; j++) {
> if (i == j){
> val = 2.0;
> }
> else{
> val = 1.0;
> }
> PetscCall(MatSetValue(A, i, j, val, INSERT_VALUES));
> }
> }
> PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
> PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
>
> // Make Identity matrix
> PetscCall(MatCreate(PETSC_COMM_WORLD, &E));
> PetscCall(MatSetType(E,MATDENSE));
> PetscCall(MatSetSizes(E, PETSC_DECIDE, PETSC_DECIDE, n, n));
> PetscCall(MatSetFromOptions(E));
> PetscCall(MatSetUp(E));
> PetscCall(MatShift(E,1.0));
> PetscCall(MatAssemblyBegin(E, MAT_FINAL_ASSEMBLY));
> PetscCall(MatAssemblyEnd(E, MAT_FINAL_ASSEMBLY));
>
> PetscCall(MatDuplicate(A, MAT_DO_NOT_COPY_VALUES, &A_temp));
> PetscCall(MatGetFactor(A, MATSOLVERPETSC, MAT_FACTOR_LU, &A_fac));
>
> IS isr, isc; MatFactorInfo info;
> MatGetOrdering(A, MATORDERINGNATURAL, &isr, &isc);
> PetscCall(MatLUFactorSymbolic(A_fac, A, isr, isc, &info));
> PetscCall(MatLUFactorNumeric(A_fac, A, &info));
> MatMatSolve(A_fac, E, A_temp);
>
> PetscCall(MatView(A_temp, PETSC_VIEWER_STDOUT_WORLD));
> MatDestroy(&A);
> MatDestroy(&A_temp);
> MatDestroy(&A_fac);
> MatDestroy(&E);
> PetscCall(PetscFinalize());
> }
>
> Best regards
> Seung Lee Kwon
>
> 2023년 5월 4일 (목) 오후 10:19, Barry Smith <bsmith at petsc.dev>님이 작성:
>
>>
>> The code in ex125.c contains
>>
>> PetscCall(MatCreate(PETSC_COMM_WORLD, &C));
>> PetscCall(MatSetOptionsPrefix(C, "rhs_"));
>> PetscCall(MatSetSizes(C, m, PETSC_DECIDE, PETSC_DECIDE, nrhs));
>> PetscCall(MatSetType(C, MATDENSE));
>> PetscCall(MatSetFromOptions(C));
>> PetscCall(MatSetUp(C));
>>
>> This dense parallel matrix is suitable for passing to MatMatSolve() as
>> the right-hand side matrix. Note it is created with PETSC_COMM_WORLD and
>> its type is set to be MATDENSE.
>>
>> You may need to make a sample code by stripping out all the excess code
>> in ex125.c to just create an MATAIJ and MATDENSE and solves with
>> MatMatSolve() to determine why you code does not work.
>>
>>
>>
>> On May 4, 2023, at 3:20 AM, 권승리 / 학생 / 항공우주공학과 <ksl7912 at snu.ac.kr>
>> wrote:
>>
>> Dear Barry Smith
>>
>> Thank you for your reply.
>>
>> I've already installed MUMPS.
>>
>> And I checked the example you said (ex125.c), I don't understand why the
>> RHS matrix becomes the SeqDense matrix.
>>
>> Could you explain in more detail?
>>
>> Best regards
>> Seung Lee Kwon
>>
>> 2023년 5월 4일 (목) 오후 12:08, Barry Smith <bsmith at petsc.dev>님이 작성:
>>
>>>
>>> You can configure with MUMPS ./configure --download-mumps
>>> --download-scalapack --download-ptscotch --download-metis
>>> --download-parmetis
>>>
>>> And then use MatMatSolve() as in src/mat/tests/ex125.c with parallel
>>> MatMatSolve() using MUMPS as the solver.
>>>
>>> Barry
>>>
>>>
>>> On May 3, 2023, at 10:29 PM, 권승리 / 학생 / 항공우주공학과 <ksl7912 at snu.ac.kr>
>>> wrote:
>>>
>>> Dear developers
>>>
>>> Thank you for your explanation.
>>>
>>> But I should use the MatCreateSeqDense because I want to use the
>>> MatMatSolve that B matrix must be a SeqDense matrix.
>>>
>>> Using MatMatSolve is an inevitable part of my code.
>>>
>>> Could you give me a comment to avoid this error?
>>>
>>> Best,
>>>
>>> Seung Lee Kwon
>>>
>>> 2023년 5월 3일 (수) 오후 7:30, Matthew Knepley <knepley at gmail.com>님이 작성:
>>>
>>>> On Wed, May 3, 2023 at 6:05 AM 권승리 / 학생 / 항공우주공학과 <ksl7912 at snu.ac.kr>
>>>> wrote:
>>>>
>>>>> Dear developers
>>>>>
>>>>> I'm trying to use parallel computing and I ran the command 'mpirun -np
>>>>> 4 ./app'
>>>>>
>>>>> In this case, there are two problems.
>>>>>
>>>>> *First,* I encountered error message
>>>>> ///
>>>>> [0]PETSC ERROR: [1]PETSC ERROR: --------------------- Error Message
>>>>> --------------------------------------------------------------
>>>>> [1]PETSC ERROR: Invalid argument
>>>>> [1]PETSC ERROR: Comm must be of size 1
>>>>> ///
>>>>> The code on the error position is
>>>>> MatCreateSeqDense(PETSC_COMM_SELF, nns, ns, NULL, &Kns));
>>>>>
>>>>
>>>> 1) "Seq" means sequential, that is "not parallel".
>>>>
>>>> 2) This line should still be fine since PETSC_COMM_SELF is a serial
>>>> communicator
>>>>
>>>> 3) You should be checking the error code for each call, maybe using the
>>>> CHKERRQ() macro
>>>>
>>>> 4) Please always send the entire error message, not a snippet
>>>>
>>>> THanks
>>>>
>>>> Matt
>>>>
>>>>
>>>>> Could "MatCreateSeqDense" not be used in parallel computing?
>>>>>
>>>>> *Second*, the same error message is repeated as many times as the
>>>>> number of cores.
>>>>> if I use command -np 4, then the error message is repeated 4 times.
>>>>> Could you recommend some advice related to this?
>>>>>
>>>>> Best,
>>>>> Seung Lee Kwon
>>>>>
>>>>> --
>>>>> Seung Lee Kwon, Ph.D.Candidate
>>>>> Aerospace Structures and Materials Laboratory
>>>>> Department of Mechanical and Aerospace Engineering
>>>>> Seoul National University
>>>>> Building 300 Rm 503, Gwanak-ro 1, Gwanak-gu, Seoul, South Korea, 08826
>>>>> E-mail : ksl7912 at snu.ac.kr
>>>>> Office : +82-2-880-7389
>>>>> C. P : +82-10-4695-1062
>>>>>
>>>>
>>>>
>>>> --
>>>> 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/>
>>>>
>>>
>>>
>>> --
>>> Seung Lee Kwon, Ph.D.Candidate
>>> Aerospace Structures and Materials Laboratory
>>> Department of Mechanical and Aerospace Engineering
>>> Seoul National University
>>> Building 300 Rm 503, Gwanak-ro 1, Gwanak-gu, Seoul, South Korea, 08826
>>> E-mail : ksl7912 at snu.ac.kr
>>> Office : +82-2-880-7389
>>> C. P : +82-10-4695-1062
>>>
>>>
>>>
>>
>> --
>> Seung Lee Kwon, Ph.D.Candidate
>> Aerospace Structures and Materials Laboratory
>> Department of Mechanical and Aerospace Engineering
>> Seoul National University
>> Building 300 Rm 503, Gwanak-ro 1, Gwanak-gu, Seoul, South Korea, 08826
>> E-mail : ksl7912 at snu.ac.kr
>> Office : +82-2-880-7389
>> C. P : +82-10-4695-1062
>>
>>
>>
>
> --
> Seung Lee Kwon, Ph.D.Candidate
> Aerospace Structures and Materials Laboratory
> Department of Mechanical and Aerospace Engineering
> Seoul National University
> Building 300 Rm 503, Gwanak-ro 1, Gwanak-gu, Seoul, South Korea, 08826
> E-mail : ksl7912 at snu.ac.kr
> Office : +82-2-880-7389
> C. P : +82-10-4695-1062
>
--
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/>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20230505/fa53bd0c/attachment-0001.html>
More information about the petsc-users
mailing list