From 2eb39750b90286b089526eaf4e652f1e272643e0 Mon Sep 17 00:00:00 2001 From: Sidafa Conde Date: Tue, 4 Jun 2019 08:21:45 -0600 Subject: [PATCH 1/9] Sidafe's Branch Squashed in Place * [WIP] Tempus: consolidate stepper observers * [WIP] Tempus: BDF2, IMEX_RK composite observer * Tempus: leapfrog stepper observer * Tempus: Observer - add RK observer for all MSRK stepper type * Tempus: observer - consolidate IMEX RK into RK observers * Tempus: observer - Partitioned IMEX RK now using RKObserver * Tempus: observer - remove redundant/speciliazed rk stpper observers * Tempus: observer - re-order the method definition * to match the order they are called by the stepper * Tempus: observer - add extra observer method to erk/dirk * Tempus: observer - add rk observer tests * Tempus: observer - add stepper rk observer logging * Tempus: getObserver - access observer * add getter to access the stepper observer * [WIP] Tempus: stepper - composite observer * guard the `setObserver` to only add obsever if the composite is empty. * This will ensure that the same observer is not added multiple times --- packages/tempus/src/Tempus_Integrator.hpp | 2 +- ...mpus_IntegratorAdjointSensitivity_decl.hpp | 2 +- ...mpus_IntegratorAdjointSensitivity_impl.hpp | 2 +- .../src/Tempus_IntegratorBasic_decl.hpp | 2 +- ...mpus_IntegratorForwardSensitivity_decl.hpp | 2 +- ...PseudoTransientAdjointSensitivity_decl.hpp | 2 +- ...PseudoTransientAdjointSensitivity_impl.hpp | 2 +- ...PseudoTransientForwardSensitivity_decl.hpp | 2 +- ...PseudoTransientForwardSensitivity_impl.hpp | 2 +- .../tempus/src/Tempus_SolutionState_decl.hpp | 2 +- .../tempus/src/Tempus_StepperBDF2_decl.hpp | 11 +- .../tempus/src/Tempus_StepperBDF2_impl.hpp | 29 ++-- .../src/Tempus_StepperBackwardEuler_decl.hpp | 3 + .../tempus/src/Tempus_StepperDIRKObserver.hpp | 82 ---------- .../tempus/src/Tempus_StepperDIRK_decl.hpp | 9 +- .../tempus/src/Tempus_StepperDIRK_impl.hpp | 59 +++---- .../src/Tempus_StepperExplicitRKObserver.hpp | 73 --------- ...tepperExplicitRKObserverComposite_decl.hpp | 79 ---------- ...tepperExplicitRKObserverComposite_impl.hpp | 84 ---------- .../src/Tempus_StepperExplicitRK_decl.hpp | 15 +- .../src/Tempus_StepperExplicitRK_impl.hpp | 57 ++++--- .../src/Tempus_StepperForwardEuler_decl.hpp | 3 + .../src/Tempus_StepperHHTAlpha_decl.hpp | 3 + .../src/Tempus_StepperIMEX_RKPartObserver.hpp | 87 ---------- .../Tempus_StepperIMEX_RK_Partition_decl.hpp | 7 +- .../Tempus_StepperIMEX_RK_Partition_impl.hpp | 43 ++--- .../tempus/src/Tempus_StepperIMEX_RK_decl.hpp | 7 +- .../tempus/src/Tempus_StepperIMEX_RK_impl.hpp | 45 ++---- .../src/Tempus_StepperLeapfrog_decl.hpp | 5 + .../src/Tempus_StepperLeapfrog_impl.hpp | 29 ++-- ...empus_StepperNewmarkExplicitAForm_decl.hpp | 3 + ...empus_StepperNewmarkImplicitAForm_decl.hpp | 3 + ...empus_StepperNewmarkImplicitDForm_decl.hpp | 3 + .../Tempus_StepperObserverComposite_decl.hpp | 4 +- .../Tempus_StepperObserverComposite_impl.hpp | 4 - .../src/Tempus_StepperOperatorSplit_decl.hpp | 8 + .../src/Tempus_StepperOperatorSplit_impl.hpp | 2 +- .../src/Tempus_StepperRKButcherTableau.hpp | 56 +++---- ...erver.hpp => Tempus_StepperRKObserver.hpp} | 47 +++--- ... => Tempus_StepperRKObserverComposite.cpp} | 6 +- ...Tempus_StepperRKObserverComposite_decl.hpp | 97 ++++++++++++ ...Tempus_StepperRKObserverComposite_impl.hpp | 106 +++++++++++++ .../src/Tempus_StepperRKObserverLogging.cpp | 19 +++ .../Tempus_StepperRKObserverLogging_decl.hpp | 112 +++++++++++++ .../Tempus_StepperRKObserverLogging_impl.hpp | 124 +++++++++++++++ ...tepperStaggeredForwardSensitivity_decl.hpp | 3 + .../src/Tempus_StepperTrapezoidal_decl.hpp | 3 + .../src/Tempus_StepperTrapezoidal_impl.hpp | 4 +- packages/tempus/src/Tempus_Stepper_decl.hpp | 4 + .../test/Observer/Tempus_ObserverTest.cpp | 149 +++++++++++++++++- .../test/Observer/Tempus_Observer_SinCos.xml | 6 + ...s_PhysicsStateTest_StepperForwardEuler.hpp | 2 + .../Tempus_UnitTest_ERK_ForwardEuler.cpp | 6 +- 53 files changed, 871 insertions(+), 650 deletions(-) delete mode 100644 packages/tempus/src/Tempus_StepperDIRKObserver.hpp delete mode 100644 packages/tempus/src/Tempus_StepperExplicitRKObserver.hpp delete mode 100644 packages/tempus/src/Tempus_StepperExplicitRKObserverComposite_decl.hpp delete mode 100644 packages/tempus/src/Tempus_StepperExplicitRKObserverComposite_impl.hpp delete mode 100644 packages/tempus/src/Tempus_StepperIMEX_RKPartObserver.hpp rename packages/tempus/src/{Tempus_StepperIMEX_RKObserver.hpp => Tempus_StepperRKObserver.hpp} (61%) rename packages/tempus/src/{Tempus_StepperExplicitRKObserverComposite.cpp => Tempus_StepperRKObserverComposite.cpp} (70%) create mode 100644 packages/tempus/src/Tempus_StepperRKObserverComposite_decl.hpp create mode 100644 packages/tempus/src/Tempus_StepperRKObserverComposite_impl.hpp create mode 100644 packages/tempus/src/Tempus_StepperRKObserverLogging.cpp create mode 100644 packages/tempus/src/Tempus_StepperRKObserverLogging_decl.hpp create mode 100644 packages/tempus/src/Tempus_StepperRKObserverLogging_impl.hpp diff --git a/packages/tempus/src/Tempus_Integrator.hpp b/packages/tempus/src/Tempus_Integrator.hpp index 863b3b0b0123..8632350c9faa 100644 --- a/packages/tempus/src/Tempus_Integrator.hpp +++ b/packages/tempus/src/Tempus_Integrator.hpp @@ -73,7 +73,7 @@ class Integrator /// Get current time virtual Scalar getTime() const = 0; /// Get current index - virtual Scalar getIndex() const = 0; + virtual int getIndex() const = 0; /// Get the Status virtual Tempus::Status getStatus() const = 0; /// Get the stepper diff --git a/packages/tempus/src/Tempus_IntegratorAdjointSensitivity_decl.hpp b/packages/tempus/src/Tempus_IntegratorAdjointSensitivity_decl.hpp index 7138f9145502..0e2648cec3f4 100644 --- a/packages/tempus/src/Tempus_IntegratorAdjointSensitivity_decl.hpp +++ b/packages/tempus/src/Tempus_IntegratorAdjointSensitivity_decl.hpp @@ -96,7 +96,7 @@ class IntegratorAdjointSensitivity : /// Get current time virtual Scalar getTime() const override; /// Get current index - virtual Scalar getIndex() const override; + virtual int getIndex() const override; /// Get Status virtual Status getStatus() const override; /// Get the Stepper diff --git a/packages/tempus/src/Tempus_IntegratorAdjointSensitivity_impl.hpp b/packages/tempus/src/Tempus_IntegratorAdjointSensitivity_impl.hpp index cc49e1d48a1a..42010e98a465 100644 --- a/packages/tempus/src/Tempus_IntegratorAdjointSensitivity_impl.hpp +++ b/packages/tempus/src/Tempus_IntegratorAdjointSensitivity_impl.hpp @@ -227,7 +227,7 @@ getTime() const } template -Scalar +int IntegratorAdjointSensitivity:: getIndex() const { diff --git a/packages/tempus/src/Tempus_IntegratorBasic_decl.hpp b/packages/tempus/src/Tempus_IntegratorBasic_decl.hpp index 572d93b74ed5..9403253626f4 100644 --- a/packages/tempus/src/Tempus_IntegratorBasic_decl.hpp +++ b/packages/tempus/src/Tempus_IntegratorBasic_decl.hpp @@ -89,7 +89,7 @@ class IntegratorBasic : virtual public Tempus::Integrator virtual Scalar getTime() const override {return solutionHistory_->getCurrentTime();} /// Get current index - virtual Scalar getIndex() const override + virtual int getIndex() const override {return solutionHistory_->getCurrentIndex();} /// Get Status virtual Status getStatus() const override diff --git a/packages/tempus/src/Tempus_IntegratorForwardSensitivity_decl.hpp b/packages/tempus/src/Tempus_IntegratorForwardSensitivity_decl.hpp index 5cc7b7a5b845..690a981f07e4 100644 --- a/packages/tempus/src/Tempus_IntegratorForwardSensitivity_decl.hpp +++ b/packages/tempus/src/Tempus_IntegratorForwardSensitivity_decl.hpp @@ -128,7 +128,7 @@ class IntegratorForwardSensitivity : virtual public Tempus::Integrator virtual Scalar getTime() const override { return integrator_->getTime(); } /// Get current index - virtual Scalar getIndex() const override + virtual int getIndex() const override { return integrator_->getIndex(); } /// Get Status virtual Status getStatus() const override diff --git a/packages/tempus/src/Tempus_IntegratorPseudoTransientAdjointSensitivity_decl.hpp b/packages/tempus/src/Tempus_IntegratorPseudoTransientAdjointSensitivity_decl.hpp index 43b0e99a055f..4936899c1acd 100644 --- a/packages/tempus/src/Tempus_IntegratorPseudoTransientAdjointSensitivity_decl.hpp +++ b/packages/tempus/src/Tempus_IntegratorPseudoTransientAdjointSensitivity_decl.hpp @@ -97,7 +97,7 @@ class IntegratorPseudoTransientAdjointSensitivity : /// Get current time virtual Scalar getTime() const override; /// Get current index - virtual Scalar getIndex() const override; + virtual int getIndex() const override; /// Get Status virtual Status getStatus() const override; /// Get the Stepper diff --git a/packages/tempus/src/Tempus_IntegratorPseudoTransientAdjointSensitivity_impl.hpp b/packages/tempus/src/Tempus_IntegratorPseudoTransientAdjointSensitivity_impl.hpp index 9d9a90752d32..afd5149b193f 100644 --- a/packages/tempus/src/Tempus_IntegratorPseudoTransientAdjointSensitivity_impl.hpp +++ b/packages/tempus/src/Tempus_IntegratorPseudoTransientAdjointSensitivity_impl.hpp @@ -110,7 +110,7 @@ getTime() const } template -Scalar +int IntegratorPseudoTransientAdjointSensitivity:: getIndex() const { diff --git a/packages/tempus/src/Tempus_IntegratorPseudoTransientForwardSensitivity_decl.hpp b/packages/tempus/src/Tempus_IntegratorPseudoTransientForwardSensitivity_decl.hpp index 68e59a926766..e54cbf694f1c 100644 --- a/packages/tempus/src/Tempus_IntegratorPseudoTransientForwardSensitivity_decl.hpp +++ b/packages/tempus/src/Tempus_IntegratorPseudoTransientForwardSensitivity_decl.hpp @@ -111,7 +111,7 @@ class IntegratorPseudoTransientForwardSensitivity : /// Get current time virtual Scalar getTime() const override; /// Get current index - virtual Scalar getIndex() const override; + virtual int getIndex() const override; /// Get Status virtual Status getStatus() const override; /// Get the Stepper diff --git a/packages/tempus/src/Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp b/packages/tempus/src/Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp index 44ed8c05ad82..704c9c5105e8 100644 --- a/packages/tempus/src/Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp +++ b/packages/tempus/src/Tempus_IntegratorPseudoTransientForwardSensitivity_impl.hpp @@ -121,7 +121,7 @@ getTime() const } template -Scalar +int IntegratorPseudoTransientForwardSensitivity:: getIndex() const { diff --git a/packages/tempus/src/Tempus_SolutionState_decl.hpp b/packages/tempus/src/Tempus_SolutionState_decl.hpp index 48892bf05653..d64022521f45 100644 --- a/packages/tempus/src/Tempus_SolutionState_decl.hpp +++ b/packages/tempus/src/Tempus_SolutionState_decl.hpp @@ -126,7 +126,7 @@ class SolutionState : return metaData_nc_; } virtual Scalar getTime() const {return metaData_->getTime();} - virtual Scalar getIndex() const {return metaData_->getIStep();} + virtual int getIndex() const {return metaData_->getIStep();} virtual Scalar getTimeStep() const {return metaData_->getDt();} virtual Scalar getOrder() const {return metaData_->getOrder();} virtual Scalar getNRunningFailures() const diff --git a/packages/tempus/src/Tempus_StepperBDF2_decl.hpp b/packages/tempus/src/Tempus_StepperBDF2_decl.hpp index a496c282ac96..d1b11c7d4a05 100644 --- a/packages/tempus/src/Tempus_StepperBDF2_decl.hpp +++ b/packages/tempus/src/Tempus_StepperBDF2_decl.hpp @@ -11,6 +11,7 @@ #include "Tempus_StepperImplicit.hpp" #include "Tempus_WrapperModelEvaluator.hpp" +#include "Tempus_StepperObserverComposite.hpp" #include "Tempus_StepperBDF2Observer.hpp" @@ -77,6 +78,9 @@ class StepperBDF2 : virtual public Tempus::StepperImplicit virtual void setObserver( Teuchos::RCP > obs = Teuchos::null); + virtual Teuchos::RCP > getObserver() const + { return this->stepperObserver_; } + /// Set the stepper to use in first step void setStartUpStepper(std::string startupStepperType = "DIRK 1 Stage Theta Method"); @@ -131,9 +135,10 @@ class StepperBDF2 : virtual public Tempus::StepperImplicit private: - Teuchos::RCP > startUpStepper_; - Teuchos::RCP > stepperBDF2Observer_; - Scalar order_; + Teuchos::RCP > startUpStepper_; + Teuchos::RCP > stepperObserver_; + Teuchos::RCP > stepperBDF2Observer_; + Scalar order_; }; /** \brief Time-derivative interface for BDF2. diff --git a/packages/tempus/src/Tempus_StepperBDF2_impl.hpp b/packages/tempus/src/Tempus_StepperBDF2_impl.hpp index 616711bbd96c..e5ae51bdf0a6 100644 --- a/packages/tempus/src/Tempus_StepperBDF2_impl.hpp +++ b/packages/tempus/src/Tempus_StepperBDF2_impl.hpp @@ -101,20 +101,17 @@ template void StepperBDF2::setObserver( Teuchos::RCP > obs) { - if (obs == Teuchos::null) { - // Create default observer, otherwise keep current observer. - if (this->stepperObserver_ == Teuchos::null) { - stepperBDF2Observer_ = - Teuchos::rcp(new StepperBDF2Observer()); - this->stepperObserver_ = - Teuchos::rcp_dynamic_cast > - (stepperBDF2Observer_,true); - } - } else { - this->stepperObserver_ = obs; - stepperBDF2Observer_ = - Teuchos::rcp_dynamic_cast >(this->stepperObserver_,true); - } + + if (this->stepperObserver_ == Teuchos::null) + this->stepperObserver_ = + Teuchos::rcp(new StepperObserverComposite()); + + if (( obs == Teuchos::null ) and (this->stepperObserver_->getSize() == 0) ) + obs = Teuchos::rcp(new StepperBDF2Observer()); + + this->stepperObserver_->addObserver( + Teuchos::rcp_dynamic_cast > (obs, true) ); + } @@ -184,7 +181,7 @@ void StepperBDF2::takeStep( //IKT, FIXME: add error checking regarding states being consecutive and //whether interpolated states are OK to use. - this->stepperObserver_->observeBeginTakeStep(solutionHistory, *this); + //this->stepperObserver_->observeBeginTakeStep(solutionHistory, *this); RCP > workingState=solutionHistory->getWorkingState(); RCP > currentState=solutionHistory->getCurrentState(); @@ -225,7 +222,7 @@ void StepperBDF2::takeStep( workingState->setSolutionStatus(sStatus); // Converged --> pass. workingState->setOrder(getOrder()); - this->stepperObserver_->observeEndTakeStep(solutionHistory, *this); + //this->stepperObserver_->observeEndTakeStep(solutionHistory, *this); } return; } diff --git a/packages/tempus/src/Tempus_StepperBackwardEuler_decl.hpp b/packages/tempus/src/Tempus_StepperBackwardEuler_decl.hpp index 4c1880fe3fba..04d25f4781fd 100644 --- a/packages/tempus/src/Tempus_StepperBackwardEuler_decl.hpp +++ b/packages/tempus/src/Tempus_StepperBackwardEuler_decl.hpp @@ -61,6 +61,9 @@ class StepperBackwardEuler : virtual void setObserver( Teuchos::RCP > obs = Teuchos::null); + virtual Teuchos::RCP > getObserver() const + { return this->stepperBEObserver_; } + /// Set the predictor void setPredictor(std::string predictorType = "None"); void setPredictor(Teuchos::RCP > predictorStepper); diff --git a/packages/tempus/src/Tempus_StepperDIRKObserver.hpp b/packages/tempus/src/Tempus_StepperDIRKObserver.hpp deleted file mode 100644 index 883478bb4b3f..000000000000 --- a/packages/tempus/src/Tempus_StepperDIRKObserver.hpp +++ /dev/null @@ -1,82 +0,0 @@ -// @HEADER -// **************************************************************************** -// Tempus: Copyright (2017) Sandia Corporation -// -// Distributed under BSD 3-clause license (See accompanying file Copyright.txt) -// **************************************************************************** -// @HEADER - -#ifndef Tempus_StepperDIRKObserver_hpp -#define Tempus_StepperDIRKObserver_hpp - -#include "Tempus_SolutionHistory.hpp" - - -namespace Tempus { - -// Forward Declaration for recursive includes (this Observer <--> Stepper) -template class StepperDIRK; - -/** \brief StepperDIRKObserver class for StepperDIRK. - * - * This is a means for application developers to perform tasks - * during the time steps, e.g., - * - Compute specific quantities - * - Output information - * - "Massage" the working solution state - * - ... - * - * Design Considerations - * - StepperDIRKObserver is not stateless! Developers may touch the - * solution state! Developers need to be careful not to break the - * restart (checkpoint) capability. - */ -template -class StepperDIRKObserver - : virtual public Tempus::StepperObserver -{ -public: - - /// Constructor - StepperDIRKObserver(){} - - /// Destructor - virtual ~StepperDIRKObserver(){} - - /// Observe Stepper at beginning of takeStep. - virtual void observeBeginTakeStep( - Teuchos::RCP > /* sh */, - Stepper & /* stepper */){} - - /// Observe Stepper at beginning of each stage. - virtual void observeBeginStage( - Teuchos::RCP > /* sh */, - StepperDIRK & /* stepperDIRK */){} - - /// Observe Stepper before Explicit evaluation of Implicit ODE ME. - virtual void observeBeforeExplicit( - Teuchos::RCP > /* sh */, - StepperDIRK & /* stepperDIRK */){} - - /// Observe Stepper before nonlinear solve. - virtual void observeBeforeSolve( - Teuchos::RCP > /* sh */, - StepperDIRK & /* stepperDIRK */){} - - /// Observe Stepper after nonlinear solve. - virtual void observeAfterSolve( - Teuchos::RCP > /* sh */, - StepperDIRK & /* stepperDIRK */){} - - /// Observe Stepper at end of each stage. - virtual void observeEndStage( - Teuchos::RCP > /* sh */, - StepperDIRK & /* stepperDIRK */){} - - /// Observe Stepper at end of takeStep. - virtual void observeEndTakeStep( - Teuchos::RCP > /* sh */, - Stepper & /* stepper */){} -}; -} // namespace Tempus -#endif // Tempus_StepperDIRKObserver_hpp diff --git a/packages/tempus/src/Tempus_StepperDIRK_decl.hpp b/packages/tempus/src/Tempus_StepperDIRK_decl.hpp index 82d817a863bf..f988c491b13a 100644 --- a/packages/tempus/src/Tempus_StepperDIRK_decl.hpp +++ b/packages/tempus/src/Tempus_StepperDIRK_decl.hpp @@ -13,7 +13,7 @@ #include "Tempus_RKButcherTableau.hpp" #include "Tempus_StepperImplicit.hpp" #include "Tempus_WrapperModelEvaluator.hpp" -#include "Tempus_StepperDIRKObserver.hpp" +#include "Tempus_StepperRKObserverComposite.hpp" namespace Tempus { @@ -96,6 +96,9 @@ class StepperDIRK : virtual public Tempus::StepperImplicit virtual Teuchos::RCP > getTableau() { return tableau_; } + virtual Teuchos::RCP > getObserver() const + { return this->stepperObserver_; } + /// Initialize during construction and after changing input parameters. virtual void initialize(); @@ -175,7 +178,7 @@ class StepperDIRK : virtual public Tempus::StepperImplicit /// Setup for constructor. virtual void setup( const Teuchos::RCP >& wrapperModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, const Teuchos::RCP >& solver, bool useFSAL, std::string ICConsistency, @@ -191,7 +194,7 @@ class StepperDIRK : virtual public Tempus::StepperImplicit Teuchos::RCP > stageX_; Teuchos::RCP > xTilde_; - Teuchos::RCP > stepperDIRKObserver_; + Teuchos::RCP > stepperObserver_; // For Embedded RK bool useEmbedded_; diff --git a/packages/tempus/src/Tempus_StepperDIRK_impl.hpp b/packages/tempus/src/Tempus_StepperDIRK_impl.hpp index adaa9bbd8fc7..9b694214a2cb 100644 --- a/packages/tempus/src/Tempus_StepperDIRK_impl.hpp +++ b/packages/tempus/src/Tempus_StepperDIRK_impl.hpp @@ -31,14 +31,14 @@ void StepperDIRK::setupDefault() this->setUseEmbedded( this->getUseEmbeddedDefault()); this->setZeroInitialGuess( false); - stepperDIRKObserver_ = Teuchos::rcp(new StepperDIRKObserver()); + stepperObserver_ = Teuchos::rcp(new StepperRKObserverComposite()); } template void StepperDIRK::setup( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, const Teuchos::RCP >& solver, bool useFSAL, std::string ICConsistency, @@ -52,7 +52,7 @@ void StepperDIRK::setup( this->setUseEmbedded( useEmbedded); this->setZeroInitialGuess( zeroInitialGuess); - stepperDIRKObserver_ = Teuchos::rcp(new StepperDIRKObserver()); + stepperObserver_ = Teuchos::rcp(new StepperRKObserverComposite()); this->setObserver(obs); if (appModel != Teuchos::null) { @@ -86,20 +86,17 @@ template void StepperDIRK::setObserver( Teuchos::RCP > obs) { - if (obs == Teuchos::null) { - // Create default observer, otherwise keep current observer. - if (this->stepperObserver_ == Teuchos::null) { - stepperDIRKObserver_ = - Teuchos::rcp(new StepperDIRKObserver()); - this->stepperObserver_ = - Teuchos::rcp_dynamic_cast > - (stepperDIRKObserver_); - } - } else { - this->stepperObserver_ = obs; - stepperDIRKObserver_ = - Teuchos::rcp_dynamic_cast >(this->stepperObserver_); - } + + if (this->stepperObserver_ == Teuchos::null) + this->stepperObserver_ = + Teuchos::rcp(new StepperRKObserverComposite()); + + if (( obs == Teuchos::null ) and (this->stepperObserver_->getSize() == 0) ) + obs = Teuchos::rcp(new StepperRKObserver()); + + this->stepperObserver_->addObserver( + Teuchos::rcp_dynamic_cast > (obs, true) ); + } @@ -122,9 +119,9 @@ void StepperDIRK::initialize() this->setObserver(); - TEUCHOS_TEST_FOR_EXCEPTION( this->stepperDIRKObserver_ == Teuchos::null, + TEUCHOS_TEST_FOR_EXCEPTION( this->stepperObserver_ == Teuchos::null, std::logic_error, - "Error - stepperDIRKObserver is null!\n"); + "Error - StepperRKObserver is null!\n"); // Initialize the stage vectors const int numStages = tableau_->numStages(); @@ -197,8 +194,11 @@ void StepperDIRK::takeStep( bool pass = true; Thyra::SolveStatus sStatus; for (int i=0; i < numStages; ++i) { - if (!Teuchos::is_null(stepperDIRKObserver_)) - stepperDIRKObserver_->observeBeginStage(solutionHistory, *this); + this->stepperObserver_->observeBeginStage(solutionHistory, *this); + + // ???: is it a good idea to leave this (no-op) here? + this->stepperObserver_ + ->observeBeforeImplicitExplicitly(solutionHistory, *this); if ( i == 0 && this->getUseFSAL() && workingState->getNConsecutiveFailures() == 0 ) { @@ -235,8 +235,7 @@ void StepperDIRK::takeStep( inArgs.set_x_dot(Teuchos::null); outArgs.set_f(stageXDot_[i]); - if (!Teuchos::is_null(stepperDIRKObserver_)) - stepperDIRKObserver_->observeBeforeExplicit(solutionHistory,*this); + this->stepperObserver_->observeBeforeExplicit(solutionHistory, *this); this->wrapperModel_->getAppModel()->evalModel(inArgs,outArgs); } } else { @@ -252,22 +251,19 @@ void StepperDIRK::takeStep( auto p = Teuchos::rcp(new ImplicitODEParameters( timeDer, dt, alpha, beta, SOLVE_FOR_X, i)); - if (!Teuchos::is_null(stepperDIRKObserver_)) - stepperDIRKObserver_->observeBeforeSolve(solutionHistory, *this); + this->stepperObserver_->observeBeforeSolve(solutionHistory, *this); sStatus = this->solveImplicitODE(stageX_, stageXDot_[i], ts, p); if (sStatus.solveStatus != Thyra::SOLVE_STATUS_CONVERGED) pass=false; - if (!Teuchos::is_null(stepperDIRKObserver_)) - stepperDIRKObserver_->observeAfterSolve(solutionHistory, *this); + this->stepperObserver_->observeAfterSolve(solutionHistory, *this); timeDer->compute(stageX_, stageXDot_[i]); } } - if (!Teuchos::is_null(stepperDIRKObserver_)) - stepperDIRKObserver_->observeEndStage(solutionHistory, *this); + this->stepperObserver_->observeEndStage(solutionHistory, *this); } // Sum for solution: x_n = x_n-1 + Sum{ dt*b(i) * f(i) } @@ -357,11 +353,6 @@ StepperDIRK::getValidParameters() const Teuchos::RCP pl = Teuchos::parameterList(); this->getValidParametersBasicDIRK(pl); - //this->getValidParametersBasic(pl); - //pl->set("Initial Condition Consistency Check", false); - //pl->set("Zero Initial Guess", false); - //pl->set("Reset Initial Guess", true); - return pl; } diff --git a/packages/tempus/src/Tempus_StepperExplicitRKObserver.hpp b/packages/tempus/src/Tempus_StepperExplicitRKObserver.hpp deleted file mode 100644 index c2be499e79cf..000000000000 --- a/packages/tempus/src/Tempus_StepperExplicitRKObserver.hpp +++ /dev/null @@ -1,73 +0,0 @@ -// @HEADER -// **************************************************************************** -// Tempus: Copyright (2017) Sandia Corporation -// -// Distributed under BSD 3-clause license (See accompanying file Copyright.txt) -// **************************************************************************** -// @HEADER - -#ifndef Tempus_StepperExplicitRKObserver_hpp -#define Tempus_StepperExplicitRKObserver_hpp - -#include "Tempus_SolutionHistory.hpp" -#include "Tempus_StepperObserver.hpp" - - -namespace Tempus { - -// Forward Declaration for recursive includes (this Observer <--> Stepper) -template class StepperExplicitRK; - -/** \brief StepperExplicitRKObserver class for StepperExplicitRK. - * - * This is a means for application developers to perform tasks - * during the time steps, e.g., - * - Compute specific quantities - * - Output information - * - "Massage" the working solution state - * - ... - * - * Design Considerations - * - StepperExplicitRKObserver is not stateless! Developers may touch the - * solution state! Developers need to be careful not to break the - * restart (checkpoint) capability. - */ -template -class StepperExplicitRKObserver - : virtual public Tempus::StepperObserver -{ -public: - - /// Constructor - StepperExplicitRKObserver(){} - - /// Destructor - virtual ~StepperExplicitRKObserver(){} - - /// Observe Stepper at beginning of takeStep. - virtual void observeBeginTakeStep( - Teuchos::RCP > /* sh */, - Stepper & /* stepper */){} - - /// Observe Stepper at beginning of each stage. - virtual void observeBeginStage( - Teuchos::RCP > /* sh */, - StepperExplicitRK & /* stepperExplicitRK */){} - - /// Observe Stepper before Explicit evaluation of Implicit ODE ME. - virtual void observeBeforeExplicit( - Teuchos::RCP > /* sh */, - StepperExplicitRK & /* stepperExplicitRK */){} - - /// Observe Stepper at end of each stage. - virtual void observeEndStage( - Teuchos::RCP > /* sh */, - StepperExplicitRK & /* stepperExplicitRK */){} - - /// Observe Stepper at end of takeStep. - virtual void observeEndTakeStep( - Teuchos::RCP > /* sh */, - Stepper & /* stepper */){} -}; -} // namespace Tempus -#endif // Tempus_StepperExplicitRKObserver_hpp diff --git a/packages/tempus/src/Tempus_StepperExplicitRKObserverComposite_decl.hpp b/packages/tempus/src/Tempus_StepperExplicitRKObserverComposite_decl.hpp deleted file mode 100644 index 9dc503665239..000000000000 --- a/packages/tempus/src/Tempus_StepperExplicitRKObserverComposite_decl.hpp +++ /dev/null @@ -1,79 +0,0 @@ -// @HEADER -// **************************************************************************** -// Tempus: Copyright (2017) Sandia Corporation -// -// Distributed under BSD 3-clause license (See accompanying file Copyright.txt) -// **************************************************************************** -// @HEADER - -#ifndef Tempus_StepperExplicitRKObserverComposite_decl_hpp -#define Tempus_StepperExplicitRKObserverComposite_decl_hpp - -#include "Tempus_StepperExplicitRKObserver.hpp" -#include "Tempus_TimeStepControl.hpp" -#include - -namespace Tempus { - -/** \brief This observer is a composite observer, - * - * which takes other StepperExplicitRKObservers and sequentially calls each - * individual observer function. - */ -template -class StepperExplicitRKObserverComposite - : virtual public Tempus::StepperExplicitRKObserver -{ -public: - - /// Default constructor - StepperExplicitRKObserverComposite(); - - /// Destructor - virtual ~StepperExplicitRKObserverComposite(); - - /// \name Override StepperExplicitRKObserver basic methods - //@{ - /// Observe Stepper at beginning of takeStep. - virtual void observeBeginTakeStep( - Teuchos::RCP > sh, - Stepper & stepper) override; - - /// Observe Stepper at beginning of each stage. - virtual void observeBeginStage( - Teuchos::RCP > sh, - StepperExplicitRK & stepperExplicitRK) override; - - /// Observe Stepper before Explicit evaluation of Implicit ODE ME. - virtual void observeBeforeExplicit( - Teuchos::RCP > sh, - StepperExplicitRK & stepperExplicitRK) override; - - /// Observe Stepper at end of each stage. - virtual void observeEndStage( - Teuchos::RCP > sh, - StepperExplicitRK & stepperExplicitRK) override; - - /// Observe Stepper at end of takeStep. - virtual void observeEndTakeStep( - Teuchos::RCP > sh, - Stepper & stepper) override; - - // add observer to the composite observer list - void addObserver(const Teuchos::RCP > &observer); - - // clear all observer from the composite observer list - void clearObservers(); - - // Size of composite observer list - bool empty(); - //@} - -private: - - std::vector > > observers_; - -}; - -} // namespace Tempus -#endif // Tempus_StepperExplicitRKObserverComposite_decl_hpp diff --git a/packages/tempus/src/Tempus_StepperExplicitRKObserverComposite_impl.hpp b/packages/tempus/src/Tempus_StepperExplicitRKObserverComposite_impl.hpp deleted file mode 100644 index 63cd8e7545b7..000000000000 --- a/packages/tempus/src/Tempus_StepperExplicitRKObserverComposite_impl.hpp +++ /dev/null @@ -1,84 +0,0 @@ -// @HEADER -// **************************************************************************** -// Tempus: Copyright (2017) Sandia Corporation -// -// Distributed under BSD 3-clause license (See accompanying file Copyright.txt) -// **************************************************************************** -// @HEADER - -#ifndef Tempus_StepperExplicitRKObserverComposite_impl_hpp -#define Tempus_StepperExplicitRKObserverComposite_impl_hpp - -#include "Tempus_StepperExplicitRKObserver.hpp" -#include "Tempus_TimeStepControl.hpp" - -namespace Tempus { - -template -StepperExplicitRKObserverComposite::StepperExplicitRKObserverComposite(){} - -template -StepperExplicitRKObserverComposite::~StepperExplicitRKObserverComposite(){} - -template -void StepperExplicitRKObserverComposite:: -observeBeginTakeStep(Teuchos::RCP > sh, - Stepper & stepper) -{ - for(auto& o : observers_) - o->observeBeginTakeStep(sh,stepper); -} - -template -void StepperExplicitRKObserverComposite::observeBeginStage( - Teuchos::RCP > sh, - StepperExplicitRK & stepperExplicitRK) -{ - for(auto& o : observers_) - o->observeBeginStage(sh,stepperExplicitRK); -} - -template -void StepperExplicitRKObserverComposite::observeBeforeExplicit( - Teuchos::RCP > sh, - StepperExplicitRK & stepperExplicitRK) -{ - for(auto& o : observers_) - o->observeBeforeExplicit(sh,stepperExplicitRK); -} - -template -void StepperExplicitRKObserverComposite::observeEndStage( - Teuchos::RCP > sh, - StepperExplicitRK & stepperExplicitRK) -{ - for(auto& o : observers_) - o->observeEndStage(sh,stepperExplicitRK); -} - -template -void StepperExplicitRKObserverComposite:: -observeEndTakeStep(Teuchos::RCP > sh, - Stepper & stepper) -{ - for(auto& o : observers_) - o->observeEndTakeStep(sh,stepper); -} - -template -void StepperExplicitRKObserverComposite:: -addObserver(const Teuchos::RCP > &observer) -{ - observers_.push_back(observer); -} - -template -void StepperExplicitRKObserverComposite:: -clearObservers() { observers_.clear();} - -template -bool StepperExplicitRKObserverComposite:: -empty() { return observers_.empty();} - -} // namespace Tempus -#endif // Tempus_StepperExplicitRKObserverComposite_impl_hpp diff --git a/packages/tempus/src/Tempus_StepperExplicitRK_decl.hpp b/packages/tempus/src/Tempus_StepperExplicitRK_decl.hpp index a8088a4bbd42..bb059417b60f 100644 --- a/packages/tempus/src/Tempus_StepperExplicitRK_decl.hpp +++ b/packages/tempus/src/Tempus_StepperExplicitRK_decl.hpp @@ -12,8 +12,7 @@ #include "Tempus_config.hpp" #include "Tempus_StepperExplicit.hpp" #include "Tempus_RKButcherTableau.hpp" -#include "Tempus_StepperObserverComposite.hpp" -#include "Tempus_StepperExplicitRKObserverComposite.hpp" +#include "Tempus_StepperRKObserverComposite.hpp" namespace Tempus { @@ -102,6 +101,9 @@ class StepperExplicitRK : virtual public Tempus::StepperExplicit virtual Teuchos::RCP > getTableau() { return tableau_; } + virtual Teuchos::RCP > getObserver() const + { return this->stepperObserver_; } + /// Initialize during construction and after changing input parameters. virtual void initialize(); @@ -161,7 +163,7 @@ class StepperExplicitRK : virtual public Tempus::StepperExplicit /// Setup for constructor. virtual void setup( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, @@ -172,11 +174,10 @@ class StepperExplicitRK : virtual public Tempus::StepperExplicit Teuchos::RCP > tableau_; - std::vector > > stageXDot_; - Teuchos::RCP > stageX_; + std::vector > > stageXDot_; + Teuchos::RCP > stageX_; - Teuchos::RCP > stepperObserver_; - Teuchos::RCP > stepperExplicitRKObserver_; + Teuchos::RCP > stepperObserver_; // For Embedded RK bool useEmbedded_; diff --git a/packages/tempus/src/Tempus_StepperExplicitRK_impl.hpp b/packages/tempus/src/Tempus_StepperExplicitRK_impl.hpp index 0593f40a6663..f51bcf997f85 100644 --- a/packages/tempus/src/Tempus_StepperExplicitRK_impl.hpp +++ b/packages/tempus/src/Tempus_StepperExplicitRK_impl.hpp @@ -26,16 +26,14 @@ void StepperExplicitRK::setupDefault() this->setUseEmbedded( this->getUseEmbeddedDefault()); this->stepperObserver_ = - Teuchos::rcp(new StepperObserverComposite()); - this->stepperExplicitRKObserver_ = - Teuchos::rcp(new StepperExplicitRKObserverComposite()); + Teuchos::rcp(new StepperRKObserverComposite()); } template void StepperExplicitRK::setup( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, @@ -47,9 +45,7 @@ void StepperExplicitRK::setup( this->setUseEmbedded( useEmbedded); this->stepperObserver_ = - Teuchos::rcp(new StepperObserverComposite()); - this->stepperExplicitRKObserver_ = - Teuchos::rcp(new StepperExplicitRKObserverComposite()); + Teuchos::rcp(new StepperRKObserverComposite()); this->setObserver(obs); if (appModel != Teuchos::null) { @@ -160,19 +156,17 @@ template void StepperExplicitRK::setObserver( Teuchos::RCP > obs) { - if (obs != Teuchos::null ) { - stepperObserver_->addObserver(obs); - auto ERKObs = - Teuchos::rcp_dynamic_cast > (obs); - if (ERKObs!=Teuchos::null) stepperExplicitRKObserver_->addObserver(ERKObs); - } else { - if (stepperExplicitRKObserver_->empty()) { - auto ERKObs = Teuchos::rcp(new StepperExplicitRKObserver()); - stepperExplicitRKObserver_->addObserver(ERKObs); - stepperObserver_->addObserver( - Teuchos::rcp_dynamic_cast > (ERKObs, true)); - } - } + + if (this->stepperObserver_ == Teuchos::null) + this->stepperObserver_ = + Teuchos::rcp(new StepperRKObserverComposite()); + + if (( obs == Teuchos::null ) and (this->stepperObserver_->getSize() == 0) ) + obs = Teuchos::rcp(new StepperRKObserver()); + + this->stepperObserver_->addObserver( + Teuchos::rcp_dynamic_cast > (obs, true) ); + } @@ -189,8 +183,8 @@ void StepperExplicitRK::initialize() this->setObserver(); - TEUCHOS_TEST_FOR_EXCEPTION( this->stepperObserver_->empty() || - this->stepperExplicitRKObserver_->empty(), std::logic_error, + TEUCHOS_TEST_FOR_EXCEPTION( this->stepperObserver_->getSize() < 1 + , std::logic_error, "Error - Composite Observer is empty!\n"); // Initialize the stage vectors @@ -256,8 +250,14 @@ void StepperExplicitRK::takeStep( // Compute stage solutions for (int i=0; i < numStages; ++i) { - if (!Teuchos::is_null(stepperExplicitRKObserver_)) - stepperExplicitRKObserver_->observeBeginStage(solutionHistory, *this); + this->stepperObserver_->observeBeginStage(solutionHistory, *this); + + // ???: is it a good idea to leave this (no-op) here? + this->stepperObserver_ + ->observeBeforeImplicitExplicitly(solutionHistory, *this); + + // ???: is it a good idea to leave this (no-op) here? + this->stepperObserver_->observeBeforeSolve(solutionHistory, *this); if ( i == 0 && this->getUseFSAL() && workingState->getNConsecutiveFailures() == 0 ) { @@ -276,18 +276,17 @@ void StepperExplicitRK::takeStep( } const Scalar ts = time + c(i)*dt; - if (!Teuchos::is_null(stepperExplicitRKObserver_)) - stepperExplicitRKObserver_->observeBeforeExplicit(solutionHistory, - *this); + // ???: is it a good idea to leave this (no-op) here? + this->stepperObserver_->observeAfterSolve(solutionHistory, *this); + this->stepperObserver_->observeBeforeExplicit(solutionHistory, *this); auto p = Teuchos::rcp(new ExplicitODEParameters(dt)); // Evaluate xDot = f(x,t). this->evaluateExplicitODE(stageXDot_[i], stageX_, ts, p); } - if (!Teuchos::is_null(stepperExplicitRKObserver_)) - stepperExplicitRKObserver_->observeEndStage(solutionHistory, *this); + this->stepperObserver_->observeEndStage(solutionHistory, *this); } // Sum for solution: x_n = x_n-1 + Sum{ b(i) * dt*f(i) } diff --git a/packages/tempus/src/Tempus_StepperForwardEuler_decl.hpp b/packages/tempus/src/Tempus_StepperForwardEuler_decl.hpp index a9ce32a6eda5..f172c479294e 100644 --- a/packages/tempus/src/Tempus_StepperForwardEuler_decl.hpp +++ b/packages/tempus/src/Tempus_StepperForwardEuler_decl.hpp @@ -74,6 +74,9 @@ class StepperForwardEuler : virtual public Tempus::StepperExplicit virtual void setObserver( Teuchos::RCP > obs = Teuchos::null); + virtual Teuchos::RCP > getObserver() const + { return this->stepperFEObserver_; } + /// Initialize during construction and after changing input parameters. virtual void initialize(); diff --git a/packages/tempus/src/Tempus_StepperHHTAlpha_decl.hpp b/packages/tempus/src/Tempus_StepperHHTAlpha_decl.hpp index 6dd21a9be78b..3336641852cf 100644 --- a/packages/tempus/src/Tempus_StepperHHTAlpha_decl.hpp +++ b/packages/tempus/src/Tempus_StepperHHTAlpha_decl.hpp @@ -75,6 +75,9 @@ class StepperHHTAlpha : virtual public Tempus::StepperImplicit virtual void setObserver( Teuchos::RCP > /* obs */ = Teuchos::null){} + virtual Teuchos::RCP > getObserver() const + { return Teuchos::null; } + /// Initialize during construction and after changing input parameters. virtual void initialize(); diff --git a/packages/tempus/src/Tempus_StepperIMEX_RKPartObserver.hpp b/packages/tempus/src/Tempus_StepperIMEX_RKPartObserver.hpp deleted file mode 100644 index 428094a5f15c..000000000000 --- a/packages/tempus/src/Tempus_StepperIMEX_RKPartObserver.hpp +++ /dev/null @@ -1,87 +0,0 @@ -// @HEADER -// **************************************************************************** -// Tempus: Copyright (2017) Sandia Corporation -// -// Distributed under BSD 3-clause license (See accompanying file Copyright.txt) -// **************************************************************************** -// @HEADER - -#ifndef Tempus_StepperIMEX_RKPartObserver_hpp -#define Tempus_StepperIMEX_RKPartObserver_hpp - -#include "Tempus_SolutionHistory.hpp" - - -namespace Tempus { - -// Forward Declaration for recursive includes (this Observer <--> Stepper) -template class StepperIMEX_RK_Partition; - -/** \brief StepperIMEX_RKPartObserver class for StepperIMEX_RK_Partition. - * - * This is a means for application developers to perform tasks - * during the time steps, e.g., - * - Compute specific quantities - * - Output information - * - "Massage" the working solution state - * - ... - * - * Design Considerations - * - StepperIMEX_RKPartObserver is not stateless! Developers may touch the - * solution state! Developers need to be careful not to break the - * restart (checkpoint) capability. - */ -template -class StepperIMEX_RKPartObserver - : virtual public Tempus::StepperObserver -{ -public: - - /// Constructor - StepperIMEX_RKPartObserver(){} - - /// Destructor - virtual ~StepperIMEX_RKPartObserver(){} - - /// Observe Stepper at beginning of takeStep. - virtual void observeBeginTakeStep( - Teuchos::RCP > /* sh */, - Stepper & /* stepper */){} - - /// Observe Stepper at beginning of each stage. - virtual void observeBeginStage( - Teuchos::RCP > /* sh */, - StepperIMEX_RK_Partition & /* stepperIMEX_RK_Part */){} - - /// Observe Stepper before Explicit evaluation of Implicit ODE ME. - virtual void observeBeforeImplicitExplicitly( - Teuchos::RCP > /* sh */, - StepperIMEX_RK_Partition & /* stepperIMEX_RK_Part */){} - - /// Observe Stepper before nonlinear solve. - virtual void observeBeforeSolve( - Teuchos::RCP > /* sh */, - StepperIMEX_RK_Partition & /* stepperIMEX_RK_Part */){} - - /// Observe Stepper after nonlinear solve. - virtual void observeAfterSolve( - Teuchos::RCP > /* sh */, - StepperIMEX_RK_Partition & /* stepperIMEX_RK_Part */){} - - /// Observe Stepper before Explicit ME evaluation. - virtual void observeBeforeExplicit( - Teuchos::RCP > /* sh */, - StepperIMEX_RK_Partition & /* stepperIMEX_RK_Part */){} - - /// Observe Stepper at end of each stage. - virtual void observeEndStage( - Teuchos::RCP > /* sh */, - StepperIMEX_RK_Partition & /* stepperIMEX_RK_Part */){} - - /// Observe Stepper at end of takeStep. - virtual void observeEndTakeStep( - Teuchos::RCP > /* sh */, - Stepper & /* stepper */){} -}; -} // namespace Tempus -#endif // Tempus_StepperIMEX_RKPartObserver_hpp diff --git a/packages/tempus/src/Tempus_StepperIMEX_RK_Partition_decl.hpp b/packages/tempus/src/Tempus_StepperIMEX_RK_Partition_decl.hpp index 3532dc9969b3..a90fcf656467 100644 --- a/packages/tempus/src/Tempus_StepperIMEX_RK_Partition_decl.hpp +++ b/packages/tempus/src/Tempus_StepperIMEX_RK_Partition_decl.hpp @@ -13,7 +13,7 @@ #include "Tempus_RKButcherTableau.hpp" #include "Tempus_StepperImplicit.hpp" #include "Tempus_WrapperModelEvaluatorPairPartIMEX_Basic.hpp" -#include "Tempus_StepperIMEX_RKPartObserver.hpp" +#include "Tempus_StepperRKObserverComposite.hpp" namespace Tempus { @@ -281,6 +281,9 @@ class StepperIMEX_RK_Partition : virtual public Tempus::StepperImplicit virtual void setObserver( Teuchos::RCP > obs = Teuchos::null); + virtual Teuchos::RCP > getObserver() const + { return this->stepperObserver_; } + /// Initialize during construction and after changing input parameters. virtual void initialize(); @@ -350,7 +353,7 @@ class StepperIMEX_RK_Partition : virtual public Tempus::StepperImplicit Teuchos::RCP > xTilde_; - Teuchos::RCP > stepperIMEX_RKPartObserver_; + Teuchos::RCP > stepperObserver_; }; diff --git a/packages/tempus/src/Tempus_StepperIMEX_RK_Partition_impl.hpp b/packages/tempus/src/Tempus_StepperIMEX_RK_Partition_impl.hpp index 87cca502afb4..3745d54af212 100644 --- a/packages/tempus/src/Tempus_StepperIMEX_RK_Partition_impl.hpp +++ b/packages/tempus/src/Tempus_StepperIMEX_RK_Partition_impl.hpp @@ -328,21 +328,15 @@ template void StepperIMEX_RK_Partition::setObserver( Teuchos::RCP > obs) { - if (obs == Teuchos::null) { - // Create default observer, otherwise keep current observer. - if (this->stepperObserver_ == Teuchos::null) { - stepperIMEX_RKPartObserver_ = - Teuchos::rcp(new StepperIMEX_RKPartObserver()); - this->stepperObserver_ = - Teuchos::rcp_dynamic_cast > - (stepperIMEX_RKPartObserver_); - } - } else { - this->stepperObserver_ = obs; - stepperIMEX_RKPartObserver_ = - Teuchos::rcp_dynamic_cast > - (this->stepperObserver_); - } + if (this->stepperObserver_ == Teuchos::null) + this->stepperObserver_ = + Teuchos::rcp(new StepperRKObserverComposite()); + + if (( obs == Teuchos::null ) and (this->stepperObserver_->getSize() == 0) ) + obs = Teuchos::rcp(new StepperRKObserver()); + + this->stepperObserver_->addObserver( + Teuchos::rcp_dynamic_cast > (obs, true) ); } @@ -550,8 +544,7 @@ void StepperIMEX_RK_Partition::takeStep( // Compute stage solutions for (int i = 0; i < numStages; ++i) { - if (!Teuchos::is_null(stepperIMEX_RKPartObserver_)) - stepperIMEX_RKPartObserver_->observeBeginStage(solutionHistory, *this); + this->stepperObserver_->observeBeginStage(solutionHistory, *this); Thyra::assign(stageY.ptr(), *(wrapperModelPairIMEX->getExplicitOnlyVector(currentState->getX()))); @@ -582,9 +575,7 @@ void StepperIMEX_RK_Partition::takeStep( assign(stageGx_[i].ptr(), Teuchos::ScalarTraits::zero()); } else { Thyra::assign(stageX.ptr(), *xTilde_); - if (!Teuchos::is_null(stepperIMEX_RKPartObserver_)) - stepperIMEX_RKPartObserver_-> - observeBeforeImplicitExplicitly(solutionHistory, *this); + this->stepperObserver_->observeBeforeImplicitExplicitly(solutionHistory, *this); evalImplicitModelExplicitly(stageX, stageY, ts, dt, i, stageGx_[i]); } } else { @@ -617,8 +608,7 @@ void StepperIMEX_RK_Partition::takeStep( wrapperModelPairIMEX->setForSolve(timeDer, inArgs, outArgs); - if (!Teuchos::is_null(stepperIMEX_RKPartObserver_)) - stepperIMEX_RKPartObserver_->observeBeforeSolve(solutionHistory, *this); + this->stepperObserver_->observeBeforeSolve(solutionHistory, *this); this->solver_->setModel(wrapperModelPairIMEX); sStatus = this->solveImplicitODE(stageX); @@ -626,18 +616,15 @@ void StepperIMEX_RK_Partition::takeStep( wrapperModelPairIMEX->setUseImplicitModel(false); - if (!Teuchos::is_null(stepperIMEX_RKPartObserver_)) - stepperIMEX_RKPartObserver_->observeAfterSolve(solutionHistory, *this); + this->stepperObserver_->observeAfterSolve(solutionHistory, *this); // Update contributions to stage values Thyra::V_StVpStV(stageGx_[i].ptr(), -alpha, *stageX, alpha, *xTilde_); } - if (!Teuchos::is_null(stepperIMEX_RKPartObserver_)) - stepperIMEX_RKPartObserver_->observeBeforeExplicit(solutionHistory,*this); + this->stepperObserver_->observeBeforeExplicit(solutionHistory, *this); evalExplicitModel(stageZ_, tHats, dt, i, stageF_[i]); - if (!Teuchos::is_null(stepperIMEX_RKPartObserver_)) - stepperIMEX_RKPartObserver_->observeEndStage(solutionHistory, *this); + this->stepperObserver_->observeEndStage(solutionHistory, *this); } // Sum for solution: y_n = y_n-1 - dt*Sum{ bHat(i)*fy(i) } diff --git a/packages/tempus/src/Tempus_StepperIMEX_RK_decl.hpp b/packages/tempus/src/Tempus_StepperIMEX_RK_decl.hpp index 3c0cc0086e4b..51a0ba79000e 100644 --- a/packages/tempus/src/Tempus_StepperIMEX_RK_decl.hpp +++ b/packages/tempus/src/Tempus_StepperIMEX_RK_decl.hpp @@ -13,7 +13,7 @@ #include "Tempus_RKButcherTableau.hpp" #include "Tempus_StepperImplicit.hpp" #include "Tempus_WrapperModelEvaluatorPairIMEX_Basic.hpp" -#include "Tempus_StepperIMEX_RKObserver.hpp" +#include "Tempus_StepperRKObserverComposite.hpp" namespace Tempus { @@ -281,6 +281,9 @@ class StepperIMEX_RK : virtual public Tempus::StepperImplicit virtual void setObserver( Teuchos::RCP > obs = Teuchos::null); + virtual Teuchos::RCP > getObserver() const + { return this->stepperObserver_; } + /// Initialize during construction and after changing input parameters. virtual void initialize(); @@ -351,7 +354,7 @@ class StepperIMEX_RK : virtual public Tempus::StepperImplicit Teuchos::RCP > xTilde_; - Teuchos::RCP > stepperIMEX_RKObserver_; + Teuchos::RCP > stepperObserver_; }; diff --git a/packages/tempus/src/Tempus_StepperIMEX_RK_impl.hpp b/packages/tempus/src/Tempus_StepperIMEX_RK_impl.hpp index 0bdf25e6d29d..0e86936524ec 100644 --- a/packages/tempus/src/Tempus_StepperIMEX_RK_impl.hpp +++ b/packages/tempus/src/Tempus_StepperIMEX_RK_impl.hpp @@ -332,21 +332,17 @@ template void StepperIMEX_RK::setObserver( Teuchos::RCP > obs) { - if (obs == Teuchos::null) { - // Create default observer, otherwise keep current observer. - if (this->stepperObserver_ == Teuchos::null) { - stepperIMEX_RKObserver_ = - Teuchos::rcp(new StepperIMEX_RKObserver()); - this->stepperObserver_ = - Teuchos::rcp_dynamic_cast > - (stepperIMEX_RKObserver_); - } - } else { - this->stepperObserver_ = obs; - stepperIMEX_RKObserver_ = - Teuchos::rcp_dynamic_cast > - (this->stepperObserver_); - } + + if (this->stepperObserver_ == Teuchos::null) + this->stepperObserver_ = + Teuchos::rcp(new StepperRKObserverComposite()); + + if (( obs == Teuchos::null ) and (this->stepperObserver_->getSize() == 0) ) + obs = Teuchos::rcp(new StepperRKObserver()); + + this->stepperObserver_->addObserver( + Teuchos::rcp_dynamic_cast > (obs, true) ); + } @@ -532,8 +528,7 @@ void StepperIMEX_RK::takeStep( // Compute stage solutions for (int i = 0; i < numStages; ++i) { - if (!Teuchos::is_null(stepperIMEX_RKObserver_)) - stepperIMEX_RKObserver_->observeBeginStage(solutionHistory, *this); + this->stepperObserver_->observeBeginStage(solutionHistory, *this); Thyra::assign(xTilde_.ptr(), *(currentState->getX())); for (int j = 0; j < i; ++j) { if (AHat(i,j) != Teuchos::ScalarTraits::zero()) @@ -554,9 +549,7 @@ void StepperIMEX_RK::takeStep( assign(stageG_[i].ptr(), Teuchos::ScalarTraits::zero()); } else { Thyra::assign(stageX_.ptr(), *xTilde_); - if (!Teuchos::is_null(stepperIMEX_RKObserver_)) - stepperIMEX_RKObserver_-> - observeBeforeImplicitExplicitly(solutionHistory, *this); + this->stepperObserver_->observeBeforeImplicitExplicitly(solutionHistory, *this); evalImplicitModelExplicitly(stageX_, ts, dt, i, stageG_[i]); } } else { @@ -572,26 +565,22 @@ void StepperIMEX_RK::takeStep( auto p = Teuchos::rcp(new ImplicitODEParameters( timeDer, dt, alpha, beta, SOLVE_FOR_X, i)); - if (!Teuchos::is_null(stepperIMEX_RKObserver_)) - stepperIMEX_RKObserver_->observeBeforeSolve(solutionHistory, *this); + this->stepperObserver_->observeBeforeSolve(solutionHistory, *this); const Thyra::SolveStatus sStatus = this->solveImplicitODE(stageX_, stageG_[i], ts, p); if (sStatus.solveStatus != Thyra::SOLVE_STATUS_CONVERGED) pass = false; - if (!Teuchos::is_null(stepperIMEX_RKObserver_)) - stepperIMEX_RKObserver_->observeAfterSolve(solutionHistory, *this); + this->stepperObserver_->observeAfterSolve(solutionHistory, *this); // Update contributions to stage values Thyra::V_StVpStV(stageG_[i].ptr(), -alpha, *stageX_, alpha, *xTilde_); } - if (!Teuchos::is_null(stepperIMEX_RKObserver_)) - stepperIMEX_RKObserver_->observeBeforeExplicit(solutionHistory, *this); + this->stepperObserver_->observeBeforeExplicit(solutionHistory, *this); evalExplicitModel(stageX_, tHats, dt, i, stageF_[i]); - if (!Teuchos::is_null(stepperIMEX_RKObserver_)) - stepperIMEX_RKObserver_->observeEndStage(solutionHistory, *this); + this->stepperObserver_->observeEndStage(solutionHistory, *this); } // Sum for solution: x_n = x_n-1 - dt*Sum{ bHat(i)*f(i) + b(i)*g(i) } diff --git a/packages/tempus/src/Tempus_StepperLeapfrog_decl.hpp b/packages/tempus/src/Tempus_StepperLeapfrog_decl.hpp index 7131c4e1aa9e..1109b2251215 100644 --- a/packages/tempus/src/Tempus_StepperLeapfrog_decl.hpp +++ b/packages/tempus/src/Tempus_StepperLeapfrog_decl.hpp @@ -11,6 +11,7 @@ #include "Tempus_config.hpp" #include "Tempus_StepperExplicit.hpp" +#include "Tempus_StepperObserverComposite.hpp" #include "Tempus_StepperLeapfrogObserver.hpp" @@ -97,6 +98,9 @@ class StepperLeapfrog : virtual public Tempus::StepperExplicit virtual void setObserver( Teuchos::RCP > obs = Teuchos::null); + virtual Teuchos::RCP > getObserver() const + { return this->stepperObserver_; } + /// Initialize during construction and after changing input parameters. virtual void initialize(); @@ -139,6 +143,7 @@ class StepperLeapfrog : virtual public Tempus::StepperExplicit protected: + Teuchos::RCP > stepperObserver_; Teuchos::RCP > stepperLFObserver_; }; diff --git a/packages/tempus/src/Tempus_StepperLeapfrog_impl.hpp b/packages/tempus/src/Tempus_StepperLeapfrog_impl.hpp index 5c1cfc309dfd..ae3c4be1aad6 100644 --- a/packages/tempus/src/Tempus_StepperLeapfrog_impl.hpp +++ b/packages/tempus/src/Tempus_StepperLeapfrog_impl.hpp @@ -55,20 +55,17 @@ template void StepperLeapfrog::setObserver( Teuchos::RCP > obs) { - if (obs == Teuchos::null) { - // Create default observer, otherwise keep current observer. - if (this->stepperObserver_ == Teuchos::null) { - stepperLFObserver_ = - Teuchos::rcp(new StepperLeapfrogObserver()); - this->stepperObserver_ = - Teuchos::rcp_dynamic_cast >(stepperLFObserver_,true); - } - } else { - this->stepperObserver_ = obs; - stepperLFObserver_ = - Teuchos::rcp_dynamic_cast > - (this->stepperObserver_,true); - } + + if (this->stepperObserver_ == Teuchos::null) + this->stepperObserver_ = + Teuchos::rcp(new StepperObserverComposite()); + + if (( obs == Teuchos::null ) and (this->stepperObserver_->getSize() == 0) ) + obs = Teuchos::rcp(new StepperLeapfrogObserver()); + + this->stepperObserver_->addObserver( + Teuchos::rcp_dynamic_cast > (obs, true) ); + } template @@ -123,7 +120,7 @@ void StepperLeapfrog::takeStep( "Try setting in \"Solution History\" \"Storage Type\" = \"Undo\"\n" " or \"Storage Type\" = \"Static\" and \"Storage Limit\" = \"2\"\n"); - this->stepperObserver_->observeBeginTakeStep(solutionHistory, *this); + //this->stepperObserver_->observeBeginTakeStep(solutionHistory, *this); RCP > currentState=solutionHistory->getCurrentState(); RCP > workingState=solutionHistory->getWorkingState(); const Scalar time = currentState->getTime(); @@ -172,7 +169,7 @@ void StepperLeapfrog::takeStep( workingState->setSolutionStatus(Status::PASSED); workingState->setOrder(this->getOrder()); - this->stepperObserver_->observeEndTakeStep(solutionHistory, *this); + //this->stepperObserver_->observeEndTakeStep(solutionHistory, *this); } return; } diff --git a/packages/tempus/src/Tempus_StepperNewmarkExplicitAForm_decl.hpp b/packages/tempus/src/Tempus_StepperNewmarkExplicitAForm_decl.hpp index 58e7d8e805c8..c978ae5a2d58 100644 --- a/packages/tempus/src/Tempus_StepperNewmarkExplicitAForm_decl.hpp +++ b/packages/tempus/src/Tempus_StepperNewmarkExplicitAForm_decl.hpp @@ -81,6 +81,9 @@ class StepperNewmarkExplicitAForm virtual void setObserver( Teuchos::RCP > /* obs */ = Teuchos::null){} + virtual Teuchos::RCP > getObserver() const + { return Teuchos::null; } + /// Initialize during construction and after changing input parameters. virtual void initialize(); diff --git a/packages/tempus/src/Tempus_StepperNewmarkImplicitAForm_decl.hpp b/packages/tempus/src/Tempus_StepperNewmarkImplicitAForm_decl.hpp index b998079f4686..c2f2ed4e4d70 100644 --- a/packages/tempus/src/Tempus_StepperNewmarkImplicitAForm_decl.hpp +++ b/packages/tempus/src/Tempus_StepperNewmarkImplicitAForm_decl.hpp @@ -103,6 +103,9 @@ class StepperNewmarkImplicitAForm virtual void setObserver( Teuchos::RCP > /* obs */ = Teuchos::null){} + virtual Teuchos::RCP > getObserver() const + { return Teuchos::null; } + /// Initialize during construction and after changing input parameters. virtual void initialize(); diff --git a/packages/tempus/src/Tempus_StepperNewmarkImplicitDForm_decl.hpp b/packages/tempus/src/Tempus_StepperNewmarkImplicitDForm_decl.hpp index 2335633a2047..942eaada87c8 100644 --- a/packages/tempus/src/Tempus_StepperNewmarkImplicitDForm_decl.hpp +++ b/packages/tempus/src/Tempus_StepperNewmarkImplicitDForm_decl.hpp @@ -68,6 +68,9 @@ class StepperNewmarkImplicitDForm : virtual public Tempus::StepperImplicit > /* obs */ = Teuchos::null){} + virtual Teuchos::RCP > getObserver() const + { return Teuchos::null; } + /// Initialize during construction and after changing input parameters. virtual void initialize(); diff --git a/packages/tempus/src/Tempus_StepperObserverComposite_decl.hpp b/packages/tempus/src/Tempus_StepperObserverComposite_decl.hpp index 177783086fd9..f66c22548032 100644 --- a/packages/tempus/src/Tempus_StepperObserverComposite_decl.hpp +++ b/packages/tempus/src/Tempus_StepperObserverComposite_decl.hpp @@ -46,8 +46,8 @@ class StepperObserverComposite // clear all observer from the composite observer list void clearObservers(); - // Size of composite observer list - bool empty(); + // get the number of RK stepper observers present in the composite + std::size_t getSize() const { return observers_.size(); } //@} private: diff --git a/packages/tempus/src/Tempus_StepperObserverComposite_impl.hpp b/packages/tempus/src/Tempus_StepperObserverComposite_impl.hpp index f385112cfb55..0e69b6b844e7 100644 --- a/packages/tempus/src/Tempus_StepperObserverComposite_impl.hpp +++ b/packages/tempus/src/Tempus_StepperObserverComposite_impl.hpp @@ -49,9 +49,5 @@ template void StepperObserverComposite:: clearObservers() { observers_.clear();} -template -bool StepperObserverComposite:: -empty() { return observers_.empty();} - } // namespace Tempus #endif // Tempus_StepperObserverComposite_impl_hpp diff --git a/packages/tempus/src/Tempus_StepperOperatorSplit_decl.hpp b/packages/tempus/src/Tempus_StepperOperatorSplit_decl.hpp index 58b4622a9ff7..06d539098df8 100644 --- a/packages/tempus/src/Tempus_StepperOperatorSplit_decl.hpp +++ b/packages/tempus/src/Tempus_StepperOperatorSplit_decl.hpp @@ -63,17 +63,25 @@ class StepperOperatorSplit : virtual public Tempus::Stepper //@{ virtual void setModel( const Teuchos::RCP >& appModel); + virtual void setNonConstModel( const Teuchos::RCP >& appModel); + virtual Teuchos::RCP > getModel(); virtual void setSolver( Teuchos::RCP > solver); + virtual Teuchos::RCP > getSolver() const { return Teuchos::null; } + virtual void setObserver( Teuchos::RCP > obs = Teuchos::null); + + virtual Teuchos::RCP > getObserver() const + { return this->stepperOSObserver_; } + virtual void setTempState(Teuchos::RCP> state) { tempState_ = state; } diff --git a/packages/tempus/src/Tempus_StepperOperatorSplit_impl.hpp b/packages/tempus/src/Tempus_StepperOperatorSplit_impl.hpp index a641e99ad65a..d904488c7f66 100644 --- a/packages/tempus/src/Tempus_StepperOperatorSplit_impl.hpp +++ b/packages/tempus/src/Tempus_StepperOperatorSplit_impl.hpp @@ -134,7 +134,7 @@ void StepperOperatorSplit::setObserver( } } else { stepperOSObserver_ = - Teuchos::rcp_dynamic_cast > (obs); + Teuchos::rcp_dynamic_cast > (obs, true); } } diff --git a/packages/tempus/src/Tempus_StepperRKButcherTableau.hpp b/packages/tempus/src/Tempus_StepperRKButcherTableau.hpp index e6891b6868be..f330cf137e93 100644 --- a/packages/tempus/src/Tempus_StepperRKButcherTableau.hpp +++ b/packages/tempus/src/Tempus_StepperRKButcherTableau.hpp @@ -50,7 +50,7 @@ class StepperERK_ForwardEuler : StepperERK_ForwardEuler( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, @@ -122,7 +122,7 @@ class StepperERK_4Stage4thOrder : StepperERK_4Stage4thOrder( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, @@ -222,7 +222,7 @@ class StepperERK_BogackiShampine32 : StepperERK_BogackiShampine32( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, @@ -335,7 +335,7 @@ class StepperERK_Merson45 : StepperERK_Merson45( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, @@ -456,7 +456,7 @@ class StepperERK_3_8Rule : StepperERK_3_8Rule( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, @@ -559,7 +559,7 @@ class StepperERK_4Stage3rdOrderRunge : StepperERK_4Stage3rdOrderRunge( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, @@ -657,7 +657,7 @@ class StepperERK_5Stage3rdOrderKandG : StepperERK_5Stage3rdOrderKandG( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, @@ -755,7 +755,7 @@ class StepperERK_3Stage3rdOrder : StepperERK_3Stage3rdOrder( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, @@ -856,7 +856,7 @@ class StepperERK_3Stage3rdOrderTVD : StepperERK_3Stage3rdOrderTVD( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, @@ -958,7 +958,7 @@ class StepperERK_3Stage3rdOrderHeun : StepperERK_3Stage3rdOrderHeun( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, @@ -1055,7 +1055,7 @@ class StepperERK_Midpoint : StepperERK_Midpoint( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, @@ -1143,7 +1143,7 @@ class StepperERK_Trapezoidal : StepperERK_Trapezoidal( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, @@ -1241,7 +1241,7 @@ class StepperERK_General : StepperERK_General( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, bool useFSAL, std::string ICConsistency, bool ICConsistencyCheck, @@ -1369,7 +1369,7 @@ class StepperDIRK_BackwardEuler : StepperDIRK_BackwardEuler( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, const Teuchos::RCP >& solver, bool useFSAL, std::string ICConsistency, @@ -1474,7 +1474,7 @@ class StepperSDIRK_2Stage2ndOrder : StepperSDIRK_2Stage2ndOrder( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, const Teuchos::RCP >& solver, bool useFSAL, std::string ICConsistency, @@ -1614,7 +1614,7 @@ class StepperSDIRK_2Stage3rdOrder : StepperSDIRK_2Stage3rdOrder( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, const Teuchos::RCP >& solver, bool useFSAL, std::string ICConsistency, @@ -1775,7 +1775,7 @@ class StepperEDIRK_2Stage3rdOrder : StepperEDIRK_2Stage3rdOrder( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, const Teuchos::RCP >& solver, bool useFSAL, std::string ICConsistency, @@ -1891,7 +1891,7 @@ class StepperDIRK_1StageTheta : StepperDIRK_1StageTheta( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, const Teuchos::RCP >& solver, bool useFSAL, std::string ICConsistency, @@ -2024,7 +2024,7 @@ class StepperEDIRK_2StageTheta : StepperEDIRK_2StageTheta( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, const Teuchos::RCP >& solver, bool useFSAL, std::string ICConsistency, @@ -2158,7 +2158,7 @@ class StepperEDIRK_TrapezoidalRule : StepperEDIRK_TrapezoidalRule( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, const Teuchos::RCP >& solver, bool useFSAL, std::string ICConsistency, @@ -2269,7 +2269,7 @@ class StepperSDIRK_ImplicitMidpoint : StepperSDIRK_ImplicitMidpoint( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, const Teuchos::RCP >& solver, bool useFSAL, std::string ICConsistency, @@ -2377,7 +2377,7 @@ class StepperDIRK_1Stage1stOrderRadauIA : StepperDIRK_1Stage1stOrderRadauIA( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, const Teuchos::RCP >& solver, bool useFSAL, std::string ICConsistency, @@ -2477,7 +2477,7 @@ class StepperDIRK_2Stage2ndOrderLobattoIIIB : StepperDIRK_2Stage2ndOrderLobattoIIIB( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, const Teuchos::RCP >& solver, bool useFSAL, std::string ICConsistency, @@ -2594,7 +2594,7 @@ class StepperSDIRK_5Stage4thOrder : StepperSDIRK_5Stage4thOrder( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, const Teuchos::RCP >& solver, bool useFSAL, std::string ICConsistency, @@ -2754,7 +2754,7 @@ class StepperSDIRK_3Stage4thOrder : StepperSDIRK_3Stage4thOrder( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, const Teuchos::RCP >& solver, bool useFSAL, std::string ICConsistency, @@ -2885,7 +2885,7 @@ class StepperSDIRK_5Stage5thOrder : StepperSDIRK_5Stage5thOrder( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, const Teuchos::RCP >& solver, bool useFSAL, std::string ICConsistency, @@ -3058,7 +3058,7 @@ class StepperSDIRK_21Pair : StepperSDIRK_21Pair( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, const Teuchos::RCP >& solver, bool useFSAL, std::string ICConsistency, @@ -3179,7 +3179,7 @@ class StepperDIRK_General : StepperDIRK_General( const Teuchos::RCP >& appModel, - const Teuchos::RCP >& obs, + const Teuchos::RCP >& obs, const Teuchos::RCP >& solver, bool useFSAL, std::string ICConsistency, diff --git a/packages/tempus/src/Tempus_StepperIMEX_RKObserver.hpp b/packages/tempus/src/Tempus_StepperRKObserver.hpp similarity index 61% rename from packages/tempus/src/Tempus_StepperIMEX_RKObserver.hpp rename to packages/tempus/src/Tempus_StepperRKObserver.hpp index 276d0834330b..6924f137de04 100644 --- a/packages/tempus/src/Tempus_StepperIMEX_RKObserver.hpp +++ b/packages/tempus/src/Tempus_StepperRKObserver.hpp @@ -6,18 +6,19 @@ // **************************************************************************** // @HEADER -#ifndef Tempus_StepperIMEX_RKObserver_hpp -#define Tempus_StepperIMEX_RKObserver_hpp +#ifndef Tempus_StepperRKObserver_hpp +#define Tempus_StepperRKObserver_hpp #include "Tempus_SolutionHistory.hpp" +#include "Tempus_StepperObserver.hpp" namespace Tempus { // Forward Declaration for recursive includes (this Observer <--> Stepper) -template class StepperIMEX_RK; +//template class Stepper; -/** \brief StepperIMEX_RKObserver class for StepperIMEX_RK. +/** \brief StepperRKObserver class for StepperRK. * * This is a means for application developers to perform tasks * during the time steps, e.g., @@ -27,61 +28,61 @@ template class StepperIMEX_RK; * - ... * * Design Considerations - * - StepperIMEX_RKObserver is not stateless! Developers may touch the + * - StepperRKObserver is not stateless! Developers may touch the * solution state! Developers need to be careful not to break the * restart (checkpoint) capability. */ template -class StepperIMEX_RKObserver +class StepperRKObserver : virtual public Tempus::StepperObserver { public: /// Constructor - StepperIMEX_RKObserver(){} + StepperRKObserver(){} /// Destructor - virtual ~StepperIMEX_RKObserver(){} + virtual ~StepperRKObserver(){} - /// Observe Stepper at beginning of takeStep. + /// 1.) Observe Stepper at beginning of takeStep. virtual void observeBeginTakeStep( Teuchos::RCP > /* sh */, Stepper & /* stepper */){} - /// Observe Stepper at beginning of each stage. + /// 2.) Observe Stepper at beginning of each stage. virtual void observeBeginStage( Teuchos::RCP > /* sh */, - StepperIMEX_RK & /* stepperIMEX_RK */){} + Stepper & /* stepper */){} - /// Observe Stepper before Explicit evaluation of Implicit ODE ME. + /// 3.) Observe Stepper before Explicit evaluation of Implicit ODE ME (IMEX). virtual void observeBeforeImplicitExplicitly( Teuchos::RCP > /* sh */, - StepperIMEX_RK & /* stepperIMEX_RK */){} + Stepper & /* stepper */){} - /// Observe Stepper before nonlinear solve. + /// 4.) Observe Stepper before nonlinear solve (DIRK/IMEX). virtual void observeBeforeSolve( Teuchos::RCP > /* sh */, - StepperIMEX_RK & /* stepperIMEX_RK */){} + Stepper & /* stepper */){} - /// Observe Stepper after nonlinear solve. + /// 5.) Observe Stepper after nonlinear solve (DIRK/IMEX). virtual void observeAfterSolve( Teuchos::RCP > /* sh */, - StepperIMEX_RK & /* stepperIMEX_RK */){} + Stepper & /* stepper */){} - /// Observe Stepper before Explicit ME evaluation. + /// 6.) Observe Stepper before Explicit evaluation of Implicit ODE ME (IMEX). virtual void observeBeforeExplicit( Teuchos::RCP > /* sh */, - StepperIMEX_RK & /* stepperIMEX_RK */){} + Stepper & /* stepper */){} - /// Observe Stepper at end of each stage. + /// 7.) Observe Stepper at end of each stage. virtual void observeEndStage( Teuchos::RCP > /* sh */, - StepperIMEX_RK & /* stepperIMEX_RK */){} + Stepper & /* stepper */){} - /// Observe Stepper at end of takeStep. + /// 8.) Observe Stepper at end of takeStep. virtual void observeEndTakeStep( Teuchos::RCP > /* sh */, Stepper & /* stepper */){} }; } // namespace Tempus -#endif // Tempus_StepperIMEX_RKObserver_hpp +#endif // Tempus_StepperRKObserver_hpp diff --git a/packages/tempus/src/Tempus_StepperExplicitRKObserverComposite.cpp b/packages/tempus/src/Tempus_StepperRKObserverComposite.cpp similarity index 70% rename from packages/tempus/src/Tempus_StepperExplicitRKObserverComposite.cpp rename to packages/tempus/src/Tempus_StepperRKObserverComposite.cpp index 98191902d763..03faea9d243d 100644 --- a/packages/tempus/src/Tempus_StepperExplicitRKObserverComposite.cpp +++ b/packages/tempus/src/Tempus_StepperRKObserverComposite.cpp @@ -9,11 +9,11 @@ #include "Tempus_ExplicitTemplateInstantiation.hpp" #ifdef HAVE_TEMPUS_EXPLICIT_INSTANTIATION -#include "Tempus_StepperExplicitRKObserverComposite.hpp" -#include "Tempus_StepperExplicitRKObserverComposite_impl.hpp" +#include "Tempus_StepperRKObserverComposite.hpp" +#include "Tempus_StepperRKObserverComposite_impl.hpp" namespace Tempus { - TEMPUS_INSTANTIATE_TEMPLATE_CLASS(StepperExplicitRKObserverComposite) + TEMPUS_INSTANTIATE_TEMPLATE_CLASS(StepperRKObserverComposite) } #endif diff --git a/packages/tempus/src/Tempus_StepperRKObserverComposite_decl.hpp b/packages/tempus/src/Tempus_StepperRKObserverComposite_decl.hpp new file mode 100644 index 000000000000..abf21f5f766b --- /dev/null +++ b/packages/tempus/src/Tempus_StepperRKObserverComposite_decl.hpp @@ -0,0 +1,97 @@ +// @HEADER +// **************************************************************************** +// Tempus: Copyright (2017) Sandia Corporation +// +// Distributed under BSD 3-clause license (See accompanying file Copyright.txt) +// **************************************************************************** +// @HEADER + +#ifndef Tempus_StepperRKObserverComposite_decl_hpp +#define Tempus_StepperRKObserverComposite_decl_hpp + +#include "Tempus_StepperRKObserver.hpp" +#include "Tempus_TimeStepControl.hpp" +#include + +namespace Tempus { + +/** \brief This observer is a composite observer, + * + * which takes other StepperRKObservers and sequentially calls each + * individual observer function. + * + * NOTE: certain RK observers (ERK,DIRK) methods execute 'back-to-back' + * without any intermediate code. + */ +template +class StepperRKObserverComposite + : virtual public Tempus::StepperRKObserver +{ +public: + + /// Default constructor + StepperRKObserverComposite(); + + /// Destructor + virtual ~StepperRKObserverComposite(); + + /// \name Override StepperRKObserver basic methods + //@{ + /// Observe Stepper at beginning of takeStep. + virtual void observeBeginTakeStep( + Teuchos::RCP > sh, + Stepper & stepper) override; + + /// Observe Stepper at beginning of each stage. + virtual void observeBeginStage( + Teuchos::RCP > sh, + Stepper & stepperRK) override; + + /// Observe Stepper before Explicit evaluation of Implicit ODE ME. + virtual void observeBeforeImplicitExplicitly( + Teuchos::RCP > /* sh */, + Stepper & /* stepper */) override; + + /// Observe Stepper before nonlinear solve. + virtual void observeBeforeSolve( + Teuchos::RCP > sh , + Stepper & stepperRK) override; + + /// Observe Stepper after nonlinear solve. + virtual void observeAfterSolve( + Teuchos::RCP > sh , + Stepper & stepperRK) override; + + /// Observe Stepper before evaluation of Implicit ODE ME. + virtual void observeBeforeExplicit( + Teuchos::RCP > sh, + Stepper & stepperRK) override; + + /// Observe Stepper at end of each stage. + virtual void observeEndStage( + Teuchos::RCP > sh, + Stepper & stepperRK) override; + + /// Observe Stepper at end of takeStep. + virtual void observeEndTakeStep( + Teuchos::RCP > sh, + Stepper & stepper) override; + + // add observer to the composite observer list + void addObserver(const Teuchos::RCP > &observer); + + // clear all observer from the composite observer list + void clearObservers(); + + // get the number of RK stepper observers present in the composite + std::size_t getSize() const { return observers_.size(); } + //@} + +private: + + std::vector > > observers_; + +}; + +} // namespace Tempus +#endif // Tempus_StepperRKObserverComposite_decl_hpp diff --git a/packages/tempus/src/Tempus_StepperRKObserverComposite_impl.hpp b/packages/tempus/src/Tempus_StepperRKObserverComposite_impl.hpp new file mode 100644 index 000000000000..4bb97c88c010 --- /dev/null +++ b/packages/tempus/src/Tempus_StepperRKObserverComposite_impl.hpp @@ -0,0 +1,106 @@ +// @HEADER +// **************************************************************************** +// Tempus: Copyright (2017) Sandia Corporation +// +// Distributed under BSD 3-clause license (See accompanying file Copyright.txt) +// **************************************************************************** +// @HEADER + +#ifndef Tempus_StepperRKObserverComposite_impl_hpp +#define Tempus_StepperRKObserverComposite_impl_hpp + +#include "Tempus_StepperRKObserver.hpp" +#include "Tempus_TimeStepControl.hpp" + +namespace Tempus { + +template +StepperRKObserverComposite::StepperRKObserverComposite(){} + +template +StepperRKObserverComposite::~StepperRKObserverComposite(){} + +template +void StepperRKObserverComposite:: +observeBeginTakeStep(Teuchos::RCP > sh, + Stepper & stepper) +{ + for(auto& o : observers_) + o->observeBeginTakeStep(sh,stepper); +} + +template +void StepperRKObserverComposite::observeBeginStage( + Teuchos::RCP > sh, + Stepper & stepper) +{ + for(auto& o : observers_) + o->observeBeginStage(sh,stepper); +} + +template +void StepperRKObserverComposite::observeBeforeImplicitExplicitly( + Teuchos::RCP > sh, + Stepper & stepper) +{ + for(auto& o : observers_) + o->observeBeforeImplicitExplicitly(sh,stepper); +} + +template +void StepperRKObserverComposite::observeBeforeSolve( + Teuchos::RCP > sh, + Stepper & stepper) +{ + for(auto& o : observers_) + o->observeBeforeSolve(sh,stepper); +} + +template +void StepperRKObserverComposite::observeAfterSolve( + Teuchos::RCP > sh, + Stepper & stepper) +{ + for(auto& o : observers_) + o->observeAfterSolve(sh,stepper); +} + +template +void StepperRKObserverComposite::observeBeforeExplicit( + Teuchos::RCP > sh, + Stepper & stepper) +{ + for(auto& o : observers_) + o->observeBeforeExplicit(sh,stepper); +} +template +void StepperRKObserverComposite::observeEndStage( + Teuchos::RCP > sh, + Stepper & stepper) +{ + for(auto& o : observers_) + o->observeEndStage(sh,stepper); +} + +template +void StepperRKObserverComposite:: +observeEndTakeStep(Teuchos::RCP > sh, + Stepper & stepper) +{ + for(auto& o : observers_) + o->observeEndTakeStep(sh,stepper); +} + +template +void StepperRKObserverComposite:: +addObserver(const Teuchos::RCP > &observer) +{ + observers_.push_back(observer); +} + +template +void StepperRKObserverComposite:: +clearObservers() { observers_.clear();} + +} // namespace Tempus +#endif // Tempus_StepperRKObserverComposite_impl_hpp diff --git a/packages/tempus/src/Tempus_StepperRKObserverLogging.cpp b/packages/tempus/src/Tempus_StepperRKObserverLogging.cpp new file mode 100644 index 000000000000..94b801c25d5c --- /dev/null +++ b/packages/tempus/src/Tempus_StepperRKObserverLogging.cpp @@ -0,0 +1,19 @@ +// @HEADER +// **************************************************************************** +// Tempus: Copyright (2017) Sandia Corporation +// +// Distributed under BSD 3-clause license (See accompanying file Copyright.txt) +// **************************************************************************** +// @HEADER + +#include "Tempus_ExplicitTemplateInstantiation.hpp" + +#ifdef HAVE_TEMPUS_EXPLICIT_INSTANTIATION +#include "Tempus_StepperRKObserverLogging.hpp" +#include "Tempus_StepperRKObserverLogging_impl.hpp" + +namespace Tempus { + TEMPUS_INSTANTIATE_TEMPLATE_CLASS(StepperRKObserverLogging) +} + +#endif diff --git a/packages/tempus/src/Tempus_StepperRKObserverLogging_decl.hpp b/packages/tempus/src/Tempus_StepperRKObserverLogging_decl.hpp new file mode 100644 index 000000000000..1576896c1eab --- /dev/null +++ b/packages/tempus/src/Tempus_StepperRKObserverLogging_decl.hpp @@ -0,0 +1,112 @@ +// @HEADER +// **************************************************************************** +// Tempus: Copyright (2017) Sandia Corporation +// +// Distributed under BSD 3-clause license (See accompanying file Copyright.txt) +// **************************************************************************** +// @HEADER + +#ifndef Tempus_StepperRKObserverLogging_decl_hpp +#define Tempus_StepperRKObserverLogging_decl_hpp + +#include "Tempus_StepperRKObserver.hpp" +#include + +namespace Tempus { + +/** \brief This observer logs calls to observer functions. + * This observer simply logs and counts the calls to each of the + * observer functions. This is useful in monirtoring and debugging + * the time integration. + */ +template +class StepperRKObserverLogging + : virtual public Tempus::StepperRKObserver +{ +public: + + /// Constructor + StepperRKObserverLogging(); + + /// Destructor + virtual ~StepperRKObserverLogging(); + + /// \name Override StepperRKObserver basic methods + //@{ + /// 1.) Observe Stepper at beginning of takeStep. + virtual void observeBeginTakeStep( + Teuchos::RCP > /* sh */, + Stepper & /* stepper */) override; + + /// 2.) Observe Stepper at beginning of each stage. + virtual void observeBeginStage( + Teuchos::RCP > /* sh */, + Stepper & /* stepper */) override; + + /// 3.) Observe Stepper before Explicit evaluation of Implicit ODE ME (IMEX). + virtual void observeBeforeImplicitExplicitly( + Teuchos::RCP > /* sh */, + Stepper & /* stepper */) override; + + /// 4.) Observe Stepper before nonlinear solve (DIRK/IMEX). + virtual void observeBeforeSolve( + Teuchos::RCP > /* sh */, + Stepper & /* stepper */) override; + + /// 5.) Observe Stepper after nonlinear solve (DIRK/IMEX). + virtual void observeAfterSolve( + Teuchos::RCP > /* sh */, + Stepper & /* stepper */) override; + + /// 6.) Observe Stepper before Explicit evaluation of Implicit ODE ME (IMEX). + virtual void observeBeforeExplicit( + Teuchos::RCP > /* sh */, + Stepper & /* stepper */) override; + + /// 7.) Observe Stepper at end of each stage. + virtual void observeEndStage( + Teuchos::RCP > /* sh */, + Stepper & /* stepper */) override; + + /// 8.) Observe Stepper at end of takeStep. + virtual void observeEndTakeStep( + Teuchos::RCP > /* sh */, + Stepper & /* stepper */) override; + //@} + + void resetLogCounters(); + + Teuchos::RCP > getCounters(); + + Teuchos::RCP > getOrder(); + + /** \name String names logged in map + Use these strings to validate a call stack with this observer + */ + //@{ + const std::string nameObserveBeginTakeStep_; + const std::string nameObserveBeginStage_; + const std::string nameObserveBeforeImplicitExplicitly_; + const std::string nameObserveBeforeSolve_; + const std::string nameObserveAfterSolve_; + const std::string nameObserveBeforeExplicit_; + const std::string nameObserveEndStage_; + const std::string nameObserveEndTakeStep_; + //@} + +private: + + /** \brief Asserts next call on the stack is correct and removes from stack + + This is a const method so that it can be called from the + derived StepperRKObserver methods that are const. + */ + void logCall(const std::string call) const; + + Teuchos::RCP< std::map > counters_; + Teuchos::RCP< std::list > order_; + +}; + +} // namespace Tempus +#endif // Tempus_StepperRKObserverLogging_decl_hpp diff --git a/packages/tempus/src/Tempus_StepperRKObserverLogging_impl.hpp b/packages/tempus/src/Tempus_StepperRKObserverLogging_impl.hpp new file mode 100644 index 000000000000..f285037494cd --- /dev/null +++ b/packages/tempus/src/Tempus_StepperRKObserverLogging_impl.hpp @@ -0,0 +1,124 @@ +// @HEADER +// **************************************************************************** +// Tempus: Copyright (2017) Sandia Corporation +// +// Distributed under BSD 3-clause license (See accompanying file Copyright.txt) +// **************************************************************************** +// @HEADER + +#ifndef Tempus_StepperRKObserverLogging_impl_hpp +#define Tempus_StepperRKObserverLogging_impl_hpp + +#include "Tempus_StepperRKObserver.hpp" +#include "Tempus_TimeStepControl.hpp" + +namespace Tempus { + +template +StepperRKObserverLogging::StepperRKObserverLogging() + : nameObserveBeginTakeStep_ ( "observeBeginTakeStep"), + nameObserveBeginStage_ ( "observeBeginStage"), + nameObserveBeforeImplicitExplicitly_ ( "observeBeforeImplicitExplicitly"), + nameObserveBeforeSolve_ ( "observeBeforeSolve"), + nameObserveAfterSolve_ ( "observeAfterSolve"), + nameObserveBeforeExplicit_ ( "observeBeforeExplicit"), + nameObserveEndStage_ ( "observeEndStage"), + nameObserveEndTakeStep_ ( "observeEndTakeStep") +{ + counters_ = Teuchos::rcp(new std::map); + order_ = Teuchos::rcp(new std::list); + this->resetLogCounters(); +} + +template +StepperRKObserverLogging::~StepperRKObserverLogging(){} + +template +void StepperRKObserverLogging:: +observeBeginTakeStep( + Teuchos::RCP > /* sh */, + Stepper & /* stepper */) +{ logCall(nameObserveBeginTakeStep_); } + +template +void StepperRKObserverLogging:: +observeBeginStage( + Teuchos::RCP > /* sh */, + Stepper & /* stepper */) +{ logCall(nameObserveBeginStage_); } + +template +void StepperRKObserverLogging:: +observeBeforeImplicitExplicitly( + Teuchos::RCP > /* sh */, + Stepper & /* stepper */) +{ logCall(nameObserveBeforeImplicitExplicitly_); } + +template +void StepperRKObserverLogging:: +observeBeforeSolve( + Teuchos::RCP > /* sh */, + Stepper & /* stepper */) +{ logCall(nameObserveBeforeSolve_); } + +template +void StepperRKObserverLogging:: +observeAfterSolve( + Teuchos::RCP > /* sh */, + Stepper & /* stepper */) +{ logCall(nameObserveAfterSolve_); } + +template +void StepperRKObserverLogging:: +observeBeforeExplicit( + Teuchos::RCP > /* sh */, + Stepper & /* stepper */) +{ logCall(nameObserveBeforeExplicit_); } + +template +void StepperRKObserverLogging:: +observeEndStage( + Teuchos::RCP > /* sh */, + Stepper & /* stepper */) +{ logCall(nameObserveEndStage_); } + +template +void StepperRKObserverLogging:: +observeEndTakeStep( + Teuchos::RCP > /* sh */, + Stepper & /* stepper */) +{ logCall(nameObserveEndTakeStep_); } + +template +void StepperRKObserverLogging::resetLogCounters() +{ + (*counters_)[nameObserveBeginTakeStep_ ] = 0; + (*counters_)[nameObserveBeginStage_ ] = 0; + (*counters_)[nameObserveBeforeImplicitExplicitly_ ] = 0; + (*counters_)[nameObserveBeforeSolve_ ] = 0; + (*counters_)[nameObserveAfterSolve_ ] = 0; + (*counters_)[nameObserveBeforeExplicit_ ] = 0; + (*counters_)[nameObserveEndStage_ ] = 0; + (*counters_)[nameObserveEndTakeStep_ ] = 0; + order_->clear(); +} + +template +Teuchos::RCP > +StepperRKObserverLogging::getCounters() +{ return counters_; } + +template +Teuchos::RCP > +StepperRKObserverLogging::getOrder() +{ return order_; } + +template +void StepperRKObserverLogging::logCall(const std::string call) const +{ + (*counters_)[call] += 1; + order_->push_back(call); +} + +} // namespace Tempus +#endif // Tempus_StepperRKObserverLogging_impl_hpp diff --git a/packages/tempus/src/Tempus_StepperStaggeredForwardSensitivity_decl.hpp b/packages/tempus/src/Tempus_StepperStaggeredForwardSensitivity_decl.hpp index d2802c51e3ba..43c71ee1f525 100644 --- a/packages/tempus/src/Tempus_StepperStaggeredForwardSensitivity_decl.hpp +++ b/packages/tempus/src/Tempus_StepperStaggeredForwardSensitivity_decl.hpp @@ -93,6 +93,9 @@ class StepperStaggeredForwardSensitivity : virtual void setObserver( Teuchos::RCP > /* obs */ = Teuchos::null){} + virtual Teuchos::RCP > getObserver() const + { return Teuchos::null; } + /// Initialize during construction and after changing input parameters. virtual void initialize(); diff --git a/packages/tempus/src/Tempus_StepperTrapezoidal_decl.hpp b/packages/tempus/src/Tempus_StepperTrapezoidal_decl.hpp index 5b008e682988..8800a116ff38 100644 --- a/packages/tempus/src/Tempus_StepperTrapezoidal_decl.hpp +++ b/packages/tempus/src/Tempus_StepperTrapezoidal_decl.hpp @@ -64,6 +64,9 @@ class StepperTrapezoidal : virtual public Tempus::StepperImplicit virtual void setObserver( Teuchos::RCP > obs = Teuchos::null); + virtual Teuchos::RCP > getObserver() const + { return this->stepperTrapObserver_; } + /// Initialize during construction and after changing input parameters. virtual void initialize(); diff --git a/packages/tempus/src/Tempus_StepperTrapezoidal_impl.hpp b/packages/tempus/src/Tempus_StepperTrapezoidal_impl.hpp index f9e063753b65..2c8e49b7748d 100644 --- a/packages/tempus/src/Tempus_StepperTrapezoidal_impl.hpp +++ b/packages/tempus/src/Tempus_StepperTrapezoidal_impl.hpp @@ -73,13 +73,13 @@ void StepperTrapezoidal::setObserver( Teuchos::rcp(new StepperTrapezoidalObserver()); this->stepperObserver_ = Teuchos::rcp_dynamic_cast > - (stepperTrapObserver_,true); + (stepperTrapObserver_, true); } } else { this->stepperObserver_ = obs; stepperTrapObserver_ = Teuchos::rcp_dynamic_cast > - (this->stepperObserver_,true); + (this->stepperObserver_, true); } } diff --git a/packages/tempus/src/Tempus_Stepper_decl.hpp b/packages/tempus/src/Tempus_Stepper_decl.hpp index 40f2fd6bc92d..791c0002f273 100644 --- a/packages/tempus/src/Tempus_Stepper_decl.hpp +++ b/packages/tempus/src/Tempus_Stepper_decl.hpp @@ -84,6 +84,10 @@ class Stepper virtual void setObserver( Teuchos::RCP > obs = Teuchos::null) = 0; + /// Get Observer + virtual Teuchos::RCP > getObserver() const = 0; + + /// Initialize during construction and after changing input parameters. virtual void initialize() = 0; diff --git a/packages/tempus/test/Observer/Tempus_ObserverTest.cpp b/packages/tempus/test/Observer/Tempus_ObserverTest.cpp index d1192d8f7473..e30b00e49bae 100644 --- a/packages/tempus/test/Observer/Tempus_ObserverTest.cpp +++ b/packages/tempus/test/Observer/Tempus_ObserverTest.cpp @@ -16,6 +16,9 @@ #include "Tempus_IntegratorObserverLogging.hpp" #include "Tempus_IntegratorObserverComposite.hpp" +#include "Tempus_StepperRKObserverLogging.hpp" +#include "Tempus_StepperRKObserverComposite.hpp" + #include "../TestModels/SinCosModel.hpp" #include "../TestModels/VanDerPolModel.hpp" #include "../TestUtils/Tempus_ConvergenceTestUtils.hpp" @@ -34,9 +37,14 @@ using Tempus::IntegratorBasic; using Tempus::SolutionHistory; using Tempus::SolutionState; +#define TEST_INTEGRATOROBSERVERLOGGING +#define TEST_INTEGRATOROBSERVERCOMPOSITE +#define TEST_STEPPERRKOBSERVERLOGGING + // ************************************************************ // ************************************************************ +#ifdef TEST_INTEGRATOROBSERVERLOGGING TEUCHOS_UNIT_TEST(Observer, IntegratorObserverLogging) { // Read params from .xml file @@ -150,8 +158,11 @@ TEUCHOS_UNIT_TEST(Observer, IntegratorObserverLogging) Teuchos::TimeMonitor::summarize(); } +#endif // TEST_INTEGRATOROBSERVERLOGGING -TEUCHOS_UNIT_TEST( Observer, IntegratorObserverComposite) { +#ifdef TEST_INTEGRATOROBSERVERCOMPOSITE +TEUCHOS_UNIT_TEST( Observer, IntegratorObserverComposite) +{ // Read params from .xml file RCP pList = @@ -228,6 +239,142 @@ TEUCHOS_UNIT_TEST( Observer, IntegratorObserverComposite) { TEST_EQUALITY( counters.find(loggingObs2->nameObserveEndIntegrator_ )->second, 1); } +#endif + + +// ************************************************************ +// ************************************************************ +#ifdef TEST_STEPPERRKOBSERVERLOGGING +TEUCHOS_UNIT_TEST(Observer, StepperRKObserverLogging) +{ + + // number of stages for the RK method + const int erkStageCount = 4; + + // Read params from .xml file + RCP pList = + getParametersFromXmlFile("Tempus_Observer_SinCos.xml"); + + // Setup the SinCosModel + RCP scm_pl = sublist(pList, "SinCosModel", true); + RCP > model = + Teuchos::rcp(new SinCosModel (scm_pl)); + + // Setup the Integrator and reset initial time step + RCP pl = sublist(pList, "Tempus", true); + pl->sublist("Demo Integrator").set("Stepper Name", "Demo RK Stepper"); + RCP > integrator = + Tempus::integratorBasic(pl, model); + + RCP > loggingObs = + Teuchos::rcp(new Tempus::StepperRKObserverLogging); + + const auto stepper_observer = + Teuchos::rcp_dynamic_cast > + (integrator->getStepper()->getObserver(), true); + + stepper_observer->clearObservers(); + stepper_observer->addObserver(loggingObs); + + // Integrate to timeMax + bool integratorStatus = integrator->advanceTime(); + TEST_ASSERT(integratorStatus) + + // Test if at 'Final Time' + double time = integrator->getTime(); + const int numberTimeStep = integrator->getIndex(); + double timeFinal = pl->sublist("Demo Integrator") + .sublist("Time Step Control").get("Final Time"); + TEST_FLOATING_EQUALITY(time, timeFinal, 1.0e-14); + + // Construct the reference counter and order for comparison. + std::map refCounters; + std::list refOrder; + + refCounters[loggingObs->nameObserveBeginTakeStep_ ] = numberTimeStep; + refCounters[loggingObs->nameObserveBeginStage_ ] = erkStageCount * numberTimeStep; + refCounters[loggingObs->nameObserveBeforeImplicitExplicitly_] = erkStageCount * numberTimeStep; + refCounters[loggingObs->nameObserveBeforeSolve_ ] = erkStageCount * numberTimeStep; + refCounters[loggingObs->nameObserveAfterSolve_ ] = erkStageCount * numberTimeStep; + refCounters[loggingObs->nameObserveBeforeExplicit_ ] = erkStageCount * numberTimeStep; + refCounters[loggingObs->nameObserveEndStage_ ] = erkStageCount * numberTimeStep; + refCounters[loggingObs->nameObserveEndTakeStep_ ] = numberTimeStep; + + for (int i=0 ; i< numberTimeStep; ++i) { + refOrder.push_back(loggingObs->nameObserveBeginTakeStep_ ); + for (int j=0; jnameObserveBeginStage_ ); + refOrder.push_back(loggingObs->nameObserveBeforeImplicitExplicitly_ ); + refOrder.push_back(loggingObs->nameObserveBeforeSolve_ ); + refOrder.push_back(loggingObs->nameObserveAfterSolve_ ); + refOrder.push_back(loggingObs->nameObserveBeforeExplicit_ ); + refOrder.push_back(loggingObs->nameObserveEndStage_ ); + } + refOrder.push_back(loggingObs->nameObserveEndTakeStep_ ); + } + + const std::map& counters = *(loggingObs->getCounters()); + const std::list& order = *(loggingObs->getOrder()); + + // Compare against reference. + TEST_EQUALITY( + counters.find(loggingObs->nameObserveBeginTakeStep_ )->second, + refCounters.find(loggingObs->nameObserveBeginTakeStep_ )->second); + TEST_EQUALITY( + counters.find(loggingObs->nameObserveBeginStage_ )->second, + refCounters.find(loggingObs->nameObserveBeginStage_ )->second); + TEST_EQUALITY( + counters.find(loggingObs->nameObserveBeforeImplicitExplicitly_ )->second, + refCounters.find(loggingObs->nameObserveBeforeImplicitExplicitly_ )->second); + TEST_EQUALITY( + counters.find(loggingObs->nameObserveBeforeSolve_ )->second, + refCounters.find(loggingObs->nameObserveBeforeSolve_ )->second); + TEST_EQUALITY( + counters.find(loggingObs->nameObserveAfterSolve_ )->second, + refCounters.find(loggingObs->nameObserveAfterSolve_ )->second); + TEST_EQUALITY( + counters.find(loggingObs->nameObserveBeforeExplicit_ )->second, + refCounters.find(loggingObs->nameObserveBeforeExplicit_ )->second); + TEST_EQUALITY( + counters.find(loggingObs->nameObserveEndStage_ )->second, + refCounters.find(loggingObs->nameObserveEndStage_ )->second); + TEST_EQUALITY( + counters.find(loggingObs->nameObserveEndTakeStep_ )->second, + refCounters.find(loggingObs->nameObserveEndTakeStep_ )->second); + + TEUCHOS_ASSERT(order.size() == refOrder.size()); + + std::list::const_iterator orderIter = order.begin(); + std::list::const_iterator refOrderIter = refOrder.begin(); + + for ( ; orderIter != order.end(); ++orderIter,++refOrderIter) { + TEST_EQUALITY(*orderIter, *refOrderIter); + } + + // Test the reset. + loggingObs->resetLogCounters(); + TEST_EQUALITY( + counters.find(loggingObs->nameObserveBeginTakeStep_)->second, 0); + TEST_EQUALITY( + counters.find(loggingObs->nameObserveBeginStage_)->second, 0); + TEST_EQUALITY( + counters.find(loggingObs->nameObserveBeforeImplicitExplicitly_)->second, 0); + TEST_EQUALITY( + counters.find(loggingObs->nameObserveBeforeSolve_)->second, 0); + TEST_EQUALITY( + counters.find(loggingObs->nameObserveAfterSolve_)->second, 0); + TEST_EQUALITY( + counters.find(loggingObs->nameObserveBeforeExplicit_)->second, 0); + TEST_EQUALITY( + counters.find(loggingObs->nameObserveEndStage_)->second, 0); + TEST_EQUALITY( + counters.find(loggingObs->nameObserveEndTakeStep_)->second, 0); + TEST_EQUALITY(order.size(), 0); + + + Teuchos::TimeMonitor::summarize(); +} +#endif } // namespace Tempus_Test diff --git a/packages/tempus/test/Observer/Tempus_Observer_SinCos.xml b/packages/tempus/test/Observer/Tempus_Observer_SinCos.xml index ccfb01e7d2c9..06683c7bebe7 100644 --- a/packages/tempus/test/Observer/Tempus_Observer_SinCos.xml +++ b/packages/tempus/test/Observer/Tempus_Observer_SinCos.xml @@ -47,8 +47,14 @@ + + + + + + diff --git a/packages/tempus/test/PhysicsState/Tempus_PhysicsStateTest_StepperForwardEuler.hpp b/packages/tempus/test/PhysicsState/Tempus_PhysicsStateTest_StepperForwardEuler.hpp index 74ce7ffdca29..f90a65b1d74b 100644 --- a/packages/tempus/test/PhysicsState/Tempus_PhysicsStateTest_StepperForwardEuler.hpp +++ b/packages/tempus/test/PhysicsState/Tempus_PhysicsStateTest_StepperForwardEuler.hpp @@ -40,6 +40,8 @@ class StepperPhysicsStateTest } void setObserver(Teuchos::RCP > /*obs*/) {} + virtual Teuchos::RCP > getObserver() const + { return Teuchos::null; } void initialize() {} Teuchos::RCP > getDefaultStepperState() { return Teuchos::null; } diff --git a/packages/tempus/unit_test/Tempus_UnitTest_ERK_ForwardEuler.cpp b/packages/tempus/unit_test/Tempus_UnitTest_ERK_ForwardEuler.cpp index 9b4d4fb2455f..b6a7589a1a86 100644 --- a/packages/tempus/unit_test/Tempus_UnitTest_ERK_ForwardEuler.cpp +++ b/packages/tempus/unit_test/Tempus_UnitTest_ERK_ForwardEuler.cpp @@ -15,7 +15,7 @@ #include "Tempus_StepperFactory.hpp" #include "Tempus_StepperRKButcherTableau.hpp" -#include "Tempus_StepperExplicitRKObserverComposite.hpp" +#include "Tempus_StepperRKObserverComposite.hpp" #include "../TestModels/SinCosModel.hpp" #include "../TestModels/VanDerPolModel.hpp" @@ -65,7 +65,7 @@ TEUCHOS_UNIT_TEST(ERK_ForwardEuler, Arglist_Construction) { auto defaultStepper = rcp(new Tempus::StepperERK_ForwardEuler()); auto model = rcp(new Tempus_Test::SinCosModel()); - auto obs = rcp(new Tempus::StepperExplicitRKObserverComposite()); + auto obs = rcp(new Tempus::StepperRKObserverComposite()); auto stepper = rcp(new Tempus::StepperERK_ForwardEuler( model, obs, @@ -118,7 +118,7 @@ TEUCHOS_UNIT_TEST(ERK_ForwardEuler, StepperFactory_Construction) auto defaultStepper = rcp(new Tempus::StepperERK_ForwardEuler()); auto model = rcp(new Tempus_Test::SinCosModel()); - auto obs = rcp(new Tempus::StepperExplicitRKObserverComposite()); + auto obs = rcp(new Tempus::StepperRKObserverComposite()); auto stepper = rcp(new Tempus::StepperERK_ForwardEuler( model, obs, From 29ef871d501c246dcffa0c533b61b2dfaa5c2f31 Mon Sep 17 00:00:00 2001 From: Matthias Mayr Date: Mon, 23 Sep 2019 16:38:13 +0200 Subject: [PATCH 2/9] MueLu: remove regionMG coordinate debug output --- .../regionMG/examples/structured/Driver_Structured_Regions.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/packages/muelu/research/regionMG/examples/structured/Driver_Structured_Regions.cpp b/packages/muelu/research/regionMG/examples/structured/Driver_Structured_Regions.cpp index c3c7c81e6334..d216b702e86c 100644 --- a/packages/muelu/research/regionMG/examples/structured/Driver_Structured_Regions.cpp +++ b/packages/muelu/research/regionMG/examples/structured/Driver_Structured_Regions.cpp @@ -1286,7 +1286,7 @@ int main_(Teuchos::CommandLineProcessor &clp, Xpetra::UnderlyingLib& lib, int ar // std::cout << "p=" << myRank << " | compositeToRegionLIDs: " << compositeToRegionLIDs() << std::endl; // std::cout << "p=" << myRank << " | quasiRegionGIDs: " << quasiRegionGIDs << std::endl; // std::cout << "p=" << myRank << " | interfaceLIDs: " << interfaceLIDs() << std::endl; - std::cout << "p=" << myRank << " | quasiRegionCoordGIDs: " << quasiRegionCoordGIDs() << std::endl; + // std::cout << "p=" << myRank << " | quasiRegionCoordGIDs: " << quasiRegionCoordGIDs() << std::endl; // In our very particular case we know that a node is at most shared by 4 regions. // Other geometries will certainly have different constrains and a parallel reduction using MAX From 9d05ce6d51e202bbcbeeb369c3c2a0d054ce468e Mon Sep 17 00:00:00 2001 From: Brian Kelley Date: Mon, 23 Sep 2019 10:55:13 -0600 Subject: [PATCH 3/9] TrilinosCouplings: fix build errors w/o deprecated Fix missing typedefs of LO, GO when Tpetra deprecated off. Fix is the same as @DrBooom used in fd78156078, just added to these two other files. --- .../scaling/HybridIntrepidPoisson3D_Pamgen_Tpetra_main.cpp | 3 +++ .../scaling/StructuredIntrepidPoisson_Pamgen_Tpetra_main.cpp | 3 +++ 2 files changed, 6 insertions(+) diff --git a/packages/trilinoscouplings/examples/scaling/HybridIntrepidPoisson3D_Pamgen_Tpetra_main.cpp b/packages/trilinoscouplings/examples/scaling/HybridIntrepidPoisson3D_Pamgen_Tpetra_main.cpp index 6b5c6c2921d4..8b5746abeb00 100644 --- a/packages/trilinoscouplings/examples/scaling/HybridIntrepidPoisson3D_Pamgen_Tpetra_main.cpp +++ b/packages/trilinoscouplings/examples/scaling/HybridIntrepidPoisson3D_Pamgen_Tpetra_main.cpp @@ -104,6 +104,9 @@ main (int argc, char *argv[]) #ifdef HAVE_TRILINOSCOUPLINGS_MUELU typedef TpetraIntrepidPoissonExample::LO LO; typedef TpetraIntrepidPoissonExample::GO GO; +#else + typedef TpetraIntrepidPoissonExample::sparse_matrix_type::local_ordinal_type LO; + typedef TpetraIntrepidPoissonExample::sparse_matrix_type::global_ordinal_type GO; #endif // HAVE_TRILINOSCOUPLINGS_MUELU typedef TpetraIntrepidPoissonExample::Node Node; typedef Teuchos::ScalarTraits STS; diff --git a/packages/trilinoscouplings/examples/scaling/StructuredIntrepidPoisson_Pamgen_Tpetra_main.cpp b/packages/trilinoscouplings/examples/scaling/StructuredIntrepidPoisson_Pamgen_Tpetra_main.cpp index 6054ee83c7ab..5fe563218068 100644 --- a/packages/trilinoscouplings/examples/scaling/StructuredIntrepidPoisson_Pamgen_Tpetra_main.cpp +++ b/packages/trilinoscouplings/examples/scaling/StructuredIntrepidPoisson_Pamgen_Tpetra_main.cpp @@ -102,6 +102,9 @@ main (int argc, char *argv[]) #ifdef HAVE_TRILINOSCOUPLINGS_MUELU typedef TpetraIntrepidPoissonExample::LO LO; typedef TpetraIntrepidPoissonExample::GO GO; +#else + typedef TpetraIntrepidPoissonExample::sparse_matrix_type::local_ordinal_type LO; + typedef TpetraIntrepidPoissonExample::sparse_matrix_type::global_ordinal_type GO; #endif // HAVE_TRILINOSCOUPLINGS_MUELU typedef TpetraIntrepidPoissonExample::Node Node; typedef TpetraIntrepidPoissonExample::MT MT; From 5cabe9e5cf9a2d8da9ba4415dc6dcf344bb3809c Mon Sep 17 00:00:00 2001 From: Jonathan Hu Date: Mon, 23 Sep 2019 10:05:20 -0700 Subject: [PATCH 4/9] Testing: Tpetra: enable ETI explicitly --- ...i_release_tpetra_deprecated_code_off_downstream_enabled.cmake | 1 + 1 file changed, 1 insertion(+) diff --git a/cmake/ctest/drivers/rocketman/mpi_release_tpetra_deprecated_code_off_downstream_enabled.cmake b/cmake/ctest/drivers/rocketman/mpi_release_tpetra_deprecated_code_off_downstream_enabled.cmake index 680c33c9412c..690f219d2929 100644 --- a/cmake/ctest/drivers/rocketman/mpi_release_tpetra_deprecated_code_off_downstream_enabled.cmake +++ b/cmake/ctest/drivers/rocketman/mpi_release_tpetra_deprecated_code_off_downstream_enabled.cmake @@ -86,6 +86,7 @@ SET(Trilinos_CTEST_DO_ALL_AT_ONCE FALSE) # below. SET(EXTRA_CONFIGURE_OPTIONS + "-DTrilinos_ENABLE_EXPLICIT_INSTANTIATION=ON" "-DTpetra_ENABLE_DEPRECATED_CODE=OFF" "-DKOKKOS_ENABLE_DEPRECATED_CODE=OFF" "-DTPL_ENABLE_Matio=OFF" From cf5a28fb033281170cac9d79ded3546c38d84ea9 Mon Sep 17 00:00:00 2001 From: Christian Glusa Date: Mon, 23 Sep 2019 14:33:08 -0600 Subject: [PATCH 5/9] MueLu: Fix issue in nullspace fix --- .../muelu/adapters/xpetra/MueLu_RefMaxwell_def.hpp | 13 ++++++++++++- .../src/Smoothers/MueLu_Amesos2Smoother_def.hpp | 2 +- 2 files changed, 13 insertions(+), 2 deletions(-) diff --git a/packages/muelu/adapters/xpetra/MueLu_RefMaxwell_def.hpp b/packages/muelu/adapters/xpetra/MueLu_RefMaxwell_def.hpp index 8e26ac0c993d..64d2d8f8e4cb 100644 --- a/packages/muelu/adapters/xpetra/MueLu_RefMaxwell_def.hpp +++ b/packages/muelu/adapters/xpetra/MueLu_RefMaxwell_def.hpp @@ -807,8 +807,19 @@ namespace MueLu { if (!reuse) { ParameterList& userParamList = precList22_.sublist("user data"); userParamList.set >("Coordinates", Coords_); - // If we detected no boundary conditions, the (2,2) problem is singular + // If we detected no boundary conditions, the (2,2) problem is singular. + // Therefore, if we want to use a direct coarse solver, we need to fix up the nullspace. + std::string coarseType = ""; + if (precList22_.isParameter("coarse: type")) { + coarseType = precList22_.get("coarse: type"); + // Transform string to "Abcde" notation + std::transform(coarseType.begin(), coarseType.end(), coarseType.begin(), ::tolower); + std::transform(coarseType.begin(), ++coarseType.begin(), coarseType.begin(), ::toupper); + } if (BCrowcount_ == 0 && + (coarseType == "" || + coarseType == "Klu" || + coarseType == "Klu2") && (!precList22_.isSublist("coarse: params") || !precList22_.sublist("coarse: params").isParameter("fix nullspace"))) precList22_.sublist("coarse: params").set("fix nullspace",true); diff --git a/packages/muelu/src/Smoothers/MueLu_Amesos2Smoother_def.hpp b/packages/muelu/src/Smoothers/MueLu_Amesos2Smoother_def.hpp index 05f428a150ba..f4ac876e7c84 100644 --- a/packages/muelu/src/Smoothers/MueLu_Amesos2Smoother_def.hpp +++ b/packages/muelu/src/Smoothers/MueLu_Amesos2Smoother_def.hpp @@ -256,7 +256,7 @@ namespace MueLu { // add A for (size_t i = 0; i < N; i++) { for (size_t jj = rowPointers[i]; jj < rowPointers[i+1]; jj++) { - LO j = A->getColMap()->getGlobalElement(colIndices[jj]); + LO j = colMap->getLocalElement(A->getColMap()->getGlobalElement(colIndices[jj])); SC v = values[jj]; newValues[i*M+j] += v; } From 9be5dba368d11bc857de2f9853570c60b0e07cbb Mon Sep 17 00:00:00 2001 From: rstumin Date: Mon, 23 Sep 2019 14:55:06 -0600 Subject: [PATCH 6/9] mods to fix a bug in the find region code and to add more flexibility to the options that we can run and the output that we can dump. --- .../src/SetupRegionHierarchy_def.hpp | 159 +++++++++--------- .../src/SetupRegionMatrix_def.hpp | 2 +- .../structured/Driver_Structured_Regions.cpp | 35 ++++ 3 files changed, 119 insertions(+), 77 deletions(-) diff --git a/packages/muelu/research/mmayr/composite_to_regions/src/SetupRegionHierarchy_def.hpp b/packages/muelu/research/mmayr/composite_to_regions/src/SetupRegionHierarchy_def.hpp index bf7697801c58..757cd58c378c 100644 --- a/packages/muelu/research/mmayr/composite_to_regions/src/SetupRegionHierarchy_def.hpp +++ b/packages/muelu/research/mmayr/composite_to_regions/src/SetupRegionHierarchy_def.hpp @@ -313,7 +313,6 @@ void MakeCoarseLevelMaps2(const int maxRegPerGID, using MT = typename Teuchos::ScalarTraits::magnitudeType; - RCP > comm = regProlong[1][0]->getRowMap()->getComm(); const GO GO_INV = Teuchos::OrdinalTraits::invalid(); const int numLevels = regProlong.size(); @@ -1048,87 +1047,95 @@ void vCycle(const int l, ///< ID of current level RCP tm = rcp(new TimeMonitor(*TimeMonitor::getNewTimer("vCycle: * - coarsest grid solve"))); - // First get the Xpetra vectors from region to composite format - // (the coarseCompMat should already exist) - RCP compX = VectorFactory::Build(coarseCompMat->getRowMap(), true); - RCP compRhs = VectorFactory::Build(coarseCompMat->getRowMap(), true); - { - for (int j = 0; j < maxRegPerProc; j++) { - RCP inverseInterfaceScaling = VectorFactory::Build(regInterfaceScalings[l][j]->getMap()); - inverseInterfaceScaling->reciprocal(*regInterfaceScalings[l][j]); - fineRegB[j]->elementWiseMultiply(SC_ONE, *fineRegB[j], *inverseInterfaceScaling, SC_ZERO); - } - - regionalToComposite(fineRegB, compRhs, maxRegPerProc, quasiRegRowMaps[l], - regRowImporters[l], Xpetra::ADD); + const bool useCoarseSmoother = false; + if (useCoarseSmoother) { + smootherApply(smootherParams[l], maxRegPerProc, fineRegX, fineRegB, regMatrices[l], + regInterfaceScalings[l], compRowMaps[l], + quasiRegRowMaps[l], regRowMaps[l], regRowImporters[l]); } + else { + // First get the Xpetra vectors from region to composite format + // (the coarseCompMat should already exist) + RCP compX = VectorFactory::Build(coarseCompMat->getRowMap(), true); + RCP compRhs = VectorFactory::Build(coarseCompMat->getRowMap(), true); + { + for (int j = 0; j < maxRegPerProc; j++) { + RCP inverseInterfaceScaling = VectorFactory::Build(regInterfaceScalings[l][j]->getMap()); + inverseInterfaceScaling->reciprocal(*regInterfaceScalings[l][j]); + fineRegB[j]->elementWiseMultiply(SC_ONE, *fineRegB[j], *inverseInterfaceScaling, SC_ZERO); + } - const bool useDirectSolver = coarseSolverData->get("use direct solver"); - if (useDirectSolver) - { -#if defined(HAVE_MUELU_TPETRA) && defined(HAVE_MUELU_AMESOS2) + regionalToComposite(fineRegB, compRhs, maxRegPerProc, quasiRegRowMaps[l], + regRowImporters[l], Xpetra::ADD); + } - using DirectCoarseSolver = Amesos2::Solver, Tpetra::MultiVector >; - RCP coarseSolver = coarseSolverData->get >("direct solver object"); - - TEUCHOS_TEST_FOR_EXCEPT_MSG(coarseCompMat->getRowMap()->lib()!=Xpetra::UseTpetra, - "Coarse solver requires Tpetra/Amesos2 stack."); - TEUCHOS_ASSERT(!coarseSolver.is_null()); - - // using Utilities = MueLu::Utilities; - - // From here on we switch to Tpetra for simplicity - // we could also implement a similar Epetra branch - using Tpetra_MultiVector = Tpetra::MultiVector; - - // *fos << "Attempting to use Amesos2 to solve the coarse grid problem" << std::endl; - RCP tX = Utilities::MV2NonConstTpetraMV2(*compX); - RCP tB = Utilities::MV2TpetraMV(compRhs); - - /* Solve! - * - * Calling solve() on the coarseSolver should just do a triangular solve, since symbolic - * and numeric factorization are supposed to have happened during hierarchy setup. - * Here, we just check if they're done and print message if not. - * - * We don't have to change the map of tX and tB since we have configured the Amesos2 solver - * during its construction to work with non-continuous maps. - */ - if (not coarseSolver->getStatus().symbolicFactorizationDone()) - *fos << "Symbolic factorization should have been done during hierarchy setup, " - "but actually is missing. Anyway ... just do it right now." << std::endl; - if (not coarseSolver->getStatus().numericFactorizationDone()) - *fos << "Numeric factorization should have been done during hierarchy setup, " - "but actually is missing. Anyway ... just do it right now." << std::endl; - coarseSolver->solve(tX.ptr(), tB.ptr()); + const bool useDirectSolver = coarseSolverData->get("use direct solver"); + if (useDirectSolver) + { +#if defined(HAVE_MUELU_TPETRA) && defined(HAVE_MUELU_AMESOS2) + + using DirectCoarseSolver = Amesos2::Solver, Tpetra::MultiVector >; + RCP coarseSolver = coarseSolverData->get >("direct solver object"); + + TEUCHOS_TEST_FOR_EXCEPT_MSG(coarseCompMat->getRowMap()->lib()!=Xpetra::UseTpetra, + "Coarse solver requires Tpetra/Amesos2 stack."); + TEUCHOS_ASSERT(!coarseSolver.is_null()); + + // using Utilities = MueLu::Utilities; + + // From here on we switch to Tpetra for simplicity + // we could also implement a similar Epetra branch + using Tpetra_MultiVector = Tpetra::MultiVector; + + // *fos << "Attempting to use Amesos2 to solve the coarse grid problem" << std::endl; + RCP tX = Utilities::MV2NonConstTpetraMV2(*compX); + RCP tB = Utilities::MV2TpetraMV(compRhs); + + /* Solve! + * + * Calling solve() on the coarseSolver should just do a triangular solve, since symbolic + * and numeric factorization are supposed to have happened during hierarchy setup. + * Here, we just check if they're done and print message if not. + * + * We don't have to change the map of tX and tB since we have configured the Amesos2 solver + * during its construction to work with non-continuous maps. + */ + if (not coarseSolver->getStatus().symbolicFactorizationDone()) + *fos << "Symbolic factorization should have been done during hierarchy setup, " + "but actually is missing. Anyway ... just do it right now." << std::endl; + if (not coarseSolver->getStatus().numericFactorizationDone()) + *fos << "Numeric factorization should have been done during hierarchy setup, " + "but actually is missing. Anyway ... just do it right now." << std::endl; + coarseSolver->solve(tX.ptr(), tB.ptr()); #else - *fos << "+++++++++++++++++++++++++++ WARNING +++++++++++++++++++++++++\n" - << "+ Coarse level direct solver requires Tpetra and Amesos2. +\n" - << "+ Skipping the coarse level solve. +\n" - << "+++++++++++++++++++++++++++ WARNING +++++++++++++++++++++++++" - << std::endl; + *fos << "+++++++++++++++++++++++++++ WARNING +++++++++++++++++++++++++\n" + << "+ Coarse level direct solver requires Tpetra and Amesos2. +\n" + << "+ Skipping the coarse level solve. +\n" + << "+++++++++++++++++++++++++++ WARNING +++++++++++++++++++++++++" + << std::endl; #endif + } + else // use AMG as coarse level solver + { + + // Extract the hierarchy from the coarseSolverData + RCP amgHierarchy = coarseSolverData->get>("amg hierarchy object"); + + // Run a single V-cycle + amgHierarchy->Iterate(*compRhs, *compX, 1); + } + + // Transform back to region format + Array > quasiRegX(maxRegPerProc); + compositeToRegional(compX, quasiRegX, fineRegX, + maxRegPerProc, + quasiRegRowMaps[l], + regRowMaps[l], + regRowImporters[l]); + + tm = Teuchos::null; } - else // use AMG as coarse level solver - { - - // Extract the hierarchy from the coarseSolverData - RCP amgHierarchy = coarseSolverData->get>("amg hierarchy object"); - - // Run a single V-cycle - amgHierarchy->Iterate(*compRhs, *compX, 1); - } - - // Transform back to region format - Array > quasiRegX(maxRegPerProc); - compositeToRegional(compX, quasiRegX, fineRegX, - maxRegPerProc, - quasiRegRowMaps[l], - regRowMaps[l], - regRowImporters[l]); - - tm = Teuchos::null; - } + } return; } // vCycle diff --git a/packages/muelu/research/mmayr/composite_to_regions/src/SetupRegionMatrix_def.hpp b/packages/muelu/research/mmayr/composite_to_regions/src/SetupRegionMatrix_def.hpp index aa2f502d55cf..2ba53c8d0532 100644 --- a/packages/muelu/research/mmayr/composite_to_regions/src/SetupRegionMatrix_def.hpp +++ b/packages/muelu/research/mmayr/composite_to_regions/src/SetupRegionMatrix_def.hpp @@ -87,7 +87,7 @@ Teuchos::Array findCommonRegions(const GlobalOrdinal nodeA, ///< GID of fir Array regionsA, regionsB; { LO nodeALID = nodesToRegionsMap->getLocalElement(nodeA); - LO nodeBLID = nodesToRegionsMap->getLocalElement(nodeA); + LO nodeBLID = nodesToRegionsMap->getLocalElement(nodeB); for (int i = 0; i < nodesToRegions.size(); ++i) { regionsA.push_back(nodesToRegions[i][nodeALID]); regionsB.push_back(nodesToRegions[i][nodeBLID]); diff --git a/packages/muelu/research/regionMG/examples/structured/Driver_Structured_Regions.cpp b/packages/muelu/research/regionMG/examples/structured/Driver_Structured_Regions.cpp index d216b702e86c..8de2f572cfb6 100644 --- a/packages/muelu/research/regionMG/examples/structured/Driver_Structured_Regions.cpp +++ b/packages/muelu/research/regionMG/examples/structured/Driver_Structured_Regions.cpp @@ -1534,6 +1534,41 @@ int main_(Teuchos::CommandLineProcessor &clp, Xpetra::UnderlyingLib& lib, int ar Array > regB(maxRegPerProc); compositeToRegional(B, quasiRegB, regB, maxRegPerProc, rowMapPerGrp, revisedRowMapPerGrp, rowImportPerGrp); +#ifdef DUMP_LOCALX_AND_A + FILE *fp; + char str[80]; + sprintf(str,"theMatrix.%d",myRank); + fp = fopen(str,"w"); + fprintf(fp, "%%%%MatrixMarket matrix coordinate real general\n"); + LO numNzs = 0; + for (size_t kkk = 0; kkk < regionGrpMats[0]->getNodeNumRows(); kkk++) { + ArrayView AAcols; + ArrayView AAvals; + regionGrpMats[0]->getLocalRowView(kkk, AAcols, AAvals); + const int *Acols = AAcols.getRawPtr(); + const SC *Avals = AAvals.getRawPtr(); + numNzs += AAvals.size(); + } + fprintf(fp, "%d %d %d\n",regionGrpMats[0]->getNodeNumRows(),regionGrpMats[0]->getNodeNumRows(),numNzs); + + for (size_t kkk = 0; kkk < regionGrpMats[0]->getNodeNumRows(); kkk++) { + ArrayView AAcols; + ArrayView AAvals; + regionGrpMats[0]->getLocalRowView(kkk, AAcols, AAvals); + const int *Acols = AAcols.getRawPtr(); + const SC *Avals = AAvals.getRawPtr(); + LO RowLeng = AAvals.size(); + for (LO kk = 0; kk < RowLeng; kk++) { + fprintf(fp, "%d %d %22.16e\n",kkk+1,Acols[kk]+1,Avals[kk]); + } + } + fclose(fp); + sprintf(str,"theX.%d",myRank); + fp = fopen(str,"w"); + ArrayRCP lX= regX[0]->getDataNonConst(0); + for (size_t kkk = 0; kkk < regionGrpMats[0]->getNodeNumRows(); kkk++) fprintf(fp, "%22.16e\n",lX[kkk]); + fclose(fp); +#endif // printRegionalObject("regB 0", regB, myRank, *fos); From 4789956c784b0a9118e3e344ee88c860c1d3f2c9 Mon Sep 17 00:00:00 2001 From: Luc Berger-Vergiat Date: Mon, 23 Sep 2019 15:26:05 -0600 Subject: [PATCH 7/9] MueLu: adding unit-test for factoryManager Adding a unit-test to check that factoryManager is set with the correct value for useKokkos_. This required adding the GetKokkosRefactor() method to the FactoryManager class. --- .../src/Interface/MueLu_HierarchyManager.hpp | 2 +- .../MueCentral/MueLu_FactoryManager_decl.hpp | 2 + .../interface/ParameterListInterpreter.cpp | 6 +- .../test/unit_tests_kokkos/CMakeLists.txt | 1 + .../unit_tests_kokkos/UseKokkos_kokkos.cpp | 105 ++++++++++++++++++ 5 files changed, 112 insertions(+), 4 deletions(-) create mode 100644 packages/muelu/test/unit_tests_kokkos/UseKokkos_kokkos.cpp diff --git a/packages/muelu/src/Interface/MueLu_HierarchyManager.hpp b/packages/muelu/src/Interface/MueLu_HierarchyManager.hpp index 5efc6bdc5ad4..52fdc7dd675f 100644 --- a/packages/muelu/src/Interface/MueLu_HierarchyManager.hpp +++ b/packages/muelu/src/Interface/MueLu_HierarchyManager.hpp @@ -353,7 +353,7 @@ namespace MueLu { if (!M.is_null()) { Xpetra::IO::Write(fileName,* M); } - } + } else if (L->IsAvailable(name)) { // Try nofactory RCP M = L->template Get< RCP >(name); diff --git a/packages/muelu/src/MueCentral/MueLu_FactoryManager_decl.hpp b/packages/muelu/src/MueCentral/MueLu_FactoryManager_decl.hpp index 9ad22688887d..68f6bf2aa540 100644 --- a/packages/muelu/src/MueCentral/MueLu_FactoryManager_decl.hpp +++ b/packages/muelu/src/MueCentral/MueLu_FactoryManager_decl.hpp @@ -176,6 +176,8 @@ namespace MueLu { useKokkos_ = useKokkos; } + bool GetKokkosRefactor() const { return useKokkos_; } + //@} void Clean() const { defaultFactoryTable_.clear(); } diff --git a/packages/muelu/test/interface/ParameterListInterpreter.cpp b/packages/muelu/test/interface/ParameterListInterpreter.cpp index d24034bd4b00..bef2a65a91c9 100644 --- a/packages/muelu/test/interface/ParameterListInterpreter.cpp +++ b/packages/muelu/test/interface/ParameterListInterpreter.cpp @@ -78,8 +78,8 @@ int main_(Teuchos::CommandLineProcessor &clp, Xpetra::UnderlyingLib& lib, int ar using Teuchos::rcp; using Teuchos::TimeMonitor; - typedef typename Teuchos::ScalarTraits::magnitudeType real_type; - typedef Xpetra::MultiVector RealValuedMultiVector; + using real_type = typename Teuchos::ScalarTraits::coordinateType; + using RealValuedMultiVector = Xpetra::MultiVector; // ========================================================================= // MPI initialization using Teuchos @@ -119,7 +119,7 @@ int main_(Teuchos::CommandLineProcessor &clp, Xpetra::UnderlyingLib& lib, int ar matrixParameters.set("nx", Teuchos::as(9999)); matrixParameters.set("matrixType", "Laplace1D"); RCP A = MueLuTests::TestHelpers::TestFactory::Build1DPoisson(matrixParameters.get("nx"), lib); - RCP coordinates = Galeri::Xpetra::Utils::CreateCartesianCoordinates("1D", A->getRowMap(), matrixParameters); + RCP coordinates = Galeri::Xpetra::Utils::CreateCartesianCoordinates("1D", A->getRowMap(), matrixParameters); std::string prefix; if (useKokkos) { diff --git a/packages/muelu/test/unit_tests_kokkos/CMakeLists.txt b/packages/muelu/test/unit_tests_kokkos/CMakeLists.txt index 202faee342dd..d1393a6563c4 100644 --- a/packages/muelu/test/unit_tests_kokkos/CMakeLists.txt +++ b/packages/muelu/test/unit_tests_kokkos/CMakeLists.txt @@ -26,6 +26,7 @@ APPEND_SET(SOURCES StructuredAggregation_kokkos.cpp #SaPFactory_kokkos.cpp TentativePFactory_kokkos.cpp + UseKokkos_kokkos.cpp ) ### Tests that require other Trilinos packages diff --git a/packages/muelu/test/unit_tests_kokkos/UseKokkos_kokkos.cpp b/packages/muelu/test/unit_tests_kokkos/UseKokkos_kokkos.cpp new file mode 100644 index 000000000000..82a9e1fd3915 --- /dev/null +++ b/packages/muelu/test/unit_tests_kokkos/UseKokkos_kokkos.cpp @@ -0,0 +1,105 @@ +// @HEADER +// +// *********************************************************************** +// +// MueLu: A package for multigrid based preconditioning +// Copyright 2012 Sandia Corporation +// +// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation, +// the U.S. Government retains certain rights in this software. +// +// Redistribution and use in source and binary forms, with or without +// modification, are permitted provided that the following conditions are +// met: +// +// 1. Redistributions of source code must retain the above copyright +// notice, this list of conditions and the following disclaimer. +// +// 2. Redistributions in binary form must reproduce the above copyright +// notice, this list of conditions and the following disclaimer in the +// documentation and/or other materials provided with the distribution. +// +// 3. Neither the name of the Corporation nor the names of the +// contributors may be used to endorse or promote products derived from +// this software without specific prior written permission. +// +// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY +// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE +// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR +// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE +// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, +// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, +// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR +// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF +// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING +// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS +// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. +// +// Questions? Contact +// Jonathan Hu (jhu@sandia.gov) +// Andrey Prokopenko (aprokop@sandia.gov) +// Ray Tuminaro (rstumin@sandia.gov) +// +// *********************************************************************** +// +// @HEADER +#include +#include + +#include +#include + +#include "MueLu_UseDefaultTypes.hpp" + +#include "MueLu_TestHelpers_kokkos.hpp" +#include "MueLu_Version.hpp" + +#include "MueLu_ParameterListInterpreter.hpp" +#include "MueLu_FactoryManager.hpp" +#include "MueLu_HierarchyManager.hpp" +#include "MueLu_Hierarchy.hpp" + +namespace MueLuTests { + + + TEUCHOS_UNIT_TEST_TEMPLATE_4_DECL(UseKokkos_kokkos, FactoryManager, Scalar, LocalOrdinal, GlobalOrdinal, Node) + { +# include "MueLu_UseShortNames.hpp" + MUELU_TESTING_SET_OSTREAM; + MUELU_TESTING_LIMIT_SCOPE(Scalar,GlobalOrdinal,Node); + out << "version: " << MueLu::Version() << std::endl; + + using Teuchos::RCP; + using Teuchos::rcp; + using real_type = typename Teuchos::ScalarTraits::coordinateType; + using RealValuedMultiVector = Xpetra::MultiVector; + + const bool useKokkos = false; + + // Build the problem + RCP A = MueLuTests::TestHelpers_kokkos::TestFactory::Build1DPoisson(2001); + RCP map = A->getMap(); + RCP coordinates = Xpetra::MultiVectorFactory::Build(map, 1); + RCP nullspace = MultiVectorFactory::Build(map, 1); + nullspace->putScalar(Teuchos::ScalarTraits::one()); + + // Setting paramters for the hierarchy + Teuchos::ParameterList paramList; + paramList.set("verbosity", "test"); + paramList.set("use kokkos refactor", useKokkos); + RCP mueluFactory = rcp(new ParameterListInterpreter(paramList)); + mueluFactory->CheckConfig(); + const RCP manager = mueluFactory->GetFactoryManager(0); + RCP factoryManager = Teuchos::rcp_dynamic_cast(manager, true); + + TEST_EQUALITY(factoryManager->GetKokkosRefactor() == useKokkos, true); + } + + +#define MUELU_ETI_GROUP(SC,LO,GO,NO) \ + TEUCHOS_UNIT_TEST_TEMPLATE_4_INSTANT(UseKokkos_kokkos, FactoryManager, SC, LO, GO, NO) + +#include + + +} //namespace MueLuTests From 261ca660c5150600d568a1c364bd74eb9351d56b Mon Sep 17 00:00:00 2001 From: Luc Berger-Vergiat Date: Mon, 23 Sep 2019 10:42:48 -0600 Subject: [PATCH 8/9] MueLu: fixing bug with kokkos refactor in FactoryManager, see issue #5961 --- .../muelu/src/Interface/MueLu_ParameterListInterpreter_def.hpp | 1 + 1 file changed, 1 insertion(+) diff --git a/packages/muelu/src/Interface/MueLu_ParameterListInterpreter_def.hpp b/packages/muelu/src/Interface/MueLu_ParameterListInterpreter_def.hpp index 5917a2c46b40..f159420c2baa 100644 --- a/packages/muelu/src/Interface/MueLu_ParameterListInterpreter_def.hpp +++ b/packages/muelu/src/Interface/MueLu_ParameterListInterpreter_def.hpp @@ -421,6 +421,7 @@ namespace MueLu { // FIXME: should it be here, or higher up RCP defaultManager = rcp(new FactoryManager()); defaultManager->SetVerbLevel(this->verbosity_); + defaultManager->SetKokkosRefactor(useKokkos_); // We will ignore keeps0 std::vector keeps0; From 73d733433a48ffac9f3ea14086a9b54c3443af38 Mon Sep 17 00:00:00 2001 From: Luc Berger-Vergiat Date: Mon, 23 Sep 2019 16:49:15 -0600 Subject: [PATCH 9/9] MueLu: tentative fix for issue #5962 There is a bad access to a Kokkos::View data in AggregationPhase2bAlgorithm_kokkos. This should be fixed with these changes. --- .../MueLu_AggregationPhase2bAlgorithm_kokkos_def.hpp | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/packages/muelu/src/Graph/UncoupledAggregation/MueLu_AggregationPhase2bAlgorithm_kokkos_def.hpp b/packages/muelu/src/Graph/UncoupledAggregation/MueLu_AggregationPhase2bAlgorithm_kokkos_def.hpp index 99f59493705e..f1e3c9d9a9e1 100644 --- a/packages/muelu/src/Graph/UncoupledAggregation/MueLu_AggregationPhase2bAlgorithm_kokkos_def.hpp +++ b/packages/muelu/src/Graph/UncoupledAggregation/MueLu_AggregationPhase2bAlgorithm_kokkos_def.hpp @@ -141,7 +141,7 @@ namespace MueLu { // (aggStat[neigh] == AGGREGATED) if (graph.isLocalNeighborVertex(neigh) && aggStat(neigh) == AGGREGATED) - Kokkos::atomic_add(&aggWeight(vertex2AggId(neigh, 0), 0), + Kokkos::atomic_add(&aggWeight(vertex2AggId(neigh, 0)), connectWeight(neigh)); } @@ -169,7 +169,7 @@ namespace MueLu { } } if (bestScore >= 0) { - aggStat(i, 0) = AGGREGATED; + aggStat(i) = AGGREGATED; vertex2AggId(i, 0) = bestAggId; procWinner(i, 0) = myRank; @@ -244,7 +244,7 @@ namespace MueLu { // (aggStat[neigh] == AGGREGATED) if (graph.isLocalNeighborVertex(neigh) && aggStat(neigh) == AGGREGATED) - Kokkos::atomic_add(&aggWeight(vertex2AggId(neigh, 0), 0), + Kokkos::atomic_add(&aggWeight(vertex2AggId(neigh, 0)), connectWeight(neigh)); } });