24 #ifndef _SHIFTEDCRANKNICOLSONProblem_H_
25 #define _SHIFTEDCRANKNICOLSONProblem_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 "shifted Crank-Nicolson";
84 if (_initial_problem == NULL)
87 VECTOR, dealdim>(*this);
89 return *_initial_problem;
101 template<
typename EDC>
104 dealii::Vector<double> &local_vector,
double scale,
double )
108 damped_cn_theta = 0.5
109 + this->
GetProblem().GetBaseProblem().GetSpaceTimeHandler()->GetStepSize();
111 dealii::Vector<double> tmp(local_vector);
121 this->
GetProblem().ElementTimeEquation(dc, tmp,
125 this->
GetProblem().ElementTimeEquationExplicit(dc, local_vector,
128 else if (this->
GetPart() ==
"Old")
130 damped_cn_theta = 0.5
131 + this->
GetProblem().GetBaseProblem().GetSpaceTimeHandler()->GetStepSize();
133 dealii::Vector<double> tmp(local_vector);
141 this->
GetProblem().ElementTimeEquation(dc, local_vector,
152 template<
typename EDC>
155 dealii::Vector<double> &local_vector,
double scale)
159 damped_cn_theta = 0.5
160 + this->
GetProblem().GetBaseProblem().GetSpaceTimeHandler()->GetStepSize();
163 else if (this->
GetPart() ==
"Old")
165 damped_cn_theta = 0.5
166 + this->
GetProblem().GetBaseProblem().GetSpaceTimeHandler()->GetStepSize();
167 this->
GetProblem().ElementRhs(edc, local_vector,
180 const std::map<std::string,
const dealii::Vector<double>*> ¶m_values,
181 const std::map<std::string, const VECTOR*> &domain_values,
182 VECTOR& rhs_vector,
double scale)
186 damped_cn_theta = 0.5
187 + this->
GetProblem().GetBaseProblem().GetSpaceTimeHandler()->GetStepSize();
188 this->
GetProblem().PointRhs(param_values, domain_values, rhs_vector,
191 else if (this->
GetPart() ==
"Old")
193 damped_cn_theta = 0.5
194 + this->
GetProblem().GetBaseProblem().GetSpaceTimeHandler()->GetStepSize();
195 this->
GetProblem().PointRhs(param_values, domain_values, rhs_vector,
206 template<
typename EDC>
209 dealii::FullMatrix<double> &local_matrix)
211 assert(this->
GetPart() ==
"New");
212 damped_cn_theta = 0.5
213 + this->
GetProblem().GetBaseProblem().GetSpaceTimeHandler()->GetStepSize();
214 dealii::FullMatrix<double> m(local_matrix);
228 this->
GetProblem().ElementTimeMatrixExplicit(edc, m);
236 template<
typename FDC>
239 dealii::Vector<double> &local_vector,
double scale,
double )
243 damped_cn_theta = 0.5
244 + this->
GetProblem().GetBaseProblem().GetSpaceTimeHandler()->GetStepSize();
247 else if (this->
GetPart() ==
"Old")
249 damped_cn_theta = 0.5
250 + this->
GetProblem().GetBaseProblem().GetSpaceTimeHandler()->GetStepSize();
262 template<
typename FDC>
265 dealii::Vector<double> &local_vector,
double scale,
double )
269 damped_cn_theta = 0.5
270 + this->
GetProblem().GetBaseProblem().GetSpaceTimeHandler()->GetStepSize();
271 this->
GetProblem().InterfaceEquation(fdc, local_vector,
274 else if (this->
GetPart() ==
"Old")
276 damped_cn_theta = 0.5
277 + this->
GetProblem().GetBaseProblem().GetSpaceTimeHandler()->GetStepSize();
289 template<
typename FDC>
292 dealii::Vector<double> &local_vector,
double scale = 1.)
299 template<
typename FDC>
302 dealii::FullMatrix<double> &local_matrix)
304 assert(this->
GetPart() ==
"New");
305 damped_cn_theta = 0.5
306 + this->
GetProblem().GetBaseProblem().GetSpaceTimeHandler()->GetStepSize();
308 dealii::FullMatrix<double> m(local_matrix);
312 this->
GetProblem().FaceMatrix(fdc, m, damped_cn_theta, 1.);
320 template<
typename FDC>
323 dealii::FullMatrix<double> &local_matrix)
325 assert(this->
GetPart() ==
"New");
326 damped_cn_theta = 0.5
327 + this->
GetProblem().GetBaseProblem().GetSpaceTimeHandler()->GetStepSize();
329 dealii::FullMatrix<double> m(local_matrix);
333 this->
GetProblem().InterfaceMatrix(fdc, m, damped_cn_theta, 1.);
341 template<
typename FDC>
344 dealii::Vector<double> &local_vector,
double scale,
double )
348 damped_cn_theta = 0.5
349 + this->
GetProblem().GetBaseProblem().GetSpaceTimeHandler()->GetStepSize();
350 this->
GetProblem().BoundaryEquation(fdc, local_vector,
353 else if (this->
GetPart() ==
"Old")
355 damped_cn_theta = 0.5
356 + this->
GetProblem().GetBaseProblem().GetSpaceTimeHandler()->GetStepSize();
357 this->
GetProblem().BoundaryEquation(fdc, local_vector,
369 template<
typename FDC>
372 dealii::Vector<double> &local_vector,
double scale)
379 template<
typename FDC>
382 dealii::FullMatrix<double> &local_matrix)
384 assert(this->
GetPart() ==
"New");
385 damped_cn_theta = 0.5
386 + this->
GetProblem().GetBaseProblem().GetSpaceTimeHandler()->GetStepSize();
387 dealii::FullMatrix<double> m(local_matrix);
391 this->
GetProblem().BoundaryMatrix(fdc, m, damped_cn_theta, 1.);
397 double damped_cn_theta;
void BoundaryEquation(const FDC &fdc, dealii::Vector< double > &local_vector, double scale, double)
Definition: shifted_crank_nicolson_problem.h:343
void BoundaryMatrix(const FDC &fdc, dealii::FullMatrix< double > &local_matrix)
Definition: shifted_crank_nicolson_problem.h:381
void FaceMatrix(const FDC &fdc, dealii::FullMatrix< double > &local_matrix)
Definition: shifted_crank_nicolson_problem.h:301
InitialProblem< ShiftedCrankNicolsonProblem, VECTOR, dealdim > & GetInitialProblem()
Definition: shifted_crank_nicolson_problem.h:82
void ElementEquation(const EDC &dc, dealii::Vector< double > &local_vector, double scale, double)
Definition: shifted_crank_nicolson_problem.h:103
Definition: primal_ts_base.h:49
void FaceEquation(const FDC &fdc, dealii::Vector< double > &local_vector, double scale, double)
Definition: shifted_crank_nicolson_problem.h:238
std::string GetName()
Definition: shifted_crank_nicolson_problem.h:75
ShiftedCrankNicolsonProblem(OPTPROBLEM &OP)
Definition: shifted_crank_nicolson_problem.h:60
Definition: shifted_crank_nicolson_problem.h:56
void InterfaceMatrix(const FDC &fdc, dealii::FullMatrix< double > &local_matrix)
Definition: shifted_crank_nicolson_problem.h:322
Definition: initialproblem.h:42
~ShiftedCrankNicolsonProblem()
Definition: shifted_crank_nicolson_problem.h:66
void BoundaryRhs(const FDC &fdc, dealii::Vector< double > &local_vector, double scale)
Definition: shifted_crank_nicolson_problem.h:371
const SpaceTimeHandler< FE, DH, SPARSITYPATTERN, VECTOR, dopedim, dealdim > * GetSpaceTimeHandler() const
Definition: ts_base.h:442
void ElementRhs(const EDC &edc, dealii::Vector< double > &local_vector, double scale)
Definition: shifted_crank_nicolson_problem.h:154
OPTPROBLEM & GetProblem()
Definition: ts_base.h:466
void FaceRhs(const FDC &fdc, dealii::Vector< double > &local_vector, double scale=1.)
Definition: shifted_crank_nicolson_problem.h:291
ShiftedCrankNicolsonProblem & GetBaseProblem()
Definition: shifted_crank_nicolson_problem.h:94
std::string GetPart() const
Definition: ts_base.h:479
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: shifted_crank_nicolson_problem.h:179
void ElementMatrix(const EDC &edc, dealii::FullMatrix< double > &local_matrix)
Definition: shifted_crank_nicolson_problem.h:208
void InterfaceEquation(const FDC &fdc, dealii::Vector< double > &local_vector, double scale, double)
Definition: shifted_crank_nicolson_problem.h:264