From ec5882ccfd2a40b986e746cf0a321f2695b01b3f Mon Sep 17 00:00:00 2001 From: maxnezddyur Date: Wed, 13 Dec 2023 08:26:38 -0700 Subject: [PATCH] add residual and jacobain trainers and transfer options update variable names fix issue with parallel implementation need to close jacobian for nonlinear update clears for other containers remove duplicate snapshot change clone to collect and add error checking update snapshot method update method names add option to variable mapping base save before refactor Fix snapshot container system. (#27101) --- framework/contrib/wasp | 2 +- framework/src/reporters/AccumulateReporter.C | 2 +- large_media | 2 +- libmesh | 2 +- .../include/reporters/JacobianContainer.h | 36 + .../include/reporters/MappingReporter.h | 1 + .../include/reporters/ResidualContainer.h | 35 + .../include/reporters/SnapshotContainerBase.h | 4 + .../transfers/SerializedSnapshotTransfer.h | 101 +++ .../include/userobjects/InverseRB.h | 165 ++++ .../stochastic_tools/include/utils/QDEIM.h | 111 +++ .../include/variablemappings/DEIMRBMapping.h | 195 +++++ .../variablemappings/VariableMappingBase.h | 9 + .../src/reporters/JacobianContainer.C | 108 +++ .../src/reporters/MappingReporter.C | 13 +- .../src/reporters/ResidualContainer.C | 51 ++ .../src/reporters/SnapshotContainerBase.C | 20 +- .../transfers/SerializedSnapshotTransfer.C | 237 ++++++ .../src/userobjects/InverseRB.C | 307 ++++++++ modules/stochastic_tools/src/utils/POD.C | 10 +- modules/stochastic_tools/src/utils/QDEIM.C | 131 ++++ .../src/variablemappings/DEIMRBMapping.C | 437 +++++++++++ .../gold/parallel_storage_main_out.json | 580 ++++++++++++++ .../gold/parallel_storage_main_out.json.1 | 580 ++++++++++++++ .../res_jac_storage/parallel_storage_main.i | 93 +++ .../tests/reporters/res_jac_storage/sub.i | 119 +++ .../tests/reporters/res_jac_storage/tests | 13 + .../gold/parallel_storage_main_out.json | 726 ++++++++++++++++++ .../gold/parallel_storage_main_out.json.1 | 726 ++++++++++++++++++ .../parallel_storage_main.i | 82 ++ .../res_jac_storage_non_linear_time/sub.i | 106 +++ .../res_jac_storage_non_linear_time/tests | 14 + .../polynomial_regression/poly_reg_vec.i | 1 + .../rb_inverse_mapping/gold/main_rb_out.e | Bin 0 -> 192040 bytes .../userobjects/rb_inverse_mapping/main_rb.i | 92 +++ .../userobjects/rb_inverse_mapping/tests | 12 + .../deim_mapping/parallel_storage_main.i | 125 +++ .../tests/variablemappings/deim_mapping/sub.i | 105 +++ .../variablemappings/deim_mapping/test_sub.i | 91 +++ .../tests/variablemappings/deim_mapping/tests | 12 + .../parallel_storage_main.i | 120 +++ .../deim_mapping_load_stepping/sub.i | 116 +++ .../deim_mapping_load_stepping/test_sub.i | 98 +++ .../deim_mapping_load_stepping/tests | 12 + .../pod_mapping/pod_mapping_main.i | 26 +- modules/stochastic_tools/unit/src/TestQDEIM.C | 84 ++ petsc | 2 +- test/tests/outputs/csv/csv_transient_vpp.i | 3 +- 48 files changed, 5887 insertions(+), 30 deletions(-) create mode 100644 modules/stochastic_tools/include/reporters/JacobianContainer.h create mode 100644 modules/stochastic_tools/include/reporters/ResidualContainer.h create mode 100644 modules/stochastic_tools/include/transfers/SerializedSnapshotTransfer.h create mode 100644 modules/stochastic_tools/include/userobjects/InverseRB.h create mode 100644 modules/stochastic_tools/include/utils/QDEIM.h create mode 100644 modules/stochastic_tools/include/variablemappings/DEIMRBMapping.h create mode 100644 modules/stochastic_tools/src/reporters/JacobianContainer.C create mode 100644 modules/stochastic_tools/src/reporters/ResidualContainer.C create mode 100644 modules/stochastic_tools/src/transfers/SerializedSnapshotTransfer.C create mode 100644 modules/stochastic_tools/src/userobjects/InverseRB.C create mode 100644 modules/stochastic_tools/src/utils/QDEIM.C create mode 100644 modules/stochastic_tools/src/variablemappings/DEIMRBMapping.C create mode 100644 modules/stochastic_tools/test/tests/reporters/res_jac_storage/gold/parallel_storage_main_out.json create mode 100644 modules/stochastic_tools/test/tests/reporters/res_jac_storage/gold/parallel_storage_main_out.json.1 create mode 100644 modules/stochastic_tools/test/tests/reporters/res_jac_storage/parallel_storage_main.i create mode 100644 modules/stochastic_tools/test/tests/reporters/res_jac_storage/sub.i create mode 100644 modules/stochastic_tools/test/tests/reporters/res_jac_storage/tests create mode 100644 modules/stochastic_tools/test/tests/reporters/res_jac_storage_non_linear_time/gold/parallel_storage_main_out.json create mode 100644 modules/stochastic_tools/test/tests/reporters/res_jac_storage_non_linear_time/gold/parallel_storage_main_out.json.1 create mode 100644 modules/stochastic_tools/test/tests/reporters/res_jac_storage_non_linear_time/parallel_storage_main.i create mode 100644 modules/stochastic_tools/test/tests/reporters/res_jac_storage_non_linear_time/sub.i create mode 100644 modules/stochastic_tools/test/tests/reporters/res_jac_storage_non_linear_time/tests create mode 100644 modules/stochastic_tools/test/tests/userobjects/rb_inverse_mapping/gold/main_rb_out.e create mode 100644 modules/stochastic_tools/test/tests/userobjects/rb_inverse_mapping/main_rb.i create mode 100644 modules/stochastic_tools/test/tests/userobjects/rb_inverse_mapping/tests create mode 100644 modules/stochastic_tools/test/tests/variablemappings/deim_mapping/parallel_storage_main.i create mode 100644 modules/stochastic_tools/test/tests/variablemappings/deim_mapping/sub.i create mode 100644 modules/stochastic_tools/test/tests/variablemappings/deim_mapping/test_sub.i create mode 100644 modules/stochastic_tools/test/tests/variablemappings/deim_mapping/tests create mode 100644 modules/stochastic_tools/test/tests/variablemappings/deim_mapping_load_stepping/parallel_storage_main.i create mode 100644 modules/stochastic_tools/test/tests/variablemappings/deim_mapping_load_stepping/sub.i create mode 100644 modules/stochastic_tools/test/tests/variablemappings/deim_mapping_load_stepping/test_sub.i create mode 100644 modules/stochastic_tools/test/tests/variablemappings/deim_mapping_load_stepping/tests create mode 100644 modules/stochastic_tools/unit/src/TestQDEIM.C diff --git a/framework/contrib/wasp b/framework/contrib/wasp index f16be318812f..28fe0d67f106 160000 --- a/framework/contrib/wasp +++ b/framework/contrib/wasp @@ -1 +1 @@ -Subproject commit f16be318812fa6b26cbf376694cf2c2beec4038d +Subproject commit 28fe0d67f10693eeb8c5c4d151af03d573962155 diff --git a/framework/src/reporters/AccumulateReporter.C b/framework/src/reporters/AccumulateReporter.C index b3b6cd245dd1..10634eeb4275 100644 --- a/framework/src/reporters/AccumulateReporter.C +++ b/framework/src/reporters/AccumulateReporter.C @@ -19,7 +19,7 @@ AccumulateReporter::validParams() "over time into a vector reporter value of the same type."); params.addRequiredParam>("reporters", "The reporters to accumulate."); - params.set("execute_on") = {EXEC_INITIAL, EXEC_TIMESTEP_END}; + params.set("execute_on") = {EXEC_MULTIAPP_FIXED_POINT_END}; params.suppressParameter("execute_on"); return params; } diff --git a/large_media b/large_media index f7f7a383fc9e..e4812f5b9245 160000 --- a/large_media +++ b/large_media @@ -1 +1 @@ -Subproject commit f7f7a383fc9eb79b2c1702beac9914423cb01d8d +Subproject commit e4812f5b9245831473dd269d4bfd5e8485f536b4 diff --git a/libmesh b/libmesh index 4e661a5ff0a3..7c03afd53c9c 160000 --- a/libmesh +++ b/libmesh @@ -1 +1 @@ -Subproject commit 4e661a5ff0a35eaf798a2817c818d10c5fe6c59c +Subproject commit 7c03afd53c9c7fd7fcb6e3e20bae37cf1b124986 diff --git a/modules/stochastic_tools/include/reporters/JacobianContainer.h b/modules/stochastic_tools/include/reporters/JacobianContainer.h new file mode 100644 index 000000000000..b992d0fba257 --- /dev/null +++ b/modules/stochastic_tools/include/reporters/JacobianContainer.h @@ -0,0 +1,36 @@ +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +#pragma once + +#include "SnapshotContainerBase.h" +#include "libmesh/id_types.h" +#include "libmesh/petsc_vector.h" +#include "libmesh/vectormap.h" + +/** + * This class is responsible for collecting jacobian row concatenated vectors in one place. The + * vectors are kept distributed with respect to the communicator of the application. + * The whole jacobian row concatenated vector is stored. + * The saving frequency can be defined using the `execute_on` parameter. + */ +class JacobianContainer : public SnapshotContainerBase +{ +public: + static InputParameters validParams(); + JacobianContainer(const InputParameters & parameters); + +protected: + virtual std::unique_ptr> collectSnapshot() override; + + std::vector> & _sparse_ind; + + NonlinearSystem & _nl_sys; + const TagID _tag_id; +}; diff --git a/modules/stochastic_tools/include/reporters/MappingReporter.h b/modules/stochastic_tools/include/reporters/MappingReporter.h index df3914b14429..06c0e2a725cb 100644 --- a/modules/stochastic_tools/include/reporters/MappingReporter.h +++ b/modules/stochastic_tools/include/reporters/MappingReporter.h @@ -57,4 +57,5 @@ class MappingReporter : public StochasticReporter, public MappingInterface std::vector> *> _vector_real_values_parallel_storage; std::vector *> _vector_real_values; ///@} + const bool _build_all_mappings_only; }; diff --git a/modules/stochastic_tools/include/reporters/ResidualContainer.h b/modules/stochastic_tools/include/reporters/ResidualContainer.h new file mode 100644 index 000000000000..b8c19378722d --- /dev/null +++ b/modules/stochastic_tools/include/reporters/ResidualContainer.h @@ -0,0 +1,35 @@ +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +#pragma once + +#include "MooseTypes.h" +#include "NonlinearSystemBase.h" +#include "SnapshotContainerBase.h" +#include "SystemBase.h" +#include "libmesh/petsc_vector.h" + +/** + * This class is responsible for collecting residual vectors in one place. The + * vectors are kept distributed with respect to the communicator of the application. + * The whole residual vector is stored. + * The saving frequency can be defined using the `execute_on` parameter. + */ +class ResidualContainer : public SnapshotContainerBase +{ +public: + static InputParameters validParams(); + ResidualContainer(const InputParameters & parameters); + +protected: + virtual std::unique_ptr> collectSnapshot() override; + + const NonlinearSystem & _nl_sys; + const TagID _tag_id; +}; diff --git a/modules/stochastic_tools/include/reporters/SnapshotContainerBase.h b/modules/stochastic_tools/include/reporters/SnapshotContainerBase.h index 3c474271fe62..24503d3152df 100644 --- a/modules/stochastic_tools/include/reporters/SnapshotContainerBase.h +++ b/modules/stochastic_tools/include/reporters/SnapshotContainerBase.h @@ -10,6 +10,7 @@ #pragma once #include "GeneralReporter.h" +#include "TaggingInterface.h" #include "libmesh/petsc_vector.h" #include "libmesh/petsc_matrix.h" #include "NonlinearSystemBase.h" @@ -59,4 +60,7 @@ class SnapshotContainerBase : public GeneralReporter /// The nonlinear system's number whose solution shall be collected const unsigned int _nonlinear_system_number; + /// Tolerance for comparing two snapshots. If comparison is below tolerance + /// snapshot is not saved. + const Real _save_tolerance; }; diff --git a/modules/stochastic_tools/include/transfers/SerializedSnapshotTransfer.h b/modules/stochastic_tools/include/transfers/SerializedSnapshotTransfer.h new file mode 100644 index 000000000000..15bc3d833be0 --- /dev/null +++ b/modules/stochastic_tools/include/transfers/SerializedSnapshotTransfer.h @@ -0,0 +1,101 @@ +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +#pragma once + +// MOOSE includes +#include "ParallelSolutionStorage.h" +#include "StochasticToolsTransfer.h" +#include "SnapshotContainerBase.h" +#include "UserObjectInterface.h" + +// Forward declarations +class ParallelSolutionStorage; + +/** + * This class is responsible for serializing solutions coming from subapps on + * specific processors. It is designed to serve as an interface between + * SolutionContainer and ParallelSolutionStorage objects. + */ +class SerializedSnapshotTransfer : public StochasticToolsTransfer +{ +public: + static InputParameters validParams(); + + SerializedSnapshotTransfer(const InputParameters & parameters); + + virtual void initialSetup() override; + virtual void execute() override; + + ///@{ + /** + * Methods used when running in batch mode (see SamplerFullSolveMultiApp) + */ + void initializeFromMultiapp() override {} + void executeFromMultiapp() override; + void finalizeFromMultiapp() override {} + + void initializeToMultiapp() override {} + void executeToMultiapp() override {} + void finalizeToMultiapp() override {} + ///@} + +protected: + /// The storage on the main application where the serialized solutions should be + /// transferred + ParallelSolutionStorage * _parallel_storage; + + /// Link to the storage spaces on the subapplications (will only hold one in batch mode) + std::vector _solution_container; + /// Link to the storage spaces on the subapplications (will only hold one in batch mode) + std::vector _residual_container; + /// Link to the storage spaces on the subapplications (will only hold one in batch mode) + std::vector _jacobian_container; + +private: + /// Serialize on the root processor of the subapplication and transfer the result to the main application + void transferToSubAppRoot(FEProblemBase & app_problem, + SnapshotContainerBase & snap_container, + const dof_id_type global_i, + const std::string snapshot_type); + + /** + * Serialize on methodically determined rank of the subapp and transfer to the main application. + * Example: Let's say we have 5 samples and 3 processors on a sub-application. + * In this case, we will serialize the first two on rank 1, the second two on rank + * 2 and the last one on rank 3. + */ + void transferInParallel(FEProblemBase & app_problem, + SnapshotContainerBase & snap_container, + const dof_id_type global_i, + const std::string snapshot_type); + + /** + * Initializes the solution container if the multiapp is run in normal mode. We need this because + * in normal mode we don't have a function for initialization besides `initialSetup()`, which + * is execute every time regardless of the multiapp settings. + */ + void initializeInNormalMode(); + + /// This routine queries the solution container addresses from the subapps. We need to redo this + /// every time initialSetup() (batch-reset) is called on the subapp because the address + /// of SolutionContainer changes. Considering that the transfer doesn't know the multiapp + /// setting, we use the same approach for batch-restore as well, which might be a little + /// wasteful if the execution of the subapps is very fast (usually not the case). + void initializeInBatchMode(); + + /// Return the system which contains the given variable. It can either be a flavor of + /// a nonlinear system or the auxiliary system. + /// @param vname The name of the variable whose system is queried + SystemBase & getSystem(FEProblemBase & app_problem, const VariableName & vname); + + /// User-selected switch that determines if we want to serialize on the root of the subapp + /// only or distribute the solutions between all the ranks of the subapp. + const bool _serialize_on_root; +}; diff --git a/modules/stochastic_tools/include/userobjects/InverseRB.h b/modules/stochastic_tools/include/userobjects/InverseRB.h new file mode 100644 index 000000000000..dfc26d478ea2 --- /dev/null +++ b/modules/stochastic_tools/include/userobjects/InverseRB.h @@ -0,0 +1,165 @@ +///* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +#pragma once + +#include "DEIMRBMapping.h" +#include "GeneralUserObject.h" +#include "MappingInterface.h" +#include "NonlinearSystemBase.h" +#include "SystemBase.h" +#include "libmesh/elem.h" +#include "libmesh/elem_range.h" +#include "libmesh/id_types.h" +#include + +/** + * InverseRB is a user object that performs inverse reduced basis mapping. + */ +class InverseRB : public GeneralUserObject, public MappingInterface +{ +public: + static InputParameters validParams(); + + /** + * Constructor for InverseRB. + * @param parameters Input parameters + */ + InverseRB(const InputParameters & parameters); + + virtual void initialize() override; + + virtual void execute() override; + + virtual void finalize() override; + + virtual void initialSetup() override; + + // only needed for ElementUserObjects and NodalUseroObjects + virtual void threadJoin(const UserObject & /*y*/) override{}; + +protected: + /// Link to the mapping object which provides the inverse mapping function + DEIMRBMapping * _mapping; + + std::vector _residual_inds; + + std::vector> _jacobian_matrix_inds; + + const DofMap & _dof_map; + + const dof_id_type _max_iter; + + const Real _tolerance; + + const Real _relax_factor; + +private: + /** + * Finds the reduced element range for the given DOFs. + * @param dofs Vector of DOF indices + * @return Vector of pointers to the reduced elements + */ + std::vector findReducedElemRange(const std::vector & dofs); + + /** + * Finds the reduced node range for the given DOFs. + * @param dofs Vector of DOF indices + * @return Vector of pointers to the reduced nodes + */ + std::vector findReducedNodeRange(const std::vector & dofs); + + /** + * Creates a range object from a vector of items. + * @tparam RangeType Type of the range object + * @tparam ItemType Type of the items in the vector + * @tparam IteratorType Type of the iterator for the range object + * @param items Vector of pointers to the items + * @return Unique pointer to the created range object + */ + template + std::unique_ptr createRangeFromVector(const std::vector & items); + + /** + * Retrieves the reduced Jacobian values from the sparse matrix. + * @param jac Sparse Jacobian matrix + * @return Vector of reduced Jacobian values + */ + std::vector getReducedJacValues(const SparseMatrix & jac); + + /** + * Retrieves the reduced residual values from the numeric vector. + * @param res Numeric residual vector + * @return Reduced residual values + */ + std::vector getReducedResValues(const NumericVector & res); + + /** + * Computes the reduced Jacobian matrix. + * @return Reduced Jacobian + */ + DenseMatrix computeReducedJacobian(); + + /** + * Computes the reduced residual vector. + * @return Reduced residual + */ + DenseVector computeReducedResidual(); + + /** + * Updates the solution vector with the reduced solution. + * @param reduced_sol Reduced solution vector + */ + void updateSolution(const DenseVector & reduced_sol); + + /** + * Computes the residual norm. + * @return Residual norm value + */ + Real computeResidual(); + + /// Vector of reduced Jacobian elements + std::vector _red_jac_elem; + + /// Vector of reduced residual elements + std::vector _red_res_elem; + + /// Vector of reduced Jacobian nodes + std::vector _red_jac_node; + + /// Vector of reduced residual nodes + std::vector _red_res_node; + + /// Local range of reduced Jacobian elements + std::unique_ptr _red_jac_elem_local_range; + + /// Local range of reduced residual elements + std::unique_ptr _red_res_elem_local_range; + + /// Local range of reduced Jacobian nodes + std::unique_ptr _red_jac_node_local_range; + + /// Local range of reduced residual nodes + std::unique_ptr _red_res_node_local_range; + + /// Nonlinear system + NonlinearSystemBase & _nl_sys; + + /// Jacobian matrix + SparseMatrix * _jac_matrix; + + /// Residual vector + NumericVector * _residual; + + /// Current solution + const NumericVector * _curr_sol; + + /// Linear solver + LinearSolver * _solver; +}; diff --git a/modules/stochastic_tools/include/utils/QDEIM.h b/modules/stochastic_tools/include/utils/QDEIM.h new file mode 100644 index 000000000000..4e5dfe7edaef --- /dev/null +++ b/modules/stochastic_tools/include/utils/QDEIM.h @@ -0,0 +1,111 @@ +///* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +#pragma once + +#include +#include +#include + +// Algorithm is derived from this resource. +// Z. Drmač and S. Gugercin, "A New Selection Operator for the Discrete Empirical Interpolation +// Method---Improved A Priori Error Bound and Extensions," SIAM Journal on Scientific Computing, +// vol. 38, no. 2, pp. A631-A648, 2016, doi: 10.1137/15M1019271. [Online]. Available: +// https://doi.org/10.1137/15M1019271 + +class QDEIM +{ +public: + /** + * Constructor for QDEIM. + * @param ortho_basis_vecs Vector of orthonormal basis vectors + */ + QDEIM(std::vector> & ortho_basis_vecs); + + /** + * Computes the selection indices using the QDEIM algorithm. + * @param selection_indices Vector to store the computed selection indices + */ + void computeSelection(std::vector & selection_indices); + + /** + * Computes both the selection indices and a new basis using the QDEIM algorithm. + * @param new_basis Pointer to a vector to store the computed new basis vectors + * @param selection_indices Vector to store the computed selection indices + */ + void computeSelectionAndBasis(std::vector> * new_basis, + std::vector & selection_indices); + +private: + /** + * Builds the orthonormal matrix from the basis vectors. + */ + void buildOrthoMatrix(); + + /** + * Performs QR decomposition on the orthonormal matrix. + */ + void performQRDecomposition(); + + /** + * Extracts the selection indices from the QR decomposition. + */ + void extractIndices(); + + /** + * Constructs the projection matrix using the selected indices. + */ + void constructProjectionMatrix(); + + /** + * Creates the inverse permutation vector. + */ + void createInversePermutationVector(); + + /** + * Applies the inverse permutation to the selection vector. + */ + void applyInversePermutation(); + + /** + * Converts an Eigen vector to a std::vector. + * @param sample_indices Vector to store the converted indices + */ + void convertEigenVectorToStdVector(std::vector & sample_indices); + + /** + * Converts an Eigen matrix to a std::vector of DenseVector. + * @param new_basis Vector to store the converted basis vectors + */ + void convertEigenMatrixToStdVector(std::vector> & new_basis); + + /// Vector of orthonormal basis vectors + std::vector> & _ortho_basis_vecs; + + /// Number of basis vectors + dof_id_type _num_basis; + + /// Length of each basis vector + dof_id_type _length_basis; + + /// Orthonormal basis matrix + RealEigenMatrix _ortho_basis; + + /// QR decomposition object + Eigen::ColPivHouseholderQR _qr; + + /// Vector of selected indices + Eigen::VectorXi _indices; + + /// Selection vector + Eigen::VectorXi _selection_vector; + + /// Projection matrix + Eigen::MatrixXd _projection_mat; +}; diff --git a/modules/stochastic_tools/include/variablemappings/DEIMRBMapping.h b/modules/stochastic_tools/include/variablemappings/DEIMRBMapping.h new file mode 100644 index 000000000000..22de87723892 --- /dev/null +++ b/modules/stochastic_tools/include/variablemappings/DEIMRBMapping.h @@ -0,0 +1,195 @@ +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +#pragma once + +#include "ReporterInterface.h" +#include "VariableMappingBase.h" +#include "ParallelSolutionStorage.h" +#include "UserObjectInterface.h" + +#include +#include "POD.h" + +/** + * Class for integrating Discrete Empirical Interpolation Method (DEIM) with + * Proper Orthogonal Decomposition (POD)-based Reduced Basis (RB) mapping. + * This class specializes in mapping between full-order model spaces and + * reduced-order model spaces, leveraging POD for dimensionality + * reduction and DEIM for efficient nonlinear approximation. This class assumes + * that the parallel storage has a solution, residual, and jacobian variable name. + */ +class DEIMRBMapping : public VariableMappingBase, + public UserObjectInterface, + public ReporterInterface +{ +public: + static InputParameters validParams(); + DEIMRBMapping(const InputParameters & parameters); + + virtual void buildMapping(const VariableName & /*vname*/) override + { + mooseError("Need to build all mappings concurrently with DEIMRBMapping."); + }; + /** + * Builds the mappings of all the components in the DEIM RB system. + */ + virtual void buildAllMappings() override; + + /** + * Method used for mapping full-order solutions for a given variable + * onto a latent space + * @param vname The name of the variable + * @param full_order_vector Serialized vector of the solution field for the given variable + * @param reduced_order_vector Storage space for the coordinates in the latent space + */ + void map(const VariableName & vname, + const DenseVector & full_order_vector, + std::vector & reduced_order_vector) const override; + + /** + * Method used for mapping full-order solutions for a given variable + * onto a latent space + * @param vname The name of the variable + * @param global_sample_i The global index of the sample whose solution should be mapped + * into the latent space + * @param reduced_order_vector Storage space for the coordinates in the latent space + */ + void map(const VariableName & vname, + const unsigned int global_sample_i, + std::vector & reduced_order_vector) const override; + + /** + * Method used for mapping reduced-order solutions for a given variable + * onto the full-order space + * @param vname The name of the variable + * @param reduced_order_vector The coordinates in the latent space + * @param full_order_vector Storage for the reconstructed solution for the given variable + */ + void inverse_map(const VariableName & vname, + const std::vector & reduced_order_vector, + DenseVector & full_order_vector) const override; + + /** + * Computes the reduced jacobian using the selection values. This function + * assumes that the jacobian values are in the same order as the selection indicies + * @param red_jac_values vector of the jacobian values + * @return RealDenseMatrix Reduced jacobian matrix + */ + RealDenseMatrix compute_reduced_jac(const std::vector & red_jac_values); + + /** + * Computes the reduced residual vector using the given residual values. This function + * assumes that the residual values are provided in a specific order corresponding to + * the problem's requirements or the arrangement of the system equations. + * + * @param red_res_values A vector of residual values used to compute the reduced residual. + * @return DenseVector The reduced residual vector. + */ + DenseVector compute_reduced_res(const std::vector & red_res_values); + + /** + * Gets the selection indices for the residual. + * @return A reference to the vector of DOF indices for the residual. + */ + const std::vector getResidualSelectionIndices() const + { + return _residual_selection_inds; + } + + /** + * Gets the reduced DOFs needed during the online phase. + * @return A reference to the set of reduced DOFs. + */ + const std::vector> getJacobianSelectionIndices() const + { + return _jac_matrix_selection_inds; + } + + dof_id_type getReducedSize() + { + checkIfReadyToUse("solution"); + return _sol_basis.size(); + } + +protected: + /// The number of modes which need to be computed + const std::vector _num_modes; + + /// The energy thresholds for truncation of the number of modes, defined by the user + const std::vector & _energy_threshold; + + /// Restartable container holding the basis functions for the solution + std::vector> & _sol_basis; + + /// Restartable vector holding the selection indices for the residual + std::vector & _residual_selection_inds; + + /// Restartable vector holding the selection indices for the jacobian + std::vector & _jacobian_selection_inds; + /// Restartable vector holding the reduced dofs that are need during online phase + std::vector> & _jac_matrix_selection_inds; + + /// Restartable vector holding the residual's reduced basis + std::vector> & _reduced_residual; + /// Restartable vector holding the reduced jacobian matrix basis + std::vector> & _reduced_jacobian; + + /// Restartable matrix holding the residual's reduced basis + DenseMatrix & _residual_selection_matrix; + /// Restartable matrix holding the jacobian's reduced basis + DenseMatrix & _jacobian_selection_matrix; + + /// Residual basis needed for computing but should not be stored + std::vector> _res_basis; + /// Jacobian basis needed for computing but should not be stored + std::vector> _jac_basis; + +private: + /** + * Build the basis for each component using POD + */ + virtual void buildComponentBases(); + /** + * Construct the jacobian matrix basis from the row-concatenated vectors + */ + virtual void constructReducedBasis(); + /** + * Get the selection indicies for the residual and jacobian + */ + virtual void computeDEIMSelectionIndices(); + /** + * Create a selection matrix from the selection indicies and basis + */ + virtual void storeSelectionMatrices(); + + std::vector> computeBasis(const VariableName & vname); + + void computeSparseMatrixConversion(); + /// Link to the parallel storage which holds the solution fields that are used for the SVD + const ParallelSolutionStorage * const _parallel_storage; + + /// The POD object which can be used to compute the basis functions/vectors + const StochasticTools::POD _pod; + + /// From the column concatenated indices, extract the unique reduced dofs + std::vector> + convertMatrixSelectionInds(dof_id_type num_dofs, const std::vector & indices); + + const std::vector> * _jacobian_indices; +}; diff --git a/modules/stochastic_tools/include/variablemappings/VariableMappingBase.h b/modules/stochastic_tools/include/variablemappings/VariableMappingBase.h index 17c88b019b16..7ffd1dbdc888 100644 --- a/modules/stochastic_tools/include/variablemappings/VariableMappingBase.h +++ b/modules/stochastic_tools/include/variablemappings/VariableMappingBase.h @@ -9,6 +9,7 @@ #pragma once +#include "MooseError.h" #include "StochasticToolsApp.h" #include "MooseObject.h" #include "RestartableModelInterface.h" @@ -34,6 +35,14 @@ class VariableMappingBase : public MooseObject, public RestartableModelInterface */ virtual void buildMapping(const VariableName & vname) = 0; + /** + * Abstract function for building all the mappings of all the variables + */ + virtual void buildAllMappings() + { + mooseError("buildAllMappings has not been implemented for this Mapping"); + }; + /** * Map a full-order solution vector (in DenseVector format) for a given variable into a * reduced-order vector (in a standard vector format) diff --git a/modules/stochastic_tools/src/reporters/JacobianContainer.C b/modules/stochastic_tools/src/reporters/JacobianContainer.C new file mode 100644 index 000000000000..38f260fde65e --- /dev/null +++ b/modules/stochastic_tools/src/reporters/JacobianContainer.C @@ -0,0 +1,108 @@ +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +#include "MooseEnum.h" +#include "ReporterName.h" +#include "JacobianContainer.h" +#include "NonlinearSystemBase.h" +#include "NonlinearSystem.h" +#include "libmesh/id_types.h" +#include "libmesh/distributed_vector.h" +#include "libmesh/nonlinear_implicit_system.h" +#include "libmesh/petsc_vector.h" +#include "libmesh/sparse_matrix.h" +#include +#include + +registerMooseObject("StochasticToolsApp", JacobianContainer); + +InputParameters +JacobianContainer::validParams() +{ + InputParameters params = SnapshotContainerBase::validParams(); + params.addClassDescription("Class responsible for collecting distributed jacobian row " + "concatenated vectors into a container. We append " + "a new distributed jacobian row concatenated at every execution."); + + params.addRequiredParam( + "tag_name", "Name of the matrix tag to collect the snapshot, defaults to the system matrix."); + params.addRequiredParam( + "jac_indices_reporter_name", "Name of the reporter that will hold the sparse indicies."); + + return params; +} + +JacobianContainer::JacobianContainer(const InputParameters & parameters) + : SnapshotContainerBase(parameters), + _sparse_ind(declareValueByName>>( + getParam("jac_indices_reporter_name"), REPORTER_MODE_REPLICATED)), + _nl_sys(_fe_problem.getNonlinearSystem(_nonlinear_system_number)), + _tag_id(_fe_problem.getMatrixTagID(getParam("tag_name"))) +{ + if (getParam("execute_on").contains(EXEC_NONLINEAR)) + paramError("execute_on", "Matrix will always be empty during NONLINEAR."); +} + +// Function to collect a snapshot of the current state of the Jacobian matrix +std::unique_ptr> +JacobianContainer::collectSnapshot() +{ + _sparse_ind.clear(); + + // Obtain a reference to the Jacobian matrix from the nonlinear system + auto & jac = _nl_sys.getMatrix(_tag_id); + + // Retrieve the dimensions of the Jacobian matrix + auto num_cols = jac.n(); + + // Close the Jacobian matrix to finalize its state for NONLINEAR systems + jac.close(); + + std::vector nnz_vals; + + // Iterate over the rows of the Jacobian matrix within the local domain + for (dof_id_type row = jac.row_start(); row < jac.row_stop(); row++) + { + // Temporary vectors to hold the non-zero values and column indices for the current row + std::vector values; + std::vector cols; + + // Extract the non-zero values and their column indices for the current row + jac.get_row(row, cols, values); + + // Adjust the column indices based on the row, for a flattened vector representation + for (auto & col : cols) + { + _sparse_ind.push_back(std::make_pair(row, col)); + col = row * num_cols + col; // Update column index + } + nnz_vals.insert(nnz_vals.end(), values.begin(), values.end()); + } + + // Now each processors has the entire vector. + _communicator.allgather(_sparse_ind); + + // Calculate the total size of the non-zero elements across all processors + dof_id_type total_size = nnz_vals.size(); + _communicator.sum(total_size); + + // Create distributed vectors for storing the non-zero values and their indices + DistributedVector sparse_vector(_communicator, total_size, nnz_vals.size()); + + // Vector to hold the distributed indices for the non-zero values + std::vector dist_ind(nnz_vals.size()); + // Initialize the distributed indices with consecutive values starting from the first local index + std::iota(dist_ind.begin(), dist_ind.end(), sparse_vector.first_local_index()); + + // Add the non-zero values to the distributed vectors + sparse_vector.add_vector(nnz_vals, dist_ind); + + // Return a clone of the flattened Jacobian vector containing the non-zero values + return sparse_vector.clone(); +} diff --git a/modules/stochastic_tools/src/reporters/MappingReporter.C b/modules/stochastic_tools/src/reporters/MappingReporter.C index 33dca0f080c2..997c97e531c8 100644 --- a/modules/stochastic_tools/src/reporters/MappingReporter.C +++ b/modules/stochastic_tools/src/reporters/MappingReporter.C @@ -37,6 +37,7 @@ MappingReporter::validParams() "Sampler be able to identify how the samples are distributed among " "the processes. Only needed if parallel storage is defined. It is important to have the " "same sampler here as the one used to prepare the snapshots in the parallel storage."); + params.addParam("build_all_mappings_only", false, "Build all the mappings only."); return params; } @@ -48,8 +49,12 @@ MappingReporter::MappingReporter(const InputParameters & parameters) : nullptr), _sampler(isParamValid("sampler") ? &getSampler("sampler") : nullptr), _mapping_name(getParam("mapping")), - _variable_names(getParam>("variables")) + _variable_names(getParam>("variables")), + _build_all_mappings_only(getParam("build_all_mappings_only")) { + if (_build_all_mappings_only) + return; + if (_parallel_storage) { if (_sampler) @@ -89,6 +94,12 @@ MappingReporter::initialSetup() void MappingReporter::execute() { + if (_build_all_mappings_only) + { + _mapping->buildAllMappings(); + return; + } + // We have two execution modes. If the parallel storage is supplied we loop over the snapshots in // the parallel storage, and project them to obtain their coefficients. if (_parallel_storage) diff --git a/modules/stochastic_tools/src/reporters/ResidualContainer.C b/modules/stochastic_tools/src/reporters/ResidualContainer.C new file mode 100644 index 000000000000..4734da7aeb04 --- /dev/null +++ b/modules/stochastic_tools/src/reporters/ResidualContainer.C @@ -0,0 +1,51 @@ +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +#include "Moose.h" +#include "MooseTypes.h" +#include "NonlinearSystem.h" +#include "ResidualContainer.h" +#include "NonlinearSystemBase.h" + +registerMooseObject("StochasticToolsApp", ResidualContainer); + +InputParameters +ResidualContainer::validParams() +{ + InputParameters params = SnapshotContainerBase::validParams(); + params.addClassDescription( + "Class responsible for collecting distributed residual vectors into a container. We append " + "a new distributed residual vector at every execution."); + params.addRequiredParam( + "tag_name", + "Name of the residual tag to collect snapshot, defaults to the total " + "residual."); + return params; +} + +ResidualContainer::ResidualContainer(const InputParameters & parameters) + : SnapshotContainerBase(parameters), + _nl_sys(_fe_problem.getNonlinearSystem(_nonlinear_system_number)), + _tag_id(_fe_problem.getVectorTagID(getParam("tag_name"))) +{ +} + +std::unique_ptr> +ResidualContainer::collectSnapshot() +{ + + auto exec_flag = _fe_problem.getCurrentExecuteOnFlag(); + // The best time to get the residual for linear problems is at the beginning + // of each timestep, but the residual is not calculated then. This forces it to run. + if (exec_flag == EXEC_TIMESTEP_BEGIN) + // This happens on timestep_begin + _fe_problem.computeResidualL2Norm(); + + return _nl_sys.getVector(_tag_id).clone(); +} diff --git a/modules/stochastic_tools/src/reporters/SnapshotContainerBase.C b/modules/stochastic_tools/src/reporters/SnapshotContainerBase.C index 63267e98bb6e..05d1b1fedfd6 100644 --- a/modules/stochastic_tools/src/reporters/SnapshotContainerBase.C +++ b/modules/stochastic_tools/src/reporters/SnapshotContainerBase.C @@ -13,11 +13,11 @@ InputParameters SnapshotContainerBase::validParams() { InputParameters params = GeneralReporter::validParams(); - params.addParam( "nonlinear_system_name", "nl0", "Option to select which nonlinear system's solution shall be stored."); + params.addParam("save_tolerance", 1e-8, "Tolerance for determining a duplicate snapshot."); return params; } @@ -26,7 +26,8 @@ SnapshotContainerBase::SnapshotContainerBase(const InputParameters & parameters) _accumulated_data(declareRestartableData>>>( "accumulated_snapshots")), _nonlinear_system_number( - _fe_problem.nlSysNum(getParam("nonlinear_system_name"))) + _fe_problem.nlSysNum(getParam("nonlinear_system_name"))), + _save_tolerance(getParam("save_tolerance")) { } @@ -49,6 +50,17 @@ SnapshotContainerBase::getSnapshot(unsigned int local_i) const void SnapshotContainerBase::execute() { - // Store the cloned snapshot. Each derived class has to implement the cloneSnapshot() method. - _accumulated_data.push_back(collectSnapshot()); + auto possible_snap = collectSnapshot(); + // Do not store zero containers. These wouldn't hold any information anyway. + if (possible_snap->linfty_norm() <= std::numeric_limits::epsilon()) + return; + + // We do this check every time we add a snap, so only need to check the last one. + if (_accumulated_data.size() > 1) + if (possible_snap->size() != _accumulated_data.back()->size()) + mooseError("Snapshot sizes changed mid-simulation."); + + // Store the cloned snapshot. Each derived class has to implement the collectSnapshot() + // method. + _accumulated_data.push_back(possible_snap->clone()); } diff --git a/modules/stochastic_tools/src/transfers/SerializedSnapshotTransfer.C b/modules/stochastic_tools/src/transfers/SerializedSnapshotTransfer.C new file mode 100644 index 000000000000..45b9fdff96b2 --- /dev/null +++ b/modules/stochastic_tools/src/transfers/SerializedSnapshotTransfer.C @@ -0,0 +1,237 @@ +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +// StochasticTools includes +#include "SerializedSnapshotTransfer.h" +#include "NonlinearSystemBase.h" +#include "Sampler.h" +#include "libmesh/id_types.h" + +registerMooseObject("StochasticToolsApp", SerializedSnapshotTransfer); + +InputParameters +SerializedSnapshotTransfer::validParams() +{ + InputParameters params = StochasticToolsTransfer::validParams(); + params.addClassDescription( + "Serializes and transfers solution vectors for given variables from sub-applications."); + params.addRequiredParam("parallel_storage", + "The name of the parallel storage reporter."); + params.addRequiredParam("solution_container", + "The name of the solution container on the subapp."); + params.addRequiredParam("residual_container", + "The name of the residual container on the subapp."); + params.addRequiredParam("jacobian_container", + "The name of the jacobian container on the subapp."); + params.addParam("jacobian_index_container", + "The name of the jacobian container on the subapp."); + params.addParam("serialize_on_root", + false, + "If we want to gather the solution fields only on the root processors of " + "the subapps before transfering to the main app."); + return params; +} + +SerializedSnapshotTransfer::SerializedSnapshotTransfer(const InputParameters & parameters) + : StochasticToolsTransfer(parameters), _serialize_on_root(getParam("serialize_on_root")) +{ + if (hasToMultiApp()) + paramError("to_multi_app", "To and between multiapp directions are not implemented"); +} + +void +SerializedSnapshotTransfer::initialSetup() +{ + // Check if we have the storage space to receive the serialized solution fields + _parallel_storage = &_fe_problem.getUserObject( + getParam("parallel_storage")); +} + +void +SerializedSnapshotTransfer::initializeInNormalMode() +{ + _solution_container.clear(); + _residual_container.clear(); + _jacobian_container.clear(); + + const auto n = getFromMultiApp()->numGlobalApps(); + const auto & serialized_solution_reporter = getParam("solution_container"); + const auto & serialized_residual_reporter = getParam("residual_container"); + const auto & serialized_jacobian_reporter = getParam("jacobian_container"); + + for (MooseIndex(n) i = 0; i < n; i++) + if (getFromMultiApp()->hasLocalApp(i)) + { + FEProblemBase & app_problem = getFromMultiApp()->appProblemBase(i); + _solution_container.push_back( + &app_problem.getUserObject(serialized_solution_reporter)); + _residual_container.push_back( + &app_problem.getUserObject(serialized_residual_reporter)); + _jacobian_container.push_back( + &app_problem.getUserObject(serialized_jacobian_reporter)); + } +} + +void +SerializedSnapshotTransfer::initializeInBatchMode() +{ + // First we fetch the solution containers from the subapps. This function is used + // in batch mode only so we will have one solution container on each rank + _solution_container.clear(); + _residual_container.clear(); + _jacobian_container.clear(); + FEProblemBase & app_problem = getFromMultiApp()->appProblemBase(_app_index); + + _solution_container.push_back(&app_problem.getUserObject( + getParam("solution_container"))); + _residual_container.push_back(&app_problem.getUserObject( + getParam("residual_container"))); + _jacobian_container.push_back(&app_problem.getUserObject( + getParam("jacobian_container"))); +} + +void +SerializedSnapshotTransfer::execute() +{ + initializeInNormalMode(); + + const auto n = getFromMultiApp()->numGlobalApps(); + + for (MooseIndex(n) i = 0; i < n; i++) + { + if (getFromMultiApp()->hasLocalApp(i)) + { + FEProblemBase & app_problem = getFromMultiApp()->appProblemBase(i); + + // Converting the local indexing to global sample indices + const unsigned int local_i = i - _sampler_ptr->getLocalRowBegin(); + + // Here we have to branch out based on if only the root processors + // need to participate in the transfer or if we would like to distribute the + // data among every processor of the subapplication + if (_serialize_on_root) + { + transferToSubAppRoot(app_problem, *_solution_container[local_i], i, "solution"); + transferToSubAppRoot(app_problem, *_residual_container[local_i], i, "residual"); + transferToSubAppRoot(app_problem, *_jacobian_container[local_i], i, "jacobian"); + } + else + { + transferInParallel(app_problem, *_solution_container[local_i], i, "solution"); + transferInParallel(app_problem, *_residual_container[local_i], i, "residual"); + transferInParallel(app_problem, *_jacobian_container[local_i], i, "jacobian"); + } + } + } +} + +void +SerializedSnapshotTransfer::executeFromMultiapp() +{ + initializeInBatchMode(); + + if (getFromMultiApp()->hasLocalApp(_app_index)) + { + FEProblemBase & app_problem = getFromMultiApp()->appProblemBase(_app_index); + + // Here we have to branch out based on if only the root processors + // need to participate in the transfer or if we would like to distribute the + // data among every processor of the subapplication + if (_serialize_on_root) + { + transferToSubAppRoot(app_problem, *_solution_container[0], _global_index, "solution"); + transferToSubAppRoot(app_problem, *_residual_container[0], _global_index, "residual"); + transferToSubAppRoot(app_problem, *_jacobian_container[0], _global_index, "jacobian"); + } + else + { + transferInParallel(app_problem, *_solution_container[0], _global_index, "solution"); + transferInParallel(app_problem, *_residual_container[0], _global_index, "residual"); + transferInParallel(app_problem, *_jacobian_container[0], _global_index, "jacobian"); + } + } +} + +void +SerializedSnapshotTransfer::transferInParallel(FEProblemBase & app_problem, + SnapshotContainerBase & solution_container, + const dof_id_type global_i, + const std::string snapshot_type) +{ + unsigned int local_app_index = global_i - _sampler_ptr->getLocalRowBegin(); + + // We need to go through this communicator because the multiapp's + // communicator is not necessarily the communicator of the underlying + // MooseObject. + const auto & comm = app_problem.comm(); + dof_id_type num_entries = _sampler_ptr->getNumberOfLocalRows(); + comm.sum(num_entries); + + // We shall distribute the samples on the given application between its processors. + // Only using a linear partitioning here for the sake of simplicity. + dof_id_type new_local_entries_begin; + dof_id_type new_local_entries_end; + dof_id_type num_new_local_entries; + + MooseUtils::linearPartitionItems(num_entries, + comm.size(), + comm.rank(), + num_new_local_entries, + new_local_entries_begin, + new_local_entries_end); + + for (unsigned int solution_i = 0; solution_i < solution_container.getContainer().size(); + ++solution_i) + { + DenseVector serialized_solution; + // Create indices of that span entire vector + auto num_vals = solution_container.getSnapshot(solution_i)->size(); + std::vector indices(num_vals); + std::iota(indices.begin(), indices.end(), 0); + + // Localize the solution and add it to the local container on the rank + // which is supposed to own it + solution_container.getSnapshot(solution_i) + ->localize( + serialized_solution.get_values(), + (local_app_index >= new_local_entries_begin && local_app_index < new_local_entries_end) + ? indices + : std::vector()); + + if (local_app_index >= new_local_entries_begin && local_app_index < new_local_entries_end) + _parallel_storage->addEntry(snapshot_type, global_i, serialized_solution); + } +} + +void +SerializedSnapshotTransfer::transferToSubAppRoot(FEProblemBase & /*app_problem*/, + SnapshotContainerBase & solution_container, + const dof_id_type global_i, + const std::string snapshot_type) +{ + + for (unsigned int solution_i = 0; solution_i < solution_container.getContainer().size(); + ++solution_i) + { + DenseVector serialized_solution; + + // In this case we always serialize on the root processor of the application. + solution_container.getSnapshot(solution_i)->localize_to_one(serialized_solution.get_values()); + + if (getFromMultiApp()->isRootProcessor()) + _parallel_storage->addEntry(snapshot_type, global_i, serialized_solution); + } +} + +SystemBase & +SerializedSnapshotTransfer::getSystem(FEProblemBase & app_problem, const VariableName & vname) +{ + auto & variable = app_problem.getVariable(_tid, vname); + return variable.sys(); +} diff --git a/modules/stochastic_tools/src/userobjects/InverseRB.C b/modules/stochastic_tools/src/userobjects/InverseRB.C new file mode 100644 index 000000000000..45a1b4065249 --- /dev/null +++ b/modules/stochastic_tools/src/userobjects/InverseRB.C @@ -0,0 +1,307 @@ +#include "DEIMRBMapping.h" +#include "InverseRB.h" +#include "NonlinearSystemBase.h" +#include "libmesh/dense_vector.h" +#include "libmesh/elem_range.h" +#include "libmesh/implicit_system.h" +#include "libmesh/numeric_vector.h" +#include "libmesh/sparse_matrix.h" + +registerMooseObject("StochasticToolsApp", InverseRB); + +InputParameters +InverseRB::validParams() +{ + InputParameters params = GeneralUserObject::validParams(); + params.addClassDescription(""); + params.addRequiredParam( + "mapping", "The name of the mapping object which provides the inverse mapping function."); + params.addParam( + "max_iter", 1e2, "Maximum number of iterations for the newton updates."); + params.addParam("tolerance", 1e-6, "Convergence tolerance for the newton loop."); + params.addParam("relaxation_factor", 1, "relaxation_factor"); + return params; +} + +InverseRB::InverseRB(const InputParameters & parameters) + : GeneralUserObject(parameters), + MappingInterface(this), + _dof_map(_sys.dofMap()), + _max_iter(getParam("max_iter")), + _tolerance(getParam("tolerance")), + _relax_factor(getParam("relaxation_factor")), + _nl_sys(_fe_problem.getNonlinearSystemBase(0)) +{ +} + +void +InverseRB::initialize() +{ + _jac_matrix = &static_cast(_nl_sys.system()).get_system_matrix(); + _residual = &_nl_sys.RHS(); + _curr_sol = &*_nl_sys.currentSolution(); +} + +void +InverseRB::execute() +{ + // Start with an initial guess of 0 for the reduced solution + // This will not work for a few things. + dof_id_type reduced_size = _mapping->getReducedSize(); + DenseVector reduced_sol(reduced_size); + + // Main loop (for instance, until convergence) + bool is_converged = false; + unsigned int iter = 0; + updateSolution(reduced_sol); + while (!is_converged && iter < _max_iter) + { + // Calculate reduced Jacobian and residual + auto reduced_jac = computeReducedJacobian(); + auto reduced_res = computeReducedResidual(); + reduced_res.scale(-1.0); + + // Perform linear solve (without updating the reduced solution yet) + DenseVector reduced_sol_update; + reduced_jac.lu_solve(reduced_res, reduced_sol_update); + reduced_sol.add(_relax_factor, reduced_sol_update); + + // Update the full solution and check for convergence + updateSolution(reduced_sol); + Real res = computeResidual(); + + is_converged = res <= _tolerance; + + iter++; + } + + if (!is_converged) + { + // Handle non-convergence scenario + } +} + +void +InverseRB::finalize() +{ +} + +void +InverseRB::initialSetup() +{ + // Retrieve the mapping object. This is used to get the indices for the residual and Jacobian. + const auto var_map = &getMapping("mapping"); + + // Dynamic cast to ensure var_map is of type DEIMRBMapping. This is crucial + // for accessing specific methods later. + _mapping = dynamic_cast(var_map); + mooseAssert(_mapping, "The DEIMRBMapping does not exist."); + + // Extract the indices for the residual and Jacobian matrices from the mapping. + _residual_inds = _mapping->getResidualSelectionIndices(); + _jacobian_matrix_inds = _mapping->getJacobianSelectionIndices(); + + // Create a set of dofs for the Jacobian. This set avoids duplicate entries. + std::set jac_dofs_set; + for (const auto & dof_pair : _jacobian_matrix_inds) + { + jac_dofs_set.insert(dof_pair.first); + jac_dofs_set.insert(dof_pair.second); + } + + std::vector jacobian_vec_inds(jac_dofs_set.begin(), jac_dofs_set.end()); + + // Determine the reduced element ranges for both the Jacobian and the residual. + // These ranges represent a subset of elements that are relevant for our calculations. + _red_jac_elem = findReducedElemRange(jacobian_vec_inds); + _red_res_elem = findReducedElemRange(_residual_inds); + + // Similarly, determine the reduced node ranges for both the Jacobian and the residual. + _red_jac_node = findReducedNodeRange(jacobian_vec_inds); + _red_res_node = findReducedNodeRange(_residual_inds); + + // Create local ranges for elements and nodes. These ranges are used to run calculations + // on a subset of elements and nodes, improving efficiency. + _red_jac_elem_local_range = + createRangeFromVector(_red_jac_elem); + _red_res_elem_local_range = + createRangeFromVector(_red_res_elem); + _red_jac_node_local_range = + createRangeFromVector(_red_jac_node); + _red_res_node_local_range = + createRangeFromVector(_red_res_node); +} + +std::vector +InverseRB::findReducedElemRange(const std::vector & dofs) +{ + std::vector reduced_elems; + + // Convert dofs to an unordered_set for faster lookup + std::unordered_set dofs_set(dofs.begin(), dofs.end()); + + for (const auto elem : *_fe_problem.mesh().getActiveLocalElementRange()) + { + std::vector elem_dofs; + _dof_map.dof_indices(elem, elem_dofs); + + // Check if any dof in elem_dofs is in dofs_set + for (const auto & dof : elem_dofs) + { + if (dofs_set.find(dof) != dofs_set.end()) + { + reduced_elems.push_back(elem); + break; // No need to check other dofs for this element + } + } + } + + return reduced_elems; +} + +std::vector +InverseRB::findReducedNodeRange(const std::vector & dofs) +{ + std::vector reduced_nodes; + + // Convert dofs to an unordered_set for faster lookup + std::unordered_set dofs_set(dofs.begin(), dofs.end()); + + // Iterate over nodes in the mesh + for (const auto node : *_fe_problem.mesh().getLocalNodeRange()) + { + + std::vector node_dofs; + _dof_map.dof_indices(node, node_dofs); + + // Check if any dof in node_dofs is in dofs_set + for (const auto & dof : node_dofs) + { + if (dofs_set.find(dof) != dofs_set.end()) + { + reduced_nodes.push_back(node); + break; // No need to check other dofs for this node + } + } + } + + return reduced_nodes; +} + +template +std::unique_ptr +InverseRB::createRangeFromVector(const std::vector & items) +{ + // Convert vector to raw pointers for range creation + ItemType * const * item_itr_begin = const_cast(items.data()); + ItemType * const * item_itr_end = item_itr_begin + items.size(); + + // Create iterators for the range + const auto range_begin = IteratorType( + item_itr_begin, item_itr_end, libMesh::Predicates::NotNull()); + const auto range_end = IteratorType( + item_itr_end, item_itr_end, libMesh::Predicates::NotNull()); + + // Create and return the unique_ptr to the range + return std::make_unique(range_begin, range_end); +} + +std::vector +InverseRB::getReducedJacValues(const SparseMatrix & jac) +{ + std::unordered_map, Real> jac_map; + + auto local_start = jac.row_start(); + auto local_end = jac.row_stop(); + + for (const auto & index_pair : _jacobian_matrix_inds) + { + dof_id_type row_index = index_pair.first; + dof_id_type col_index = index_pair.second; + + if (row_index >= local_start && row_index < local_end) + { + Real value = jac(row_index, col_index); + jac_map[index_pair] = value; + } + } + // Root has the correct map now + _communicator.set_union(jac_map); + std::vector red_jac_values; + for (const auto & index_pair : _jacobian_matrix_inds) + red_jac_values.push_back(jac_map.at(index_pair)); + + _communicator.broadcast(red_jac_values); + + return red_jac_values; +} + +std::vector +InverseRB::getReducedResValues(const NumericVector & res) +{ + std::unordered_map res_map; + + auto local_start = res.first_local_index(); + auto local_end = res.last_local_index(); + + for (const auto & index : _residual_inds) + { + if (index >= local_start && index < local_end) + { + Real value = res(index); + res_map[index] = value; + } + } + + // Root has the correct map now + _communicator.set_union(res_map); + std::vector red_res_values; + for (const auto & index : _residual_inds) + red_res_values.push_back(res_map.at(index)); + + _communicator.broadcast(red_res_values); + + return red_res_values; +} +DenseMatrix +InverseRB::computeReducedJacobian() +{ + // Compute the full Jacobian. + // TODO: We want to use the ranges later to only compute a subset of the elements + _fe_problem.computeJacobian(*_curr_sol, *_jac_matrix, 0); + // Extract reduced Jacobian values and compute reduced Jacobian + auto red_jac_values = getReducedJacValues(*_jac_matrix); + return _mapping->compute_reduced_jac(red_jac_values); +} + +DenseVector +InverseRB::computeReducedResidual() +{ + // Compute the full residual + _fe_problem.computeResidual(*_curr_sol, *_residual, 0); + // Extract reduced residual values and compute reduced residual + auto red_res_values = getReducedResValues(*_residual); + return _mapping->compute_reduced_res(red_res_values); +} + +void +InverseRB::updateSolution(const DenseVector & reduced_sol) +{ + DenseVector full_solution(_curr_sol->size()); + _mapping->inverse_map("solution", reduced_sol.get_values(), full_solution); + + auto & local_sol = _nl_sys.solution(); + for (const auto & i : make_range(local_sol.first_local_index(), local_sol.last_local_index())) + local_sol.set(i, full_solution(i)); + + local_sol.close(); + _nl_sys.update(); +} + +Real +InverseRB::computeResidual() +{ + _fe_problem.computeResidual(*_curr_sol, *_residual, 0); + auto error = _residual->l2_norm(); + return error; +} diff --git a/modules/stochastic_tools/src/utils/POD.C b/modules/stochastic_tools/src/utils/POD.C index 5c45a9143578..b7ef852b95ff 100644 --- a/modules/stochastic_tools/src/utils/POD.C +++ b/modules/stochastic_tools/src/utils/POD.C @@ -10,6 +10,7 @@ #include "MooseError.h" #include "POD.h" +#include namespace StochasticTools { @@ -126,10 +127,11 @@ POD::computePOD(const VariableName & vname, // Set the subspace size for the Lanczos method, we take twice as many // basis vectors as the requested number of POD modes. This guarantees in most of the case the // convergence of the singular triplets. - SVDSetFromOptions(svd); - SVDSetDimensions( + ierr = SVDSetFromOptions(svd); + LIBMESH_CHKERR(ierr); + ierr = SVDSetDimensions( svd, num_modes, std::min(2 * num_modes, global_rows), std::min(2 * num_modes, global_rows)); - + LIBMESH_CHKERR(ierr); // Compute the singular value triplets ierr = SVDSolve(svd); LIBMESH_CHKERR(ierr); @@ -207,6 +209,6 @@ POD::determineNumberOfModes(const std::vector & singular_values, if (ev_sum[num_modes] / ev_sum.back() > 1 - threshold) break; - return num_modes + 1; + return num_modes_compute; } } diff --git a/modules/stochastic_tools/src/utils/QDEIM.C b/modules/stochastic_tools/src/utils/QDEIM.C new file mode 100644 index 000000000000..5b44caf1181c --- /dev/null +++ b/modules/stochastic_tools/src/utils/QDEIM.C @@ -0,0 +1,131 @@ +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +#include "QDEIM.h" +#include "libmesh/int_range.h" +#include + +QDEIM::QDEIM(std::vector> & ortho_basis_vecs) + : _ortho_basis_vecs(ortho_basis_vecs) +{ +} + +void +QDEIM::computeSelection(std::vector & selection_indices) +{ + // Determine the number and lengths of the basis + _num_basis = _ortho_basis_vecs.size(); + _length_basis = (_num_basis > 0) ? _ortho_basis_vecs[0].size() : 0; + + buildOrthoMatrix(); + performQRDecomposition(); + extractIndices(); + + convertEigenVectorToStdVector(selection_indices); +} + +void +QDEIM::computeSelectionAndBasis(std::vector> * new_basis, + std::vector & selection_indices) +{ + + computeSelection(selection_indices); + + constructProjectionMatrix(); + applyInversePermutation(); + + // convert eigen vector to std vector + // Convert _M to new_basis + convertEigenMatrixToStdVector(*new_basis); +} + +void +QDEIM::buildOrthoMatrix() +{ + _ortho_basis.resize(_num_basis, _length_basis); + + // We expect the we got a vector of DenseVectors where each vector made up from + // rows of the ortho-normal basis matrix. For the QDEIM algorithm, we need the + // transpose of that matrix. + for (const auto i : make_range(_num_basis)) + { + _ortho_basis.row(i) = + Eigen::Map(_ortho_basis_vecs[i].get_values().data(), _length_basis); + } +} + +void +QDEIM::performQRDecomposition() +{ + _qr = _ortho_basis.colPivHouseholderQr(); +} + +void +QDEIM::extractIndices() +{ + // Get the column permutation matrix from QR decomposition + auto perm = _qr.colsPermutation(); + _indices = perm.indices(); + + // Extract the first _num_cols indices (pivotal indices) + _selection_vector = _indices.head(_num_basis); +} + +void +QDEIM::constructProjectionMatrix() +{ + // Construct the matrix M using the QR decomposition results + + auto R = _qr.matrixR().topLeftCorner(_num_basis, _num_basis); + auto rightSide = _qr.matrixR().topRightCorner(_num_basis, _length_basis - _num_basis); + _projection_mat = Eigen::MatrixXd::Identity(_length_basis, _num_basis); + _projection_mat.bottomRows(_length_basis - _num_basis) = + R.triangularView().solve(rightSide).transpose(); +} + +void +QDEIM::applyInversePermutation() +{ + // Create the inverse permutation vector + Eigen::VectorXi Pinverse(_length_basis); + for (const auto i : make_range(_length_basis)) + Pinverse(_indices(i)) = i; + + Eigen::MatrixXd mat_permuted(_length_basis, _num_basis); + for (const auto i : make_range(_length_basis)) + mat_permuted.row(Pinverse(i)) = _projection_mat.row(i); + + _projection_mat = mat_permuted; +} + +void +QDEIM::convertEigenMatrixToStdVector(std::vector> & new_basis) +{ + + new_basis.resize(_projection_mat.rows()); // Reserve space for rows + + for (const auto i : make_range(_projection_mat.rows())) + { + DenseVector rowVector(_projection_mat.cols()); // Create a dense vector for each row + + for (const auto j : make_range(_projection_mat.cols())) + rowVector(j) = _projection_mat(i, j); // Copy each element + + new_basis[i] = rowVector; // Add the row vector to new_basis + } +} + +void +QDEIM::convertEigenVectorToStdVector(std::vector & vector) +{ + vector.resize(_selection_vector.size()); // Reserve space for elements + + for (const auto i : make_range(_selection_vector.size())) + vector[i] = static_cast(_selection_vector(i)); +} diff --git a/modules/stochastic_tools/src/variablemappings/DEIMRBMapping.C b/modules/stochastic_tools/src/variablemappings/DEIMRBMapping.C new file mode 100644 index 000000000000..96b73a4f9e0d --- /dev/null +++ b/modules/stochastic_tools/src/variablemappings/DEIMRBMapping.C @@ -0,0 +1,437 @@ +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html + +#include "DEIMRBMapping.h" +#include +#include +#include "InputParameters.h" + +#include "MooseError.h" +#include "MooseTypes.h" +#include "QDEIM.h" +#include "ReporterInterface.h" +#include "VariableMappingBase.h" +#include "libmesh/dense_matrix.h" +#include "libmesh/dense_vector.h" + +#include "libmesh/int_range.h" +#include +#include + +registerMooseObject("StochasticToolsApp", DEIMRBMapping); + +InputParameters +DEIMRBMapping::validParams() +{ + InputParameters params = VariableMappingBase::validParams(); + params.addClassDescription("Class for integrating Discrete Empirical Interpolation Method (DEIM) " + "with Proper Orthogonal Decomposition (POD)-based Reduced Basis (RB) " + "mapping. This class specializes in mapping between full-order model " + "spaces and reduced-order model spaces."); + params.addParam( + "solution_storage", "The name of the storage reporter where the snapshots are located."); + params.addParam>( + "num_modes_to_compute", + "The number of modes that this object should compute. " + "Modes with 0 eigenvalues are filtered out, so the real number of modes " + "might be lower than this. This is also used for setting the " + "subspace sizes for distributed singular value solves. By default, the subspace used for the " + "SVD is twice as big as the number of requested vectors. For more information see the SLEPc " + "manual. If not specified, only one mode is computed per variable."); + params.addParam>( + "energy_threshold", + std::vector(), + "The energy threshold for the automatic truncation of the number of modes. In general, the " + "lower this number is the more information is retained about the system by keeping more POD " + "modes."); + params.addParam("extra_slepc_options", + "", + "Additional options for the singular/eigenvalue solvers in SLEPc."); + params.addParam("jac_index_name", + "Name of the reporter that holds the jacobian indicies."); + return params; +} + +DEIMRBMapping::DEIMRBMapping(const InputParameters & parameters) + : VariableMappingBase(parameters), + UserObjectInterface(this), + ReporterInterface(this), + _num_modes(isParamValid("num_modes_to_compute") + ? getParam>("num_modes_to_compute") + : std::vector(_variable_names.size(), 1)), + _energy_threshold(getParam>("energy_threshold")), + _sol_basis(declareModelData>>("solution_basis")), + _residual_selection_inds( + declareModelData>("residual_selection_indices")), + _jacobian_selection_inds( + declareModelData>("jacobian_selection_indices")), + _jac_matrix_selection_inds(declareModelData>>( + "jacobian_matrix_selection_indices")), + _reduced_residual(declareModelData>>("reduced_residual_basis")), + _reduced_jacobian(declareModelData>>("reduced_jacobian_basis")), + _residual_selection_matrix(declareModelData>("residual_selection_matrix")), + _jacobian_selection_matrix(declareModelData>(("jacobian_selection_matrix"))), + _parallel_storage(isParamValid("solution_storage") + ? &getUserObject("solution_storage") + : nullptr), + _pod(StochasticTools::POD( + _parallel_storage, getParam("extra_slepc_options"), _communicator)), + _jacobian_indices( + isParamValid("jac_index_name") + ? &getReporterValue>>("jac_index_name") + : nullptr) +{ + if (!isParamValid("filename")) + { + if (_num_modes.size() != _variable_names.size()) + paramError("num_modes_to_compute", + "The number of modes should be defined for each variable!"); + + for (const auto & mode : _num_modes) + if (!mode) + paramError("num_modes_to_compute", + "The number of modes should always be a positive integer!"); + + if (_energy_threshold.size()) + { + if (_energy_threshold.size() != _variable_names.size()) + paramError("energy_threshold", + "The energy thresholds should be defined for each variable!"); + + for (const auto & threshold : _energy_threshold) + if (threshold < 0 || threshold >= 1) + paramError("energy_threshold", + "The energy thresholds should always be in the [0,1) range!"); + } + +#if PETSC_VERSION_LESS_THAN(3, 14, 0) + mooseError("DEIMRBMapping is not supported with PETSc version below 3.14!"); +#else + // This mapping will only map these three "variables" + std::vector rom_components({"solution", "residual", "jacobian"}); + // Make a copy of _variable_names to sort. Don't want to change the original vector. + std::vector variable_names_sorted = _variable_names; + + // Sort both vectors + std::sort(variable_names_sorted.begin(), variable_names_sorted.end()); + std::sort(rom_components.begin(), rom_components.end()); + + if (rom_components != variable_names_sorted) + paramError("variables", + "The variable name should only be in the set {solution, residual, jacobian}."); + + for (const auto & vname : _variable_names) + { + _mapping_ready_to_use.emplace(vname, false); + } +#endif + } +} +void +DEIMRBMapping::buildAllMappings() +{ + + if (!_jacobian_indices) + mooseError("Jacobian indices reporter was not provided."); + + // Build the basis for each component + buildComponentBases(); + + // The rest of these are currently only done in serial. No need for all ranks to + // do the work. + if (_communicator.rank() == 0) + { + // Construct jacobian matrix basis from vector basis + constructReducedBasis(); + // Compute DEIM on the on residual and jacobian to get the selection + // indicies + computeDEIMSelectionIndices(); + // store the selection matrices + storeSelectionMatrices(); + + // Now we need to fix the selection indicies. They are currently the selection + // indices for the dense nonzero data. Need to convert back to sparse. + computeSparseMatrixConversion(); + + // All mappings are now ready to use + for (auto & comp : _mapping_ready_to_use) + comp.second = true; + } +} + +void +DEIMRBMapping::map(const VariableName & /*vname*/, + const DenseVector & /*full_order_vector*/, + std::vector & /*reduced_order_vector*/) const +{ + mooseError("DEIMRBMapping does not provide this functionalty."); +} + +void +DEIMRBMapping::map(const VariableName & /*vname*/, + const unsigned int /*global_sample_i*/, + std::vector & /*reduced_order_vector*/) const +{ + mooseError("DEIMRBMapping does not provide this functionalty."); +} + +void +DEIMRBMapping::inverse_map(const VariableName & vname, + const std::vector & reduced_order_vector, + DenseVector & full_order_vector) const +{ + if (vname != "solution") + mooseError("We are only inverse mapping the solution currently"); + + if (reduced_order_vector.size() != _sol_basis.size()) + mooseError("The number of supplied reduced-order coefficients (", + reduced_order_vector.size(), + ") is not the same as the number of basis functions (", + _sol_basis.size(), + ") for the solution!"); + + // This zeros the DenseVector too + full_order_vector.resize(_sol_basis[0].size()); + + for (auto base_i : index_range(reduced_order_vector)) + for (unsigned int dof_i = 0; dof_i < _sol_basis[base_i].size(); ++dof_i) + full_order_vector(dof_i) += reduced_order_vector[base_i] * _sol_basis[base_i](dof_i); +} + +void +DEIMRBMapping::buildComponentBases() +{ + _sol_basis = computeBasis("solution"); + + _res_basis = computeBasis("residual"); + + _jac_basis = computeBasis("jacobian"); + + for (const auto & basis : {_sol_basis, _res_basis, _jac_basis}) + if (basis.size() == 0) + mooseError("Basis is empty. This will not work for DEIMRBMapping."); +} + +void +DEIMRBMapping::constructReducedBasis() +{ + // Need to generate a basis matrix for multiplication + // _sol_basis is a vector of columns of the matrix + + RealEigenMatrix sol_basis(_sol_basis[0].size(), _sol_basis.size()); + for (MooseIndex(_sol_basis) col = 0; col < _sol_basis.size(); ++col) + { + const auto & column_vector = _sol_basis[col]; + for (MooseIndex(column_vector) row = 0; row < column_vector.size(); ++row) + { + sol_basis(row, col) = column_vector(row); + } + } + + RealEigenMatrix sol_basis_transpose = sol_basis.transpose(); + + // Calculate the dimensions of the square matrix + const auto matrix_size = _sol_basis[0].size(); + + for (const auto & basis : _jac_basis) + { + // Create a new matrix for each basis + Eigen::SparseMatrix matrix(matrix_size, matrix_size); + + // Prepare to fill the matrix using triplets + std::vector> triplets; + triplets.reserve(_jacobian_indices->size()); // Reserve space to avoid reallocations + + for (const auto i : make_range(_jacobian_indices->size())) + { + auto index_pair = _jacobian_indices->at(i); + Real value = basis(i); // Corresponding value for this entry + + auto row = index_pair.first; // division to get the row + auto col = index_pair.second; // Modulus to get the column + + // Add a triplet for the current non-zero entry + triplets.emplace_back(row, col, value); + } + + // Bulk insert the triplets into the matrix + matrix.setFromTriplets(triplets.begin(), triplets.end()); + + // Optionally, make the matrix compressed to optimize operations + // matrix.makeCompressed(); + + // Compute reduced matrix: V^T * J_i * V + RealEigenMatrix reduced_matrix = sol_basis_transpose * matrix * sol_basis; + + // Fill the real dense matrix with the values from the reduced matrix + RealDenseMatrix real_dense_mat(reduced_matrix.rows(), reduced_matrix.cols()); + for (const auto i : make_range(reduced_matrix.rows())) + { + for (const auto j : make_range(reduced_matrix.cols())) + { + real_dense_mat(i, j) = reduced_matrix(i, j); + } + } + + // Add the matrix to the vector of matrices + _reduced_jacobian.push_back(real_dense_mat); + } + + RealDenseMatrix real_dense_solution_basis(sol_basis.rows(), sol_basis.cols()); + for (const auto i : make_range(sol_basis.rows())) + { + for (const auto j : make_range(sol_basis.cols())) + { + real_dense_solution_basis(i, j) = sol_basis(i, j); + } + } + + // * Residual Basis; + for (const auto & basis : _res_basis) + { + DenseVector red_resid; + real_dense_solution_basis.vector_mult_transpose(red_resid, basis); + + _reduced_residual.push_back(red_resid); + } +} + +void +DEIMRBMapping::computeDEIMSelectionIndices() +{ + + // QDEIM for Jacobian + QDEIM jac_qd(_jac_basis); + jac_qd.computeSelection(_jacobian_selection_inds); + + // QDEIM for residual + QDEIM res_qd(_res_basis); + res_qd.computeSelection(_residual_selection_inds); +} + +void +DEIMRBMapping::storeSelectionMatrices() +{ + + // Save jacobian selection matrix + const auto num_jac_sel = _jacobian_selection_inds.size(); + _jacobian_selection_matrix.resize(num_jac_sel, num_jac_sel); + for (MooseIndex(num_jac_sel) i = 0; i < num_jac_sel; ++i) + { + dof_id_type selected_row = _jacobian_selection_inds[i]; // Get the index of the selected row + for (MooseIndex(_jac_basis) col = 0; col < _jac_basis.size(); ++col) + // Assign the value from the selected row of each column in _jac_basis + _jacobian_selection_matrix(i, col) = _jac_basis[col](selected_row); + } + + // save residual selection matrix + const auto num_res_sel = _residual_selection_inds.size(); + _residual_selection_matrix.resize(num_res_sel, + _res_basis.size()); // Resize to the correct dimensions + + for (MooseIndex(num_res_sel) selection = 0; selection < num_res_sel; ++selection) + { + dof_id_type selected_row = + _residual_selection_inds[selection]; // Get the index of the selected row + for (MooseIndex(_res_basis) col = 0; col < _res_basis.size(); ++col) + { + // Assign the value from the selected row of each column in _res_basis + _residual_selection_matrix(selection, col) = _res_basis[col](selected_row); + } + } +} + +std::vector> +DEIMRBMapping::computeBasis(const VariableName & vname) +{ + + std::vector> basis; + // TODO: These are throwaway vectors that are not needed. Change POD class so + // it doesn't need it + std::vector> right_basis; + std::vector singular_values; + + // First find and create the basis for the solution + auto it = std::find(_variable_names.begin(), _variable_names.end(), vname); + unsigned int index = std::distance(_variable_names.begin(), it); + // Find the number of modes that we would want to compute + std::size_t num_modes_compute = (std::size_t)_num_modes[index]; + const Real threshold = (_energy_threshold.empty() ? 0.0 : _energy_threshold[index]) + + std::numeric_limits::epsilon(); + + _pod.computePOD(vname, basis, right_basis, singular_values, num_modes_compute, threshold); + // Again we don't need these so clearing them right away + + return basis; +} + +std::vector> +DEIMRBMapping::convertMatrixSelectionInds(dof_id_type num_dofs, + const std::vector & indices) +{ + std::vector> dof_pairs; + + for (const auto & index : indices) + { + dof_id_type row = index / num_dofs; // Determine the row in the Jacobian matrix + dof_id_type col = index % num_dofs; // Determine the column in the Jacobian matrix + + // Add the row and column as a pair + dof_pairs.emplace_back(row, col); + } + + return dof_pairs; +} + +RealDenseMatrix +DEIMRBMapping::compute_reduced_jac(const std::vector & red_jac_values) +{ + // We assume the reduced values come in the same order as + // _jacobian_selection_inds + DenseVector jac_theta; + + _jacobian_selection_matrix.lu_solve(red_jac_values, jac_theta); + + // Now we can construct the reduced jacobian + RealDenseMatrix reduce_jac(_reduced_jacobian[0].m(), _reduced_jacobian[0].n()); + + for (MooseIndex(jac_theta) i = 0; i < jac_theta.size(); i++) + reduce_jac.add(jac_theta(i), _reduced_jacobian[i]); + + return reduce_jac; +} + +DenseVector +DEIMRBMapping::compute_reduced_res(const std::vector & red_res_values) +{ + // We assume the reduced values come in the same order as + // _residual_selection_inds + DenseVector res_theta; + + _residual_selection_matrix.lu_solve(red_res_values, res_theta); + + // Now we can construct the reduced residual + DenseVector reduce_res(_reduced_residual[0].size()); + + for (MooseIndex(res_theta) i = 0; i < res_theta.size(); i++) + reduce_res.add(res_theta(i), _reduced_residual[i]); + + return reduce_res; +} + +void +DEIMRBMapping::computeSparseMatrixConversion() +{ + + _jac_matrix_selection_inds.clear(); + + // Convert from the selection indices + for (auto & sel_index : _jacobian_selection_inds) + _jac_matrix_selection_inds.push_back(_jacobian_indices->at(sel_index)); +} diff --git a/modules/stochastic_tools/test/tests/reporters/res_jac_storage/gold/parallel_storage_main_out.json b/modules/stochastic_tools/test/tests/reporters/res_jac_storage/gold/parallel_storage_main_out.json new file mode 100644 index 000000000000..5187f0dfeb77 --- /dev/null +++ b/modules/stochastic_tools/test/tests/reporters/res_jac_storage/gold/parallel_storage_main_out.json @@ -0,0 +1,580 @@ +{ + "number_of_parts": 2, + "part": 0, + "reporters": { + "jac_indices": { + "type": "ConstantReporter", + "values": { + "jac_indices:jacobian_storage:indices": { + "type": "std::vector>>" + } + } + }, + "parallel_storage": { + "type": "ParallelSolutionStorage", + "values": { + "parallel_solution_storage": { + "type": "std::map>, std::hash, std::equal_to, std::allocator>>>>, std::less, std::allocator>, std::hash, std::equal_to, std::allocator>>>>>>>" + } + } + } + }, + "time_steps": [ + { + "jac_indices": { + "jac_indices:jacobian_storage:indices": [ + [ + [ + 0, + 0 + ], + [ + 0, + 1 + ], + [ + 0, + 2 + ], + [ + 0, + 3 + ], + [ + 1, + 0 + ], + [ + 1, + 1 + ], + [ + 1, + 2 + ], + [ + 1, + 3 + ], + [ + 2, + 0 + ], + [ + 2, + 1 + ], + [ + 2, + 2 + ], + [ + 2, + 3 + ], + [ + 2, + 4 + ], + [ + 2, + 5 + ], + [ + 3, + 0 + ], + [ + 3, + 1 + ], + [ + 3, + 2 + ], + [ + 3, + 3 + ], + [ + 3, + 4 + ], + [ + 3, + 5 + ], + [ + 4, + 2 + ], + [ + 4, + 3 + ], + [ + 4, + 4 + ], + [ + 4, + 5 + ], + [ + 4, + 6 + ], + [ + 4, + 7 + ], + [ + 5, + 2 + ], + [ + 5, + 3 + ], + [ + 5, + 4 + ], + [ + 5, + 5 + ], + [ + 5, + 6 + ], + [ + 5, + 7 + ], + [ + 6, + 4 + ], + [ + 6, + 5 + ], + [ + 6, + 6 + ], + [ + 6, + 7 + ], + [ + 6, + 8 + ], + [ + 6, + 9 + ], + [ + 7, + 4 + ], + [ + 7, + 5 + ], + [ + 7, + 6 + ], + [ + 7, + 7 + ], + [ + 7, + 8 + ], + [ + 7, + 9 + ], + [ + 8, + 6 + ], + [ + 8, + 7 + ], + [ + 8, + 8 + ], + [ + 8, + 9 + ], + [ + 8, + 10 + ], + [ + 8, + 11 + ], + [ + 9, + 6 + ], + [ + 9, + 7 + ], + [ + 9, + 8 + ], + [ + 9, + 9 + ], + [ + 9, + 10 + ], + [ + 9, + 11 + ], + [ + 10, + 8 + ], + [ + 10, + 9 + ], + [ + 10, + 10 + ], + [ + 10, + 11 + ], + [ + 10, + 12 + ], + [ + 10, + 13 + ], + [ + 11, + 8 + ], + [ + 11, + 9 + ], + [ + 11, + 10 + ], + [ + 11, + 11 + ], + [ + 11, + 12 + ], + [ + 11, + 13 + ], + [ + 12, + 10 + ], + [ + 12, + 11 + ], + [ + 12, + 12 + ], + [ + 12, + 13 + ], + [ + 13, + 10 + ], + [ + 13, + 11 + ], + [ + 13, + 12 + ], + [ + 13, + 13 + ] + ] + ] + }, + "parallel_storage": { + "parallel_solution_storage": { + "jacobian": { + "0": [ + [ + 1.0, + 0.0, + 0.0, + 0.0, + 0.0, + 1.0, + 0.0, + 0.0, + -2.0, + 0.0, + 4.0, + 0.0, + -2.0, + 0.0, + 0.0, + -4.0, + 0.0, + 8.0, + 0.0, + -4.0, + -2.0, + 0.0, + 4.0, + 0.0, + -2.0, + 0.0, + 0.0, + -4.0, + 0.0, + 8.0, + 0.0, + -4.0, + -2.0, + 0.0, + 4.0, + 0.0, + -2.0, + 0.0, + 0.0, + -4.0, + 0.0, + 8.0, + 0.0, + -4.0, + -2.0, + 0.0, + 4.0, + 0.0, + -2.0, + 0.0, + 0.0, + -4.0, + 0.0, + 8.0, + 0.0, + -4.0, + -2.0, + 0.0, + 4.0, + 0.0, + -2.0, + 0.0, + 0.0, + -4.0, + 0.0, + 8.0, + 0.0, + -4.0, + 0.0, + 0.0, + 1.0, + 0.0, + 0.0, + 0.0, + 0.0, + 1.0 + ] + ], + "1": [ + [ + 1.0, + 0.0, + 0.0, + 0.0, + 0.0, + 1.0, + 0.0, + 0.0, + -2.0, + 0.0, + 4.0, + 0.0, + -2.0, + 0.0, + 0.0, + -4.0, + 0.0, + 8.0, + 0.0, + -4.0, + -2.0, + 0.0, + 4.0, + 0.0, + -2.0, + 0.0, + 0.0, + -4.0, + 0.0, + 8.0, + 0.0, + -4.0, + -2.0, + 0.0, + 4.0, + 0.0, + -2.0, + 0.0, + 0.0, + -4.0, + 0.0, + 8.0, + 0.0, + -4.0, + -2.0, + 0.0, + 4.0, + 0.0, + -2.0, + 0.0, + 0.0, + -4.0, + 0.0, + 8.0, + 0.0, + -4.0, + -2.0, + 0.0, + 4.0, + 0.0, + -2.0, + 0.0, + 0.0, + -4.0, + 0.0, + 8.0, + 0.0, + -4.0, + 0.0, + 0.0, + 1.0, + 0.0, + 0.0, + 0.0, + 0.0, + 1.0 + ] + ] + }, + "residual": { + "0": [ + [ + 0.0, + 0.0, + 0.0, + 0.0, + 0.0, + 0.0, + 8.881784197001252e-16, + 0.0, + -7.105427357601002e-15, + -8.881784197001252e-16, + 3.552713678800501e-15, + 0.0, + 0.0, + 0.0 + ] + ], + "1": [ + [ + 0.0, + 0.0, + -1.7763568394002505e-15, + 0.0, + 0.0, + 0.0, + 0.0, + -3.552713678800501e-15, + 0.0, + 0.0, + -1.7763568394002505e-15, + -3.552713678800501e-15, + 0.0, + 0.0 + ] + ] + }, + "solution": { + "0": [ + [ + 0.0, + 0.0, + 9.92699954624491, + 0.8913743705221333, + 15.883199273991854, + 1.5327487410442666, + 17.868599183240836, + 1.9241231115664, + 15.883199273991853, + 2.0654974820885332, + 9.92699954624491, + 1.9568718526106668, + 0.0, + 1.5982462231328003 + ] + ], + "1": [ + [ + 0.0, + 0.0, + 4.884840824098017, + 2.1858085930557767, + 7.815745318556828, + 4.121617186111553, + 8.792713483376431, + 5.80742577916733, + 7.815745318556828, + 7.243234372223108, + 4.884840824098017, + 8.429042965278885, + 0.0, + 9.364851558334664 + ] + ] + } + } + }, + "time": 2.0, + "time_step": 2 + } + ] +} diff --git a/modules/stochastic_tools/test/tests/reporters/res_jac_storage/gold/parallel_storage_main_out.json.1 b/modules/stochastic_tools/test/tests/reporters/res_jac_storage/gold/parallel_storage_main_out.json.1 new file mode 100644 index 000000000000..2c9d04ca0e43 --- /dev/null +++ b/modules/stochastic_tools/test/tests/reporters/res_jac_storage/gold/parallel_storage_main_out.json.1 @@ -0,0 +1,580 @@ +{ + "number_of_parts": 2, + "part": 1, + "reporters": { + "jac_indices": { + "type": "ConstantReporter", + "values": { + "jac_indices:jacobian_storage:indices": { + "type": "std::vector>>" + } + } + }, + "parallel_storage": { + "type": "ParallelSolutionStorage", + "values": { + "parallel_solution_storage": { + "type": "std::map>, std::hash, std::equal_to, std::allocator>>>>, std::less, std::allocator>, std::hash, std::equal_to, std::allocator>>>>>>>" + } + } + } + }, + "time_steps": [ + { + "jac_indices": { + "jac_indices:jacobian_storage:indices": [ + [ + [ + 0, + 0 + ], + [ + 0, + 1 + ], + [ + 0, + 2 + ], + [ + 0, + 3 + ], + [ + 1, + 0 + ], + [ + 1, + 1 + ], + [ + 1, + 2 + ], + [ + 1, + 3 + ], + [ + 2, + 0 + ], + [ + 2, + 1 + ], + [ + 2, + 2 + ], + [ + 2, + 3 + ], + [ + 2, + 4 + ], + [ + 2, + 5 + ], + [ + 3, + 0 + ], + [ + 3, + 1 + ], + [ + 3, + 2 + ], + [ + 3, + 3 + ], + [ + 3, + 4 + ], + [ + 3, + 5 + ], + [ + 4, + 2 + ], + [ + 4, + 3 + ], + [ + 4, + 4 + ], + [ + 4, + 5 + ], + [ + 4, + 6 + ], + [ + 4, + 7 + ], + [ + 5, + 2 + ], + [ + 5, + 3 + ], + [ + 5, + 4 + ], + [ + 5, + 5 + ], + [ + 5, + 6 + ], + [ + 5, + 7 + ], + [ + 6, + 4 + ], + [ + 6, + 5 + ], + [ + 6, + 6 + ], + [ + 6, + 7 + ], + [ + 6, + 8 + ], + [ + 6, + 9 + ], + [ + 7, + 4 + ], + [ + 7, + 5 + ], + [ + 7, + 6 + ], + [ + 7, + 7 + ], + [ + 7, + 8 + ], + [ + 7, + 9 + ], + [ + 8, + 6 + ], + [ + 8, + 7 + ], + [ + 8, + 8 + ], + [ + 8, + 9 + ], + [ + 8, + 10 + ], + [ + 8, + 11 + ], + [ + 9, + 6 + ], + [ + 9, + 7 + ], + [ + 9, + 8 + ], + [ + 9, + 9 + ], + [ + 9, + 10 + ], + [ + 9, + 11 + ], + [ + 10, + 8 + ], + [ + 10, + 9 + ], + [ + 10, + 10 + ], + [ + 10, + 11 + ], + [ + 10, + 12 + ], + [ + 10, + 13 + ], + [ + 11, + 8 + ], + [ + 11, + 9 + ], + [ + 11, + 10 + ], + [ + 11, + 11 + ], + [ + 11, + 12 + ], + [ + 11, + 13 + ], + [ + 12, + 10 + ], + [ + 12, + 11 + ], + [ + 12, + 12 + ], + [ + 12, + 13 + ], + [ + 13, + 10 + ], + [ + 13, + 11 + ], + [ + 13, + 12 + ], + [ + 13, + 13 + ] + ] + ] + }, + "parallel_storage": { + "parallel_solution_storage": { + "jacobian": { + "2": [ + [ + 1.0, + 0.0, + 0.0, + 0.0, + 0.0, + 1.0, + 0.0, + 0.0, + -2.0, + 0.0, + 4.0, + 0.0, + -2.0, + 0.0, + 0.0, + -4.0, + 0.0, + 8.0, + 0.0, + -4.0, + -2.0, + 0.0, + 4.0, + 0.0, + -2.0, + 0.0, + 0.0, + -4.0, + 0.0, + 8.0, + 0.0, + -4.0, + -2.0, + 0.0, + 4.0, + 0.0, + -2.0, + 0.0, + 0.0, + -4.0, + 0.0, + 8.0, + 0.0, + -4.0, + -2.0, + 0.0, + 4.0, + 0.0, + -2.0, + 0.0, + 0.0, + -4.0, + 0.0, + 8.0, + 0.0, + -4.0, + -2.0, + 0.0, + 4.0, + 0.0, + -2.0, + 0.0, + 0.0, + -4.0, + 0.0, + 8.0, + 0.0, + -4.0, + 0.0, + 0.0, + 1.0, + 0.0, + 0.0, + 0.0, + 0.0, + 1.0 + ] + ], + "3": [ + [ + 1.0, + 0.0, + 0.0, + 0.0, + 0.0, + 1.0, + 0.0, + 0.0, + -2.0, + 0.0, + 4.0, + 0.0, + -2.0, + 0.0, + 0.0, + -4.0, + 0.0, + 8.0, + 0.0, + -4.0, + -2.0, + 0.0, + 4.0, + 0.0, + -2.0, + 0.0, + 0.0, + -4.0, + 0.0, + 8.0, + 0.0, + -4.0, + -2.0, + 0.0, + 4.0, + 0.0, + -2.0, + 0.0, + 0.0, + -4.0, + 0.0, + 8.0, + 0.0, + -4.0, + -2.0, + 0.0, + 4.0, + 0.0, + -2.0, + 0.0, + 0.0, + -4.0, + 0.0, + 8.0, + 0.0, + -4.0, + -2.0, + 0.0, + 4.0, + 0.0, + -2.0, + 0.0, + 0.0, + -4.0, + 0.0, + 8.0, + 0.0, + -4.0, + 0.0, + 0.0, + 1.0, + 0.0, + 0.0, + 0.0, + 0.0, + 1.0 + ] + ] + }, + "residual": { + "2": [ + [ + 0.0, + 0.0, + -3.552713678800501e-15, + 0.0, + 0.0, + 8.881784197001252e-16, + -3.552713678800501e-15, + -1.7763568394002505e-15, + 0.0, + 0.0, + -3.552713678800501e-15, + -3.552713678800501e-15, + 0.0, + 0.0 + ] + ], + "3": [ + [ + 0.0, + 0.0, + 3.552713678800501e-15, + 0.0, + 0.0, + -8.881784197001252e-16, + 0.0, + 8.881784197001252e-16, + 0.0, + 0.0, + 3.552713678800501e-15, + 0.0, + 0.0, + 0.0 + ] + ] + }, + "solution": { + "2": [ + [ + 0.0, + 0.0, + 8.206340405803251, + 1.5569165456566985, + 13.130144649285203, + 2.863833091313397, + 14.771412730445853, + 3.920749636970095, + 13.130144649285203, + 4.727666182626794, + 8.206340405803251, + 5.2845827282834925, + 0.0, + 5.591499273940192 + ] + ], + "3": [ + [ + 0.0, + 0.0, + 8.696932561755883, + 0.8105539617363589, + 13.915092098809412, + 1.3711079234727177, + 15.654478611160588, + 1.6816618852090768, + 13.915092098809412, + 1.7422158469454356, + 8.696932561755883, + 1.5527698086817945, + 0.0, + 1.1133237704181533 + ] + ] + } + } + }, + "time": 2.0, + "time_step": 2 + } + ] +} diff --git a/modules/stochastic_tools/test/tests/reporters/res_jac_storage/parallel_storage_main.i b/modules/stochastic_tools/test/tests/reporters/res_jac_storage/parallel_storage_main.i new file mode 100644 index 000000000000..38a5ddf9c3e8 --- /dev/null +++ b/modules/stochastic_tools/test/tests/reporters/res_jac_storage/parallel_storage_main.i @@ -0,0 +1,93 @@ +[Mesh] + [gmg] + type = GeneratedMeshGenerator + dim = 1 + [] +[] + +[Problem] + solve = false +[] +[Executioner] + type = Steady +[] + +[Distributions] + [S_dist] + type = Uniform + lower_bound = 0 + upper_bound = 10 + [] + [D_dist] + type = Uniform + lower_bound = 0 + upper_bound = 10 + [] +[] + +[Samplers] + [sample] + type = MonteCarlo + num_rows = 4 + distributions = 'S_dist D_dist' + execute_on = PRE_MULTIAPP_SETUP + seed = 0 + [] +[] + +[MultiApps] + [worker] + type = SamplerFullSolveMultiApp + input_files = sub.i + sampler = sample + mode = batch-reset + [] +[] + +[Transfers] + [param_transfer] + type = SamplerParameterTransfer + to_multi_app = worker + sampler = sample + parameters = 'Kernels/source_u/value BCs/right_v/value' + [] + [solution_transfer] + type = SerializedSnapshotTransfer + parallel_storage = parallel_storage + from_multi_app = worker + sampler = sample + solution_container = solution_storage + residual_container = residual_storage + jacobian_container = jacobian_storage + serialize_on_root = true + [] + # We are only transferring the one copy of the indices. We assume that the + # order, number, anything does not change from one run to another. Breaking + # this assumption would not allow M/DEIM to work in the current configuration. + [jac_indices] + type = MultiAppCloneReporterTransfer + from_multi_app = worker + from_reporters = 'jacobian_storage/indices' + to_reporter = jac_indices + [] +[] + +[Reporters] + [parallel_storage] + type = ParallelSolutionStorage + variables = 'solution residual jacobian' + outputs = out + [] + [jac_indices] + type = ConstantReporter + outputs = out + [] +[] + +[Outputs] + [out] + type = JSON + execute_on = FINAL + execute_system_information_on = none + [] +[] diff --git a/modules/stochastic_tools/test/tests/reporters/res_jac_storage/sub.i b/modules/stochastic_tools/test/tests/reporters/res_jac_storage/sub.i new file mode 100644 index 000000000000..d94130ccf68b --- /dev/null +++ b/modules/stochastic_tools/test/tests/reporters/res_jac_storage/sub.i @@ -0,0 +1,119 @@ +[Mesh] + type = GeneratedMesh + dim = 1 + nx = 6 + xmax = 6 +[] + +[Variables] + [u] + [] + [v] + [] +[] + +[Kernels] + [diffusion_u] + type = MatDiffusion + variable = u + diffusivity = D_u + [] + [source_u] + type = BodyForce + variable = u + value = 1.0 + [] + [diffusion_v] + type = MatDiffusion + variable = v + diffusivity = D_v + [] + [source_v] + type = BodyForce + variable = v + value = 1.0 + [] +[] + +[Materials] + [diffusivity_u] + type = GenericConstantMaterial + prop_names = D_u + prop_values = 2.0 + [] + [diffusivity_v] + type = GenericConstantMaterial + prop_names = D_v + prop_values = 4.0 + [] +[] + +[BCs] + [left_u] + type = DirichletBC + variable = u + boundary = left + value = 0 + [] + [right_u] + type = DirichletBC + variable = u + boundary = right + value = 0 + [] + [left_v] + type = DirichletBC + variable = v + boundary = left + value = 0 + [] + [right_v] + type = DirichletBC + variable = v + boundary = right + value = 0 + [] +[] + +[Executioner] + type = Transient + dt = 0.01 + num_steps = 1 + solve_type = NEWTON + petsc_options_iname = '-pc_type -pc_hypre_type' + petsc_options_value = 'hypre boomeramg' +[] + +[Controls] + [stochastic] + type = SamplerReceiver + [] +[] + +[Reporters] + [solution_storage] + type = SolutionContainer + execute_on = 'TIMESTEP_END' + [] + [residual_storage] + type = ResidualContainer + tag_name = total + execute_on = 'TIMESTEP_END' + [] + [jacobian_storage] + type = JacobianContainer + tag_name = total + jac_indices_reporter_name = indices + execute_on = 'TIMESTEP_END' + [] +[] + +[Problem] + extra_tag_vectors = 'total' + extra_tag_matrices = 'total' +[] + +[GlobalParams] + extra_matrix_tags = 'total' + extra_vector_tags = 'total' +[] diff --git a/modules/stochastic_tools/test/tests/reporters/res_jac_storage/tests b/modules/stochastic_tools/test/tests/reporters/res_jac_storage/tests new file mode 100644 index 000000000000..7e3b41f93e8e --- /dev/null +++ b/modules/stochastic_tools/test/tests/reporters/res_jac_storage/tests @@ -0,0 +1,13 @@ +[Tests] + issues = '#' + design = '' + + [p] + type = JSONDiff + requirement = 'The system should be able to print serialized snapshots of the residual and jacobian in a distributed fashion in a json format.' + input = parallel_storage_main.i + jsondiff = 'parallel_storage_main_out.json parallel_storage_main_out.json.1' + min_parallel = 2 + max_parallel = 2 + [] +[] diff --git a/modules/stochastic_tools/test/tests/reporters/res_jac_storage_non_linear_time/gold/parallel_storage_main_out.json b/modules/stochastic_tools/test/tests/reporters/res_jac_storage_non_linear_time/gold/parallel_storage_main_out.json new file mode 100644 index 000000000000..6058b58a2357 --- /dev/null +++ b/modules/stochastic_tools/test/tests/reporters/res_jac_storage_non_linear_time/gold/parallel_storage_main_out.json @@ -0,0 +1,726 @@ +{ + "number_of_parts": 2, + "part": 0, + "reporters": { + "jac_indices": { + "type": "ConstantReporter", + "values": { + "jac_indices:jacobian_storage:indices": { + "type": "std::vector>>" + } + } + }, + "parallel_storage": { + "type": "ParallelSolutionStorage", + "values": { + "parallel_solution_storage": { + "type": "std::map>, std::hash, std::equal_to, std::allocator>>>>, std::less, std::allocator>, std::hash, std::equal_to, std::allocator>>>>>>>" + } + } + } + }, + "time_steps": [ + { + "jac_indices": { + "jac_indices:jacobian_storage:indices": [ + [ + [ + 0, + 0 + ], + [ + 0, + 1 + ], + [ + 1, + 0 + ], + [ + 1, + 1 + ], + [ + 1, + 2 + ], + [ + 2, + 1 + ], + [ + 2, + 2 + ], + [ + 2, + 3 + ], + [ + 3, + 2 + ], + [ + 3, + 3 + ], + [ + 3, + 4 + ], + [ + 4, + 3 + ], + [ + 4, + 4 + ], + [ + 4, + 5 + ], + [ + 5, + 4 + ], + [ + 5, + 5 + ], + [ + 5, + 6 + ], + [ + 6, + 5 + ], + [ + 6, + 6 + ] + ] + ] + }, + "parallel_storage": { + "parallel_solution_storage": { + "jacobian": { + "0": [ + [ + 2.3333333333333335, + 1.1666666666666667, + 1.1666666666666667, + 4.666666666666667, + 1.1666666666666667, + 1.1666666666666667, + 4.666666666666667, + 1.1666666666666667, + 1.1666666666666667, + 4.666666666666667, + 1.1666666666666667, + 1.1666666666666667, + 4.666666666666667, + 1.1666666666666667, + 1.1666666666666667, + 5.0840577305704, + 1.633513801893624, + 1.633513801893624, + 3.7833308103374295 + ], + [ + 2.3544168928218334, + 1.1878772285707466, + 1.1878772285707466, + 4.812014043012618, + 1.2055826025438514, + 1.2055826025438514, + 4.80730004762024, + 1.1980659118805859, + 1.1980659118805859, + 4.802678342190141, + 1.2161276062003816, + 1.2161276062003816, + 4.830927377570498, + 1.1411886904026562, + 1.1411886904026562, + 4.733248783254279, + 1.558598484894702, + 1.558598484894702, + 3.6891766115020865 + ], + [ + 2.3543952787351587, + 1.1878553542258512, + 1.1878553542258512, + 4.811865934425313, + 1.2055446381889112, + 1.2055446381889112, + 4.807162703001351, + 1.1980298892878685, + 1.1980298892878685, + 4.802520210995715, + 1.2160861905834386, + 1.2160861905834386, + 4.830645923890426, + 1.1410397796775196, + 1.1410397796775196, + 4.732418691104339, + 1.5584397335151732, + 1.5584397335151732, + 3.6889735997627464 + ], + [ + 2.3543952787351587, + 1.1878553542258512, + 1.1878553542258512, + 4.811865934425313, + 1.2055446381889112, + 1.2055446381889112, + 4.807162703001351, + 1.1980298892878685, + 1.1980298892878685, + 4.802520210995715, + 1.2160861905834386, + 1.2160861905834386, + 4.830645923890426, + 1.1410397796775196, + 1.1410397796775196, + 4.732418691104339, + 1.5584397335151732, + 1.5584397335151732, + 3.6889735997627464 + ], + [ + 2.3543952787351587, + 1.1878553542258512, + 1.1878553542258512, + 4.811865934425313, + 1.2055446381889112, + 1.2055446381889112, + 4.807162703001351, + 1.1980298892878685, + 1.1980298892878685, + 4.802520210995715, + 1.2160861905834386, + 1.2160861905834386, + 4.830645923890426, + 1.1410397796775196, + 1.1410397796775196, + 4.732418691104339, + 1.5584397335151732, + 1.5584397335151732, + 3.6889735997627464 + ], + [ + 2.3543952787351587, + 1.1878553542258512, + 1.1878553542258512, + 4.811865934425313, + 1.2055446381889112, + 1.2055446381889112, + 4.807162703001351, + 1.1980298892878685, + 1.1980298892878685, + 4.802520210995715, + 1.2160861905834386, + 1.2160861905834386, + 4.830645923890426, + 1.1410397796775196, + 1.1410397796775196, + 4.732418691104339, + 1.5584397335151732, + 1.5584397335151732, + 3.6889735997627464 + ], + [ + 2.354395278739353, + 1.1878553542300962, + 1.1878553542300962, + 4.811865934770823, + 1.2055446385123894, + 1.2055446385123894, + 4.807162704225325, + 1.19802988891571, + 1.19802988891571, + 4.802520207280105, + 1.2160861900061197, + 1.2160861900061197, + 4.830645923162296, + 1.1410397790984006, + 1.1410397790984006, + 4.7324186868437925, + 1.558439732659978, + 1.558439732659978, + 3.688973598669079 + ], + [ + 2.373890148333171, + 1.2076934090490368, + 1.2076934090490368, + 4.947917352689835, + 1.2413813790690353, + 1.2413813790690353, + 4.940523413455807, + 1.232903174382027, + 1.232903174382027, + 4.943322796715252, + 1.245105573194195, + 1.245105573194195, + 4.941548816441129, + 1.1827971129347343, + 1.1827971129347343, + 4.946138155837161, + 1.5983345307934678, + 1.5983345307934678, + 3.7395200342176262 + ], + [ + 2.3738721110201246, + 1.2076749538638332, + 1.2076749538638332, + 4.9477895989003216, + 1.2413471541657235, + 1.2413471541657235, + 4.940393427980927, + 1.2328688380305373, + 1.2328688380305373, + 4.943187524850628, + 1.2450806307860733, + 1.2450806307860733, + 4.9414449349410905, + 1.182743552668974, + 1.182743552668974, + 4.945830083124296, + 1.5982758981762584, + 1.5982758981762584, + 3.73944642994821 + ], + [ + 2.3738721110201246, + 1.2076749538638332, + 1.2076749538638332, + 4.9477895989003216, + 1.2413471541657235, + 1.2413471541657235, + 4.940393427980927, + 1.2328688380305373, + 1.2328688380305373, + 4.943187524850628, + 1.2450806307860733, + 1.2450806307860733, + 4.9414449349410905, + 1.182743552668974, + 1.182743552668974, + 4.945830083124296, + 1.5982758981762584, + 1.5982758981762584, + 3.73944642994821 + ] + ], + "1": [ + [ + 2.3333333333333335, + 1.1666666666666667, + 1.1666666666666667, + 4.666666666666667, + 1.1666666666666667, + 1.1666666666666667, + 4.666666666666667, + 1.1666666666666667, + 1.1666666666666667, + 4.666666666666667, + 1.1666666666666667, + 1.1666666666666667, + 4.666666666666667, + 1.1666666666666667, + 1.1666666666666667, + 10.137214434011078, + 12.480891399536972, + 12.480891399536972, + 42.11968449747014 + ], + [ + 2.3425465602332274, + 1.175904145925543, + 1.175904145925543, + 4.736216241555631, + 1.1899623481988801, + 1.1899623481988801, + 4.742241183070809, + 1.1630756836269476, + 1.1630756836269476, + 4.7056034180327515, + 1.2852187619971283, + 1.2852187619971283, + 4.9434398940592015, + 0.869482700607402, + 0.869482700607402, + 6.080460822704838, + 9.860563680844733, + 9.860563680844733, + 34.626321165146074 + ], + [ + 2.3425159184375057, + 1.175873343079298, + 1.175873343079298, + 4.736142599309302, + 1.1900435934569322, + 1.1900435934569322, + 4.742344699509202, + 1.1626596524845256, + 1.1626596524845256, + 4.70458125332083, + 1.286916761956245, + 1.286916761956245, + 4.941767474781399, + 0.8617371491496417, + 0.8617371491496417, + 5.9944690159466045, + 9.795210595438295, + 9.795210595438295, + 34.427943420954904 + ], + [ + 2.342515904202989, + 1.175873328770091, + 1.175873328770091, + 4.7361425851105245, + 1.1900436512674153, + 1.1900436512674153, + 4.742344842264675, + 1.1626594542544806, + 1.1626594542544806, + 4.7045808600573675, + 1.286917756729955, + 1.286917756729955, + 4.941766527556073, + 0.8617327591506659, + 0.8617327591506659, + 5.9944203897980755, + 9.795173465562753, + 9.795173465562753, + 34.42783055395174 + ], + [ + 2.342515904202989, + 1.175873328770091, + 1.175873328770091, + 4.7361425851105245, + 1.1900436512674153, + 1.1900436512674153, + 4.742344842264675, + 1.1626594542544806, + 1.1626594542544806, + 4.7045808600573675, + 1.286917756729955, + 1.286917756729955, + 4.941766527556073, + 0.8617327591506659, + 0.8617327591506659, + 5.9944203897980755, + 9.795173465562753, + 9.795173465562753, + 34.42783055395174 + ], + [ + 2.342515904202989, + 1.175873328770091, + 1.175873328770091, + 4.7361425851105245, + 1.1900436512674153, + 1.1900436512674153, + 4.742344842264675, + 1.1626594542544806, + 1.1626594542544806, + 4.7045808600573675, + 1.286917756729955, + 1.286917756729955, + 4.941766527556073, + 0.8617327591506659, + 0.8617327591506659, + 5.9944203897980755, + 9.795173465562753, + 9.795173465562753, + 34.42783055395174 + ], + [ + 2.342515904202989, + 1.175873328770091, + 1.175873328770091, + 4.7361425851105245, + 1.1900436512674153, + 1.1900436512674153, + 4.742344842264675, + 1.1626594542544806, + 1.1626594542544806, + 4.7045808600573675, + 1.286917756729955, + 1.286917756729955, + 4.941766527556073, + 0.8617327591506659, + 0.8617327591506659, + 5.9944203897980755, + 9.795173465562753, + 9.795173465562753, + 34.42783055395174 + ], + [ + 2.342515904201445, + 1.1758733287685388, + 1.1758733287685388, + 4.736142585108242, + 1.190043651272942, + 1.190043651272942, + 4.742344842306106, + 1.162659454262084, + 1.162659454262084, + 4.704580860068304, + 1.286917756730937, + 1.286917756730937, + 4.941766527556064, + 0.8617327591491883, + 0.8617327591491883, + 5.994420389782454, + 9.795173465550862, + 9.795173465550862, + 34.427830553915584 + ], + [ + 2.352912562065284, + 1.1863554211776122, + 1.1863554211776122, + 4.803676127442335, + 1.2049030647161327, + 1.2049030647161327, + 4.7980231399261335, + 1.1861200752737562, + 1.1861200752737562, + 4.803627235271355, + 1.2833698519766013, + 1.2833698519766013, + 4.857216977975279, + 0.8626092574042399, + 0.8626092574042399, + 6.159426373024409, + 9.927880761312405, + 9.927880761312405, + 34.830071311738195 + ], + [ + 2.3529078970292874, + 1.1863507039564247, + 1.1863507039564247, + 4.803645596249034, + 1.2048962749947316, + 1.2048962749947316, + 4.797993223488914, + 1.1861051412243835, + 1.1861051412243835, + 4.803544588565185, + 1.283354128950204, + 1.283354128950204, + 4.857161933724973, + 0.8625872935294466, + 0.8625872935294466, + 6.159234928202109, + 9.927739884366117, + 9.927739884366117, + 34.82964554082572 + ], + [ + 2.3529078970292874, + 1.1863507039564247, + 1.1863507039564247, + 4.803645596249034, + 1.2048962749947316, + 1.2048962749947316, + 4.797993223488914, + 1.1861051412243835, + 1.1861051412243835, + 4.803544588565185, + 1.283354128950204, + 1.283354128950204, + 4.857161933724973, + 0.8625872935294466, + 0.8625872935294466, + 6.159234928202109, + 9.927739884366117, + 9.927739884366117, + 34.82964554082572 + ] + ] + }, + "residual": { + "0": [ + [ + 0.0, + 0.0, + 0.0, + 0.0, + 0.0, + 2.210595497826724, + 4.792111530577633 + ], + [ + 0.1057353034814488, + 0.5136846853180381, + 0.5273714051116921, + 0.5421796506777693, + 0.47060918543383745, + 1.0825898963790923, + 4.36943774099197 + ], + [ + 0.10562658240252507, + 0.5131648171001897, + 0.5268430761953242, + 0.5415907271675529, + 0.46942973537845956, + 1.0797455107425784, + 4.368533333194797 + ], + [ + 0.10562658242362297, + 0.5131648187832683, + 0.526843079133558, + 0.5415907155048361, + 0.46942973066703914, + 1.0797454965054272, + 4.368533328322642 + ], + [ + 0.20395889345551788, + 0.992480352019768, + 1.0370199172671712, + 1.0533288607286826, + 0.9236287564251432, + 1.8181744989134063, + 4.594636412527734 + ], + [ + 0.20386766220989427, + 0.9920292673246955, + 1.0365235504429684, + 1.0528424841680946, + 0.9231727959903445, + 1.8171238349238175, + 4.59430582031117 + ] + ], + "1": [ + [ + 0.0, + 0.0, + 0.0, + 0.0, + 0.0, + 41.961931250536786, + 127.75143974251777 + ], + [ + 0.046126765396925465, + 0.25520683920013565, + 0.23819803724159075, + 0.3847446591420674, + 0.24535339165932957, + 24.52626801039243, + 102.46721211497702 + ], + [ + 0.04597315379200863, + 0.25514883961382906, + 0.23761986362664958, + 0.3853941694040008, + 0.2260534647182142, + 24.12854190133635, + 101.807885040983 + ], + [ + 0.045973082432700485, + 0.2551489128700742, + 0.23761986946642757, + 0.38539517760450703, + 0.22604260859173286, + 24.12831653627873, + 101.80751004878623 + ], + [ + 0.04597308242495893, + 0.25514891287430436, + 0.2376198696028292, + 0.3853951776533107, + 0.22604260859047343, + 24.128316536206256, + 101.8075100486661 + ], + [ + 0.09816995810723955, + 0.4873365333401991, + 0.4726156997900545, + 0.6827929063042816, + 0.007990218390300452, + 24.87479097935264, + 103.1448801826265 + ], + [ + 0.09814650246427967, + 0.4872314380004758, + 0.47248659927007375, + 0.6825096468494267, + 0.0077583905115601315, + 24.873905265244172, + 103.14346356297958 + ] + ] + }, + "solution": { + "0": [ + [ + 0.0, + 0.08972553906448123, + 0.0741702312038282, + 0.05846038498952314, + 0.14870306758471305, + -0.26547267155707605, + 1.5982462231328003 + ], + [ + 0.0, + 0.17174834221430837, + 0.13844731666450286, + 0.13751846608738863, + 0.18764583332945178, + -0.12241930589790913, + 1.5982462231328003 + ] + ], + "1": [ + [ + 0.0, + 0.03925089682303096, + 0.0599313889256762, + -0.0777667168961355, + 0.5549317116846882, + -2.344961299583184, + 9.364851558334664 + ], + [ + 0.0, + 0.08342441687944108, + 0.07778854759746291, + 0.004654511931757624, + 0.4649192743751917, + -2.2094863720435023, + 9.364851558334664 + ] + ] + } + } + }, + "time": 2.0, + "time_step": 2 + } + ] +} diff --git a/modules/stochastic_tools/test/tests/reporters/res_jac_storage_non_linear_time/gold/parallel_storage_main_out.json.1 b/modules/stochastic_tools/test/tests/reporters/res_jac_storage_non_linear_time/gold/parallel_storage_main_out.json.1 new file mode 100644 index 000000000000..bcee343e5371 --- /dev/null +++ b/modules/stochastic_tools/test/tests/reporters/res_jac_storage_non_linear_time/gold/parallel_storage_main_out.json.1 @@ -0,0 +1,726 @@ +{ + "number_of_parts": 2, + "part": 1, + "reporters": { + "jac_indices": { + "type": "ConstantReporter", + "values": { + "jac_indices:jacobian_storage:indices": { + "type": "std::vector>>" + } + } + }, + "parallel_storage": { + "type": "ParallelSolutionStorage", + "values": { + "parallel_solution_storage": { + "type": "std::map>, std::hash, std::equal_to, std::allocator>>>>, std::less, std::allocator>, std::hash, std::equal_to, std::allocator>>>>>>>" + } + } + } + }, + "time_steps": [ + { + "jac_indices": { + "jac_indices:jacobian_storage:indices": [ + [ + [ + 0, + 0 + ], + [ + 0, + 1 + ], + [ + 1, + 0 + ], + [ + 1, + 1 + ], + [ + 1, + 2 + ], + [ + 2, + 1 + ], + [ + 2, + 2 + ], + [ + 2, + 3 + ], + [ + 3, + 2 + ], + [ + 3, + 3 + ], + [ + 3, + 4 + ], + [ + 4, + 3 + ], + [ + 4, + 4 + ], + [ + 4, + 5 + ], + [ + 5, + 4 + ], + [ + 5, + 5 + ], + [ + 5, + 6 + ], + [ + 6, + 5 + ], + [ + 6, + 6 + ] + ] + ] + }, + "parallel_storage": { + "parallel_solution_storage": { + "jacobian": { + "2": [ + [ + 2.3333333333333335, + 1.1666666666666667, + 1.1666666666666667, + 4.666666666666667, + 1.1666666666666667, + 1.1666666666666667, + 4.666666666666667, + 1.1666666666666667, + 1.1666666666666667, + 4.666666666666667, + 1.1666666666666667, + 1.1666666666666667, + 4.666666666666667, + 1.1666666666666667, + 1.1666666666666667, + 6.737909560968637, + 4.339815493040649, + 4.339815493040649, + 12.954685744527293 + ], + [ + 2.350218643817411, + 1.1836334372664614, + 1.1836334372664614, + 4.78583679647457, + 1.2007577771579823, + 1.2007577771579823, + 4.7858239453883815, + 1.1836449057501925, + 1.1836449057501925, + 4.765387825846486, + 1.250328778862542, + 1.250328778862542, + 4.884622710404396, + 0.995868641006004, + 0.995868641006004, + 4.889394584965988, + 3.710044180763676, + 3.710044180763676, + 11.6358263763302 + ], + [ + 2.350198249900927, + 1.183612846697581, + 1.183612846697581, + 4.785725436428283, + 1.2007499629591258, + 1.2007499629591258, + 4.7857601698372765, + 1.183538326139575, + 1.183538326139575, + 4.765022551429736, + 1.2505660360446837, + 1.2505660360446837, + 4.883151675175163, + 0.9933915432112888, + 0.9933915432112888, + 4.8691343785539685, + 3.7028534261384136, + 3.7028534261384136, + 11.61985610343732 + ], + [ + 2.3501982506131434, + 1.1836128474166614, + 1.1836128474166614, + 4.785725446452288, + 1.2007499693671146, + 1.2007499693671146, + 4.785760192812833, + 1.1835383201466856, + 1.1835383201466856, + 4.76502253722391, + 1.2505660769321731, + 1.2505660769321731, + 4.883151504283087, + 0.9933912117104229, + 0.9933912117104229, + 4.8691316428277105, + 3.702852453323218, + 3.702852453323218, + 11.619853941372133 + ], + [ + 2.3501982506131434, + 1.1836128474166614, + 1.1836128474166614, + 4.785725446452288, + 1.2007499693671146, + 1.2007499693671146, + 4.785760192812833, + 1.1835383201466856, + 1.1835383201466856, + 4.76502253722391, + 1.2505660769321731, + 1.2505660769321731, + 4.883151504283087, + 0.9933912117104229, + 0.9933912117104229, + 4.8691316428277105, + 3.702852453323218, + 3.702852453323218, + 11.619853941372133 + ], + [ + 2.3501982506131434, + 1.1836128474166614, + 1.1836128474166614, + 4.785725446452288, + 1.2007499693671146, + 1.2007499693671146, + 4.785760192812833, + 1.1835383201466856, + 1.1835383201466856, + 4.76502253722391, + 1.2505660769321731, + 1.2505660769321731, + 4.883151504283087, + 0.9933912117104229, + 0.9933912117104229, + 4.8691316428277105, + 3.702852453323218, + 3.702852453323218, + 11.619853941372133 + ], + [ + 2.3501982506131434, + 1.1836128474166614, + 1.1836128474166614, + 4.785725446452288, + 1.2007499693671146, + 1.2007499693671146, + 4.785760192812833, + 1.1835383201466856, + 1.1835383201466856, + 4.76502253722391, + 1.2505660769321731, + 1.2505660769321731, + 4.883151504283087, + 0.9933912117104229, + 0.9933912117104229, + 4.8691316428277105, + 3.702852453323218, + 3.702852453323218, + 11.619853941372133 + ], + [ + 2.350198250612853, + 1.1836128474163679, + 1.1836128474163679, + 4.78572544645094, + 1.2007499693672434, + 1.2007499693672434, + 4.785760192815085, + 1.1835383201471106, + 1.1835383201471106, + 4.765022537224349, + 1.250566076932175, + 1.250566076932175, + 4.883151504283074, + 0.9933912117104162, + 0.9933912117104162, + 4.8691316428276705, + 3.702852453323205, + 3.702852453323205, + 11.619853941372101 + ], + [ + 2.366716356679297, + 1.2003680831394088, + 1.2003680831394088, + 4.898355874102501, + 1.228896698526583, + 1.228896698526583, + 4.891014961629189, + 1.2159427504140377, + 1.2159427504140377, + 4.896388118081415, + 1.2636355474876253, + 1.2636355474876253, + 4.90564973811094, + 1.0280186865329402, + 1.0280186865329402, + 5.176645335093257, + 3.811582015014376, + 3.811582015014376, + 11.859004740763284 + ], + [ + 2.366703732496522, + 1.200355218204775, + 1.200355218204775, + 4.898269227398156, + 1.228875091291298, + 1.228875091291298, + 4.890929308639586, + 1.215912791048344, + 1.215912791048344, + 4.8962651442911085, + 1.263629251765192, + 1.263629251765192, + 4.905613284785698, + 1.027944315579063, + 1.027944315579063, + 5.175998943909416, + 3.8113569381850168, + 3.8113569381850168, + 11.858514782774266 + ], + [ + 2.366703732496522, + 1.200355218204775, + 1.200355218204775, + 4.898269227398156, + 1.228875091291298, + 1.228875091291298, + 4.890929308639586, + 1.215912791048344, + 1.215912791048344, + 4.8962651442911085, + 1.263629251765192, + 1.263629251765192, + 4.905613284785698, + 1.027944315579063, + 1.027944315579063, + 5.175998943909416, + 3.8113569381850168, + 3.8113569381850168, + 11.858514782774266 + ] + ], + "3": [ + [ + 2.3333333333333335, + 1.1666666666666667, + 1.1666666666666667, + 4.666666666666667, + 1.1666666666666667, + 1.1666666666666667, + 4.666666666666667, + 1.1666666666666667, + 1.1666666666666667, + 4.666666666666667, + 1.1666666666666667, + 1.1666666666666667, + 4.666666666666667, + 1.1666666666666667, + 1.1666666666666667, + 4.947267970585088, + 1.4696966789525714, + 1.4696966789525714, + 3.264852078558531 + ], + [ + 2.3518247482738155, + 1.1852557752284463, + 1.1852557752284463, + 4.793796515552001, + 1.2005602208839088, + 1.2005602208839088, + 4.789439816589232, + 1.1947128115782348, + 1.1947128115782348, + 4.786516460856358, + 1.2071804440006706, + 1.2071804440006706, + 4.805487621626911, + 1.155350510204914, + 1.155350510204914, + 4.73284240392091, + 1.4250450740567844, + 1.4250450740567844, + 3.2125117769926104 + ], + [ + 2.351808291977864, + 1.1852391451288642, + 1.1852391451288642, + 4.793683611716044, + 1.200530983356131, + 1.200530983356131, + 4.789334266277203, + 1.1946862550754886, + 1.1946862550754886, + 4.7864008996736835, + 1.2071478619757228, + 1.2071478619757228, + 4.805313319850939, + 1.1552839462212054, + 1.1552839462212054, + 4.732523036310873, + 1.424988990360586, + 1.424988990360586, + 3.2124453026653863 + ], + [ + 2.351808291977864, + 1.1852391451288642, + 1.1852391451288642, + 4.793683611716044, + 1.200530983356131, + 1.200530983356131, + 4.789334266277203, + 1.1946862550754886, + 1.1946862550754886, + 4.7864008996736835, + 1.2071478619757228, + 1.2071478619757228, + 4.805313319850939, + 1.1552839462212054, + 1.1552839462212054, + 4.732523036310873, + 1.424988990360586, + 1.424988990360586, + 3.2124453026653863 + ], + [ + 2.351808291977864, + 1.1852391451288642, + 1.1852391451288642, + 4.793683611716044, + 1.200530983356131, + 1.200530983356131, + 4.789334266277203, + 1.1946862550754886, + 1.1946862550754886, + 4.7864008996736835, + 1.2071478619757228, + 1.2071478619757228, + 4.805313319850939, + 1.1552839462212054, + 1.1552839462212054, + 4.732523036310873, + 1.424988990360586, + 1.424988990360586, + 3.2124453026653863 + ], + [ + 2.351808291977864, + 1.1852391451288642, + 1.1852391451288642, + 4.793683611716044, + 1.200530983356131, + 1.200530983356131, + 4.789334266277203, + 1.1946862550754886, + 1.1946862550754886, + 4.7864008996736835, + 1.2071478619757228, + 1.2071478619757228, + 4.805313319850939, + 1.1552839462212054, + 1.1552839462212054, + 4.732523036310873, + 1.424988990360586, + 1.424988990360586, + 3.2124453026653863 + ], + [ + 2.3518082919539824, + 1.1852391451047308, + 1.1852391451047308, + 4.793683611656197, + 1.2005309834175, + 1.2005309834175, + 4.789334266537436, + 1.1946862549323807, + 1.1946862549323807, + 4.786400898453074, + 1.2071478618140197, + 1.2071478618140197, + 4.805313319920181, + 1.1552839461852336, + 1.1552839461852336, + 4.732523035739838, + 1.4249889902427384, + 1.4249889902427384, + 3.212445302525703 + ], + [ + 2.368826752545427, + 1.2025200054663303, + 1.2025200054663303, + 4.912179798275157, + 1.2318402813551292, + 1.2318402813551292, + 4.9058672526088545, + 1.224824173992519, + 1.224824173992519, + 4.908067608086599, + 1.2331816473153738, + 1.2331816473153738, + 4.90638194038827, + 1.1911940392180826, + 1.1911940392180826, + 4.9065414832260155, + 1.4554858015053802, + 1.4554858015053802, + 3.2483189602809377 + ], + [ + 2.3688130304224595, + 1.202506005114011, + 1.202506005114011, + 4.912082887755403, + 1.2318142062925515, + 1.2318142062925515, + 4.9057683216936, + 1.2247985152086216, + 1.2247985152086216, + 4.9079665117648945, + 1.233161673416547, + 1.233161673416547, + 4.906299004876681, + 1.1911561690235724, + 1.1911561690235724, + 4.906337286695724, + 1.4554492692928522, + 1.4554492692928522, + 3.2482763101838605 + ], + [ + 2.3688130304224595, + 1.202506005114011, + 1.202506005114011, + 4.912082887755403, + 1.2318142062925515, + 1.2318142062925515, + 4.9057683216936, + 1.2247985152086216, + 1.2247985152086216, + 4.9079665117648945, + 1.233161673416547, + 1.233161673416547, + 4.906299004876681, + 1.1911561690235724, + 1.1911561690235724, + 4.906337286695724, + 1.4554492692928522, + 1.4554492692928522, + 3.2482763101838605 + ] + ] + }, + "residual": { + "2": [ + [ + 0.0, + 0.0, + 0.0, + 0.0, + 0.0, + 13.110979301689879, + 34.48625309391985 + ], + [ + 0.08463020270967941, + 0.4255700272475326, + 0.4255665707413909, + 0.49840377614804926, + 0.3270503256823529, + 6.488268516839167, + 29.614676392734683 + ], + [ + 0.08452774149626902, + 0.4252206152124731, + 0.42512114733994166, + 0.4978172840349864, + 0.31777313607783586, + 6.413448369759177, + 29.556773823939324 + ], + [ + 0.08452774507451127, + 0.4252206580901563, + 0.4251212058165829, + 0.4978173357569222, + 0.3177719823142059, + 6.413438269653375, + 29.556765986738377 + ], + [ + 0.08452774507305144, + 0.4252206580863758, + 0.4251212058235949, + 0.4978173357590854, + 0.31777198231415893, + 6.413438269653226, + 29.556765986738263 + ], + [ + 0.16771109954676394, + 0.819051639421229, + 0.8396360264245233, + 0.9399160399576953, + 0.493259930328764, + 7.540615091601431, + 30.426466889444146 + ], + [ + 0.16764737675324112, + 0.818748842235572, + 0.8392929774480691, + 0.939517967761611, + 0.4929671303248813, + 7.538250494183736, + 30.424679302398204 + ] + ], + "3": [ + [ + 0.0, + 0.0, + 0.0, + 0.0, + 0.0, + 1.4590782905108135, + 3.0863718937777547 + ], + [ + 0.09270130875565305, + 0.449031279160888, + 0.4617821226284367, + 0.47102429108815697, + 0.42004643958123916, + 0.7830949704565189, + 2.8438921276234854 + ], + [ + 0.09261859276681904, + 0.44863435050259454, + 0.4613787617720564, + 0.4705875418122355, + 0.4193628201196647, + 0.7819899322316597, + 2.84358573256493 + ], + [ + 0.09261859264678349, + 0.44863435044607125, + 0.46137876221828944, + 0.4705875379986856, + 0.4193628197985848, + 0.7819899304195223, + 2.8435857319211024 + ], + [ + 0.17836689502939346, + 0.8663502127415394, + 0.9063292698912562, + 0.9151835734862286, + 0.8268940673043158, + 1.383053309873694, + 3.009511904465793 + ], + [ + 0.17829758884117575, + 0.8660077479049137, + 0.9059526079869302, + 0.914816750975155, + 0.826542118292001, + 1.382356812530369, + 3.0093139486917813 + ] + ] + }, + "solution": { + "2": [ + [ + 0.0, + 0.07193161101050198, + 0.07204666966914233, + -0.000432534042346283, + 0.34371383450422044, + -1.2296909179555944, + 5.591499273940192 + ], + [ + 0.0, + 0.1416653031303965, + 0.11805867864511127, + 0.08863437725879098, + 0.3089343558950541, + -0.9986610488804344, + 5.591499273940192 + ] + ], + "3": [ + [ + 0.0, + 0.07876264080036548, + 0.06429671016359721, + 0.05436461617410629, + 0.11605731041157243, + -0.1677646188658123, + 1.1133237704181533 + ], + [ + 0.0, + 0.1505299144563643, + 0.12112810701140952, + 0.12200023363432676, + 0.1551164983810106, + -0.05268714078933466, + 1.1133237704181533 + ] + ] + } + } + }, + "time": 2.0, + "time_step": 2 + } + ] +} diff --git a/modules/stochastic_tools/test/tests/reporters/res_jac_storage_non_linear_time/parallel_storage_main.i b/modules/stochastic_tools/test/tests/reporters/res_jac_storage_non_linear_time/parallel_storage_main.i new file mode 100644 index 000000000000..c4569bcd774b --- /dev/null +++ b/modules/stochastic_tools/test/tests/reporters/res_jac_storage_non_linear_time/parallel_storage_main.i @@ -0,0 +1,82 @@ +[StochasticTools] +[] + +[Distributions] + [S_dist] + type = Uniform + lower_bound = 0 + upper_bound = 10 + [] + [D_dist] + type = Uniform + lower_bound = 0 + upper_bound = 10 + [] +[] + +[Samplers] + [sample] + type = MonteCarlo + num_rows = 4 + distributions = 'S_dist D_dist' + execute_on = PRE_MULTIAPP_SETUP + seed = 0 + [] +[] + +[MultiApps] + [worker] + type = SamplerFullSolveMultiApp + input_files = sub.i + sampler = sample + mode = batch-reset + [] +[] + +[Transfers] + [param_transfer] + type = SamplerParameterTransfer + to_multi_app = worker + sampler = sample + parameters = 'Kernels/source_u/value BCs/right_u/value' + [] + [solution_transfer] + type = SerializedSnapshotTransfer + parallel_storage = parallel_storage + from_multi_app = worker + sampler = sample + solution_container = solution_storage + residual_container = residual_storage + jacobian_container = jacobian_storage + serialize_on_root = true + [] + # We are only transferring the one copy of the indices. We assume that the + # order, number, anything does not change from one run to another. Breaking + # this assumption would not allow M/DEIM to work in the current configuration. + [jac_indices] + type = MultiAppCloneReporterTransfer + from_multi_app = worker + from_reporters = 'jacobian_storage/indices' + to_reporter = jac_indices + [] +[] + +[Reporters] + [parallel_storage] + type = ParallelSolutionStorage + variables = 'solution residual jacobian' + outputs = out + [] + [jac_indices] + type = ConstantReporter + outputs = out + [] +[] + +[Outputs] + [out] + type = JSON + execute_on = FINAL + execute_system_information_on = none + [] +[] diff --git a/modules/stochastic_tools/test/tests/reporters/res_jac_storage_non_linear_time/sub.i b/modules/stochastic_tools/test/tests/reporters/res_jac_storage_non_linear_time/sub.i new file mode 100644 index 000000000000..12e2a0a76ee7 --- /dev/null +++ b/modules/stochastic_tools/test/tests/reporters/res_jac_storage_non_linear_time/sub.i @@ -0,0 +1,106 @@ +[Mesh] + type = GeneratedMesh + dim = 1 + nx = 6 + xmax = 6 +[] + +[Variables] + [u] + [] +[] + +[Kernels] + [diffusion_u] + type = MatDiffusion + variable = u + diffusivity = D_u + [] + [source_u] + type = BodyForce + variable = u + value = 1.0 + [] + [non_linear_u] + type = ExponentialReaction + variable = u + mu1 = 7 + mu2 = 0.4 + extra_matrix_tags = 'nonlin_tag' + extra_vector_tags = 'nonlin_tag' + [] + [time_u] + type = TimeDerivative + variable = u + [] + +[] + +[Materials] + [diffusivity_u] + type = GenericConstantMaterial + prop_names = D_u + prop_values = 2.0 + [] + [diffusivity_v] + type = GenericConstantMaterial + prop_names = D_v + prop_values = 4.0 + [] +[] + +[BCs] + [left_u] + type = DirichletBC + variable = u + boundary = left + value = 0 + [] + [right_u] + type = DirichletBC + variable = u + boundary = right + value = 0 + [] +[] + +[Executioner] + type = Transient + dt = 0.01 + num_steps = 2 + solve_type = NEWTON + petsc_options_iname = '-pc_type -pc_hypre_type' + petsc_options_value = 'hypre boomeramg' + nl_rel_tol = 1e-16 + nl_abs_tol = 1e-12 +[] + +[Controls] + [stochastic] + type = SamplerReceiver + [] +[] + +[Reporters] + [solution_storage] + type = SolutionContainer + execute_on = 'TIMESTEP_END' + [] + [residual_storage] + type = ResidualContainer + tag_name = nonlin_tag + execute_on = 'NONLINEAR' + [] + [jacobian_storage] + type = JacobianContainer + tag_name = nonlin_tag + jac_indices_reporter_name = indices + execute_on = 'LINEAR TIMESTEP_END' + [] +[] + +[Problem] + extra_tag_vectors = 'nonlin_tag' + extra_tag_matrices = 'nonlin_tag' +[] + diff --git a/modules/stochastic_tools/test/tests/reporters/res_jac_storage_non_linear_time/tests b/modules/stochastic_tools/test/tests/reporters/res_jac_storage_non_linear_time/tests new file mode 100644 index 000000000000..70db0b0b8856 --- /dev/null +++ b/modules/stochastic_tools/test/tests/reporters/res_jac_storage_non_linear_time/tests @@ -0,0 +1,14 @@ +[Tests] + issues = '#' + design = '' + + [p] + type = JSONDiff + requirement = 'The system should be able to print serialized snapshots of a nonlinear unsteady residual and jacobian in a distributed fashion in a json format.' + input = parallel_storage_main.i + jsondiff = 'parallel_storage_main_out.json parallel_storage_main_out.json.1' + allow_test_objects = true + min_parallel = 2 + max_parallel = 2 + [] +[] diff --git a/modules/stochastic_tools/test/tests/surrogates/polynomial_regression/poly_reg_vec.i b/modules/stochastic_tools/test/tests/surrogates/polynomial_regression/poly_reg_vec.i index ace8ce08f546..87bd4662a96b 100644 --- a/modules/stochastic_tools/test/tests/surrogates/polynomial_regression/poly_reg_vec.i +++ b/modules/stochastic_tools/test/tests/surrogates/polynomial_regression/poly_reg_vec.i @@ -92,4 +92,5 @@ type = JSON execute_on = timestep_end [] + [] diff --git a/modules/stochastic_tools/test/tests/userobjects/rb_inverse_mapping/gold/main_rb_out.e b/modules/stochastic_tools/test/tests/userobjects/rb_inverse_mapping/gold/main_rb_out.e new file mode 100644 index 0000000000000000000000000000000000000000..465fb75aad78596ee8531838b471279c9fb67a03 GIT binary patch literal 192040 zcmeF)b(odq-?shf84wdYu`yA>ju}GI0_pDVumc0TTWl%2y9paC>_h>(yM6clJip`P zjnDJm_xsnwyZzQ|+sEsegT*y7*Q{e^O{mRVv@9tqDpFmnfvS;%M)w)rZ`7coqGFyc zsxxues6l-uOdK@U=kinK=t-mc3>-Gf=lr{Me17zpfrBQbKd-wLhbp^4BL{J#I`eqG zZhpQ`|B)lojXe2^wSTYAgh3PY-z_Su?I#Qy$jkjN{Yr6;VWU0Gj-Dr-d*wJ!pRt3+ z^D@&Pc%FLp`8*RWcUq?N`DMGt|C_E!bSu=poYCoYMc(!wou^v1^Cru7p6Jh#^?y?l;!1`iv_wak3aLz>yPSp z*sxK)zqya};guUQaVYiQf4>2_9w*7W1B`2Kx{r89GVIBtZx>Hd4@&@u5seH&E!R6(dp&y=(~ zg=3vT{r}f1&hy^ZYg_x|brMf?e8a1_J+Jjx?_agI^?hV6%R5SMXK9-jH2Lq={{4|W zUuX1$f4w|kr`jiPyGnkhjuWin_Wyi&-uAzHd9{CkGPkWW;omRM`%~@HKVRNClvUjR zpD)kb{!d@t)!epXz!>gjhmrC%DNQN=gV()%>m2Lc#xdXa@YHzy{tvh9{xg5(b@%dm z|MGrR`#!3D-RkyfIc(&}&izJC%6B$z<2ijI*T$P_c;0LJ&vVsY|Nr5+;@w$vb2=B- z>!0VUy|(`!&vm=!s_V7;`)v8`^zHEJe|Wyyo8MjI-|=H~A3unV>gRoF*=x)j(r`LG zPxsdFkE#JX^Zp`F5TBJui^GQy|;va z)JwPjuP?v$;&~(TmifL~eAH$`cYb=l`1*+)Yj0b79ac_% zmjCgO+ul&~kK5M1ZTWrruiM_}^{G8??d!X0zFGR`KhAsO>ooNomJc-gRMfB&(kq}@6E8K2(&yJ-pZ)Xi)c(1*FRH%(QlF;Q9#d)mbo8~y zg^klwb?xzu`19K1qWm|jYmbZL_O-_k2m3Uny7u@<{QcVF=kwE3b?x!%IA87Yo49@L z@w+&G?eY8c*B4aR9)FCNS9|m$Jx|rO$6w-nwa4Gm`I=PM9)FM5r}p?~I^W>xIM$@I z6&3xjV}AKHMR6?7=dYb*gkCk!0 z+M^#~J>TEQI8Rj^>%{qLk9FgGwZ|3m^{A;mx(d(t_c5JkeoehNt{CU5J+73_H@K$u zxN@AY_83d9`TH2RU#0dKujeXpTs7XlzmM_qR*hq^~Espi$e19L~ zJoV!k3#^Ib|8XA+^LeV@$bau&_x}&`%dcMce!lGe`nx#KviIL*@2|_=KmT=qjO(@R z{qA4)v$)-|_osi|f8ut_-cSB@zliI)?EC*;@9%Uy2G{)S{TkQ(U+=#-Z~c6IYLNuT|Fyi?cv8a+e1WQ$eY{)vyq_*i&9aZ{mVF#o_`L717St3TkNJRR>c74YT=wg}!sFz;t7X5w zDLhV2k7NGzb#LM8m-O+adfD%X{Nmi$TjA@-*1fKl#`9nY|lbI<@fmrSSR1JI=p;PF46mf4XDT)K9Z& zT=+hJ`Y_AS?+SY2Y^q zFNM7n_QKa4Z#W0O-@7P3zAXH_GImk;J*9MqDExds_OUR39a`8+`ucNzP2u;)(uFAe z{#f2e{_OFO-xDwFC4c;^SJ(?5ISao(R`|Vfhp?cg@N?vJ9SXm9QTVx-Uw741#C0hA z965cyEc_fff1OcN_&IWZ{Iu-9?^f7LdO%rS`1i%~1BT5^Gmu83u+3VUkaaJcqCf*xmfx-b#P7L z=VJNmlET0Dnjc>+OdlQ!KS!>YKW`R(E*6g;3jZCFw2$<`!p|R1`z&9N6!wxoeiruP z2hPIJ#R@+c^8;q#=VFC_e?Na*Ed2YG>GS0OUwT^&@MK%l3um z^A#vupL{?0PrE-x^P^+uESlim%=gcWqG*{_xrGg zukQ=T&n7$`MM6Tr-z3fstZRS04+b);SImc<@@}F(-%a9qBw!@KeQvBEJu0IFVIuEV1EJ+;CyU5LN`dy)Uw zuSa~}a#7)!zJI&CT!;K{taWwa=)=gV# zmKKhD;-6c+yj+J5^MCIo{`*LS3&-W;dVCZ=uPhw%&nuUk>)@BUt!}CQsBm1auEWRi zbL7G?{~UR_xeop~sNs**9~X|x)phv9A8YgT{lanix(=VFPj`D)7mn!%x699U$lu_v z{Y#N z5x=KYIOg9|T3)Whm+8y>)vF7~^aGFO=Q_mik9}D<`p~C(dASZ>r5|hetS%hW1EA&R zI@JDubPC7(dwt8vbx5zwfABy0kC`}e%&7m|3ow4zkf9U*bKm@;iDSn82fy`2MI#3d z{!f1Q7Zvtj*#Cd>2Y~8pdMDlN>~}I zU{$P!`tbj)Jra4aU{I2?}?a3ZGQB%F*>a4Js2={N&t;w+qvb8s%s z!}+)X7vdsZj7xASrs6VOjw^5_uEN#02G?R5reg+XVivB$^|%2y;wIdTTW~9G!|k{O zcj7MGjeBq}?!*0f01x6JJd8*1C?3OX%)wmD!{c}YPvR*&jc4#I=3@b#!}E9nFXAOE z#LIXEui`bljyLco-oo2>2k+uNEW-Qv0E@8%OYtE-!pHaopW-uojxX>fzQWh|2H)a4 ze2*XSBYwiq_yxb>H~fx2@F)H%Dk|ns%s+}tP>M2?qXLzvLLJn_3aEz_u@Y9sDp(b( zp*~i}8dwu+VQs8~b+I1S#|GFC8=(O<#wOSln_+Wofi2MxTVZQ#gKe=LwnroEfE}?D z8lwqz#xB?uyJ2_ifjzMo_QpP_Mh%*x8JeR7_C-syLTj`^TeQP|*dGVrKpcej=zxRK z5uMN(UCcO{6TQ$Ieb5*E&>sUZ5Q8unLogJ>FdQQ=5~DC0V=xxuFdh>y5tDET z4#iz;6hx4 zi*X4q#Z+8|%W(y+#8tQ&*Wg-A!*tBROw7V{xE?p)M%;v(aSLw6ZMYqG;7;6yyKxWh z#eKLR58y#Ogop769>rsrjX9W$d3YR8;7L4%r|}G)#e6Klb9f#v;6=QIg?JgS;8nba z*YO74#9Me9@8Dg$hedcFA7C+-U@1PtNB9_@;8T2t&+!Gm#8>zl-{4z(hwt$Ne#B4s z8Nc9H{D$B02mZug`O|($F~<^=q73DzKqabB2X(Om>S0B!gq5)hR>f+lkJYgT*2G#^ z8|z?QtcUfn0XD=&Xn>8e2{y%M*c@A6OEkn**c#hlTWp8z(Fi+WN9=^gXo8)w3wFhB z*d2RdPwa)gu@9Y6EOuR;bfeGQ*jzj#~C;iXW?v|gL82n&c_9~5EtQMT!Kq66_??1 zT!AZb6|TlLxE9ke9WyW!vv3`*#|^j)`S5>Mf2JcDO39}Dmtp2rJ#5iemOUdAhU6|doSyn#3I7T(4? zco*+s5#GlKSd1lDiVyJ-KE@~b6rbU9e1R|V6~4wd_!i&cd;EYO@e_W=FZdO|;dlIj zKk-+7$X{B_u>_?kLpdr?i7M1VU95n5SP?5>Wvqf#u^Q@Qb*zCku@=_GI#?I$VSQ|X z4Y3g#U}J28O|cm^#}?QU4Y3ur#x~d%+hKb&!VcIGJE1X}U}x-tU9lT>#~#=ddtqJ(o8pq&ROvZ6I z9w*>LOuY>oQBhJ2F}D;I2-5ST%3pVaRDyGMYtH3;8IM*Ww;zy;7VMDt8opk z#WYOE49vtVT!-s%18&4kxEZ(LR@{c$aR=_iUAPM;B92KZU73!caRzN+hh?TH1R>7)R4fU})*1(!r3u|K?tc&%qJ~qIH*a!`KM+Y2?j_8EW=z^~3hVJNrp6G?%=!3rKhyECVff$6r z7=ob~hT#~2kr;*17=y7Ghw+$ziI{{#a3~JL;Wz?E;wT)AV{j}c<2W3T6L2D?;3S-k zQ*bIy!|6B!XW}fJjdO4=&cpe*02ksST#QR_DW>8wT#hSnC9cBNxCYl^8m40gW?~ku z!}YiUH{vGTj9YLkZo}=k19##s+>Lv1FYd$rcmNOLAv}yn@F*U`Y|O!2%){e&0#D*8 zJdJ1YEaqbYp2PEa0WabuEX2!r1+U^YypA{UCf>r^cn9y|JuJfe_yCKs1WWNDKElWN z1fSwFe2y>hCBDMf_y*tNJA98H@FRZ0&-ewu;y3(`Kkz61Dk>^3LNQ8EiZW1^SD+GA zsDrv#0rjvVR>I0y1*>8;)W_;r18ZU}tc`WBF4n{P*Z>=1BQ(Iq*aVwmGi;76uq7H| zD{PHzur0R3_GpA1up@RtV>H3e*af>{H|&l*uqXDy-q;7#s6kUSLvysizG#V7XpJ^# zi+0!#`{Mu{h=b4`9dIx@q7yo!3%a5kx}yhrq8ECj5Bj1X`eOhFVh{#n2!>)9hGPUq zViZPW48~#{#$y5|ViFF)p*ReO;|Lsyqi{5i!LgW(<8VAqz=@cGlW;Ol!KpY6r{fHq ziL-Dv&cV4j59i|oT!@QsF)qQSn2O7AIj+E!xC&R}8eEHMn2s5kiCMS~*W(7I0y1*>8;)W_;r18ZU}tc`WBF4n{P*Z>=1BQ(Iq*aVwmGi;76uq7H|D{PHzur0R3 z_GpA1up@RtV>H3e*af>{H|&l*uqXDy-q;7#s6kUSLvysizG#V7XpJ^#i+0!#`{Mu{ zh=b4`9dIx@q7yo!3%a5kx}yhrq8ECj5Bj1X`eOhFVh{#n2!>)9hGPUqViZPW48~#{ z#$y5|ViFF)p*ReO;|Lsyqi{5i!LgW(<8VAqz=@cGlW;Ol!KpY6r{fHqiL-Dv&cV4j z59i|oT!@QsF)qQSn2O7AIj+E!xC&R}8eEHMn2s5kiCMS~*W(7I0y1*>8; z)W_;r18ZU}tc`WBF4n{P*Z>=1BQ(Iq*aVwmGi;76uq7H|D{PHzur0R3_GpA1up@Rt zV>H3e*af>{H|&l*uqXDy-q;7#s6kUSLvysizG#V7XpJ^#i+0!#`{Mu{h=b4`9dIx@ zq7yo!3%a5kx}yhrq8ECj5Bj1X`eOhFVh{#n2!>)9hGPUqViZPW48~#{#$y5|ViFF) zp*ReO;|Lsyqi{5i!LgW(<8VAqz=@cGlW;Ol!KpY6r{fHqiL-Dv&cV4j59i|oT!@Qs zF)qQSn2O7AIj+E!xC&R}8eEHMn2s5kiCMS~*W(7p5^>SJ}Rfi1WMLX<={c!*e#6f6}4mcPc(FvW=1zph% z-O&R*(F?uN2Yt~G{V@OoF$jY(1Vb?l!!ZIQF$$wG24gV}<1qmfF$ss@P#lKCaRiRU zQ8*gM;8;w?aX20);6zNpNjMp&;8dK3({TpQ#925S=ipqNhx2g(F2qH+7?v02a#7(#vx8PRXhTCxm?!;ZV8~5N|+=u(|03O6c zco>i1Q9Op(n1i{PhsW^*p2Sml8qeTa%*O&ehv)GEUc^gSh?nsSUd3y89dF=GyoI;% z4&KFkScLcS0TyEkmf}NvgpctFKE-GF9ADr|e1)&^4Zg*9_#QvtNBo4J@e6*%Z}=U5 z;7|NjR8*%3#VA23%2199RH6#ps7_t1fO=RFD`91p5^>SJ}Rfi1WMLX<={c!*e#6f6}4mcPc(FvW=1zph%-O&R*(F?uN z2Yt~G{V@OoF$jY(1Vb?l!!ZIQF$$wG24gV}<1qmfF$ss@P#lKCaRiRUQ8*gM;8;w? zaX20);6zNpNjMp&;8dK3({TpQ#925S=ipqNhx2g(F2qH+7?v02a#7(#vx8PRXhTCxm?!;ZV8~5N|+=u(|03O6cco>i1Q9Op( zn1i{PhsW^*p2Sml8qeTa%*O&ehv)GEUc^gSh?nsSUd3y89dF=GyoI;%4&KFkScLcS z0TyEkmf}NvgpctFKE-GF9ADr|e1)&^4Zg*9_#QvtNBo4J@e6*%Z}=U5;7|NjR8+SJ z#VA23%2199RH6!XzzyrJfO=RFD`91p5^>SJ}RfiXIrOwRtD&l&j z&g^yDIy&<<*fu)zcGx~Tb0h2!oq0#>6rH&-nnY*b8M{Ph-W9t=XWkuqL}%U;dqrp7 z8~a3Ou0{>)%uUe@HIbWhYyms-zGw+Mb1Srlow*I#!p__d`@znXYPjXurv2SPuQ7zp*QTzeb5(n=6>jpzL5uT90)t}APk0G z8A>q(LnBvk92T9q3d5r_*TsnF%x+_3bY{0PDmt^<7#*G2ZH$S|>^8z+!g1+&fFd6!_M3j7r@Tk8yCXP+!q(Y&fFgt!_GVqm%z?E7?;9sDBQ+WTo&1F zTppcyIIf7!>^829&O8!VMQ3&!SHsRc8rR_J$Zq3W*qO&-8tlw&V>;~2<1qtvX16gD zcIJtg1v|4lxDIyaLvTIp%3uft~p%+zLChJGc#Y z=3{U>?9A@q4%nF|<4)L_-N9Y3Garw;VP`%8_rT74BJPEqc?#}>o%tl(4?D9vcmQ_h zQ}7_{%%|cZ*qKkm!>}`-jz?f;J_C=!&U_{wgPqwO%*JDp&*nG>cII<17k1`zF%Nd; z^YA$A%;)0?*qJZDldv;ih^Js@z6ejl&U`VRft~phJPSMXrI-)9%P39ux<{5Y$cIKIQ19s+F zcoTN!>+lxr%-7>>*qLv@JFqj~h<9OUb_egl&g>2r!OrXs-p8WI?%)I1nccx+*qPnI z64;sVz*5+m@5G0&Gv9@eU}wG?AH&Xk4?cmN`CfbqJM(?`40dLB@Hy7+=HA{0P2*o%vCG3p?{;_zrgF+4vrI<~jHQcILVG5q9Qz_z8CA$MG}l z%unDK*qNWiudp*eh2LOjej2~S&ioAifSvhS{0TereEbDFH?bgp8B$!Fv)d>x$z9sX z3*^$=rJcM;F3Vlo$xGz&+@(z{Bv<4vZQ^BeW$w}@ULjZIE^Xpfa-H0zO}s{~o4d4$ z*U2m7E^Xosa=qN8O}t57F?VSbZ;@BZUE0Lk*P`q{S(k9qJ@iw_jn_vgU+vYB9f*llZm%FqHc2K;1?$RdML2;wp zrA@Gd;vI6AHo*>xcg$Vd1luUyDLOXcEQ%YWNo2RNb9836u}gGjx3OzEu}5@fx3On*X1B3dbY{1)cXVdAu}^enw^1FP*=^KBXLcJ+qcgjWX3?45M)T;* zZlgtXX1B3#bY{2FGCH%{Xce8=ZM2Tg>^9m&XLcKHqcgjWcF~#L#(vS6-NydWncc<# z(V5-Gfzg@W#zE1U-A4Q9%x7v)kwoJG0vu z06VkW7zjJF+ZY79CTNVo7!rABjzgm}?}A~`nRmtT=*+ueM0DoeF)}*y9vBs!c~6Xv z&b$}KL}%U`W1}2921?{Z5$h&*=^4q}&g?d(L}zvfCq-v=2Pa2ob_b_KXLbjtMrU>hr$uLW2d778b_ZufXLbi?MrU>h zXGLdr2WLlTb_eG~XLbkYMrU>h=S63B2j@p;b_W+kXYP*+qcab{MbVkv!Nt*;2jP6jLL+jmx4lyN%1EGrNr|qBFaVE2A^JjjN(FyN#=(GrNszqBFaVYojx}jcL)D z-Ny9j%x+^wbY{0PGdi=|m=&GbZCn?f*=<}Oo!M>N5S`g=+!&qNZQK-{*=^h$o!M>N z5}ny?+!~$PZQK@}*=^h&o!M>N5uMp>+!>wOZQK=|`FPwNo!M>N6P@`)+#8+QZQK`~ z*=^h(o!M>%}-Nq}?ncc>#(V3^?wdl+<@OpG+ zxA8`F=2>_%I}2=csn|?+ju8Bv)gz#I9fo!MjyNw^CGta|M(V5-G z&(WEmz%S96pTw`xncc>3(V5-G@6nmv#vjp{-Nv8Mxe2%NSN_pp%50;gICp6$Y@?(k zIC%3azC+bF4%yR;Lw zQBpT|X(w!>WQE+Nov@9PdbvwGVH+hY<}U4oZIrB(yR;LwQL=LG(oWb$$tt-^J7F6o ztL854gl&|pmb*X%(gl&|ppS!dZwo$S{?$S=!M#+Y`OFLm3B^%`~ z?SyTVG{{}r3EL>yICp6$Y@=k8+@+ncjgn1smv+K7N;b<~+6mhz**te?Cv2l+i`=E1 zu#J)}bC-6)HcA@iF71SElx&r|v=g>bvUTp#PS{4tHn~eXVH+je<}U4oZIo=6yR;Lw zQL=sR(oWb$Nu%7Qov@9P9deg;!Zu2FjIJoQu@lcXj?6Ynn!wKNHg<-c*=_6sJG0x^ z6?SH~u^a5nZew@Yncc=7urs@jJz-~d8+*ad>^Am>o!M>d13R^7Rh z&g?dt!OrY9n#0cQHd?^W>^Am=o!M=)gq_)Kw1S=4ZM24+*=@9eo!M=)g`L@Lw1b`5 zZR`g-v)kAoc4oJ60PM_e<3QM%-Nr$%GrNuUurs@j4zM#fz`?LHZ;XzxGjD=UurqIp z&ag9YhAyx(Z;r08GjD-zurqIo?yxgAL=V`Rw?a?YnYTtS*qOIMZ`hf)MIYFiw?kjp znYTwj*qIxlKkUpqU;ymQJ7OU08e=C6!r;hFI1Y)<*5Ao!K262|KepI0|-VcW^Z9%9L+juTIv)gz+It5}o-;{2HCvZTuFU*=_tDo!M>t5uKZG8-L=j$ZVsu zC_1v+C@qeT>^4eEa+h|(HcI&mc4;SUqqHn{X(w!>v^;ldCv2m%B6n#gY@@U?cWEbV zqqHh_X(w!>v`+5QPS{3i-Q1;}u#M6ca+h|(HcIQ|F71SEl&+Y&v=g>bx>D}aPS{52 z%DGEBVH>5ZyR;LwQMy|0(oWb$Y5m-lkljY<>bXlh`GCAe?$S;clh@2$ z+Q|~~TDeO*SxR0zcWEbVqja6zrJb;i(sgr}cEUDF*UMen3EL=LKX+*-Y@>98+@+nc zjnWNsmv+K7N;k?~+6mhzZIHXP6Sh&haqiMi*hc9lxl21?8>O4(F71SElx~*0v=g>b zx_R!>PS{527P(71VH>4e<}U4oZIm|5UD^rTDBUV|X(w!>bnD!uov@A4ZE}}(!Zu2` z&0X3F+bG>GcWEbVqjdY+rJb;i(nh&UJ7F88JLE2Hf^C%Uh@E)fn%PEabY{1)M|5Vlv1fEL(V5-GLD8ArM*HZ@ZlgnVX18&0bmooG zF*@@m=oFoKQ*@5b>^8bYXLcK1qcgjWZqb?DM)&BfMB2tcj)Nk*jlr-pyNw~RGrNtUurs@jVX!m1 zjp48}yNwaBGrNtEurs@jQLr<+jnS|(yNxlhGrNtkurs@jaj-MHjq$KE*I)wd%uO*7 zcIIZ71Us|aI0Sa)7C01kX18$|?9453IPA==a0Kkkt#Ksm%x!QK?96R(H0;dna189s z`{7vFnfJ$J*qIN&aj-KVh~r^rb{i+a&g?c$gq_)KOo5&GV4MUyb4Q#EJ98(T0y}eO zoC-U07n}w=b61=WJ99Uj0XuVdoC!N~51a)%b5EQNJ996b13PnXoC`a1ADjm}b6=bf zJ99r=06TMkTnIb!09*vSK^TaOaY^LC950Q|JOoptGY`dO(V2(g^61RNaYc0I5x6os z^GIA3op}_lj?O$9*F^2reXLcLUMQ3&!&qrr=8!tp>b{j87XLcJeMQ3&!3!^i; zjhCY{yNy?(GrNsfqcgjW*P=7e!0XYOXX1_M%(L)jbY{2lR&-{!@pg1(xA9JNX1DQf zbY{2lUUX)+u_!vT+ju`Zv)lL}Iu#K|Sa+h|(Hp=Sfu9WOH%2v-^+6mhz zTO)UACyU8z<}U4I33;vDrJXD#ubsQJ6Sh&dPVUl9*hbm9xl21?8)fU|F71SEl&zn; zv=g>bwn6UFPS{4-hPg{SVH;%|Y?Iujov@9vO>>ua z!Zyk_%U#+D+bG*ScWEbVqil=ZrJb;ivMqC$cEUEw8s;wTgl&{W6Sh&deeTju*hX2S+@(#hjj|oEBkxOh!Zylw%3aEA zqpWdsX1CEKI^5qmGrNtZ(V5*wv*^rjqj_{@x6vXxv)kA=Id==*({8!060wZZFGvx>^3?_XLcK1qBC!fuF;vdK)2}3Zlim2X1CEJIAb{l=7GrNtx(V5*wzv#?vqknX4B5h*;$AOXE#vs_4-Ns#xU5K-Nta(ncc<+*qPnNNZ6U(#wggC-NtCxncc=1*qPnNSlF4}#yHrS-Nty> zncc<&*qPnNMA(_##w6I8o8u7Jncc>ruru$A!(eA_iNj%MZiORYXKsxnVP|fGqhM!l zi=$y@Zii!FXWkFT!p^)uCd1Br0FHy5*=-yTJG0w30d{7$aU$%@9WVuU=7Vt(?93f; zGVIKqa0=|qopCDc%w2FA?95$pI_%8da0cwm-Ek)D%sp@x?94rJHtfv3a1QLuy>Tw= z%zbbk?96>}KJ3i>Z~^Sh{c$1e24Vm%!o`sXal9lt^I%*Wop}hRMrR(1%c3(6!{yPL zhvSOq%p-7Rbmoz`DmwEhTpgWxG_HxxJO<%7;o!K2c z1Us`kco=qOckl@8%4Jgq_)KJO#UpaS@)zGm$Uh_-u6MOEEt>^HeN|&U_i3i_UyGo{!Fa1zw2Gd?j9t z&U_VKiq3pB7Di{j1}{fvz80@UXP$;vqccy(Ytfl!;PvRtGx0`r=2>_%I}2= zcsn|?+ju8Bv)gz#I9fo!M;D_kUbMa$zW_R#YbY^$(b982R@Jn=NckpX;X1DQM zbZ)|J{Ej~&yNy4iGrNtyqBFCN@}lU-Zlk<7IydrmLCv2m^t$!Zyn5 z>`Rch#JNbaTM()y17L(V^UE0YK@>;n|J6TFzJ9lX(Y@>Xg+@+ncjq-JKmv+K7 z%Gb+X+6mhzUq5$gCv2m9gWRQ^u#NH!bC-6)Hp(~3UD^rTC~uIvv=g>bzH#o-PS{5I zCb>&HVH@R}<}U4oZIo}8yR;LwQNDTZ(oWb$`4+iLJ7F8;Tjnn9gl&{J%w5_E+bG{E zcWEbVqkQY!rJb;i@@;aLcEUEwx6NJJ3EL>&E_Z1sY@>Yp+@(#hjq*m=f%ly|VH@Q; z<}U4oZIthnJ9pwX${R;#b{kEiGrNtQqcgjWU7|CygYsRYGrNu5qBFaV-J>(RjXk0> zyNx}gGrNtwqBFaVy`wX`jeVjsyN&AT%xgvcHkwCgb{j3CGrNs_ zqcgjWmeHBrMyu$|ZliT{X1CEMI^Am~&g?c0h|cUb4vfz1 zHV%r;>^9m*XLcJMqBFaVgQGLMjgHZo-A1SA%x0Z_ZsR!Encc?murs@j6JTdgi~N=?u=7mXYPX2U}x@%(_v@shBIJi?v68IXYPTsU}x@$ zvteiMg>ztM?u~O{XYPaZU}x@&^I>Q1hYMgg0R3?xE{Z&mCu_T zV+QQZ6EG8YW_K_Pc4l{Q9qi2R;Ck4Z-N6m8GrNNuVP|#+H^I*A4sM2>*&W;hJF`2u z6?SHKa2xE*?%;OVnccx1urp7_ovhr3{Bb_aLE&g>5Eft~q8+zUJN6x;_p^GUcL zcIK1u0PM`C;6d1#PsKy9GoOZsVP`%akHF4+1|Ef-`Aj?pJM&qX4LkGMm;*cWIhYGO z^SPJ@JM(#X9CqgO@dWJ57vM?QU4#qq6rPTJF~?`3Ghc#dqcdNM`O%rDVnKA~%kW%u z=F9PXbmlAYLUiUU@nUr5tMF2E=Bu$VI`cJnIXd&TcqKaXG`t#}c{*N;&O8IJ!_GVt zZ@|v%4&H>F*&VzEJF`1@8+K-Q@DA+E?%-Y6nccyAurs@ZMX)owgZE))b_XB8&g>2r z!_MpumcY*J4wk~s><&JJo!K3H1Us`k_!xF(ckl`9%=h9`*qQIcXRtHhkI!LeegI#< z&io+0gq`^zd<8r6!}uC@=11@i?97kiTiBT&!*{SV&&Kz#Gta>furtrakFYb(!%wg? zKaQVaXMO^|z|Q<6eudrBcnZJa_sGw1{3AN^v-mSQ^L+djoq0h~QAJU7WVcaK939zh zRFvc{?c@b=Y3|ZaUL=?0F74zca(V93P8O1R0d{F8FOw^Cmv-_Bxhi*QC$Ez0*P_a_((k9qJ#mc!$n_veOtK=?if*n+> zn!B_Kc2Kce?$RdMK}G%Cm6P2;#p=0BoA`jdM()xk7L(V^UD^aYs8}m^X%kDyYv(R) z;zRN}xl5a12Nmn)E^UGxRIHb~v* zP|+ZFX%p&nu!44`m&0X3AJE+(!cWD#spknjfrA@GdiY;=NHo*=m zw#;4H1Usl`n7gzIc2KcZ?$RdMLB-a&OPgQ^72D)4ZGs(CY@55Z3ARzOUGCB**ha)%w5_E+o;$nckaY(R5Xsx>^7Q2XLcJqM`v~$yF_Pp8@onl zb{o4zXLcLAM`v~$dqih;8+%4)b{l&|XLcKVM`v~$`$T7U8`aU7-9}AxX1CEaI^{gHd;hyb{qRfXLcJcqcgjWR?(T=M(gOzZlg_fX1CEcI=&Kc zZR{VN*=-yUo!M<17@gT|92A||ZM2Wh>^3??XLcJ0M`v~$9iub5jZV>--A3o=%$uQ0 zbY{2FH9GSa=oX!MOLULU>^6EtXLcJsqcgjWUeTG|M(^m%Zlh0hX1CEdIyRBE(T`*Q z$Zlf*?96UsAneR;V-W1jZeuX)%x+@{?96UsDD2E`V;JnrZeuv?%x+@@?96UsB<##? zV-)PnZeuj;%x+^0?96UsEbPo~V;t?95Fu5q4&`F$s3&<~Rg)<`y^< zcIJI?80^d~aX9SEt#AbF%&l=G?96R&6zt4xaWw4A?Qjh2%=_V3*qPnNWZ0SA#&NJS zyN%;vXFdogz|PzrC&JF$0aIXSJ{Tv#&fF0v!_M3Zr@+qK8K=U|+y$q>&fFEJ!_M3d zXTZ+f9cRMM+yiI9&fF7c!_M3b=fKX~8|T8#+z02u&fFL0!>&L2;R0M3c>u?YqB9S~ z#nG7u;gaaggK=qe<{_9Gop~rOi_Sa@mq%wFjw_-wkHD4DnMdNP=***Vb#&&@xF$OD z7+f2jc`T+yXC8;?(V53%Ms#L(Ff%%{JD3%n*&SRLo!K2+AD!79+z_4F9o!h5*&W;z zo!K4S9G%%6+!CGH9o!n7*&W;#o!K4S9-Y}8+!3AG9o!k6*&W;!o!K4S9i90E+!LMI z9o!q8c?#}}&g>5EkIw849*EBD4jzoo><%7^&g>2zj?U~39*NHE4jzrp><%7_&g>3m zM`v~ibD}f5gSpX}-NC%*%^9~{ zXLcJ4qBFaV=b|&ajpw5?yNws3GrNr!qcgjWm!dOYjfK&f-Nwt&nXknw(V3^=)#%J_ z^7E0XLcJOMrU>#A4O+&8y`n!b{n5WXLcK(MrU>#pG9YO z8=psKb{k(rXLcK3MrU>#Uqxqj8(&9fb{pSBXLcLkMrVEu-$iGhjqjr~yNw^BGtb43 z(V5-GPtlp(!Ozi|-NrAla}#dkSNs;)ZTudc*=_t0o!M>t8J*c}{1u&3e+@+ncjmnjCmv+K7 zDp$!}+6mjJTs3!TCv2l~wcMqhEF#y>T?N@~RIZ-8v=g>bxkm2NPS{4}nz>6mVH=fe z*X%(gl$x=pS!dZwo$o3?$S=!M&*XN zOFLm3l^f+Q?SySqHppGt3EQaLICp6$Y@>3M+@+ncjmk}Pmv+K7DmTkr+6mjJ+&p(_ zCv2l~i`=E1u#L(sbC-6)HYywDF71SERBn~Kv=g>bxpnT+PS`=^Hn~fiU>lX&<}Pi5 zZB%ZD?Rnq06Sh&=D0gWmY@>3A+@+ncjmjN!mv+K7DtF4AJ8>J8jiWQWjV95V-Nw$* zncc=N(V5-GuF;v@#%|G>*+J#*(V5-G9?_ZI#-7oc-Ns(gncc?T(V5-GKGB)oMs;*% zw^0+F*=;n9&g?dtMQ3&!&7(8BjTX_F-NwGrncYUq=*(`TRdi;z(K#2SsOg8||YryNwRfncc?0(V5*w$LP#% zqf>Ndx6wH|v)kwroq2O~jn2FUxvncYUu=*(`TS9E5#(K|Y`+vpP= zn@HQ}%duZ%x6vPVX16f_c4oIR5O!v_F$i{Mw=o!YX16f}c4oIR6n18}F${KQw=o=c zX16f{c4oIR5_V>{F$#8Ow=o)aX16g0c4oIR7ItR0F%EWSw=o`e<{C_Zo!MZ7hn;>^9zy&g?cmh|cUb7Ds1x8%v@yyN#vMncc>R(V5-GN3b)yjgMhxb{n6- z&g?cmg`L@Lddn*gc6S@C$y8{1nIEqBB2@-=i}>gFm7(KZ`%3 zGtb9g(U}(%6;%~QM|K-k#nF-7Mpa4f(oS9=m*y_*?$Sddt z!49f6$X(h5JE+<)cWD#splYMsrA@Gdss_1Bn_vf38|N-P;Bg0!0ztu?#9;t`?}BPUYu>fGxz=fzVYnk+jH&Oc5}}CIiI^uJbMQ@ z=Lvk!sCAHYp1=o<_6c&%6ZoLfzCq4;0v|NmFUUDh;EhK62RY{nywT_YXoD8!6TH!= zZIE+5!5fX*1v%#vywRwAkaIr48;v>yIp-6+(WqmPGoP3@8g)|6V)I64P+s&SLY% zVai!--bhi-V)I51UVB4o^65q_RCqv3Vm6EN8KKBONSfv3Vl{EN8KKBNHrVv3Vm3 zmb2Krkp-5s*u0Sqmb2KrF&Qjpv3Vm0EN8KKV+vT#V)MpSu$;x_ja;yt#paDXu$;x_ zjeM}2#paCyu$;x_jY6=T#paD^U^$BqgCek;#paD-u$;v`pad*uaZe}(%URqDri0}y z?hP}*au)Z2nP54K`@$@+oW=cMHdxN${%|x{&f)=Z3|P+Mfp9EX&f-BZ2P|i?dE+>+ zoW=S)2}MD`#;AoTHq@nXpVbi=%L^au#R7a^)<}hVzuOcru)? zoW(h?LOF}4zy-=#JQXff&f;9SNI8r1;9})0&WDxCSzG{@C}(jYT&kSK=7UwrS!_PI zOgW3q2bU{nvH9Q%b?t4{lJ-V)Ma`%2{kaxJfyS%?CFtXR-O<7Ue89AKa>(#pZ+Cl(X1;aJzCA zn-A_#&SLYyoyu8kKDbLci_HgjD`&C!V2yGXn-A71XYnbpPC1Jg!#&Dbyaeu5&f-(y zKIJSnAKb5;#pZ(tl(X2p@t|_%3G>E7@UUX@#v{sEY~FZOIg8C3k11!ddE;^AEH-aE zp`69$jVG0}*u3$Sau%C6o>tCc^TspES!~{TRym8!8_y|cv3cWp%EJZ@i(L#paDSm9yBq@s@HH zn>XH8&f?YZj&c@X3-2mt@pbT?au#0??<;4qdE*1+EH-a^sGP;-jgOSG*u3$vau%C6 zK2gqM^TwyjS!~|;OgW3q8=os@v3cVQk+@@!b3VZvi8}>3=M%h6TFewEXX;Z z;ElxQLC*OEZzM*7ob$=o#4UoH^9kNaY#HR7Pw+-!t03omf;STP406sVcq4JIAm@C7 zHxl;_a?U4sBe8Xmb3VZviTeaO=M%h51SkB@>Fb6DW zv3cV-u$;w1;CQf{#Y15(SkB^MFb^ze@o+c+ENAiIFdr;u@d!8(ENAgZSOAu@coZxI z%UL`cP6EqWJO&nl7DwS+n{8feVzgcq&||oW;3tk#ZL2!Ntm1oDVCNv$y~*QO;uX!KKPsY(7|} zoWo-l9x48JHgZ~Us9 z#paFQl(X2p@w;*sn>YSY&SLY%pUPQm-uO#7i=TwQm9rRcBqbmcWRf;W=333ARScq3`sAm@C7HE{JVHgaD!(jxB zgi$aW#=uw@2jgJ^OoSs~5*!Jsa1^9LI%Gg5L?H{ZVKU^v6qpLRkO%ot0EI9Oil7)u zpcJOV444VCU^W~L$H1{L2abc|VJ^&r6JS1^2n%2#oCJ&DWH<#D!xA_ZPJ`27DVzak z!dY-OoCC|?Tv!h0!TGQPE`ST+BDfe!P# za1-1Nx4^A%8{7_ez@2ax+zo4BEv$ok;9j^7?uQ59L3jurhDYF0cnltgC*VnV3Z8~% z;8}PMo`)CUMR*BbhF9QKcnw~MH{eZp3*Lrz;9YnR-iHt1L-+_jhEL#A_zXUWFW^h~ z3ciMK;9K|(zK0*+NB9ZW!_V*w{0hIp@9+ow34g)g3C{ZU>and44WJ=x2Ajhc&DVzz>csJ>;cW7IYgiZw1if$C+r1#Lu=Rv z_J#dme>eczKwD@B?V$s7gig>Iy1;>O5F89$p&J|m-60tcg~K2PdO%O;1-+pU^o4%V z9|pic7zBf12n>Z`FdPnt5ik-)!Dtu*V__VOhY2tdj(|yUB&5PokOt|H0hthmEXan* zkONa-D&#^QC7UVsq7%*2%Ev?umv=NL`Z@yVJp}gwt;P7J7^5s!w#?`>;yYQ6W9fIg{H6@><)WC zGiVMGXaOyu73>Ln!QRjs_JMt2KiD4*fHu$;+Ch8h03D$dbcQZ)ARGh-Ls#eqhd_5o zhC|^nNP!;E6M8{!=mUMBAM}R-Fc1d8U>E{JVHgaD!(jxBgi$aW#=uw@2jgJ^OoSs~ z5*!Jsa1^9LI%Gg5L?H{ZVKU^v6qpLRkO%ot0EI9Oil7)upcJOV444VCU^W~L$H1{L z2abc|VJ^&r6JS1^2n%2#oCJ&DWH<#D!xA_ZPJ`27DVzak!dY-OoCC|?Tv!h0!TGQP zE`ST+BDfe!P#a1-1Nx4^A%8{7_ez@2ax z+zo4BEv$ok;9j^7?uQ59L3jurhDYF0cnltgC*VnV3Z8~%;8}PMo`)CUMR*BbhF9QK zcnw~MH{eZp3*Lrz;9YnR-iHt1L-+_jhEL#A_zXUWFW^h~3ciMK;9K|(zK0*+NB9ZW z!_V*w{0hIp@9+ow34g)g!H~Z}J+}3s0W^foU~||48bKl?!IrQUYz^DMwy+&EhV5Yo z*b#PuouLWr0=q&}*bR1vJ)jvhhX}NQme30JguP&IXbtJM z&GG>97>e zfHUDNI2+D^WpFMmhx6ckSOFKng>VsE3@hOhxD-~wWpFuM0awCRa5Y>5tKnL>4z7nA z;6}I!ZiZXnR=5prhdba-xC`!vHLw=e!98#<+z0o=1Mna`1P{X_@F+Y6kHZu2Bs>LA z!!z(KJO|Ii3-BVm1TVuY@G86pufrSgCcFi2!#nUUya(^Y2k;?$1RujE@F{!-pTigM zC42>6!#D6Pd@Fr`p^Ix!e+2JYypiR z5t3j_*b26WZD3p24jRMuumkJ}JHgJ-1a^U4p(*SJyTcyP44Oj(T0l!^1$)9?us5`Z zePCbM5B7%xpbfNzcF-O=Ku72VouLaH2nWHz&=tDDA2-57zV@Pa2Nq2VHAvpF)$X!!FZSe6X6J$1V=(D90h5R4jGUMQOJU9 zm<%~E1*SqSn6F>ox*f#cwKm<#jZ1egye!U9+b zC&3~(8BT%4umnzp)8KSi3TMEXa2A{m=fE;J7nZ|$a6YVn3*bVy2rh<|a0y%rtKc%Y z9Ik*X;VQTqu7TBXEnElJ!wqmF+ypnnEpRK`2DifxEJn&`{4n2 z5FUbu;SqQg9)ri>33w8of~Vmbcov?6=ivo-5nh6q;T3olUW3=+4R{mYg16xvco*J- z_u&Kh5I%yB;S=~2K7-HU3-}Vgg0JBl_!ho{@8Jjd5q^U8@H6}Zzrt_uJNyBE!e8)r zLc(SVP!H-u184{&H`^SxfJTrANw6hs1zW>5uq|u{jbVG(0d|C)U}tCoyTGo{6n2B% zVGn2q%^?CUpe3||Jz+1{8(PCYurKTf`@;dy2HHY9Xb&BrBXok!&;<^JgWzE33f#bU&;xoxFX#<@pfB`;{xARr!XOw7LtrQjgW+&EjDV3a3P!^i7z^WIJWPOz za0E<(BOw)zf;32n49J8iWI;Aeh8&mzQy~}fARh{#5T-#96hjG=!gQDcGhr6YhNIya zI2Pu>ad14$g?VrS%!dTeX%i%mY zA6CExa3Nd-7sE=p1TKYDa2Z?b#M>d3-`hO@Blmr55dFm2s{dp!Q=1*JPA+1)9?&D3(vvx@B+LDFTu<3 z3cL!h!Rzn_ya{i?+wcy&3-7`E@Bw@XAHm1)3498l!RPP=dFZeqlVeHrEePKV?9}a*v&=%T3 zd*}cip%Zk5E^r_m1P4P`=mv*CcSwdq;V?*n9?%ndL2u{-eW4%phXF7U2Ekw$0z+XK z42Q#E1dN1HFdD|dSQrQ6VFFBqBVZC738`=tq(M4lKqf>X3$kG{(U?H3Yi{NB91s1~+I2BHV z(_tx`0cXNla5kI+%ivsC4(Gx7umUcC3*jQT7*@h1a4D>U%iwaj0QS$ z9b6ALz>RPd+zhwCt#BLM4tKzva2MPSYhW#`gL~j!xDW1!2jD??2p)z<;8A!C9)~C3 zNq7pLhG*becn+S27vM#B30{U*;8l1HUWYf}O?V65hIimycn{u(58y-i2tI~S;8XYv zK8G*hOZW=DhHv0o_zu2@AK*v$3D(2U@C*D3zrpYD2mA?t!QTl9TO>d|s1FUGA#4Vl z13%oN5hOwqYzbSz*02q13)?|s*dBI(9bqTf8JfT@uq!l$-C%dv1DZi|h(HTy39VpH z*bDZC*02xk3;V(TZ~(M{w$KjRLkH*xouD&xfdksX;1{kPy(ee9cI8xm<6-pXgCIrg*k8>91nA09-IL4 z;Y3&f3*jVK1Si8Ouo#xWsc;&c4ol$-I1|o-v*8?A2IstGu2p7S{uo5nT zOJNmU2A9JXa3x#?SHm^18m@)w;Ci?LZiJiQX1E1zh1=kExC8ElyWnnE18ZR&+ynQ* zeQ-ZK01v`L@Gv|AkHTZ{I6MJQ!c*`xJOj_dbMQR8058Hz@G`stufl8aI=lgI!dvh* zyaVsTd+pTcMGIeY2~c7xqv4`>F>Ap$L+CA5M)VK3MlTEjlDFYE{V!vW9++Cn>M4;`Q*bb`*%1rCIR z;9%$q-QW=D4#{vR90n=S1A0O)=nZ|KFZ6@{FaQR^AQ%ioU?>cO;cz&NfRQi?M#C5w z3*%rsOn`}S1WbY>Ar+2-G)RXG$b=|lK{iZ=9GC)AAs6x>9}1ulra=)DLkX0^beI7% zVHV7Wqv04h7UsZla6HU~d2j;EhZA7|EQFI_5u6OCz+zYer^0D)IxK}V;7m9R&W3Yf z8Jr8t;XF7WR=@>tAzTC(!%DaWE`?Qa8C(umz?Efvu;SRVH?t;5v4XlNAa1Y!I_rd+}06Yi}!Nc$fJPMD&!OQRpyb7+lA=32(vM@D98S@4@@<0elD_!N>3kdJM&GG>97>efHUDNI2+D^WpFMmhx6ckSOFKng>VsE3@hOhxD-~wWpFuM0awCR za5Y>5tKnL>4z7nA;6}I!ZiZXnR=5prhdba-xC`!vHLw=e!98#<+z0o=1Mpx%LP9;L z4-KFpYzCXd7SIS1AqlpGtzc`|2DXLmpfPL@JHU>x6YLC4U>DdGn!;|dJM00?pgBaK z1+;`#uqW&VdqZp32lj>iV1GCO+CW=q2koH)bObw&{|4+>GJ$g{;ZgXF!{WH>9SNL& z_0X=~k4Q!x)gK3(tM$>j{?|wtfe|ntI42r>ilo%NTgRWnmh;EfUfsH7U$!%WbIsOi z+!{v;+ag#6FLMf@zne?svb_+V<&>xoogoX(g~#C!P7^dqJPOW$2ib3a@4tcn8lRiL z{~IumcY-L;<|aG_zvBtog?gz#`&{n<_=(exHl+S|pe?UYJ6Qi4s(|*U!ALj}Xfqpp zh9_v7+;y(w*wVh*wF$Sb-;eDqpe?p_Ot%KugKaTf2Cv`=+DA9n&SQHKJcl2sXE!d- zX1g4ofIsmCZIl3U$u)e_h0qXu14P#+E zOoU003Tco5Q9-+wfMy9#!rz=nw8!=MSv}gidgxg1FAkS|s-FqVfcCFG+BTpsYCt&+ zXtx_Y2DFI{+5d*=a5m7E+chb7Er);8T|0P$E$zEo<2{Nk=Y?B~`<1Gs&39`iw5M*4 z{tM?2_2TAIOWD%K+L*WO*;@Jr{~Ut6Q!2W7cuJ@~iXGac-!dv=3kOX3jSa~bI8 zp1MF{ZPGZI1ipe z2YytabEQ5S)_M?Uz#a@PGNMRW)D# zs6JXndKE>}N}`#O^jVQ1!-kFQ6X8Vq2Vd3f;6E-f(yJgZFD*YaGB_tc8tGA#U78oo zFDZ`jZ`EQ?BW>FimlR}7PAe|S$w)0JD99~tQ&3nEX_FIamzh?Q)-FFfJ2Nl8U13qd zlxRjtal5>Ng5oIu&n(T27Ki?}ZQHg9xy9{DqQxcdf5q)ei=#ya>A_zXrKjfPPmdNA zM^p3C3JY`cv)kpR<>aRprMJyt=ObhIQBFaAWI%paL4~8~Ta=a;ol#IU)%~C-FRg@M z>Q;+G^?wowtADLqr2H#8rzDa=<+3>_J(`t~bzoXn=Lmb$v3#Giqc;)#4(NPHj83Z%dmL>(|J^eTLUuAdOe7U%g6;ia6^c zqjK`1<;P1^jw+2thNKlmI&_S5=+v#_!QDC>=&EmH|3IIa(TviPwDjDlE3>#HEw9iO zT(Q88-8yyd)}{R>6}Yi~@c;hjcFqsTFDxyI^v&V0Yv0wNqRL8|fBB`bMInzp` zsc9KzOI#T83Z zoE#Z2bWD%I1A5m~LMR`JDNy?8$U``bv(()`oG~l+x|xrKaZ= zWK1n9#usPjVXP{vz=EPoy0P+td=z^fMsuUqv8TEceKF-*oHBZpO-VNKD-x-pAYWwc z1s3OIMr+zVN{zk1l39f{+0OPdInuAk=#e8g@qn3_)phf=$^Y`E{^G=*y^8<6cQO*0 z*tT76G^-@Fv~G6Pgxp{!;JSdGy?VF~A}upBH9dn`(mqsEu~y&)JAv4(w~jmYH=Fve@=dRNs4MCN@h-R zVQyN6YiFyp*|j<3Xucc4R@ph1sEx|enI%POOqWWEa%Orh%1W=*9=ZZ6i~ElX%*iXv z&B z*5*7cj*!oEp879Eojd z%l4**Q##nuqG_44Dpvf4VnoTvh!(PZn9uZa zChN{KqnW9N1v#uJ=9EN>+yXnLZgK+=F0e2sBdsX2c4fI8*J+cMIXx{uqq=3k*j2Xk z(5??=q!wqSrr0vWE-IEh)$? z^Kz`#)518*OI%*|#rUSdIVIYr%O;BWAEi6ebEDYAt!n+uQ>*N&X@-|tn5`)Ipz~%=hVrI_tYUMVS_m8R9Q-OBdKy@X?-rA_ZnodS{9>%(dH63oCGKU<) zgw)lcYgcWQS=OvtQNdTO4%v-`OY<49RsC`*$Dp&@!)hApvXYY{F;pNYs4LLg)jCzV znufZ<9?~Ah)MDqC=A~vAe#UJHSy^1GNXstBPs=T4r0tsj@&SD8 zDi@TNFlet~|7;Iq9l^ThE|t3?nb}rT>_04!L#}cY)IZw8s)Np2>|rpwtEn%Iy>ZY$ zvzl>|I~QYDSq0W)zu3>ivY~MG9ke^-n4E&|(q_7Oc4`(gs)7o)XxTBxzK5gjmf7NCWs}|^nj2cN;*Tv*}lLj{K zkp0p?c2Pm8O&~Qg_WKySmy4(76sAUJ+G>#9av97NeQ|ak#wtp2ezZ6*E43(UM_!#M zv33F6c*0E`HICyXTrYWbFmLaD9}By!TrBgZ^g65 zVeo*VeR_=8*q-^HyFHv~tLg3uAp5rR7x#Vjj_emPW=8gpblBLPU-s9P>zV(gR0|H+%=Z8?PdYobMk z1w6pY#GbXQva4OOPH(OHr&rLRN?4oL#xD$1Da;c8g1hdypk^y4l~1=Czxxll>%O?4 z)Gqed$*kGx{htn&{QkIlwF2XVC98$>}TtM9qMeL#A zu)6+x?sKTver`0jzxy8>iu?yGRk5NF{^KE-orfDAeffF`7OVVs zT&)O}UxGae&Uq^^KKim+DF5QE9Z{ZL$;t4SD|`*s&`j80*7;y*Bkv-JePLWT@|PcG z#i%yD9;Xe{MSI^%H5XL40zPbLpBgHxnzZdqW;7!=jeBdW%;xQfwVBb^`MOP=k8OfG zzqq@jpfFlgl2iGzPi+o)gSV0;N0{Dj+Kj?UWc6O>V;fgs^_RFz<=B_|T!GbJ;<5s( zxkOvd`R7x1$p8K?H>l&m(*NE?E_?f7UEbobkzU!(!{8Cc-g<#c=5XFJ$@{X@MYB8p z3=Fm#>6L5K9&T`1+zO1fPpR|w&e>`0_p-`MKx_|Vy-h18KO?u4SAOx%G+t-Nox0hR znfsO9Wfl7hr`qev#aNpWe8pa%8#kop)buiy*hL8n4Bw#Zcd!;J`+-OeI~jW?RaRj2 zI~jX{ZXi-r%JXaP{cgEA#oT^ruP)3lFzxx4!`}2;9U>7YFRS6&(@?ieVX@wsR<%G* z#;o9-kw_D}${Tz+S8be9>>h=#@!fc_)SoqAbfR_-XC|)%bhmI8MY&10#64-9nwDjk zY>J{;_BP;6Em#GXl~U$me`arLLB5%#K)dDA-KXozsI-k;zMQ3h=C@xLU|H=wZsFv#)S?0|Y;m_pmE(?m57VP18z`?9UgnPT;?#;4 zFDvh4?1$_MtmayJt?SSgIJ@$0RrszJduT5)tx!TJtc570Ik}nk7`eItdz)*85k%}x z$>2q%MFly$fWWs~v{}JZ(smlfzK7+VR9UYSbW&w~7KLv5_{PSEydE`p9c-xj+_$ya zLo2YxMl3lpZv2EA8>^~6u=_A7Jgr=HA2;y<|FV*~gOdSs%{y2NQ8u_hP@8TvJB>Hp z*JW#yomN_0%=k4ut%Nt0*kf#!Dl+yH+Yeniu}gUKLUGRQN{i&NpEkMC{Opp+wkTQc z8w+DE&>oP=jizN)*|*@xaY;AU3#_@d z4hr1#md5uP^CaQhoG3z#c`t}$-cvIgA?g|XXgO#3UvYm5Zu|ES|{rT+6 zLzK=jo;E>&dZ?tB8@5WK*%fqRF)hjiB4y3f#sVGWs-zZvU+`q>ZHuwaR3S z%`C;e{G__kVQuzsChb?*SywfP{mSFaylP$^TwZYO&DX4r?P08q!>qdA!`XGchqLQ; z5C0Rcjn?rBXZ(;WHNO8I-%K19yg8uGma;cx9h=V^D##E1;b8SzJ5gRzK_jacl>SA_!YaD%6#E`#6y%1Vov%(U zK5?gpuUZwj;m@Fnb%m2Xv;xB&yx-wk)M4>P?j??CG?&0^v>=NQ@$x%VWgGh*Zt!_R zJF#027&>6ofF6S*qXrD&<8t~8Pwg|bca6^h+5?JwMp;cYXa&Z4+PIJIsrtCPKla$2 zlATw_qGacn7UmXIHn874wy}#cnNN7BP(V10eNng0FLJs2?YUeX-Zgz6yUJ6Ft6A49 zFFp1;WUQ2v$4hhUM)=_IfJ(a?dx3R60oP@?oGkY;2WD0IC9FS}Sts@ro4YTvSY1gg ztoX45|D+Df3*-V+)sNIti zeehaOUPp!`Opdd;h{oqrx)Zb(hg`8h|Du*Mk!#T&27Ab2&}7~zl9q3K$7Q$zuXF4t zc2HpXOW3L$b!-Q1ciUJUmKPX0>ez~ND!T&R$&*`9@fFAZw28GFyrOW6s|EFk#wUe#NT!+%xiMm{FZyxh0X1~jOmi@h8* zy@(GrtLa(#*qbk3APeNx4VPjsF!;RX8s2+ihaBr^Q^Xsec%HZ@RiEVSuC$hEQi}qE z+YJgPYozeoxwYq8REEizHNFd#;S<>cu+oFnzz1l3x&OvpyC2_><#P&Z$P&? zjCDp|cAJ(R9XG?<_<{o?H%;j=c_pt25PK;iAl1TmcV0bT@R{VcCfo`yST$1?y;05$~^NuvGT$Lk6F- z`3Q6yw{wnW_>W!KGz;1yhP~Q=&o&M|1CUq5l^w1Hv7d+K1$viMD)FBb7<>X-rQ&L{ zhb%(!bhU2;DjjZc$gz&Gv!>ezG`OjEMWgwWW1j`ED8ggyS+mMs&d6MscQv|~*XHLU z%BDEDyldfoWWig|?L={f__Ig582gExlgE1nbGXEsJBz0wY!XzZZH%SRNF;bY5A#r_ zp}Y<}{1VC751Dml_pn(FMgD^VGo!Xb9?j42Evl=<*b8*;0nabBo1l0Jy4_068&Jw$ zPZmr4nbz*z&OBYafqpqZHW zS1>ntwC;y)4_$%o5Q{SDnB9QWeqdj*wl?-R!6=(XyM*bUB`SXXg94f7R4wnH>|s$f zre z?YQpC1%5lRdD@$6w*}mp)~ejofL=Q&GFHlXvpX}fId#edVbD6^6} zwkP)CGdZl-%&hL3N^OS3_7PR=d(=K~l=UU|iK`V3IrcGj+0unw0S|o|Uabc0zR3v9 zLxY`+{gD56JE?!a-btso>!j-B6`^GX#?nazy_4^SuxKUMxo_-btMfi3D=jZ4H++wX z_a2`L(WBpp9z*+W;+u@@?>2bDU~SsNDpwu--q|*G>#>I$`Xun$jInucw1zo_+rz#C zMvUCl2g8DSivKLtibq-&*XEFGyryTR@@jhPJ^Z=XlMn5xdh%7s*xz_j+J@HZO_t1ply*M4AMwK-%K#M5dzuPHFr^}~O`E6%K=<9$^5)qe6} zZO+4JPIk2FBmO}Z*6NV!eDsyIq5G$6BTBH$u0srOctkPA3p9MV<%T|svNq?Tn=hBv zaKCfxhx|V`a@b(Y5&vPl2I4gkuYq_C#A_g41MwP&*Fd}m;x!Pjfp`tXYam_&@fwKN zK)eRxH4v|Xcn!pBAYKFU8i?0GyawVm5U+uF4a938UIXzOh}S^82I4gkuYq_C#A_g4 z1MwP&*Fd}m;x!Pjfp`tXYam_&@fwKNK)eRxH4v|Xcn!pBAYKFU8i?0GyawVm5U+uF z4a938UIXzOh}S^82I4gkuYq_C#A_g41MwP&*Fd}m;x!Pjfp`tXYam_&@fxUu8aV9R z%7vsPHhAsJnI~sz8(p*P?)_G@`u3NU#_ftON?qMb+tw*NZZlfuMmq?BzR_$XhUUE4iM#*W!zVT%;MZMIX>5$k(jof2$Ur!>3r=7zH; z z`)Eqbxli3Wq^MNel9ZM!?pyI#t2x>pm(udK7gx<3c5(Uk*&Uu)*M#!j_O+CjpPsVW zvWM@`_Rf@634J>Cc>Z^7|4M1qWN_AFSGUr3uas78uKKv`l-b%IozkjTRDlwfpYWlUi%bzyD<42CrN(U)zN#t={@-t9~i3X!~kPt8eyQamDtd%C>tt zpY%=%w*BOfZ0{cVV!yx-s}e#!xLoJ`_d4I(={zsc`8`_abwi!cEp;AWsq;5o=k0u* zuNUe(E!FwiTIc1?Iv?kR&cpRO?{)~CZ~N)|`djDK4LYBul%G%TJo4wyP@OLi={%XC z^W$ip7b|o=yjgZWl${4o|Gh(cua>^Er0c!XbF*?iosMbJjlcK0dA*!Y!=y`Z>5{kR z;;-hvF!Hp!q)SWb(oni2sa?7>q)S`rQYu}>N|#2`WgqEsg>=b~E+>AQ#uG$5!%l6V`p!~n5boo=d+$>$9(&bY5`c?V4f4QHXkI}>V_n43_Z%LQiC#NP&@Q|7 z_$KL+DP2~Ebm=W!9+zLdF0`wzJ>~aZ`}wYPxmmibmoClamlLGRHPYpH>9S6K*(%iT z(f+yiZXfw&u5@`#?b#sd(qFn9z2?Tv&L}){x%W#e>5>>~$8M1>ze<;brOUC>Wt?=` zLVnprx?C(>#z~jcL%PhCU-p(Rjirn4Ps6(GB){~LE={Ei{jKjOUB8$mT`rO?FUl_i zrOSiTWv2Xcf^@k+y1XM@ZjvtFm-l0?|DykJ{npjeh4bI}WsP)64Ebe(bfF(|?btEW zUFu1ft)$Cx z=`vQjERim^Nta{fmp!D*w$f#RbeW}o|ErKL$4M*!F7!ukJoct^84zm6=5@bp^80O) z+kYuteE;Hg*-n1BN4jh+T?T~wk`d^#_E+iBNV=RA^2@2x<#y>pKjQlRh;-RTx-69V7A85+_hQ@X5|E}Ki2 z#nNTC{Ni=FCFGagLhaZo(q)GH^11x7rTYCY(j`&4d@5a5N|#X~U6`l3`R`!0WA{pz z+49Rt(&cLD;^(vBcFgxLr%D$;Po-UQehKR`R=Nz4F3sc@#??-jccn}JkY5&sblF$B z?5cLGvHaq7*-g3(mM$k&?qAma9O$y(Pw8^1bQvaHPLW^ks6iJ$PCiw-lt`B^qzgK@ z^^A7XrGa#LN4lIRU4}`QOEey%o}DfOqzmf@ZvJw7s2#gTy1b$B*aM;Q*zWSn(vU8t zp?<%C{6hQX{L)vt;5XNf@%v7fx1~#XJhrg=IWwkAy?m?nYRC4MF1t&Y-KEQ;@=HXz z43{qRrAxGOJGNlGbU8)343;h@NSB+Wi=S79{c>^*eo2rnZ3A7F|0G>rmM&*X7uqYg ze&pAUXqTKWz8$+wy3Ca>z8$+Zq|0L=U3QQz%cRTo>R)C^m%G%Cd0p^{)1{|$*;Trf zg#0o>e(~!t3%h^*{#nJJE}1D^+DVt)q|0v7ZiGE|gywUpZa;I?OH7#jn#X zRy%f;bn)Y{$EC}!>R)!0F6T-Y-@i$-Q!FGcG2ytf zwPX9pFJ6~Nqzmi1ZeB@y<^1Ax@#`=nLhYFEUzUXW{dd%kJtJLymoAN^%Q@1;_b<~b z`{jOq->r8JmM+~xe#x)FFTJG8;x$j5_UZ8p9!!y6I!l*b(XUE>5?H`PLeKK zhgtiTbop4i%#bdJNta_o?Vr^^cY zh5o|%<#6fZ+p%!}a-;m>$73&sblF+DoGo4aeBZZYcL#pCe|<=o14Dk9B3<5B|1vVv z?=KGd#p~k7$$QH$FH4sLLb`qFaG+IU*|tbe(~3v-j^=u;B;YK%&sdfe_OiTDqZrV z3+t0kmvckovEkCi_b=gg%#V{lkS;F=y1YWUEY!|O{B}~(uH*~H{YKuT@I5j%R;&^esh3X*=Ehu=(*TvR;AzfaPU&c$9?jc>gU#6*l@#`?&FX4IRo6?2xu=5M+FU~J~@A{X^ zrAwxCNs%sRhsI;Rf4Q@UcFd2*-jiRvF3jKEbumAGIY)lEUb@gORb02c|4ZpIAf!uU z>5?X0UX?C>Jm%NOPG56bv-cl;V%4M4rKfahu6C@GboofSw39B=2hYM9rJ#GHn(!a9lFPw+ybbYUIE_4~_1{R{m?#dXX3 zeZS9q(#W&Xw6AC%f(e(}2O6Y|Sf(xr|3QXJC7 zw_`uaFE2=!Nz%oSll{E%2K6ub(#4O*c2WQ0`+dI-gC2H&$*WHUy8OZQNv8|rW7m!? zl`f+~ehKU1`xoDi`FZ8(@(bfMx6c2ebYZ;YblFn6>?vI?sG%LZO@8t1nAc_hkS-HK zehK&c=X8Jn^y}8%z1acMrJwxLM!IyDUw)J>?WIfDFB9dL??U~{xR75Ks^7m}x_H0% z`M$rte}cwi-Y-u|mxrax_tIspbXg!>hD(?B@(Vt2_g(mP{_wixF&ZcPc8qbF)8!NC z@{V-bO1iX^E-R(WRcgnWzdOI&DqX%-JJv^jX)9e=2XKBF7xIhmU;KLVdEFmMd!^p1 zAEil`KGLN_s2$UNU>!P2m!m^|@w)hR%QxkhiPEKy+A+Uw>H8P23+i)9_jM6bh%BsoG4wmFU!sMJ4hF0G%wI?LekoKt#`=rX#rH1@rHfzZe^h>1BVAbc zckS4<(uMUJH%>k*qzl)*oGwpG7uIL|{cUbN*|%d}m)qr+&!h|2M_fC`ILh_=z8$+v zy7+lzSQqb?-tvpuv4!iU%TVd^f^-Sb_m{7EZ>Lx1|8?{6(xt!r(p|dr2y{7jGwE`W zbeSVv=7!c`bp2()>(ZqtQFnzJm2Z!*I~G?(cRDJ`>&d{4xY54f?=4-tUs9w?qIBshU1my` zMbc%w`WLUuJJMyEbQz#_Y?=BOe| z_nGHAzjT)_dq@}j;B!b_wdgmABF;18Bq)U!;=^|a4N|zbZ0r^>*xEo%P-ug;ntIV|H6HicK_JIkEIKFr^^%4 z#rOMuoP6$@m5uwKy7fb2)Q$}ZwPS;%%U05*N2q@}OS)vqFTZI#_Of*G>&Zjpm-D2H zU+14HU6Mk&_C#@h>?U0br3?LrTTkvTUDisM6Ez<5 z^OvwLew^(6(lw;Zjnc)>Uv82vetpdQrKNPif36*)eRH}zDqTj0`u%gdzdPyY);(W& zNxFEy945aEmoD2$mwuuC!ur1(k2R4lxzgoU>5?4Mh5Py3eM*bej`{U5-|v5m2T7NU z-o!G_UmS82Y&W#m_70H{5*R>#{<+Fg|m-oFHA& zrAt%k(pb9WNEiAI=NHxiTsyW#epx787DyMr&i}RAF|SKw^)H7?m&>G!_e)rp&!meV zkL?liON#vRfOPTe{C<7xtTijs^JXsCX`poR^GffRBc#iY(q*W0nXCR~iTuKS{cb(^ z2eo5fmlNfeku~_m_bGFVdVSeNK{qfSJiF8SnF6q*RaaJW= zXxCi-;;-+|mo8VSfARCm59ODgq>FFIE|o5R-Ey-0!nnE8`dB~th4oI?zc5dAe(}0I zEM5FMj319J3-$YcUg`ZZHl)jlQ2*lFu@dRx=ldVZFJ6}s(q)zW@=&PX-z;dy*8VNO z`0JHk7sg?(e__4F_4~{(onN^A?{pa{UA7PDk|AA~FFIWgl3!SFak`*~Yscn=+OaR> z7y4;8e`ylZ3Lw=bjzxZ~nP`dc*m0wGjW91ir-?)EX>k;WPS^a(k^)Csb{v}uK*p1TVbm?-J zbfKMdx(t*qTL-$_{ddSO=-_nev=LprU;O+ftPA%~xcQ53#~6RP{v|`Y`1w9QbNvgx zaO;-79b;X;k}m7y7e5~Jx-9R$GIGnef8M)V?N~33$9k!MnIv6ykuD=c?bvy0$7VRual?~;W-=)h$>2iZ~=`3CNzEV5JdX3Y?kH@|VwPRkFUee`K zwPRiv-@kajd=+ZPc9kyA%P)th9lIyg?|Z*24e8>qS9-r>NtY(l2fEDqTe|!z zUB*inu6Mb1tYgS8^edJ6m+RE;d%tjf%k?i_7y9u^>oC3@%a&jK{Ke~nj;VCuMznA=c&6{e+dW8HkIiyRfbXh1}d^=VwU4EA?z8za6zhtW&yE)YFXGs@) z;&k!-zQ3QJ>-DZ3yGpunz2Eu8>vDv2VLoYonZfiar z&tF)Nsnm}7{)O?b(zu(uHw=8;|*Q%l;u`}&y8I?x{*qtjYke&2 zm*b_&YWc<;xWetpchW8rm}rfSEYmM)C9-F)A-V~1bz`@l?Ax*Zr3>}%bn)x_ zT<5RU@B4Mjch!!uj^f5+14Fv34Xwk3+cAFM^)FtRr=-h)(uMn)+*B8?XG@n$ zrOPbo!uQTEUKiG1T>tWwboo%a43RDu%P(G+%cRS!kY9Yi|7mDE=KbPz@$;88(k0x# z`0Fpef7wOtnD2jBJ85&x*3_oA#?;H2J6iSyCAzeJxgykBU)Dvif{|Ki&*uCqA5aJ}F8<-(2V^1gH#Bwa2F^)D->%XI1D`+dKj z{F!v|?N}@MWsv-GMQEIy7n)c4e&4rayGfVF<(Cd2UBc~H_<1p}OIP`2s&w(~*r#d%s*2>i6@d3*%B4pX zO1k*xl6*UsD_vGgmnWr*Ux#7+#p$v>G#>Nw{fk3>SrO8uxpeuW23>qR=5<*mzZA*D?5b@A<(?_a!M{B^Mt)?EDb6{~mp zcD{7!B3%xYE}0>}_;nbs%TnpW{e4asKY#g4epwo7$1V%$;@dGlucV)G>-_$GMO;^L z^B3B!O6}N1(xtEb!gUwdzkDWL-jH8X$uEo>T)*#iDf3IH-|r+{4w5c8HR$5UV}6|M`+fgBas%mdx^&42 zwPW5detoQ-{PI)CFJWCu)s8KfE`6j+B-D<*C0%+-mle`wU$tYGN|$NU#kXUu&p5w) zss827kY9$Vf4NkC@%_HaUA=zxaMXPk!;b_ zE_;Uh{e0=dxXktY*UK-g7r6J3ydl5rC0#~`bfJH6?U?Uh{JhfZ5*{Z%5z@u4C;RtZ zc)!e9^U}iNMX$fTJk;;|`;`26%tlXA=5_J? zi=V%EUA$k;3bkW?Ug`UnuwO1vJJw3Npp)~9*Tv8G{dml;C%-OTykAC!+A*(7Zpbfw V{^I?zTWJ2$UVibq`1dLG{(r)D_zD03 literal 0 HcmV?d00001 diff --git a/modules/stochastic_tools/test/tests/userobjects/rb_inverse_mapping/main_rb.i b/modules/stochastic_tools/test/tests/userobjects/rb_inverse_mapping/main_rb.i new file mode 100644 index 000000000000..789e2f122a27 --- /dev/null +++ b/modules/stochastic_tools/test/tests/userobjects/rb_inverse_mapping/main_rb.i @@ -0,0 +1,92 @@ +S = 10 +D = 10 +L = 5 + +[Mesh] + type = GeneratedMesh + dim = 2 + nx = 50 + ny = 50 + xmax = ${L} +[] + +[Problem] + solve = false +[] + +[Variables] + [u] + [] +[] + +[Kernels] + [diffusion_u] + type = MatDiffusion + variable = u + diffusivity = D_u + [] + [source_u] + type = BodyForce + variable = u + value = 5.0 + [] +[] + +[Functions] + [du] + type = ParsedFunction + expression = 'D * D * x + 1' + symbol_names = D + symbol_values = ${D} + [] +[] + +[Materials] + [diffusivity_u] + type = GenericFunctionMaterial + prop_names = D_u + prop_values = du + [] +[] + +[BCs] + [left_u] + type = DirichletBC + variable = u + boundary = left + value = 0 + preset = true + [] + [right_u] + type = DirichletBC + variable = u + boundary = right + value = ${S} + preset = true + [] +[] + +[Executioner] + type = Steady + solve_type = NEWTON +[] + +[VariableMappings] + [rb_mapping] + type = DEIMRBMapping + filename = 'gold_mapping.rd' + [] +[] + +[UserObjects] + [im] + type = InverseRB + mapping = rb_mapping + execute_on = TIMESTEP_END + max_iter = 1 + [] +[] + +[Outputs] + exodus = true +[] diff --git a/modules/stochastic_tools/test/tests/userobjects/rb_inverse_mapping/tests b/modules/stochastic_tools/test/tests/userobjects/rb_inverse_mapping/tests new file mode 100644 index 000000000000..5e481bb87664 --- /dev/null +++ b/modules/stochastic_tools/test/tests/userobjects/rb_inverse_mapping/tests @@ -0,0 +1,12 @@ +[Tests] + issues = '#' + design = '' + + [p] + type = JSONDiff + requirement = 'The system should be able to print serialized snapshots in a distributed fashion in a json format.' + input = parallel_storage_main.i + jsondiff = 'parallel_storage_main_out.json' + + [] +[] diff --git a/modules/stochastic_tools/test/tests/variablemappings/deim_mapping/parallel_storage_main.i b/modules/stochastic_tools/test/tests/variablemappings/deim_mapping/parallel_storage_main.i new file mode 100644 index 000000000000..48ed490c4eeb --- /dev/null +++ b/modules/stochastic_tools/test/tests/variablemappings/deim_mapping/parallel_storage_main.i @@ -0,0 +1,125 @@ +[StochasticTools] +[] + +[Distributions] + [S_dist] + type = Uniform + lower_bound = 0 + upper_bound = 20 + [] + [D_dist] + type = Uniform + lower_bound = 0 + upper_bound = 20 + [] + [L_dist] + type = Uniform + lower_bound = 1 + upper_bound = 10 + [] +[] + +[Samplers] + [sample] + type = MonteCarlo + num_rows = 100 + distributions = 'S_dist D_dist L_dist' + execute_on = PRE_MULTIAPP_SETUP + # min_procs_per_row = 4 + [] +[] + +[VariableMappings] + [rb_mapping] + type = DEIMRBMapping + solution_storage = parallel_storage + variables = 'solution residual jacobian' + num_modes_to_compute = '30 10 10' + extra_slepc_options = "-svd_monitor_all" + jac_index_name = 'dummy_storage/indices' + [] +[] + +[MultiApps] + [worker] + type = SamplerFullSolveMultiApp + input_files = sub.i + sampler = sample + mode = batch-reset + [] +[] + +[Transfers] + + [solution_transfer] + type = SerializedSnapshotTransfer + parallel_storage = parallel_storage + from_multi_app = worker + sampler = sample + solution_container = solution_storage + residual_container = residual_storage + jacobian_container = jacobian_storage + serialize_on_root = true + [] + # We are only transferring the one copy of the indices. We assume that the + # order, number, anything does not change from one run to another. Breaking + # this assumption would not allow M/DEIM to work in the current configuration. + [jac_indices] + type = MultiAppReporterTransfer + from_multi_app = worker + from_reporters = 'jacobian_storage/indices' + to_reporters = 'dummy_storage/indices' + subapp_index = 0 + [] +[] + +[Controls] + [cmd_line] + type = MultiAppSamplerControl + multi_app = worker + sampler = sample + param_names = 'S D L' + [] +[] + +[Reporters] + [parallel_storage] + type = ParallelSolutionStorage + variables = 'solution residual jacobian' + outputs = none + [] + [build_mapping] + type = MappingReporter + mapping = rb_mapping + variables = "solution residual jacobian" + build_all_mappings_only = true + execute_on = TIMESTEP_END # To make sure the trainer sees the results on + outputs = out + [] + # Allow for the transfer of the indices. Creates the reporter type needed. + [dummy_storage] + type = JacobianContainer + tag_name = dummy + jac_indices_reporter_name = indices + execute_on = 'NONE' + outputs = none + [] +[] + +[Outputs] + [mapping] + type = MappingOutput + mappings = rb_mapping + execute_on = FINAL + [] + [out] + type = JSON + execute_on = TIMESTEP_END + execute_system_information_on = none + [] +[] + +[Problem] + extra_tag_matrices = 'dummy' +[] + diff --git a/modules/stochastic_tools/test/tests/variablemappings/deim_mapping/sub.i b/modules/stochastic_tools/test/tests/variablemappings/deim_mapping/sub.i new file mode 100644 index 000000000000..c470b44a6ef0 --- /dev/null +++ b/modules/stochastic_tools/test/tests/variablemappings/deim_mapping/sub.i @@ -0,0 +1,105 @@ +S = 10 +D = 10 +L = 5 + +[Mesh] + type = GeneratedMesh + dim = 1 + nx = 50 + xmax = ${L} +[] + +[Variables] + [u] + [] +[] + +[Kernels] + [diffusion_u] + type = MatDiffusion + variable = u + diffusivity = D_u + [] + [source_u] + type = BodyForce + variable = u + value = 1.0 + [] +[] + +[Functions] + [du] + type = ParsedFunction + expression = 'D * D * x + 1' + symbol_names = D + symbol_values = ${D} + [] +[] + +[Materials] + [diffusivity_u] + type = GenericFunctionMaterial + prop_names = D_u + prop_values = du + [] +[] + +[BCs] + [left_u] + type = DirichletBC + variable = u + boundary = left + value = 0 + preset = true + [] + [right_u] + type = DirichletBC + variable = u + boundary = right + value = ${S} + preset = true + [] +[] + +[Executioner] + type = Steady + solve_type = NEWTON + petsc_options_iname = '-pc_type -pc_factor_shift_type -pc_factor_mat_solver_type' + petsc_options_value = 'lu NONZERO strumpack' + nl_abs_tol = 1e-8 + nl_rel_tol = 1e-18 +[] + +[Controls] + [stochastic] + type = SamplerReceiver + [] +[] + +[Reporters] + [solution_storage] + type = SolutionContainer + execute_on = 'FINAL' + [] + [residual_storage] + type = ResidualContainer + tag_name = total + execute_on = 'TIMESTEP_BEGIN' + [] + [jacobian_storage] + type = JacobianContainer + tag_name = total + jac_indices_reporter_name = indices + execute_on = 'FINAL' + [] +[] + +[Problem] + extra_tag_vectors = 'total' + extra_tag_matrices = 'total' +[] + +[GlobalParams] + extra_matrix_tags = 'total' + extra_vector_tags = 'total' +[] diff --git a/modules/stochastic_tools/test/tests/variablemappings/deim_mapping/test_sub.i b/modules/stochastic_tools/test/tests/variablemappings/deim_mapping/test_sub.i new file mode 100644 index 000000000000..d1199ae461f9 --- /dev/null +++ b/modules/stochastic_tools/test/tests/variablemappings/deim_mapping/test_sub.i @@ -0,0 +1,91 @@ +S = 10 +D = 10 +L = 5 + +[Mesh] + type = GeneratedMesh + dim = 1 + nx = 50 + xmax = ${L} +[] + +[Problem] + solve = false +[] + +[Variables] + [u] + [] +[] + +[Kernels] + [diffusion_u] + type = MatDiffusion + variable = u + diffusivity = D_u + [] + [source_u] + type = BodyForce + variable = u + value = 1.0 + [] +[] + +[Functions] + [du] + type = ParsedFunction + expression = 'D * D * x + 1' + symbol_names = D + symbol_values = ${D} + [] +[] + +[Materials] + [diffusivity_u] + type = GenericFunctionMaterial + prop_names = D_u + prop_values = du + [] +[] + +[BCs] + [left_u] + type = DirichletBC + variable = u + boundary = left + value = 0 + preset = true + [] + [right_u] + type = DirichletBC + variable = u + boundary = right + value = ${S} + preset = true + [] +[] + +[Executioner] + type = Steady + solve_type = NEWTON +[] + +[VariableMappings] + [rb_mapping] + type = DEIMRBMapping + filename = 'parallel_storage_main_mapping_rb_mapping.rd' + [] +[] + +[UserObjects] + [im] + type = InverseRB + mapping = rb_mapping + execute_on = TIMESTEP_END + max_iter = 3 + [] +[] + +[Outputs] + exodus = true +[] diff --git a/modules/stochastic_tools/test/tests/variablemappings/deim_mapping/tests b/modules/stochastic_tools/test/tests/variablemappings/deim_mapping/tests new file mode 100644 index 000000000000..5e481bb87664 --- /dev/null +++ b/modules/stochastic_tools/test/tests/variablemappings/deim_mapping/tests @@ -0,0 +1,12 @@ +[Tests] + issues = '#' + design = '' + + [p] + type = JSONDiff + requirement = 'The system should be able to print serialized snapshots in a distributed fashion in a json format.' + input = parallel_storage_main.i + jsondiff = 'parallel_storage_main_out.json' + + [] +[] diff --git a/modules/stochastic_tools/test/tests/variablemappings/deim_mapping_load_stepping/parallel_storage_main.i b/modules/stochastic_tools/test/tests/variablemappings/deim_mapping_load_stepping/parallel_storage_main.i new file mode 100644 index 000000000000..0ca1e41f9335 --- /dev/null +++ b/modules/stochastic_tools/test/tests/variablemappings/deim_mapping_load_stepping/parallel_storage_main.i @@ -0,0 +1,120 @@ +[StochasticTools] +[] + +[Distributions] + [S_dist] + type = Uniform + lower_bound = 0 + upper_bound = 20 + [] + [D_dist] + type = Uniform + lower_bound = 0 + upper_bound = 20 + [] + [L_dist] + type = Uniform + lower_bound = 1 + upper_bound = 10 + [] +[] + +[Samplers] + [sample] + type = MonteCarlo + num_rows = 50 + distributions = 'S_dist D_dist L_dist' + execute_on = PRE_MULTIAPP_SETUP + # min_procs_per_row = 4 + [] +[] + +[VariableMappings] + [rb_mapping] + type = DEIMRBMapping + solution_storage = parallel_storage + variables = 'solution residual jacobian' + num_modes_to_compute = '30 4 10' + extra_slepc_options = "-svd_monitor_all" + jac_index_name = 'jacobian_storage/indices' + [] +[] + +[MultiApps] + [worker] + type = SamplerFullSolveMultiApp + input_files = sub.i + sampler = sample + mode = batch-reset + # min_procs_per_app = 4 + [] +[] + +[Transfers] + + [solution_transfer] + type = SerializedSnapshotTransfer + parallel_storage = parallel_storage + from_multi_app = worker + sampler = sample + solution_container = solution_storage + residual_container = residual_storage + jacobian_container = jacobian_storage + serialize_on_root = true + [] + # We are only transferring the one copy of the indices. We assume that the + # order, number, anything does not change from one run to another. Breaking + # this assumption would not allow M/DEIM to work in the current configuration. + [jac_indices] + type = MultiAppReporterTransfer + from_multi_app = worker + from_reporters = 'jacobian_storage/indices' + to_reporters = 'jacobian_storage/indices' + subapp_index = 0 + [] +[] + +[Controls] + [cmd_line] + type = MultiAppSamplerControl + multi_app = worker + sampler = sample + param_names = 'S D L' + [] +[] + +[Reporters] + [parallel_storage] + type = ParallelSolutionStorage + variables = 'solution residual jacobian' + outputs = none + [] + [build_mapping] + type = MappingReporter + mapping = rb_mapping + variables = "solution residual jacobian" + build_all_mappings_only = true + execute_on = TIMESTEP_END # To make sure the trainer sees the results on + outputs = out + [] + # Allow for the transfer of the indices. Create the reporter type needed. + [jacobian_storage] + type = JacobianContainer + jac_indices_reporter_name = indices + execute_on = 'NONE' + outputs = none + [] +[] + +[Outputs] + [mapping] + type = MappingOutput + mappings = rb_mapping + execute_on = FINAL + [] + [out] + type = JSON + execute_on = TIMESTEP_END + execute_system_information_on = none + [] +[] diff --git a/modules/stochastic_tools/test/tests/variablemappings/deim_mapping_load_stepping/sub.i b/modules/stochastic_tools/test/tests/variablemappings/deim_mapping_load_stepping/sub.i new file mode 100644 index 000000000000..61cb8d629939 --- /dev/null +++ b/modules/stochastic_tools/test/tests/variablemappings/deim_mapping_load_stepping/sub.i @@ -0,0 +1,116 @@ +S = 10 +D = 10 +L = 5 + +[Mesh] + [gmg] + type = GeneratedMeshGenerator + dim = 2 + nx = 50 + ny = 10 + xmax = ${L} + ymin = -1 + ymax = 1 + + [] + +[] + +[Variables] + [u] + [] +[] + +[Kernels] + [diffusion_u] + type = ADMatDiffusion + variable = u + diffusivity = D_u + [] + [source_u] + type = ADBodyForce + variable = u + value = 1.0 + [] +[] + +[Functions] + [du] + type = ParsedFunction + expression = 'D * x^2 + 1' + symbol_names = D + symbol_values = ${D} + [] +[] + +[Materials] + [diffusivity_u] + type = ADGenericFunctionMaterial + prop_names = D_u + prop_values = du + [] +[] + +[BCs] + [left_u] + type = DirichletBC + variable = u + boundary = 'left top bottom' + value = 0 + [] + [right_u] + type = FunctionDirichletBC + variable = u + boundary = right + function = '${S} * t^2 * -(y^2-1)' + preset = false + [] +[] + +[Executioner] + type = Transient + dt = 0.25 + num_steps = 4 + solve_type = NEWTON + petsc_options_iname = '-pc_type -pc_factor_shift_type -pc_factor_mat_solver_type' + petsc_options_value = 'lu NONZERO strumpack' + +[] + +[Controls] + [stochastic] + type = SamplerReceiver + [] +[] + +[Reporters] + [solution_storage] + type = SolutionContainer + execute_on = 'TIMESTEP_END' + [] + [residual_storage] + type = ResidualContainer + tag_name = 'total' + execute_on = 'TIMESTEP_BEGIN' + [] + [jacobian_storage] + type = JacobianContainer + tag_name = 'total' + jac_indices_reporter_name = indices + execute_on = 'TIMESTEP_END' + [] +[] + +[Outputs] + console = true +[] + +[Problem] + extra_tag_vectors = 'total' + extra_tag_matrices = 'total' +[] + +[GlobalParams] + extra_matrix_tags = 'total' + extra_vector_tags = 'total' +[] diff --git a/modules/stochastic_tools/test/tests/variablemappings/deim_mapping_load_stepping/test_sub.i b/modules/stochastic_tools/test/tests/variablemappings/deim_mapping_load_stepping/test_sub.i new file mode 100644 index 000000000000..40ce2d1b8f0c --- /dev/null +++ b/modules/stochastic_tools/test/tests/variablemappings/deim_mapping_load_stepping/test_sub.i @@ -0,0 +1,98 @@ +S = 10 +D = 10 +L = 5 + +[Mesh] + [gmg] + type = GeneratedMeshGenerator + dim = 2 + nx = 50 + ny = 10 + xmax = ${L} + ymin = -1 + ymax = 1 + [] + +[] + +[Problem] + solve = false +[] + +[Variables] + [u] + [] +[] + +[Kernels] + [diffusion_u] + type = ADMatDiffusion + variable = u + diffusivity = D_u + [] + [source_u] + type = ADBodyForce + variable = u + value = 1.0 + [] +[] + +[Functions] + [du] + type = ParsedFunction + expression = 'D * x^2 + 1' + symbol_names = D + symbol_values = ${D} + [] +[] + +[Materials] + [diffusivity_u] + type = ADGenericFunctionMaterial + prop_names = D_u + prop_values = du + [] +[] + +[BCs] + [left_u] + type = DirichletBC + variable = u + boundary = 'left top bottom' + value = 0 + [] + [right_u] + type = FunctionDirichletBC + variable = u + boundary = right + function = '${S} * t^2 * -(y^2-1)' + preset = false + [] +[] + +[Executioner] + type = Transient + dt = 0.25 + num_steps = 4 + solve_type = NEWTON +[] + +[VariableMappings] + [rb_mapping] + type = DEIMRBMapping + filename = 'parallel_storage_main_mapping_rb_mapping.rd' + [] +[] + +[UserObjects] + [im] + type = InverseRB + mapping = rb_mapping + execute_on = TIMESTEP_END + max_iter = 2 + [] +[] + +[Outputs] + exodus = true +[] diff --git a/modules/stochastic_tools/test/tests/variablemappings/deim_mapping_load_stepping/tests b/modules/stochastic_tools/test/tests/variablemappings/deim_mapping_load_stepping/tests new file mode 100644 index 000000000000..5e481bb87664 --- /dev/null +++ b/modules/stochastic_tools/test/tests/variablemappings/deim_mapping_load_stepping/tests @@ -0,0 +1,12 @@ +[Tests] + issues = '#' + design = '' + + [p] + type = JSONDiff + requirement = 'The system should be able to print serialized snapshots in a distributed fashion in a json format.' + input = parallel_storage_main.i + jsondiff = 'parallel_storage_main_out.json' + + [] +[] diff --git a/modules/stochastic_tools/test/tests/variablemappings/pod_mapping/pod_mapping_main.i b/modules/stochastic_tools/test/tests/variablemappings/pod_mapping/pod_mapping_main.i index 4096635821f7..295d6c44ea73 100644 --- a/modules/stochastic_tools/test/tests/variablemappings/pod_mapping/pod_mapping_main.i +++ b/modules/stochastic_tools/test/tests/variablemappings/pod_mapping/pod_mapping_main.i @@ -4,21 +4,17 @@ [Distributions] [S_dist] type = Uniform - lower_bound = 0 - upper_bound = 10 - [] - [D_dist] - type = Uniform - lower_bound = 0 + lower_bound = 9 upper_bound = 10 [] + [] [Samplers] [sample] type = MonteCarlo - num_rows = 4 - distributions = 'S_dist D_dist' + num_rows = 8 + distributions = 'S_dist' execute_on = initial min_procs_per_row = 2 [] @@ -39,18 +35,18 @@ type = PODMapping solution_storage = parallel_storage variables = "v" - num_modes_to_compute = '2' + num_modes_to_compute = '1' extra_slepc_options = "-svd_monitor_all" [] [] [Transfers] - [param_transfer] - type = SamplerParameterTransfer - to_multi_app = worker - sampler = sample - parameters = 'Kernels/source_v/value BCs/right_v/value' - [] + # [param_transfer] + # type = SamplerParameterTransfer + # to_multi_app = worker + # sampler = sample + # parameters = 'Kernels/source_v/value BCs/right_v/value' + # [] [solution_transfer] type = SerializedSolutionTransfer parallel_storage = parallel_storage diff --git a/modules/stochastic_tools/unit/src/TestQDEIM.C b/modules/stochastic_tools/unit/src/TestQDEIM.C new file mode 100644 index 000000000000..f044694c66f7 --- /dev/null +++ b/modules/stochastic_tools/unit/src/TestQDEIM.C @@ -0,0 +1,84 @@ +//* This file is part of the MOOSE framework +//* https://www.mooseframework.org +//* +//* All rights reserved, see COPYRIGHT for full restrictions +//* https://github.com/idaholab/moose/blob/master/COPYRIGHT +//* +//* Licensed under LGPL 2.1, please see LICENSE for details +//* https://www.gnu.org/licenses/lgpl-2.1.html +#include "gtest/gtest.h" + +#include "QDEIM.h" +#include "libmesh/dense_vector.h" +#include "libmesh/id_types.h" + +TEST(StochasticTools, QDEIM) +{ + // Here is the gold U file + DenseVector U1 = {-0.064096, -0.40913, -0.16302, -0.02448, -0.89517}; + DenseVector U2 = {-0.39344, -0.53986, -0.63356, 0.023065, 0.38965}; + DenseVector U3 = {-0.80489, 0.35679, 0.1099, -0.44717, -0.11322}; + DenseVector U4 = {0.43825, 0.15973, -0.43841, -0.76825, -0.0035379}; + std::vector> gold_U = {U1, U2, U3, U4}; + + // Compute QDEIM + std::vector> new_basis; + std::vector sample_indices; + QDEIM qd(gold_U); + qd.computeSelectionAndBasis(&new_basis, sample_indices); + + // Gold results file + std::vector gold_s = {0, 4, 3, 2}; + DenseVector M1 = {1, 0, 0, 0}; + DenseVector M2 = {0.0552399996, 0.2959201626, -0.7330802260, 0.9731093669}; + DenseVector M3 = {0, 0, 0, 1}; + DenseVector M4 = {0, 0, 1, 0}; + DenseVector M5 = {0, 1, 0, 0}; + std::vector> gold_M = {M1, M2, M3, M4, M5}; + + // Test Results + ASSERT_EQ(new_basis.size(), gold_M.size()); + + for (size_t i = 0; i < new_basis.size(); ++i) + { + ASSERT_EQ(new_basis[i].size(), gold_M[i].size()); + for (size_t j = 0; j < new_basis[i].size(); ++j) + { + EXPECT_NEAR(new_basis[i](j), gold_M[i](j), 1e-6); + } + } + + // Check if sample_indices and gold_s have the same size + ASSERT_EQ(sample_indices.size(), gold_s.size()); + + for (size_t i = 0; i < sample_indices.size(); ++i) + { + EXPECT_EQ(sample_indices[i], gold_s[i]); + } +} + +TEST(StochasticTools, QDEIM_Selection_Only) +{ + // Here is the gold U file + DenseVector U1 = {-0.064096, -0.40913, -0.16302, -0.02448, -0.89517}; + DenseVector U2 = {-0.39344, -0.53986, -0.63356, 0.023065, 0.38965}; + DenseVector U3 = {-0.80489, 0.35679, 0.1099, -0.44717, -0.11322}; + DenseVector U4 = {0.43825, 0.15973, -0.43841, -0.76825, -0.0035379}; + std::vector> gold_U = {U1, U2, U3, U4}; + + // Compute QDEIM + std::vector sample_indices; + QDEIM qd(gold_U); + qd.computeSelection(sample_indices); + + // Gold results file + std::vector gold_s = {0, 4, 3, 2}; + + // Check if sample_indices and gold_s have the same size + ASSERT_EQ(sample_indices.size(), gold_s.size()); + + for (size_t i = 0; i < sample_indices.size(); ++i) + { + EXPECT_EQ(sample_indices[i], gold_s[i]); + } +} diff --git a/petsc b/petsc index 1580af9fcba3..9a61b8cbbacd 160000 --- a/petsc +++ b/petsc @@ -1 +1 @@ -Subproject commit 1580af9fcba34756b20ac647dbc5659db3fdd05b +Subproject commit 9a61b8cbbacd10b14c2a18c9bc314b0ecb8ec018 diff --git a/test/tests/outputs/csv/csv_transient_vpp.i b/test/tests/outputs/csv/csv_transient_vpp.i index 0d145cf42e76..e53165042d8b 100644 --- a/test/tests/outputs/csv/csv_transient_vpp.i +++ b/test/tests/outputs/csv/csv_transient_vpp.i @@ -72,6 +72,7 @@ [Outputs] [out] type = CSV - execute_on = 'LINEAR' + # execute_on = 'LINEAR' + execute_on = 'final' [] []