1*fe1fc275SAlexander static char help[] = "Tests assembly of a matrix from another matrix's hash table.\n\n"; 2*fe1fc275SAlexander 3*fe1fc275SAlexander #include <petscmat.h> 4*fe1fc275SAlexander 5*fe1fc275SAlexander int main(int argc, char **argv) 6*fe1fc275SAlexander { 7*fe1fc275SAlexander Mat A, B; 8*fe1fc275SAlexander PetscInt m, n, i, j; 9*fe1fc275SAlexander PetscScalar v; 10*fe1fc275SAlexander 11*fe1fc275SAlexander PetscFunctionBeginUser; 12*fe1fc275SAlexander PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 13*fe1fc275SAlexander 14*fe1fc275SAlexander /* ------- Set values in A --------- */ 15*fe1fc275SAlexander PetscCall(MatCreate(PETSC_COMM_WORLD, &A)); 16*fe1fc275SAlexander PetscCall(MatSetSizes(A, 1, 1, PETSC_DETERMINE, PETSC_DETERMINE)); 17*fe1fc275SAlexander PetscCall(MatSetFromOptions(A)); 18*fe1fc275SAlexander PetscCall(MatSetUp(A)); 19*fe1fc275SAlexander PetscCall(MatGetSize(A, &m, &n)); 20*fe1fc275SAlexander for (i = 0; i < m; i++) { 21*fe1fc275SAlexander for (j = 0; j < n; j++) { 22*fe1fc275SAlexander v = 10.0 * i + j + 1; 23*fe1fc275SAlexander PetscCall(MatSetValues(A, 1, &i, 1, &j, &v, ADD_VALUES)); 24*fe1fc275SAlexander } 25*fe1fc275SAlexander } 26*fe1fc275SAlexander 27*fe1fc275SAlexander /* Create B */ 28*fe1fc275SAlexander PetscCall(MatDuplicate(A, MAT_DO_NOT_COPY_VALUES, &B)); 29*fe1fc275SAlexander PetscCall(MatCopyHashToXAIJ(A, B)); 30*fe1fc275SAlexander PetscCall(MatView(B, PETSC_VIEWER_STDOUT_WORLD)); 31*fe1fc275SAlexander 32*fe1fc275SAlexander PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY)); 33*fe1fc275SAlexander PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY)); 34*fe1fc275SAlexander PetscCall(MatView(A, PETSC_VIEWER_STDOUT_WORLD)); 35*fe1fc275SAlexander 36*fe1fc275SAlexander PetscCall(MatDestroy(&A)); 37*fe1fc275SAlexander PetscCall(MatDestroy(&B)); 38*fe1fc275SAlexander PetscCall(PetscFinalize()); 39*fe1fc275SAlexander return 0; 40*fe1fc275SAlexander } 41*fe1fc275SAlexander 42*fe1fc275SAlexander /*TEST 43*fe1fc275SAlexander 44*fe1fc275SAlexander test: 45*fe1fc275SAlexander suffix: seq 46*fe1fc275SAlexander args: -mat_type seqaij 47*fe1fc275SAlexander filter: grep -v "Mat Object" 48*fe1fc275SAlexander 49*fe1fc275SAlexander test: 50*fe1fc275SAlexander suffix: mpi 51*fe1fc275SAlexander args: -mat_type mpiaij 52*fe1fc275SAlexander nsize: 4 53*fe1fc275SAlexander filter: grep -v "Mat Object" 54*fe1fc275SAlexander 55*fe1fc275SAlexander TEST*/ 56