24 #ifndef PDEProblemContainer_H_
25 #define PDEProblemContainer_H_
43 #include <deal.II/lac/vector.h>
44 #include <deal.II/lac/full_matrix.h>
45 #include <deal.II/grid/tria_iterator.h>
46 #include <deal.II/dofs/dof_handler.h>
47 #include <deal.II/dofs/dof_accessor.h>
48 #include <deal.II/dofs/dof_tools.h>
49 #include <deal.II/fe/fe_system.h>
50 #include <deal.II/fe/fe_values.h>
51 #include <deal.II/fe/mapping.h>
52 #include <deal.II/base/quadrature_lib.h>
53 #include <deal.II/lac/block_sparsity_pattern.h>
54 #include <deal.II/lac/compressed_simple_sparsity_pattern.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_;
764 unsigned int time_dof_number,
808 template<
typename INTEGRATOR>
817 template<
typename INTEGRATOR>
829 return this->
GetPDE().GetStateNBlocks();
834 std::vector<unsigned int>&
837 return this->
GetPDE().GetStateBlockComponent();
848 const std::map<std::string, unsigned int>&
851 return functional_position_;
854 template<
typename ELEMENTITERATOR>
855 bool AtInterface(ELEMENTITERATOR& element,
unsigned int face)
const
857 return this->
GetPDE().AtInterface(element,face);
865 std::string algo_type_;
869 VECTOR, dealdim>*> aux_functionals_;
870 std::map<std::string, unsigned int> functional_position_;
872 unsigned int functional_for_ee_num_;
875 std::vector<unsigned int> dirichlet_colors_;
876 std::vector<std::vector<bool> > dirichlet_comps_;
877 std::vector<PrimalDirichletData<DD, VECTOR, dealdim>*> primal_dirichlet_values_;
878 const dealii::Function<dealdim>* zero_dirichlet_values_;
880 const dealii::Function<dealdim>* initial_values_;
882 std::vector<unsigned int> state_boundary_equation_colors_;
883 std::vector<unsigned int> adjoint_boundary_equation_colors_;
885 std::vector<unsigned int> boundary_functional_colors_;
889 DH>, PDE, DD, SPARSITYPATTERN, VECTOR, dealdim>* state_problem_;
893 DH>, PDE, DD, SPARSITYPATTERN, VECTOR, dealdim> ;
897 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
898 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
904 ExceptionHandler_ = NULL;
905 OutputHandler_ = NULL;
906 zero_dirichlet_values_ =
new ZeroFunction<dealdim>(
907 this->
GetPDE().GetStateNComponents());
909 functional_for_ee_num_ = dealii::numbers::invalid_unsigned_int;
914 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
915 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
918 if (zero_dirichlet_values_ != NULL)
920 delete zero_dirichlet_values_;
922 for (
unsigned int i = 0; i < primal_dirichlet_values_.size(); i++)
924 if (primal_dirichlet_values_[i] != NULL)
925 delete primal_dirichlet_values_[i];
927 if (state_problem_ != NULL)
929 delete state_problem_;
935 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
936 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
939 std::string algo_type)
941 if (state_problem_ != NULL)
943 delete state_problem_;
944 state_problem_ = NULL;
947 if (algo_type_ != algo_type && algo_type_ !=
"")
950 "PDEProblemContainer::ReInit");
954 algo_type_ = algo_type;
955 this->SetTypeInternal(
"");
957 if (algo_type_ ==
"reduced")
959 GetSpaceTimeHandler()->ReInit(this->GetPDE().GetStateNBlocks(),
960 this->GetPDE().GetStateBlockComponent());
965 "PDEProblemContainer::ReInit");
972 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
973 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
976 std::string type,
unsigned int num)
978 if (this->GetType() != type || this->GetTypeNum() != num)
980 this->SetTypeInternal(type);
981 this->SetTypeNumInternal(num);
982 this->GetPDE().SetProblemType(type);
983 if (functional_for_ee_num_ != dealii::numbers::invalid_unsigned_int)
984 aux_functionals_[functional_for_ee_num_]->SetProblemType(type);
991 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
992 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
993 template<
typename DATACONTAINER>
996 const DATACONTAINER& edc)
999 if (this->GetType() ==
"cost_functional")
1003 else if (this->GetType() ==
"aux_functional")
1006 return aux_functionals_[this->GetTypeNum()]->ElementValue(edc);
1008 else if (this->GetType() ==
"error_evaluation")
1010 return aux_functionals_[functional_for_ee_num_]->ElementValue(edc);
1015 "PDEProblemContainer::ElementFunctional");
1020 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1021 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1024 const std::map<std::string,
const dealii::Vector<double>*> ¶m_values,
1025 const std::map<std::string, const VECTOR*> &domain_values)
1027 if (this->GetType() ==
"cost_functional")
1031 else if (this->GetType() ==
"aux_functional")
1034 return aux_functionals_[this->GetTypeNum()]->PointValue(
1035 this->GetSpaceTimeHandler()->GetStateDoFHandler(),
1036 this->GetSpaceTimeHandler()->GetStateDoFHandler(), param_values,
1040 else if (this->GetType() ==
"error_evaluation")
1042 return aux_functionals_[functional_for_ee_num_]->PointValue(
1043 this->GetSpaceTimeHandler()->GetStateDoFHandler(),
1044 this->GetSpaceTimeHandler()->GetStateDoFHandler(), param_values,
1050 "PDEProblemContainer::PointFunctional");
1056 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1057 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1058 template<
typename FACEDATACONTAINER>
1061 const FACEDATACONTAINER& fdc)
1063 if (this->GetType() ==
"cost_functional")
1068 else if (this->GetType() ==
"aux_functional")
1071 return aux_functionals_[this->GetTypeNum()]->BoundaryValue(fdc);
1073 else if (this->GetType() ==
"error_evaluation")
1076 return aux_functionals_[functional_for_ee_num_]->BoundaryValue(fdc);
1081 "PDEProblemContainer::BoundaryFunctional");
1087 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1088 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1089 template<
typename FACEDATACONTAINER>
1092 const FACEDATACONTAINER& fdc)
1094 if (this->GetType() ==
"cost_functional")
1099 else if (this->GetType() ==
"aux_functional")
1102 return aux_functionals_[this->GetTypeNum()]->FaceValue(fdc);
1104 else if (this->GetType() ==
"error_evaluation")
1107 return aux_functionals_[functional_for_ee_num_]->FaceValue(fdc);
1112 "PDEProblemContainer::FaceFunctional");
1118 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1119 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1122 const std::map<std::string,
const dealii::Vector<double>*> ¶m_values,
1123 const std::map<std::string, const VECTOR*> &domain_values)
1125 if (this->GetType() ==
"cost_functional")
1130 else if (this->GetType() ==
"aux_functional")
1133 return aux_functionals_[this->GetTypeNum()]->AlgebraicValue(
1134 param_values, domain_values);
1136 else if (this->GetType() ==
"error_evaluation")
1139 return aux_functionals_[functional_for_ee_num_]->AlgebraicValue(
1140 param_values, domain_values);
1145 "PDEProblemContainer::AlgebraicFunctional");
1151 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1152 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1153 template<
typename DATACONTAINER>
1156 const DATACONTAINER& edc, dealii::Vector<double> &local_vector,
1157 double scale,
double scale_ico)
1159 if (this->GetType() ==
"state")
1162 this->GetPDE().ElementEquation(edc, local_vector, scale, scale_ico);
1164 else if (this->GetType() ==
"adjoint_for_ee")
1167 this->GetPDE().ElementEquation_U(edc, local_vector, scale,
1173 "PDEProblemContainer::ElementEquation");
1179 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1180 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1181 template<
typename DATACONTAINER>
1184 const DATACONTAINER& edc, dealii::Vector<double> &local_vector,
1188 if (this->GetType() ==
"state")
1191 this->GetPDE().ElementTimeEquation(edc, local_vector, scale);
1193 else if (this->GetType() ==
"adjoint_for_ee")
1196 "OptProblem::ElementTimeEquation");
1201 "PDEProblemContainer::ElementTimeEquation");
1207 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1208 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1209 template<
typename DATACONTAINER>
1212 const DATACONTAINER& edc, dealii::Vector<double> &local_vector,
1216 if (this->GetType() ==
"state")
1219 this->GetPDE().ElementTimeEquationExplicit(edc, local_vector,
1222 else if (this->GetType() ==
"adjoint_for_ee")
1225 "OptProblem::ElementTimeEquation");
1230 "PDEProblemContainer::ElementTimeEquation");
1236 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1237 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1238 template<
typename FACEDATACONTAINER>
1241 const FACEDATACONTAINER& fdc,
1242 dealii::Vector<double> &local_vector,
double scale,
1245 if (this->GetType() ==
"state")
1248 this->GetPDE().FaceEquation(fdc, local_vector, scale, scale_ico);
1250 else if (this->GetType() ==
"adjoint_for_ee")
1253 this->GetPDE().FaceEquation_U(fdc, local_vector, scale,
1259 "PDEProblemContainer::FaceEquation");
1265 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1266 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1267 template<
typename FACEDATACONTAINER>
1270 const FACEDATACONTAINER& fdc,
1271 dealii::Vector<double> &local_vector,
double scale,
1274 if (this->GetType() ==
"state")
1277 this->GetPDE().InterfaceEquation(fdc, local_vector, scale,
1280 else if (this->GetType() ==
"adjoint_for_ee")
1283 this->GetPDE().InterfaceEquation_U(fdc, local_vector, scale,
1289 "PDEProblemContainer::InterfaceEquation");
1294 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1295 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1296 template<
typename FACEDATACONTAINER>
1299 const FACEDATACONTAINER& fdc,
1300 dealii::Vector<double> &local_vector,
double scale,
1303 if (this->GetType() ==
"state")
1306 this->GetPDE().BoundaryEquation(fdc, local_vector, scale,
1309 else if (this->GetType() ==
"adjoint_for_ee")
1312 this->GetPDE().BoundaryEquation_U(fdc, local_vector, scale,
1318 "PDEProblemContainer::ElementBoundaryEquation");
1324 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1325 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1326 template<
typename DATACONTAINER>
1329 const DATACONTAINER& edc, dealii::Vector<double> &local_vector,
1333 if (this->GetType() ==
"state")
1336 this->GetPDE().ElementRightHandSide(edc, local_vector, scale);
1338 else if (this->GetType() ==
"adjoint_for_ee")
1343 if (aux_functionals_[functional_for_ee_num_]->GetType().find(
"domain")
1344 != std::string::npos)
1345 aux_functionals_[functional_for_ee_num_]->ElementValue_U(edc,
1346 local_vector, scale);
1351 "PDEProblemContainer::ElementRhs");
1357 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1358 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1361 const std::map<std::string,
const dealii::Vector<double>*> ¶m_values,
1362 const std::map<std::string, const VECTOR*> &domain_values,
1363 VECTOR& rhs_vector,
double scale)
1365 if (this->GetType() ==
"adjoint_for_ee")
1368 if (aux_functionals_[functional_for_ee_num_]->GetType().find(
"point")
1369 != std::string::npos)
1370 aux_functionals_[functional_for_ee_num_]->PointValue_U(
1371 this->GetSpaceTimeHandler()->GetStateDoFHandler(),
1372 this->GetSpaceTimeHandler()->GetStateDoFHandler(), param_values,
1373 domain_values, rhs_vector, scale);
1377 throw DOpEException(
"Not implemented",
"OptProblem::ElementRhs");
1383 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1384 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1385 template<
typename FACEDATACONTAINER>
1388 const FACEDATACONTAINER& fdc,
1389 dealii::Vector<double> &local_vector,
double scale)
1391 if (this->GetType() ==
"state")
1394 this->GetPDE().FaceRightHandSide(fdc, local_vector, scale);
1396 else if (this->GetType() ==
"adjoint_for_ee")
1399 if (aux_functionals_[functional_for_ee_num_]->GetType().find(
"face")
1400 != std::string::npos)
1401 aux_functionals_[functional_for_ee_num_]->FaceValue_U(fdc,
1402 local_vector, scale);
1407 "PDEProblemContainer::ElementFaceRhs");
1413 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1414 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1415 template<
typename FACEDATACONTAINER>
1418 const FACEDATACONTAINER& fdc,
1419 dealii::Vector<double> &local_vector,
double scale)
1421 if (this->GetType() ==
"state")
1424 this->GetPDE().BoundaryRightHandSide(fdc, local_vector, scale);
1426 else if (this->GetType() ==
"adjoint_for_ee")
1429 if (aux_functionals_[functional_for_ee_num_]->GetType().find(
1430 "boundary") != std::string::npos)
1431 aux_functionals_[functional_for_ee_num_]->BoundaryValue_U(fdc,
1432 local_vector, scale);
1437 "PDEProblemContainer::ElementBoundaryRhs");
1443 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1444 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1445 template<
typename DATACONTAINER>
1448 const DATACONTAINER& edc,
1449 dealii::FullMatrix<double> &local_entry_matrix,
double scale,
1453 if (this->GetType() ==
"state")
1456 this->GetPDE().ElementMatrix(edc, local_entry_matrix, scale, scale_ico);
1458 else if (this->GetType() ==
"adjoint_for_ee")
1461 this->GetPDE().ElementMatrix_T(edc, local_entry_matrix, scale,
1467 "PDEProblemContainer::ElementMatrix");
1474 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1475 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1476 template<
typename DATACONTAINER>
1479 const DATACONTAINER& edc, FullMatrix<double> &local_entry_matrix)
1482 if (this->GetType() ==
"state")
1485 this->GetPDE().ElementTimeMatrix(edc, local_entry_matrix);
1487 else if (this->GetType() ==
"dual_for_ee")
1490 "PDEProblemContainer::ElementTimeMatrix");
1495 "PDEProblemContainer::ElementTimeMatrix");
1502 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1503 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1504 template<
typename DATACONTAINER>
1507 const DATACONTAINER& edc,
1508 dealii::FullMatrix<double> &local_entry_matrix)
1511 if (this->GetType() ==
"state")
1514 this->GetPDE().ElementTimeMatrixExplicit(edc, local_entry_matrix);
1516 else if (this->GetType() ==
"dual_for_ee")
1519 "PDEProblemContainer::ElementTimeMatrix");
1524 "PDEProblemContainer::ElementTimeMatrix");
1531 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1532 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1533 template<
typename FACEDATACONTAINER>
1536 const FACEDATACONTAINER& fdc, FullMatrix<double> &local_entry_matrix,
1537 double scale,
double scale_ico)
1539 if (this->GetType() ==
"state")
1542 this->GetPDE().FaceMatrix(fdc, local_entry_matrix, scale, scale_ico);
1544 else if (this->GetType() ==
"adjoint_for_ee")
1546 this->GetPDE().FaceMatrix_T(fdc, local_entry_matrix, scale,
1552 "PDEProblemContainer::NewtonFaceMatrix");
1559 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1560 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1561 template<
typename FACEDATACONTAINER>
1564 const FACEDATACONTAINER& fdc, FullMatrix<double> &local_entry_matrix,
1565 double scale,
double scale_ico)
1567 if (this->GetType() ==
"state")
1570 this->GetPDE().InterfaceMatrix(fdc, local_entry_matrix, scale,
1573 else if (this->GetType() ==
"adjoint_for_ee")
1575 this->GetPDE().InterfaceMatrix_T(fdc, local_entry_matrix, scale,
1581 "PDEProblemContainer::NewtonInterfaceMatrix");
1587 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1588 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1589 template<
typename FACEDATACONTAINER>
1592 const FACEDATACONTAINER& fdc, FullMatrix<double> &local_matrix,
1593 double scale,
double scale_ico)
1595 if (this->GetType() ==
"state")
1598 this->GetPDE().BoundaryMatrix(fdc, local_matrix, scale,
1601 else if (this->GetType() ==
"adjoint_for_ee")
1604 this->GetPDE().BoundaryMatrix_T(fdc, local_matrix, scale,
1610 "PDEProblemContainer::ElementBoundaryMatrix");
1617 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1618 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1622 if (this->GetType() ==
"state" || this->GetType() ==
"adjoint_for_ee"
1623 || this->GetType() ==
"error_evaluation")
1630 "PDEProblemContainer::GetDoFType");
1636 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1637 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1638 const FE<dealdim, dealdim>&
1641 if ((this->GetType() ==
"state") || this->GetType() ==
"adjoint_for_ee")
1643 return this->GetSpaceTimeHandler()->GetFESystem(
"state");
1648 "PDEProblemContainer::GetFESystem");
1654 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1655 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1661 if (this->GetType().find(
"aux_functional") != std::string::npos)
1663 r = aux_functionals_[this->GetTypeNum()]->GetUpdateFlags();
1667 r = this->GetPDE().GetUpdateFlags();
1668 if (this->GetType() ==
"adjoint_for_ee"
1669 || this->GetType() ==
"error_evaluation")
1671 if (functional_for_ee_num_ != dealii::numbers::invalid_unsigned_int)
1672 r = r | aux_functionals_[functional_for_ee_num_]->GetUpdateFlags();
1675 return r | update_JxW_values;
1680 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1681 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1686 if (this->GetType().find(
"aux_functional") != std::string::npos)
1688 r = aux_functionals_[this->GetTypeNum()]->GetFaceUpdateFlags();
1692 r = this->GetPDE().GetFaceUpdateFlags();
1693 if (this->GetType() ==
"adjoint_for_ee"
1694 || this->GetType() ==
"error_evaluation")
1696 if (functional_for_ee_num_ != dealii::numbers::invalid_unsigned_int)
1699 | aux_functionals_[functional_for_ee_num_]->GetFaceUpdateFlags();
1702 return r | update_JxW_values;
1707 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1708 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1712 if (this->GetType() ==
"aux_functional")
1714 return aux_functionals_[this->GetTypeNum()]->GetType();
1716 else if (this->GetType() ==
"error_evaluation")
1718 return aux_functionals_[functional_for_ee_num_]->GetType();
1725 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1726 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1730 if (this->GetType() ==
"aux_functional")
1732 return aux_functionals_[this->GetTypeNum()]->GetName();
1734 else if (this->GetType() ==
"error_evaluation")
1736 return aux_functionals_[functional_for_ee_num_]->GetName();
1743 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1744 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1748 unsigned int time_dof_number,
1751 GetSpaceTimeHandler()->SetInterval(interval);
1754 for (
unsigned int i = 0; i < primal_dirichlet_values_.size(); i++)
1755 primal_dirichlet_values_[i]->SetTime(time);
1756 for (
unsigned int i = 0; i < aux_functionals_.size(); i++)
1757 aux_functionals_[i]->SetTime(time);
1759 this->GetPDE().SetTime(time);
1765 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1766 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1769 SPARSITYPATTERN & sparsity)
const
1771 if (this->GetType() ==
"state" || this->GetType() ==
"adjoint_for_ee")
1773 this->GetSpaceTimeHandler()->ComputeStateSparsityPattern(sparsity);
1778 "PDEProblemContainer::ComputeSparsityPattern");
1827 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1828 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1832 if (this->GetType().find(
"aux_functional") != std::string::npos)
1834 return aux_functionals_[this->GetTypeNum()]->HasFaces();
1838 if ((this->GetType() ==
"state"))
1840 return this->GetPDE().HasFaces();
1842 else if (this->GetType() ==
"adjoint_for_ee")
1844 return this->GetPDE().HasFaces()
1845 || aux_functionals_[functional_for_ee_num_]->HasFaces();
1849 throw DOpEException(
"Unknown Type: '" + this->GetType() +
"'!",
1850 "PDEProblemContainer::HasFaces");
1857 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1858 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1862 if ((this->GetType() ==
"state") || this->GetType() ==
"aux_functional")
1867 else if (this->GetType() ==
"adjoint_for_ee")
1869 return aux_functionals_[functional_for_ee_num_]->HasPoints();
1873 throw DOpEException(
"Unknown Type: '" + this->GetType() +
"'!",
1874 "PDEProblemContainer::HasFaces");
1881 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1882 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1886 if (this->GetType().find(
"aux_functional") != std::string::npos)
1892 if ((this->GetType() ==
"state") || this->GetType() ==
"adjoint_for_ee")
1894 return this->GetPDE().HasInterfaces();
1898 throw DOpEException(
"Unknown Type: '" + this->GetType() +
"'!",
1899 "PDEProblemContainer::HasFaces");
1906 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1907 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1910 unsigned int color,
const std::vector<bool>& comp_mask,
1914 assert(values->n_components() == comp_mask.size());
1916 unsigned int comp = dirichlet_colors_.size();
1917 for (
unsigned int i = 0; i < dirichlet_colors_.size(); ++i)
1919 if (dirichlet_colors_[i] == color)
1925 if (comp != dirichlet_colors_.size())
1927 std::stringstream s;
1928 s <<
"DirichletColor" << color <<
" has multiple occurrences !";
1930 "PDEProblemContainer::SetDirichletBoundary");
1932 dirichlet_colors_.push_back(color);
1933 dirichlet_comps_.push_back(comp_mask);
1935 DD, VECTOR, dealdim>(*values);
1936 primal_dirichlet_values_.push_back(data);
1941 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1942 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1943 const std::vector<unsigned int>&
1946 if ((this->GetType() ==
"state") || this->GetType() ==
"adjoint_for_ee")
1948 return dirichlet_colors_;
1953 "PDEProblemContainer::GetDirichletColors");
1959 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1960 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1961 const std::vector<bool>&
1963 unsigned int color)
const
1965 if ((this->GetType() ==
"state" || this->GetType() ==
"adjoint_for_ee"))
1967 unsigned int comp = dirichlet_colors_.size();
1968 for (
unsigned int i = 0; i < dirichlet_colors_.size(); ++i)
1970 if (dirichlet_colors_[i] == color)
1976 if (comp == dirichlet_colors_.size())
1978 std::stringstream s;
1979 s <<
"DirichletColor" << color <<
" has not been found !";
1981 "PDEProblemContainer::GetDirichletCompMask");
1983 return dirichlet_comps_[comp];
1988 "PDEProblemContainer::GetDirichletCompMask");
1994 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
1995 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
1996 const Function<dealdim>&
1999 const std::map<std::string,
const dealii::Vector<double>*> ¶m_values,
2000 const std::map<std::string, const VECTOR*> &domain_values)
const
2003 unsigned int col = dirichlet_colors_.size();
2004 if ((this->GetType() ==
"state") || this->GetType() ==
"adjoint_for_ee")
2006 for (
unsigned int i = 0; i < dirichlet_colors_.size(); ++i)
2008 if (dirichlet_colors_[i] == color)
2014 if (col == dirichlet_colors_.size())
2016 std::stringstream s;
2017 s <<
"DirichletColor" << color <<
" has not been found !";
2019 "PDEProblemContainer::GetDirichletValues");
2025 "PDEProblemContainer::GetDirichletValues");
2028 if (this->GetType() ==
"state")
2030 primal_dirichlet_values_[col]->ReInit(param_values, domain_values,
2032 return *(primal_dirichlet_values_[col]);
2034 else if (this->GetType() ==
"adjoint_for_ee"
2035 || (this->GetType() ==
"adjoint_hessian"))
2037 return *(zero_dirichlet_values_);
2042 "PDEProblemContainer::GetDirichletValues");
2048 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
2049 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
2050 const std::vector<unsigned int>&
2053 if (this->GetType() ==
"state")
2055 return state_boundary_equation_colors_;
2057 else if (this->GetType() ==
"adjoint_for_ee")
2059 return adjoint_boundary_equation_colors_;
2064 "PDEProblemContainer::GetBoundaryEquationColors");
2070 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
2071 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
2077 unsigned int comp = state_boundary_equation_colors_.size();
2078 for (
unsigned int i = 0; i < state_boundary_equation_colors_.size();
2081 if (state_boundary_equation_colors_[i] == color)
2087 if (comp != state_boundary_equation_colors_.size())
2089 std::stringstream s;
2090 s <<
"Boundary Equation Color" << color
2091 <<
" has multiple occurences !";
2093 "PDEProblemContainer::SetBoundaryEquationColors");
2095 state_boundary_equation_colors_.push_back(color);
2098 unsigned int comp = adjoint_boundary_equation_colors_.size();
2099 for (
unsigned int i = 0; i < adjoint_boundary_equation_colors_.size();
2102 if (adjoint_boundary_equation_colors_[i] == color)
2108 if (comp != adjoint_boundary_equation_colors_.size())
2115 adjoint_boundary_equation_colors_.push_back(color);
2123 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
2124 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
2125 const std::vector<unsigned int>&
2129 if (this->GetType() ==
"cost_functional"
2130 || this->GetType() ==
"aux_functional"
2131 || this->GetType() ==
"error_evaluation")
2133 return boundary_functional_colors_;
2138 "PDEProblemContainer::GetBoundaryFunctionalColors");
2144 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
2145 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
2151 unsigned int comp = boundary_functional_colors_.size();
2152 for (
unsigned int i = 0; i < boundary_functional_colors_.size(); ++i)
2154 if (boundary_functional_colors_[i] == color)
2160 if (comp != boundary_functional_colors_.size())
2162 std::stringstream s;
2163 s <<
"Boundary Functional Color" << color
2164 <<
" has multiple occurences !";
2166 "PDEProblemContainer::SetBoundaryFunctionalColors");
2168 boundary_functional_colors_.push_back(color);
2171 unsigned int comp = adjoint_boundary_equation_colors_.size();
2172 for (
unsigned int i = 0; i < adjoint_boundary_equation_colors_.size();
2175 if (adjoint_boundary_equation_colors_[i] == color)
2181 if (comp != adjoint_boundary_equation_colors_.size())
2188 adjoint_boundary_equation_colors_.push_back(color);
2195 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
2196 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
2200 return this->GetPDE().GetStateNBlocks();
2205 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
2206 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
2210 if ((this->GetType() ==
"state") || (this->GetType() ==
"adjoint_for_ee"))
2212 return this->GetStateNBlocks();
2217 "PDEProblemContainer::GetNBlocks");
2223 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
2224 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
2227 unsigned int b)
const
2229 if ((this->GetType() ==
"state") || (this->GetType() ==
"adjoint_for_ee"))
2231 return GetSpaceTimeHandler()->GetStateDoFsPerBlock(b);
2236 "PDEProblemContainer::GetDoFsPerBlock");
2242 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
2243 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
2244 const std::vector<unsigned int>&
2247 if ((this->GetType() ==
"state") || (this->GetType() ==
"adjoint_for_ee"))
2249 return GetSpaceTimeHandler()->GetStateDoFsPerBlock();
2254 "PDEProblemContainer::GetDoFsPerBlock");
2260 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
2261 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
2262 const dealii::ConstraintMatrix&
2268 if ((this->GetType() ==
"state") || (this->GetType() ==
"adjoint_for_ee"))
2270 return GetSpaceTimeHandler()->GetStateDoFConstraints();
2275 "PDEProblemContainer::GetDoFConstraints");
2281 template<
typename PDE,
typename DD,
typename SPARSITYPATTERN,
typename VECTOR,
2282 int dealdim,
template<
int,
int>
class FE,
template<
int,
int>
class DH>
2286 if (this->GetType() ==
"cost_functional")
2288 else if (this->GetType() ==
"aux_functional")
2289 return aux_functionals_[this->GetTypeNum()]->NeedTime();
2290 else if (this->GetType() ==
"error_evaluation")
2291 return aux_functionals_[functional_for_ee_num_]->NeedTime();
2294 "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:1060
void ElementMatrix(const DATACONTAINER &dc, dealii::FullMatrix< double > &local_entry_matrix, double scale=1., double scale_ico=1.)
Definition: pdeproblemcontainer.h:1447
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:1884
PDEProblemContainer< PDE, DD, SPARSITYPATTERN, VECTOR, dealdim, FE, DH > & GetBaseProblem()
Definition: pdeproblemcontainer.h:141
void ReInit(std::string algo_type)
Definition: pdeproblemcontainer.h:938
unsigned int GetStateNBlocks()
Definition: pdeproblemcontainer.h:827
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:1997
void FaceRhs(const FACEDATACONTAINER &dc, dealii::Vector< double > &local_vector, double scale=1.)
Definition: pdeproblemcontainer.h:1387
dealii::UpdateFlags GetFaceUpdateFlags() const
Definition: pdeproblemcontainer.h:1683
virtual std::string GetName() const
Definition: pdeproblemcontainer.h:114
unsigned int GetStateNBlocks() const
Definition: pdeproblemcontainer.h:2198
std::string GetDoFType() const
Definition: pdeproblemcontainer.h:1620
void SetBoundaryFunctionalColors(unsigned int color)
Definition: pdeproblemcontainer.h:2147
Definition: elementdatacontainer.h:58
unsigned int GetNFunctionals() const
Definition: pdeproblemcontainer.h:698
void ComputeSparsityPattern(SPARSITYPATTERN &sparsity) const
Definition: pdeproblemcontainer.h:1768
void FaceEquation(const FACEDATACONTAINER &dc, dealii::Vector< double > &local_vector, double scale, double scale_ico)
Definition: pdeproblemcontainer.h:1240
Definition: problemcontainer_internal.h:38
Definition: optproblemcontainer.h:70
const std::vector< unsigned int > & GetDirichletColors() const
Definition: pdeproblemcontainer.h:1944
std::vector< unsigned int > & GetStateBlockComponent()
Definition: pdeproblemcontainer.h:835
const dealii::ConstraintMatrix & GetDoFConstraints() const
Definition: pdeproblemcontainer.h:2263
void AddFunctional(FunctionalInterface< ElementDataContainer, FaceDataContainer, DH, VECTOR, dealdim > *F)
Definition: pdeproblemcontainer.h:641
std::string GetFunctionalType() const
Definition: pdeproblemcontainer.h:1710
const StateSpaceTimeHandler< FE, DH, SPARSITYPATTERN, VECTOR, dealdim > * GetSpaceTimeHandler() const
Definition: pdeproblemcontainer.h:770
Definition: timeiterator.h:62
bool NeedTimeFunctional() const
Definition: pdeproblemcontainer.h:2284
const PDE & GetPDE() const
Definition: problemcontainer_internal.h:103
unsigned int GetNBlocks() const
Definition: pdeproblemcontainer.h:2208
const dealii::Function< dealdim > & GetInitialValues() const
Definition: pdeproblemcontainer.h:614
const std::vector< unsigned int > & GetDoFsPerBlock() const
Definition: pdeproblemcontainer.h:2245
StateSpaceTimeHandler< FE, DH, SPARSITYPATTERN, VECTOR, dealdim > * GetSpaceTimeHandler()
Definition: pdeproblemcontainer.h:778
Definition: facedatacontainer.h:50
void ElementTimeEquationExplicit(const DATACONTAINER &dc, dealii::Vector< double > &local_vector, double scale=1.)
Definition: pdeproblemcontainer.h:1211
void SetBoundaryEquationColors(unsigned int color)
Definition: pdeproblemcontainer.h:2073
dealii::UpdateFlags GetUpdateFlags() const
Definition: pdeproblemcontainer.h:1657
bool HasPoints() const
Definition: pdeproblemcontainer.h:1860
const std::map< std::string, unsigned int > & GetFunctionalPosition() const
Definition: pdeproblemcontainer.h:849
std::string GetFunctionalName() const
Definition: pdeproblemcontainer.h:1728
Definition: statespacetimehandler.h:59
double ElementFunctional(const DATACONTAINER &edc)
Definition: pdeproblemcontainer.h:995
const std::vector< unsigned int > & GetBoundaryEquationColors() const
Definition: pdeproblemcontainer.h:2051
Definition: primaldirichletdata.h:41
void ElementRhs(const DATACONTAINER &dc, dealii::Vector< double > &local_vector, double scale=1.)
Definition: pdeproblemcontainer.h:1328
void ElementEquation(const DATACONTAINER &edc, dealii::Vector< double > &local_vector, double scale, double scale_ico)
Definition: pdeproblemcontainer.h:1155
Definition: functionalinterface.h:53
void SetDirichletBoundaryColors(unsigned int color, const std::vector< bool > &comp_mask, const DD *values)
Definition: pdeproblemcontainer.h:1909
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:1417
bool HasFaces() const
Definition: pdeproblemcontainer.h:1830
void AddAuxiliaryToIntegrator(INTEGRATOR &)
Definition: pdeproblemcontainer.h:810
void ElementTimeMatrixExplicit(const DATACONTAINER &dc, dealii::FullMatrix< double > &local_entry_matrix)
Definition: pdeproblemcontainer.h:1506
DOpEExceptionHandler< VECTOR > * GetExceptionHandler()
Definition: pdeproblemcontainer.h:745
PDEProblemContainer(PDE &pde, StateSpaceTimeHandler< FE, DH, SPARSITYPATTERN, VECTOR, dealdim > &STH)
Definition: pdeproblemcontainer.h:899
double FaceFunctional(const FACEDATACONTAINER &fdc)
Definition: pdeproblemcontainer.h:1091
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:916
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:1639
void BoundaryEquation(const FACEDATACONTAINER &dc, dealii::Vector< double > &local_vector, double scale, double scale_ico)
Definition: pdeproblemcontainer.h:1298
void ElementTimeEquation(const DATACONTAINER &dc, dealii::Vector< double > &local_vector, double scale=1.)
Definition: pdeproblemcontainer.h:1183
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:1360
Definition: pdeproblemcontainer.h:100
bool AtInterface(ELEMENTITERATOR &element, unsigned int face) const
Definition: pdeproblemcontainer.h:855
const std::vector< bool > & GetDirichletCompMask(unsigned int color) const
Definition: pdeproblemcontainer.h:1962
double AlgebraicFunctional(const std::map< std::string, const dealii::Vector< double > * > &values, const std::map< std::string, const VECTOR * > &block_values)
Definition: pdeproblemcontainer.h:1121
Definition: dopeexception.h:35
void SetTime(double time, unsigned int time_dof_number, const TimeIterator &interval)
Definition: pdeproblemcontainer.h:1746
void InterfaceEquation(const FACEDATACONTAINER &dc, dealii::Vector< double > &local_vector, double scale, double scale_ico)
Definition: pdeproblemcontainer.h:1269
void DeleteAuxiliaryFromIntegrator(INTEGRATOR &)
Definition: pdeproblemcontainer.h:819
void SetType(std::string type, unsigned int num=0)
Definition: pdeproblemcontainer.h:975
const std::vector< unsigned int > & GetBoundaryFunctionalColors() const
Definition: pdeproblemcontainer.h:2126
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:1023
Definition: optproblemcontainer.h:72
void FaceMatrix(const FACEDATACONTAINER &dc, dealii::FullMatrix< double > &local_entry_matrix, double scale=1., double scale_ico=1.)