1! 2! Program to test recently added F90 features for Mat 3! 4 program main 5 6#include <petsc/finclude/petscmat.h> 7 use petscmat 8 implicit none 9 10 PetscErrorCode ierr 11 Mat A,B 12 Mat C,SC 13 MatNullSpace sp,sp1 14 PetscInt one,zero,rend 15 PetscScalar sone 16 Vec x,y 17 18 zero = 0 19 one = 1 20 sone = 1 21 call PetscInitialize(PETSC_NULL_CHARACTER,ierr) 22 if (ierr .ne. 0) then 23 print*, 'Unable to begin PETSc program' 24 endif 25 26 call MatCreate(PETSC_COMM_WORLD,A,ierr) 27 call MatCreate(PETSC_COMM_WORLD,B,ierr) 28 29 call MatGetNullSpace(A,sp,ierr) 30 if (sp .ne. PETSC_NULL_MATNULLSPACE) then; SETERRA(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Matrix null space should not exist"); endif 31 32 call MatSetNullSpace(A,PETSC_NULL_MATNULLSPACE,ierr) 33 call MatGetNullSpace(A,sp,ierr) 34 if (sp .ne. PETSC_NULL_MATNULLSPACE) then; SETERRA(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Matrix null space should not exist"); endif 35 36 call MatNullSpaceCreate(PETSC_COMM_WORLD,PETSC_TRUE,zero,PETSC_NULL_VEC,sp,ierr) 37 call MatSetNullSpace(A,sp,ierr) 38 call MatGetNullSpace(A,sp1,ierr) 39 if (sp1 .eq. PETSC_NULL_MATNULLSPACE) then; SETERRA(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Matrix null space should not be null"); endif 40 call MatNullSpaceDestroy(sp,ierr) 41 42 call MatCreateSeqDense(PETSC_COMM_WORLD,one,one,PETSC_NULL_SCALAR,C,ierr) 43 call MatSetValues(C,one,zero,one,zero,sone,INSERT_VALUES,ierr) 44 call MatAssemblyBegin(C,MAT_FINAL_ASSEMBLY,ierr) 45 call MatAssemblyEnd(C,MAT_FINAL_ASSEMBLY,ierr) 46 call MatCreateSchurComplement(C,C,C,C,PETSC_NULL_MAT,SC,ierr) 47 call MatGetOwnershipRange(SC,PETSC_NULL_INTEGER,rend,ierr) 48 call VecCreateSeq(PETSC_COMM_SELF,one,x,ierr) 49 call VecDuplicate(x,y,ierr) 50 call VecSetValues(x,one,zero,sone,INSERT_VALUES,ierr) 51 call VecAssemblyBegin(x,ierr) 52 call VecAssemblyEnd(x,ierr) 53 call MatMult(SC,x,y,ierr) 54 call VecView(y,PETSC_VIEWER_STDOUT_SELF,ierr) 55 call VecSetRandom(x,PETSC_NULL_RANDOM,ierr) 56 call VecView(x,PETSC_VIEWER_STDOUT_SELF,ierr) 57 58 call MatDestroy(SC,ierr) 59 call MatDestroy(C,ierr) 60 call VecDestroy(x,ierr) 61 call VecDestroy(y,ierr) 62 call MatDestroy(A,ierr) 63 call MatDestroy(B,ierr) 64 call PetscFinalize(ierr) 65 end 66 67!/*TEST 68! 69! test: 70! requires: !complex 71! 72!TEST*/ 73