Actual source code: ex205.c

  1: static char help[] = "Tests MatCopy() for SHELL matrices\n\n";

  3: #include <petscmat.h>

  5: typedef struct _n_User *User;
  6: struct _n_User {
  7:   Mat A;
  8: };

 10: static PetscErrorCode MatMult_User(Mat A, Vec X, Vec Y)
 11: {
 12:   User user;

 14:   PetscFunctionBegin;
 15:   PetscCall(MatShellGetContext(A, &user));
 16:   PetscCall(MatMult(user->A, X, Y));
 17:   PetscFunctionReturn(PETSC_SUCCESS);
 18: }

 20: static PetscErrorCode MatCopy_User(Mat A, Mat B, MatStructure str)
 21: {
 22:   User userA, userB;

 24:   PetscFunctionBegin;
 25:   PetscCall(MatShellGetContext(A, &userA));
 26:   if (userA) {
 27:     PetscCall(PetscNew(&userB));
 28:     PetscCall(MatDuplicate(userA->A, MAT_COPY_VALUES, &userB->A));
 29:     PetscCall(MatShellSetContext(B, userB));
 30:   }
 31:   PetscFunctionReturn(PETSC_SUCCESS);
 32: }

 34: static PetscErrorCode MatDestroy_User(Mat A)
 35: {
 36:   User user;

 38:   PetscFunctionBegin;
 39:   PetscCall(MatShellGetContext(A, &user));
 40:   if (user) {
 41:     PetscCall(MatDestroy(&user->A));
 42:     PetscCall(PetscFree(user));
 43:   }
 44:   PetscFunctionReturn(PETSC_SUCCESS);
 45: }

 47: int main(int argc, char **args)
 48: {
 49:   const PetscScalar xvals[] = {11, 13}, yvals[] = {17, 19};
 50:   const PetscInt    inds[]  = {0, 1};
 51:   PetscScalar       avals[] = {2, 3, 5, 7};
 52:   Mat               S1, S2;
 53:   Vec               X, Y;
 54:   User              user;

 56:   PetscFunctionBeginUser;
 57:   PetscCall(PetscInitialize(&argc, &args, NULL, help));

 59:   PetscCall(PetscNew(&user));
 60:   PetscCall(MatCreateSeqAIJ(PETSC_COMM_WORLD, 2, 2, 2, NULL, &user->A));
 61:   PetscCall(MatSetUp(user->A));
 62:   PetscCall(MatSetValues(user->A, 2, inds, 2, inds, avals, INSERT_VALUES));
 63:   PetscCall(MatAssemblyBegin(user->A, MAT_FINAL_ASSEMBLY));
 64:   PetscCall(MatAssemblyEnd(user->A, MAT_FINAL_ASSEMBLY));
 65:   PetscCall(VecCreateSeq(PETSC_COMM_WORLD, 2, &X));
 66:   PetscCall(VecSetValues(X, 2, inds, xvals, INSERT_VALUES));
 67:   PetscCall(VecAssemblyBegin(X));
 68:   PetscCall(VecAssemblyEnd(X));
 69:   PetscCall(VecDuplicate(X, &Y));
 70:   PetscCall(VecSetValues(Y, 2, inds, yvals, INSERT_VALUES));
 71:   PetscCall(VecAssemblyBegin(Y));
 72:   PetscCall(VecAssemblyEnd(Y));

 74:   PetscCall(MatCreateShell(PETSC_COMM_WORLD, 2, 2, 2, 2, user, &S1));
 75:   PetscCall(MatSetUp(S1));
 76:   PetscCall(MatShellSetOperation(S1, MATOP_MULT, (void (*)(void))MatMult_User));
 77:   PetscCall(MatShellSetOperation(S1, MATOP_COPY, (void (*)(void))MatCopy_User));
 78:   PetscCall(MatShellSetOperation(S1, MATOP_DESTROY, (void (*)(void))MatDestroy_User));
 79:   PetscCall(MatCreateShell(PETSC_COMM_WORLD, 2, 2, 2, 2, NULL, &S2));
 80:   PetscCall(MatSetUp(S2));
 81:   PetscCall(MatShellSetOperation(S2, MATOP_MULT, (void (*)(void))MatMult_User));
 82:   PetscCall(MatShellSetOperation(S2, MATOP_COPY, (void (*)(void))MatCopy_User));
 83:   PetscCall(MatShellSetOperation(S2, MATOP_DESTROY, (void (*)(void))MatDestroy_User));

 85:   PetscCall(MatScale(S1, 31));
 86:   PetscCall(MatShift(S1, 37));
 87:   PetscCall(MatDiagonalScale(S1, X, Y));
 88:   PetscCall(MatCopy(S1, S2, SAME_NONZERO_PATTERN));
 89:   PetscCall(MatMult(S1, X, Y));
 90:   PetscCall(VecView(Y, PETSC_VIEWER_STDOUT_WORLD));
 91:   PetscCall(MatMult(S2, X, Y));
 92:   PetscCall(VecView(Y, PETSC_VIEWER_STDOUT_WORLD));

 94:   PetscCall(MatDestroy(&S1));
 95:   PetscCall(MatDestroy(&S2));
 96:   PetscCall(VecDestroy(&X));
 97:   PetscCall(VecDestroy(&Y));
 98:   PetscCall(PetscFinalize());
 99:   return 0;
100: }

102: /*TEST

104:    test:
105:       args: -malloc_dump

107: TEST*/