Actual source code: ex76.c

petsc-3.12.2 2019-11-22
Report Typos and Errors
  1:  #include <petscksp.h>

  3: static char help[] = "Solves a linear system using PCHPDDM.\n\n";

  5: int main(int argc,char **args)
  6: {
  7:   Vec            x,b;        /* computed solution and RHS */
  8:   Mat            A,aux;      /* linear system matrix */
  9:   KSP            ksp;        /* linear solver context */
 10:   PC             pc;
 11:   IS             is,sizes;
 12:   const PetscInt *idx;
 13:   PetscMPIInt    rank,size;
 14:   char           dir[PETSC_MAX_PATH_LEN],name[PETSC_MAX_PATH_LEN];
 15:   PetscViewer    viewer;

 18:   PetscInitialize(&argc,&args,NULL,help);if (ierr) return ierr;
 19:   MPI_Comm_size(PETSC_COMM_WORLD,&size);
 20:   if (size != 4) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_USER,"This example requires 4 processes");
 21:   MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
 22:   MatCreate(PETSC_COMM_WORLD,&A);
 23:   MatCreate(PETSC_COMM_SELF,&aux);
 24:   ISCreate(PETSC_COMM_SELF,&is);
 25:   PetscStrcpy(dir,".");
 26:   PetscOptionsGetString(NULL,NULL,"-load_dir",dir,sizeof(dir),NULL);
 27:   /* loading matrices */
 28:   PetscSNPrintf(name,sizeof(name),"%s/sizes_%d_%d.dat",dir,rank,size);
 29:   PetscViewerBinaryOpen(PETSC_COMM_SELF,name,FILE_MODE_READ,&viewer);
 30:   ISCreate(PETSC_COMM_SELF,&sizes);
 31:   ISLoad(sizes,viewer);
 32:   ISGetIndices(sizes,&idx);
 33:   MatSetSizes(A,idx[0],idx[1],idx[2],idx[3]);
 34:   ISRestoreIndices(sizes,&idx);
 35:   ISDestroy(&sizes);
 36:   PetscViewerDestroy(&viewer);
 37:   MatSetUp(A);
 38:   PetscSNPrintf(name,sizeof(name),"%s/A.dat",dir);
 39:   PetscViewerBinaryOpen(PETSC_COMM_WORLD,name,FILE_MODE_READ,&viewer);
 40:   MatLoad(A,viewer);
 41:   PetscViewerDestroy(&viewer);
 42:   PetscSNPrintf(name,sizeof(name),"%s/is_%d_%d.dat",dir,rank,size);
 43:   PetscViewerBinaryOpen(PETSC_COMM_SELF,name,FILE_MODE_READ,&viewer);
 44:   ISLoad(is,viewer);
 45:   PetscViewerDestroy(&viewer);
 46:   PetscSNPrintf(name,sizeof(name),"%s/Neumann_%d_%d.dat",dir,rank,size);
 47:   PetscViewerBinaryOpen(PETSC_COMM_SELF,name,FILE_MODE_READ,&viewer);
 48:   MatLoad(aux,viewer);
 49:   PetscViewerDestroy(&viewer);
 50:   /* ready for testing */
 51:   KSPCreate(PETSC_COMM_WORLD,&ksp);
 52:   KSPSetOperators(ksp,A,A);
 53:   KSPSetFromOptions(ksp);
 54:   KSPGetPC(ksp,&pc);
 55:   PetscObjectCompose((PetscObject)pc,"_PCHPDDM_Neumann_IS",(PetscObject)is);
 56:   PetscObjectCompose((PetscObject)pc,"_PCHPDDM_Neumann_Mat",(PetscObject)aux);
 57:   ISDestroy(&is);
 58:   MatDestroy(&aux);
 59:   MatCreateVecs(A,&x,&b);
 60:   VecSet(b,1.0);
 61:   KSPSolve(ksp,b,x);
 62:   VecDestroy(&x);
 63:   VecDestroy(&b);
 64:   KSPDestroy(&ksp);
 65:   MatDestroy(&A);
 66:   PetscFinalize();
 67:   return ierr;
 68: }

 70: /*TEST

 72:    test:
 73:       requires: hpddm datafilespath double !complex !define(PETSC_USE_64BIT_INDICES)
 74:       nsize: 4
 75:       args: -ksp_rtol 1e-3 -ksp_converged_reason -pc_type {{bjacobi hpddm}shared output} -pc_hpddm_coarse_sub_pc_type lu -sub_pc_type lu -options_left no -load_dir ${DATAFILESPATH}/matrices/hpddm/GENEO

 77:    test:
 78:       requires: hpddm datafilespath double !complex !define(PETSC_USE_64BIT_INDICES)
 79:       suffix: geneo
 80:       nsize: 4
 81:       args: -ksp_converged_reason -pc_type hpddm -pc_hpddm_levels_1_sub_pc_type lu -pc_hpddm_levels_1_eps_nev {{5 10 15}separate output} -pc_hpddm_coarse_p {{1 2}shared output} -pc_hpddm_coarse_pc_type redundant -load_dir ${DATAFILESPATH}/matrices/hpddm/GENEO

 83:    test:
 84:       requires: hpddm datafilespath double !complex !define(PETSC_USE_64BIT_INDICES)
 85:       suffix: fgmres_geneo_20_p_2
 86:       nsize: 4
 87:       args: -ksp_converged_reason -pc_type hpddm -pc_hpddm_levels_1_sub_pc_type lu -pc_hpddm_levels_1_eps_nev 20 -pc_hpddm_coarse_p 2 -pc_hpddm_coarse_pc_type redundant -ksp_type fgmres -pc_hpddm_coarse_mat_type {{baij sbaij}shared output} -load_dir ${DATAFILESPATH}/matrices/hpddm/GENEO

 89:    test:
 90:       requires: hpddm datafilespath double !complex !define(PETSC_USE_64BIT_INDICES)
 91:       suffix: fgmres_geneo_20_p_2_geneo
 92:       nsize: 4
 93:       args: -ksp_converged_reason -pc_type hpddm -pc_hpddm_levels_1_sub_pc_type lu -pc_hpddm_levels_1_eps_nev 20 -pc_hpddm_levels_2_p 2 -pc_hpddm_levels_2_mat_type {{baij sbaij}shared output} -pc_hpddm_levels_2_eps_nev {{5 20}separate output} -pc_hpddm_levels_2_sub_pc_type cholesky -pc_hpddm_levels_2_ksp_type gmres -pc_hpddm_levels_2_ksp_converged_reason -ksp_type fgmres -pc_hpddm_coarse_mat_type {{baij sbaij}shared output} -load_dir ${DATAFILESPATH}/matrices/hpddm/GENEO

 95: TEST*/