Lines Matching refs:factors
147 /* We do not expect one needs factors on host */
1678 Mat_SeqAIJKokkosTriFactors *factors = (Mat_SeqAIJKokkosTriFactors *)A->spptr;
1679 const PetscBool has_lower = factors->iL_d.extent(0) ? PETSC_TRUE : PETSC_FALSE; // false with Choleksy
1680 const PetscBool has_upper = factors->iU_d.extent(0) ? PETSC_TRUE : PETSC_FALSE; // true with LU and Choleksy
1683 if (!factors->sptrsv_symbolic_completed) { // If sptrsv_symbolic was not called yet
1684 if (has_upper) PetscCallCXX(sptrsv_symbolic(&factors->khU, factors->iU_d, factors->jU_d, factors->aU_d));
1685 if (has_lower) PetscCallCXX(sptrsv_symbolic(&factors->khL, factors->iL_d, factors->jL_d, factors->aL_d));
1686 factors->sptrsv_symbolic_completed = PETSC_TRUE;
1694 Mat_SeqAIJKokkosTriFactors *factors = (Mat_SeqAIJKokkosTriFactors *)A->spptr;
1695 const PetscBool has_lower = factors->iL_d.extent(0) ? PETSC_TRUE : PETSC_FALSE; // false with Choleksy
1696 const PetscBool has_upper = factors->iU_d.extent(0) ? PETSC_TRUE : PETSC_FALSE; // true with LU or Choleksy
1699 if (!factors->transpose_updated) {
1701 if (!factors->iUt_d.extent(0)) { // Allocate Ut on device if not yet
1702 factors->iUt_d = MatRowMapKokkosView("factors->iUt_d", n + 1); // KK requires this view to be initialized to 0 to call transpose_matrix
1703 factors->jUt_d = MatColIdxKokkosView(NoInit("factors->jUt_d"), factors->jU_d.extent(0));
1704 factors->aUt_d = MatScalarKokkosView(NoInit("factors->aUt_d"), factors->aU_d.extent(0));
1707 if (factors->iU_h.extent(0)) { // If U is on host (factorization was done on host), we also compute the transpose on host
1708 if (!factors->U) {
1711 PetscCall(MatCreateSeqAIJWithArrays(PETSC_COMM_SELF, n, n, factors->iU_h.data(), factors->jU_h.data(), factors->aU_h.data(), &factors->U));
1712 PetscCall(MatTranspose(factors->U, MAT_INITIAL_MATRIX, &factors->Ut));
1714 seq = static_cast<Mat_SeqAIJ *>(factors->Ut->data);
1715 factors->iUt_h = MatRowMapKokkosViewHost(seq->i, n + 1);
1716 factors->jUt_h = MatColIdxKokkosViewHost(seq->j, seq->nz);
1717 factors->aUt_h = MatScalarKokkosViewHost(seq->a, seq->nz);
1719 PetscCall(MatTranspose(factors->U, MAT_REUSE_MATRIX, &factors->Ut)); // Matrix Ut' data is aliased with {i, j, a}Ut_h
1722 PetscCallCXX(Kokkos::deep_copy(factors->iUt_d, factors->iUt_h));
1723 PetscCallCXX(Kokkos::deep_copy(factors->jUt_d, factors->jUt_h));
1724 PetscCallCXX(Kokkos::deep_copy(factors->aUt_d, factors->aUt_h));
1727 PetscCallCXX(transpose_matrix<ConstMatRowMapKokkosView, ConstMatColIdxKokkosView, ConstMatScalarKokkosView, MatRowMapKokkosView, MatColIdxKokkosView, MatScalarKokkosView, MatRowMapKokkosView, DefaultExecutionSpace>(n, n, factors->iU_d,
1728 factors->jU_d, factors->aU_d,
1729 factors->iUt_d, factors->jUt_d,
1730 factors->aUt_d));
1731 PetscCallCXX(sort_crs_matrix<DefaultExecutionSpace, MatRowMapKokkosView, MatColIdxKokkosView, MatScalarKokkosView>(factors->iUt_d, factors->jUt_d, factors->aUt_d));
1733 PetscCallCXX(sptrsv_symbolic(&factors->khUt, factors->iUt_d, factors->jUt_d, factors->aUt_d));
1738 if (!factors->iLt_d.extent(0)) { // Allocate Lt on device if not yet
1739 factors->iLt_d = MatRowMapKokkosView("factors->iLt_d", n + 1); // KK requires this view to be initialized to 0 to call transpose_matrix
1740 factors->jLt_d = MatColIdxKokkosView(NoInit("factors->jLt_d"), factors->jL_d.extent(0));
1741 factors->aLt_d = MatScalarKokkosView(NoInit("factors->aLt_d"), factors->aL_d.extent(0));
1744 if (factors->iL_h.extent(0)) { // If L is on host, we also compute the transpose on host
1745 if (!factors->L) {
1748 PetscCall(MatCreateSeqAIJWithArrays(PETSC_COMM_SELF, n, n, factors->iL_h.data(), factors->jL_h.data(), factors->aL_h.data(), &factors->L));
1749 PetscCall(MatTranspose(factors->L, MAT_INITIAL_MATRIX, &factors->Lt));
1751 seq = static_cast<Mat_SeqAIJ *>(factors->Lt->data);
1752 factors->iLt_h = MatRowMapKokkosViewHost(seq->i, n + 1);
1753 factors->jLt_h = MatColIdxKokkosViewHost(seq->j, seq->nz);
1754 factors->aLt_h = MatScalarKokkosViewHost(seq->a, seq->nz);
1756 PetscCall(MatTranspose(factors->L, MAT_REUSE_MATRIX, &factors->Lt)); // Matrix Lt' data is aliased with {i, j, a}Lt_h
1759 PetscCallCXX(Kokkos::deep_copy(factors->iLt_d, factors->iLt_h));
1760 PetscCallCXX(Kokkos::deep_copy(factors->jLt_d, factors->jLt_h));
1761 PetscCallCXX(Kokkos::deep_copy(factors->aLt_d, factors->aLt_h));
1764 PetscCallCXX(transpose_matrix<ConstMatRowMapKokkosView, ConstMatColIdxKokkosView, ConstMatScalarKokkosView, MatRowMapKokkosView, MatColIdxKokkosView, MatScalarKokkosView, MatRowMapKokkosView, DefaultExecutionSpace>(n, n, factors->iL_d,
1765 factors->jL_d, factors->aL_d,
1766 factors->iLt_d, factors->jLt_d,
1767 factors->aLt_d));
1768 PetscCallCXX(sort_crs_matrix<DefaultExecutionSpace, MatRowMapKokkosView, MatColIdxKokkosView, MatScalarKokkosView>(factors->iLt_d, factors->jLt_d, factors->aLt_d));
1770 PetscCallCXX(sptrsv_symbolic(&factors->khLt, factors->iLt_d, factors->jLt_d, factors->aLt_d));
1773 factors->transpose_updated = PETSC_TRUE;
1779 // R is represented by rowperm in factors. If R is identity (i.e, no reordering), then rowperm is empty.
1783 Mat_SeqAIJKokkosTriFactors *factors = (Mat_SeqAIJKokkosTriFactors *)A->spptr;
1785 PetscScalarKokkosView D = factors->D_d;
1789 PetscIntKokkosView &rowperm = factors->rowperm;
1802 Y = factors->workVector;
1804 B = factors->workVector;
1808 PetscCallCXX(sptrsv_solve(exec, &factors->khUt, factors->iUt_d, factors->jUt_d, factors->aUt_d, B, Y));
1818 X = factors->workVector; // B is not needed anymore
1820 PetscCallCXX(sptrsv_solve(exec, &factors->khU, factors->iU_d, factors->jU_d, factors->aU_d, Y, X));
1832 // R and C are represented by rowperm and colperm in factors.
1837 Mat_SeqAIJKokkosTriFactors *factors = (Mat_SeqAIJKokkosTriFactors *)A->spptr;
1842 PetscIntKokkosView &rowperm = factors->rowperm;
1843 PetscIntKokkosView &colperm = factors->colperm;
1856 Y = factors->workVector;
1858 B = factors->workVector;
1862 PetscCallCXX(sptrsv_solve(exec, &factors->khL, factors->iL_d, factors->jL_d, factors->aL_d, B, Y));
1868 X = factors->workVector;
1870 PetscCallCXX(sptrsv_solve(exec, &factors->khU, factors->iU_d, factors->jU_d, factors->aU_d, Y, X));
1882 // R and C are represented by rowperm and colperm in factors.
1888 Mat_SeqAIJKokkosTriFactors *factors = (Mat_SeqAIJKokkosTriFactors *)A->spptr;
1893 PetscIntKokkosView &rowperm = factors->rowperm;
1894 PetscIntKokkosView &colperm = factors->colperm;
1907 Y = factors->workVector;
1909 B = factors->workVector;
1913 PetscCallCXX(sptrsv_solve(exec, &factors->khUt, factors->iUt_d, factors->jUt_d, factors->aUt_d, B, Y));
1919 X = factors->workVector;
1921 PetscCallCXX(sptrsv_solve(exec, &factors->khLt, factors->iLt_d, factors->jLt_d, factors->aLt_d, Y, X));
1939 Mat_SeqAIJKokkosTriFactors *factors = (Mat_SeqAIJKokkosTriFactors *)B->spptr;
1945 if (factors->iL_h.extent(0) == 0) { // Allocate memory and copy the L, U structure for the first time
1947 factors->iL_h = MatRowMapKokkosViewHost(NoInit("iL_h"), m + 1);
1948 factors->jL_h = MatColIdxKokkosViewHost(NoInit("jL_h"), (Bi[m] - Bi[0]) + m); // + the diagonal entries
1949 factors->aL_h = MatScalarKokkosViewHost(NoInit("aL_h"), (Bi[m] - Bi[0]) + m);
1950 factors->iU_h = MatRowMapKokkosViewHost(NoInit("iU_h"), m + 1);
1951 factors->jU_h = MatColIdxKokkosViewHost(NoInit("jU_h"), (Bdiag[0] - Bdiag[m]));
1952 factors->aU_h = MatScalarKokkosViewHost(NoInit("aU_h"), (Bdiag[0] - Bdiag[m]));
1954 PetscInt *Li = factors->iL_h.data();
1955 PetscInt *Lj = factors->jL_h.data();
1956 PetscInt *Ui = factors->iU_h.data();
1957 PetscInt *Uj = factors->jU_h.data();
1974 factors->iL_d = Kokkos::create_mirror_view_and_copy(DefaultMemorySpace(), factors->iL_h);
1975 factors->jL_d = Kokkos::create_mirror_view_and_copy(DefaultMemorySpace(), factors->jL_h);
1976 factors->iU_d = Kokkos::create_mirror_view_and_copy(DefaultMemorySpace(), factors->iU_h);
1977 factors->jU_d = Kokkos::create_mirror_view_and_copy(DefaultMemorySpace(), factors->jU_h);
1978 factors->aL_d = Kokkos::create_mirror_view(DefaultMemorySpace(), factors->aL_h);
1979 factors->aU_d = Kokkos::create_mirror_view(DefaultMemorySpace(), factors->aU_h);
1989 factors->rowperm = PetscIntKokkosView(NoInit("rowperm"), m);
1990 PetscCallCXX(Kokkos::deep_copy(factors->rowperm, PetscIntKokkosViewHost(const_cast<PetscInt *>(ip), m)));
2002 factors->colperm = PetscIntKokkosView(NoInit("colperm"), n);
2003 PetscCallCXX(Kokkos::deep_copy(factors->colperm, PetscIntKokkosViewHost(const_cast<PetscInt *>(ip), n)));
2014 factors->khL.create_sptrsv_handle(sptrsv_alg, m, true /* L is lower tri */);
2015 factors->khU.create_sptrsv_handle(sptrsv_alg, m, false /* U is not lower tri */);
2016 factors->khLt.create_sptrsv_handle(sptrsv_alg, m, false /* L^T is not lower tri */);
2017 factors->khUt.create_sptrsv_handle(sptrsv_alg, m, true /* U^T is lower tri */);
2024 const PetscInt *Li = factors->iL_h.data();
2025 const PetscInt *Ui = factors->iU_h.data();
2027 PetscScalar *La = factors->aL_h.data();
2028 PetscScalar *Ua = factors->aU_h.data();
2037 PetscCallCXX(Kokkos::deep_copy(factors->aL_d, factors->aL_h));
2038 PetscCallCXX(Kokkos::deep_copy(factors->aU_d, factors->aU_h));
2039 // Once the factors' values have changed, we need to update their transpose and redo sptrsv symbolic
2040 factors->transpose_updated = PETSC_FALSE;
2041 factors->sptrsv_symbolic_completed = PETSC_FALSE;
2055 Mat_SeqAIJKokkosTriFactors *factors = (Mat_SeqAIJKokkosTriFactors *)B->spptr;
2067 PetscCallCXX(spiluk_numeric(&factors->kh, fill_lev, i_d, j_d, a_d, factors->iL_d, factors->jL_d, factors->aL_d, factors->iU_d, factors->jU_d, factors->aU_d));
2076 /* Once the factors' value changed, we need to update their transpose and sptrsv handle */
2077 factors->transpose_updated = PETSC_FALSE;
2078 factors->sptrsv_symbolic_completed = PETSC_FALSE;
2089 Mat_SeqAIJKokkosTriFactors *factors = (Mat_SeqAIJKokkosTriFactors *)B->spptr;
2100 factors->kh.create_spiluk_handle(SPILUKAlgorithm::SEQLVLSCHD_TP1, n, nnzL, nnzU);
2102 auto spiluk_handle = factors->kh.get_spiluk_handle();
2104 Kokkos::realloc(factors->iL_d, n + 1); /* Free old arrays and realloc */
2105 Kokkos::realloc(factors->jL_d, spiluk_handle->get_nnzL());
2106 Kokkos::realloc(factors->iU_d, n + 1);
2107 Kokkos::realloc(factors->jU_d, spiluk_handle->get_nnzU());
2112 PetscCallCXX(spiluk_symbolic(&factors->kh, fill_lev, i_d, j_d, factors->iL_d, factors->jL_d, factors->iU_d, factors->jU_d));
2115 Kokkos::resize(factors->jL_d, spiluk_handle->get_nnzL()); /* Shrink or expand, and retain old value */
2116 Kokkos::resize(factors->jU_d, spiluk_handle->get_nnzU());
2117 Kokkos::realloc(factors->aL_d, spiluk_handle->get_nnzL()); /* No need to retain old value */
2118 Kokkos::realloc(factors->aU_d, spiluk_handle->get_nnzU());
2128 factors->khL.create_sptrsv_handle(sptrsv_alg, n, true /* L is lower tri */);
2129 factors->khU.create_sptrsv_handle(sptrsv_alg, n, false /* U is not lower tri */);
2130 factors->khLt.create_sptrsv_handle(sptrsv_alg, n, false /* L^T is not lower tri */);
2131 factors->khUt.create_sptrsv_handle(sptrsv_alg, n, true /* U^T is lower tri */);
2183 Mat_SeqAIJKokkosTriFactors *factors = (Mat_SeqAIJKokkosTriFactors *)B->spptr;
2189 if (factors->iU_h.extent(0) == 0) { // First time of numeric factorization
2191 factors->iU_h = PetscIntKokkosViewHost(const_cast<PetscInt *>(Bi), m + 1); // wrap Bi as iU_h
2192 factors->jU_h = MatColIdxKokkosViewHost(NoInit("jU_h"), Bi[m]);
2193 factors->aU_h = MatScalarKokkosViewHost(NoInit("aU_h"), Bi[m]);
2194 factors->D_h = MatScalarKokkosViewHost(NoInit("D_h"), m);
2195 factors->aU_d = Kokkos::create_mirror_view(DefaultMemorySpace(), factors->aU_h);
2196 factors->D_d = Kokkos::create_mirror_view(DefaultMemorySpace(), factors->D_h);
2199 PetscInt *Uj = factors->jU_h.data();
2207 PetscCallCXX(factors->iU_d = Kokkos::create_mirror_view_and_copy(DefaultMemorySpace(), factors->iU_h));
2208 PetscCallCXX(factors->jU_d = Kokkos::create_mirror_view_and_copy(DefaultMemorySpace(), factors->jU_h));
2218 PetscCallCXX(factors->rowperm = PetscIntKokkosView(NoInit("rowperm"), m));
2219 PetscCallCXX(Kokkos::deep_copy(factors->rowperm, PetscIntKokkosViewHost(const_cast<PetscInt *>(ip), m)));
2230 factors->khU.create_sptrsv_handle(sptrsv_alg, m, false /* U is not lower tri */);
2231 factors->khUt.create_sptrsv_handle(sptrsv_alg, m, true /* U^T is lower tri */);
2238 PetscScalar *Ua = factors->aU_h.data();
2239 PetscScalar *D = factors->D_h.data();
2245 PetscCallCXX(Kokkos::deep_copy(factors->aU_d, factors->aU_h));
2246 PetscCallCXX(Kokkos::deep_copy(factors->D_d, factors->D_h));
2248 factors->sptrsv_symbolic_completed = PETSC_FALSE; // When numeric value changed, we must do these again
2249 factors->transpose_updated = PETSC_FALSE;