xref: /petsc/src/ksp/ksp/tests/ex16f.F90 (revision bcd4bb4a4158aa96f212e9537e87b40407faf83e)
1!
2#include <petsc/finclude/petscksp.h>
3program main
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  PetscCallA(PetscInitialize(ierr))
38
39!     Read in matrix and RHS
40  PetscCallA(PetscOptionsGetString(PETSC_NULL_OPTIONS, PETSC_NULL_CHARACTER, '-f', f, flg, ierr))
41  PetscCallA(PetscViewerBinaryOpen(PETSC_COMM_WORLD, f, FILE_MODE_READ, fd, ierr))
42
43  PetscCallA(MatCreate(PETSC_COMM_WORLD, A, ierr))
44  PetscCallA(MatSetType(A, MATSEQAIJ, ierr))
45  PetscCallA(MatLoad(A, fd, ierr))
46
47  PetscCallA(VecCreate(PETSC_COMM_WORLD, b, ierr))
48  PetscCallA(VecLoad(b, fd, ierr))
49  PetscCallA(PetscViewerDestroy(fd, ierr))
50
51! Set up solution
52  PetscCallA(VecDuplicate(b, x, ierr))
53  PetscCallA(VecDuplicate(b, u, ierr))
54
55! Solve system-1
56  PetscCallA(KSPCreate(PETSC_COMM_WORLD, ksp1, ierr))
57  PetscCallA(KSPSetOptionsPrefix(ksp1, 'a', ierr))
58  PetscCallA(KSPAppendOptionsPrefix(ksp1, '_', ierr))
59  PetscCallA(KSPSetOperators(ksp1, A, A, ierr))
60  PetscCallA(KSPSetFromOptions(ksp1, ierr))
61  PetscCallA(KSPSolve(ksp1, b, x, ierr))
62
63! Show result
64  PetscCallA(MatMult(A, x, u, ierr))
65  PetscCallA(VecAXPY(u, none, b, ierr))
66  PetscCallA(VecNorm(u, NORM_2, norm, ierr))
67  PetscCallA(KSPGetIterationNumber(ksp1, its, ierr))
68
69  write (6, 100) norm, its
70100 format('Residual norm ', e11.4, ' iterations ', i5)
71
72! Create system 2 by striping off some rows of the matrix
73  PetscCallA(ISCreateStride(PETSC_COMM_SELF, ifive, izero, ione, isrow, ierr))
74  PetscCallA(MatZeroRowsIS(A, isrow, five, PETSC_NULL_VEC, PETSC_NULL_VEC, ierr))
75
76! Solve system-2
77  PetscCallA(KSPCreate(PETSC_COMM_WORLD, ksp2, ierr))
78  PetscCallA(KSPSetOptionsPrefix(ksp2, 'b', ierr))
79  PetscCallA(KSPAppendOptionsPrefix(ksp2, '_', ierr))
80  PetscCallA(KSPSetOperators(ksp2, A, A, ierr))
81  PetscCallA(KSPSetFromOptions(ksp2, ierr))
82  PetscCallA(KSPSolve(ksp2, b, x, ierr))
83
84! Show result
85  PetscCallA(MatMult(A, x, u, ierr))
86  PetscCallA(VecAXPY(u, none, b, ierr))
87  PetscCallA(VecNorm(u, NORM_2, norm, ierr))
88  PetscCallA(KSPGetIterationNumber(ksp2, its, ierr))
89  write (6, 100) norm, its
90
91! Cleanup
92  PetscCallA(KSPDestroy(ksp1, ierr))
93  PetscCallA(KSPDestroy(ksp2, ierr))
94  PetscCallA(VecDestroy(b, ierr))
95  PetscCallA(VecDestroy(x, ierr))
96  PetscCallA(VecDestroy(u, ierr))
97  PetscCallA(MatDestroy(A, ierr))
98  PetscCallA(ISDestroy(isrow, ierr))
99
100  PetscCallA(PetscFinalize(ierr))
101end
102
103!/*TEST
104!
105!    test:
106!      args: -f ${DATAFILESPATH}/matrices/arco1 -options_left no
107!      requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
108!
109!TEST*/
110