Actual source code: ex71f.F90
1: ! Contributed by leonardo.mutti01@universitadipavia.it
2: program main
3: #include <petsc/finclude/petscsys.h>
4: #include <petsc/finclude/petscmat.h>
5: #include <petsc/finclude/petscpc.h>
6: USE petscmat
7: USE petscpc
8: implicit none
10: Mat :: A
11: PetscInt :: M, M2, NSubx, dof, overlap, NSub,i1
12: PetscInt :: I, J
13: PetscMPIInt :: size
14: PetscErrorCode :: ierr
15: PetscScalar :: v
16: PC :: pc
17: IS :: subdomains_IS(20), inflated_IS(20)
18: PetscViewer :: singleton
20: PetscCallA(PetscInitialize(PETSC_NULL_CHARACTER, ierr))
21: PetscCallMPIA(MPI_Comm_size(PETSC_COMM_WORLD, size, ierr))
22: M = 16
23: M2 = M*M
24: i1 = 1
25: PetscCallA(MatCreateFromOptions(PETSC_COMM_WORLD, PETSC_NULL_CHARACTER,i1, PETSC_DECIDE, PETSC_DECIDE, M2, M2,A, ierr))
26: DO I=1,M2
27: DO J=1,M2
28: v = I*J
29: PetscCallA(MatSetValue(A, I-1, J-1, v, INSERT_VALUES, ierr))
30: END DO
31: END DO
33: PetscCallA(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY, ierr))
34: PetscCallA(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY, ierr))
35: PetscCallA(PCCreate(PETSC_COMM_WORLD, pc, ierr))
36: PetscCallA(PCSetOperators(pc, A, A, ierr))
37: PetscCallA(PCSetType(pc, PCGASM, ierr))
39: NSubx = 4
40: dof = 1
41: overlap = 0
43: PetscCallA(PCGASMCreateSubdomains2D(pc, M, M,NSubx, NSubx, dof, overlap, NSub, subdomains_IS, inflated_IS, ierr))
44: PetscCallA(PCGASMSetSubdomains(pc,NSub,subdomains_IS, inflated_IS, ierr))
45: PetscCallA(PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, singleton, ierr))
46: PetscCallA(PetscViewerASCIIPrintf(singleton, 'GASM index sets from this MPI process\n', ierr))
47: do i=1,Nsub
48: PetscCallA(ISView(subdomains_IS(i), singleton, ierr))
49: end do
50: PetscCallA(PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, singleton, ierr))
51: PetscCallA(PCGASMDestroySubdomains(NSub, subdomains_IS, inflated_IS, ierr))
53: if (size == 1) then
54: ! this routine only works on one rank
55: PetscCallA(PCASMCreateSubdomains2D(M, M, NSubx, NSubx, dof, overlap, NSub, subdomains_IS, inflated_IS, ierr))
56: do i=1,Nsub
57: PetscCallA(ISView(subdomains_IS(i), PETSC_VIEWER_STDOUT_SELF, ierr))
58: end do
59: PetscCallA(PCASMDestroySubdomains(NSub, subdomains_IS, inflated_IS, ierr))
60: endif
62: PetscCallA(MatDestroy(A, ierr))
63: PetscCallA(PCDestroy(pc, ierr))
64: PetscCallA(PetscFinalize(ierr))
65: end
67: !/*TEST
68: !
69: ! test:
70: ! suffix: 1
71: !
72: ! test:
73: ! suffix: 2
74: ! nsize: 2
75: !
76: !TEST*/