1! 2 program main 3#include <petsc/finclude/petscksp.h> 4 use petscksp 5 implicit none 6 7! 8! This example is a modified Fortran version of ex6.c. It tests the use of 9! options prefixes in PETSc. Two linear problems are solved in this program. 10! The first problem is read from a file. The second problem is constructed 11! from the first, by eliminating some of the entries of the linear matrix 'A'. 12 13! Each solve is distinguished by a unique prefix - 'a' for the first, 'b' 14! for the second. With the prefix the user can distinguish between the various 15! options (command line, from .petscrc file, etc.) for each of the solvers. 16! Input arguments are: 17! -f <input_file> : file to load. For a 5X5 example of the 5-pt. stencil 18! use the file petsc/src/mat/examples/mat.ex.binary 19 20 PetscErrorCode ierr 21 PetscInt its,ione,ifive,izero 22 PetscBool flg 23 PetscScalar none,five 24 PetscReal norm 25 Vec x,b,u 26 Mat A 27 KSP ksp1,ksp2 28 character*(PETSC_MAX_PATH_LEN) f 29 PetscViewer fd 30 IS isrow 31 none = -1.0 32 five = 5.0 33 ifive = 5 34 ione = 1 35 izero = 0 36 37 call PetscInitialize(PETSC_NULL_CHARACTER,ierr) 38 if (ierr .ne. 0) then 39 print*,'Unable to initialize PETSc' 40 stop 41 endif 42 43! Read in matrix and RHS 44 call PetscOptionsGetString(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-f',f,flg,ierr);CHKERRA(ierr) 45 call PetscViewerBinaryOpen(PETSC_COMM_WORLD,f,FILE_MODE_READ,fd,ierr);CHKERRA(ierr) 46 47 call MatCreate(PETSC_COMM_WORLD,A,ierr) 48 call MatSetType(A, MATSEQAIJ,ierr) 49 call MatLoad(A,fd,ierr) 50 51 call VecCreate(PETSC_COMM_WORLD,b,ierr) 52 call VecLoad(b,fd,ierr) 53 call PetscViewerDestroy(fd,ierr) 54 55! Set up solution 56 call VecDuplicate(b,x,ierr) 57 call VecDuplicate(b,u,ierr) 58 59! Solve system-1 60 call KSPCreate(PETSC_COMM_WORLD,ksp1,ierr) 61 call KSPSetOptionsPrefix(ksp1,'a',ierr) 62 call KSPAppendOptionsPrefix(ksp1,'_',ierr) 63 call KSPSetOperators(ksp1,A,A,ierr) 64 call KSPSetFromOptions(ksp1,ierr) 65 call KSPSolve(ksp1,b,x,ierr) 66 67! Show result 68 call MatMult(A,x,u,ierr) 69 call VecAXPY(u,none,b,ierr) 70 call VecNorm(u,NORM_2,norm,ierr) 71 call KSPGetIterationNumber(ksp1,its,ierr) 72 73 write(6,100) norm,its 74 100 format('Residual norm ',e11.4,' iterations ',i5) 75 76! Create system 2 by striping off some rows of the matrix 77 call ISCreateStride(PETSC_COMM_SELF,ifive,izero,ione,isrow,ierr) 78 call MatZeroRowsIS(A,isrow,five,PETSC_NULL_VEC, & 79 & PETSC_NULL_VEC,ierr) 80 81! Solve system-2 82 call KSPCreate(PETSC_COMM_WORLD,ksp2,ierr) 83 call KSPSetOptionsPrefix(ksp2,'b',ierr) 84 call KSPAppendOptionsPrefix(ksp2,'_',ierr) 85 call KSPSetOperators(ksp2,A,A,ierr) 86 call KSPSetFromOptions(ksp2,ierr) 87 call KSPSolve(ksp2,b,x,ierr) 88 89! Show result 90 call MatMult(A,x,u,ierr) 91 call VecAXPY(u,none,b,ierr) 92 call VecNorm(u,NORM_2,norm,ierr) 93 call KSPGetIterationNumber(ksp2,its,ierr) 94 write(6,100) norm,its 95 96! Cleanup 97 call KSPDestroy(ksp1,ierr) 98 call KSPDestroy(ksp2,ierr) 99 call VecDestroy(b,ierr) 100 call VecDestroy(x,ierr) 101 call VecDestroy(u,ierr) 102 call MatDestroy(A,ierr) 103 call ISDestroy(isrow,ierr) 104 105 call PetscFinalize(ierr) 106 end 107 108!/*TEST 109! 110! test: 111! args: -f ${DATAFILESPATH}/matrices/arco1 -options_left no 112! requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) 113! 114!TEST*/ 115