24 #ifndef _PDEProblemContainer_H_
25 #define _PDEProblemContainer_H_
39 #include <deal.II/multigrid/mg_dof_handler.h>
43 #include <lac/vector.h>
44 #include <lac/full_matrix.h>
45 #include <grid/tria_iterator.h>
46 #include <dofs/dof_handler.h>
47 #include <dofs/dof_accessor.h>
48 #include <dofs/dof_tools.h>
49 #include <fe/fe_system.h>
50 #include <fe/fe_values.h>
51 #include <fe/mapping.h>
52 #include <base/quadrature_lib.h>
53 #include <lac/block_sparsity_pattern.h>
54 #include <lac/compressed_simple_sparsity_pattern.h>
57 #include <deal.II/multigrid/mg_dof_handler.h>
58 #include <deal.II/multigrid/mg_constrained_dofs.h>
59 #include <deal.II/multigrid/multigrid.h>
60 #include <deal.II/multigrid/mg_transfer.h>
61 #include <deal.II/multigrid/mg_tools.h>
62 #include <deal.II/multigrid/mg_coarse.h>
63 #include <deal.II/multigrid/mg_smoother.h>
64 #include <deal.II/multigrid/mg_matrix.h>
75 template<
typename VECTOR>
76 class DOpEOutputHandler;
77 template<
typename VECTOR>
78 class DOpEExceptionHandler;
97 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
98 int dealdim,
template<
int,
int>
class FE = dealii::FESystem,
99 template<
int,
int>
class DH = dealii::DoFHandler>
116 return "PDEProblemContainer";
125 DH>, PDE, DD, SPARSITYPATTERN, VECTOR, dealdim>&
128 if (_state_problem == NULL)
132 FE, DH>, PDE, DD, SPARSITYPATTERN, VECTOR, dealdim>(*
this,
135 return *_state_problem;
155 ReInit(std::string algo_type);
170 _ExceptionHandler = OH;
183 SetType(std::string type,
unsigned int num = 0);
206 template<
typename DATACONTAINER>
226 const std::map<std::string,
const dealii::Vector<double>*> ¶m_values,
227 const std::map<std::string, const VECTOR*> &domain_values);
250 template<
typename FACEDATACONTAINER>
262 template<
typename FACEDATACONTAINER>
275 const std::map<std::string,
const dealii::Vector<double>*> &values,
276 const std::map<std::string, const VECTOR*> &block_values);
303 template<
typename DATACONTAINER>
306 dealii::Vector<double> &local_vector,
double scale,
318 template<
typename DATACONTAINER>
321 dealii::Vector<double> &local_vector,
double scale = 1.);
331 template<
typename DATACONTAINER>
334 dealii::Vector<double> &local_vector,
double scale = 1.);
349 template<
typename DATACONTAINER>
352 dealii::Vector<double> &local_vector,
double scale = 1.);
370 const std::map<std::string,
const dealii::Vector<double>*> ¶m_values,
371 const std::map<std::string, const VECTOR*> &domain_values,
372 VECTOR& rhs_vector,
double scale = 1.);
400 template<
typename DATACONTAINER>
403 dealii::FullMatrix<double> &local_entry_matrix,
double scale = 1.,
404 double scale_ico = 1.);
419 template<
typename DATACONTAINER>
422 dealii::FullMatrix<double> &local_entry_matrix);
435 template<
typename DATACONTAINER>
438 dealii::FullMatrix<double> &local_entry_matrix);
448 template<
typename FACEDATACONTAINER>
451 dealii::Vector<double> &local_vector,
double scale,
461 template<
typename FACEDATACONTAINER>
464 dealii::Vector<double> &local_vector,
double scale,
474 template<
typename FACEDATACONTAINER>
476 FaceRhs(
const FACEDATACONTAINER& dc,
477 dealii::Vector<double> &local_vector,
double scale = 1.);
487 template<
typename FACEDATACONTAINER>
490 dealii::FullMatrix<double> &local_entry_matrix,
double scale = 1.,
491 double scale_ico = 1.);
500 template<
typename FACEDATACONTAINER>
503 dealii::FullMatrix<double> &local_entry_matrix,
double scale = 1.,
504 double scale_ico = 1.);
513 template<
typename FACEDATACONTAINER>
516 dealii::Vector<double> &local_vector,
double scale,
526 template<
typename FACEDATACONTAINER>
529 dealii::Vector<double> &local_vector,
double scale = 1.);
539 template<
typename FACEDATACONTAINER>
542 dealii::FullMatrix<double> &local_matrix,
double scale = 1.,
543 double scale_ico = 1.);
547 const FE<dealdim, dealdim>&
590 const std::vector<bool>& comp_mask,
const DD* values);
594 const std::vector<unsigned int>&
596 const std::vector<bool>&
601 const dealii::Function<dealdim> &
603 const std::map<std::string,
const dealii::Vector<double>*> ¶m_values,
604 const std::map<std::string, const VECTOR*> &domain_values)
const;
611 _initial_values = values;
613 const dealii::Function<dealdim>&
616 return *_initial_values;
623 const std::vector<unsigned int>&
627 const std::vector<unsigned int>&
645 _aux_functionals.push_back(F);
646 if (_functional_position.find(F->GetName())
647 != _functional_position.end())
650 "You cant add two functionals with the same name.",
651 "PDEProblemContainer::AddFunctional");
653 _functional_position[F->GetName()] = _aux_functionals.size() - 1;
676 if (_aux_functionals[i]->
GetName() == functional_name)
680 _functional_for_ee_num = i;
689 "Can't find functional " + functional_name
690 +
" in _aux_functionals",
691 "PDEProblemContainer::SetFunctionalForErrorEstimation");
700 return _aux_functionals.size();
720 const std::vector<unsigned int>&
725 const dealii::ConstraintMatrix&
747 assert(_ExceptionHandler);
748 return _ExceptionHandler;
756 assert(_OutputHandler);
757 return _OutputHandler;
793 unsigned int n_levels)
const;
802 unsigned int n_levels)
const;
806 template<
typename INTEGRATOR>
815 template<
typename INTEGRATOR>
827 return this->
GetPDE().GetStateNBlocks();
832 std::vector<unsigned int>&
835 return this->
GetPDE().GetStateBlockComponent();
846 const std::map<std::string, unsigned int>&
849 return _functional_position;
857 std::string _algo_type;
861 VECTOR, dealdim>*> _aux_functionals;
862 std::map<std::string, unsigned int> _functional_position;
864 unsigned int _functional_for_ee_num;
867 std::vector<unsigned int> _dirichlet_colors;
868 std::vector<std::vector<bool> > _dirichlet_comps;
869 std::vector<PrimalDirichletData<DD, VECTOR, dealdim>*> _primal_dirichlet_values;
870 const dealii::Function<dealdim>* _zero_dirichlet_values;
872 const dealii::Function<dealdim>* _initial_values;
874 std::vector<unsigned int> _state_boundary_equation_colors;
875 std::vector<unsigned int> _adjoint_boundary_equation_colors;
877 std::vector<unsigned int> _boundary_functional_colors;
881 DH>, PDE, DD, SPARSITYPATTERN, VECTOR, dealdim>* _state_problem;
885 DH>, PDE, DD, SPARSITYPATTERN, VECTOR, dealdim> ;
889 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
890 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
896 _ExceptionHandler = NULL;
897 _OutputHandler = NULL;
898 _zero_dirichlet_values =
new ZeroFunction<dealdim>(
899 this->
GetPDE().GetStateNComponents());
901 _functional_for_ee_num = dealii::numbers::invalid_unsigned_int;
906 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
907 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
910 if (_zero_dirichlet_values != NULL)
912 delete _zero_dirichlet_values;
914 for (
unsigned int i = 0; i < _primal_dirichlet_values.size(); i++)
916 if (_primal_dirichlet_values[i] != NULL)
917 delete _primal_dirichlet_values[i];
919 if (_state_problem != NULL)
921 delete _state_problem;
927 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
928 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
931 std::string algo_type)
933 if (_state_problem != NULL)
935 delete _state_problem;
936 _state_problem = NULL;
939 if (_algo_type != algo_type && _algo_type !=
"")
942 "PDEProblemContainer::ReInit");
946 _algo_type = algo_type;
947 this->SetTypeInternal(
"");
949 if (_algo_type ==
"reduced")
951 GetSpaceTimeHandler()->ReInit(this->GetPDE().GetStateNBlocks(),
952 this->GetPDE().GetStateBlockComponent());
957 "PDEProblemContainer::ReInit");
964 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
965 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
968 std::string type,
unsigned int num)
970 if (this->GetType() != type || this->GetTypeNum() != num)
972 this->SetTypeInternal(type);
973 this->SetTypeNumInternal(num);
974 this->GetPDE().SetProblemType(type);
975 if (_functional_for_ee_num != dealii::numbers::invalid_unsigned_int)
976 _aux_functionals[_functional_for_ee_num]->SetProblemType(type);
983 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
984 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
985 template<
typename DATACONTAINER>
988 const DATACONTAINER& edc)
991 if (this->GetType() ==
"cost_functional")
995 else if (this->GetType() ==
"aux_functional")
998 return _aux_functionals[this->GetTypeNum()]->ElementValue(edc);
1000 else if (this->GetType() ==
"error_evaluation")
1002 return _aux_functionals[_functional_for_ee_num]->ElementValue(edc);
1007 "PDEProblemContainer::ElementFunctional");
1012 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1013 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1016 const std::map<std::string,
const dealii::Vector<double>*> ¶m_values,
1017 const std::map<std::string, const VECTOR*> &domain_values)
1019 if (this->GetType() ==
"cost_functional")
1023 else if (this->GetType() ==
"aux_functional")
1026 return _aux_functionals[this->GetTypeNum()]->PointValue(
1027 this->GetSpaceTimeHandler()->GetStateDoFHandler(),
1028 this->GetSpaceTimeHandler()->GetStateDoFHandler(), param_values,
1032 else if (this->GetType() ==
"error_evaluation")
1034 return _aux_functionals[_functional_for_ee_num]->PointValue(
1035 this->GetSpaceTimeHandler()->GetStateDoFHandler(),
1036 this->GetSpaceTimeHandler()->GetStateDoFHandler(), param_values,
1042 "PDEProblemContainer::PointFunctional");
1048 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1049 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1050 template<
typename FACEDATACONTAINER>
1053 const FACEDATACONTAINER& fdc)
1055 if (this->GetType() ==
"cost_functional")
1060 else if (this->GetType() ==
"aux_functional")
1063 return _aux_functionals[this->GetTypeNum()]->BoundaryValue(fdc);
1065 else if (this->GetType() ==
"error_evaluation")
1068 return _aux_functionals[_functional_for_ee_num]->BoundaryValue(fdc);
1073 "PDEProblemContainer::BoundaryFunctional");
1079 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1080 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1081 template<
typename FACEDATACONTAINER>
1084 const FACEDATACONTAINER& fdc)
1086 if (this->GetType() ==
"cost_functional")
1091 else if (this->GetType() ==
"aux_functional")
1094 return _aux_functionals[this->GetTypeNum()]->FaceValue(fdc);
1096 else if (this->GetType() ==
"error_evaluation")
1099 return _aux_functionals[_functional_for_ee_num]->FaceValue(fdc);
1104 "PDEProblemContainer::FaceFunctional");
1110 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1111 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1114 const std::map<std::string,
const dealii::Vector<double>*> ¶m_values,
1115 const std::map<std::string, const VECTOR*> &domain_values)
1117 if (this->GetType() ==
"cost_functional")
1122 else if (this->GetType() ==
"aux_functional")
1125 return _aux_functionals[this->GetTypeNum()]->AlgebraicValue(
1126 param_values, domain_values);
1128 else if (this->GetType() ==
"error_evaluation")
1131 return _aux_functionals[_functional_for_ee_num]->AlgebraicValue(
1132 param_values, domain_values);
1137 "PDEProblemContainer::AlgebraicFunctional");
1143 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1144 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1145 template<
typename DATACONTAINER>
1148 const DATACONTAINER& edc, dealii::Vector<double> &local_vector,
1149 double scale,
double scale_ico)
1151 if (this->GetType() ==
"state")
1154 this->GetPDE().ElementEquation(edc, local_vector, scale, scale_ico);
1156 else if (this->GetType() ==
"adjoint_for_ee")
1159 this->GetPDE().ElementEquation_U(edc, local_vector, scale,
1165 "PDEProblemContainer::ElementEquation");
1171 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1172 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1173 template<
typename DATACONTAINER>
1176 const DATACONTAINER& edc, dealii::Vector<double> &local_vector,
1180 if (this->GetType() ==
"state")
1183 this->GetPDE().ElementTimeEquation(edc, local_vector, scale);
1185 else if (this->GetType() ==
"adjoint_for_ee")
1188 "OptProblem::ElementTimeEquation");
1193 "PDEProblemContainer::ElementTimeEquation");
1199 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1200 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1201 template<
typename DATACONTAINER>
1204 const DATACONTAINER& edc, dealii::Vector<double> &local_vector,
1208 if (this->GetType() ==
"state")
1211 this->GetPDE().ElementTimeEquationExplicit(edc, local_vector,
1214 else if (this->GetType() ==
"adjoint_for_ee")
1217 "OptProblem::ElementTimeEquation");
1222 "PDEProblemContainer::ElementTimeEquation");
1228 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1229 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1230 template<
typename FACEDATACONTAINER>
1233 const FACEDATACONTAINER& fdc,
1234 dealii::Vector<double> &local_vector,
double scale,
1237 if (this->GetType() ==
"state")
1240 this->GetPDE().FaceEquation(fdc, local_vector, scale, scale_ico);
1242 else if (this->GetType() ==
"adjoint_for_ee")
1245 this->GetPDE().FaceEquation_U(fdc, local_vector, scale,
1251 "PDEProblemContainer::FaceEquation");
1257 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1258 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1259 template<
typename FACEDATACONTAINER>
1262 const FACEDATACONTAINER& fdc,
1263 dealii::Vector<double> &local_vector,
double scale,
1266 if (this->GetType() ==
"state")
1269 this->GetPDE().InterfaceEquation(fdc, local_vector, scale,
1272 else if (this->GetType() ==
"adjoint_for_ee")
1275 this->GetPDE().InterfaceEquation_U(fdc, local_vector, scale,
1281 "PDEProblemContainer::InterfaceEquation");
1286 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1287 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1288 template<
typename FACEDATACONTAINER>
1291 const FACEDATACONTAINER& fdc,
1292 dealii::Vector<double> &local_vector,
double scale,
1295 if (this->GetType() ==
"state")
1298 this->GetPDE().BoundaryEquation(fdc, local_vector, scale,
1301 else if (this->GetType() ==
"adjoint_for_ee")
1304 this->GetPDE().BoundaryEquation_U(fdc, local_vector, scale,
1310 "PDEProblemContainer::ElementBoundaryEquation");
1316 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1317 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1318 template<
typename DATACONTAINER>
1321 const DATACONTAINER& edc, dealii::Vector<double> &local_vector,
1325 if (this->GetType() ==
"state")
1328 this->GetPDE().ElementRightHandSide(edc, local_vector, scale);
1330 else if (this->GetType() ==
"adjoint_for_ee")
1335 if (_aux_functionals[_functional_for_ee_num]->GetType().find(
"domain")
1336 != std::string::npos)
1337 _aux_functionals[_functional_for_ee_num]->ElementValue_U(edc,
1338 local_vector, scale);
1343 "PDEProblemContainer::ElementRhs");
1349 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1350 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1353 const std::map<std::string,
const dealii::Vector<double>*> ¶m_values,
1354 const std::map<std::string, const VECTOR*> &domain_values,
1355 VECTOR& rhs_vector,
double scale)
1357 if (this->GetType() ==
"adjoint_for_ee")
1360 if (_aux_functionals[_functional_for_ee_num]->GetType().find(
"point")
1361 != std::string::npos)
1362 _aux_functionals[_functional_for_ee_num]->PointValue_U(
1363 this->GetSpaceTimeHandler()->GetStateDoFHandler(),
1364 this->GetSpaceTimeHandler()->GetStateDoFHandler(), param_values,
1365 domain_values, rhs_vector, scale);
1369 throw DOpEException(
"Not implemented",
"OptProblem::ElementRhs");
1375 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1376 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1377 template<
typename FACEDATACONTAINER>
1380 const FACEDATACONTAINER& fdc,
1381 dealii::Vector<double> &local_vector,
double scale)
1383 if (this->GetType() ==
"state")
1386 this->GetPDE().FaceRightHandSide(fdc, local_vector, scale);
1388 else if (this->GetType() ==
"adjoint_for_ee")
1391 if (_aux_functionals[_functional_for_ee_num]->GetType().find(
"face")
1392 != std::string::npos)
1393 _aux_functionals[_functional_for_ee_num]->FaceValue_U(fdc,
1394 local_vector, scale);
1399 "PDEProblemContainer::ElementFaceRhs");
1405 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1406 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1407 template<
typename FACEDATACONTAINER>
1410 const FACEDATACONTAINER& fdc,
1411 dealii::Vector<double> &local_vector,
double scale)
1413 if (this->GetType() ==
"state")
1416 this->GetPDE().BoundaryRightHandSide(fdc, local_vector, scale);
1418 else if (this->GetType() ==
"adjoint_for_ee")
1421 if (_aux_functionals[_functional_for_ee_num]->GetType().find(
1422 "boundary") != std::string::npos)
1423 _aux_functionals[_functional_for_ee_num]->BoundaryValue_U(fdc,
1424 local_vector, scale);
1429 "PDEProblemContainer::ElementBoundaryRhs");
1435 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1436 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1437 template<
typename DATACONTAINER>
1440 const DATACONTAINER& edc,
1441 dealii::FullMatrix<double> &local_entry_matrix,
double scale,
1445 if (this->GetType() ==
"state")
1448 this->GetPDE().ElementMatrix(edc, local_entry_matrix, scale, scale_ico);
1450 else if (this->GetType() ==
"adjoint_for_ee")
1453 this->GetPDE().ElementMatrix_T(edc, local_entry_matrix, scale,
1459 "PDEProblemContainer::ElementMatrix");
1466 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1467 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1468 template<
typename DATACONTAINER>
1471 const DATACONTAINER& edc, FullMatrix<double> &local_entry_matrix)
1474 if (this->GetType() ==
"state")
1477 this->GetPDE().ElementTimeMatrix(edc, local_entry_matrix);
1479 else if (this->GetType() ==
"dual_for_ee")
1482 "PDEProblemContainer::ElementTimeMatrix");
1487 "PDEProblemContainer::ElementTimeMatrix");
1494 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1495 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1496 template<
typename DATACONTAINER>
1499 const DATACONTAINER& edc,
1500 dealii::FullMatrix<double> &local_entry_matrix)
1503 if (this->GetType() ==
"state")
1506 this->GetPDE().ElementTimeMatrixExplicit(edc, local_entry_matrix);
1508 else if (this->GetType() ==
"dual_for_ee")
1511 "PDEProblemContainer::ElementTimeMatrix");
1516 "PDEProblemContainer::ElementTimeMatrix");
1523 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1524 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1525 template<
typename FACEDATACONTAINER>
1528 const FACEDATACONTAINER& fdc, FullMatrix<double> &local_entry_matrix,
1529 double scale,
double scale_ico)
1531 if (this->GetType() ==
"state")
1534 this->GetPDE().FaceMatrix(fdc, local_entry_matrix, scale, scale_ico);
1536 else if (this->GetType() ==
"adjoint_for_ee")
1538 this->GetPDE().FaceMatrix_T(fdc, local_entry_matrix, scale,
1544 "PDEProblemContainer::NewtonFaceMatrix");
1551 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1552 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1553 template<
typename FACEDATACONTAINER>
1556 const FACEDATACONTAINER& fdc, FullMatrix<double> &local_entry_matrix,
1557 double scale,
double scale_ico)
1559 if (this->GetType() ==
"state")
1562 this->GetPDE().InterfaceMatrix(fdc, local_entry_matrix, scale,
1565 else if (this->GetType() ==
"adjoint_for_ee")
1567 this->GetPDE().InterfaceMatrix_T(fdc, local_entry_matrix, scale,
1573 "PDEProblemContainer::NewtonInterfaceMatrix");
1579 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1580 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1581 template<
typename FACEDATACONTAINER>
1584 const FACEDATACONTAINER& fdc, FullMatrix<double> &local_matrix,
1585 double scale,
double scale_ico)
1587 if (this->GetType() ==
"state")
1590 this->GetPDE().BoundaryMatrix(fdc, local_matrix, scale,
1593 else if (this->GetType() ==
"adjoint_for_ee")
1596 this->GetPDE().BoundaryMatrix_T(fdc, local_matrix, scale,
1602 "PDEProblemContainer::ElementBoundaryMatrix");
1609 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1610 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1614 if (this->GetType() ==
"state" || this->GetType() ==
"adjoint_for_ee"
1615 || this->GetType() ==
"error_evaluation")
1622 "PDEProblemContainer::GetDoFType");
1628 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1629 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1630 const FE<dealdim, dealdim>&
1633 if ((this->GetType() ==
"state") || this->GetType() ==
"adjoint_for_ee")
1635 return this->GetSpaceTimeHandler()->GetFESystem(
"state");
1640 "PDEProblemContainer::GetFESystem");
1646 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1647 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1653 if (this->GetType().find(
"aux_functional") != std::string::npos)
1655 r = _aux_functionals[this->GetTypeNum()]->GetUpdateFlags();
1659 r = this->GetPDE().GetUpdateFlags();
1660 if (this->GetType() ==
"adjoint_for_ee"
1661 || this->GetType() ==
"error_evaluation")
1663 if (_functional_for_ee_num != dealii::numbers::invalid_unsigned_int)
1664 r = r | _aux_functionals[_functional_for_ee_num]->GetUpdateFlags();
1667 return r | update_JxW_values;
1672 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1673 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1678 if (this->GetType().find(
"aux_functional") != std::string::npos)
1680 r = _aux_functionals[this->GetTypeNum()]->GetFaceUpdateFlags();
1684 r = this->GetPDE().GetFaceUpdateFlags();
1685 if (this->GetType() ==
"adjoint_for_ee"
1686 || this->GetType() ==
"error_evaluation")
1688 if (_functional_for_ee_num != dealii::numbers::invalid_unsigned_int)
1691 | _aux_functionals[_functional_for_ee_num]->GetFaceUpdateFlags();
1694 return r | update_JxW_values;
1699 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1700 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1704 if (this->GetType() ==
"aux_functional")
1706 return _aux_functionals[this->GetTypeNum()]->GetType();
1708 else if (this->GetType() ==
"error_evaluation")
1710 return _aux_functionals[_functional_for_ee_num]->GetType();
1717 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1718 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1722 if (this->GetType() ==
"aux_functional")
1724 return _aux_functionals[this->GetTypeNum()]->GetName();
1726 else if (this->GetType() ==
"error_evaluation")
1728 return _aux_functionals[_functional_for_ee_num]->GetName();
1735 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1736 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1741 GetSpaceTimeHandler()->SetInterval(interval);
1744 for (
unsigned int i = 0; i < _primal_dirichlet_values.size(); i++)
1745 _primal_dirichlet_values[i]->SetTime(time);
1746 for (
unsigned int i = 0; i < _aux_functionals.size(); i++)
1747 _aux_functionals[i]->SetTime(time);
1749 this->GetPDE().SetTime(time);
1755 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1756 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1759 SPARSITYPATTERN & sparsity)
const
1761 if (this->GetType() ==
"state" || this->GetType() ==
"adjoint_for_ee")
1763 this->GetSpaceTimeHandler()->ComputeStateSparsityPattern(sparsity);
1768 "PDEProblemContainer::ComputeSparsityPattern");
1774 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1775 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1778 dealii::MGLevelObject<dealii::BlockSparsityPattern> & mg_sparsity_patterns,
1779 unsigned int n_levels)
const
1781 if (this->GetType() ==
"state")
1783 this->GetSpaceTimeHandler()->ComputeMGStateSparsityPattern(mg_sparsity_patterns, n_levels);
1788 "PDEProblemContainer::ComputeMGSparsityPattern");
1794 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1795 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1798 dealii::MGLevelObject<dealii::SparsityPattern> & mg_sparsity_patterns,
1799 unsigned int n_levels)
const
1801 if (this->GetType() ==
"state")
1803 this->GetSpaceTimeHandler()->ComputeMGStateSparsityPattern(mg_sparsity_patterns, n_levels);
1808 "PDEProblemContainer::ComputeMGSparsityPattern");
1817 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1818 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1822 if (this->GetType().find(
"aux_functional") != std::string::npos)
1824 return _aux_functionals[this->GetTypeNum()]->HasFaces();
1828 if ((this->GetType() ==
"state"))
1830 return this->GetPDE().HasFaces();
1832 else if (this->GetType() ==
"adjoint_for_ee")
1834 return this->GetPDE().HasFaces()
1835 || _aux_functionals[_functional_for_ee_num]->HasFaces();
1839 throw DOpEException(
"Unknown Type: '" + this->GetType() +
"'!",
1840 "PDEProblemContainer::HasFaces");
1847 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1848 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1852 if ((this->GetType() ==
"state") || this->GetType() ==
"aux_functional")
1857 else if (this->GetType() ==
"adjoint_for_ee")
1859 return _aux_functionals[_functional_for_ee_num]->HasPoints();
1863 throw DOpEException(
"Unknown Type: '" + this->GetType() +
"'!",
1864 "PDEProblemContainer::HasFaces");
1871 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1872 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1876 if (this->GetType().find(
"aux_functional") != std::string::npos)
1882 if ((this->GetType() ==
"state") || this->GetType() ==
"adjoint_for_ee")
1884 return this->GetPDE().HasInterfaces();
1888 throw DOpEException(
"Unknown Type: '" + this->GetType() +
"'!",
1889 "PDEProblemContainer::HasFaces");
1896 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1897 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1900 unsigned int color,
const std::vector<bool>& comp_mask,
1904 assert(values->n_components() == comp_mask.size());
1906 unsigned int comp = _dirichlet_colors.size();
1907 for (
unsigned int i = 0; i < _dirichlet_colors.size(); ++i)
1909 if (_dirichlet_colors[i] == color)
1915 if (comp != _dirichlet_colors.size())
1917 std::stringstream s;
1918 s <<
"DirichletColor" << color <<
" has multiple occurrences !";
1920 "PDEProblemContainer::SetDirichletBoundary");
1922 _dirichlet_colors.push_back(color);
1923 _dirichlet_comps.push_back(comp_mask);
1925 DD, VECTOR, dealdim>(*values);
1926 _primal_dirichlet_values.push_back(data);
1931 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1932 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1933 const std::vector<unsigned int>&
1936 if ((this->GetType() ==
"state") || this->GetType() ==
"adjoint_for_ee")
1938 return _dirichlet_colors;
1943 "PDEProblemContainer::GetDirichletColors");
1949 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1950 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1951 const std::vector<bool>&
1953 unsigned int color)
const
1955 if ((this->GetType() ==
"state" || this->GetType() ==
"adjoint_for_ee"))
1957 unsigned int comp = _dirichlet_colors.size();
1958 for (
unsigned int i = 0; i < _dirichlet_colors.size(); ++i)
1960 if (_dirichlet_colors[i] == color)
1966 if (comp == _dirichlet_colors.size())
1968 std::stringstream s;
1969 s <<
"DirichletColor" << color <<
" has not been found !";
1971 "PDEProblemContainer::GetDirichletCompMask");
1973 return _dirichlet_comps[comp];
1978 "PDEProblemContainer::GetDirichletCompMask");
1984 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1985 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1986 const Function<dealdim>&
1989 const std::map<std::string,
const dealii::Vector<double>*> ¶m_values,
1990 const std::map<std::string, const VECTOR*> &domain_values)
const
1993 unsigned int col = _dirichlet_colors.size();
1994 if ((this->GetType() ==
"state") || this->GetType() ==
"adjoint_for_ee")
1996 for (
unsigned int i = 0; i < _dirichlet_colors.size(); ++i)
1998 if (_dirichlet_colors[i] == color)
2004 if (col == _dirichlet_colors.size())
2006 std::stringstream s;
2007 s <<
"DirichletColor" << color <<
" has not been found !";
2009 "PDEProblemContainer::GetDirichletValues");
2015 "PDEProblemContainer::GetDirichletValues");
2018 if (this->GetType() ==
"state")
2020 _primal_dirichlet_values[col]->ReInit(param_values, domain_values,
2022 return *(_primal_dirichlet_values[col]);
2024 else if (this->GetType() ==
"adjoint_for_ee"
2025 || (this->GetType() ==
"adjoint_hessian"))
2027 return *(_zero_dirichlet_values);
2032 "PDEProblemContainer::GetDirichletValues");
2038 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
2039 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
2040 const std::vector<unsigned int>&
2043 if (this->GetType() ==
"state")
2045 return _state_boundary_equation_colors;
2047 else if (this->GetType() ==
"adjoint_for_ee")
2049 return _adjoint_boundary_equation_colors;
2054 "PDEProblemContainer::GetBoundaryEquationColors");
2060 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
2061 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
2067 unsigned int comp = _state_boundary_equation_colors.size();
2068 for (
unsigned int i = 0; i < _state_boundary_equation_colors.size();
2071 if (_state_boundary_equation_colors[i] == color)
2077 if (comp != _state_boundary_equation_colors.size())
2079 std::stringstream s;
2080 s <<
"Boundary Equation Color" << color
2081 <<
" has multiple occurences !";
2083 "PDEProblemContainer::SetBoundaryEquationColors");
2085 _state_boundary_equation_colors.push_back(color);
2088 unsigned int comp = _adjoint_boundary_equation_colors.size();
2089 for (
unsigned int i = 0; i < _adjoint_boundary_equation_colors.size();
2092 if (_adjoint_boundary_equation_colors[i] == color)
2098 if (comp != _adjoint_boundary_equation_colors.size())
2105 _adjoint_boundary_equation_colors.push_back(color);
2113 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
2114 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
2115 const std::vector<unsigned int>&
2119 if (this->GetType() ==
"cost_functional"
2120 || this->GetType() ==
"aux_functional"
2121 || this->GetType() ==
"error_evaluation")
2123 return _boundary_functional_colors;
2128 "PDEProblemContainer::GetBoundaryFunctionalColors");
2134 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
2135 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
2141 unsigned int comp = _boundary_functional_colors.size();
2142 for (
unsigned int i = 0; i < _boundary_functional_colors.size(); ++i)
2144 if (_boundary_functional_colors[i] == color)
2150 if (comp != _boundary_functional_colors.size())
2152 std::stringstream s;
2153 s <<
"Boundary Functional Color" << color
2154 <<
" has multiple occurences !";
2156 "PDEProblemContainer::SetBoundaryFunctionalColors");
2158 _boundary_functional_colors.push_back(color);
2161 unsigned int comp = _adjoint_boundary_equation_colors.size();
2162 for (
unsigned int i = 0; i < _adjoint_boundary_equation_colors.size();
2165 if (_adjoint_boundary_equation_colors[i] == color)
2171 if (comp != _adjoint_boundary_equation_colors.size())
2178 _adjoint_boundary_equation_colors.push_back(color);
2185 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
2186 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
2190 return this->GetPDE().GetStateNBlocks();
2195 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
2196 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
2200 if ((this->GetType() ==
"state") || (this->GetType() ==
"adjoint_for_ee"))
2202 return this->GetStateNBlocks();
2207 "PDEProblemContainer::GetNBlocks");
2213 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
2214 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
2217 unsigned int b)
const
2219 if ((this->GetType() ==
"state") || (this->GetType() ==
"adjoint_for_ee"))
2221 return GetSpaceTimeHandler()->GetStateDoFsPerBlock(b);
2226 "PDEProblemContainer::GetDoFsPerBlock");
2232 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
2233 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
2234 const std::vector<unsigned int>&
2237 if ((this->GetType() ==
"state") || (this->GetType() ==
"adjoint_for_ee"))
2239 return GetSpaceTimeHandler()->GetStateDoFsPerBlock();
2244 "PDEProblemContainer::GetDoFsPerBlock");
2250 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
2251 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
2252 const dealii::ConstraintMatrix&
2258 if ((this->GetType() ==
"state") || (this->GetType() ==
"adjoint_for_ee"))
2260 return GetSpaceTimeHandler()->GetStateDoFConstraints();
2265 "PDEProblemContainer::GetDoFConstraints");
2271 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
2272 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
2276 if (this->GetType() ==
"cost_functional")
2278 else if (this->GetType() ==
"aux_functional")
2279 return _aux_functionals[this->GetTypeNum()]->NeedTime();
2280 else if (this->GetType() ==
"error_evaluation")
2281 return _aux_functionals[_functional_for_ee_num]->NeedTime();
2284 "PDEProblemContainer::NeedTimeFunctional");
Definition: stateproblem.h:48
StateProblem< PDEProblemContainer< PDE, DD, SPARSITYPATTERN, VECTOR, dealdim, FE, DH >, PDE, DD, SPARSITYPATTERN, VECTOR, dealdim > & GetStateProblem()
Definition: pdeproblemcontainer.h:126
double BoundaryFunctional(const FACEDATACONTAINER &fdc)
Definition: pdeproblemcontainer.h:1052
void ElementMatrix(const DATACONTAINER &dc, dealii::FullMatrix< double > &local_entry_matrix, double scale=1., double scale_ico=1.)
Definition: pdeproblemcontainer.h:1439
void RegisterExceptionHandler(DOpEExceptionHandler< VECTOR > *OH)
Definition: pdeproblemcontainer.h:168
void SetInitialValues(const dealii::Function< dealdim > *values)
Definition: pdeproblemcontainer.h:609
bool HasInterfaces() const
Definition: pdeproblemcontainer.h:1874
PDEProblemContainer< PDE, DD, SPARSITYPATTERN, VECTOR, dealdim, FE, DH > & GetBaseProblem()
Definition: pdeproblemcontainer.h:141
void ReInit(std::string algo_type)
Definition: pdeproblemcontainer.h:930
unsigned int GetStateNBlocks()
Definition: pdeproblemcontainer.h:825
void RegisterOutputHandler(DOpEOutputHandler< VECTOR > *OH)
Definition: pdeproblemcontainer.h:160
const dealii::Function< dealdim > & GetDirichletValues(unsigned int color, const std::map< std::string, const dealii::Vector< double > * > ¶m_values, const std::map< std::string, const VECTOR * > &domain_values) const
Definition: pdeproblemcontainer.h:1987
void FaceRhs(const FACEDATACONTAINER &dc, dealii::Vector< double > &local_vector, double scale=1.)
Definition: pdeproblemcontainer.h:1379
dealii::UpdateFlags GetFaceUpdateFlags() const
Definition: pdeproblemcontainer.h:1675
virtual std::string GetName() const
Definition: pdeproblemcontainer.h:114
unsigned int GetStateNBlocks() const
Definition: pdeproblemcontainer.h:2188
std::string GetDoFType() const
Definition: pdeproblemcontainer.h:1612
void SetBoundaryFunctionalColors(unsigned int color)
Definition: pdeproblemcontainer.h:2137
Definition: elementdatacontainer.h:56
unsigned int GetNFunctionals() const
Definition: pdeproblemcontainer.h:698
void ComputeSparsityPattern(SPARSITYPATTERN &sparsity) const
Definition: pdeproblemcontainer.h:1758
void FaceEquation(const FACEDATACONTAINER &dc, dealii::Vector< double > &local_vector, double scale, double scale_ico)
Definition: pdeproblemcontainer.h:1232
Definition: problemcontainer_internal.h:38
Definition: optproblemcontainer.h:70
const std::vector< unsigned int > & GetDirichletColors() const
Definition: pdeproblemcontainer.h:1934
std::vector< unsigned int > & GetStateBlockComponent()
Definition: pdeproblemcontainer.h:833
void ComputeMGSparsityPattern(dealii::MGLevelObject< dealii::BlockSparsityPattern > &mg_sparsity_patterns, unsigned int n_levels) const
Definition: pdeproblemcontainer.h:1777
const dealii::ConstraintMatrix & GetDoFConstraints() const
Definition: pdeproblemcontainer.h:2253
void AddFunctional(FunctionalInterface< ElementDataContainer, FaceDataContainer, DH, VECTOR, dealdim > *F)
Definition: pdeproblemcontainer.h:641
std::string GetFunctionalType() const
Definition: pdeproblemcontainer.h:1702
void SetTime(double time, const TimeIterator &interval)
Definition: pdeproblemcontainer.h:1738
const StateSpaceTimeHandler< FE, DH, SPARSITYPATTERN, VECTOR, dealdim > * GetSpaceTimeHandler() const
Definition: pdeproblemcontainer.h:768
Definition: timeiterator.h:63
bool NeedTimeFunctional() const
Definition: pdeproblemcontainer.h:2274
const PDE & GetPDE() const
Definition: problemcontainer_internal.h:103
unsigned int GetNBlocks() const
Definition: pdeproblemcontainer.h:2198
const dealii::Function< dealdim > & GetInitialValues() const
Definition: pdeproblemcontainer.h:614
const std::vector< unsigned int > & GetDoFsPerBlock() const
Definition: pdeproblemcontainer.h:2235
StateSpaceTimeHandler< FE, DH, SPARSITYPATTERN, VECTOR, dealdim > * GetSpaceTimeHandler()
Definition: pdeproblemcontainer.h:776
Definition: facedatacontainer.h:48
void ElementTimeEquationExplicit(const DATACONTAINER &dc, dealii::Vector< double > &local_vector, double scale=1.)
Definition: pdeproblemcontainer.h:1203
void SetBoundaryEquationColors(unsigned int color)
Definition: pdeproblemcontainer.h:2063
dealii::UpdateFlags GetUpdateFlags() const
Definition: pdeproblemcontainer.h:1649
bool HasPoints() const
Definition: pdeproblemcontainer.h:1850
const std::map< std::string, unsigned int > & GetFunctionalPosition() const
Definition: pdeproblemcontainer.h:847
std::string GetFunctionalName() const
Definition: pdeproblemcontainer.h:1720
Definition: statespacetimehandler.h:59
double ElementFunctional(const DATACONTAINER &edc)
Definition: pdeproblemcontainer.h:987
const std::vector< unsigned int > & GetBoundaryEquationColors() const
Definition: pdeproblemcontainer.h:2041
Definition: primaldirichletdata.h:41
void ElementRhs(const DATACONTAINER &dc, dealii::Vector< double > &local_vector, double scale=1.)
Definition: pdeproblemcontainer.h:1320
void ElementEquation(const DATACONTAINER &edc, dealii::Vector< double > &local_vector, double scale, double scale_ico)
Definition: pdeproblemcontainer.h:1147
Definition: functionalinterface.h:53
void SetDirichletBoundaryColors(unsigned int color, const std::vector< bool > &comp_mask, const DD *values)
Definition: pdeproblemcontainer.h:1899
void SetFunctionalForErrorEstimation(std::string functional_name)
Definition: pdeproblemcontainer.h:670
void BoundaryRhs(const FACEDATACONTAINER &dc, dealii::Vector< double > &local_vector, double scale=1.)
Definition: pdeproblemcontainer.h:1409
bool HasFaces() const
Definition: pdeproblemcontainer.h:1820
void AddAuxiliaryToIntegrator(INTEGRATOR &)
Definition: pdeproblemcontainer.h:808
void ElementTimeMatrixExplicit(const DATACONTAINER &dc, dealii::FullMatrix< double > &local_entry_matrix)
Definition: pdeproblemcontainer.h:1498
DOpEExceptionHandler< VECTOR > * GetExceptionHandler()
Definition: pdeproblemcontainer.h:745
PDEProblemContainer(PDE &pde, StateSpaceTimeHandler< FE, DH, SPARSITYPATTERN, VECTOR, dealdim > &STH)
Definition: pdeproblemcontainer.h:891
double FaceFunctional(const FACEDATACONTAINER &fdc)
Definition: pdeproblemcontainer.h:1083
void InterfaceMatrix(const FACEDATACONTAINER &dc, dealii::FullMatrix< double > &local_entry_matrix, double scale=1., double scale_ico=1.)
void ElementTimeMatrix(const DATACONTAINER &dc, dealii::FullMatrix< double > &local_entry_matrix)
virtual ~PDEProblemContainer()
Definition: pdeproblemcontainer.h:908
void BoundaryMatrix(const FACEDATACONTAINER &dc, dealii::FullMatrix< double > &local_matrix, double scale=1., double scale_ico=1.)
const FE< dealdim, dealdim > & GetFESystem() const
Definition: pdeproblemcontainer.h:1631
void BoundaryEquation(const FACEDATACONTAINER &dc, dealii::Vector< double > &local_vector, double scale, double scale_ico)
Definition: pdeproblemcontainer.h:1290
void ElementTimeEquation(const DATACONTAINER &dc, dealii::Vector< double > &local_vector, double scale=1.)
Definition: pdeproblemcontainer.h:1175
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=1.)
Definition: pdeproblemcontainer.h:1352
Definition: pdeproblemcontainer.h:100
const std::vector< bool > & GetDirichletCompMask(unsigned int color) const
Definition: pdeproblemcontainer.h:1952
double AlgebraicFunctional(const std::map< std::string, const dealii::Vector< double > * > &values, const std::map< std::string, const VECTOR * > &block_values)
Definition: pdeproblemcontainer.h:1113
Definition: dopeexception.h:35
void InterfaceEquation(const FACEDATACONTAINER &dc, dealii::Vector< double > &local_vector, double scale, double scale_ico)
Definition: pdeproblemcontainer.h:1261
void DeleteAuxiliaryFromIntegrator(INTEGRATOR &)
Definition: pdeproblemcontainer.h:817
void SetType(std::string type, unsigned int num=0)
Definition: pdeproblemcontainer.h:967
const std::vector< unsigned int > & GetBoundaryFunctionalColors() const
Definition: pdeproblemcontainer.h:2116
DOpEOutputHandler< VECTOR > * GetOutputHandler()
Definition: pdeproblemcontainer.h:754
double PointFunctional(const std::map< std::string, const dealii::Vector< double > * > ¶m_values, const std::map< std::string, const VECTOR * > &domain_values)
Definition: pdeproblemcontainer.h:1015
Definition: optproblemcontainer.h:72
void FaceMatrix(const FACEDATACONTAINER &dc, dealii::FullMatrix< double > &local_entry_matrix, double scale=1., double scale_ico=1.)