24 #ifndef _ForwardEulerProblem_H_
25 #define _ForwardEulerProblem_H_
52 template<
typename OPTPROBLEM,
typename SPARSITYPATTERN,
typename VECTOR,
53 int dopedim,
int dealdim,
54 template <
int,
int>
class FE = dealii::FESystem,
55 template <
int,
int>
class DH = dealii::DoFHandler>
57 SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>
61 PrimalTSBase<OPTPROBLEM, SPARSITYPATTERN, VECTOR, dopedim, dealdim,
64 _initial_problem = NULL;
68 if (_initial_problem != NULL)
69 delete _initial_problem;
77 return "forward Euler";
85 if (_initial_problem == NULL)
90 return *_initial_problem;
102 template<
typename EDC>
105 dealii::Vector<double> &local_vector,
double scale,
double )
109 dealii::Vector<double> tmp(local_vector);
115 this->
GetProblem().ElementTimeEquation(edc, tmp,
119 this->
GetProblem().ElementTimeEquationExplicit(edc, local_vector,
123 else if (this->
GetPart() ==
"Old")
125 dealii::Vector<double> tmp(local_vector);
143 template<
typename EDC>
146 dealii::Vector<double> &local_vector,
double scale)
152 else if (this->
GetPart() ==
"Old")
165 const std::map<std::string,
const dealii::Vector<double>*> ¶m_values,
166 const std::map<std::string, const VECTOR*> &domain_values,
167 VECTOR& rhs_vector,
double scale)
172 else if (this->
GetPart() ==
"Old")
184 template<
typename EDC>
187 dealii::FullMatrix<double> &local_matrix)
189 assert(this->
GetPart() ==
"New");
190 dealii::FullMatrix<double> m(local_matrix);
200 this->
GetProblem().ElementTimeMatrixExplicit(edc, m);
208 template<
typename FDC>
211 dealii::Vector<double> &local_vector,
double scale,
double )
217 else if (this->
GetPart() ==
"Old")
229 template<
typename FDC>
232 dealii::Vector<double> &local_vector,
double scale,
double )
238 else if (this->
GetPart() ==
"Old")
251 template<
typename FDC>
254 dealii::Vector<double> &local_vector,
double scale = 1.)
261 template<
typename FDC>
264 dealii::FullMatrix<double> &local_matrix)
266 assert(this->
GetPart() ==
"New");
273 template<
typename FDC>
276 dealii::FullMatrix<double> &local_matrix)
278 assert(this->
GetPart() ==
"New");
285 template<
typename FDC>
288 dealii::Vector<double> &local_vector,
double scale,
double )
294 else if (this->
GetPart() ==
"Old")
307 template<
typename FDC>
310 dealii::Vector<double> &local_vector,
double scale)
317 template<
typename FDC>
320 dealii::FullMatrix<double> &local_matrix)
322 assert(this->
GetPart() ==
"New");
void BoundaryRhs(const FDC &fdc, dealii::Vector< double > &local_vector, double scale)
Definition: forward_euler_problem.h:309
void PointRhs(const std::map< std::string, const dealii::Vector< double > * > ¶m_values, const std::map< std::string, const VECTOR * > &domain_values, VECTOR &rhs_vector, double scale)
Definition: forward_euler_problem.h:164
void FaceEquation(const FDC &fdc, dealii::Vector< double > &local_vector, double scale, double)
Definition: forward_euler_problem.h:210
void ElementMatrix(const EDC &edc, dealii::FullMatrix< double > &local_matrix)
Definition: forward_euler_problem.h:186
ForwardEulerProblem< OPTPROBLEM, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH > & GetBaseProblem()
Definition: forward_euler_problem.h:96
Definition: primal_ts_base.h:49
void BoundaryEquation(const FDC &fdc, dealii::Vector< double > &local_vector, double scale, double)
Definition: forward_euler_problem.h:287
void FaceMatrix(const FDC &fdc, dealii::FullMatrix< double > &local_matrix)
Definition: forward_euler_problem.h:263
void InterfaceEquation(const FDC &fdc, dealii::Vector< double > &local_vector, double scale, double)
Definition: forward_euler_problem.h:231
void FaceRhs(const FDC &fdc, dealii::Vector< double > &local_vector, double scale=1.)
Definition: forward_euler_problem.h:253
~ForwardEulerProblem()
Definition: forward_euler_problem.h:66
ForwardEulerProblem(OPTPROBLEM &OP)
Definition: forward_euler_problem.h:60
Definition: initialproblem.h:42
const SpaceTimeHandler< FE, DH, SPARSITYPATTERN, VECTOR, dopedim, dealdim > * GetSpaceTimeHandler() const
Definition: ts_base.h:442
OPTPROBLEM & GetProblem()
Definition: ts_base.h:466
void ElementEquation(const EDC &edc, dealii::Vector< double > &local_vector, double scale, double)
Definition: forward_euler_problem.h:104
Definition: forward_euler_problem.h:56
std::string GetPart() const
Definition: ts_base.h:479
std::string GetName()
Definition: forward_euler_problem.h:75
InitialProblem< ForwardEulerProblem< OPTPROBLEM, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH >, VECTOR, dealdim > & GetInitialProblem()
Definition: forward_euler_problem.h:83
void InterfaceMatrix(const FDC &fdc, dealii::FullMatrix< double > &local_matrix)
Definition: forward_euler_problem.h:275
void BoundaryMatrix(const FDC &fdc, dealii::FullMatrix< double > &local_matrix)
Definition: forward_euler_problem.h:319
void ElementRhs(const EDC &edc, dealii::Vector< double > &local_vector, double scale)
Definition: forward_euler_problem.h:145