1c080761bSJose E. Roman! Test MatCreateMPIAdj() with NULL argument 'values' 2c080761bSJose E. Roman 3c080761bSJose E. Roman#include <petsc/finclude/petscmat.h> 4*c5e229c2SMartin Diehlprogram main 5c080761bSJose E. Roman use petscmat 6c080761bSJose E. Roman implicit none 7c080761bSJose E. Roman 8c080761bSJose E. Roman Mat :: mesh, dual 9c080761bSJose E. Roman MatPartitioning :: part 10c080761bSJose E. Roman IS :: is 11c080761bSJose E. Roman PetscInt, parameter :: Nvertices = 6, ncells = 2, two = 2 12c080761bSJose E. Roman PetscInt :: ii(3), jj(6) 13c080761bSJose E. Roman PetscMPIInt :: sz, rnk 14c080761bSJose E. Roman PetscErrorCode :: ierr 15c080761bSJose E. Roman 16c080761bSJose E. Roman PetscCallA(PetscInitialize(ierr)) 17c080761bSJose E. Roman 18c080761bSJose E. Roman PetscCallMPIA(MPI_Comm_size(PETSC_COMM_WORLD, sz, ierr)) 194820e4eaSBarry Smith PetscCheckA(sz == 2, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, 'This example is for exactly two processes') 20c080761bSJose E. Roman PetscCallMPIA(MPI_Comm_rank(PETSC_COMM_WORLD, rnk, ierr)) 21c080761bSJose E. Roman ii(1) = 0 22c080761bSJose E. Roman ii(2) = 3 23c080761bSJose E. Roman ii(3) = 6 244820e4eaSBarry Smith if (rnk == 0) then 25c080761bSJose E. Roman jj(1) = 0 26c080761bSJose E. Roman jj(2) = 1 27c080761bSJose E. Roman jj(3) = 2 28c080761bSJose E. Roman jj(4) = 1 29c080761bSJose E. Roman jj(5) = 2 30c080761bSJose E. Roman jj(6) = 3 31c080761bSJose E. Roman else 32c080761bSJose E. Roman jj(1) = 1 33c080761bSJose E. Roman jj(2) = 4 34c080761bSJose E. Roman jj(3) = 5 35c080761bSJose E. Roman jj(4) = 1 36c080761bSJose E. Roman jj(5) = 3 37c080761bSJose E. Roman jj(6) = 5 38c080761bSJose E. Roman end if 39c080761bSJose E. Roman 40c080761bSJose E. Roman PetscCallA(MatCreateMPIAdj(PETSC_COMM_WORLD, ncells, Nvertices, ii, jj, PETSC_NULL_INTEGER_ARRAY, mesh, ierr)) 41c080761bSJose E. Roman PetscCallA(MatMeshToCellGraph(mesh, two, dual, ierr)) 42c080761bSJose E. Roman PetscCallA(MatView(dual, PETSC_VIEWER_STDOUT_WORLD, ierr)) 43c080761bSJose E. Roman 44c080761bSJose E. Roman PetscCallA(MatPartitioningCreate(PETSC_COMM_WORLD, part, ierr)) 45c080761bSJose E. Roman PetscCallA(MatPartitioningSetAdjacency(part, dual, ierr)) 46c080761bSJose E. Roman PetscCallA(MatPartitioningSetFromOptions(part, ierr)) 47c080761bSJose E. Roman PetscCallA(MatPartitioningApply(part, is, ierr)) 48c080761bSJose E. Roman PetscCallA(ISView(is, PETSC_VIEWER_STDOUT_WORLD, ierr)) 49c080761bSJose E. Roman PetscCallA(ISDestroy(is, ierr)) 50c080761bSJose E. Roman PetscCallA(MatPartitioningDestroy(part, ierr)) 51c080761bSJose E. Roman 52c080761bSJose E. Roman PetscCallA(MatDestroy(mesh, ierr)) 53c080761bSJose E. Roman PetscCallA(MatDestroy(dual, ierr)) 54c080761bSJose E. Roman PetscCallA(PetscFinalize(ierr)) 55c080761bSJose E. Roman 56c080761bSJose E. Romanend program 57c080761bSJose E. Roman 58c080761bSJose E. Roman!/*TEST 59c080761bSJose E. Roman! 60c080761bSJose E. Roman! build: 61c080761bSJose E. Roman! requires: parmetis 62c080761bSJose E. Roman! 63c080761bSJose E. Roman! test: 64c080761bSJose E. Roman! nsize: 2 65c080761bSJose E. Roman! 66c080761bSJose E. Roman!TEST*/ 67