1c4762a1bSJed Brown static char help[] = "Tests MatComputeOperator() and MatComputeOperatorTranspose()\n\n"; 2c4762a1bSJed Brown 3c4762a1bSJed Brown #include <petscmat.h> 4c4762a1bSJed Brown 5*d71ae5a4SJacob Faibussowitsch int main(int argc, char **argv) 6*d71ae5a4SJacob Faibussowitsch { 7c4762a1bSJed Brown Mat A, Ae, Aet; 8c4762a1bSJed Brown char filename[PETSC_MAX_PATH_LEN]; 9c4762a1bSJed Brown char expltype[128], *etype = NULL; 10c4762a1bSJed Brown PetscInt bs = 1; 11c4762a1bSJed Brown PetscBool flg, check = PETSC_TRUE; 12c4762a1bSJed Brown 13327415f7SBarry Smith PetscFunctionBeginUser; 149566063dSJacob Faibussowitsch PetscCall(PetscInitialize(&argc, &argv, (char *)0, help)); 15c4762a1bSJed Brown 169566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetString(NULL, NULL, "-expl_type", expltype, sizeof(expltype), &flg)); 1748a46eb9SPierre Jolivet if (flg) PetscCall(PetscStrallocpy(expltype, &etype)); 189566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetString(NULL, NULL, "-f", filename, sizeof(filename), &flg)); 199566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetInt(NULL, NULL, "-bs", &bs, NULL)); 20c4762a1bSJed Brown if (!flg) { 21c4762a1bSJed Brown PetscInt M = 13, N = 6; 22c4762a1bSJed Brown 239566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetInt(NULL, NULL, "-M", &M, NULL)); 249566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetInt(NULL, NULL, "-N", &N, NULL)); 259566063dSJacob Faibussowitsch PetscCall(MatCreateDense(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE, M, N, NULL, &A)); 269566063dSJacob Faibussowitsch PetscCall(MatSetBlockSize(A, bs)); 279566063dSJacob Faibussowitsch PetscCall(MatSetRandom(A, NULL)); 28c4762a1bSJed Brown } else { 29c4762a1bSJed Brown PetscViewer viewer; 30c4762a1bSJed Brown 319566063dSJacob Faibussowitsch PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, filename, FILE_MODE_READ, &viewer)); 329566063dSJacob Faibussowitsch PetscCall(MatCreate(PETSC_COMM_WORLD, &A)); 339566063dSJacob Faibussowitsch PetscCall(MatSetBlockSize(A, bs)); 349566063dSJacob Faibussowitsch PetscCall(MatSetFromOptions(A)); 359566063dSJacob Faibussowitsch PetscCall(MatLoad(A, viewer)); 369566063dSJacob Faibussowitsch PetscCall(PetscViewerDestroy(&viewer)); 37c4762a1bSJed Brown } 389566063dSJacob Faibussowitsch PetscCall(PetscObjectSetName((PetscObject)A, "Matrix")); 399566063dSJacob Faibussowitsch PetscCall(MatViewFromOptions(A, NULL, "-view_expl")); 40c4762a1bSJed Brown 419566063dSJacob Faibussowitsch PetscCall(MatComputeOperator(A, etype, &Ae)); 429566063dSJacob Faibussowitsch PetscCall(PetscObjectSetName((PetscObject)Ae, "Explicit matrix")); 439566063dSJacob Faibussowitsch PetscCall(MatViewFromOptions(Ae, NULL, "-view_expl")); 44c4762a1bSJed Brown 459566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetBool(NULL, NULL, "-check", &check, NULL)); 46c4762a1bSJed Brown if (check) { 47c4762a1bSJed Brown Mat A2; 48c4762a1bSJed Brown PetscReal err, tol = PETSC_SMALL; 49c4762a1bSJed Brown 509566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetReal(NULL, NULL, "-tol", &tol, NULL)); 519566063dSJacob Faibussowitsch PetscCall(MatConvert(A, etype, MAT_INITIAL_MATRIX, &A2)); 529566063dSJacob Faibussowitsch PetscCall(MatAXPY(A2, -1.0, Ae, DIFFERENT_NONZERO_PATTERN)); 539566063dSJacob Faibussowitsch PetscCall(MatNorm(A2, NORM_FROBENIUS, &err)); 5448a46eb9SPierre Jolivet if (err > tol) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Error %g > %g (type %s)\n", (double)err, (double)tol, etype)); 559566063dSJacob Faibussowitsch PetscCall(MatDestroy(&A2)); 56c4762a1bSJed Brown } 57c4762a1bSJed Brown 589566063dSJacob Faibussowitsch PetscCall(MatComputeOperatorTranspose(A, etype, &Aet)); 599566063dSJacob Faibussowitsch PetscCall(PetscObjectSetName((PetscObject)Aet, "Explicit matrix transpose")); 609566063dSJacob Faibussowitsch PetscCall(MatViewFromOptions(Aet, NULL, "-view_expl")); 61c4762a1bSJed Brown 629566063dSJacob Faibussowitsch PetscCall(PetscFree(etype)); 639566063dSJacob Faibussowitsch PetscCall(MatDestroy(&Ae)); 649566063dSJacob Faibussowitsch PetscCall(MatDestroy(&Aet)); 659566063dSJacob Faibussowitsch PetscCall(MatDestroy(&A)); 669566063dSJacob Faibussowitsch PetscCall(PetscFinalize()); 67b122ec5aSJacob Faibussowitsch return 0; 68c4762a1bSJed Brown } 69c4762a1bSJed Brown 70c4762a1bSJed Brown /*TEST 71c4762a1bSJed Brown 72c4762a1bSJed Brown test: 73c4762a1bSJed Brown output_file: output/ex222_null.out 74c4762a1bSJed Brown 75c4762a1bSJed Brown testset: 76c4762a1bSJed Brown suffix: matexpl_rect 77c4762a1bSJed Brown output_file: output/ex222_null.out 78c4762a1bSJed Brown nsize: {{1 3}} 79c4762a1bSJed Brown args: -expl_type {{dense aij baij}} 80c4762a1bSJed Brown 81c4762a1bSJed Brown testset: 82c4762a1bSJed Brown suffix: matexpl_square 83c4762a1bSJed Brown output_file: output/ex222_null.out 84c4762a1bSJed Brown nsize: {{1 3}} 85c4762a1bSJed Brown args: -bs {{1 2 3}} -M 36 -N 36 -expl_type {{dense aij baij sbaij}} 86c4762a1bSJed Brown 87c4762a1bSJed Brown TEST*/ 88