9 #ifndef Tempus_RKButcherTableau_hpp
10 #define Tempus_RKButcherTableau_hpp
14 #pragma clang system_header
20 #include "Tempus_config.hpp"
44 template<
class Scalar>
52 std::string stepperType,
80 Scalar sumb = ST::zero();
81 for (
size_t i = 0; i < this->
numStages(); i++) sumb +=
b_(i);
84 "Error - Butcher Tableau b fails to satisfy Sum(b_i) = 1.\n"
85 <<
" Sum(b_i) = " << sumb <<
"\n");
89 for (
size_t i = 0; i < this->
numStages(); i++) {
90 Scalar sumai = ST::zero();
91 for (
size_t j = 0; j < this->
numStages(); j++) sumai +=
A_(i,j);
93 if (std::abs(sumai) > 1.0e-08)
94 failed = (std::abs((sumai-
c_(i))/sumai) > 1.0e-08);
96 failed = (std::abs(
c_(i)) > 1.0e-08);
99 "Error - Butcher Tableau fails to satisfy c_i = Sum_j(a_ij).\n"
100 <<
" Stage i = " << i <<
"\n"
101 <<
" c_i = " <<
c_(i) <<
"\n"
102 <<
" Sum_j(a_ij) = " << sumai <<
"\n"
103 <<
" This may be OK as some valid tableaus do not satisfy\n"
104 <<
" this condition. If so, construct this RKButcherTableau\n"
105 <<
" with checkC = false.\n");
165 out <<
"number of Stages = " << this->
numStages() << std::endl;
166 out <<
"A = " <<
printMat(this->
A()) << std::endl;
167 out <<
"b = " <<
printMat(this->
b()) << std::endl;
168 out <<
"c = " <<
printMat(this->
c()) << std::endl;
170 out <<
"order = " << this->
order() << std::endl;
171 out <<
"orderMin = " << this->
orderMin() << std::endl;
172 out <<
"orderMax = " << this->
orderMax() << std::endl;
173 out <<
"isImplicit = " << this->
isImplicit() << std::endl;
174 out <<
"isDIRK = " << this->
isDIRK() << std::endl;
175 out <<
"isEmbedded = " << this->
isEmbedded() << std::endl;
177 out <<
"TVD Coeff = " << this->
getTVDCoeff() << std::endl;
183 const Scalar relTol = 1.0e-15;
191 if(std::abs((t.
A_(i,j) -
A_(i,j))/
A_(i,j)) > relTol)
return false;
203 if(std::abs((t.
b_(i) -
b_(i))/
b_(i)) > relTol)
return false;
204 if(std::abs((t.
c_(i) -
c_(i))/
c_(i)) > relTol)
return false;
212 return !((*this) == t);
223 for (
size_t i = 0; i < this->
numStages(); i++)
224 for (
size_t j = i; j < this->
numStages(); j++)
231 bool nonZero =
false;
232 for (
size_t i = 0; i < this->
numStages(); i++) {
233 if (
A_(i,i) != 0.0) nonZero =
true;
234 for (
size_t j = i+1; j < this->
numStages(); j++)
237 if (nonZero ==
false)
isDIRK_ =
false;
260 #endif // Tempus_RKButcherTableau_hpp
Teuchos::SerialDenseMatrix< int, Scalar > A_
virtual void setEmbedded(bool a)
virtual const Teuchos::SerialDenseVector< int, Scalar > & b() const
Return the vector of quadrature weights.
virtual const Teuchos::SerialDenseVector< int, Scalar > & c() const
Return the vector of stage positions.
RKButcherTableau(std::string stepperType, const Teuchos::SerialDenseMatrix< int, Scalar > &A, const Teuchos::SerialDenseVector< int, Scalar > &b, const Teuchos::SerialDenseVector< int, Scalar > &c, const int order, const int orderMin, const int orderMax, const Teuchos::SerialDenseVector< int, Scalar > &bstar=Teuchos::SerialDenseVector< int, Scalar >(), bool checkC=true)
bool operator!=(const RKButcherTableau &t) const
SerialBandDenseMatrixPrinter< OrdinalType, ScalarType > printMat(const SerialBandDenseMatrix< OrdinalType, ScalarType > &obj)
virtual const Teuchos::SerialDenseVector< int, Scalar > & bstar() const
Return the vector of quadrature weights for embedded methods.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
virtual void setTVDCoeff(const Scalar a)
set TVD coefficient of RK method
virtual bool isTVD() const
Return true if the RK method is TVD.
virtual void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const
Teuchos::SerialDenseVector< int, Scalar > bstar_
virtual void setTVD(bool a)
virtual bool isImplicit() const
Return true if the RK method is implicit.
virtual bool isEmbedded() const
Return true if the RK method has embedded capabilities.
virtual bool isDIRK() const
Return true if the RK method is Diagonally Implicit.
virtual int orderMax() const
Return the maximum order.
virtual Scalar getTVDCoeff() const
Return TVD coefficient of RK method.
void set_isDIRK()
DIRK is defined as if a_ij = 0 for j>i and a_ii != 0 for at least one i.
Teuchos::SerialDenseVector< int, Scalar > b_
virtual std::string description() const
OrdinalType length() const
basic_FancyOStream & setOutputToRootOnly(const int rootRank)
OrdinalType numCols() const
virtual std::size_t numStages() const
Return the number of stages.
Teuchos::SerialDenseVector< int, Scalar > c_
#define TEUCHOS_ASSERT(assertion_test)
#define TEUCHOS_ASSERT_EQUALITY(val1, val2)
virtual const Teuchos::SerialDenseMatrix< int, Scalar > & A() const
Return the matrix coefficients.
virtual int order() const
Return the order.
virtual int orderMin() const
Return the minimum order.
bool operator==(const RKButcherTableau &t) const
OrdinalType numRows() const