From: Luca Heltai Date: Sat, 15 May 2021 13:59:40 +0000 (+0200) Subject: Moved to new interface in sundials. X-Git-Tag: v9.3.0-rc1~53^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=b108fe75d06f08db71b33b031f2ea57daa7f4251;p=dealii.git Moved to new interface in sundials. --- diff --git a/include/deal.II/sundials/ida.h b/include/deal.II/sundials/ida.h index 4a6fd13415..72dc4e87fd 100644 --- a/include/deal.II/sundials/ida.h +++ b/include/deal.II/sundials/ida.h @@ -828,6 +828,12 @@ namespace SUNDIALS * contains algebraic constraints, or Lagrange multipliers), you should * overwrite this function in order to return only the differential * components of your system. + * + * When running in parallel, every process will call this function + * independently, and syncronization will happen at the end of the + * initialization setup to communicate what components are local. Make sure + * you only return the locally owned (or locally relevant) components, in + * order to minimize communication between processes. */ std::function differential_components; @@ -877,32 +883,11 @@ namespace SUNDIALS */ void *ida_mem; - /** - * IDA solution vector. - */ - N_Vector yy; - - /** - * IDA solution derivative vector. - */ - N_Vector yp; - - /** - * IDA absolute tolerances vector. - */ - N_Vector abs_tolls; - - /** - * IDA differential components vector. - */ - N_Vector diff_id; - /** * Number of iteration required to solve the Jacobian system */ int n_iter; - /** * MPI communicator. SUNDIALS solver runs happily in * parallel. Note that if the library is compiled without MPI diff --git a/include/deal.II/sundials/n_vector.templates.h b/include/deal.II/sundials/n_vector.templates.h index 49ea702018..e289f97933 100644 --- a/include/deal.II/sundials/n_vector.templates.h +++ b/include/deal.II/sundials/n_vector.templates.h @@ -220,6 +220,10 @@ namespace SUNDIALS realtype weighted_rms_norm(N_Vector x, N_Vector w); + template + realtype + weighted_rms_norm_mask(N_Vector x, N_Vector w, N_Vector mask); + template realtype max_norm(N_Vector x); @@ -665,6 +669,24 @@ SUNDIALS::internal::NVectorOperations::weighted_rms_norm(N_Vector x, N_Vector w) +template +realtype +SUNDIALS::internal::NVectorOperations::weighted_rms_norm_mask(N_Vector x, + N_Vector w, + N_Vector mask) +{ + // TODO copy can be avoided by a custom kernel + VectorType tmp = *unwrap_nvector_const(x); + auto * w_dealii = unwrap_nvector_const(w); + auto * mask_dealii = unwrap_nvector_const(mask); + const auto n = tmp.size(); + tmp.scale(*w_dealii); + tmp.scale(*mask_dealii); + return tmp.l2_norm() / std::sqrt(n); +} + + + template realtype SUNDIALS::internal::NVectorOperations::max_norm(N_Vector x) @@ -940,11 +962,11 @@ SUNDIALS::internal::create_empty_nvector() v->ops->nvdotprod = NVectorOperations::dot_product; v->ops->nvmaxnorm = NVectorOperations::max_norm; v->ops->nvwrmsnorm = NVectorOperations::weighted_rms_norm; - // v->ops->nvwrmsnormmask = undef; v->ops->nvmin = NVectorOperations::min_element; v->ops->nvwl2norm = NVectorOperations::weighted_l2_norm; v->ops->nvl1norm = NVectorOperations::l1_norm; - ; + v->ops->nvwrmsnormmask = + NVectorOperations::weighted_rms_norm_mask; // v->ops->nvcompare = undef; // v->ops->nvinvtest = undef; // v->ops->nvconstrmask = undef; diff --git a/source/sundials/ida.cc b/source/sundials/ida.cc index 6788c62fca..bd207108cf 100644 --- a/source/sundials/ida.cc +++ b/source/sundials/ida.cc @@ -14,12 +14,17 @@ //----------------------------------------------------------- #include +#include + #include #ifdef DEAL_II_WITH_SUNDIALS # include # include + +# include +# include # ifdef DEAL_II_WITH_TRILINOS # include # include @@ -29,7 +34,6 @@ # include # endif -# include # include # include @@ -63,24 +67,13 @@ namespace SUNDIALS void * user_data) { IDA &solver = *static_cast *>(user_data); - GrowingVectorMemory mem; - - typename VectorMemory::Pointer src_yy(mem); - solver.reinit_vector(*src_yy); - typename VectorMemory::Pointer src_yp(mem); - solver.reinit_vector(*src_yp); - - typename VectorMemory::Pointer residual(mem); - solver.reinit_vector(*residual); - - copy(*src_yy, yy); - copy(*src_yp, yp); + auto *src_yy = internal::unwrap_nvector_const(yy); + auto *src_yp = internal::unwrap_nvector_const(yp); + auto *residual = internal::unwrap_nvector(rr); int err = solver.residual(tt, *src_yy, *src_yp, *residual); - copy(rr, *residual); - return err; } @@ -103,16 +96,9 @@ namespace SUNDIALS (void)resp; IDA &solver = *static_cast *>(IDA_mem->ida_user_data); - GrowingVectorMemory mem; - typename VectorMemory::Pointer src_yy(mem); - solver.reinit_vector(*src_yy); - - typename VectorMemory::Pointer src_yp(mem); - solver.reinit_vector(*src_yp); - - copy(*src_yy, yy); - copy(*src_yp, yp); + auto *src_yy = internal::unwrap_nvector_const(yy); + auto *src_yp = internal::unwrap_nvector_const(yp); int err = solver.setup_jacobian(IDA_mem->ida_tn, *src_yy, @@ -141,16 +127,13 @@ namespace SUNDIALS *static_cast *>(IDA_mem->ida_user_data); GrowingVectorMemory mem; - typename VectorMemory::Pointer src(mem); - solver.reinit_vector(*src); - typename VectorMemory::Pointer dst(mem); solver.reinit_vector(*dst); - copy(*src, b); + auto *src = internal::unwrap_nvector(b); int err = solver.solve_jacobian_system(*src, *dst); - copy(b, *dst); + *src = *dst; return err; } @@ -173,21 +156,12 @@ namespace SUNDIALS { Assert(user_data != nullptr, ExcInternalError()); IDA &solver = *static_cast *>(user_data); - GrowingVectorMemory mem; - - typename VectorMemory::Pointer src_yy(mem); - solver.reinit_vector(*src_yy); - - typename VectorMemory::Pointer src_yp(mem); - solver.reinit_vector(*src_yp); - - copy(*src_yy, yy); - copy(*src_yp, yp); + auto *src_yy = internal::unwrap_nvector_const(yy); + auto *src_yp = internal::unwrap_nvector_const(yp); int err = solver.setup_jacobian(tt, *src_yy, *src_yp, cj); - return err; } @@ -204,21 +178,11 @@ namespace SUNDIALS const IDA &solver = *static_cast *>(LS->content); - // Allocate temporary (deal.II-type) vectors into which to copy the - // N_vectors - GrowingVectorMemory mem; - typename VectorMemory::Pointer src_b(mem); - typename VectorMemory::Pointer dst_x(mem); - - solver.reinit_vector(*src_b); - solver.reinit_vector(*dst_x); - - copy(*src_b, b); + auto *src_b = internal::unwrap_nvector_const(b); + auto *dst_x = internal::unwrap_nvector(x); const int err = solver.solve_jacobian_system(*src_b, *dst_x); - copy(x, *dst_x); - return err; } @@ -233,20 +197,12 @@ namespace SUNDIALS { IDA &solver = *static_cast *>(LS->content); - // Allocate temporary (deal.II-type) vectors into which to copy the - // N_vectors - GrowingVectorMemory mem; - typename VectorMemory::Pointer src_b(mem); - typename VectorMemory::Pointer dst_x(mem); - - solver.reinit_vector(*src_b); - solver.reinit_vector(*dst_x); + auto *src_b = internal::unwrap_nvector_const(b); + auto *dst_x = internal::unwrap_nvector(x); - copy(*src_b, b); int n_iter; const int err = solver.solve_with_jacobian(*src_b, *dst_x, n_iter, tol); solver.set_n_iter(n_iter > 0 ? n_iter : 1); - copy(x, *dst_x); return err; } @@ -259,10 +215,6 @@ namespace SUNDIALS IDA::IDA(const AdditionalData &data, const MPI_Comm mpi_comm) : data(data) , ida_mem(nullptr) - , yy(nullptr) - , yp(nullptr) - , abs_tolls(nullptr) - , diff_id(nullptr) , communicator(is_serial_vector::value ? MPI_COMM_SELF : Utilities::MPI::duplicate_communicator(mpi_comm)) @@ -291,8 +243,6 @@ namespace SUNDIALS unsigned int IDA::solve_dae(VectorType &solution, VectorType &solution_dot) { - unsigned int system_size = solution.size(); - double t = data.initial_time; double h = data.initial_step_size; unsigned int step_number = 0; @@ -301,36 +251,11 @@ namespace SUNDIALS int status; (void)status; + reset(data.initial_time, data.initial_step_size, solution, solution_dot); + // The solution is stored in // solution. Here we take only a // view of it. -# ifdef DEAL_II_WITH_MPI - if (is_serial_vector::value == false) - { - const IndexSet is = solution.locally_owned_elements(); - const std::size_t local_system_size = is.n_elements(); - - yy = N_VNew_Parallel(communicator, local_system_size, system_size); - - yp = N_VNew_Parallel(communicator, local_system_size, system_size); - - diff_id = N_VNew_Parallel(communicator, local_system_size, system_size); - - abs_tolls = - N_VNew_Parallel(communicator, local_system_size, system_size); - } - else -# endif - { - Assert(is_serial_vector::value, - ExcInternalError( - "Trying to use a serial code with a parallel vector.")); - yy = N_VNew_Serial(system_size); - yp = N_VNew_Serial(system_size); - diff_id = N_VNew_Serial(system_size); - abs_tolls = N_VNew_Serial(system_size); - } - reset(data.initial_time, data.initial_step_size, solution, solution_dot); double next_time = data.initial_time; @@ -340,15 +265,15 @@ namespace SUNDIALS { next_time += data.output_period; + auto yy = internal::make_nvector_view(solution); + auto yp = internal::make_nvector_view(solution_dot); + status = IDASolve(ida_mem, next_time, &t, yy, yp, IDA_NORMAL); AssertIDA(status); status = IDAGetLastStep(ida_mem, &h); AssertIDA(status); - copy(solution, yy); - copy(solution_dot, yp); - while (solver_should_restart(t, solution, solution_dot)) reset(t, h, solution, solution_dot); @@ -357,24 +282,6 @@ namespace SUNDIALS output_step(t, solution, solution_dot, step_number); } - // Free the vectors which are no longer used. -# ifdef DEAL_II_WITH_MPI - if (is_serial_vector::value == false) - { - N_VDestroy_Parallel(yy); - N_VDestroy_Parallel(yp); - N_VDestroy_Parallel(abs_tolls); - N_VDestroy_Parallel(diff_id); - } - else -# endif - { - N_VDestroy_Serial(yy); - N_VDestroy_Serial(yp); - N_VDestroy_Serial(abs_tolls); - N_VDestroy_Serial(diff_id); - } - return step_number; } @@ -387,71 +294,26 @@ namespace SUNDIALS VectorType & solution, VectorType & solution_dot) { - unsigned int system_size; - bool first_step = (current_time == data.initial_time); + bool first_step = (current_time == data.initial_time); if (ida_mem) IDAFree(&ida_mem); ida_mem = IDACreate(); - // Free the vectors which are no longer used. - if (yy) - { -# ifdef DEAL_II_WITH_MPI - if (is_serial_vector::value == false) - { - N_VDestroy_Parallel(yy); - N_VDestroy_Parallel(yp); - N_VDestroy_Parallel(abs_tolls); - N_VDestroy_Parallel(diff_id); - } - else -# endif - { - N_VDestroy_Serial(yy); - N_VDestroy_Serial(yp); - N_VDestroy_Serial(abs_tolls); - N_VDestroy_Serial(diff_id); - } - } - int status; (void)status; - system_size = solution.size(); -# ifdef DEAL_II_WITH_MPI - if (is_serial_vector::value == false) - { - const IndexSet is = solution.locally_owned_elements(); - const std::size_t local_system_size = is.n_elements(); - - yy = N_VNew_Parallel(communicator, local_system_size, system_size); - yp = N_VNew_Parallel(communicator, local_system_size, system_size); - - diff_id = N_VNew_Parallel(communicator, local_system_size, system_size); - - abs_tolls = - N_VNew_Parallel(communicator, local_system_size, system_size); - } - else -# endif - { - yy = N_VNew_Serial(system_size); - yp = N_VNew_Serial(system_size); - diff_id = N_VNew_Serial(system_size); - abs_tolls = N_VNew_Serial(system_size); - } - - copy(yy, solution); - copy(yp, solution_dot); + auto yy = internal::make_nvector_view(solution); + auto yp = internal::make_nvector_view(solution_dot); status = IDAInit(ida_mem, t_dae_residual, current_time, yy, yp); AssertIDA(status); if (get_local_tolerances) { - copy(abs_tolls, get_local_tolerances()); - status = IDASVtolerances(ida_mem, data.relative_tolerance, abs_tolls); + const auto abs_tols = + internal::make_nvector_view(get_local_tolerances()); + status = IDASVtolerances(ida_mem, data.relative_tolerance, abs_tols); AssertIDA(status); } else @@ -477,9 +339,10 @@ namespace SUNDIALS auto dc = differential_components(); for (auto i = dc.begin(); i != dc.end(); ++i) diff_comp_vector[*i] = 1.0; + diff_comp_vector.compress(VectorOperation::insert); - copy(diff_id, diff_comp_vector); - status = IDASetId(ida_mem, diff_id); + const auto diff_id = internal::make_nvector_view(diff_comp_vector); + status = IDASetId(ida_mem, diff_id); AssertIDA(status); } @@ -625,9 +488,6 @@ namespace SUNDIALS status = IDAGetConsistentIC(ida_mem, yy, yp); AssertIDA(status); - - copy(solution, yy); - copy(solution_dot, yp); } else if (type == AdditionalData::use_y_diff) { @@ -637,9 +497,6 @@ namespace SUNDIALS status = IDAGetConsistentIC(ida_mem, yy, yp); AssertIDA(status); - - copy(solution, yy); - copy(solution_dot, yp); } } diff --git a/tests/sundials/n_vector.cc b/tests/sundials/n_vector.cc index 4968fa1938..1a9a2697c3 100644 --- a/tests/sundials/n_vector.cc +++ b/tests/sundials/n_vector.cc @@ -585,6 +585,33 @@ test_weighted_rms_norm() +template +void +test_weighted_rms_norm_mask() +{ + const auto vector_a = create_test_vector(2.0); + const auto vector_b = create_test_vector(3.0); + const auto vector_c = create_test_vector(1.0); + const auto vector_d = create_test_vector(0.0); + + auto nv_a = make_nvector_view(vector_a); + auto nv_b = make_nvector_view(vector_b); + auto nv_c = make_nvector_view(vector_c); + auto nv_d = make_nvector_view(vector_d); + { + const auto result = N_VWrmsNormMask(nv_a, nv_b, nv_c); + Assert(std::fabs(result - 6.0) < 1e-12, NVectorTestError()); + deallog << "test_weighted_rms_norm_mask 1 OK" << std::endl; + } + { + const auto result = N_VWrmsNormMask(nv_a, nv_b, nv_d); + Assert(std::fabs(result) < 1e-12, NVectorTestError()); + deallog << "test_weighted_rms_norm_mask 0 OK" << std::endl; + } +} + + + template void test_max_norm() @@ -674,6 +701,7 @@ run_all_tests(const std::string &prefix) test_elementwise_inv(); test_elementwise_abs(); test_weighted_rms_norm(); + test_weighted_rms_norm_mask(); test_max_norm(); test_min_element(); test_scale(); diff --git a/tests/sundials/n_vector.with_sundials.geq.5.0.0.with_mpi=true.mpirun=2.with_trilinos=true.with_petsc=true.output b/tests/sundials/n_vector.with_sundials.geq.5.0.0.with_mpi=true.mpirun=2.with_trilinos=true.with_petsc=true.output index 137c7ab77c..1ff8003c1c 100644 --- a/tests/sundials/n_vector.with_sundials.geq.5.0.0.with_mpi=true.mpirun=2.with_trilinos=true.with_petsc=true.output +++ b/tests/sundials/n_vector.with_sundials.geq.5.0.0.with_mpi=true.mpirun=2.with_trilinos=true.with_petsc=true.output @@ -14,6 +14,8 @@ DEAL:0:Vector::test_elementwise_div OK DEAL:0:Vector::test_elementwise_inv OK DEAL:0:Vector::test_elementwise_abs OK DEAL:0:Vector::test_weighted_rms_norm OK +DEAL:0:Vector::test_weighted_rms_norm_mask 1 OK +DEAL:0:Vector::test_weighted_rms_norm_mask 0 OK DEAL:0:Vector::test_max_norm OK DEAL:0:Vector::test_min_element OK DEAL:0:Vector::test_scale OK @@ -32,6 +34,8 @@ DEAL:0:BlockVector::test_elementwise_div OK DEAL:0:BlockVector::test_elementwise_inv OK DEAL:0:BlockVector::test_elementwise_abs OK DEAL:0:BlockVector::test_weighted_rms_norm OK +DEAL:0:BlockVector::test_weighted_rms_norm_mask 1 OK +DEAL:0:BlockVector::test_weighted_rms_norm_mask 0 OK DEAL:0:BlockVector::test_max_norm OK DEAL:0:BlockVector::test_min_element OK DEAL:0:BlockVector::test_scale OK @@ -50,6 +54,8 @@ DEAL:0:LinearAlgebra::distributed::Vector::test_elementwise_div OK DEAL:0:LinearAlgebra::distributed::Vector::test_elementwise_inv OK DEAL:0:LinearAlgebra::distributed::Vector::test_elementwise_abs OK DEAL:0:LinearAlgebra::distributed::Vector::test_weighted_rms_norm OK +DEAL:0:LinearAlgebra::distributed::Vector::test_weighted_rms_norm_mask 1 OK +DEAL:0:LinearAlgebra::distributed::Vector::test_weighted_rms_norm_mask 0 OK DEAL:0:LinearAlgebra::distributed::Vector::test_max_norm OK DEAL:0:LinearAlgebra::distributed::Vector::test_min_element OK DEAL:0:LinearAlgebra::distributed::Vector::test_scale OK @@ -68,6 +74,8 @@ DEAL:0:LinearAlgebra::distributed::BlockVector::test_elementwise_div OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_elementwise_inv OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_elementwise_abs OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_weighted_rms_norm OK +DEAL:0:LinearAlgebra::distributed::BlockVector::test_weighted_rms_norm_mask 1 OK +DEAL:0:LinearAlgebra::distributed::BlockVector::test_weighted_rms_norm_mask 0 OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_max_norm OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_min_element OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_scale OK @@ -86,6 +94,8 @@ DEAL:0:TrilinosWrappers::MPI::Vector::test_elementwise_div OK DEAL:0:TrilinosWrappers::MPI::Vector::test_elementwise_inv OK DEAL:0:TrilinosWrappers::MPI::Vector::test_elementwise_abs OK DEAL:0:TrilinosWrappers::MPI::Vector::test_weighted_rms_norm OK +DEAL:0:TrilinosWrappers::MPI::Vector::test_weighted_rms_norm_mask 1 OK +DEAL:0:TrilinosWrappers::MPI::Vector::test_weighted_rms_norm_mask 0 OK DEAL:0:TrilinosWrappers::MPI::Vector::test_max_norm OK DEAL:0:TrilinosWrappers::MPI::Vector::test_min_element OK DEAL:0:TrilinosWrappers::MPI::Vector::test_scale OK @@ -104,6 +114,8 @@ DEAL:0:TrilinosWrappers::MPI::BlockVector::test_elementwise_div OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_elementwise_inv OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_elementwise_abs OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_weighted_rms_norm OK +DEAL:0:TrilinosWrappers::MPI::BlockVector::test_weighted_rms_norm_mask 1 OK +DEAL:0:TrilinosWrappers::MPI::BlockVector::test_weighted_rms_norm_mask 0 OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_max_norm OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_min_element OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_scale OK @@ -122,6 +134,8 @@ DEAL:0:PETScWrappers::MPI::Vector::test_elementwise_div OK DEAL:0:PETScWrappers::MPI::Vector::test_elementwise_inv OK DEAL:0:PETScWrappers::MPI::Vector::test_elementwise_abs OK DEAL:0:PETScWrappers::MPI::Vector::test_weighted_rms_norm OK +DEAL:0:PETScWrappers::MPI::Vector::test_weighted_rms_norm_mask 1 OK +DEAL:0:PETScWrappers::MPI::Vector::test_weighted_rms_norm_mask 0 OK DEAL:0:PETScWrappers::MPI::Vector::test_max_norm OK DEAL:0:PETScWrappers::MPI::Vector::test_min_element OK DEAL:0:PETScWrappers::MPI::Vector::test_scale OK @@ -140,6 +154,8 @@ DEAL:0:PETScWrappers::MPI::BlockVector::test_elementwise_div OK DEAL:0:PETScWrappers::MPI::BlockVector::test_elementwise_inv OK DEAL:0:PETScWrappers::MPI::BlockVector::test_elementwise_abs OK DEAL:0:PETScWrappers::MPI::BlockVector::test_weighted_rms_norm OK +DEAL:0:PETScWrappers::MPI::BlockVector::test_weighted_rms_norm_mask 1 OK +DEAL:0:PETScWrappers::MPI::BlockVector::test_weighted_rms_norm_mask 0 OK DEAL:0:PETScWrappers::MPI::BlockVector::test_max_norm OK DEAL:0:PETScWrappers::MPI::BlockVector::test_min_element OK DEAL:0:PETScWrappers::MPI::BlockVector::test_scale OK @@ -159,6 +175,8 @@ DEAL:1:Vector::test_elementwise_div OK DEAL:1:Vector::test_elementwise_inv OK DEAL:1:Vector::test_elementwise_abs OK DEAL:1:Vector::test_weighted_rms_norm OK +DEAL:1:Vector::test_weighted_rms_norm_mask 1 OK +DEAL:1:Vector::test_weighted_rms_norm_mask 0 OK DEAL:1:Vector::test_max_norm OK DEAL:1:Vector::test_min_element OK DEAL:1:Vector::test_scale OK @@ -177,6 +195,8 @@ DEAL:1:BlockVector::test_elementwise_div OK DEAL:1:BlockVector::test_elementwise_inv OK DEAL:1:BlockVector::test_elementwise_abs OK DEAL:1:BlockVector::test_weighted_rms_norm OK +DEAL:1:BlockVector::test_weighted_rms_norm_mask 1 OK +DEAL:1:BlockVector::test_weighted_rms_norm_mask 0 OK DEAL:1:BlockVector::test_max_norm OK DEAL:1:BlockVector::test_min_element OK DEAL:1:BlockVector::test_scale OK @@ -195,6 +215,8 @@ DEAL:1:LinearAlgebra::distributed::Vector::test_elementwise_div OK DEAL:1:LinearAlgebra::distributed::Vector::test_elementwise_inv OK DEAL:1:LinearAlgebra::distributed::Vector::test_elementwise_abs OK DEAL:1:LinearAlgebra::distributed::Vector::test_weighted_rms_norm OK +DEAL:1:LinearAlgebra::distributed::Vector::test_weighted_rms_norm_mask 1 OK +DEAL:1:LinearAlgebra::distributed::Vector::test_weighted_rms_norm_mask 0 OK DEAL:1:LinearAlgebra::distributed::Vector::test_max_norm OK DEAL:1:LinearAlgebra::distributed::Vector::test_min_element OK DEAL:1:LinearAlgebra::distributed::Vector::test_scale OK @@ -213,6 +235,8 @@ DEAL:1:LinearAlgebra::distributed::BlockVector::test_elementwise_div OK DEAL:1:LinearAlgebra::distributed::BlockVector::test_elementwise_inv OK DEAL:1:LinearAlgebra::distributed::BlockVector::test_elementwise_abs OK DEAL:1:LinearAlgebra::distributed::BlockVector::test_weighted_rms_norm OK +DEAL:1:LinearAlgebra::distributed::BlockVector::test_weighted_rms_norm_mask 1 OK +DEAL:1:LinearAlgebra::distributed::BlockVector::test_weighted_rms_norm_mask 0 OK DEAL:1:LinearAlgebra::distributed::BlockVector::test_max_norm OK DEAL:1:LinearAlgebra::distributed::BlockVector::test_min_element OK DEAL:1:LinearAlgebra::distributed::BlockVector::test_scale OK @@ -231,6 +255,8 @@ DEAL:1:TrilinosWrappers::MPI::Vector::test_elementwise_div OK DEAL:1:TrilinosWrappers::MPI::Vector::test_elementwise_inv OK DEAL:1:TrilinosWrappers::MPI::Vector::test_elementwise_abs OK DEAL:1:TrilinosWrappers::MPI::Vector::test_weighted_rms_norm OK +DEAL:1:TrilinosWrappers::MPI::Vector::test_weighted_rms_norm_mask 1 OK +DEAL:1:TrilinosWrappers::MPI::Vector::test_weighted_rms_norm_mask 0 OK DEAL:1:TrilinosWrappers::MPI::Vector::test_max_norm OK DEAL:1:TrilinosWrappers::MPI::Vector::test_min_element OK DEAL:1:TrilinosWrappers::MPI::Vector::test_scale OK @@ -249,6 +275,8 @@ DEAL:1:TrilinosWrappers::MPI::BlockVector::test_elementwise_div OK DEAL:1:TrilinosWrappers::MPI::BlockVector::test_elementwise_inv OK DEAL:1:TrilinosWrappers::MPI::BlockVector::test_elementwise_abs OK DEAL:1:TrilinosWrappers::MPI::BlockVector::test_weighted_rms_norm OK +DEAL:1:TrilinosWrappers::MPI::BlockVector::test_weighted_rms_norm_mask 1 OK +DEAL:1:TrilinosWrappers::MPI::BlockVector::test_weighted_rms_norm_mask 0 OK DEAL:1:TrilinosWrappers::MPI::BlockVector::test_max_norm OK DEAL:1:TrilinosWrappers::MPI::BlockVector::test_min_element OK DEAL:1:TrilinosWrappers::MPI::BlockVector::test_scale OK @@ -267,6 +295,8 @@ DEAL:1:PETScWrappers::MPI::Vector::test_elementwise_div OK DEAL:1:PETScWrappers::MPI::Vector::test_elementwise_inv OK DEAL:1:PETScWrappers::MPI::Vector::test_elementwise_abs OK DEAL:1:PETScWrappers::MPI::Vector::test_weighted_rms_norm OK +DEAL:1:PETScWrappers::MPI::Vector::test_weighted_rms_norm_mask 1 OK +DEAL:1:PETScWrappers::MPI::Vector::test_weighted_rms_norm_mask 0 OK DEAL:1:PETScWrappers::MPI::Vector::test_max_norm OK DEAL:1:PETScWrappers::MPI::Vector::test_min_element OK DEAL:1:PETScWrappers::MPI::Vector::test_scale OK @@ -285,6 +315,8 @@ DEAL:1:PETScWrappers::MPI::BlockVector::test_elementwise_div OK DEAL:1:PETScWrappers::MPI::BlockVector::test_elementwise_inv OK DEAL:1:PETScWrappers::MPI::BlockVector::test_elementwise_abs OK DEAL:1:PETScWrappers::MPI::BlockVector::test_weighted_rms_norm OK +DEAL:1:PETScWrappers::MPI::BlockVector::test_weighted_rms_norm_mask 1 OK +DEAL:1:PETScWrappers::MPI::BlockVector::test_weighted_rms_norm_mask 0 OK DEAL:1:PETScWrappers::MPI::BlockVector::test_max_norm OK DEAL:1:PETScWrappers::MPI::BlockVector::test_min_element OK DEAL:1:PETScWrappers::MPI::BlockVector::test_scale OK diff --git a/tests/sundials/n_vector.with_sundials.geq.5.0.0.with_trilinos=true.with_petsc=true.output b/tests/sundials/n_vector.with_sundials.geq.5.0.0.with_trilinos=true.with_petsc=true.output index 7cb9f49ebd..ec2624d748 100644 --- a/tests/sundials/n_vector.with_sundials.geq.5.0.0.with_trilinos=true.with_petsc=true.output +++ b/tests/sundials/n_vector.with_sundials.geq.5.0.0.with_trilinos=true.with_petsc=true.output @@ -14,6 +14,8 @@ DEAL:0:Vector::test_elementwise_div OK DEAL:0:Vector::test_elementwise_inv OK DEAL:0:Vector::test_elementwise_abs OK DEAL:0:Vector::test_weighted_rms_norm OK +DEAL:0:Vector::test_weighted_rms_norm_mask 1 OK +DEAL:0:Vector::test_weighted_rms_norm_mask 0 OK DEAL:0:Vector::test_max_norm OK DEAL:0:Vector::test_min_element OK DEAL:0:Vector::test_scale OK @@ -32,6 +34,8 @@ DEAL:0:BlockVector::test_elementwise_div OK DEAL:0:BlockVector::test_elementwise_inv OK DEAL:0:BlockVector::test_elementwise_abs OK DEAL:0:BlockVector::test_weighted_rms_norm OK +DEAL:0:BlockVector::test_weighted_rms_norm_mask 1 OK +DEAL:0:BlockVector::test_weighted_rms_norm_mask 0 OK DEAL:0:BlockVector::test_max_norm OK DEAL:0:BlockVector::test_min_element OK DEAL:0:BlockVector::test_scale OK @@ -50,6 +54,8 @@ DEAL:0:LinearAlgebra::distributed::Vector::test_elementwise_div OK DEAL:0:LinearAlgebra::distributed::Vector::test_elementwise_inv OK DEAL:0:LinearAlgebra::distributed::Vector::test_elementwise_abs OK DEAL:0:LinearAlgebra::distributed::Vector::test_weighted_rms_norm OK +DEAL:0:LinearAlgebra::distributed::Vector::test_weighted_rms_norm_mask 1 OK +DEAL:0:LinearAlgebra::distributed::Vector::test_weighted_rms_norm_mask 0 OK DEAL:0:LinearAlgebra::distributed::Vector::test_max_norm OK DEAL:0:LinearAlgebra::distributed::Vector::test_min_element OK DEAL:0:LinearAlgebra::distributed::Vector::test_scale OK @@ -68,6 +74,8 @@ DEAL:0:LinearAlgebra::distributed::BlockVector::test_elementwise_div OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_elementwise_inv OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_elementwise_abs OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_weighted_rms_norm OK +DEAL:0:LinearAlgebra::distributed::BlockVector::test_weighted_rms_norm_mask 1 OK +DEAL:0:LinearAlgebra::distributed::BlockVector::test_weighted_rms_norm_mask 0 OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_max_norm OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_min_element OK DEAL:0:LinearAlgebra::distributed::BlockVector::test_scale OK @@ -86,6 +94,8 @@ DEAL:0:TrilinosWrappers::MPI::Vector::test_elementwise_div OK DEAL:0:TrilinosWrappers::MPI::Vector::test_elementwise_inv OK DEAL:0:TrilinosWrappers::MPI::Vector::test_elementwise_abs OK DEAL:0:TrilinosWrappers::MPI::Vector::test_weighted_rms_norm OK +DEAL:0:TrilinosWrappers::MPI::Vector::test_weighted_rms_norm_mask 1 OK +DEAL:0:TrilinosWrappers::MPI::Vector::test_weighted_rms_norm_mask 0 OK DEAL:0:TrilinosWrappers::MPI::Vector::test_max_norm OK DEAL:0:TrilinosWrappers::MPI::Vector::test_min_element OK DEAL:0:TrilinosWrappers::MPI::Vector::test_scale OK @@ -104,6 +114,8 @@ DEAL:0:TrilinosWrappers::MPI::BlockVector::test_elementwise_div OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_elementwise_inv OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_elementwise_abs OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_weighted_rms_norm OK +DEAL:0:TrilinosWrappers::MPI::BlockVector::test_weighted_rms_norm_mask 1 OK +DEAL:0:TrilinosWrappers::MPI::BlockVector::test_weighted_rms_norm_mask 0 OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_max_norm OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_min_element OK DEAL:0:TrilinosWrappers::MPI::BlockVector::test_scale OK @@ -122,6 +134,8 @@ DEAL:0:PETScWrappers::MPI::Vector::test_elementwise_div OK DEAL:0:PETScWrappers::MPI::Vector::test_elementwise_inv OK DEAL:0:PETScWrappers::MPI::Vector::test_elementwise_abs OK DEAL:0:PETScWrappers::MPI::Vector::test_weighted_rms_norm OK +DEAL:0:PETScWrappers::MPI::Vector::test_weighted_rms_norm_mask 1 OK +DEAL:0:PETScWrappers::MPI::Vector::test_weighted_rms_norm_mask 0 OK DEAL:0:PETScWrappers::MPI::Vector::test_max_norm OK DEAL:0:PETScWrappers::MPI::Vector::test_min_element OK DEAL:0:PETScWrappers::MPI::Vector::test_scale OK @@ -140,6 +154,8 @@ DEAL:0:PETScWrappers::MPI::BlockVector::test_elementwise_div OK DEAL:0:PETScWrappers::MPI::BlockVector::test_elementwise_inv OK DEAL:0:PETScWrappers::MPI::BlockVector::test_elementwise_abs OK DEAL:0:PETScWrappers::MPI::BlockVector::test_weighted_rms_norm OK +DEAL:0:PETScWrappers::MPI::BlockVector::test_weighted_rms_norm_mask 1 OK +DEAL:0:PETScWrappers::MPI::BlockVector::test_weighted_rms_norm_mask 0 OK DEAL:0:PETScWrappers::MPI::BlockVector::test_max_norm OK DEAL:0:PETScWrappers::MPI::BlockVector::test_min_element OK DEAL:0:PETScWrappers::MPI::BlockVector::test_scale OK