1c4762a1bSJed Brown! 2a5b23f4aSJose E. Roman! Description: Tests PCFieldSplitGetIS and PCFieldSplitSetIS from Fortran. 3c4762a1bSJed Brown! 4c4762a1bSJed Brown program main 5c4762a1bSJed Brown#include <petsc/finclude/petscksp.h> 6c4762a1bSJed Brown use petscksp 7c4762a1bSJed Brown implicit none 8c4762a1bSJed Brown 9c4762a1bSJed Brown Vec x,b,u 10c4762a1bSJed Brown Mat A 11c4762a1bSJed Brown KSP ksp 12c4762a1bSJed Brown PC pc 13*ccfd86f1SBarry Smith PetscReal norm 14c4762a1bSJed Brown PetscErrorCode ierr 15c4762a1bSJed Brown PetscInt i,n,col(3),its,i1,i2,i3 16c4762a1bSJed Brown PetscInt ione,izero 17c4762a1bSJed Brown PetscBool flg 18c4762a1bSJed Brown PetscMPIInt size 19c4762a1bSJed Brown PetscScalar none,one,value(3) 20c4762a1bSJed Brown IS isin,isout 21c4762a1bSJed Brown 22c4762a1bSJed Brown! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 23c4762a1bSJed Brown! Beginning of program 24c4762a1bSJed Brown! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 25c4762a1bSJed Brown 26d8606c27SBarry Smith PetscCallA(PetscInitialize(ierr)) 27d8606c27SBarry Smith PetscCallMPIA(MPI_Comm_size(PETSC_COMM_WORLD,size,ierr)) 28dcb3e689SBarry Smith PetscCheckA(size .eq. 1,PETSC_COMM_WORLD,PETSC_ERR_WRONG_MPI_SIZE,'This is a uniprocessor example only') 29c4762a1bSJed Brown none = -1.0 30c4762a1bSJed Brown one = 1.0 31c4762a1bSJed Brown n = 10 32c4762a1bSJed Brown i1 = 1 33c4762a1bSJed Brown i2 = 2 34c4762a1bSJed Brown i3 = 3 35d8606c27SBarry Smith PetscCallA(PetscOptionsGetInt(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-n',n,flg,ierr)) 36c4762a1bSJed Brown 37c4762a1bSJed Brown! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 38c4762a1bSJed Brown! Compute the matrix and right-hand-side vector that define 39c4762a1bSJed Brown! the linear system, Ax = b. 40c4762a1bSJed Brown! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 41c4762a1bSJed Brown 42c4762a1bSJed Brown! Create matrix. When using MatCreate(), the matrix format can 43c4762a1bSJed Brown! be specified at runtime. 44c4762a1bSJed Brown 45d8606c27SBarry Smith PetscCallA(MatCreate(PETSC_COMM_WORLD,A,ierr)) 46d8606c27SBarry Smith PetscCallA(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n,ierr)) 47d8606c27SBarry Smith PetscCallA(MatSetFromOptions(A,ierr)) 48d8606c27SBarry Smith PetscCallA(MatSetUp(A,ierr)) 49c4762a1bSJed Brown 50c4762a1bSJed Brown! Assemble matrix. 51c4762a1bSJed Brown! - Note that MatSetValues() uses 0-based row and column numbers 52c4762a1bSJed Brown! in Fortran as well as in C (as set here in the array "col"). 53c4762a1bSJed Brown 54c4762a1bSJed Brown value(1) = -1.0 55c4762a1bSJed Brown value(2) = 2.0 56c4762a1bSJed Brown value(3) = -1.0 57c4762a1bSJed Brown do 50 i=1,n-2 58c4762a1bSJed Brown col(1) = i-1 59c4762a1bSJed Brown col(2) = i 60c4762a1bSJed Brown col(3) = i+1 615d83a8b1SBarry Smith PetscCallA(MatSetValues(A,i1,[i],i3,col,value,INSERT_VALUES,ierr)) 62c4762a1bSJed Brown 50 continue 63c4762a1bSJed Brown i = n - 1 64c4762a1bSJed Brown col(1) = n - 2 65c4762a1bSJed Brown col(2) = n - 1 665d83a8b1SBarry Smith PetscCallA(MatSetValues(A,i1,[i],i2,col,value,INSERT_VALUES,ierr)) 67c4762a1bSJed Brown i = 0 68c4762a1bSJed Brown col(1) = 0 69c4762a1bSJed Brown col(2) = 1 70c4762a1bSJed Brown value(1) = 2.0 71c4762a1bSJed Brown value(2) = -1.0 725d83a8b1SBarry Smith PetscCallA(MatSetValues(A,i1,[i],i2,col,value,INSERT_VALUES,ierr)) 73d8606c27SBarry Smith PetscCallA(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY,ierr)) 74d8606c27SBarry Smith PetscCallA(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY,ierr)) 75c4762a1bSJed Brown 76c4762a1bSJed Brown! Create vectors. Note that we form 1 vector from scratch and 77c4762a1bSJed Brown! then duplicate as needed. 78c4762a1bSJed Brown 79d8606c27SBarry Smith PetscCallA(VecCreate(PETSC_COMM_WORLD,x,ierr)) 80d8606c27SBarry Smith PetscCallA(VecSetSizes(x,PETSC_DECIDE,n,ierr)) 81d8606c27SBarry Smith PetscCallA(VecSetFromOptions(x,ierr)) 82d8606c27SBarry Smith PetscCallA(VecDuplicate(x,b,ierr)) 83d8606c27SBarry Smith PetscCallA(VecDuplicate(x,u,ierr)) 84c4762a1bSJed Brown 85c4762a1bSJed Brown! Set exact solution; then compute right-hand-side vector. 86c4762a1bSJed Brown 87d8606c27SBarry Smith PetscCallA(VecSet(u,one,ierr)) 88d8606c27SBarry Smith PetscCallA(MatMult(A,u,b,ierr)) 89c4762a1bSJed Brown 90c4762a1bSJed Brown! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 91c4762a1bSJed Brown! Create the linear solver and set various options 92c4762a1bSJed Brown! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 93c4762a1bSJed Brown 94c4762a1bSJed Brown! Create linear solver context 95c4762a1bSJed Brown 96d8606c27SBarry Smith PetscCallA(KSPCreate(PETSC_COMM_WORLD,ksp,ierr)) 97c4762a1bSJed Brown 98c4762a1bSJed Brown! Set operators. Here the matrix that defines the linear system 997addb90fSBarry Smith! also serves as the matrix used to construct the preconditioner. 100c4762a1bSJed Brown 101d8606c27SBarry Smith PetscCallA(KSPSetOperators(ksp,A,A,ierr)) 102c4762a1bSJed Brown 103c4762a1bSJed Brown! Set linear solver defaults for this problem (optional). 104c4762a1bSJed Brown! - By extracting the KSP and PC contexts from the KSP context, 105d8606c27SBarry Smith! we can then directly call any KSP and PC routines 106c4762a1bSJed Brown! to set various options. 107c4762a1bSJed Brown! - The following four statements are optional; all of these 108c4762a1bSJed Brown! parameters could alternatively be specified at runtime via 109*ccfd86f1SBarry Smith! KSPSetFromOptions() 110c4762a1bSJed Brown 111d8606c27SBarry Smith PetscCallA(KSPGetPC(ksp,pc,ierr)) 112d8606c27SBarry Smith PetscCallA(PCSetType(pc,PCFIELDSPLIT,ierr)) 113c4762a1bSJed Brown izero = 0 114c4762a1bSJed Brown ione = 1 115d8606c27SBarry Smith PetscCallA(ISCreateStride(PETSC_COMM_SELF,n,izero,ione,isin,ierr)) 116dcb3e689SBarry Smith PetscCallA(PCFieldSplitSetIS(pc,'splitname',isin,ierr)) 117dcb3e689SBarry Smith PetscCallA(PCFieldSplitGetIS(pc,'splitname',isout,ierr)) 118dcb3e689SBarry Smith PetscCheckA(isin .eq. isout,PETSC_COMM_SELF,PETSC_ERR_PLIB,'PCFieldSplitGetIS() failed') 119c4762a1bSJed Brown 120c4762a1bSJed Brown! Set runtime options, e.g., 121c4762a1bSJed Brown! -ksp_type <type> -pc_type <type> -ksp_monitor -ksp_rtol <rtol> 122c4762a1bSJed Brown! These options will override those specified above as long as 123c4762a1bSJed Brown! KSPSetFromOptions() is called _after_ any other customization 124c4762a1bSJed Brown! routines. 125c4762a1bSJed Brown 126d8606c27SBarry Smith PetscCallA(KSPSetFromOptions(ksp,ierr)) 127c4762a1bSJed Brown 128c4762a1bSJed Brown! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 129c4762a1bSJed Brown! Solve the linear system 130c4762a1bSJed Brown! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 131c4762a1bSJed Brown 132d8606c27SBarry Smith PetscCallA(KSPSolve(ksp,b,x,ierr)) 133c4762a1bSJed Brown 134c4762a1bSJed Brown! View solver info; we could instead use the option -ksp_view 135c4762a1bSJed Brown 136d8606c27SBarry Smith PetscCallA(KSPView(ksp,PETSC_VIEWER_STDOUT_WORLD,ierr)) 137c4762a1bSJed Brown 138c4762a1bSJed Brown! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 139c4762a1bSJed Brown! Check solution and clean up 140c4762a1bSJed Brown! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 141c4762a1bSJed Brown 142c4762a1bSJed Brown! Check the error 143c4762a1bSJed Brown 144d8606c27SBarry Smith PetscCallA(VecAXPY(x,none,u,ierr)) 145d8606c27SBarry Smith PetscCallA(VecNorm(x,NORM_2,norm,ierr)) 146d8606c27SBarry Smith PetscCallA(KSPGetIterationNumber(ksp,its,ierr)) 147c4762a1bSJed Brown if (norm .gt. 1.e-12) then 148c4762a1bSJed Brown write(6,100) norm,its 149c4762a1bSJed Brown else 150c4762a1bSJed Brown write(6,200) its 151c4762a1bSJed Brown endif 152c4762a1bSJed Brown 100 format('Norm of error ',e11.4,', Iterations = ',i5) 153c4762a1bSJed Brown 200 format('Norm of error < 1.e-12, Iterations = ',i5) 154c4762a1bSJed Brown 155c4762a1bSJed Brown! Free work space. All PETSc objects should be destroyed when they 156c4762a1bSJed Brown! are no longer needed. 157c4762a1bSJed Brown 158d8606c27SBarry Smith PetscCallA(ISDestroy(isin,ierr)) 159d8606c27SBarry Smith PetscCallA(VecDestroy(x,ierr)) 160d8606c27SBarry Smith PetscCallA(VecDestroy(u,ierr)) 161d8606c27SBarry Smith PetscCallA(VecDestroy(b,ierr)) 162d8606c27SBarry Smith PetscCallA(MatDestroy(A,ierr)) 163d8606c27SBarry Smith PetscCallA(KSPDestroy(ksp,ierr)) 164d8606c27SBarry Smith PetscCallA(PetscFinalize(ierr)) 165c4762a1bSJed Brown 166c4762a1bSJed Brown end 167c4762a1bSJed Brown 168c4762a1bSJed Brown!/*TEST 169c4762a1bSJed Brown! 170c4762a1bSJed Brown! test: 171c4762a1bSJed Brown! args: -ksp_monitor 172c4762a1bSJed Brown! 173c4762a1bSJed Brown!TEST*/ 174