| 02
+ | Use ModelEvaluator
+ | Use ModelEvaluator for problem description (under construction).
-
Also, one can look at the tests and unit tests for examples of usage,
- but these may not show standard usage or context.
+ but these may not show standard usage or context as they are meant for
+ testing.
*/
diff --git a/packages/tempus/examples/00_Basic_Problem/00_Basic_Problem.cpp b/packages/tempus/examples/00_Basic_Problem/00_Basic_Problem.cpp
index 39774629687e..0c9f0f3439e4 100644
--- a/packages/tempus/examples/00_Basic_Problem/00_Basic_Problem.cpp
+++ b/packages/tempus/examples/00_Basic_Problem/00_Basic_Problem.cpp
@@ -246,7 +246,7 @@ int main(int argc, char *argv[])
if ( std::isnan(x_n[0]) || std::isnan(x_n[1]) ) {
passing = false;
} else {
- // Promote to next step (n -> n+1).
+ // Promote to next step (n <- n+1).
n++;
x_n[0] = x_np1[0];
x_n[1] = x_np1[1];
diff --git a/packages/tempus/examples/01_Utilize_Thyra/01_Utilize_Thyra.cpp b/packages/tempus/examples/01_Utilize_Thyra/01_Utilize_Thyra.cpp
new file mode 100644
index 000000000000..02ff6b19199c
--- /dev/null
+++ b/packages/tempus/examples/01_Utilize_Thyra/01_Utilize_Thyra.cpp
@@ -0,0 +1,292 @@
+// @HEADER
+// ****************************************************************************
+// Tempus: Copyright (2017) Sandia Corporation
+//
+// Distributed under BSD 3-clause license (See accompanying file Copyright.txt)
+// ****************************************************************************
+// @HEADER
+
+#include
+#include
+#include
+#include
+#include "Teuchos_StandardCatchMacros.hpp"
+
+#include "Thyra_VectorStdOps.hpp"
+#include "Thyra_DefaultSpmdVectorSpace.hpp"
+#include "Thyra_DetachedVectorView.hpp"
+
+
+using namespace std;
+using Teuchos::RCP;
+
+/** @file */
+
+/** \brief Description:
+ *
+ * \section example-01 Example 1: Utilize Thyra
+ *
+ * This problem takes \ref 00_Basic_Problem.cpp "Basic Problem"
+ * and utilizes Thyra vectors, replacing the standard double vectors,
+ * 01_Utilize_Thyra.cpp. Tempus uses Thyra for its Abstract Numerical
+ * Algorithms (ANAs), which are the mathematical concepts of
+ * vectors, vector spaces, and linear operators. All other ANA
+ * interfaces and support software are built on these fundamental
+ * operator/vector interfaces (including ModelEvaluators).
+ *
+ * \subsection declaration Declaration and Initialization
+ *
+ * We first need to create a vector space to construct the Thyra::Vector
+ * we need. Thus we replace
+ * @code
+ * // Solution and its time-derivative.
+ * double x_n[2]; // at time index n
+ * double xDot_n[2]; // at time index n
+ * @endcode
+ * with
+ * @code
+ * // Solution and its time-derivative.
+ * int vectorLength = 2; // number state unknowns
+ * RCP > xSpace =
+ * Thyra::defaultSpmdVectorSpace(vectorLength);
+ *
+ * RCP > x_n = Thyra::createMember(xSpace);
+ * RCP > xDot_n = Thyra::createMember(xSpace);
+ * @endcode
+ *
+ * The initialization can be achieved by replacing
+ * @code
+ * // Initial Conditions
+ * double time = 0.0;
+ * double epsilon = 1.0e-1;
+ * x_n [0] = 2.0;
+ * x_n [1] = 0.0;
+ * xDot_n[0] = 0.0;
+ * xDot_n[1] = -2.0/epsilon;
+ * @endcode
+ * with
+ * @code
+ * // Initial Conditions
+ * double time = 0.0;
+ * double epsilon = 1.0e-1;
+ * { // Scope to delete DetachedVectorViews
+ * Thyra::DetachedVectorView x_n_view(*x_n);
+ * x_n_view[0] = 2.0;
+ * x_n_view[1] = 0.0;
+ * Thyra::DetachedVectorView xDot_n_view(*xDot_n);
+ * xDot_n_view[0] = 0.0;
+ * xDot_n_view[1] = -2.0/epsilon;
+ * }
+ * @endcode
+ * and the scoped Thyra::DetachedVectorView's are used to delete them
+ * after use.
+ *
+ * Elements of the thyra::Vector can be quickly accessed with get_ele by
+ * replacing
+ * @code
+ * cout << n << " " << time << " " << x_n[0] << " " << x_n[1] << endl;
+ * @endcode
+ * with
+ * @code
+ * cout << n << " " << time << " " << get_ele(*(x_n), 0)
+ * << " " << get_ele(*(x_n), 1) << endl;
+ * @endcode
+ *
+ * To initialize the solution at the next time step, we can simply
+ * replace
+ * @code
+ * // Initialize next time step
+ * double x_np1[2]; // at time index n+1
+ * x_np1[0] = x_n[0];
+ * x_np1[1] = x_n[1];
+ * @endcode
+ * with
+ * @code
+ * // Initialize next time step
+ * RCP > x_np1 = x_n->clone_v(); // at time index n+1
+ * @endcode
+ *
+ * The model evaluation is again achieved through Thyra::DetachedVectorView
+ * by replacing
+ * @code
+ * // Righthand side evaluation and time-derivative at n.
+ * xDot_n[0] = x_n[1];
+ * xDot_n[1] = ((1.0 - x_n[0]*x_n[0])*x_n[1] - x_n[0])/epsilon;
+ * @endcode
+ * with
+ * @code
+ * // Righthand side evaluation and time-derivative at n.
+ * {
+ * Thyra::ConstDetachedVectorView x_n_view(*x_n);
+ * Thyra::DetachedVectorView xDot_n_view(*xDot_n);
+ * xDot_n_view[0] = x_n_view[1];
+ * xDot_n_view[1] =
+ * ((1.0-x_n_view[0]*x_n_view[0])*x_n_view[1]-x_n_view[0])/epsilon;
+ * }
+ * @endcode
+ *
+ * The Forward Euler time step
+ * @code
+ * // Take the timestep - Forward Euler
+ * x_np1[0] = x_n[0] + dt*xDot_n[0];
+ * x_np1[1] = x_n[1] + dt*xDot_n[1];
+ * @endcode
+ * can be replaced with
+ * @code
+ * // Take the timestep - Forward Euler
+ * Thyra::V_VpStV(x_np1.ptr(), *x_n, dt, *xDot_n);
+ * @endcode
+ * where Thyra::V_VpStV performs an axpy.
+ *
+ * We can also take advantage other Thyra features, Thyra::norm, to check
+ * if the solution is passing, i.e.,
+ * @code
+ * // Test if solution is passing.
+ * if ( std::isnan(x_n[0]) || std::isnan(x_n[1]) ) {
+ * @endcode
+ * is replaced with
+ * @code
+ * // Test if solution is passing.
+ * if ( std::isnan(Thyra::norm(*x_np1) ) {
+ * @endcode
+ *
+ * The solution update is changed from
+ * @code
+ * // Promote to next step (n <- n+1).
+ * n++;
+ * x_n[0] = x_np1[0];
+ * x_n[1] = x_np1[1];
+ * @endcode
+ * to
+ * @code
+ * // Promote to next step (n <- n+1).
+ * n++;
+ * Thyra::V_V(x_n.ptr(), *x_np1);
+ * @endcode
+ *
+ * The remainder of 01_Utilize_Thyra.cpp (e.g., regression testing) has
+ * similar changes to those from above.
+ *
+ * \subsection example-01_Next Links
+ *
+ * - Back to: \ref tutorials
+ * - Previous: \ref 00_Basic_Problem.cpp "Basic Problem"
+ * - Next: Use ModelEvaluator
+ */
+int main(int argc, char *argv[])
+{
+ bool verbose = true;
+ bool success = false;
+ try {
+ // Solution and its time-derivative.
+ int vectorLength = 2; // number state unknowns
+ RCP > xSpace =
+ Thyra::defaultSpmdVectorSpace(vectorLength);
+
+ RCP > x_n = Thyra::createMember(xSpace);
+ RCP > xDot_n = Thyra::createMember(xSpace);
+
+ // Initial Conditions
+ double time = 0.0;
+ double epsilon = 1.0e-1;
+ { // Scope to delete DetachedVectorViews
+ Thyra::DetachedVectorView x_n_view(*x_n);
+ x_n_view[0] = 2.0;
+ x_n_view[1] = 0.0;
+ Thyra::DetachedVectorView xDot_n_view(*xDot_n);
+ xDot_n_view[0] = 0.0;
+ xDot_n_view[1] = -2.0/epsilon;
+ }
+
+ // Timestep size
+ double finalTime = 2.0;
+ int nTimeSteps = 2001;
+ const double constDT = finalTime/(nTimeSteps-1);
+
+ // Advance the solution to the next timestep.
+ int n = 0;
+ bool passing = true;
+ cout << n << " " << time << " " << get_ele(*(x_n), 0)
+ << " " << get_ele(*(x_n), 1) << endl;
+ while (passing && time < finalTime && n < nTimeSteps) {
+
+ // Initialize next time step
+ RCP > x_np1 = x_n->clone_v(); // at time index n+1
+
+ // Set the timestep and time.
+ double dt = constDT;
+ time = (n+1)*dt;
+
+ // Righthand side evaluation and time-derivative at n.
+ {
+ Thyra::ConstDetachedVectorView x_n_view(*x_n);
+ Thyra::DetachedVectorView xDot_n_view(*xDot_n);
+ xDot_n_view[0] = x_n_view[1];
+ xDot_n_view[1] =
+ ((1.0-x_n_view[0]*x_n_view[0])*x_n_view[1]-x_n_view[0])/epsilon;
+ }
+
+ // Take the timestep - Forward Euler
+ Thyra::V_VpStV(x_np1.ptr(), *x_n, dt, *xDot_n);
+
+ // Test if solution is passing.
+ if ( std::isnan(Thyra::norm(*x_np1)) ) {
+ passing = false;
+ } else {
+ // Promote to next step (n <- n+1).
+ n++;
+ Thyra::V_V(x_n.ptr(), *x_np1);
+ }
+
+ // Output
+ if ( n%100 == 0 )
+ cout << n << " " << time << " " << get_ele(*(x_n), 0)
+ << " " << get_ele(*(x_n), 1) << endl;
+ }
+
+ // Test for regression.
+ RCP > x_regress = x_n->clone_v();
+ {
+ Thyra::DetachedVectorView x_regress_view(*x_regress);
+ x_regress_view[0] = -1.59496108218721311;
+ x_regress_view[1] = 0.96359412806611255;
+ }
+
+ RCP > x_error = x_n->clone_v();
+ Thyra::V_VmV(x_error.ptr(), *x_n, *x_regress);
+ double x_L2norm_error = Thyra::norm_2(*x_error );
+ double x_L2norm_regress = Thyra::norm_2(*x_regress);
+
+ cout << "Relative L2 Norm of the error (regression) = "
+ << x_L2norm_error/x_L2norm_regress << endl;
+ if ( x_L2norm_error > 1.0e-08*x_L2norm_regress) {
+ passing = false;
+ cout << "FAILED regression constraint!" << endl;
+ }
+
+ RCP > x_best = x_n->clone_v();
+ {
+ Thyra::DetachedVectorView x_best_view(*x_best);
+ x_best_view[0] = -1.59496108218721311;
+ x_best_view[1] = 0.96359412806611255;
+ }
+
+ Thyra::V_VmV(x_error.ptr(), *x_n, *x_best);
+ x_L2norm_error = Thyra::norm_2(*x_error);
+ double x_L2norm_best = Thyra::norm_2(*x_best );
+
+ cout << "Relative L2 Norm of the error (best) = "
+ << x_L2norm_error/x_L2norm_best << endl;
+ if ( x_L2norm_error > 0.02*x_L2norm_best) {
+ passing = false;
+ cout << "FAILED best constraint!" << endl;
+ }
+ if (passing) success = true;
+ }
+ TEUCHOS_STANDARD_CATCH_STATEMENTS(verbose, std::cerr, success);
+
+ if(success)
+ cout << "\nEnd Result: Test Passed!" << std::endl;
+
+ return ( success ? EXIT_SUCCESS : EXIT_FAILURE );
+}
diff --git a/packages/tempus/examples/01_Utilize_Thyra/CMakeLists.txt b/packages/tempus/examples/01_Utilize_Thyra/CMakeLists.txt
new file mode 100644
index 000000000000..a1eb75b337b0
--- /dev/null
+++ b/packages/tempus/examples/01_Utilize_Thyra/CMakeLists.txt
@@ -0,0 +1,7 @@
+INCLUDE_DIRECTORIES(${CMAKE_CURRENT_SOURCE_DIR})
+
+TRIBITS_ADD_EXECUTABLE_AND_TEST(
+ 01_Utilize_Thyra
+ SOURCES 01_Utilize_Thyra.cpp
+ NUM_MPI_PROCS 1
+ )
diff --git a/packages/tempus/examples/CMakeLists.txt b/packages/tempus/examples/CMakeLists.txt
index 21162fbe91aa..55d1823e7ba8 100644
--- a/packages/tempus/examples/CMakeLists.txt
+++ b/packages/tempus/examples/CMakeLists.txt
@@ -1 +1,2 @@
ADD_SUBDIRECTORY(00_Basic_Problem)
+ADD_SUBDIRECTORY(01_Utilize_Thyra)
|