1*d0de2241SAndrew Spott 2*d0de2241SAndrew Spott #include <petsc-private/matimpl.h> /*I "petscmat.h" I*/ 3*d0de2241SAndrew Spott 4*d0de2241SAndrew Spott typedef struct { 5*d0de2241SAndrew Spott Mat A; 6*d0de2241SAndrew Spott } Mat_HTranspose; 7*d0de2241SAndrew Spott 8*d0de2241SAndrew Spott #undef __FUNCT__ 9*d0de2241SAndrew Spott #define __FUNCT__ "MatMult_HTranspose" 10*d0de2241SAndrew Spott PetscErrorCode MatMult_HTranspose(Mat N,Vec x,Vec y) 11*d0de2241SAndrew Spott { 12*d0de2241SAndrew Spott Mat_HTranspose *Na = (Mat_HTranspose*)N->data; 13*d0de2241SAndrew Spott PetscErrorCode ierr; 14*d0de2241SAndrew Spott 15*d0de2241SAndrew Spott PetscFunctionBegin; 16*d0de2241SAndrew Spott ierr = MatMultHermitianTranspose(Na->A,x,y);CHKERRQ(ierr); 17*d0de2241SAndrew Spott PetscFunctionReturn(0); 18*d0de2241SAndrew Spott } 19*d0de2241SAndrew Spott 20*d0de2241SAndrew Spott #undef __FUNCT__ 21*d0de2241SAndrew Spott #define __FUNCT__ "MatMultAdd_HTranspose" 22*d0de2241SAndrew Spott PetscErrorCode MatMultAdd_HTranspose(Mat N,Vec v1,Vec v2,Vec v3) 23*d0de2241SAndrew Spott { 24*d0de2241SAndrew Spott Mat_HTranspose *Na = (Mat_HTranspose*)N->data; 25*d0de2241SAndrew Spott PetscErrorCode ierr; 26*d0de2241SAndrew Spott 27*d0de2241SAndrew Spott PetscFunctionBegin; 28*d0de2241SAndrew Spott ierr = MatMultHermitianTransposeAdd(Na->A,v1,v2,v3);CHKERRQ(ierr); 29*d0de2241SAndrew Spott PetscFunctionReturn(0); 30*d0de2241SAndrew Spott } 31*d0de2241SAndrew Spott 32*d0de2241SAndrew Spott #undef __FUNCT__ 33*d0de2241SAndrew Spott #define __FUNCT__ "MatMultHTranspose_HTranspose" 34*d0de2241SAndrew Spott PetscErrorCode MatMultHTranspose_HTranspose(Mat N,Vec x,Vec y) 35*d0de2241SAndrew Spott { 36*d0de2241SAndrew Spott Mat_HTranspose *Na = (Mat_HTranspose*)N->data; 37*d0de2241SAndrew Spott PetscErrorCode ierr; 38*d0de2241SAndrew Spott 39*d0de2241SAndrew Spott PetscFunctionBegin; 40*d0de2241SAndrew Spott ierr = MatMult(Na->A,x,y);CHKERRQ(ierr); 41*d0de2241SAndrew Spott PetscFunctionReturn(0); 42*d0de2241SAndrew Spott } 43*d0de2241SAndrew Spott 44*d0de2241SAndrew Spott #undef __FUNCT__ 45*d0de2241SAndrew Spott #define __FUNCT__ "MatMultHTransposeAdd_HTranspose" 46*d0de2241SAndrew Spott PetscErrorCode MatMultHTransposeAdd_HTranspose(Mat N,Vec v1,Vec v2,Vec v3) 47*d0de2241SAndrew Spott { 48*d0de2241SAndrew Spott Mat_HTranspose *Na = (Mat_HTranspose*)N->data; 49*d0de2241SAndrew Spott PetscErrorCode ierr; 50*d0de2241SAndrew Spott 51*d0de2241SAndrew Spott PetscFunctionBegin; 52*d0de2241SAndrew Spott ierr = MatMultAdd(Na->A,v1,v2,v3);CHKERRQ(ierr); 53*d0de2241SAndrew Spott PetscFunctionReturn(0); 54*d0de2241SAndrew Spott } 55*d0de2241SAndrew Spott 56*d0de2241SAndrew Spott #undef __FUNCT__ 57*d0de2241SAndrew Spott #define __FUNCT__ "MatDestroy_HTranspose" 58*d0de2241SAndrew Spott PetscErrorCode MatDestroy_HTranspose(Mat N) 59*d0de2241SAndrew Spott { 60*d0de2241SAndrew Spott Mat_HTranspose *Na = (Mat_HTranspose*)N->data; 61*d0de2241SAndrew Spott PetscErrorCode ierr; 62*d0de2241SAndrew Spott 63*d0de2241SAndrew Spott PetscFunctionBegin; 64*d0de2241SAndrew Spott ierr = MatDestroy(&Na->A);CHKERRQ(ierr); 65*d0de2241SAndrew Spott ierr = PetscFree(N->data);CHKERRQ(ierr); 66*d0de2241SAndrew Spott PetscFunctionReturn(0); 67*d0de2241SAndrew Spott } 68*d0de2241SAndrew Spott 69*d0de2241SAndrew Spott #undef __FUNCT__ 70*d0de2241SAndrew Spott #define __FUNCT__ "MatDuplicate_HTranspose" 71*d0de2241SAndrew Spott PetscErrorCode MatDuplicate_HTranspose(Mat N, MatDuplicateOption op, Mat* m) 72*d0de2241SAndrew Spott { 73*d0de2241SAndrew Spott Mat_HTranspose *Na = (Mat_HTranspose*)N->data; 74*d0de2241SAndrew Spott PetscErrorCode ierr; 75*d0de2241SAndrew Spott 76*d0de2241SAndrew Spott PetscFunctionBegin; 77*d0de2241SAndrew Spott if (op == MAT_COPY_VALUES) { 78*d0de2241SAndrew Spott ierr = MatHermitianTranspose(Na->A,MAT_INITIAL_MATRIX,m);CHKERRQ(ierr); 79*d0de2241SAndrew Spott } else if (op == MAT_DO_NOT_COPY_VALUES) { 80*d0de2241SAndrew Spott ierr = MatDuplicate(Na->A,MAT_DO_NOT_COPY_VALUES,m);CHKERRQ(ierr); 81*d0de2241SAndrew Spott ierr = MatTranspose(*m,MAT_REUSE_MATRIX,m);CHKERRQ(ierr); 82*d0de2241SAndrew Spott } else { 83*d0de2241SAndrew Spott SETERRQ(PetscObjectComm((PetscObject)N),PETSC_ERR_SUP,"MAT_SHARE_NONZERO_PATTERN not supported for this matrix type"); 84*d0de2241SAndrew Spott } 85*d0de2241SAndrew Spott PetscFunctionReturn(0); 86*d0de2241SAndrew Spott } 87*d0de2241SAndrew Spott 88*d0de2241SAndrew Spott #undef __FUNCT__ 89*d0de2241SAndrew Spott #define __FUNCT__ "MatCreateHermitianTranspose" 90*d0de2241SAndrew Spott /*@ 91*d0de2241SAndrew Spott MatCreateHermitianTranspose - Creates a new matrix object that behaves like A'* 92*d0de2241SAndrew Spott 93*d0de2241SAndrew Spott Collective on Mat 94*d0de2241SAndrew Spott 95*d0de2241SAndrew Spott Input Parameter: 96*d0de2241SAndrew Spott . A - the (possibly rectangular) matrix 97*d0de2241SAndrew Spott 98*d0de2241SAndrew Spott Output Parameter: 99*d0de2241SAndrew Spott . N - the matrix that represents A'* 100*d0de2241SAndrew Spott 101*d0de2241SAndrew Spott Level: intermediate 102*d0de2241SAndrew Spott 103*d0de2241SAndrew Spott Notes: The hermitian transpose A' is NOT actually formed! Rather the new matrix 104*d0de2241SAndrew Spott object performs the matrix-vector product by using the MatMultHermitianTranspose() on 105*d0de2241SAndrew Spott the original matrix 106*d0de2241SAndrew Spott 107*d0de2241SAndrew Spott .seealso: MatCreateNormal(), MatMult(), MatMultHermitianTranspose(), MatCreate() 108*d0de2241SAndrew Spott 109*d0de2241SAndrew Spott @*/ 110*d0de2241SAndrew Spott PetscErrorCode MatCreateHermitianTranspose(Mat A,Mat *N) 111*d0de2241SAndrew Spott { 112*d0de2241SAndrew Spott PetscErrorCode ierr; 113*d0de2241SAndrew Spott PetscInt m,n; 114*d0de2241SAndrew Spott Mat_HTranspose *Na; 115*d0de2241SAndrew Spott 116*d0de2241SAndrew Spott PetscFunctionBegin; 117*d0de2241SAndrew Spott ierr = MatGetLocalSize(A,&m,&n);CHKERRQ(ierr); 118*d0de2241SAndrew Spott ierr = MatCreate(PetscObjectComm((PetscObject)A),N);CHKERRQ(ierr); 119*d0de2241SAndrew Spott ierr = MatSetSizes(*N,n,m,PETSC_DECIDE,PETSC_DECIDE);CHKERRQ(ierr); 120*d0de2241SAndrew Spott ierr = PetscLayoutSetUp((*N)->rmap);CHKERRQ(ierr); 121*d0de2241SAndrew Spott ierr = PetscLayoutSetUp((*N)->cmap);CHKERRQ(ierr); 122*d0de2241SAndrew Spott ierr = PetscObjectChangeTypeName((PetscObject)*N,MATTRANSPOSEMAT);CHKERRQ(ierr); 123*d0de2241SAndrew Spott 124*d0de2241SAndrew Spott ierr = PetscNewLog(*N,&Na);CHKERRQ(ierr); 125*d0de2241SAndrew Spott (*N)->data = (void*) Na; 126*d0de2241SAndrew Spott ierr = PetscObjectReference((PetscObject)A);CHKERRQ(ierr); 127*d0de2241SAndrew Spott Na->A = A; 128*d0de2241SAndrew Spott 129*d0de2241SAndrew Spott (*N)->ops->destroy = MatDestroy_HTranspose; 130*d0de2241SAndrew Spott (*N)->ops->mult = MatMult_HTranspose; 131*d0de2241SAndrew Spott (*N)->ops->multadd = MatMultAdd_HTranspose; 132*d0de2241SAndrew Spott (*N)->ops->multhermitiantranspose = MatMultHTranspose_HTranspose; 133*d0de2241SAndrew Spott (*N)->ops->multhermitiantransposeadd = MatMultHTransposeAdd_HTranspose; 134*d0de2241SAndrew Spott (*N)->ops->duplicate = MatDuplicate_HTranspose; 135*d0de2241SAndrew Spott (*N)->assembled = PETSC_TRUE; 136*d0de2241SAndrew Spott 137*d0de2241SAndrew Spott ierr = MatSetBlockSizes(*N,PetscAbs(A->cmap->bs),PetscAbs(A->rmap->bs));CHKERRQ(ierr); 138*d0de2241SAndrew Spott ierr = MatSetUp(*N);CHKERRQ(ierr); 139*d0de2241SAndrew Spott PetscFunctionReturn(0); 140*d0de2241SAndrew Spott } 141