1c4762a1bSJed Brown static char help[] = "Tests Cholesky factorization for a SBAIJ matrix, (bs=2).\n"; 2c4762a1bSJed Brown /* 3c4762a1bSJed Brown This code is modified from the code contributed by JUNWANG@uwm.edu on Apr 13, 2007 4c4762a1bSJed Brown */ 5c4762a1bSJed Brown 6c4762a1bSJed Brown #include <petscmat.h> 7c4762a1bSJed Brown 8*d71ae5a4SJacob Faibussowitsch int main(int argc, char **args) 9*d71ae5a4SJacob Faibussowitsch { 10c4762a1bSJed Brown Mat mat, fact, B; 11c4762a1bSJed Brown PetscInt ind1[2], ind2[2]; 12c4762a1bSJed Brown PetscScalar temp[4]; 13c4762a1bSJed Brown PetscInt nnz[3]; 14c4762a1bSJed Brown IS perm, colp; 15c4762a1bSJed Brown MatFactorInfo info; 16c4762a1bSJed Brown PetscMPIInt size; 17c4762a1bSJed Brown 18327415f7SBarry Smith PetscFunctionBeginUser; 199566063dSJacob Faibussowitsch PetscCall(PetscInitialize(&argc, &args, 0, help)); 209566063dSJacob Faibussowitsch PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 21be096a46SBarry Smith PetscCheck(size == 1, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "This is a uniprocessor example only!"); 22c4762a1bSJed Brown 239371c9d4SSatish Balay nnz[0] = 2; 249371c9d4SSatish Balay nnz[1] = 1; 259371c9d4SSatish Balay nnz[2] = 1; 269566063dSJacob Faibussowitsch PetscCall(MatCreateSeqSBAIJ(PETSC_COMM_SELF, 2, 6, 6, 0, nnz, &mat)); 27c4762a1bSJed Brown 289371c9d4SSatish Balay ind1[0] = 0; 299371c9d4SSatish Balay ind1[1] = 1; 309371c9d4SSatish Balay temp[0] = 3; 319371c9d4SSatish Balay temp[1] = 2; 329371c9d4SSatish Balay temp[2] = 0; 339371c9d4SSatish Balay temp[3] = 3; 349566063dSJacob Faibussowitsch PetscCall(MatSetValues(mat, 2, ind1, 2, ind1, temp, INSERT_VALUES)); 359371c9d4SSatish Balay ind2[0] = 4; 369371c9d4SSatish Balay ind2[1] = 5; 379371c9d4SSatish Balay temp[0] = 1; 389371c9d4SSatish Balay temp[1] = 1; 399371c9d4SSatish Balay temp[2] = 2; 409371c9d4SSatish Balay temp[3] = 1; 419566063dSJacob Faibussowitsch PetscCall(MatSetValues(mat, 2, ind1, 2, ind2, temp, INSERT_VALUES)); 429371c9d4SSatish Balay ind1[0] = 2; 439371c9d4SSatish Balay ind1[1] = 3; 449371c9d4SSatish Balay temp[0] = 4; 459371c9d4SSatish Balay temp[1] = 1; 469371c9d4SSatish Balay temp[2] = 1; 479371c9d4SSatish Balay temp[3] = 5; 489566063dSJacob Faibussowitsch PetscCall(MatSetValues(mat, 2, ind1, 2, ind1, temp, INSERT_VALUES)); 499371c9d4SSatish Balay ind1[0] = 4; 509371c9d4SSatish Balay ind1[1] = 5; 519371c9d4SSatish Balay temp[0] = 5; 529371c9d4SSatish Balay temp[1] = 1; 539371c9d4SSatish Balay temp[2] = 1; 549371c9d4SSatish Balay temp[3] = 6; 559566063dSJacob Faibussowitsch PetscCall(MatSetValues(mat, 2, ind1, 2, ind1, temp, INSERT_VALUES)); 56c4762a1bSJed Brown 579566063dSJacob Faibussowitsch PetscCall(MatAssemblyBegin(mat, MAT_FINAL_ASSEMBLY)); 589566063dSJacob Faibussowitsch PetscCall(MatAssemblyEnd(mat, MAT_FINAL_ASSEMBLY)); 59c4762a1bSJed Brown 609566063dSJacob Faibussowitsch PetscCall(MatDuplicate(mat, MAT_SHARE_NONZERO_PATTERN, &B)); 619371c9d4SSatish Balay ind1[0] = 0; 629371c9d4SSatish Balay ind1[1] = 1; 639371c9d4SSatish Balay temp[0] = 3; 649371c9d4SSatish Balay temp[1] = 2; 659371c9d4SSatish Balay temp[2] = 0; 669371c9d4SSatish Balay temp[3] = 3; 679566063dSJacob Faibussowitsch PetscCall(MatSetValues(mat, 2, ind1, 2, ind1, temp, INSERT_VALUES)); 689371c9d4SSatish Balay ind2[0] = 4; 699371c9d4SSatish Balay ind2[1] = 5; 709371c9d4SSatish Balay temp[0] = 1; 719371c9d4SSatish Balay temp[1] = 1; 729371c9d4SSatish Balay temp[2] = 2; 739371c9d4SSatish Balay temp[3] = 1; 749566063dSJacob Faibussowitsch PetscCall(MatSetValues(mat, 2, ind1, 2, ind2, temp, INSERT_VALUES)); 759371c9d4SSatish Balay ind1[0] = 2; 769371c9d4SSatish Balay ind1[1] = 3; 779371c9d4SSatish Balay temp[0] = 4; 789371c9d4SSatish Balay temp[1] = 1; 799371c9d4SSatish Balay temp[2] = 1; 809371c9d4SSatish Balay temp[3] = 5; 819566063dSJacob Faibussowitsch PetscCall(MatSetValues(mat, 2, ind1, 2, ind1, temp, INSERT_VALUES)); 829371c9d4SSatish Balay ind1[0] = 4; 839371c9d4SSatish Balay ind1[1] = 5; 849371c9d4SSatish Balay temp[0] = 5; 859371c9d4SSatish Balay temp[1] = 1; 869371c9d4SSatish Balay temp[2] = 1; 879371c9d4SSatish Balay temp[3] = 6; 889566063dSJacob Faibussowitsch PetscCall(MatSetValues(mat, 2, ind1, 2, ind1, temp, INSERT_VALUES)); 89c4762a1bSJed Brown 909566063dSJacob Faibussowitsch PetscCall(MatAssemblyBegin(mat, MAT_FINAL_ASSEMBLY)); 919566063dSJacob Faibussowitsch PetscCall(MatAssemblyEnd(mat, MAT_FINAL_ASSEMBLY)); 92c4762a1bSJed Brown 939566063dSJacob Faibussowitsch PetscCall(PetscPrintf(PETSC_COMM_WORLD, "mat: \n")); 949566063dSJacob Faibussowitsch PetscCall(MatView(mat, PETSC_VIEWER_STDOUT_SELF)); 95c4762a1bSJed Brown 96c4762a1bSJed Brown /* begin cholesky factorization */ 979566063dSJacob Faibussowitsch PetscCall(MatGetOrdering(mat, MATORDERINGNATURAL, &perm, &colp)); 989566063dSJacob Faibussowitsch PetscCall(ISDestroy(&colp)); 99c4762a1bSJed Brown 100c4762a1bSJed Brown info.fill = 1.0; 1019566063dSJacob Faibussowitsch PetscCall(MatGetFactor(mat, MATSOLVERPETSC, MAT_FACTOR_CHOLESKY, &fact)); 1029566063dSJacob Faibussowitsch PetscCall(MatCholeskyFactorSymbolic(fact, mat, perm, &info)); 1039566063dSJacob Faibussowitsch PetscCall(MatCholeskyFactorNumeric(fact, mat, &info)); 104c4762a1bSJed Brown 1059566063dSJacob Faibussowitsch PetscCall(ISDestroy(&perm)); 1069566063dSJacob Faibussowitsch PetscCall(MatDestroy(&mat)); 1079566063dSJacob Faibussowitsch PetscCall(MatDestroy(&fact)); 1089566063dSJacob Faibussowitsch PetscCall(MatDestroy(&B)); 1099566063dSJacob Faibussowitsch PetscCall(PetscFinalize()); 110b122ec5aSJacob Faibussowitsch return 0; 111c4762a1bSJed Brown } 112