Lines Matching defs:A
63 static PetscErrorCode MatAssemblyEnd_SeqAIJKokkos(Mat A, MatAssemblyType mode)
70 PetscCall(MatAssemblyEnd_SeqAIJ(A, mode));
72 aijseq = static_cast<Mat_SeqAIJ *>(A->data);
73 aijkok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
79 if (!aijkok || aijkok->nonzerostate != A->nonzerostate) { /* aijkok might not exist yet or nonzero pattern has changed */
96 aijkok = new Mat_SeqAIJKokkos(A, A->rmap->n, A->cmap->n, aijseq, A->nonzerostate, PETSC_FALSE /* don't copy mat values to device */);
97 A->spptr = aijkok;
98 } else if (A->rmap->n && aijkok->diag_dual.extent(0) == 0) { // MatProduct might directly produce AIJ on device, but not the diag.
101 PetscCall(MatGetDiagonalMarkers_SeqAIJ(A, &adiag, NULL));
102 MatRowMapKokkosViewHost diag_h(aijseq->diag, A->rmap->n);
110 PETSC_INTERN PetscErrorCode MatSeqAIJKokkosSyncDevice(Mat A)
112 Mat_SeqAIJKokkos *aijkok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
115 PetscCheck(A->factortype == MAT_FACTOR_NONE, PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "Can't sync factorized matrix from host to device");
116 PetscCheck(aijkok, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Unexpected NULL (Mat_SeqAIJKokkos*)A->spptr");
126 PETSC_INTERN PetscErrorCode MatSeqAIJKokkosModifyDevice(Mat A)
128 Mat_SeqAIJKokkos *aijkok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
131 PetscCheck(A->factortype == MAT_FACTOR_NONE, PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "Not supported for factorized matries");
136 PetscCall(PetscObjectStateIncrease((PetscObject)A));
140 static PetscErrorCode MatSeqAIJKokkosSyncHost(Mat A)
142 Mat_SeqAIJKokkos *aijkok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
146 PetscCheckTypeName(A, MATSEQAIJKOKKOS);
148 PetscCheck(A->factortype == MAT_FACTOR_NONE, PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "Can't sync factorized matrix from device to host");
149 PetscCheck(aijkok, PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "Missing AIJKOK");
154 static PetscErrorCode MatSeqAIJGetArray_SeqAIJKokkos(Mat A, PetscScalar *array[])
156 Mat_SeqAIJKokkos *aijkok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
164 if (aijkok && A->nonzerostate == aijkok->nonzerostate) {
168 *array = static_cast<Mat_SeqAIJ *>(A->data)->a;
173 static PetscErrorCode MatSeqAIJRestoreArray_SeqAIJKokkos(Mat A, PetscScalar *array[])
175 Mat_SeqAIJKokkos *aijkok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
178 if (aijkok && A->nonzerostate == aijkok->nonzerostate) aijkok->a_dual.modify_host();
182 static PetscErrorCode MatSeqAIJGetArrayRead_SeqAIJKokkos(Mat A, const PetscScalar *array[])
184 Mat_SeqAIJKokkos *aijkok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
187 if (aijkok && A->nonzerostate == aijkok->nonzerostate) {
191 *array = static_cast<Mat_SeqAIJ *>(A->data)->a;
196 static PetscErrorCode MatSeqAIJRestoreArrayRead_SeqAIJKokkos(Mat A, const PetscScalar *array[])
203 static PetscErrorCode MatSeqAIJGetArrayWrite_SeqAIJKokkos(Mat A, PetscScalar *array[])
205 Mat_SeqAIJKokkos *aijkok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
208 if (aijkok && A->nonzerostate == aijkok->nonzerostate) {
211 *array = static_cast<Mat_SeqAIJ *>(A->data)->a;
216 static PetscErrorCode MatSeqAIJRestoreArrayWrite_SeqAIJKokkos(Mat A, PetscScalar *array[])
218 Mat_SeqAIJKokkos *aijkok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
221 if (aijkok && A->nonzerostate == aijkok->nonzerostate) {
228 static PetscErrorCode MatSeqAIJGetCSRAndMemType_SeqAIJKokkos(Mat A, const PetscInt **i, const PetscInt **j, PetscScalar **a, PetscMemType *mtype)
230 Mat_SeqAIJKokkos *aijkok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
233 PetscCheck(aijkok != NULL, PetscObjectComm((PetscObject)A), PETSC_ERR_ARG_WRONGSTATE, "aijkok is NULL");
245 static PetscErrorCode MatGetCurrentMemType_SeqAIJKokkos(PETSC_UNUSED Mat A, PetscMemType *mtype)
256 . A - the MATSEQAIJKOKKOS matrix
262 static PetscErrorCode MatSeqAIJKokkosGenerateTransposeStructure(Mat A, MatRowMapKokkosView &perm_d, KokkosCsrMatrix &T_d)
264 Mat_SeqAIJ *aseq = static_cast<Mat_SeqAIJ *>(A->data);
265 PetscInt nz = aseq->nz, m = A->rmap->N, n = A->cmap->n;
286 for (PetscInt j = Ai[i]; j < Ai[i + 1]; j++) { // A's (i,j) is T's (j,i)
310 PETSC_INTERN PetscErrorCode MatSeqAIJKokkosGenerateTranspose_Private(Mat A, KokkosCsrMatrix *csrmatT)
312 Mat_SeqAIJ *aseq = static_cast<Mat_SeqAIJ *>(A->data);
313 Mat_SeqAIJKokkos *akok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
314 PetscInt nz = aseq->nz, m = A->rmap->N, n = A->cmap->n;
318 PetscCheck(akok, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Unexpected NULL (Mat_SeqAIJKokkos*)A->spptr");
319 PetscCall(KokkosDualViewSyncDevice(akok->a_dual, PetscGetKokkosExecutionSpace())); // Sync A's values since we are going to access them on device
323 if (A->symmetric == PETSC_BOOL3_TRUE) {
337 PetscCall(MatSeqAIJKokkosGenerateTransposeStructure(A, perm, T));
348 static PetscErrorCode MatSeqAIJKokkosGenerateHermitian_Private(Mat A, KokkosCsrMatrix *csrmatH)
350 Mat_SeqAIJ *aseq = static_cast<Mat_SeqAIJ *>(A->data);
351 Mat_SeqAIJKokkos *akok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
352 PetscInt nz = aseq->nz, m = A->rmap->N, n = A->cmap->n;
356 PetscCheck(akok, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Unexpected NULL (Mat_SeqAIJKokkos*)A->spptr");
357 PetscCall(KokkosDualViewSyncDevice(akok->a_dual, PetscGetKokkosExecutionSpace())); // Sync A's values since we are going to access them on device
361 if (A->hermitian == PETSC_BOOL3_TRUE) {
375 PetscCall(MatSeqAIJKokkosGenerateTransposeStructure(A, perm, T));
385 /* y = A x */
386 static PetscErrorCode MatMult_SeqAIJKokkos(Mat A, Vec xx, Vec yy)
394 PetscCall(MatSeqAIJKokkosSyncDevice(A));
397 aijkok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
398 PetscCallCXX(KokkosSparse::spmv(PetscGetKokkosExecutionSpace(), "N", 1.0 /*alpha*/, aijkok->csrmat, xv, 0.0 /*beta*/, yv)); /* y = alpha A x + beta y */
407 /* y = A^T x */
408 static PetscErrorCode MatMultTranspose_SeqAIJKokkos(Mat A, Vec xx, Vec yy)
418 PetscCall(MatSeqAIJKokkosSyncDevice(A));
421 if (A->form_explicit_transpose) {
422 PetscCall(MatSeqAIJKokkosGenerateTranspose_Private(A, &csrmat));
425 aijkok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
429 PetscCallCXX(KokkosSparse::spmv(PetscGetKokkosExecutionSpace(), mode, 1.0 /*alpha*/, csrmat, xv, 0.0 /*beta*/, yv)); /* y = alpha A^T x + beta y */
437 /* y = A^H x */
438 static PetscErrorCode MatMultHermitianTranspose_SeqAIJKokkos(Mat A, Vec xx, Vec yy)
448 PetscCall(MatSeqAIJKokkosSyncDevice(A));
451 if (A->form_explicit_transpose) {
452 PetscCall(MatSeqAIJKokkosGenerateHermitian_Private(A, &csrmat));
455 aijkok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
459 PetscCallCXX(KokkosSparse::spmv(PetscGetKokkosExecutionSpace(), mode, 1.0 /*alpha*/, csrmat, xv, 0.0 /*beta*/, yv)); /* y = alpha A^H x + beta y */
467 /* z = A x + y */
468 static PetscErrorCode MatMultAdd_SeqAIJKokkos(Mat A, Vec xx, Vec yy, Vec zz)
476 PetscCall(MatSeqAIJKokkosSyncDevice(A));
480 aijkok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
481 PetscCallCXX(KokkosSparse::spmv(PetscGetKokkosExecutionSpace(), "N", 1.0 /*alpha*/, aijkok->csrmat, xv, 1.0 /*beta*/, zv)); /* z = alpha A x + beta z */
489 /* z = A^T x + y */
490 static PetscErrorCode MatMultTransposeAdd_SeqAIJKokkos(Mat A, Vec xx, Vec yy, Vec zz)
500 PetscCall(MatSeqAIJKokkosSyncDevice(A));
504 if (A->form_explicit_transpose) {
505 PetscCall(MatSeqAIJKokkosGenerateTranspose_Private(A, &csrmat));
508 aijkok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
512 PetscCallCXX(KokkosSparse::spmv(PetscGetKokkosExecutionSpace(), mode, 1.0 /*alpha*/, csrmat, xv, 1.0 /*beta*/, zv)); /* z = alpha A^T x + beta z */
520 /* z = A^H x + y */
521 static PetscErrorCode MatMultHermitianTransposeAdd_SeqAIJKokkos(Mat A, Vec xx, Vec yy, Vec zz)
531 PetscCall(MatSeqAIJKokkosSyncDevice(A));
535 if (A->form_explicit_transpose) {
536 PetscCall(MatSeqAIJKokkosGenerateHermitian_Private(A, &csrmat));
539 aijkok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
543 PetscCallCXX(KokkosSparse::spmv(PetscGetKokkosExecutionSpace(), mode, 1.0 /*alpha*/, csrmat, xv, 1.0 /*beta*/, zv)); /* z = alpha A^H x + beta z */
551 static PetscErrorCode MatSetOption_SeqAIJKokkos(Mat A, MatOption op, PetscBool flg)
553 Mat_SeqAIJKokkos *aijkok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
559 if (A->form_explicit_transpose && !flg && aijkok) PetscCall(aijkok->DestroyMatTranspose());
560 A->form_explicit_transpose = flg;
563 PetscCall(MatSetOption_SeqAIJ(A, op, flg));
570 PETSC_INTERN PetscErrorCode MatConvert_SeqAIJ_SeqAIJKokkos(Mat A, MatType mtype, MatReuse reuse, Mat *newmat)
577 PetscCall(MatDuplicate(A, MAT_COPY_VALUES, newmat));
580 PetscCall(MatCopy(A, *newmat, SAME_NONZERO_PATTERN)); /* newmat is already a SeqAIJKokkos */
581 } else if (reuse == MAT_INPLACE_MATRIX) { /* newmat is A */
582 PetscCheck(A == *newmat, PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "A != *newmat with MAT_INPLACE_MATRIX");
583 PetscCall(PetscFree(A->defaultvectype));
584 PetscCall(PetscStrallocpy(VECKOKKOS, &A->defaultvectype)); /* Allocate and copy the string */
585 PetscCall(PetscObjectChangeTypeName((PetscObject)A, MATSEQAIJKOKKOS));
586 PetscCall(MatSetOps_SeqAIJKokkos(A));
587 aseq = static_cast<Mat_SeqAIJ *>(A->data);
588 if (A->assembled) { /* Copy i, j (but not values) to device for an assembled matrix if not yet */
589 PetscCheck(!A->spptr, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Expect NULL (Mat_SeqAIJKokkos*)A->spptr");
590 A->spptr = new Mat_SeqAIJKokkos(A, A->rmap->n, A->cmap->n, aseq, A->nonzerostate, PETSC_FALSE);
599 static PetscErrorCode MatDuplicate_SeqAIJKokkos(Mat A, MatDuplicateOption dupOption, Mat *B)
602 Mat_SeqAIJKokkos *akok = static_cast<Mat_SeqAIJKokkos *>(A->spptr), *bkok;
606 /* Do not copy values on host as A's latest values might be on device. We don't want to do sync blindly */
607 PetscCall(MatDuplicate_SeqAIJ(A, MAT_DO_NOT_COPY_VALUES, B));
609 if (A->assembled) {
636 static PetscErrorCode MatTranspose_SeqAIJKokkos(Mat A, MatReuse reuse, Mat *B)
643 if (reuse == MAT_REUSE_MATRIX) PetscCall(MatTransposeCheckNonzeroState_Private(A, *B));
644 PetscCall(MatSeqAIJKokkosGenerateTranspose_Private(A, &internT)); /* Generate a transpose internally */
648 PetscCall(MatCreateSeqAIJKokkosWithCSRMatrix(PetscObjectComm((PetscObject)A), atkok, &At));
650 else PetscCall(MatHeaderReplace(A, &At)); /* Replace A with At inplace */
662 } else SETERRQ(PetscObjectComm((PetscObject)A), PETSC_ERR_ARG_WRONGSTATE, "B must be assembled or preallocated");
667 static PetscErrorCode MatDestroy_SeqAIJKokkos(Mat A)
672 if (A->factortype == MAT_FACTOR_NONE) {
673 aijkok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
676 delete static_cast<Mat_SeqAIJKokkosTriFactors *>(A->spptr);
678 A->spptr = NULL;
679 PetscCall(PetscObjectComposeFunction((PetscObject)A, "MatFactorGetSolverType_C", NULL));
680 PetscCall(PetscObjectComposeFunction((PetscObject)A, "MatSetPreallocationCOO_C", NULL));
681 PetscCall(PetscObjectComposeFunction((PetscObject)A, "MatSetValuesCOO_C", NULL));
683 PetscCall(PetscObjectComposeFunction((PetscObject)A, "MatConvert_seqaijkokkos_hypre_C", NULL));
685 PetscCall(MatDestroy_SeqAIJ(A));
690 MATSEQAIJKOKKOS - MATAIJKOKKOS = "(seq)aijkokkos" - A matrix type to be used for sparse matrices with Kokkos
692 A matrix type using Kokkos-Kernels CrsMatrix type for portability across different device types
701 PETSC_EXTERN PetscErrorCode MatCreate_SeqAIJKokkos(Mat A)
705 PetscCall(MatCreate_SeqAIJ(A));
706 PetscCall(MatConvert_SeqAIJ_SeqAIJKokkos(A, MATSEQAIJKOKKOS, MAT_INPLACE_MATRIX, &A));
710 /* Merge A, B into a matrix C. A is put before B. C's size would be A->rmap->n by (A->cmap->n + B->cmap->n) */
711 PetscErrorCode MatSeqAIJKokkosMergeMats(Mat A, Mat B, MatReuse reuse, Mat *C)
721 PetscValidHeaderSpecific(A, MAT_CLASSID, 1);
724 PetscCheckTypeName(A, MATSEQAIJKOKKOS);
726 PetscCheck(A->rmap->n == B->rmap->n, PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Invalid number or rows %" PetscInt_FMT " != %" PetscInt_FMT, A->rmap->n, B->rmap->n);
729 PetscCall(MatSeqAIJKokkosSyncDevice(A));
731 a = static_cast<Mat_SeqAIJ *>(A->data);
733 akok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
739 m = A->rmap->n; /* M, N and nnz of C */
740 n = A->cmap->n + B->cmap->n;
742 aN = A->cmap->n; /* N of A */
750 /* Concatenate A and B in parallel using Kokkos hierarchical parallelism */
804 Mat A, B;
816 // See if numeric has already been done in symbolic (e.g., user calls MatMatMult(A,B,MAT_INITIAL_MATRIX,..,C)).
841 A = product->A;
843 PetscCall(MatSeqAIJKokkosSyncDevice(A));
845 akok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
856 PetscCall(MatSeqAIJKokkosGenerateTranspose_Private(A, &csrmatA));
890 Mat A, B;
901 A = product->A;
903 PetscCall(MatSeqAIJKokkosSyncDevice(A));
905 akok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
912 if (A->symmetric == PETSC_BOOL3_TRUE && ptype == MATPRODUCT_AtB) {
925 PetscCall(MatSetBlockSizesFromMats(C, A, B));
930 if (A->cmap->bs > 0) PetscCall(PetscLayoutSetBlockSize(C->rmap, A->cmap->bs));
936 if (A->rmap->bs > 0) PetscCall(PetscLayoutSetBlockSize(C->rmap, A->rmap->bs));
959 PetscCall(MatSeqAIJKokkosGenerateTranspose_Private(A, &csrmatA));
1019 static PetscErrorCode MatScale_SeqAIJKokkos(Mat A, PetscScalar a)
1025 PetscCall(MatSeqAIJKokkosSyncDevice(A));
1026 aijkok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
1028 PetscCall(MatSeqAIJKokkosModifyDevice(A));
1034 // add a to A's diagonal (if A is square) or main diagonal (if A is rectangular)
1035 static PetscErrorCode MatShift_SeqAIJKokkos(Mat A, PetscScalar a)
1037 Mat_SeqAIJ *aijseq = static_cast<Mat_SeqAIJ *>(A->data);
1040 if (A->assembled && aijseq->diagDense) { // no missing diagonals
1041 PetscInt n = PetscMin(A->rmap->n, A->cmap->n);
1044 PetscCall(MatSeqAIJKokkosSyncDevice(A));
1045 const auto aijkok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
1049 PetscCall(MatSeqAIJKokkosModifyDevice(A));
1053 PetscCall(MatShift_Basic(A, a));
1092 static PetscErrorCode MatDiagonalScale_SeqAIJKokkos(Mat A, Vec ll, Vec rr)
1094 Mat_SeqAIJ *aijseq = static_cast<Mat_SeqAIJ *>(A->data);
1095 PetscInt m = A->rmap->n, n = A->cmap->n, nz = aijseq->nz;
1100 PetscCall(MatSeqAIJKokkosSyncDevice(A));
1101 const auto aijkok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
1107 PetscCheck(m == A->rmap->n, PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Left scaling vector wrong length");
1121 PetscCheck(n == A->cmap->n, PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Right scaling vector wrong length");
1128 PetscCall(MatSeqAIJKokkosModifyDevice(A));
1133 static PetscErrorCode MatZeroEntries_SeqAIJKokkos(Mat A)
1138 aijkok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
1141 PetscCall(MatSeqAIJKokkosModifyDevice(A));
1143 PetscCall(MatZeroEntries_SeqAIJ(A));
1148 static PetscErrorCode MatGetDiagonal_SeqAIJKokkos(Mat A, Vec x)
1156 PetscCheck(n == A->rmap->n, PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Nonconforming matrix and vector");
1157 PetscCheck(A->factortype == MAT_FACTOR_NONE, PETSC_COMM_SELF, PETSC_ERR_SUP, "MatGetDiagonal_SeqAIJKokkos not supported on factored matrices");
1159 PetscCall(MatSeqAIJKokkosSyncDevice(A));
1160 aijkok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
1177 PetscErrorCode MatSeqAIJGetKokkosView(Mat A, ConstMatScalarKokkosView *kv)
1182 PetscValidHeaderSpecific(A, MAT_CLASSID, 1);
1184 PetscCheckTypeName(A, MATSEQAIJKOKKOS);
1185 PetscCall(MatSeqAIJKokkosSyncDevice(A));
1186 aijkok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
1191 PetscErrorCode MatSeqAIJRestoreKokkosView(Mat A, ConstMatScalarKokkosView *kv)
1194 PetscValidHeaderSpecific(A, MAT_CLASSID, 1);
1196 PetscCheckTypeName(A, MATSEQAIJKOKKOS);
1200 PetscErrorCode MatSeqAIJGetKokkosView(Mat A, MatScalarKokkosView *kv)
1205 PetscValidHeaderSpecific(A, MAT_CLASSID, 1);
1207 PetscCheckTypeName(A, MATSEQAIJKOKKOS);
1208 PetscCall(MatSeqAIJKokkosSyncDevice(A));
1209 aijkok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
1214 PetscErrorCode MatSeqAIJRestoreKokkosView(Mat A, MatScalarKokkosView *kv)
1217 PetscValidHeaderSpecific(A, MAT_CLASSID, 1);
1219 PetscCheckTypeName(A, MATSEQAIJKOKKOS);
1220 PetscCall(MatSeqAIJKokkosModifyDevice(A));
1224 PetscErrorCode MatSeqAIJGetKokkosViewWrite(Mat A, MatScalarKokkosView *kv)
1229 PetscValidHeaderSpecific(A, MAT_CLASSID, 1);
1231 PetscCheckTypeName(A, MATSEQAIJKOKKOS);
1232 aijkok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
1237 PetscErrorCode MatSeqAIJRestoreKokkosViewWrite(Mat A, MatScalarKokkosView *kv)
1240 PetscValidHeaderSpecific(A, MAT_CLASSID, 1);
1242 PetscCheckTypeName(A, MATSEQAIJKOKKOS);
1243 PetscCall(MatSeqAIJKokkosModifyDevice(A));
1247 PetscErrorCode MatCreateSeqAIJKokkosWithKokkosViews(MPI_Comm comm, PetscInt m, PetscInt n, Kokkos::View<PetscInt *> &i_d, Kokkos::View<PetscInt *> &j_d, Kokkos::View<PetscScalar *> &a_d, Mat *A)
1253 PetscCall(MatCreate(comm, A));
1254 PetscCall(MatSetSeqAIJKokkosWithCSRMatrix(*A, akok));
1283 /* cusparseDcsrgeam2() computes C = alpha A + beta B. If one knew sparsity pattern of C, one can skip
1286 KokkosSparse::spadd(alpha,A,beta,B,C) has symbolic and numeric phases, MatAXPY does not.
1393 static PetscErrorCode MatSetValuesCOO_SeqAIJKokkos(Mat A, const PetscScalar v[], InsertMode imode)
1402 PetscCall(PetscObjectQuery((PetscObject)A, "__PETSc_MatCOOStruct_Device", (PetscObject *)&container));
1417 if (imode == INSERT_VALUES) PetscCall(MatSeqAIJGetKokkosViewWrite(A, &Aa)); /* write matrix values */
1418 else PetscCall(MatSeqAIJGetKokkosView(A, &Aa)); /* read & write matrix values */
1429 if (imode == INSERT_VALUES) PetscCall(MatSeqAIJRestoreKokkosViewWrite(A, &Aa));
1430 else PetscCall(MatSeqAIJRestoreKokkosView(A, &Aa));
1434 static PetscErrorCode MatSetOps_SeqAIJKokkos(Mat A)
1436 Mat_SeqAIJ *a = (Mat_SeqAIJ *)A->data;
1439 A->offloadmask = PETSC_OFFLOAD_KOKKOS; /* We do not really use this flag */
1440 A->boundtocpu = PETSC_FALSE;
1442 A->ops->assemblyend = MatAssemblyEnd_SeqAIJKokkos;
1443 A->ops->destroy = MatDestroy_SeqAIJKokkos;
1444 A->ops->duplicate = MatDuplicate_SeqAIJKokkos;
1445 A->ops->axpy = MatAXPY_SeqAIJKokkos;
1446 A->ops->scale = MatScale_SeqAIJKokkos;
1447 A->ops->zeroentries = MatZeroEntries_SeqAIJKokkos;
1448 A->ops->productsetfromoptions = MatProductSetFromOptions_SeqAIJKokkos;
1449 A->ops->mult = MatMult_SeqAIJKokkos;
1450 A->ops->multadd = MatMultAdd_SeqAIJKokkos;
1451 A->ops->multtranspose = MatMultTranspose_SeqAIJKokkos;
1452 A->ops->multtransposeadd = MatMultTransposeAdd_SeqAIJKokkos;
1453 A->ops->multhermitiantranspose = MatMultHermitianTranspose_SeqAIJKokkos;
1454 A->ops->multhermitiantransposeadd = MatMultHermitianTransposeAdd_SeqAIJKokkos;
1455 A->ops->productnumeric = MatProductNumeric_SeqAIJKokkos_SeqAIJKokkos;
1456 A->ops->transpose = MatTranspose_SeqAIJKokkos;
1457 A->ops->setoption = MatSetOption_SeqAIJKokkos;
1458 A->ops->getdiagonal = MatGetDiagonal_SeqAIJKokkos;
1459 A->ops->shift = MatShift_SeqAIJKokkos;
1460 A->ops->diagonalset = MatDiagonalSet_SeqAIJKokkos;
1461 A->ops->diagonalscale = MatDiagonalScale_SeqAIJKokkos;
1462 A->ops->getcurrentmemtype = MatGetCurrentMemType_SeqAIJKokkos;
1471 PetscCall(PetscObjectComposeFunction((PetscObject)A, "MatSetPreallocationCOO_C", MatSetPreallocationCOO_SeqAIJKokkos));
1472 PetscCall(PetscObjectComposeFunction((PetscObject)A, "MatSetValuesCOO_C", MatSetValuesCOO_SeqAIJKokkos));
1474 PetscCall(PetscObjectComposeFunction((PetscObject)A, "MatConvert_seqaijkokkos_hypre_C", MatConvert_AIJ_HYPRE));
1483 + A - the MATSEQAIJKOKKOS matrix
1492 PETSC_INTERN PetscErrorCode MatInvertVariableBlockDiagonal_SeqAIJKokkos(Mat A, const PetscIntKokkosView &bs, const PetscIntKokkosView &bs2, const PetscIntKokkosView &blkMap, PetscScalarKokkosView &work, PetscScalarKokkosView &diagVal)
1494 Mat_SeqAIJKokkos *akok = static_cast<Mat_SeqAIJKokkos *>(A->spptr);
1498 PetscCall(MatSeqAIJKokkosSyncDevice(A)); // Since we'll access A's value on device
1520 PetscInt i = rstart + r; // i-th row in A
1553 PETSC_INTERN PetscErrorCode MatSetSeqAIJKokkosWithCSRMatrix(Mat A, Mat_SeqAIJKokkos *akok)
1560 PetscCheck(!A->spptr, PETSC_COMM_SELF, PETSC_ERR_PLIB, "A->spptr is supposed to be empty");
1564 PetscCall(MatSetSizes(A, m, n, m, n));
1565 PetscCall(MatSetType(A, MATSEQAIJKOKKOS));
1567 /* Set up data structures of A as a MATSEQAIJ */
1568 PetscCall(MatSeqAIJSetPreallocation_SeqAIJ(A, MAT_SKIP_ALLOCATION, NULL));
1569 aseq = (Mat_SeqAIJ *)A->data;
1588 akok->nonzerostate = A->nonzerostate;
1589 A->spptr = akok; /* Set A->spptr before MatAssembly so that A->spptr won't be allocated again there */
1590 PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
1591 PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
1595 PETSC_INTERN PetscErrorCode MatSeqAIJKokkosGetKokkosCsrMatrix(Mat A, KokkosCsrMatrix *csr)
1598 PetscCall(MatSeqAIJKokkosSyncDevice(A));
1599 *csr = static_cast<Mat_SeqAIJKokkos *>(A->spptr)->csrmat;
1603 PETSC_INTERN PetscErrorCode MatCreateSeqAIJKokkosWithKokkosCsrMatrix(MPI_Comm comm, KokkosCsrMatrix csr, Mat *A)
1609 PetscCall(MatCreate(comm, A));
1610 PetscCall(MatSetSeqAIJKokkosWithCSRMatrix(*A, akok));
1618 PETSC_INTERN PetscErrorCode MatCreateSeqAIJKokkosWithCSRMatrix(MPI_Comm comm, Mat_SeqAIJKokkos *akok, Mat *A)
1621 PetscCall(MatCreate(comm, A));
1622 PetscCall(MatSetSeqAIJKokkosWithCSRMatrix(*A, akok));
1641 . A - the matrix
1661 PetscErrorCode MatCreateSeqAIJKokkos(MPI_Comm comm, PetscInt m, PetscInt n, PetscInt nz, const PetscInt nnz[], Mat *A)
1665 PetscCall(MatCreate(comm, A));
1666 PetscCall(MatSetSizes(*A, m, n, m, n));
1667 PetscCall(MatSetType(*A, MATSEQAIJKOKKOS));
1668 PetscCall(MatSeqAIJSetPreallocation_SeqAIJ(*A, nz, (PetscInt *)nnz));
1676 static PetscErrorCode MatSeqAIJKokkosSolveCheck(Mat A)
1678 Mat_SeqAIJKokkosTriFactors *factors = (Mat_SeqAIJKokkosTriFactors *)A->spptr;
1691 static PetscErrorCode MatSeqAIJKokkosTransposeSolveCheck(Mat A)
1693 const PetscInt n = A->rmap->n;
1694 Mat_SeqAIJKokkosTriFactors *factors = (Mat_SeqAIJKokkosTriFactors *)A->spptr;
1778 // Solve Ax = b, with RAR = U^T D U, where R is the row (and col) permutation matrix on A.
1780 static PetscErrorCode MatSolve_SeqAIJKokkos_Cholesky(Mat A, Vec bb, Vec xx)
1783 Mat_SeqAIJKokkosTriFactors *factors = (Mat_SeqAIJKokkosTriFactors *)A->spptr;
1784 PetscInt m = A->rmap->n;
1794 PetscCall(MatSeqAIJKokkosSolveCheck(A)); // for UX = T
1795 PetscCall(MatSeqAIJKokkosTransposeSolveCheck(A)); // for U^T Y = B
1831 // Solve Ax = b, with RAC = LU, where R and C are row and col permutation matrices on A respectively.
1834 static PetscErrorCode MatSolve_SeqAIJKokkos_LU(Mat A, Vec bb, Vec xx)
1837 Mat_SeqAIJKokkosTriFactors *factors = (Mat_SeqAIJKokkosTriFactors *)A->spptr;
1838 PetscInt m = A->rmap->n;
1849 PetscCall(MatSeqAIJKokkosSolveCheck(A));
1881 // Solve A^T x = b, with RAC = LU, where R and C are row and col permutation matrices on A respectively.
1884 // A = R^-1 L U C^-1, so A^T = C^-T U^T L^T R^-T. But since C^- = C^T, R^- = R^T, we have A^T = C U^T L^T R.
1885 static PetscErrorCode MatSolveTranspose_SeqAIJKokkos_LU(Mat A, Vec bb, Vec xx)
1888 Mat_SeqAIJKokkosTriFactors *factors = (Mat_SeqAIJKokkosTriFactors *)A->spptr;
1889 PetscInt m = A->rmap->n;
1900 PetscCall(MatSeqAIJKokkosTransposeSolveCheck(A)); // Update L^T, U^T if needed, and do sptrsv symbolic for L^T, U^T
1932 static PetscErrorCode MatLUFactorNumeric_SeqAIJKokkos(Mat B, Mat A, const MatFactorInfo *info)
1935 PetscCall(MatSeqAIJKokkosSyncHost(A));
1936 PetscCall(MatLUFactorNumeric_SeqAIJ(B, A, info));
2052 static PetscErrorCode MatILUFactorNumeric_SeqAIJKokkos_ILU0(Mat B, Mat A, const MatFactorInfo *info)
2054 Mat_SeqAIJKokkos *aijkok = (Mat_SeqAIJKokkos *)A->spptr;
2060 PetscCheck(!info->factoronhost, PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "MatFactorInfo.factoronhost should be false");
2061 PetscCall(MatSeqAIJKokkosSyncDevice(A));
2085 static PetscErrorCode MatILUFactorSymbolic_SeqAIJKokkos_ILU0(Mat B, Mat A, IS, IS, const MatFactorInfo *info)
2091 PetscInt nnzA = ((Mat_SeqAIJ *)A->data)->nz, nnzL, nnzU;
2092 PetscInt n = A->rmap->n;
2095 PetscCheck(!info->factoronhost, PetscObjectComm((PetscObject)A), PETSC_ERR_PLIB, "MatFactorInfo's factoronhost should be false as we are doing it on device right now");
2096 PetscCall(MatSeqAIJKokkosSyncDevice(A));
2109 aijkok = (Mat_SeqAIJKokkos *)A->spptr;
2144 static PetscErrorCode MatLUFactorSymbolic_SeqAIJKokkos(Mat B, Mat A, IS isrow, IS iscol, const MatFactorInfo *info)
2147 PetscCall(MatLUFactorSymbolic_SeqAIJ(B, A, isrow, iscol, info));
2154 static PetscErrorCode MatILUFactorSymbolic_SeqAIJKokkos(Mat B, Mat A, IS isrow, IS iscol, const MatFactorInfo *info)
2168 PetscCall(MatILUFactorSymbolic_SeqAIJKokkos_ILU0(B, A, isrow, iscol, info));
2170 PetscCall(MatILUFactorSymbolic_SeqAIJ(B, A, isrow, iscol, info)); // otherwise, use PETSc's ILU on host
2176 static PetscErrorCode MatCholeskyFactorNumeric_SeqAIJKokkos(Mat B, Mat A, const MatFactorInfo *info)
2179 PetscCall(MatSeqAIJKokkosSyncHost(A));
2180 PetscCall(MatCholeskyFactorNumeric_SeqAIJ(B, A, info));
2257 static PetscErrorCode MatICCFactorSymbolic_SeqAIJKokkos(Mat B, Mat A, IS perm, const MatFactorInfo *info)
2266 PetscCall(MatICCFactorSymbolic_SeqAIJ(B, A, perm, info));
2277 static PetscErrorCode MatCholeskyFactorSymbolic_SeqAIJKokkos(Mat B, Mat A, IS perm, const MatFactorInfo *info)
2286 PetscCall(MatCholeskyFactorSymbolic_SeqAIJ(B, A, perm, info)); // it sets B's two ISes ((Mat_SeqAIJ*)B->data)->{row, col} to perm
2298 static PetscErrorCode MatFactorGetSolverType_SeqAIJKokkos_Kokkos(Mat A, MatSolverType *type)
2306 MATSOLVERKOKKOS = "Kokkos" - A matrix solver type providing triangular solvers for sequential matrices
2313 PETSC_EXTERN PetscErrorCode MatGetFactor_SeqAIJKokkos_Kokkos(Mat A, MatFactorType ftype, Mat *B) /* MatGetFactor_<MatType>_<MatSolverType> */
2315 PetscInt n = A->rmap->n;
2319 PetscCall(PetscObjectGetComm((PetscObject)A, &comm));
2322 PetscCall(MatSetBlockSizesFromMats(*B, A, A));