Actual source code: ex12f.F90
1: !
2: program main
3: #include <petsc/finclude/petscksp.h>
4: use petscksp
5: implicit none
7: !
8: ! This example is the Fortran version of ex6.c. The program reads a PETSc matrix
9: ! and vector from a file and solves a linear system. Input arguments are:
10: ! -f <input_file> : file to load. For example see $PETSC_DIR/share/petsc/datafiles/matrices
11: !
13: PetscErrorCode ierr
14: PetscInt its,m,n,mlocal,nlocal
15: PetscBool flg
16: PetscScalar none
17: PetscReal norm
18: Vec x,b,u
19: Mat A
20: character*(128) f
21: PetscViewer fd
22: MatInfo info(MAT_INFO_SIZE)
23: KSP ksp
25: none = -1.0
26: PetscCallA(PetscInitialize(ierr))
28: ! Read in matrix and RHS
29: PetscCallA(PetscOptionsGetString(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-f',f,flg,ierr))
30: PetscCallA(PetscViewerBinaryOpen(PETSC_COMM_WORLD,f,FILE_MODE_READ,fd,ierr))
32: PetscCallA(MatCreate(PETSC_COMM_WORLD,A,ierr))
33: PetscCallA(MatSetType(A, MATSEQAIJ,ierr))
34: PetscCallA(MatLoad(A,fd,ierr))
36: ! Get information about matrix
37: PetscCallA(MatGetSize(A,m,n,ierr))
38: PetscCallA(MatGetLocalSize(A,mlocal,nlocal,ierr))
39: PetscCallA(MatGetInfo(A,MAT_GLOBAL_SUM,info,ierr))
40: write(*,100) m, &
41: & n, &
42: & mlocal,nlocal, &
43: & info(MAT_INFO_BLOCK_SIZE),info(MAT_INFO_NZ_ALLOCATED), &
44: & info(MAT_INFO_NZ_USED),info(MAT_INFO_NZ_UNNEEDED), &
45: & info(MAT_INFO_MEMORY),info(MAT_INFO_ASSEMBLIES), &
46: & info(MAT_INFO_MALLOCS)
48: 100 format(4(i4,1x),7(1pe9.2,1x))
49: PetscCallA(VecCreate(PETSC_COMM_WORLD,b,ierr))
50: PetscCallA(VecLoad(b,fd,ierr))
51: PetscCallA(PetscViewerDestroy(fd,ierr))
53: ! Set up solution
54: PetscCallA(VecDuplicate(b,x,ierr))
55: PetscCallA(VecDuplicate(b,u,ierr))
57: ! Solve system
58: PetscCallA(KSPCreate(PETSC_COMM_WORLD,ksp,ierr))
59: PetscCallA(KSPSetOperators(ksp,A,A,ierr))
60: PetscCallA(KSPSetFromOptions(ksp,ierr))
61: PetscCallA(KSPSolve(ksp,b,x,ierr))
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(ksp,its,ierr))
68: write(6,101) norm,its
69: 101 format('Residual norm ',1pe9.2,' iterations ',i5)
71: ! Cleanup
72: PetscCallA(KSPDestroy(ksp,ierr))
73: PetscCallA(VecDestroy(b,ierr))
74: PetscCallA(VecDestroy(x,ierr))
75: PetscCallA(VecDestroy(u,ierr))
76: PetscCallA(MatDestroy(A,ierr))
78: PetscCallA(PetscFinalize(ierr))
79: end
81: !/*TEST
82: !
83: ! test:
84: ! args: -f ${DATAFILESPATH}/matrices/arco1 -options_left no
85: ! requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
86: !
87: !TEST*/