Actual source code: ex12f.F
petsc-3.5.4 2015-05-23
1: !
2: program main
3: implicit none
5: #include <finclude/petscsys.h>
6: #include <finclude/petscvec.h>
7: #include <finclude/petscmat.h>
8: #include <finclude/petscpc.h>
9: #include <finclude/petscksp.h>
10: #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 example see $PETSC_DIR/share/petsc/datafiles/matrices
15: !
17: PetscErrorCode ierr
18: PetscInt its,m,n,mlocal,nlocal
19: PetscBool flg
20: PetscScalar norm,none
21: Vec x,b,u
22: Mat A
23: character*(128) f
24: PetscViewer fd
25: MatInfo info(MAT_INFO_SIZE)
26: KSP ksp
28: none = -1.0
29: call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
31: ! Read in matrix and RHS
32: call PetscOptionsGetString(PETSC_NULL_CHARACTER,'-f',f,flg,ierr)
33: call PetscViewerBinaryOpen(PETSC_COMM_WORLD,f,FILE_MODE_READ, &
34: & fd,ierr)
36: call MatCreate(PETSC_COMM_WORLD,A,ierr)
37: call MatSetType(A, MATSEQAIJ,ierr)
38: call MatLoad(A,fd,ierr)
40: ! Get information about matrix
41: call MatGetSize(A,m,n,ierr)
42: call MatGetLocalSize(A,mlocal,nlocal,ierr)
43: call MatGetInfo(A,MAT_GLOBAL_SUM,info,ierr)
44: write(*,100) m, &
45: & n, &
46: & mlocal,nlocal, &
47: & info(MAT_INFO_BLOCK_SIZE),info(MAT_INFO_NZ_ALLOCATED), &
48: & info(MAT_INFO_NZ_USED),info(MAT_INFO_NZ_UNNEEDED), &
49: & info(MAT_INFO_MEMORY),info(MAT_INFO_ASSEMBLIES), &
50: & info(MAT_INFO_MALLOCS)
52: 100 format(4(i4,1x),7(g7.1,1x))
53: call VecCreate(PETSC_COMM_WORLD,b,ierr)
54: call VecLoad(b,fd,ierr)
55: call PetscViewerDestroy(fd,ierr)
57: ! Set up solution
58: call VecDuplicate(b,x,ierr)
59: call VecDuplicate(b,u,ierr)
61: ! Solve system
62: call KSPCreate(PETSC_COMM_WORLD,ksp,ierr)
63: call KSPSetOperators(ksp,A,A,ierr)
64: call KSPSetFromOptions(ksp,ierr)
65: call KSPSolve(ksp,b,x,ierr)
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(ksp,its,ierr)
72: write(6,101) norm,its
73: 101 format('Residual norm ',e10.4,' iterations ',i5)
75: ! Cleanup
76: call KSPDestroy(ksp,ierr)
77: call VecDestroy(b,ierr)
78: call VecDestroy(x,ierr)
79: call VecDestroy(u,ierr)
80: call MatDestroy(A,ierr)
82: call PetscFinalize(ierr)
83: end