Actual source code: ex12f.F
1: !
2: program main
3: implicit none
5: #include include/finclude/petsc.h
6: #include include/finclude/petscvec.h
7: #include include/finclude/petscmat.h
8: #include include/finclude/petscpc.h
9: #include include/finclude/petscksp.h
10: #include include/finclude/petscviewer.h
11: !
12: ! This example is the Fortran version of ex6.c. The program reads a PETSc matrix
13: ! and vector from a file and solves a linear system. Input arguments are:
14: ! -f <input_file> : file to load. For a 5X5 example of the 5-pt. stencil
15: ! use the file petsc/src/mat/examples/matbinary.ex
16: !
18: PetscErrorCode ierr
19: PetscInt its
20: PetscTruth flg
21: PetscScalar norm,none
22: Vec x,b,u
23: Mat A
24: character*(128) f
25: PetscViewer fd
26: MatInfo info(MAT_INFO_SIZE)
27: KSP ksp
29: none = -1.0
30: call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
32: ! Read in matrix and RHS
33: call PetscOptionsGetString(PETSC_NULL_CHARACTER,'-f',f,flg,ierr)
34: print *,f
35: call PetscViewerBinaryOpen(PETSC_COMM_WORLD,f,PETSC_FILE_RDONLY, &
36: & fd,ierr)
38: call MatLoad(fd,MATSEQAIJ,A,ierr)
40: ! Get information about matrix
41: call MatGetInfo(A,MAT_GLOBAL_SUM,info,ierr)
42: write(*,100) info(MAT_INFO_ROWS_GLOBAL), &
43: & info(MAT_INFO_COLUMNS_GLOBAL), &
44: & info(MAT_INFO_ROWS_LOCAL),info(MAT_INFO_COLUMNS_LOCAL), &
45: & info(MAT_INFO_BLOCK_SIZE),info(MAT_INFO_NZ_ALLOCATED), &
46: & info(MAT_INFO_NZ_USED),info(MAT_INFO_NZ_UNNEEDED), &
47: & info(MAT_INFO_MEMORY),info(MAT_INFO_ASSEMBLIES), &
48: & info(MAT_INFO_MALLOCS)
50: 100 format(11(g7.1,1x))
51: call VecLoad(fd,PETSC_NULL_CHARACTER,b,ierr)
52: call PetscViewerDestroy(fd,ierr)
54: ! Set up solution
55: call VecDuplicate(b,x,ierr)
56: call VecDuplicate(b,u,ierr)
58: ! Solve system
59: call KSPCreate(PETSC_COMM_WORLD,ksp,ierr)
60: call KSPSetOperators(ksp,A,A,DIFFERENT_NONZERO_PATTERN,ierr)
61: call KSPSetFromOptions(ksp,ierr)
62: call KSPSolve(ksp,b,x,ierr)
64: ! Show result
65: call MatMult(A,x,u,ierr)
66: call VecAXPY(u,none,b,ierr)
67: call VecNorm(u,NORM_2,norm,ierr)
68: call KSPGetIterationNumber(ksp,its,ierr)
69: print*, 'Number of iterations = ',its
70: print*, 'Residual norm = ',norm
72: ! Cleanup
73: call KSPDestroy(ksp,ierr)
74: call VecDestroy(b,ierr)
75: call VecDestroy(x,ierr)
76: call VecDestroy(u,ierr)
77: call MatDestroy(A,ierr)
79: call PetscFinalize(ierr)
80: end