24 #ifndef _BackwardEulerProblem_H_
25 #define _BackwardEulerProblem_H_
50 template<
typename OPTPROBLEM,
typename SPARSITYPATTERN,
typename VECTOR,
51 int dopedim,
int dealdim,
template<
int,
int>
class FE = dealii::FESystem,
52 template<
int,
int>
class DH = dealii::DoFHandler>
54 SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>
58 PrimalTSBase<OPTPROBLEM, SPARSITYPATTERN, VECTOR, dopedim, dealdim,
61 _initial_problem = NULL;
65 if (_initial_problem != NULL)
66 delete _initial_problem;
79 return "backward Euler";
89 dealdim, FE, DH>, VECTOR, dealdim>&
92 if (_initial_problem == NULL)
96 dopedim, dealdim, FE, DH>, VECTOR, dealdim>(*this);
98 return *_initial_problem;
149 template<
typename EDC>
152 dealii::Vector<double> &local_vector,
double scale,
double )
156 dealii::Vector<double> tmp(local_vector);
166 this->
GetProblem().ElementTimeEquation(edc, tmp, scale);
169 this->
GetProblem().ElementTimeEquationExplicit(edc, local_vector,
172 else if (this->
GetPart() ==
"Old")
174 this->
GetProblem().ElementTimeEquation(edc, local_vector,
204 template<
typename EDC>
207 dealii::Vector<double> &local_vector,
double scale)
211 this->
GetProblem().ElementRhs(edc, local_vector,
215 else if (this->
GetPart() ==
"Old")
252 const std::map<std::string,
const dealii::Vector<double>*> ¶m_values,
253 const std::map<std::string, const VECTOR*> &domain_values,
254 VECTOR& rhs_vector,
double scale)
258 this->
GetProblem().PointRhs(param_values, domain_values, rhs_vector,
262 else if (this->
GetPart() ==
"Old")
304 template<
typename EDC>
307 dealii::FullMatrix<double> &local_matrix)
309 assert(this->
GetPart() ==
"New");
310 dealii::FullMatrix<double> m(local_matrix);
312 this->
GetProblem().ElementMatrix(edc, local_matrix,
320 local_matrix.add(1.0, m);
323 this->
GetProblem().ElementTimeMatrixExplicit(edc, m);
324 local_matrix.add(1.0, m);
345 template<
typename FDC>
348 dealii::Vector<double> &local_vector,
double scale,
353 this->
GetProblem().FaceEquation(fdc, local_vector,
359 else if (this->
GetPart() ==
"Old")
386 template<
typename FDC>
389 dealii::Vector<double> &local_vector,
double scale,
394 this->
GetProblem().InterfaceEquation(fdc, local_vector,
400 else if (this->
GetPart() ==
"Old")
425 template<
typename FDC>
428 dealii::Vector<double> &local_vector,
double scale = 1.)
448 template<
typename FDC>
451 dealii::FullMatrix<double> &local_matrix)
453 assert(this->
GetPart() ==
"New");
454 this->
GetProblem().FaceMatrix(fdc, local_matrix,
475 template<
typename FDC>
478 dealii::FullMatrix<double> &local_matrix)
480 assert(this->
GetPart() ==
"New");
481 this->
GetProblem().InterfaceMatrix(fdc, local_matrix,
506 template<
typename FDC>
509 dealii::Vector<double> &local_vector,
double scale,
514 this->
GetProblem().BoundaryEquation(fdc, local_vector,
520 else if (this->
GetPart() ==
"Old")
545 template<
typename FDC>
548 dealii::Vector<double> &local_vector,
double scale)
550 this->
GetProblem().BoundaryRhs(fdc, local_vector,
569 template<
typename FDC>
572 dealii::FullMatrix<double> &local_matrix)
574 assert(this->
GetPart() ==
"New");
575 this->
GetProblem().BoundaryMatrix(fdc, local_matrix,
585 dealdim, FE, DH>, VECTOR, dealdim> * _initial_problem;
void FaceRhs(const FDC &fdc, dealii::Vector< double > &local_vector, double scale=1.)
Definition: backward_euler_problem.h:427
void ElementEquation(const EDC &edc, dealii::Vector< double > &local_vector, double scale, double)
Definition: backward_euler_problem.h:151
void InterfaceEquation(const FDC &fdc, dealii::Vector< double > &local_vector, double scale, double)
Definition: backward_euler_problem.h:388
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: backward_euler_problem.h:251
void BoundaryMatrix(const FDC &fdc, dealii::FullMatrix< double > &local_matrix)
Definition: backward_euler_problem.h:571
void ElementRhs(const EDC &edc, dealii::Vector< double > &local_vector, double scale)
Definition: backward_euler_problem.h:206
BackwardEulerProblem(OPTPROBLEM &OP)
Definition: backward_euler_problem.h:57
Definition: primal_ts_base.h:49
Definition: backward_euler_problem.h:53
void BoundaryEquation(const FDC &fdc, dealii::Vector< double > &local_vector, double scale, double)
Definition: backward_euler_problem.h:508
InitialProblem< BackwardEulerProblem< OPTPROBLEM, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH >, VECTOR, dealdim > & GetInitialProblem()
Definition: backward_euler_problem.h:90
BackwardEulerProblem< OPTPROBLEM, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH > & GetBaseProblem()
Definition: backward_euler_problem.h:111
void BoundaryRhs(const FDC &fdc, dealii::Vector< double > &local_vector, double scale)
Definition: backward_euler_problem.h:547
~BackwardEulerProblem()
Definition: backward_euler_problem.h:63
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
std::string GetName()
Definition: backward_euler_problem.h:77
std::string GetPart() const
Definition: ts_base.h:479
void FaceEquation(const FDC &fdc, dealii::Vector< double > &local_vector, double scale, double)
Definition: backward_euler_problem.h:347
void FaceMatrix(const FDC &fdc, dealii::FullMatrix< double > &local_matrix)
Definition: backward_euler_problem.h:450
void InterfaceMatrix(const FDC &fdc, dealii::FullMatrix< double > &local_matrix)
Definition: backward_euler_problem.h:477
void ElementMatrix(const EDC &edc, dealii::FullMatrix< double > &local_matrix)
Definition: backward_euler_problem.h:306