1: static char help[] = "Tests that MatSetValues_MPISELL() invalidates the device copy when it inserts a new nonzero.\n\n";

  3: /*
  4:   The matrix is assembled and multiplied once so that both local submatrices hold current device
  5:   data, then new nonzeros are inserted into the diagonal block and into the off-diagonal block. If
  6:   MatSetValues_MPISELL() fails to mark the host copy authoritative, MatSeqSELLCUDACopyToGPU() skips
  7:   the transfer and the second MatMult() returns the pre-insertion result.

  9:   The off-diagonal insertion deliberately reuses a global column that another local row already
 10:   references. An unseen column would trigger MatDisAssemble_MPISELL(), which installs a fresh
 11:   submatrix whose offload mask starts out unallocated, and the stale-data path would not be taken.
 12: */

 14: #include <petscmat.h>

 16: int main(int argc, char **args)
 17: {
 18:   Mat         A, B;
 19:   PetscInt    i, rstart, rend, col;
 20:   PetscMPIInt rank, size;
 21:   PetscScalar value = 1.0;
 22:   PetscBool   flg;

 24:   PetscFunctionBeginUser;
 25:   PetscCall(PetscInitialize(&argc, &args, NULL, help));
 26:   PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
 27:   PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
 28:   PetscCheck(size == 2, PETSC_COMM_WORLD, PETSC_ERR_USER, "This test requires 2 processes");

 30:   /* A is the matrix under test, B an MATMPIAIJ reference that receives the same values */
 31:   PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
 32:   PetscCall(MatSetSizes(A, 4, 4, PETSC_DETERMINE, PETSC_DETERMINE));
 33:   PetscCall(MatSetType(A, MATSELL));
 34:   PetscCall(MatSetFromOptions(A));
 35:   PetscCall(MatMPISELLSetPreallocation(A, 4, NULL, 4, NULL));

 37:   PetscCall(MatCreate(PETSC_COMM_WORLD, &B));
 38:   PetscCall(MatSetSizes(B, 4, 4, PETSC_DETERMINE, PETSC_DETERMINE));
 39:   PetscCall(MatSetType(B, MATAIJ));
 40:   PetscCall(MatMPIAIJSetPreallocation(B, 4, NULL, 4, NULL));

 42:   /* both matrices gain nonzeros after their first assembly */
 43:   PetscCall(MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE));
 44:   PetscCall(MatSetOption(B, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE));

 46:   PetscCall(MatGetOwnershipRange(A, &rstart, &rend));
 47:   for (i = rstart; i < rend; i++) {
 48:     PetscCall(MatSetValues(A, 1, &i, 1, &i, &value, INSERT_VALUES));
 49:     PetscCall(MatSetValues(B, 1, &i, 1, &i, &value, INSERT_VALUES));
 50:   }
 51:   /* one off-diagonal entry, so that its global column enters garray */
 52:   col = rank ? 0 : 4;
 53:   PetscCall(MatSetValues(A, 1, &rstart, 1, &col, &value, INSERT_VALUES));
 54:   PetscCall(MatSetValues(B, 1, &rstart, 1, &col, &value, INSERT_VALUES));
 55:   PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
 56:   PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
 57:   PetscCall(MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY));
 58:   PetscCall(MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY));

 60:   /* pull both submatrices onto the device */
 61:   PetscCall(MatMultEqual(A, B, 4, &flg));
 62:   PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "MatMult() differs before insertion");

 64:   /* new nonzero in the diagonal block */
 65:   i     = rstart + 1;
 66:   col   = rstart + 2;
 67:   value = 2.0;
 68:   PetscCall(MatSetValues(A, 1, &i, 1, &col, &value, INSERT_VALUES));
 69:   PetscCall(MatSetValues(B, 1, &i, 1, &col, &value, INSERT_VALUES));
 70:   PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
 71:   PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
 72:   PetscCall(MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY));
 73:   PetscCall(MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY));
 74:   PetscCall(MatMultEqual(A, B, 4, &flg));
 75:   PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "MatMult() differs after inserting into the diagonal block");

 77:   /* new nonzero in the off-diagonal block, in a column garray already knows */
 78:   i     = rstart + 1;
 79:   col   = rank ? 0 : 4;
 80:   value = 3.0;
 81:   PetscCall(MatSetValues(A, 1, &i, 1, &col, &value, INSERT_VALUES));
 82:   PetscCall(MatSetValues(B, 1, &i, 1, &col, &value, INSERT_VALUES));
 83:   PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
 84:   PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
 85:   PetscCall(MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY));
 86:   PetscCall(MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY));
 87:   PetscCall(MatMultEqual(A, B, 4, &flg));
 88:   PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "MatMult() differs after inserting into the off-diagonal block");

 90:   PetscCall(MatDestroy(&A));
 91:   PetscCall(MatDestroy(&B));
 92:   PetscCall(PetscFinalize());
 93:   return 0;
 94: }

 96: /*TEST

 98:    test:
 99:       suffix: 1
100:       nsize: 2
101:       args: -mat_type sell
102:       output_file: output/empty.out

104:    test:
105:       suffix: cuda
106:       nsize: 2
107:       requires: cuda !complex
108:       args: -mat_type sellcuda
109:       output_file: output/empty.out

111:    test:
112:       suffix: hip
113:       nsize: 2
114:       requires: hip !complex
115:       args: -mat_type sellhip
116:       output_file: output/empty.out

118: TEST*/