diff --git a/cmake/CTestScript.cmake b/cmake/CTestScript.cmake
index 0b018ff8a58..5ba277eea80 100644
--- a/cmake/CTestScript.cmake
+++ b/cmake/CTestScript.cmake
@@ -161,8 +161,9 @@ ctest_submit(PARTS Start)
 
 if(CTEST_MEMORYCHECK_TYPE STREQUAL "CudaMemcheck")
     # generate line number information for CUDA
+    # we disable bfloat16 because our job is executed on the gpu does not support bfloat16 natively
     set(GINKGO_CONFIGURE_OPTIONS
-        "-DGINKGO_DEVEL_TOOLS=OFF;-DGINKGO_BUILD_REFERENCE=ON;-DGINKGO_BUILD_OMP=OFF;-DGINKGO_BUILD_CUDA=ON;-DGINKGO_BUILD_HIP=OFF;-DGINKGO_BUILD_SYCL=OFF;-DCMAKE_BUILD_TYPE=${CTEST_BUILD_CONFIGURATION};-DCMAKE_CUDA_FLAGS=-lineinfo"
+        "-DGINKGO_ENABLE_BFLOAT16=OFF;-DGINKGO_DEVEL_TOOLS=OFF;-DGINKGO_BUILD_REFERENCE=ON;-DGINKGO_BUILD_OMP=OFF;-DGINKGO_BUILD_CUDA=ON;-DGINKGO_BUILD_HIP=OFF;-DGINKGO_BUILD_SYCL=OFF;-DCMAKE_BUILD_TYPE=${CTEST_BUILD_CONFIGURATION};-DCMAKE_CUDA_FLAGS=-lineinfo"
     )
 elseif(
     (
diff --git a/cmake/generate_ginkgo_hpp.cmake b/cmake/generate_ginkgo_hpp.cmake
index 7121e0cb87d..2d749678779 100644
--- a/cmake/generate_ginkgo_hpp.cmake
+++ b/cmake/generate_ginkgo_hpp.cmake
@@ -13,6 +13,7 @@ function(ginkgo_generate_ginkgo_hpp)
             (file MATCHES "^ginkgo/extensions/.*$")
             OR (file MATCHES "^ginkgo/core/stop/residual_norm_reduction.hpp$")
             OR (file MATCHES "^ginkgo/core/solver/.*_trs.hpp$")
+            OR (file MATCHES "^ginkgo/core/preconditioner/utils.hpp$")
         )
             continue()
         endif()
diff --git a/core/config/preconditioner_ic_config.cpp b/core/config/preconditioner_ic_config.cpp
index d90792c2609..ff95535f597 100644
--- a/core/config/preconditioner_ic_config.cpp
+++ b/core/config/preconditioner_ic_config.cpp
@@ -27,9 +27,9 @@ parse<gko::config::LinOpFactoryType::Ic>(const gko::config::pnode& config,
     if (config.get("l_solver_type_or_value_type")) {
         GKO_INVALID_STATE(
             "preconditioner::Ic only allows value_type from "
-            "l_solver_type_or_value_type. To avoid type confusion between "
-            "these types and value_type, l_solver_type_or_value_type uses "
-            "the value_type directly.");
+            "l_solver_type_or_value_type. Please use value_type key to set the "
+            "value type used by the preconditioner and the l_lover  key to set "
+            "the solvers used for the lower triangular systems.");
     }
     return gko::config::dispatch<gko::LinOpFactory, gko::preconditioner::Ic>(
         config, context, updated,
diff --git a/core/config/preconditioner_ilu_config.cpp b/core/config/preconditioner_ilu_config.cpp
index 3714d356e71..0eb11d094e2 100644
--- a/core/config/preconditioner_ilu_config.cpp
+++ b/core/config/preconditioner_ilu_config.cpp
@@ -1,16 +1,11 @@
-// SPDX-FileCopyrightText: 2017 - 2024 The Ginkgo authors
+// SPDX-FileCopyrightText: 2017 - 2025 The Ginkgo authors
 //
 // SPDX-License-Identifier: BSD-3-Clause
 
 #include <ginkgo/core/base/exception_helpers.hpp>
 #include <ginkgo/core/config/config.hpp>
 #include <ginkgo/core/config/registry.hpp>
-#include <ginkgo/core/preconditioner/ic.hpp>
 #include <ginkgo/core/preconditioner/ilu.hpp>
-#include <ginkgo/core/preconditioner/isai.hpp>
-#include <ginkgo/core/solver/gmres.hpp>
-#include <ginkgo/core/solver/ir.hpp>
-#include <ginkgo/core/solver/triangular.hpp>
 
 #include "core/config/config_helper.hpp"
 #include "core/config/dispatch.hpp"
@@ -22,19 +17,18 @@ namespace gko {
 namespace config {
 
 
-// For Ic and Ilu, we use additional ValueType to help Solver type decision
-template <typename LSolver, typename USolver, bool ReverseApply>
+template <bool ReverseApply>
 class IluSolverHelper {
 public:
     template <typename ValueType, typename IndexType>
     class Configurator {
     public:
-        static typename preconditioner::Ilu<LSolver, USolver, ReverseApply,
+        static typename preconditioner::Ilu<ValueType, ValueType, ReverseApply,
                                             IndexType>::parameters_type
         parse(const pnode& config, const registry& context,
               const type_descriptor& td_for_child)
         {
-            return preconditioner::Ilu<LSolver, USolver, ReverseApply,
+            return preconditioner::Ilu<ValueType, ValueType, ReverseApply,
                                        IndexType>::parse(config, context,
                                                          td_for_child);
         }
@@ -42,103 +36,35 @@ class IluSolverHelper {
 };
 
 
-template <template <typename V> class LSolverBase,
-          template <typename V> class USolverBase, bool ReverseApply>
-class IluHelper1 {
-public:
-    template <typename ValueType, typename IndexType>
-    class Configurator
-        : public IluSolverHelper<
-              LSolverBase<ValueType>, USolverBase<ValueType>,
-              ReverseApply>::template Configurator<ValueType, IndexType> {};
-};
-
-
-template <template <typename V, typename I> class LSolverBase,
-          template <typename V, typename I> class USolverBase,
-          bool ReverseApply>
-class IluHelper2 {
-public:
-    template <typename ValueType, typename IndexType>
-    class Configurator
-        : public IluSolverHelper<
-              LSolverBase<ValueType, IndexType>,
-              USolverBase<ValueType, IndexType>,
-              ReverseApply>::template Configurator<ValueType, IndexType> {};
-};
-
-
 template <>
 deferred_factory_parameter<gko::LinOpFactory> parse<LinOpFactoryType::Ilu>(
     const pnode& config, const registry& context, const type_descriptor& td)
 {
     auto updated = update_type(config, td);
-    auto dispatch_solver = [&](auto reverse_apply)
-        -> deferred_factory_parameter<gko::LinOpFactory> {
-        using ReverseApply = decltype(reverse_apply);
-        // always use symmetric solver for USolverType
-        if (config.get("u_solver_type")) {
-            GKO_INVALID_STATE(
-                "preconditioner::Ilu only allows l_solver_type. The "
-                "u_solver_type automatically uses the transposed type of "
-                "l_solver_type.");
-        }
-        std::string str("solver::LowerTrs");
-        if (auto& obj = config.get("l_solver_type")) {
-            str = obj.get_string();
-        }
-        if (str == "solver::LowerTrs") {
-            return dispatch<
-                gko::LinOpFactory,
-                IluHelper2<solver::LowerTrs, solver::UpperTrs,
-                           ReverseApply::value>::template Configurator>(
-                config, context, updated,
-                make_type_selector(updated.get_value_typestr(),
-                                   value_type_list()),
-                make_type_selector(updated.get_index_typestr(),
-                                   index_type_list()));
-        } else if (str == "solver::Ir") {
-            return dispatch<
-                gko::LinOpFactory,
-                IluHelper1<solver::Ir, solver::Ir,
-                           ReverseApply::value>::template Configurator>(
-                config, context, updated,
-                make_type_selector(updated.get_value_typestr(),
-                                   value_type_list()),
-                make_type_selector(updated.get_index_typestr(),
-                                   index_type_list()));
-        } else if (str == "preconditioner::LowerIsai") {
-            return dispatch<
-                gko::LinOpFactory,
-                IluHelper2<preconditioner::LowerIsai, preconditioner::UpperIsai,
-                           ReverseApply::value>::template Configurator>(
-                config, context, updated,
-                make_type_selector(updated.get_value_typestr(),
-                                   value_type_list()),
-                make_type_selector(updated.get_index_typestr(),
-                                   index_type_list()));
-        } else if (str == "solver::Gmres") {
-            return dispatch<
-                gko::LinOpFactory,
-                IluHelper1<solver::Gmres, solver::Gmres,
-                           ReverseApply::value>::template Configurator>(
-                config, context, updated,
-                make_type_selector(updated.get_value_typestr(),
-                                   value_type_list()),
-                make_type_selector(updated.get_index_typestr(),
-                                   index_type_list()));
-        } else {
-            GKO_INVALID_CONFIG_VALUE("l_solver_type", str);
-        }
-    };
+    if (config.get("l_solver_type_or_value_type") ||
+        config.get("u_solver_type_or_value_type")) {
+        GKO_INVALID_STATE(
+            "preconditioner::Ilu only allows value_type from "
+            "l_solver_type_or_value_type/u_solver_type_or_value_type. Please "
+            "use value_type key to set the value type used by the "
+            "preconditioner and the l_lover or u_solver key to set the solvers "
+            "used for the lower and upper triangular systems.");
+    }
     bool reverse_apply = false;
     if (auto& obj = config.get("reverse_apply")) {
         reverse_apply = obj.get_boolean();
     }
     if (reverse_apply) {
-        return dispatch_solver(std::true_type{});
+        return dispatch<gko::LinOpFactory, IluSolverHelper<true>::Configurator>(
+            config, context, updated,
+            make_type_selector(updated.get_value_typestr(), value_type_list()),
+            make_type_selector(updated.get_index_typestr(), index_type_list()));
     } else {
-        return dispatch_solver(std::false_type{});
+        return dispatch<gko::LinOpFactory,
+                        IluSolverHelper<false>::Configurator>(
+            config, context, updated,
+            make_type_selector(updated.get_value_typestr(), value_type_list()),
+            make_type_selector(updated.get_index_typestr(), index_type_list()));
     }
 }
 
diff --git a/core/preconditioner/ilu.cpp b/core/preconditioner/ilu.cpp
index d6f49e49588..1c495b6ccaf 100644
--- a/core/preconditioner/ilu.cpp
+++ b/core/preconditioner/ilu.cpp
@@ -1,17 +1,14 @@
-// SPDX-FileCopyrightText: 2017 - 2024 The Ginkgo authors
+// SPDX-FileCopyrightText: 2017 - 2025 The Ginkgo authors
 //
 // SPDX-License-Identifier: BSD-3-Clause
 
 #include "ginkgo/core/preconditioner/ilu.hpp"
 
 #include <ginkgo/core/base/types.hpp>
+#include <ginkgo/core/base/utils.hpp>
 #include <ginkgo/core/config/config.hpp>
 #include <ginkgo/core/config/registry.hpp>
 #include <ginkgo/core/config/type_descriptor.hpp>
-#include <ginkgo/core/preconditioner/isai.hpp>
-#include <ginkgo/core/preconditioner/utils.hpp>
-#include <ginkgo/core/solver/gmres.hpp>
-#include <ginkgo/core/solver/ir.hpp>
 
 #include "core/config/config_helper.hpp"
 #include "core/config/dispatch.hpp"
@@ -22,24 +19,27 @@ namespace preconditioner {
 namespace detail {
 
 
-template <typename Ilu,
-          std::enable_if_t<support_ilu_parse<typename Ilu::l_solver_type,
-                                             typename Ilu::u_solver_type>>*>
+template <typename Ilu, std::enable_if_t<support_ilu_parse<Ilu>>*>
 typename Ilu::parameters_type ilu_parse(
     const config::pnode& config, const config::registry& context,
     const config::type_descriptor& td_for_child)
 {
     auto params = Ilu::build();
-
+    using l_solver_type = typename Ilu::l_solver_type;
+    using u_solver_type = typename Ilu::u_solver_type;
+    static_assert(std::is_same_v<l_solver_type, LinOp>,
+                  "only support ILU parse when l_solver_type is LinOp.");
+    static_assert(std::is_same_v<u_solver_type, LinOp>,
+                  "only support ILU parse when u_solver_type is LinOp.");
     if (auto& obj = config.get("l_solver")) {
         params.with_l_solver(
-            gko::config::parse_or_get_specific_factory<
-                const typename Ilu::l_solver_type>(obj, context, td_for_child));
+            gko::config::parse_or_get_factory<const LinOpFactory>(
+                obj, context, td_for_child));
     }
     if (auto& obj = config.get("u_solver")) {
         params.with_u_solver(
-            gko::config::parse_or_get_specific_factory<
-                const typename Ilu::u_solver_type>(obj, context, td_for_child));
+            gko::config::parse_or_get_factory<const LinOpFactory>(
+                obj, context, td_for_child));
     }
     if (auto& obj = config.get("factorization")) {
         params.with_factorization(
@@ -51,82 +51,33 @@ typename Ilu::parameters_type ilu_parse(
 }
 
 
-#define GKO_DECLARE_TRS_ILU_FALSE_PARSE(ValueType, IndexType)                 \
-    typename Ilu<solver::LowerTrs<ValueType, IndexType>,                      \
-                 solver::UpperTrs<ValueType, IndexType>, false,               \
-                 IndexType>::parameters_type                                  \
-    ilu_parse<Ilu<solver::LowerTrs<ValueType, IndexType>,                     \
-                  solver::UpperTrs<ValueType, IndexType>, false, IndexType>>( \
-        const config::pnode&, const config::registry&,                        \
-        const config::type_descriptor&)
-GKO_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE(GKO_DECLARE_TRS_ILU_FALSE_PARSE);
-
-#define GKO_DECLARE_TRS_ILU_TRUE_PARSE(ValueType, IndexType)                 \
-    typename Ilu<solver::LowerTrs<ValueType, IndexType>,                     \
-                 solver::UpperTrs<ValueType, IndexType>, true,               \
-                 IndexType>::parameters_type                                 \
-    ilu_parse<Ilu<solver::LowerTrs<ValueType, IndexType>,                    \
-                  solver::UpperTrs<ValueType, IndexType>, true, IndexType>>( \
-        const config::pnode&, const config::registry&,                       \
-        const config::type_descriptor&)
-GKO_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE(GKO_DECLARE_TRS_ILU_TRUE_PARSE);
-
-#define GKO_DECLARE_GMRES_ILU_FALSE_PARSE(ValueType, IndexType)              \
-    typename Ilu<solver::Gmres<ValueType>, solver::Gmres<ValueType>, false,  \
-                 IndexType>::parameters_type                                 \
-    ilu_parse<Ilu<solver::Gmres<ValueType>, solver::Gmres<ValueType>, false, \
-                  IndexType>>(const config::pnode&, const config::registry&, \
-                              const config::type_descriptor&)
-GKO_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE(
-    GKO_DECLARE_GMRES_ILU_FALSE_PARSE);
-
-#define GKO_DECLARE_GMRES_ILU_TRUE_PARSE(ValueType, IndexType)               \
-    typename Ilu<solver::Gmres<ValueType>, solver::Gmres<ValueType>, true,   \
-                 IndexType>::parameters_type                                 \
-    ilu_parse<Ilu<solver::Gmres<ValueType>, solver::Gmres<ValueType>, true,  \
-                  IndexType>>(const config::pnode&, const config::registry&, \
-                              const config::type_descriptor&)
-GKO_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE(GKO_DECLARE_GMRES_ILU_TRUE_PARSE);
-
-#define GKO_DECLARE_IR_ILU_FALSE_PARSE(ValueType, IndexType)                  \
-    typename Ilu<solver::Ir<ValueType>, solver::Ir<ValueType>, false,         \
-                 IndexType>::parameters_type                                  \
-    ilu_parse<                                                                \
-        Ilu<solver::Ir<ValueType>, solver::Ir<ValueType>, false, IndexType>>( \
-        const config::pnode&, const config::registry&,                        \
-        const config::type_descriptor&)
-GKO_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE(GKO_DECLARE_IR_ILU_FALSE_PARSE);
-
-#define GKO_DECLARE_IR_ILU_TRUE_PARSE(ValueType, IndexType)                  \
-    typename Ilu<solver::Ir<ValueType>, solver::Ir<ValueType>, true,         \
-                 IndexType>::parameters_type                                 \
-    ilu_parse<                                                               \
-        Ilu<solver::Ir<ValueType>, solver::Ir<ValueType>, true, IndexType>>( \
-        const config::pnode&, const config::registry&,                       \
+#define GKO_DECLARE_ILU_PARSE_FALSE(ValueType, IndexType)                 \
+    typename Ilu<ValueType, ValueType, false, IndexType>::parameters_type \
+    ilu_parse<Ilu<ValueType, ValueType, false, IndexType>>(               \
+        const config::pnode&, const config::registry&,                    \
         const config::type_descriptor&)
-GKO_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE(GKO_DECLARE_IR_ILU_TRUE_PARSE);
-
-#define GKO_DECLARE_ISAI_ILU_FALSE_PARSE(ValueType, IndexType)         \
-    typename Ilu<LowerIsai<ValueType, IndexType>,                      \
-                 UpperIsai<ValueType, IndexType>, false,               \
-                 IndexType>::parameters_type                           \
-    ilu_parse<Ilu<LowerIsai<ValueType, IndexType>,                     \
-                  UpperIsai<ValueType, IndexType>, false, IndexType>>( \
-        const config::pnode&, const config::registry&,                 \
+#define GKO_DECLARE_ILU_PARSE_TRUE(ValueType, IndexType)                 \
+    typename Ilu<ValueType, ValueType, true, IndexType>::parameters_type \
+    ilu_parse<Ilu<ValueType, ValueType, true, IndexType>>(               \
+        const config::pnode&, const config::registry&,                   \
         const config::type_descriptor&)
-GKO_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE(GKO_DECLARE_ISAI_ILU_FALSE_PARSE);
-
-#define GKO_DECLARE_ISAI_ILU_TRUE_PARSE(ValueType, IndexType)         \
-    typename Ilu<LowerIsai<ValueType, IndexType>,                     \
-                 UpperIsai<ValueType, IndexType>, true,               \
-                 IndexType>::parameters_type                          \
-    ilu_parse<Ilu<LowerIsai<ValueType, IndexType>,                    \
-                  UpperIsai<ValueType, IndexType>, true, IndexType>>( \
-        const config::pnode&, const config::registry&,                \
-        const config::type_descriptor&)
-GKO_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE(GKO_DECLARE_ISAI_ILU_TRUE_PARSE);
+
+GKO_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE(GKO_DECLARE_ILU_PARSE_FALSE);
+GKO_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE(GKO_DECLARE_ILU_PARSE_TRUE);
 
 
 }  // namespace detail
+
+
+// only instantiate the value type variants of ILU, whose solver is LinOp.
+#define GKO_DECLARE_ILU_FALSE(ValueType, IndexType) \
+    class Ilu<ValueType, ValueType, false, IndexType>
+#define GKO_DECLARE_ILU_TRUE(ValueType, IndexType) \
+    class Ilu<ValueType, ValueType, true, IndexType>
+
+GKO_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE(GKO_DECLARE_ILU_FALSE);
+GKO_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE(GKO_DECLARE_ILU_TRUE);
+
+
 }  // namespace preconditioner
 }  // namespace gko
diff --git a/core/test/config/preconditioner.cpp b/core/test/config/preconditioner.cpp
index 153eafd057f..b831fe4e491 100644
--- a/core/test/config/preconditioner.cpp
+++ b/core/test/config/preconditioner.cpp
@@ -106,11 +106,9 @@ struct Ic : PreconditionerConfigTest<::gko::preconditioner::Ic<float, int>,
 };
 
 
-struct Ilu
-    : PreconditionerConfigTest<
-          ::gko::preconditioner::Ilu<DummyIr, DummyIr, true, int>,
-          ::gko::preconditioner::Ilu<gko::solver::LowerTrs<>,
-                                     gko::solver::UpperTrs<>, false, int>> {
+struct Ilu : PreconditionerConfigTest<
+                 ::gko::preconditioner::Ilu<float, float, true, int>,
+                 ::gko::preconditioner::Ilu<double, double, false, int>> {
     static pnode::map_type setup_base()
     {
         return {{"type", pnode{"preconditioner::Ilu"}}};
@@ -119,7 +117,6 @@ struct Ilu
     static void change_template(pnode::map_type& config_map)
     {
         config_map["value_type"] = pnode{"float32"};
-        config_map["l_solver_type"] = pnode{"solver::Ir"};
         config_map["reverse_apply"] = pnode{true};
     }
 
@@ -129,13 +126,13 @@ struct Ilu
     {
         if (from_reg) {
             config_map["l_solver"] = pnode{"l_solver"};
-            param.with_l_solver(detail::registry_accessor::get_data<
-                                typename changed_type::l_solver_type::Factory>(
-                reg, "l_solver"));
+            param.with_l_solver(
+                detail::registry_accessor::get_data<gko::LinOpFactory>(
+                    reg, "l_solver"));
             config_map["u_solver"] = pnode{"u_solver"};
-            param.with_u_solver(detail::registry_accessor::get_data<
-                                typename changed_type::u_solver_type::Factory>(
-                reg, "u_solver"));
+            param.with_u_solver(
+                detail::registry_accessor::get_data<gko::LinOpFactory>(
+                    reg, "u_solver"));
             config_map["factorization"] = pnode{"factorization"};
             param.with_factorization(
                 detail::registry_accessor::get_data<gko::LinOpFactory>(
@@ -143,10 +140,10 @@ struct Ilu
         } else {
             config_map["l_solver"] = pnode{{{"type", pnode{"solver::Ir"}},
                                             {"value_type", pnode{"float32"}}}};
-            param.with_l_solver(changed_type::l_solver_type::build().on(exec));
+            param.with_l_solver(DummyIr::build().on(exec));
             config_map["u_solver"] = pnode{{{"type", pnode{"solver::Ir"}},
                                             {"value_type", pnode{"float32"}}}};
-            param.with_u_solver(changed_type::u_solver_type::build().on(exec));
+            param.with_u_solver(DummyIr::build().on(exec));
             config_map["factorization"] =
                 pnode{{{"type", pnode{"solver::Ir"}},
                        {"value_type", pnode{"float32"}}}};
diff --git a/examples/file-config-solver/config/parilu-gmres.json b/examples/file-config-solver/config/parilu-gmres.json
index 57ebb8fea52..4806f32f864 100644
--- a/examples/file-config-solver/config/parilu-gmres.json
+++ b/examples/file-config-solver/config/parilu-gmres.json
@@ -2,7 +2,6 @@
     "type": "solver::Gmres",
     "preconditioner": {
         "type": "preconditioner::Ilu",
-        "l_solver_type": "solver::LowerTrs",
         "reverse_apply": false,
         "factorization": {
             "type": "factorization::ParIlu"
diff --git a/examples/file-config-solver/file-config-solver.cpp b/examples/file-config-solver/file-config-solver.cpp
index 9f5a2c1f3e1..17be4c2dfbf 100644
--- a/examples/file-config-solver/file-config-solver.cpp
+++ b/examples/file-config-solver/file-config-solver.cpp
@@ -1,4 +1,4 @@
-// SPDX-FileCopyrightText: 2017 - 2024 The Ginkgo authors
+// SPDX-FileCopyrightText: 2017 - 2025 The Ginkgo authors
 //
 // SPDX-License-Identifier: BSD-3-Clause
 
@@ -85,8 +85,9 @@ int main(int argc, char* argv[])
     // cg.json: simple-solver
     // blockjacobi-cg.json: preconditioned-solver
     // ir.json: iterative-refinement
-    // parilu.json: ilu-preconditioned-solver (by using factoization parameter
-    //              directly)
+    // parilu.json: ilu-preconditioned-solver (by using factorization parameter
+    //              directly) Ilu only supports the value_type variant for
+    //              parse.
     // pgm-multigrid-cg.json: multigrid-preconditioned-solver (set
     //                        min_coarse_rows additionally due to this small
     //                        example matrix)
diff --git a/include/ginkgo/core/base/type_traits.hpp b/include/ginkgo/core/base/type_traits.hpp
index e3977aa61b0..5d846eb3cd3 100644
--- a/include/ginkgo/core/base/type_traits.hpp
+++ b/include/ginkgo/core/base/type_traits.hpp
@@ -35,17 +35,24 @@ using factory_type = typename factory_type_impl<Type>::type;
 
 
 // helper for handle the transposed type of concrete type and LinOp
-template <typename Type>
+template <typename Type, typename = void>
 struct transposed_type_impl {
     using type = typename Type::transposed_type;
 };
 
 // It requires LinOp to be complete type
 template <>
-struct transposed_type_impl<LinOp> {
+struct transposed_type_impl<LinOp, void> {
     using type = LinOp;
 };
 
+// return the same type when Type is the precision format.
+// it is used in ILU.
+template <typename Type>
+struct transposed_type_impl<Type, std::enable_if_t<!is_ginkgo_linop<Type>>> {
+    using type = Type;
+};
+
 
 template <typename Type>
 using transposed_type = typename transposed_type_impl<Type>::type;
diff --git a/include/ginkgo/core/preconditioner/ic.hpp b/include/ginkgo/core/preconditioner/ic.hpp
index 6c0b4be4573..3976d01a89c 100644
--- a/include/ginkgo/core/preconditioner/ic.hpp
+++ b/include/ginkgo/core/preconditioner/ic.hpp
@@ -401,24 +401,21 @@ class Ic : public EnableLinOp<Ic<LSolverTypeOrValueType, IndexType>>,
 
         // If no factories are provided, generate default ones
         if (!parameters_.l_solver_factory) {
-            if constexpr (std::is_same_v<l_solver_type, LinOp>) {
-                // linop can not have a default factory generation
-                GKO_NOT_IMPLEMENTED;
-            } else {
-                l_solver_ =
-                    generate_default_solver<l_solver_type>(exec, l_factor);
-                // If comp contains both factors: use the transposed factor to
-                // avoid transposing twice
-                if (comp->get_operators().size() == 2) {
-                    auto lh_factor = comp->get_operators()[1];
-                    GKO_ASSERT_EQUAL_DIMENSIONS(l_factor, lh_factor);
-                    lh_solver_ = as<lh_solver_type>(
-                        as<Transposable>(l_solver_)->conj_transpose());
-                } else {
-                    lh_solver_ = as<lh_solver_type>(
-                        as<Transposable>(l_solver_)->conj_transpose());
-                }
+            // when l_solver_type is LinOp, use LowerTrs as the default one
+            l_solver_ = generate_default_solver<std::conditional_t<
+                std::is_same_v<l_solver_type, LinOp>,
+                solver::LowerTrs<value_type, index_type>, l_solver_type>>(
+                exec, l_factor);
+            // If comp contains both factors: We only check the dimension from
+            // the second factor. However, we still use the l_solver^H not
+            // generate the solver on L^H to preserve the Hermitian property of
+            // this preconditioner. LSolver(L)^H is not always LSolver^H(L^H).
+            if (comp->get_operators().size() == 2) {
+                auto lh_factor = comp->get_operators()[1];
+                GKO_ASSERT_EQUAL_DIMENSIONS(l_factor, lh_factor);
             }
+            lh_solver_ = as<lh_solver_type>(
+                as<Transposable>(l_solver_)->conj_transpose());
         } else {
             l_solver_ = parameters_.l_solver_factory->generate(l_factor);
             lh_solver_ = as<lh_solver_type>(
diff --git a/include/ginkgo/core/preconditioner/ilu.hpp b/include/ginkgo/core/preconditioner/ilu.hpp
index 98aa3ce70c1..c08ef6dd3f1 100644
--- a/include/ginkgo/core/preconditioner/ilu.hpp
+++ b/include/ginkgo/core/preconditioner/ilu.hpp
@@ -1,4 +1,4 @@
-// SPDX-FileCopyrightText: 2017 - 2024 The Ginkgo authors
+// SPDX-FileCopyrightText: 2017 - 2025 The Ginkgo authors
 //
 // SPDX-License-Identifier: BSD-3-Clause
 
@@ -15,14 +15,11 @@
 #include <ginkgo/core/base/exception_helpers.hpp>
 #include <ginkgo/core/base/lin_op.hpp>
 #include <ginkgo/core/base/precision_dispatch.hpp>
+#include <ginkgo/core/base/type_traits.hpp>
 #include <ginkgo/core/config/config.hpp>
 #include <ginkgo/core/config/registry.hpp>
 #include <ginkgo/core/factorization/par_ilu.hpp>
 #include <ginkgo/core/matrix/dense.hpp>
-#include <ginkgo/core/preconditioner/isai.hpp>
-#include <ginkgo/core/preconditioner/utils.hpp>
-#include <ginkgo/core/solver/gmres.hpp>
-#include <ginkgo/core/solver/ir.hpp>
 #include <ginkgo/core/solver/solver_traits.hpp>
 #include <ginkgo/core/solver/triangular.hpp>
 #include <ginkgo/core/stop/combined.hpp>
@@ -35,19 +32,13 @@ namespace preconditioner {
 namespace detail {
 
 
-template <typename LSolverType, typename USolverType>
+template <typename Type>
 constexpr bool support_ilu_parse =
-    std::is_same<typename USolverType::transposed_type, LSolverType>::value &&
-    (is_instantiation_of<LSolverType, solver::LowerTrs>::value ||
-     is_instantiation_of<LSolverType, solver::Ir>::value ||
-     is_instantiation_of<LSolverType, solver::Gmres>::value ||
-     is_instantiation_of<LSolverType, preconditioner::LowerIsai>::value);
+    std::is_same_v<typename Type::l_solver_type, LinOp>&&
+        std::is_same_v<typename Type::u_solver_type, LinOp>;
 
 
-template <typename Ilu,
-          std::enable_if_t<!support_ilu_parse<typename Ilu::l_solver_type,
-                                              typename Ilu::u_solver_type>>* =
-              nullptr>
+template <typename Ilu, std::enable_if_t<!support_ilu_parse<Ilu>>* = nullptr>
 typename Ilu::parameters_type ilu_parse(
     const config::pnode& config, const config::registry& context,
     const config::type_descriptor& td_for_child)
@@ -56,10 +47,7 @@ typename Ilu::parameters_type ilu_parse(
         "preconditioner::Ilu only supports limited type for parse.");
 }
 
-template <
-    typename Ilu,
-    std::enable_if_t<support_ilu_parse<typename Ilu::l_solver_type,
-                                       typename Ilu::u_solver_type>>* = nullptr>
+template <typename Ilu, std::enable_if_t<support_ilu_parse<Ilu>>* = nullptr>
 typename Ilu::parameters_type ilu_parse(
     const config::pnode& config, const config::registry& context,
     const config::type_descriptor& td_for_child);
@@ -101,11 +89,14 @@ typename Ilu::parameters_type ilu_parse(
  *       uses an internal cache to accelerate multiple (sequential) applies.
  *       Using it in parallel can lead to segmentation faults, wrong results
  *       and other unwanted behavior.
+ * @note The default template during parse is <ValueType, IndexType> not
+ *       <LowerTrs, IndexType>. Only the variants with ValueType are supported
+ *       in parse.
  *
- * @tparam LSolverType  type of the solver used for the L matrix.
- *                      Defaults to solver::LowerTrs
- * @tparam USolverType  type of the solver used for the U matrix
- *                      Defaults to solver::UpperTrs
+ * @tparam LSolverTypeOrValueType  type of the solver or the value type used for
+ *                                 the L matrix. Defaults to solver::LowerTrs
+ * @tparam USolverTypeOrValueType  type of the solver or the value type used for
+ *                                 the U matrix Defaults to solver::UpperTrs
  * @tparam ReverseApply  default behavior (ReverseApply = false) is first to
  *                       solve with L (Ly = b) and then with U (Ux = y).
  *                       When set to true, it will solve first with U, and then
@@ -116,28 +107,35 @@ typename Ilu::parameters_type ilu_parse(
  * @ingroup precond
  * @ingroup LinOp
  */
-template <typename LSolverType = solver::LowerTrs<>,
-          typename USolverType = solver::UpperTrs<>, bool ReverseApply = false,
-          typename IndexType = int32>
-class Ilu : public EnableLinOp<
-                Ilu<LSolverType, USolverType, ReverseApply, IndexType>>,
-            public Transposable {
+template <typename LSolverTypeOrValueType = solver::LowerTrs<>,
+          typename USolverTypeOrValueType =
+              gko::detail::transposed_type<LSolverTypeOrValueType>,
+          bool ReverseApply = false, typename IndexType = int32>
+class Ilu
+    : public EnableLinOp<Ilu<LSolverTypeOrValueType, USolverTypeOrValueType,
+                             ReverseApply, IndexType>>,
+      public Transposable {
     friend class EnableLinOp<Ilu>;
     friend class EnablePolymorphicObject<Ilu, LinOp>;
 
 public:
     static_assert(
-        std::is_same<typename LSolverType::value_type,
-                     typename USolverType::value_type>::value,
+        std::is_same_v<gko::detail::get_value_type<LSolverTypeOrValueType>,
+                       gko::detail::get_value_type<USolverTypeOrValueType>>,
         "Both the L- and the U-solver must use the same `value_type`!");
-    using value_type = typename LSolverType::value_type;
-    using l_solver_type = LSolverType;
-    using u_solver_type = USolverType;
+    using value_type = gko::detail::get_value_type<LSolverTypeOrValueType>;
+    using l_solver_type =
+        std::conditional_t<gko::detail::is_ginkgo_linop<LSolverTypeOrValueType>,
+                           LSolverTypeOrValueType, LinOp>;
+    using u_solver_type =
+        std::conditional_t<gko::detail::is_ginkgo_linop<USolverTypeOrValueType>,
+                           USolverTypeOrValueType, LinOp>;
     static constexpr bool performs_reverse_apply = ReverseApply;
     using index_type = IndexType;
     using transposed_type =
-        Ilu<typename USolverType::transposed_type,
-            typename LSolverType::transposed_type, ReverseApply, IndexType>;
+        Ilu<gko::detail::transposed_type<USolverTypeOrValueType>,
+            gko::detail::transposed_type<LSolverTypeOrValueType>, ReverseApply,
+            IndexType>;
 
     class Factory;
 
@@ -146,13 +144,13 @@ class Ilu : public EnableLinOp<
         /**
          * Factory for the L solver
          */
-        std::shared_ptr<const typename l_solver_type::Factory>
+        std::shared_ptr<const gko::detail::factory_type<l_solver_type>>
             l_solver_factory{};
 
         /**
          * Factory for the U solver
          */
-        std::shared_ptr<const typename u_solver_type::Factory>
+        std::shared_ptr<const gko::detail::factory_type<u_solver_type>>
             u_solver_factory{};
 
         /**
@@ -162,14 +160,21 @@ class Ilu : public EnableLinOp<
 
         GKO_DEPRECATED("use with_l_solver instead")
         parameters_type& with_l_solver_factory(
-            deferred_factory_parameter<const typename l_solver_type::Factory>
+            deferred_factory_parameter<
+                const gko::detail::factory_type<l_solver_type>>
                 solver)
         {
             return with_l_solver(std::move(solver));
         }
 
+        /**
+         * When LSolverTypeOrValueType is a concrete solver type, this only
+         * accepts the factory from the same concrete solver type. When
+         * LSolverTypeOrValueType is a value type, it accepts any LinOpFactory.
+         */
         parameters_type& with_l_solver(
-            deferred_factory_parameter<const typename l_solver_type::Factory>
+            deferred_factory_parameter<
+                const gko::detail::factory_type<l_solver_type>>
                 solver)
         {
             this->l_solver_generator = std::move(solver);
@@ -185,14 +190,21 @@ class Ilu : public EnableLinOp<
 
         GKO_DEPRECATED("use with_u_solver instead")
         parameters_type& with_u_solver_factory(
-            deferred_factory_parameter<const typename u_solver_type::Factory>
+            deferred_factory_parameter<
+                const gko::detail::factory_type<u_solver_type>>
                 solver)
         {
             return with_u_solver(std::move(solver));
         }
 
+        /**
+         * When USolverTypeOrValueType is a concrete solver type, this only
+         * accepts the factory from the same concrete solver type. When
+         * USolverTypeOrValueType is a value type, it accepts any LinOpFactory.
+         */
         parameters_type& with_u_solver(
-            deferred_factory_parameter<const typename u_solver_type::Factory>
+            deferred_factory_parameter<
+                const gko::detail::factory_type<u_solver_type>>
                 solver)
         {
             this->u_solver_generator = std::move(solver);
@@ -228,10 +240,12 @@ class Ilu : public EnableLinOp<
         }
 
     private:
-        deferred_factory_parameter<const typename l_solver_type::Factory>
+        deferred_factory_parameter<
+            const gko::detail::factory_type<l_solver_type>>
             l_solver_generator;
 
-        deferred_factory_parameter<const typename u_solver_type::Factory>
+        deferred_factory_parameter<
+            const gko::detail::factory_type<u_solver_type>>
             u_solver_generator;
 
         deferred_factory_parameter<const LinOpFactory> factorization_generator;
@@ -253,15 +267,16 @@ class Ilu : public EnableLinOp<
      *
      * @return parameters
      *
-     * @note only support the following pairs for <l_solver, u_solver>:
-     *       <Ir, Ir>, <Gmres, Gmres>, <LowerTrs, UpperTrs>,
-     *       and <LowerIsai, UpperIsai>
+     * @note only support the following when using <ValueType, ValueType,
+     *       ReverseApply, IndexType> not <LSolverType, USolverType,
+     *       ReverseApply, IndexType> variants
      */
     static parameters_type parse(
         const config::pnode& config, const config::registry& context,
         const config::type_descriptor& td_for_child =
             config::make_type_descriptor<value_type, index_type>())
     {
+        // parse is not templated, so we can only use SFINAE later
         return detail::ilu_parse<Ilu>(config, context, td_for_child);
     }
 
@@ -291,11 +306,11 @@ class Ilu : public EnableLinOp<
             new transposed_type{this->get_executor()}};
         transposed->set_size(gko::transpose(this->get_size()));
         transposed->l_solver_ =
-            share(as<typename u_solver_type::transposed_type>(
-                this->get_u_solver()->transpose()));
+            share(as<gko::detail::transposed_type<u_solver_type>>(
+                as<Transposable>(this->get_u_solver())->transpose()));
         transposed->u_solver_ =
-            share(as<typename l_solver_type::transposed_type>(
-                this->get_l_solver()->transpose()));
+            share(as<gko::detail::transposed_type<l_solver_type>>(
+                as<Transposable>(this->get_l_solver())->transpose()));
 
         return std::move(transposed);
     }
@@ -306,11 +321,11 @@ class Ilu : public EnableLinOp<
             new transposed_type{this->get_executor()}};
         transposed->set_size(gko::transpose(this->get_size()));
         transposed->l_solver_ =
-            share(as<typename u_solver_type::transposed_type>(
-                this->get_u_solver()->conj_transpose()));
+            share(as<gko::detail::transposed_type<u_solver_type>>(
+                as<Transposable>(this->get_u_solver())->conj_transpose()));
         transposed->u_solver_ =
-            share(as<typename l_solver_type::transposed_type>(
-                this->get_l_solver()->conj_transpose()));
+            share(as<gko::detail::transposed_type<l_solver_type>>(
+                as<Transposable>(this->get_l_solver())->conj_transpose()));
 
         return std::move(transposed);
     }
@@ -438,11 +453,7 @@ class Ilu : public EnableLinOp<
             auto fact = std::shared_ptr<const LinOp>(
                 parameters_.factorization_factory->generate(lin_op));
             // ensure that the result is a composition
-            comp =
-                std::dynamic_pointer_cast<const Composition<value_type>>(fact);
-            if (!comp) {
-                GKO_NOT_SUPPORTED(comp);
-            }
+            comp = as<const Composition<value_type>>(fact);
         }
         if (comp->get_operators().size() == 2) {
             l_factor = comp->get_operators()[0];
@@ -456,12 +467,20 @@ class Ilu : public EnableLinOp<
 
         // If no factories are provided, generate default ones
         if (!parameters_.l_solver_factory) {
-            l_solver_ = generate_default_solver<l_solver_type>(exec, l_factor);
+            // when l_solver_type is LinOp, use LowerTrs as the default one
+            l_solver_ = generate_default_solver<std::conditional_t<
+                std::is_same_v<l_solver_type, LinOp>,
+                solver::LowerTrs<value_type, index_type>, l_solver_type>>(
+                exec, l_factor);
         } else {
             l_solver_ = parameters_.l_solver_factory->generate(l_factor);
         }
         if (!parameters_.u_solver_factory) {
-            u_solver_ = generate_default_solver<u_solver_type>(exec, u_factor);
+            // when u_solver_type is LinOp, use UpperTrs as the default one
+            u_solver_ = generate_default_solver<std::conditional_t<
+                std::is_same_v<u_solver_type, LinOp>,
+                solver::UpperTrs<value_type, index_type>, u_solver_type>>(
+                exec, u_factor);
         } else {
             u_solver_ = parameters_.u_solver_factory->generate(u_factor);
         }
diff --git a/include/ginkgo/core/preconditioner/utils.hpp b/include/ginkgo/core/preconditioner/utils.hpp
index 1e3f35c8ada..f8cc3da2293 100644
--- a/include/ginkgo/core/preconditioner/utils.hpp
+++ b/include/ginkgo/core/preconditioner/utils.hpp
@@ -1,4 +1,4 @@
-// SPDX-FileCopyrightText: 2017 - 2024 The Ginkgo authors
+// SPDX-FileCopyrightText: 2017 - 2025 The Ginkgo authors
 //
 // SPDX-License-Identifier: BSD-3-Clause
 
@@ -6,31 +6,17 @@
 #define GKO_PUBLIC_CORE_PRECONDITIONER_UTILS_HPP_
 
 
-#include <type_traits>
+#ifdef __GNUC__
+#pragma message \
+    "This file is deprecated and will be removed in a later major release."
+#elif defined(_MSC_VER)
+#pragma message WARN( \
+    "This file is deprecated and will be removed in a later major release.")
+#endif
 
-#include <ginkgo/core/preconditioner/isai.hpp>
 
+// We do not have anything here, but still keep the file for non-breaking
+// changes if they include this file
 
-namespace gko {
-namespace preconditioner {
-namespace detail {
-
-
-// true_type if ExplicitType is an instantiation of TemplateType.
-template <typename ExplicitType, template <typename...> class TemplateType>
-struct is_instantiation_of : std::false_type {};
-
-template <template <typename...> class Type, typename... Param>
-struct is_instantiation_of<Type<Param...>, Type> : std::true_type {};
-
-// LowerIsai will be treated as Isai<...> so it does not match LowerIsai<...>
-template <typename ValueType, typename IndexType>
-struct is_instantiation_of<LowerIsai<ValueType, IndexType>, LowerIsai>
-    : std::true_type {};
-
-
-}  // namespace detail
-}  // namespace preconditioner
-}  // namespace gko
 
 #endif  // GKO_PUBLIC_CORE_PRECONDITIONER_UTILS_HPP_
diff --git a/reference/test/preconditioner/ilu.cpp b/reference/test/preconditioner/ilu.cpp
index 39fd67459dd..03507c38bc1 100644
--- a/reference/test/preconditioner/ilu.cpp
+++ b/reference/test/preconditioner/ilu.cpp
@@ -34,6 +34,8 @@ class Ilu : public ::testing::Test {
     using default_ilu_prec_type = gko::preconditioner::Ilu<>;
     using ilu_prec_type =
         gko::preconditioner::Ilu<l_solver_type, u_solver_type, false>;
+    using ilu_linop_prec_type =
+        gko::preconditioner::Ilu<value_type, value_type, false>;
     using ilu_rev_prec_type =
         gko::preconditioner::Ilu<l_solver_type, u_solver_type, true>;
     using Composition = gko::Composition<value_type>;
@@ -70,7 +72,14 @@ class Ilu : public ::testing::Test {
           ilu_rev_pre_factory(ilu_rev_prec_type::build()
                                   .with_l_solver(l_factory)
                                   .with_u_solver(u_factory)
-                                  .on(exec))
+                                  .on(exec)),
+          ilu_linop_pre_factory(
+              ilu_linop_prec_type::build()
+                  .with_l_solver(l_factory)
+                  .with_u_solver(u_factory)
+                  .with_factorization(
+                      gko::factorization::ParIlu<value_type>::build())
+                  .on(exec))
     {}
 
     std::shared_ptr<const gko::Executor> exec;
@@ -82,6 +91,8 @@ class Ilu : public ::testing::Test {
     std::shared_ptr<typename u_solver_type::Factory> u_factory;
     std::shared_ptr<typename ilu_prec_type::Factory> ilu_pre_factory;
     std::shared_ptr<typename ilu_rev_prec_type::Factory> ilu_rev_pre_factory;
+    std::shared_ptr<typename ilu_linop_prec_type::Factory>
+        ilu_linop_pre_factory;
 };
 
 TYPED_TEST_SUITE(Ilu, gko::test::ValueTypes, TypenameNameGenerator);
@@ -108,6 +119,12 @@ TYPED_TEST(Ilu, BuildsCustomWithoutThrowing2)
 }
 
 
+TYPED_TEST(Ilu, BuildsWithLinOpWithoutThrowing)
+{
+    ASSERT_NO_THROW(this->ilu_linop_pre_factory->generate(this->mtx));
+}
+
+
 TYPED_TEST(Ilu, ThrowOnWrongCompositionInput)
 {
     using Composition = typename TestFixture::Composition;
@@ -248,6 +265,64 @@ TYPED_TEST(Ilu, CanBeConjTransposed)
 }
 
 
+TYPED_TEST(Ilu, CanBeTransposedWithLinOp)
+{
+    using Ilu = typename TestFixture::ilu_linop_prec_type;
+    using Mtx = typename TestFixture::Mtx;
+    using l_solver_type = typename TestFixture::l_solver_type;
+    using u_solver_type = typename TestFixture::u_solver_type;
+    auto ilu = this->ilu_linop_pre_factory->generate(this->l_u_composition);
+    auto l_ref = gko::as<Mtx>(
+        gko::as<l_solver_type>(ilu->get_l_solver())->get_system_matrix());
+    auto u_ref = gko::as<Mtx>(
+        gko::as<u_solver_type>(ilu->get_u_solver())->get_system_matrix());
+
+    auto transp = gko::as<Ilu>(ilu->transpose());
+
+    auto l_transp = gko::as<Mtx>(
+        gko::as<Mtx>(
+            gko::as<u_solver_type>(transp->get_u_solver())->get_system_matrix())
+            ->transpose());
+    auto u_transp = gko::as<Mtx>(
+        gko::as<Mtx>(
+            gko::as<l_solver_type>(transp->get_l_solver())->get_system_matrix())
+            ->transpose());
+    GKO_ASSERT_MTX_EQ_SPARSITY(l_ref, l_transp);
+    GKO_ASSERT_MTX_EQ_SPARSITY(u_ref, u_transp);
+    GKO_ASSERT_MTX_NEAR(l_ref, l_transp, 0);
+    GKO_ASSERT_MTX_NEAR(u_ref, u_transp, 0);
+}
+
+
+TYPED_TEST(Ilu, CanBeConjTransposedWithLinOp)
+{
+    using Ilu = typename TestFixture::ilu_linop_prec_type;
+    using Mtx = typename TestFixture::Mtx;
+    using l_solver_type = typename TestFixture::l_solver_type;
+    using u_solver_type = typename TestFixture::u_solver_type;
+    auto ilu = this->ilu_linop_pre_factory->generate(this->l_u_composition);
+    auto l_ref = gko::as<Mtx>(
+        gko::as<l_solver_type>(ilu->get_l_solver())->get_system_matrix());
+    auto u_ref = gko::as<Mtx>(
+        gko::as<u_solver_type>(ilu->get_u_solver())->get_system_matrix());
+
+    auto transp = gko::as<Ilu>(ilu->conj_transpose());
+
+    auto l_transp = gko::as<Mtx>(
+        gko::as<Mtx>(
+            gko::as<u_solver_type>(transp->get_u_solver())->get_system_matrix())
+            ->conj_transpose());
+    auto u_transp = gko::as<Mtx>(
+        gko::as<Mtx>(
+            gko::as<l_solver_type>(transp->get_l_solver())->get_system_matrix())
+            ->conj_transpose());
+    GKO_ASSERT_MTX_EQ_SPARSITY(l_ref, l_transp);
+    GKO_ASSERT_MTX_EQ_SPARSITY(u_ref, u_transp);
+    GKO_ASSERT_MTX_NEAR(l_ref, l_transp, 0);
+    GKO_ASSERT_MTX_NEAR(u_ref, u_transp, 0);
+}
+
+
 TYPED_TEST(Ilu, SolvesCustomTypeDefaultFactorySingleRhs)
 {
     using ilu_prec_type = typename TestFixture::ilu_prec_type;
@@ -519,6 +594,24 @@ TYPED_TEST(Ilu, SolvesMultipleRhs)
 }
 
 
+TYPED_TEST(Ilu, SolvesMultipleRhsWithLinOp)
+{
+    using Mtx = typename TestFixture::Mtx;
+    using T = typename TestFixture::value_type;
+    const auto b = gko::initialize<Mtx>(
+        {I<T>{1.0, 8.0}, I<T>{3.0, 21.0}, I<T>{6.0, 24.0}}, this->exec);
+    auto x = Mtx::create(this->exec, gko::dim<2>{3, 2});
+    x->copy_from(b);
+    auto preconditioner =
+        this->ilu_linop_pre_factory->generate(this->l_u_composition);
+
+    preconditioner->apply(b, x);
+
+    GKO_ASSERT_MTX_NEAR(x, l({{-0.125, 2.0}, {0.25, 3.0}, {1.0, 1.0}}),
+                        r<TypeParam>::value * 1e+1);
+}
+
+
 TYPED_TEST(Ilu, SolvesDifferentNumberOfRhs)
 {
     using Mtx = typename TestFixture::Mtx;