xref: /petsc/src/ksp/ksp/tests/ex16f.F90 (revision d9acb416d05abeed0a33bde3a81aeb2ea0364f6a)
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      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
70  100 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))
101      end
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