1 #include <petscmat.h> /*I "mat.h" I*/ 2 3 typedef struct{ 4 Mat A; 5 Vec D1; 6 Vec D2; 7 Vec W; 8 Vec W2; 9 Vec ADADiag; 10 PetscInt GotDiag; 11 } _p_TaoMatADACtx; 12 typedef _p_TaoMatADACtx* TaoMatADACtx; 13 14 extern PetscErrorCode MatCreateADA(Mat,Vec, Vec, Mat*); 15 16 #undef __FUNCT__ 17 #define __FUNCT__ "MatMult_ADA" 18 static PetscErrorCode MatMult_ADA(Mat mat,Vec a,Vec y) 19 { 20 TaoMatADACtx ctx; 21 PetscReal one = 1.0; 22 PetscErrorCode ierr; 23 24 PetscFunctionBegin; 25 ierr = MatShellGetContext(mat,(void **)&ctx);CHKERRQ(ierr); 26 ierr = MatMult(ctx->A,a,ctx->W);CHKERRQ(ierr); 27 if (ctx->D1){ 28 ierr = VecPointwiseMult(ctx->W,ctx->D1,ctx->W);CHKERRQ(ierr); 29 } 30 ierr = MatMultTranspose(ctx->A,ctx->W,y);CHKERRQ(ierr); 31 if (ctx->D2){ 32 ierr = VecPointwiseMult(ctx->W2, ctx->D2, a);CHKERRQ(ierr); 33 ierr = VecAXPY(y, one, ctx->W2);CHKERRQ(ierr); 34 } 35 PetscFunctionReturn(0); 36 } 37 38 #undef __FUNCT__ 39 #define __FUNCT__ "MatMultTranspose_ADA" 40 static PetscErrorCode MatMultTranspose_ADA(Mat mat,Vec a,Vec y) 41 { 42 PetscErrorCode ierr; 43 44 PetscFunctionBegin; 45 ierr = MatMult_ADA(mat,a,y);CHKERRQ(ierr); 46 PetscFunctionReturn(0); 47 } 48 49 #undef __FUNCT__ 50 #define __FUNCT__ "MatDiagonalSet_ADA" 51 static PetscErrorCode MatDiagonalSet_ADA(Vec D, Mat M) 52 { 53 TaoMatADACtx ctx; 54 PetscReal zero=0.0,one = 1.0; 55 PetscErrorCode ierr; 56 57 PetscFunctionBegin; 58 ierr = MatShellGetContext(M,(void **)&ctx);CHKERRQ(ierr); 59 if (!ctx->D2){ 60 ierr = VecDuplicate(D,&ctx->D2);CHKERRQ(ierr); 61 ierr = VecSet(ctx->D2, zero);CHKERRQ(ierr); 62 } 63 ierr = VecAXPY(ctx->D2, one, D);CHKERRQ(ierr); 64 PetscFunctionReturn(0); 65 } 66 67 #undef __FUNCT__ 68 #define __FUNCT__ "MatDestroy_ADA" 69 static PetscErrorCode MatDestroy_ADA(Mat mat) 70 { 71 PetscErrorCode ierr; 72 TaoMatADACtx ctx; 73 74 PetscFunctionBegin; 75 ierr=MatShellGetContext(mat,(void **)&ctx);CHKERRQ(ierr); 76 ierr=VecDestroy(&ctx->W);CHKERRQ(ierr); 77 ierr=VecDestroy(&ctx->W2);CHKERRQ(ierr); 78 ierr=VecDestroy(&ctx->ADADiag);CHKERRQ(ierr); 79 ierr=MatDestroy(&ctx->A);CHKERRQ(ierr); 80 ierr=VecDestroy(&ctx->D1);CHKERRQ(ierr); 81 ierr=VecDestroy(&ctx->D2);CHKERRQ(ierr); 82 ierr = PetscFree(ctx);CHKERRQ(ierr); 83 PetscFunctionReturn(0); 84 } 85 86 #undef __FUNCT__ 87 #define __FUNCT__ "MatView_ADA" 88 static PetscErrorCode MatView_ADA(Mat mat,PetscViewer viewer) 89 { 90 PetscFunctionBegin; 91 PetscFunctionReturn(0); 92 } 93 94 #undef __FUNCT__ 95 #define __FUNCT__ "MatShift_ADA" 96 static PetscErrorCode MatShift_ADA(Mat Y, PetscReal a) 97 { 98 PetscErrorCode ierr; 99 TaoMatADACtx ctx; 100 101 PetscFunctionBegin; 102 ierr = MatShellGetContext(Y,(void **)&ctx);CHKERRQ(ierr); 103 ierr = VecShift(ctx->D2,a);CHKERRQ(ierr); 104 PetscFunctionReturn(0); 105 } 106 107 #undef __FUNCT__ 108 #define __FUNCT__ "MatDuplicate_ADA" 109 static PetscErrorCode MatDuplicate_ADA(Mat mat,MatDuplicateOption op,Mat *M) 110 { 111 PetscErrorCode ierr; 112 TaoMatADACtx ctx; 113 Mat A2; 114 Vec D1b=NULL,D2b; 115 116 PetscFunctionBegin; 117 ierr = MatShellGetContext(mat,(void **)&ctx);CHKERRQ(ierr); 118 ierr = MatDuplicate(ctx->A,op,&A2);CHKERRQ(ierr); 119 if (ctx->D1){ 120 ierr = VecDuplicate(ctx->D1,&D1b);CHKERRQ(ierr); 121 ierr = VecCopy(ctx->D1,D1b);CHKERRQ(ierr); 122 } 123 ierr = VecDuplicate(ctx->D2,&D2b);CHKERRQ(ierr); 124 ierr = VecCopy(ctx->D2,D2b);CHKERRQ(ierr); 125 ierr = MatCreateADA(A2,D1b,D2b,M);CHKERRQ(ierr); 126 if (ctx->D1){ 127 ierr = PetscObjectDereference((PetscObject)D1b);CHKERRQ(ierr); 128 } 129 ierr = PetscObjectDereference((PetscObject)D2b);CHKERRQ(ierr); 130 ierr = PetscObjectDereference((PetscObject)A2);CHKERRQ(ierr); 131 PetscFunctionReturn(0); 132 } 133 134 #undef __FUNCT__ 135 #define __FUNCT__ "MatEqual_ADA" 136 static PetscErrorCode MatEqual_ADA(Mat A,Mat B,PetscBool *flg) 137 { 138 PetscErrorCode ierr; 139 TaoMatADACtx ctx1,ctx2; 140 141 PetscFunctionBegin; 142 ierr = MatShellGetContext(A,(void **)&ctx1);CHKERRQ(ierr); 143 ierr = MatShellGetContext(B,(void **)&ctx2);CHKERRQ(ierr); 144 ierr = VecEqual(ctx1->D2,ctx2->D2,flg);CHKERRQ(ierr); 145 if (*flg==PETSC_TRUE){ 146 ierr = VecEqual(ctx1->D1,ctx2->D1,flg);CHKERRQ(ierr); 147 } 148 if (*flg==PETSC_TRUE){ 149 ierr = MatEqual(ctx1->A,ctx2->A,flg);CHKERRQ(ierr); 150 } 151 PetscFunctionReturn(0); 152 } 153 154 #undef __FUNCT__ 155 #define __FUNCT__ "MatScale_ADA" 156 static PetscErrorCode MatScale_ADA(Mat mat, PetscReal a) 157 { 158 PetscErrorCode ierr; 159 TaoMatADACtx ctx; 160 161 PetscFunctionBegin; 162 ierr = MatShellGetContext(mat,(void **)&ctx);CHKERRQ(ierr); 163 ierr = VecScale(ctx->D1,a);CHKERRQ(ierr); 164 if (ctx->D2){ 165 ierr = VecScale(ctx->D2,a);CHKERRQ(ierr); 166 } 167 PetscFunctionReturn(0); 168 } 169 170 #undef __FUNCT__ 171 #define __FUNCT__ "MatTranspose_ADA" 172 static PetscErrorCode MatTranspose_ADA(Mat mat,MatReuse reuse,Mat *B) 173 { 174 PetscErrorCode ierr; 175 TaoMatADACtx ctx; 176 177 PetscFunctionBegin; 178 ierr = MatShellGetContext(mat,(void **)&ctx);CHKERRQ(ierr); 179 if (reuse == MAT_INITIAL_MATRIX){ 180 ierr = MatDuplicate(mat,MAT_COPY_VALUES,B);CHKERRQ(ierr); 181 } else if (reuse == MAT_REUSE_MATRIX){ 182 ierr = MatCopy(mat,*B,SAME_NONZERO_PATTERN);CHKERRQ(ierr); 183 } else SETERRQ(PetscObjectComm((PetscObject)mat),PETSC_ERR_SUP,"Does not support inplace transpose"); 184 PetscFunctionReturn(0); 185 } 186 187 #undef __FUNCT__ 188 #define __FUNCT__ "MatADAComputeDiagonal" 189 PetscErrorCode MatADAComputeDiagonal(Mat mat) 190 { 191 PetscErrorCode ierr; 192 PetscInt i,m,n,low,high; 193 PetscScalar *dtemp,*dptr; 194 TaoMatADACtx ctx; 195 196 PetscFunctionBegin; 197 ierr = MatShellGetContext(mat,(void **)&ctx);CHKERRQ(ierr); 198 ierr = MatGetOwnershipRange(mat, &low, &high);CHKERRQ(ierr); 199 ierr = MatGetSize(mat,&m,&n);CHKERRQ(ierr); 200 201 ierr = PetscMalloc1(n,&dtemp );CHKERRQ(ierr); 202 203 for (i=0; i<n; i++){ 204 ierr = MatGetColumnVector(ctx->A, ctx->W, i);CHKERRQ(ierr); 205 ierr = VecPointwiseMult(ctx->W,ctx->W,ctx->W);CHKERRQ(ierr); 206 ierr = VecDotBegin(ctx->D1, ctx->W,dtemp+i);CHKERRQ(ierr); 207 } 208 for (i=0; i<n; i++){ 209 ierr = VecDotEnd(ctx->D1, ctx->W,dtemp+i);CHKERRQ(ierr); 210 } 211 212 ierr = VecGetArray(ctx->ADADiag,&dptr);CHKERRQ(ierr); 213 for (i=low; i<high; i++){ 214 dptr[i-low]= dtemp[i]; 215 } 216 ierr = VecRestoreArray(ctx->ADADiag,&dptr);CHKERRQ(ierr); 217 ierr = PetscFree(dtemp);CHKERRQ(ierr); 218 PetscFunctionReturn(0); 219 } 220 221 #undef __FUNCT__ 222 #define __FUNCT__ "MatGetDiagonal_ADA" 223 static PetscErrorCode MatGetDiagonal_ADA(Mat mat,Vec v) 224 { 225 PetscErrorCode ierr; 226 PetscReal one=1.0; 227 TaoMatADACtx ctx; 228 229 PetscFunctionBegin; 230 ierr = MatShellGetContext(mat,(void **)&ctx);CHKERRQ(ierr); 231 ierr = MatADAComputeDiagonal(mat);CHKERRQ(ierr); 232 ierr=VecCopy(ctx->ADADiag,v);CHKERRQ(ierr); 233 if (ctx->D2){ 234 ierr=VecAXPY(v, one, ctx->D2);CHKERRQ(ierr); 235 } 236 PetscFunctionReturn(0); 237 } 238 239 #undef __FUNCT__ 240 #define __FUNCT__ "MatGetSubMatrix_ADA" 241 static PetscErrorCode MatGetSubMatrix_ADA(Mat mat,IS isrow,IS iscol,MatReuse cll, Mat *newmat) 242 { 243 PetscErrorCode ierr; 244 PetscInt low,high; 245 IS ISrow; 246 Vec D1,D2; 247 Mat Atemp; 248 TaoMatADACtx ctx; 249 PetscBool isequal; 250 251 PetscFunctionBegin; 252 ierr = ISEqual(isrow,iscol,&isequal);CHKERRQ(ierr); 253 if (!isequal) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Only for idential column and row indices"); 254 ierr = MatShellGetContext(mat,(void **)&ctx);CHKERRQ(ierr); 255 256 ierr = MatGetOwnershipRange(ctx->A,&low,&high);CHKERRQ(ierr); 257 ierr = ISCreateStride(PetscObjectComm((PetscObject)mat),high-low,low,1,&ISrow);CHKERRQ(ierr); 258 ierr = MatGetSubMatrix(ctx->A,ISrow,iscol,cll,&Atemp);CHKERRQ(ierr); 259 ierr = ISDestroy(&ISrow);CHKERRQ(ierr); 260 261 if (ctx->D1){ 262 ierr=VecDuplicate(ctx->D1,&D1);CHKERRQ(ierr); 263 ierr=VecCopy(ctx->D1,D1);CHKERRQ(ierr); 264 } else { 265 D1 = NULL; 266 } 267 268 if (ctx->D2){ 269 Vec D2sub; 270 271 ierr=VecGetSubVector(ctx->D2,isrow,&D2sub);CHKERRQ(ierr); 272 ierr=VecDuplicate(D2sub,&D2);CHKERRQ(ierr); 273 ierr=VecCopy(D2sub,D2);CHKERRQ(ierr); 274 ierr=VecRestoreSubVector(ctx->D2,isrow,&D2sub);CHKERRQ(ierr); 275 } else { 276 D2 = NULL; 277 } 278 279 ierr = MatCreateADA(Atemp,D1,D2,newmat);CHKERRQ(ierr); 280 ierr = MatShellGetContext(*newmat,(void **)&ctx);CHKERRQ(ierr); 281 ierr = PetscObjectDereference((PetscObject)Atemp);CHKERRQ(ierr); 282 if (ctx->D1){ 283 ierr = PetscObjectDereference((PetscObject)D1);CHKERRQ(ierr); 284 } 285 if (ctx->D2){ 286 ierr = PetscObjectDereference((PetscObject)D2);CHKERRQ(ierr); 287 } 288 PetscFunctionReturn(0); 289 } 290 291 #undef __FUNCT__ 292 #define __FUNCT__ "MatGetSubMatrices_ADA" 293 static PetscErrorCode MatGetSubMatrices_ADA(Mat A,PetscInt n, IS *irow,IS *icol,MatReuse scall,Mat **B) 294 { 295 PetscErrorCode ierr; 296 PetscInt i; 297 298 PetscFunctionBegin; 299 if (scall == MAT_INITIAL_MATRIX) { 300 ierr = PetscMalloc1(n+1,B );CHKERRQ(ierr); 301 } 302 for ( i=0; i<n; i++ ) { 303 ierr = MatGetSubMatrix_ADA(A,irow[i],icol[i],scall,&(*B)[i]);CHKERRQ(ierr); 304 } 305 PetscFunctionReturn(0); 306 } 307 308 #undef __FUNCT__ 309 #define __FUNCT__ "MatGetColumnVector_ADA" 310 static PetscErrorCode MatGetColumnVector_ADA(Mat mat,Vec Y, PetscInt col) 311 { 312 PetscErrorCode ierr; 313 PetscInt low,high; 314 PetscScalar zero=0.0,one=1.0; 315 316 PetscFunctionBegin; 317 ierr=VecSet(Y, zero);CHKERRQ(ierr); 318 ierr=VecGetOwnershipRange(Y,&low,&high);CHKERRQ(ierr); 319 if (col>=low && col<high){ 320 ierr=VecSetValue(Y,col,one,INSERT_VALUES);CHKERRQ(ierr); 321 } 322 ierr=VecAssemblyBegin(Y);CHKERRQ(ierr); 323 ierr=VecAssemblyEnd(Y);CHKERRQ(ierr); 324 ierr=MatMult_ADA(mat,Y,Y);CHKERRQ(ierr); 325 PetscFunctionReturn(0); 326 } 327 328 #undef __FUNCT__ 329 #define __FUNCT__ "MatConvert_ADA" 330 PETSC_INTERN PetscErrorCode MatConvert_ADA(Mat mat,MatType newtype,Mat *NewMat) 331 { 332 PetscErrorCode ierr; 333 PetscMPIInt size; 334 PetscBool sametype, issame, isdense, isseqdense; 335 TaoMatADACtx ctx; 336 337 PetscFunctionBegin; 338 ierr = MatShellGetContext(mat,(void **)&ctx);CHKERRQ(ierr); 339 ierr = MPI_Comm_size(PetscObjectComm((PetscObject)mat),&size);CHKERRQ(ierr); 340 341 ierr = PetscObjectTypeCompare((PetscObject)mat,newtype,&sametype);CHKERRQ(ierr); 342 ierr = PetscObjectTypeCompare((PetscObject)mat,MATSAME,&issame);CHKERRQ(ierr); 343 ierr = PetscObjectTypeCompare((PetscObject)mat,MATMPIDENSE,&isdense);CHKERRQ(ierr); 344 ierr = PetscObjectTypeCompare((PetscObject)mat,MATSEQDENSE,&isseqdense);CHKERRQ(ierr); 345 346 if (sametype || issame) { 347 ierr=MatDuplicate(mat,MAT_COPY_VALUES,NewMat);CHKERRQ(ierr); 348 } else if (isdense) { 349 PetscInt i,j,low,high,m,n,M,N; 350 const PetscScalar *dptr; 351 Vec X; 352 353 ierr = VecDuplicate(ctx->D2,&X);CHKERRQ(ierr); 354 ierr=MatGetSize(mat,&M,&N);CHKERRQ(ierr); 355 ierr=MatGetLocalSize(mat,&m,&n);CHKERRQ(ierr); 356 ierr = MatCreateDense(PetscObjectComm((PetscObject)mat),m,m,N,N,NULL,NewMat);CHKERRQ(ierr); 357 ierr = MatGetOwnershipRange(*NewMat,&low,&high);CHKERRQ(ierr); 358 for (i=0;i<M;i++){ 359 ierr = MatGetColumnVector_ADA(mat,X,i);CHKERRQ(ierr); 360 ierr = VecGetArrayRead(X,&dptr);CHKERRQ(ierr); 361 for (j=0; j<high-low; j++){ 362 ierr = MatSetValue(*NewMat,low+j,i,dptr[j],INSERT_VALUES);CHKERRQ(ierr); 363 } 364 ierr = VecRestoreArrayRead(X,&dptr);CHKERRQ(ierr); 365 } 366 ierr = MatAssemblyBegin(*NewMat,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 367 ierr = MatAssemblyEnd(*NewMat,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 368 ierr = VecDestroy(&X);CHKERRQ(ierr); 369 } else if (isseqdense && size==1){ 370 PetscInt i,j,low,high,m,n,M,N; 371 const PetscScalar *dptr; 372 Vec X; 373 374 ierr = VecDuplicate(ctx->D2,&X);CHKERRQ(ierr); 375 ierr = MatGetSize(mat,&M,&N);CHKERRQ(ierr); 376 ierr = MatGetLocalSize(mat,&m,&n);CHKERRQ(ierr); 377 ierr = MatCreateSeqDense(PetscObjectComm((PetscObject)mat),N,N,NULL,NewMat);CHKERRQ(ierr); 378 ierr = MatGetOwnershipRange(*NewMat,&low,&high);CHKERRQ(ierr); 379 for (i=0;i<M;i++){ 380 ierr = MatGetColumnVector_ADA(mat,X,i);CHKERRQ(ierr); 381 ierr = VecGetArrayRead(X,&dptr);CHKERRQ(ierr); 382 for (j=0; j<high-low; j++){ 383 ierr = MatSetValue(*NewMat,low+j,i,dptr[j],INSERT_VALUES);CHKERRQ(ierr); 384 } 385 ierr = VecRestoreArrayRead(X,&dptr);CHKERRQ(ierr); 386 } 387 ierr = MatAssemblyBegin(*NewMat,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 388 ierr = MatAssemblyEnd(*NewMat,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 389 ierr = VecDestroy(&X);CHKERRQ(ierr); 390 } else SETERRQ(PETSC_COMM_SELF,1,"No support to convert objects to that type"); 391 PetscFunctionReturn(0); 392 } 393 394 #undef __FUNCT__ 395 #define __FUNCT__ "MatNorm_ADA" 396 static PetscErrorCode MatNorm_ADA(Mat mat,NormType type,PetscReal *norm) 397 { 398 PetscErrorCode ierr; 399 TaoMatADACtx ctx; 400 401 PetscFunctionBegin; 402 ierr = MatShellGetContext(mat,(void **)&ctx);CHKERRQ(ierr); 403 if (type == NORM_FROBENIUS) { 404 *norm = 1.0; 405 } else if (type == NORM_1 || type == NORM_INFINITY) { 406 *norm = 1.0; 407 } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"No two norm"); 408 PetscFunctionReturn(0); 409 } 410 411 #undef __FUNCT__ 412 #define __FUNCT__ "MatCreateADA" 413 /*@C 414 MatCreateADA - Creates a matrix M=A^T D1 A + D2 where D1, D2 are diagonal 415 416 Collective on matrix 417 418 Input Parameters: 419 + mat - matrix of arbitrary type 420 . d1 - A vector with diagonal elements of D1 421 - d2 - A vector 422 423 Output Parameters: 424 . J - New matrix whose operations are defined in terms of mat, D1, and D2. 425 426 Notes: 427 The user provides the input data and is responsible for destroying 428 this data after matrix J has been destroyed. 429 The operation MatMult(A,D2,D1) must be well defined. 430 Before calling the operation MatGetDiagonal(), the function 431 MatADAComputeDiagonal() must be called. The matrices A and D1 must 432 be the same during calls to MatADAComputeDiagonal() and 433 MatGetDiagonal(). 434 435 Level: developer 436 437 .seealso: MatCreate() 438 @*/ 439 PetscErrorCode MatCreateADA(Mat mat,Vec d1, Vec d2, Mat *J) 440 { 441 MPI_Comm comm=PetscObjectComm((PetscObject)mat); 442 TaoMatADACtx ctx; 443 PetscErrorCode ierr; 444 PetscInt nloc,n; 445 446 PetscFunctionBegin; 447 ierr = PetscNew(&ctx);CHKERRQ(ierr); 448 ctx->A=mat; 449 ctx->D1=d1; 450 ctx->D2=d2; 451 if (d1){ 452 ierr = VecDuplicate(d1,&ctx->W);CHKERRQ(ierr); 453 ierr = PetscObjectReference((PetscObject)d1);CHKERRQ(ierr); 454 } else { 455 ctx->W = NULL; 456 } 457 if (d2){ 458 ierr = VecDuplicate(d2,&ctx->W2);CHKERRQ(ierr); 459 ierr = VecDuplicate(d2,&ctx->ADADiag);CHKERRQ(ierr); 460 ierr = PetscObjectReference((PetscObject)d2);CHKERRQ(ierr); 461 } else { 462 ctx->W2 = NULL; 463 ctx->ADADiag = NULL; 464 } 465 466 ctx->GotDiag=0; 467 ierr = PetscObjectReference((PetscObject)mat);CHKERRQ(ierr); 468 469 ierr=VecGetLocalSize(d2,&nloc);CHKERRQ(ierr); 470 ierr=VecGetSize(d2,&n);CHKERRQ(ierr); 471 472 ierr = MatCreateShell(comm,nloc,nloc,n,n,ctx,J);CHKERRQ(ierr); 473 474 ierr = MatShellSetOperation(*J,MATOP_MULT,(void(*)(void))MatMult_ADA);CHKERRQ(ierr); 475 ierr = MatShellSetOperation(*J,MATOP_DESTROY,(void(*)(void))MatDestroy_ADA);CHKERRQ(ierr); 476 ierr = MatShellSetOperation(*J,MATOP_VIEW,(void(*)(void))MatView_ADA);CHKERRQ(ierr); 477 ierr = MatShellSetOperation(*J,MATOP_MULT_TRANSPOSE,(void(*)(void))MatMultTranspose_ADA);CHKERRQ(ierr); 478 ierr = MatShellSetOperation(*J,MATOP_DIAGONAL_SET,(void(*)(void))MatDiagonalSet_ADA);CHKERRQ(ierr); 479 ierr = MatShellSetOperation(*J,MATOP_SHIFT,(void(*)(void))MatShift_ADA);CHKERRQ(ierr); 480 ierr = MatShellSetOperation(*J,MATOP_EQUAL,(void(*)(void))MatEqual_ADA);CHKERRQ(ierr); 481 ierr = MatShellSetOperation(*J,MATOP_SCALE,(void(*)(void))MatScale_ADA);CHKERRQ(ierr); 482 ierr = MatShellSetOperation(*J,MATOP_TRANSPOSE,(void(*)(void))MatTranspose_ADA);CHKERRQ(ierr); 483 ierr = MatShellSetOperation(*J,MATOP_GET_DIAGONAL,(void(*)(void))MatGetDiagonal_ADA);CHKERRQ(ierr); 484 ierr = MatShellSetOperation(*J,MATOP_GET_SUBMATRICES,(void(*)(void))MatGetSubMatrices_ADA);CHKERRQ(ierr); 485 ierr = MatShellSetOperation(*J,MATOP_NORM,(void(*)(void))MatNorm_ADA);CHKERRQ(ierr); 486 ierr = MatShellSetOperation(*J,MATOP_DUPLICATE,(void(*)(void))MatDuplicate_ADA);CHKERRQ(ierr); 487 ierr = MatShellSetOperation(*J,MATOP_GET_SUBMATRIX,(void(*)(void))MatGetSubMatrix_ADA);CHKERRQ(ierr); 488 489 ierr = PetscLogObjectParent((PetscObject)(*J),(PetscObject)ctx->W);CHKERRQ(ierr); 490 ierr = PetscLogObjectParent((PetscObject)mat,(PetscObject)(*J));CHKERRQ(ierr); 491 492 ierr = MatSetOption(*J,MAT_SYMMETRIC,PETSC_TRUE);CHKERRQ(ierr); 493 PetscFunctionReturn(0); 494 } 495