Actual source code: ex15f.F90
1: program main
3: #include <petsc/finclude/petscvec.h>
4: #include <petsc/finclude/petscmat.h>
6: use petscvec
7: use petscmat
9: implicit none
11: Mat :: A
12: MatPartitioning :: part
13: IS :: is
14: PetscInt :: r,myStart,myEnd
15: PetscInt :: N = 10
16: PetscErrorCode :: ierr
17: PetscScalar,pointer,dimension(:) :: vals
18: PetscInt,pointer,dimension(:) :: cols
19: PetscBool :: flg
20: PetscInt,parameter :: one = 1, two = 2, three = 3
22: PetscCallA(PetscInitialize(ierr))
24: PetscCallA(PetscOptionsGetInt(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-N',N,flg,ierr))
25: PetscCallA(MatCreate(PETSC_COMM_WORLD, A,ierr))
26: PetscCallA(MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, N, N,ierr))
27: PetscCallA(MatSetFromOptions(A,ierr))
28: PetscCallA(MatSeqAIJSetPreallocation(A, three, PETSC_NULL_INTEGER,ierr))
29: PetscCallA(MatMPIAIJSetPreallocation(A, three, PETSC_NULL_INTEGER, two, PETSC_NULL_INTEGER,ierr))
31: !/* Create a linear mesh */
32: PetscCallA(MatGetOwnershipRange(A, myStart, myEnd,ierr))
34: do r=myStart,myEnd-1
35: if (r == 0) then
36: allocate(vals(2))
37: vals = 1.0
38: allocate(cols(2),source=[r,r+1])
39: PetscCallA(MatSetValues(A, one, r, two, cols, vals, INSERT_VALUES,ierr))
40: deallocate(cols)
41: deallocate(vals)
42: else if (r == N-1) then
43: allocate(vals(2))
44: vals = 1.0
45: allocate(cols(2),source=[r-1,r])
46: PetscCallA(MatSetValues(A, one, r, two, cols, vals, INSERT_VALUES,ierr))
47: deallocate(cols)
48: deallocate(vals)
49: else
50: allocate(vals(3))
51: vals = 1.0
52: allocate(cols(3),source=[r-1,r,r+1])
53: PetscCallA(MatSetValues(A, one, r, three, cols, vals, INSERT_VALUES,ierr))
54: deallocate(cols)
55: deallocate(vals)
56: end if
57: end do
58: PetscCallA(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY,ierr))
59: PetscCallA(MatAssemblyend(A, MAT_FINAL_ASSEMBLY,ierr))
60: PetscCallA(MatPartitioningCreate(PETSC_COMM_WORLD, part,ierr))
61: PetscCallA(MatPartitioningSetAdjacency(part, A,ierr))
62: PetscCallA(MatPartitioningSetFromOptions(part,ierr))
63: PetscCallA(MatPartitioningApply(part, is,ierr))
64: PetscCallA(ISView(is, PETSC_VIEWER_STDOUT_WORLD,ierr))
65: PetscCallA(ISDestroy(is,ierr))
66: PetscCallA(MatPartitioningDestroy(part,ierr))
67: PetscCallA(MatDestroy(A,ierr))
68: PetscCallA(PetscFinalize(ierr))
70: end program
72: !/*TEST
73: !
74: ! test:
75: ! nsize: 3
76: ! requires: parmetis
77: ! args: -mat_partitioning_type parmetis
78: ! output_file: output/ex15_1.out
79: !
80: ! test:
81: ! suffix: 2
82: ! nsize: 3
83: ! requires: ptscotch
84: ! args: -mat_partitioning_type ptscotch
85: ! output_file: output/ex15_2.out
86: !
87: ! test:
88: ! suffix: 3
89: ! nsize: 4
90: ! requires: party
91: ! args: -mat_partitioning_type party
92: ! output_file: output/ex15_3.out
93: !
94: ! test:
95: ! suffix: 4
96: ! nsize: 3
97: ! requires: chaco
98: ! args: -mat_partitioning_type chaco
99: ! output_file: output/ex15_4.out
100: !
101: ! test:
102: ! suffix: 5
103: ! nsize: 3
104: ! requires: parmetis
105: ! args: -mat_partitioning_type hierarch -mat_partitioning_hierarchical_nfineparts 3 -mat_partitioning_nparts 10 -N 100
106: ! output_file: output/ex15_5.out
107: !
108: !TEST*/