#include <iomanip>
#include <iostream>
#include <stdlib.h>
#include <math.h>
#include "Teuchos_StandardCatchMacros.hpp"
#include "Thyra_VectorStdOps.hpp"
#include "Thyra_DefaultSpmdVectorSpace.hpp"
#include "Thyra_DetachedVectorView.hpp"
#include "../02_Use_ModelEvaluator/VanDerPol_ModelEvaluator_02.hpp"
#include "Tempus_SolutionState.hpp"

Functions
int	main (int argc, char *argv[])
	Description: More...

Function Documentation

int main	(	int	argc,
		char *	argv[]
	)

Description:

Example 3: Introduce SolutionState

The basic idea behind the SolutionState object is that it contains all the information about the solution at a particular time, i.e., $x(t)$ . The primary requirement for SolutionState is that it contains all the information needed to restart the solution from scratch, e.g., check-pointing and recover from a failed timestep. The four primary components of the Tempus::SolutionState are

The state - , $\dot{x}(t)$ (optional), and $\ddot{x}(t)$ (optional)
The metastate - timestep, time, order, error, etc.
The StepperState - any data needed by the stepper
The PhysicsState - any application-specific data

For additional details about the Tempus::SolutionState, please see SolutionState Description!

For this example, we will only utilize and demonstrate the state and the metastate (Tempus::SolutionStateMetaData).

Changes to Use SolutionState

In the following table are code snippets of the changes from the 02_Use_ModelEvaluator.cpp to 03_Intro_SolutionState.cpp. This is similar taking a difference between them for the main() function.

Comments	"Use ModelEvaluator" Code Snippet	"Intro SolutionState" Code Snippet
Initial the SolutionState. Create a SolutionState object with the initial-condition solution, , from the ModelEvaluator, using the one of the non-member constructors. These non-member constructors set the remaining member data to default values (see Tempus::createSolutionStateX). Member data can also be set through the basic set functions, e.g., setIndex(0). Note: we did not include its time derivative, $\dot{x}$ , as part of the SolutionState because it is not needed for output or diagnostics. However we still require memory for it in the Forward Euler timestepping. Foreshadowing: the time derivative can be handled and maintained by all the time steppers internally.	// Get initial conditions from ModelEvaluator::getNominalValues. int n = 0; double time = 0.0; RCP<Thyra::VectorBase<double> > x_n = model->getNominalValues().get_x()->clone_v(); RCP<Thyra::VectorBase<double> > xDot_n = model->getNominalValues().get_x_dot()->clone_v();	// Setup initial condition SolutionState -------------------- auto solState = Tempus::createSolutionStateX( model->getNominalValues().get_x()->clone_v()); solState->setIndex (0); solState->setTime (0.0); solState->setTimeStep(0.0); // By convention, the IC has dt=0. solState->setSolutionStatus(Tempus::Status::PASSED); // ICs are considered passed. RCP<Thyra::VectorBase<double> > xDot_n = model->getNominalValues().get_x_dot()->clone_v();
Access member data of the SolutionState. All the information about the solution and state are within the SolutionState object, and can be accessed through get functions, e.g., the index, time, solution, and solution status.	cout << n << " " << time << " " << get_ele((x_n), 0) << " " << get_ele((x_n), 1) << endl; while (passed && time < finalTime && n < nTimeSteps) {	cout << solState->getIndex() << " " << solState->getTime() << " " << get_ele((solState->getX()), 0) << " " << get_ele((solState->getX()), 1) << endl; while (solState->getSolutionStatus() == Tempus::Status::PASSED && solState->getTime() < finalTime && solState->getIndex() < nTimeSteps) {
Setup the next timestep. To initialize the next time step, we set some Referenced Counted Pointers (RCPs), and create a scratch vector. Additionally, we set the index, time and timestep for the next ("working") time step solution. Although, not evident at this point, we will use the label "working" to help manage solutions that have not met the "passing" criteria, and we can try again with new settings, e.g., timestep size.	// Initialize next time step RCP<Thyra::VectorBase<double> > x_np1 = x_n->clone_v(); // at time index n+1 // Set the timestep and time. double dt = constDT; time = (n+1)*dt;	// Initialize next time step RCP<Thyra::VectorBase<double> > x_n = solState->getX(); RCP<Thyra::VectorBase<double> > x_np1 = solState->getX()->clone_v(); // at time index n+1 solState->setSolutionStatus(Tempus::Status::WORKING); // Set the timestep and time for the working solution i.e., n+1. int index = solState->getIndex()+1; double dt = constDT; double time = index*dt;
Set status of SolutionState. Need to set SolutionState status to PASSED or FAILED, increment the index and time, and set the timestep.	// Test if solution has passed. if ( std::isnan(Thyra::norm(x_np1)) ) { passed = false; } else { // Promote to next step (n <- n+1). Thyra::V_V(x_n.ptr(), x_np1); n++; }	// Test if solution has passed. if ( std::isnan(Thyra::norm(x_np1)) ) { solState->setSolutionStatus(Tempus::Status::FAILED); } else { // Promote to next step (n <- n+1). Thyra::V_V(x_n.ptr(), x_np1); solState->setIndex (index); solState->setTime (time); solState->setTimeStep(constDT); solState->setSolutionStatus(Tempus::Status::PASSED); }
Access the SolutionState. Need to access the SolutionState for the solution index. Note the RCP x_n still has access to the solution vector.	// Output if ( n%100 == 0 ) cout << n << " " << time << " " << get_ele((x_n), 0) << " " << get_ele((x_n), 1) << endl;	// Output if ( solState->getIndex()%100 == 0 ) cout << solState->getIndex() << " " << time << " " << get_ele((x_n), 0) << " " << get_ele((x_n), 1) << endl;

The remainder of 03_Intro_SolutionState.cpp is regression testing, etc.

Links

Back to: Tutorials
Previous: Example 2: Use ModelEvaluator
Next: Adding SolutionHistory

Definition at line 213 of file 03_Intro_SolutionState.cpp.

Functions

Function Documentation

Example 3: Introduce SolutionState

Changes to Use SolutionState

Links