1 #include <petscmat.h> 2 3 /*@C 4 MatCreateDenseFromVecType - Create a matrix that matches the type of a Vec. 5 6 Collective 7 8 Input Parameters: 9 + comm - the communicator 10 . vtype - the vector type 11 . m - number of local rows (or `PETSC_DECIDE` to have calculated if `M` is given) 12 . n - number of local columns (or `PETSC_DECIDE` to have calculated if `N` is given) 13 . M - number of global rows (or `PETSC_DECIDE` to have calculated if `m` is given) 14 . N - number of global columns (or `PETSC_DECIDE` to have calculated if `n` is given) 15 . lda - optional leading dimension. Pass any non-positive number to use the default. 16 - data - optional location of matrix data, which should have the same memory type as the vector. Pass `NULL` to have PETSc take care of matrix memory allocation. 17 18 Output Parameter: 19 . A - the dense matrix 20 21 Level: advanced 22 23 .seealso: [](ch_matrices), `Mat`, `MatCreateDense()', `MatCreateDenseCUDA()`, `MatCreateDenseHIP()`, `PetscMemType` 24 @*/ 25 PetscErrorCode MatCreateDenseFromVecType(MPI_Comm comm, VecType vtype, PetscInt m, PetscInt n, PetscInt M, PetscInt N, PetscInt lda, PetscScalar *data, Mat *A) 26 { 27 VecType root_type = VECSTANDARD; 28 PetscBool isstd, iscuda, iship, iskokkos; 29 30 PetscFunctionBegin; 31 PetscCall(PetscStrcmpAny(vtype, &isstd, VECSTANDARD, VECMPI, VECSEQ, "")); 32 PetscCall(PetscStrcmpAny(vtype, &iscuda, VECCUDA, VECMPICUDA, VECSEQCUDA, "")); 33 PetscCall(PetscStrcmpAny(vtype, &iship, VECHIP, VECMPIHIP, VECSEQHIP, "")); 34 PetscCall(PetscStrcmpAny(vtype, &iskokkos, VECKOKKOS, VECMPIKOKKOS, VECSEQKOKKOS, "")); 35 PetscCheck(isstd || iscuda || iship || iskokkos, comm, PETSC_ERR_SUP, "Not for type %s\n", vtype); 36 if (iscuda) root_type = VECCUDA; 37 else if (iship) root_type = VECHIP; 38 else if (iskokkos) { 39 /* We support only one type of kokkos device */ 40 PetscCheck(!PetscDefined(HAVE_MACRO_KOKKOS_ENABLE_SYCL), comm, PETSC_ERR_SUP, "Not for sycl backend"); 41 if (PetscDefined(HAVE_MACRO_KOKKOS_ENABLE_CUDA)) iscuda = PETSC_TRUE; 42 else if (PetscDefined(HAVE_MACRO_KOKKOS_ENABLE_HIP)) iship = PETSC_TRUE; 43 else isstd = PETSC_TRUE; 44 root_type = VECKOKKOS; 45 } 46 PetscCall(MatCreate(comm, A)); 47 PetscCall(MatSetSizes(*A, m, n, M, N)); 48 if (isstd) { 49 PetscCall(MatSetType(*A, MATDENSE)); 50 if (lda > 0) PetscCall(MatDenseSetLDA(*A, lda)); 51 PetscCall(MatSeqDenseSetPreallocation(*A, data)); 52 PetscCall(MatMPIDenseSetPreallocation(*A, data)); 53 } else if (iscuda) { 54 PetscCheck(PetscDefined(HAVE_CUDA), comm, PETSC_ERR_SUP, "PETSc not compiled with CUDA support"); 55 #if defined(PETSC_HAVE_CUDA) 56 PetscCall(MatSetType(*A, MATDENSECUDA)); 57 if (lda > 0) PetscCall(MatDenseSetLDA(*A, lda)); 58 PetscCall(MatDenseCUDASetPreallocation(*A, data)); 59 #endif 60 } else if (iship) { 61 PetscCheck(PetscDefined(HAVE_HIP), comm, PETSC_ERR_SUP, "PETSc not compiled with HIP support"); 62 #if defined(PETSC_HAVE_HIP) 63 PetscCall(MatSetType(*A, MATDENSEHIP)); 64 if (lda > 0) PetscCall(MatDenseSetLDA(*A, lda)); 65 PetscCall(MatDenseHIPSetPreallocation(*A, data)); 66 #endif 67 } 68 PetscCall(MatSetVecType(*A, root_type)); 69 PetscFunctionReturn(PETSC_SUCCESS); 70 } 71