[petsc-users] Solving complex linear sparse matrix in parallel + external library

[Sal Am]

Hi,

I have a few questions:

1. The following issue/misunderstanding:
My code reads in two files one PETSc vector and one PETSc matrix (b and A
from Ax=b, size ~65000x65000).
and then calls KSP solver to solve it by running the following in the
terminal:

mpiexec - n 2 ./SolveSys -ksp_type preonly -pc_type lu
-pc_factor_mat_solver mumps

Now mumps is supposed to work in parallel and complex, but the code is not
solved in parallel it seems. It just prints the result twice. Adding
-log_view gives me

"./SolveSys on a linux-opt named F8434 with 1 processor..." printed twice.

2. Using iterative solvers, I am having difficulty getting convergence. I
found that there is a way to set the maximum number of iterations, but is
there a minimum I can increase?

3. The residual is not computed when using direct external solvers. What is
proper PETSc way of doing this?
-----------------------------------------------------The
code---------------------------------------
#include <petscksp.h>
#include <iostream>
int main(int argc,char **args)
{
  Vec            x,b;      /* approx solution, RHS */
  Mat            A;            /* linear system matrix */
  KSP            ksp;         /* linear solver context */
  PetscReal      norm;         /* norm of solution error */
  PC             pc;
  PetscMPIInt    rank, size;
  PetscViewer    viewer;
  PetscInt       its, i;
  PetscErrorCode ierr;
  PetscScalar    *xa;
  PetscBool      flg = PETSC_FALSE;

  ierr = PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr;
  MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
  MPI_Comm_size(PETSC_COMM_WORLD,&size);

  #if !defined(PETSC_USE_COMPLEX)
    SETERRQ(PETSC_COMM_WORLD,1,"This example requires complex numbers");
  #endif
  /*
   - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
         Compute the matrix and right-hand-side vector that define
         the linear system, Ax = b.
     - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  */
  ierr = PetscPrintf(PETSC_COMM_WORLD,"reading vector in binary from
Vector_b.dat ...\n");CHKERRQ(ierr);
  ierr =
PetscViewerBinaryOpen(PETSC_COMM_WORLD,"../../python/petscpy/Vector_b.dat",FILE_MODE_READ,&viewer);CHKERRQ(ierr);
  ierr = VecCreate(PETSC_COMM_WORLD, &b);CHKERRQ(ierr);
  ierr = VecLoad(b,viewer); CHKERRQ(ierr);

  ierr = PetscPrintf(PETSC_COMM_WORLD,"reading matrix in binary from
Matrix_A.dat ...\n");CHKERRQ(ierr);
  ierr =
PetscViewerBinaryOpen(PETSC_COMM_WORLD,"../../python/petscpy/Matrix_A.dat",FILE_MODE_READ,&viewer);CHKERRQ(ierr);
  ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr);
  ierr = MatLoad(A,viewer);CHKERRQ(ierr);
  ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr);

  ierr = VecDuplicate(b,&x);CHKERRQ(ierr);

  /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
                Create the linear solver and set various options
     - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  */
  PetscPrintf(PETSC_COMM_WORLD, "Creating KSP\n");
  ierr = KSPCreate(PETSC_COMM_WORLD,&ksp);CHKERRQ(ierr);


  PetscPrintf(PETSC_COMM_WORLD, "KSP Operators\n");
  ierr = KSPSetOperators(ksp,A,A);CHKERRQ(ierr);

  ierr = KSPGetPC(ksp, &pc);CHKERRQ(ierr);
  ierr = PCSetType(pc, PCLU);CHKERRQ(ierr);
  ierr = KSPSetFromOptions(ksp);CHKERRQ(ierr);

  /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
                      Solve the linear system
     - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
  ierr = KSPSetUp(ksp);CHKERRQ(ierr);
  ierr = KSPSetUpOnBlocks(ksp);CHKERRQ(ierr);
  ierr = KSPSolve(ksp,b,x);CHKERRQ(ierr);
  PetscPrintf(PETSC_COMM_WORLD, "Solved");

  /*
     Free work space.  All PETSc objects should be destroyed when they
     are no longer needed.
  */
  ierr = KSPDestroy(&ksp);CHKERRQ(ierr);
  ierr = VecDestroy(&x);CHKERRQ(ierr);
  ierr = VecDestroy(&b);CHKERRQ(ierr);
  ierr = MatDestroy(&A);CHKERRQ(ierr);
  ierr = PetscFinalize();
  return ierr;
}


Kind regards,
Sal
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20181116/abda2b19/attachment.html>