24 #ifndef _FRACTIONALSTEPTHETAPROBLEM_H_
25 #define _FRACTIONALSTEPTHETAPROBLEM_H_
57 template<
typename OPTPROBLEM,
typename SPARSITYPATTERN,
typename VECTOR,
58 int dopedim,
int dealdim,
template<
int,
int>
class FE = dealii::FESystem,
59 template<
int,
int>
class DH = dealii::DoFHandler>
61 SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>
65 PrimalTSBase<OPTPROBLEM, SPARSITYPATTERN, VECTOR, dopedim, dealdim,
68 _fs_theta = 1.0 - std::sqrt(2.0) / 2.0;
69 _fs_theta_prime = 1.0 - 2.0 * _fs_theta;
70 _fs_alpha = (1.0 - 2.0 * _fs_theta) / (1.0 - _fs_theta);
71 _fs_beta = 1.0 - _fs_alpha;
72 _initial_problem = NULL;
77 if (_initial_problem != NULL)
78 delete _initial_problem;
86 return "Fractional-Step-Theta";
92 dopedim, dealdim, FE, DH>, VECTOR, dealdim>&
95 if (_initial_problem == NULL)
99 dopedim, dealdim, FE, DH>, VECTOR, dealdim>(*this);
101 return *_initial_problem;
113 template<
typename EDC>
116 dealii::Vector<double> &local_vector,
double scale,
double )
118 if (this->
GetPart() ==
"New_for_1st_and_3rd_cycle")
120 dealii::Vector<double> tmp(local_vector);
131 this->
GetProblem().ElementTimeEquation(edc, tmp, scale / (_fs_theta));
134 this->
GetProblem().ElementTimeEquationExplicit(edc, local_vector,
135 scale / (_fs_theta));
137 else if (this->
GetPart() ==
"Old_for_1st_cycle")
139 dealii::Vector<double> tmp(local_vector);
149 this->
GetProblem().ElementTimeEquation(edc, local_vector,
150 (-1) * scale / (_fs_theta));
152 else if (this->
GetPart() ==
"Old_for_3rd_cycle")
154 dealii::Vector<double> tmp(local_vector);
164 this->
GetProblem().ElementTimeEquation(edc, local_vector,
165 (-1) * scale / (_fs_theta));
167 else if (this->
GetPart() ==
"New_for_2nd_cycle")
169 dealii::Vector<double> tmp(local_vector);
180 this->
GetProblem().ElementTimeEquation(edc, tmp,
181 scale / (_fs_theta_prime));
184 this->
GetProblem().ElementTimeEquationExplicit(edc, local_vector,
185 scale / (_fs_theta_prime));
187 else if (this->
GetPart() ==
"Old_for_2nd_cycle")
189 dealii::Vector<double> tmp(local_vector);
199 this->
GetProblem().ElementTimeEquation(edc, local_vector,
200 (-1) * scale / (_fs_theta_prime));
210 template<
typename EDC>
213 dealii::Vector<double> &local_vector,
double scale)
215 if (this->
GetPart() ==
"New_for_1st_and_3rd_cycle")
219 else if (this->
GetPart() ==
"Old_for_1st_cycle"
220 || this->
GetPart() ==
"Old_for_3rd_cycle")
222 this->
GetProblem().ElementRhs(edc, local_vector,
226 else if (this->
GetPart() ==
"New_for_2nd_cycle")
228 this->
GetProblem().ElementRhs(edc, local_vector,
232 else if (this->
GetPart() ==
"Old_for_2nd_cycle")
244 const std::map<std::string,
const dealii::Vector<double>*> ¶m_values,
245 const std::map<std::string, const VECTOR*> &domain_values,
246 VECTOR& rhs_vector,
double scale)
248 if (this->
GetPart() ==
"New_for_1st_and_3rd_cycle")
252 else if (this->
GetPart() ==
"Old_for_1st_cycle"
253 || this->
GetPart() ==
"Old_for_3rd_cycle")
255 this->
GetProblem().PointRhs(param_values, domain_values, rhs_vector,
259 else if (this->
GetPart() ==
"New_for_2nd_cycle")
261 this->
GetProblem().PointRhs(param_values, domain_values, rhs_vector,
265 else if (this->
GetPart() ==
"Old_for_2nd_cycle")
276 template<
typename EDC>
279 dealii::FullMatrix<double> &local_matrix)
281 if (this->
GetPart() ==
"New_for_1st_and_3rd_cycle")
283 dealii::FullMatrix<double> m(local_matrix);
290 local_matrix.add(1.0, m);
294 local_matrix.add(1.0 / (_fs_theta), m);
297 this->
GetProblem().ElementTimeMatrixExplicit(edc, m);
298 local_matrix.add(1.0 / (_fs_theta), m);
300 else if (this->
GetPart() ==
"New_for_2nd_cycle")
302 dealii::FullMatrix<double> m(local_matrix);
304 this->
GetProblem().ElementMatrix(edc, local_matrix,
309 local_matrix.add(1.0, m);
313 local_matrix.add(1.0 / (_fs_theta_prime), m);
316 this->
GetProblem().ElementTimeMatrixExplicit(edc, m);
317 local_matrix.add(1.0 / (_fs_theta_prime), m);
323 template<
typename FDC>
326 dealii::Vector<double> &local_vector,
double scale,
double)
328 if (this->
GetPart() ==
"New_for_1st_and_3rd_cycle")
330 this->
GetProblem().FaceEquation(fdc, local_vector,
337 else if ((this->
GetPart() ==
"Old_for_1st_cycle")
338 || (this->
GetPart() ==
"Old_for_3rd_cycle"))
340 this->
GetProblem().FaceEquation(fdc, local_vector,
345 else if (this->
GetPart() ==
"New_for_2nd_cycle")
347 this->
GetProblem().FaceEquation(fdc, local_vector,
353 else if (this->
GetPart() ==
"Old_for_2nd_cycle")
355 this->
GetProblem().FaceEquation(fdc, local_vector,
369 template<
typename FDC>
372 dealii::Vector<double> &local_vector,
double scale,
375 if (this->
GetPart() ==
"New_for_1st_and_3rd_cycle")
377 this->
GetProblem().InterfaceEquation(fdc, local_vector,
384 else if ((this->
GetPart() ==
"Old_for_1st_cycle")
385 || (this->
GetPart() ==
"Old_for_3rd_cycle"))
387 this->
GetProblem().InterfaceEquation(fdc, local_vector,
392 else if (this->
GetPart() ==
"New_for_2nd_cycle")
394 this->
GetProblem().InterfaceEquation(fdc, local_vector,
400 else if (this->
GetPart() ==
"Old_for_2nd_cycle")
402 this->
GetProblem().InterfaceEquation(fdc, local_vector,
416 template<
typename FDC>
419 dealii::Vector<double> &local_vector,
double scale = 1.)
428 template<
typename FDC>
431 dealii::FullMatrix<double> &local_matrix)
433 if (this->
GetPart() ==
"New_for_1st_and_3rd_cycle")
435 dealii::FullMatrix<double> m(local_matrix);
442 local_matrix.add(1., m);
444 else if (this->
GetPart() ==
"New_for_2nd_cycle")
446 dealii::FullMatrix<double> m(local_matrix);
453 local_matrix.add(1.0, m);
460 template<
typename FDC>
463 dealii::FullMatrix<double> &local_matrix)
465 if (this->
GetPart() ==
"New_for_1st_and_3rd_cycle")
467 dealii::FullMatrix<double> m(local_matrix);
474 local_matrix.add(1.0, m);
476 else if (this->
GetPart() ==
"New_for_2nd_cycle")
478 dealii::FullMatrix<double> m(local_matrix);
485 local_matrix.add(1.0, m);
492 template<
typename FDC>
495 dealii::Vector<double> &local_vector,
double scale,
double)
497 if (this->
GetPart() ==
"New_for_1st_and_3rd_cycle")
499 this->
GetProblem().BoundaryEquation(fdc, local_vector,
505 else if ((this->
GetPart() ==
"Old_for_1st_cycle")
506 || (this->
GetPart() ==
"Old_for_3rd_cycle"))
508 this->
GetProblem().BoundaryEquation(fdc, local_vector,
513 else if (this->
GetPart() ==
"New_for_2nd_cycle")
515 this->
GetProblem().BoundaryEquation(fdc, local_vector,
521 else if (this->
GetPart() ==
"Old_for_2nd_cycle")
523 this->
GetProblem().BoundaryEquation(fdc, local_vector,
537 template<
typename FDC>
540 dealii::Vector<double> &local_vector,
double scale)
542 this->
GetProblem().BoundaryRhs(fdc, local_vector,
549 template<
typename FDC>
552 dealii::FullMatrix<double> &local_matrix)
554 if (this->
GetPart() ==
"New_for_1st_and_3rd_cycle")
556 dealii::FullMatrix<double> m(local_matrix);
563 local_matrix.add(1.0, m);
565 else if (this->
GetPart() ==
"New_for_2nd_cycle")
567 dealii::FullMatrix<double> m(local_matrix);
574 local_matrix.add(1.0, m);
580 double _fs_theta_prime;
586 dopedim, dealdim, FE, DH>, VECTOR, dealdim> * _initial_problem;
std::string GetName()
Definition: fractional_step_theta_problem.h:84
~FractionalStepThetaProblem()
Definition: fractional_step_theta_problem.h:75
void FaceRhs(const FDC &fdc, dealii::Vector< double > &local_vector, double scale=1.)
Definition: fractional_step_theta_problem.h:418
FractionalStepThetaProblem< OPTPROBLEM, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH > & GetBaseProblem()
Definition: fractional_step_theta_problem.h:107
void BoundaryMatrix(const FDC &fdc, dealii::FullMatrix< double > &local_matrix)
Definition: fractional_step_theta_problem.h:551
void InterfaceMatrix(const FDC &fdc, dealii::FullMatrix< double > &local_matrix)
Definition: fractional_step_theta_problem.h:462
Definition: fractional_step_theta_problem.h:60
Definition: primal_ts_base.h:49
void FaceEquation(const FDC &fdc, dealii::Vector< double > &local_vector, double scale, double)
Definition: fractional_step_theta_problem.h:325
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: fractional_step_theta_problem.h:243
void ElementEquation(const EDC &edc, dealii::Vector< double > &local_vector, double scale, double)
Definition: fractional_step_theta_problem.h:115
void ElementRhs(const EDC &edc, dealii::Vector< double > &local_vector, double scale)
Definition: fractional_step_theta_problem.h:212
Definition: initialproblem.h:42
InitialProblem< FractionalStepThetaProblem< OPTPROBLEM, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH >, VECTOR, dealdim > & GetInitialProblem()
Definition: fractional_step_theta_problem.h:93
void InterfaceEquation(const FDC &fdc, dealii::Vector< double > &local_vector, double scale, double)
Definition: fractional_step_theta_problem.h:371
const SpaceTimeHandler< FE, DH, SPARSITYPATTERN, VECTOR, dopedim, dealdim > * GetSpaceTimeHandler() const
Definition: ts_base.h:442
OPTPROBLEM & GetProblem()
Definition: ts_base.h:466
void BoundaryEquation(const FDC &fdc, dealii::Vector< double > &local_vector, double scale, double)
Definition: fractional_step_theta_problem.h:494
void FaceMatrix(const FDC &fdc, dealii::FullMatrix< double > &local_matrix)
Definition: fractional_step_theta_problem.h:430
FractionalStepThetaProblem(OPTPROBLEM &OP)
Definition: fractional_step_theta_problem.h:64
std::string GetPart() const
Definition: ts_base.h:479
void BoundaryRhs(const FDC &fdc, dealii::Vector< double > &local_vector, double scale)
Definition: fractional_step_theta_problem.h:539
void ElementMatrix(const EDC &edc, dealii::FullMatrix< double > &local_matrix)
Definition: fractional_step_theta_problem.h:278