DOpE
optproblemcontainer.h
Go to the documentation of this file.
1 
24 #ifndef OptProblemContainer_H_
25 #define OptProblemContainer_H_
26 
27 #include "dopeexceptionhandler.h"
28 #include "outputhandler.h"
29 #include "functionalinterface.h"
30 #include "dofhandler_wrapper.h"
31 #include "fevalues_wrapper.h"
32 #include "function_wrapper.h"
33 #include "spacetimehandler.h"
34 #include "primaldirichletdata.h"
35 #include "tangentdirichletdata.h"
36 #include "transposeddirichletdatainterface.h"
37 #include "transposedgradientdirichletdata.h"
38 #include "transposedhessiandirichletdata.h"
39 #include "constraintvector.h"
40 #include "controlvector.h"
41 #include "statevector.h"
42 #include "elementdatacontainer.h"
43 #include "facedatacontainer.h"
44 #include "stateproblem.h"
45 #include "dopetypes.h"
46 #include "dwrdatacontainer.h"
47 #include "problemcontainer_internal.h"
48 
49 #include <lac/vector.h>
50 #include <lac/full_matrix.h>
51 #include <grid/tria_iterator.h>
52 #include <dofs/dof_handler.h>
53 #include <dofs/dof_accessor.h>
54 #include <dofs/dof_tools.h>
55 #include <fe/fe_system.h>
56 #include <fe/fe_values.h>
57 #include <fe/mapping.h>
58 #include <base/quadrature_lib.h>
59 #include <lac/block_sparsity_pattern.h>
60 #include <lac/compressed_simple_sparsity_pattern.h>
61 
62 #include <assert.h>
63 #include <string>
64 #include <vector>
65 
66 namespace DOpE
67 {
68  //Predeclaration necessary
69  template<typename VECTOR>
70  class DOpEOutputHandler;
71  template<typename VECTOR>
72  class DOpEExceptionHandler;
74 
96  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
97  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
98  typename VECTOR, int dopedim, int dealdim,
99  template<int, int> class FE = dealii::FESystem,
100  template<int, int> class DH = dealii::DoFHandler>
101  class OptProblemContainer : public ProblemContainerInternal<PDE>
102  {
103  public:
104  OptProblemContainer(FUNCTIONAL& functional, PDE& pde,
105  CONSTRAINTS& constraints,
106  SpaceTimeHandler<FE, DH, SPARSITYPATTERN, VECTOR, dopedim, dealdim>& STH);
107 
108  /******************************************************/
109 
110  virtual ~OptProblemContainer();
111 
112  /******************************************************/
113 
114  virtual std::string
115  GetName() const
116  {
117  return "OptProblem";
118  }
119 
120  /******************************************************/
124  StateProblem<
125  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD,
126  CONSTRAINTS, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>,
127  PDE, DD, SPARSITYPATTERN, VECTOR, dealdim>&
128  GetStateProblem()
129  {
130  if (state_problem_ == NULL)
131  {
132  state_problem_ = new StateProblem<
133  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD,
134  CONSTRAINTS, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE,
135  DH>, PDE, DD, SPARSITYPATTERN, VECTOR, dealdim>(
136  *this, this->GetPDE());
137  }
138  return *state_problem_;
139  }
140 
141  //TODO This is Pfush needed to split into different subproblems and allow optproblem to
142  //be substituted as any of these problems. Can be removed once the splitting is complete.
143  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD,
144  CONSTRAINTS, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>&
145  GetBaseProblem()
146  {
147  return *this;
148  }
149  /******************************************************/
150 
158  void
159  ReInit(std::string algo_type);
160 
161  /******************************************************/
162 
163  void
164  RegisterOutputHandler(DOpEOutputHandler<VECTOR>* OH)
165  {
166  OutputHandler_ = OH;
167  }
168 
169  /******************************************************/
170 
171  void
172  RegisterExceptionHandler(DOpEExceptionHandler<VECTOR>* OH)
173  {
174  ExceptionHandler_ = OH;
175  }
176 
177  /******************************************************/
178 
186  void
187  SetType(std::string type, unsigned int num = 0);
188 
189  /******************************************************/
190 
206  template<typename DATACONTAINER>
207  double
208  ElementFunctional(const DATACONTAINER& edc);
209 
210  /******************************************************/
211 
221  double
222  PointFunctional(
223  const std::map<std::string, const dealii::Vector<double>*> &param_values,
224  const std::map<std::string, const VECTOR*> &domain_values);
225 
226  /******************************************************/
227 
244  template<typename FACEDATACONTAINER>
245  double
246  BoundaryFunctional(const FACEDATACONTAINER& fdc);
247 
248  /******************************************************/
249 
256  template<typename FACEDATACONTAINER>
257  double
258  FaceFunctional(const FACEDATACONTAINER& fdc);
259 
260  /******************************************************/
267  double
268  AlgebraicFunctional(
269  const std::map<std::string, const dealii::Vector<double>*> &values,
270  const std::map<std::string, const VECTOR*> &block_values);
271 
278  void
279  AlgebraicResidual(VECTOR& residual,
280  const std::map<std::string, const dealii::Vector<double>*> &values,
281  const std::map<std::string, const VECTOR*> &block_values);
282 
283  /******************************************************/
284 
308  template<typename DATACONTAINER>
309  void
310  ElementEquation(const DATACONTAINER& edc,
311  dealii::Vector<double> &local_vector, double scale,
312  double scale_ico);
313 
314  /******************************************************/
315 
321  template<typename DATACONTAINER>
322  void
323  ElementTimeEquation(const DATACONTAINER& dc,
324  dealii::Vector<double> &local_vector, double scale = 1.);
325 
326  /******************************************************/
327 
337  template<typename DATACONTAINER>
338  void
339  ElementTimeEquationExplicit(const DATACONTAINER& dc,
340  dealii::Vector<double> &local_vector, double scale = 1.);
341 
342  /******************************************************/
343 
355  template<typename DATACONTAINER>
356  void
357  ElementRhs(const DATACONTAINER& dc,
358  dealii::Vector<double> &local_vector, double scale = 1.);
359 
360  /******************************************************/
361 
374  void
375  PointRhs(
376  const std::map<std::string, const dealii::Vector<double>*> &param_values,
377  const std::map<std::string, const VECTOR*> &domain_values,
378  VECTOR& rhs_vector, double scale = 1.);
379 
380  /******************************************************/
381 
405  template<typename DATACONTAINER>
406  void
407  ElementMatrix(const DATACONTAINER& dc,
408  dealii::FullMatrix<double> &local_entry_matrix, double scale = 1.,
409  double scale_ico = 1.);
410 
411  /******************************************************/
412 
424  template<typename DATACONTAINER>
425  void
426  ElementTimeMatrix(const DATACONTAINER& dc,
427  dealii::FullMatrix<double> &local_entry_matrix);
428 
429  /******************************************************/
430 
440  template<typename DATACONTAINER>
441  void
442  ElementTimeMatrixExplicit(const DATACONTAINER& dc,
443  dealii::FullMatrix<double> &local_entry_matrix);
444 
445  /******************************************************/
446 
453  template<typename FACEDATACONTAINER>
454  void
455  FaceEquation(const FACEDATACONTAINER& dc,
456  dealii::Vector<double> &local_vector, double scale,
457  double scale_ico);
458 
459  /******************************************************/
466  //FIXME maybe InterfaceEquation and InterfaceMatrix could get
467  //integrated into FaceEquation and FaceMatrix?
468  template<typename FACEDATACONTAINER>
469  void
470  InterfaceEquation(const FACEDATACONTAINER& dc,
471  dealii::Vector<double> &local_vector, double scale,
472  double scale_ico);
473 
474  /******************************************************/
481  template<typename FACEDATACONTAINER>
482  void
483  FaceRhs(const FACEDATACONTAINER& dc,
484  dealii::Vector<double> &local_vector, double scale = 1.);
485 
486  /******************************************************/
487 
494  template<typename FACEDATACONTAINER>
495  void
496  FaceMatrix(const FACEDATACONTAINER& dc,
497  dealii::FullMatrix<double> &local_entry_matrix, double scale = 1.,
498  double scale_ico = 1.);
499 
500  /******************************************************/
507  template<typename FACEDATACONTAINER>
508  void
509  InterfaceMatrix(const FACEDATACONTAINER& dc,
510  dealii::FullMatrix<double> &local_entry_matrix, double scale = 1.,
511  double scale_ico = 1.);
512 
513  /******************************************************/
514 
521  template<typename FACEDATACONTAINER>
522  void
523  BoundaryEquation(const FACEDATACONTAINER& dc,
524  dealii::Vector<double> &local_vector, double scale,
525  double scale_ico);
526 
527  /******************************************************/
528 
535  template<typename FACEDATACONTAINER>
536  void
537  BoundaryRhs(const FACEDATACONTAINER& dc,
538  dealii::Vector<double> &local_vector, double scale = 1.);
539 
540  /******************************************************/
541 
548  template<typename FACEDATACONTAINER>
549  void
550  BoundaryMatrix(const FACEDATACONTAINER& dc,
551  dealii::FullMatrix<double> &local_matrix, double scale = 1.,
552  double scale_ico = 1.);
553  /******************************************************/
554  void
555  ComputeLocalControlConstraints(VECTOR& constraints,
556  const std::map<std::string, const dealii::Vector<double>*> &values,
557  const std::map<std::string, const VECTOR*> &block_values);
558  /******************************************************/
559 
565  void
566  GetControlBoxConstraints(VECTOR& lb, VECTOR& ub) const
567  {
568  this->GetConstraints()->GetControlBoxConstraints(lb, ub);
569  }
570 
571  /******************************************************/
572 
573  const FE<dealdim, dealdim>&
574  GetFESystem() const;
575 
576  /******************************************************/
583  bool
584  HasFaces() const;
585 
586  /******************************************************/
593  bool
594  HasPoints() const;
595 
596  /******************************************************/
603  bool
604  HasInterfaces() const;
605 
606  /******************************************************/
607 
608  dealii::UpdateFlags
609  GetUpdateFlags() const;
610 
611  /******************************************************/
612 
613  dealii::UpdateFlags
614  GetFaceUpdateFlags() const;
615 
616  /******************************************************/
617 
618  void
619  SetControlDirichletBoundaryColors(unsigned int color,
620  const std::vector<bool>& comp_mask,
621  const DOpEWrapper::Function<dealdim>* values);
622 
623  /******************************************************/
624 
625  void
626  SetDirichletBoundaryColors(unsigned int color,
627  const std::vector<bool>& comp_mask, const DD* values);
628 
629  /******************************************************/
630 
631  const std::vector<unsigned int>&
632  GetDirichletColors() const;
633  const std::vector<unsigned int>&
634  GetTransposedDirichletColors() const;
635  const std::vector<bool>&
636  GetDirichletCompMask(unsigned int color) const;
637  const std::vector<bool>&
638  GetTransposedDirichletCompMask(unsigned int color) const;
639 
640  /******************************************************/
641 
642  const dealii::Function<dealdim> &
643  GetDirichletValues(unsigned int color,
644  const std::map<std::string, const dealii::Vector<double>*> &param_values,
645  const std::map<std::string, const VECTOR*> &domain_values) const;
646 
647  /******************************************************/
648 
649  const TransposedDirichletDataInterface<dealdim> &
650  GetTransposedDirichletValues(unsigned int color,
651  const std::map<std::string, const dealii::Vector<double>*> &param_values,
652  const std::map<std::string, const VECTOR*> &domain_values) const;
653 
654  /******************************************************/
655 
656  void
657  SetInitialValues(const dealii::Function<dealdim>* values)
658  {
659  initial_values_ = values;
660  }
661  const dealii::Function<dealdim>&
662  GetInitialValues() const
663  {
664  return *initial_values_;
665  }
666 
667  /******************************************************/
668 
669  void
670  SetControlBoundaryEquationColors(unsigned int color);
671  void
672  SetBoundaryEquationColors(unsigned int color);
673  const std::vector<unsigned int>&
674  GetBoundaryEquationColors() const;
675  void
676  SetBoundaryFunctionalColors(unsigned int color);
677  const std::vector<unsigned int>&
678  GetBoundaryFunctionalColors() const;
679 
680  /******************************************************/
681 
691  void
692  AddFunctional(FUNCTIONAL_INTERFACE* F)
693  {
694  aux_functionals_.push_back(F);
695  if (functional_position_.find(F->GetName())
696  != functional_position_.end())
697  {
698  throw DOpEException(
699  "You cant add two functionals with the same name.",
700  "OPTProblemContainer::AddFunctional");
701  }
702  functional_position_[F->GetName()] = aux_functionals_.size();
703  //remember! At functional_values_[0] we store always the cost functional!
704  }
705 
706  /******************************************************/
707 
717  void
718  SetFunctionalForErrorEstimation(std::string functional_name)
719  {
720  if (GetFunctional()->GetName() == functional_name)
721  {
722  functional_for_ee_is_cost_ = true;
723  functional_for_ee_num_ = dealii::numbers::invalid_unsigned_int;
724  }
725  else
726  {
727  functional_for_ee_is_cost_ = false;
728  bool found = false;
729  //we go through all aux functionals.
730  for (unsigned int i = 0; i < this->GetNFunctionals(); i++)
731  {
732  if (aux_functionals_[i]->GetName() == functional_name)
733  {
734  //if the names match, we have found our functional.
735  found = true;
736  functional_for_ee_num_ = i;
737  }
738  }
739  //If we have not found a functional with the given name,
740  //we throw an error.
741  if (!found)
742  {
743  throw DOpEException(
744  "Can't find functional " + functional_name
745  + " in aux_functionals_",
746  "Optproblem::SetFunctionalForErrorEstimation");
747  }
748  }
749  }
750 
751  /******************************************************/
752 
753  unsigned int
754  GetNFunctionals() const
755  {
756  return aux_functionals_.size();
757  }
758 
759  /******************************************************/
760 
761  unsigned int
762  GetControlNBlocks() const;
763 
764  /******************************************************/
765 
766  unsigned int
767  GetStateNBlocks() const;
768 
769  /******************************************************/
770 
771  unsigned int
772  GetNBlocks() const;
773 
774  /******************************************************/
775 
776  unsigned int
777  GetDoFsPerBlock(unsigned int b) const;
778 
779  /******************************************************/
780 
781  const std::vector<unsigned int>&
782  GetDoFsPerBlock() const;
783 
784  /******************************************************/
785 
786  const dealii::ConstraintMatrix&
787  GetDoFConstraints() const;
788 
789  /******************************************************/
790 
791  std::string
792  GetDoFType() const;
793  std::string
794  GetFunctionalType() const;
795  std::string
796  GetFunctionalName() const;
797  std::string
798  GetConstraintType() const;
799 
800  /******************************************************/
801 
802  bool
803  NeedTimeFunctional() const;
804 
805  /******************************************************/
806 
807  bool
808  HasControlInDirichletData() const;
809 
810  /******************************************************/
811 
812  DOpEExceptionHandler<VECTOR>*
813  GetExceptionHandler()
814  {
815  assert(ExceptionHandler_);
816  return ExceptionHandler_;
817  }
818 
819  /******************************************************/
820 
821  DOpEOutputHandler<VECTOR>*
822  GetOutputHandler()
823  {
824  assert(OutputHandler_);
825  return OutputHandler_;
826  }
827 
828  /******************************************************/
829 
838  void
839  SetTime(double time, const TimeIterator& interval, bool initial = false);
840 
841  /******************************************************/
842 
843  const SpaceTimeHandler<FE, DH, SPARSITYPATTERN, VECTOR, dopedim, dealdim>*
844  GetSpaceTimeHandler() const
845  {
846  return STH_;
847  }
848 
849  /******************************************************/
850 
851  SpaceTimeHandler<FE, DH, SPARSITYPATTERN, VECTOR, dopedim, dealdim>*
852  GetSpaceTimeHandler()
853  {
854  return STH_;
855  }
856 
857  /******************************************************/
858 
859  void
860  ComputeSparsityPattern(SPARSITYPATTERN & sparsity) const;
861 
862  /******************************************************/
863 
864  void
865  PostProcessConstraints(ConstraintVector<VECTOR>& g) const;
866 
867  /******************************************************/
868 
869  void
870  AddAuxiliaryControl(const ControlVector<VECTOR>* c, std::string name);
871 
872  /******************************************************/
873 
874  void
875  AddAuxiliaryState(const StateVector<VECTOR>* c, std::string name);
876 
877  /******************************************************/
878 
879  void
880  AddAuxiliaryConstraint(const ConstraintVector<VECTOR>* c,
881  std::string name);
882 
883  /******************************************************/
884 
885  const ControlVector<VECTOR>*
886  GetAuxiliaryControl(std::string name) const;
887 
888  /******************************************************/
889 
890  const StateVector<VECTOR>*
891  GetAuxiliaryState(std::string name) const;
892 
893  /******************************************************/
894 
895  void
896  DeleteAuxiliaryControl(std::string name);
897 
898  /******************************************************/
899 
900  void
901  DeleteAuxiliaryState(std::string name);
902 
903  /******************************************************/
904 
905  void
906  DeleteAuxiliaryConstraint(std::string name);
907 
908  /*****************************************************************/
909 
910  const ConstraintVector<VECTOR>*
911  GetAuxiliaryConstraint(std::string name)
912  {
913  typename std::map<std::string, const ConstraintVector<VECTOR> *>::iterator it =
914  auxiliary_constraints_.find(name);
915  if (it == auxiliary_constraints_.end())
916  {
917  throw DOpEException("Could not find data" + name,
918  "OptProblemContainer::GetAuxiliaryConstraint");
919  }
920  return it->second;
921  }
922  /*****************************************************************/
931  template<typename INTEGRATOR>
932  void
933  AddAuxiliaryToIntegrator(INTEGRATOR& integrator)
934  {
935  {
936  //Only add control vector if the vecor is distributed in time, or
937  //we are calculating the initial value.
938  if((GetSpaceTimeHandler()->GetControlType()==DOpEtypes::ControlType::initial && initial_)
939  || (GetSpaceTimeHandler()->GetControlType()!=DOpEtypes::ControlType::initial))
940  {
941  typename std::map<std::string, const ControlVector<VECTOR> *>::iterator it =
942  auxiliary_controls_.begin();
943  for (; it != auxiliary_controls_.end(); it++)
944  {
945  if (dopedim == dealdim)
946  {
947  integrator.AddDomainData(it->first,
948  &(it->second->GetSpacialVector()));
949  }
950  else if (dopedim == 0)
951  {
952  integrator.AddParamData(it->first,
953  &(it->second->GetSpacialVectorCopy()));
954  }
955  else
956  {
957  throw DOpEException("dopedim not implemented",
958  "OptProblemContainer::AddAuxiliaryToIntegrator");
959  }
960  }
961  }
962  }
963  {
964  typename std::map<std::string, const StateVector<VECTOR> *>::iterator it =
965  auxiliary_state_.begin();
966  for (; it != auxiliary_state_.end(); it++)
967  {
968  integrator.AddDomainData(it->first,
969  &(it->second->GetSpacialVector()));
970  }
971  }
972  {
973  typename std::map<std::string, const ConstraintVector<VECTOR> *>::iterator it =
974  auxiliary_constraints_.begin();
975  for (; it != auxiliary_constraints_.end(); it++)
976  {
977  integrator.AddDomainData(it->first + "_local",
978  &(it->second->GetSpacialVector("local")));
979  integrator.AddParamData(it->first + "_global",
980  &(it->second->GetGlobalConstraints()));
981  }
982  }
983  }
984 
985  /******************************************************/
994  template<typename INTEGRATOR>
995  void
996  DeleteAuxiliaryFromIntegrator(INTEGRATOR& integrator)
997  {
998  {
999  //Only delete control vector if the vecor is distributed in time, or
1000  //we are calculating the initial value. Otherwise it would not have been added.
1001  if((GetSpaceTimeHandler()->GetControlType()==DOpEtypes::ControlType::initial && initial_)
1002  || (GetSpaceTimeHandler()->GetControlType()!=DOpEtypes::ControlType::initial))
1003  {
1004  typename std::map<std::string, const ControlVector<VECTOR> *>::iterator it =
1005  auxiliary_controls_.begin();
1006  for (; it != auxiliary_controls_.end(); it++)
1007  {
1008  if (dopedim == dealdim)
1009  {
1010  integrator.DeleteDomainData(it->first);
1011  }
1012  else if (dopedim == 0)
1013  {
1014  integrator.DeleteParamData(it->first);
1015  it->second->UnLockCopy();
1016  }
1017  else
1018  {
1019  throw DOpEException("dopedim not implemented",
1020  "OptProblemContainer::AddAuxiliaryToIntegrator");
1021  }
1022  }
1023  }
1024  }
1025  {
1026  typename std::map<std::string, const StateVector<VECTOR> *>::iterator it =
1027  auxiliary_state_.begin();
1028  for (; it != auxiliary_state_.end(); it++)
1029  {
1030  integrator.DeleteDomainData(it->first);
1031  }
1032  }
1033  {
1034  typename std::map<std::string, const ConstraintVector<VECTOR> *>::iterator it =
1035  auxiliary_constraints_.begin();
1036  for (; it != auxiliary_constraints_.end(); it++)
1037  {
1038  integrator.DeleteDomainData(it->first + "_local");
1039  integrator.DeleteParamData(it->first + "_global");
1040  }
1041  }
1042  }
1043 
1044  /******************************************************/
1045  const std::map<std::string, unsigned int>&
1046  GetFunctionalPosition() const
1047  {
1048  return functional_position_;
1049  }
1050 
1051  unsigned int
1052  GetStateNBlocks()
1053  {
1054  return this->GetPDE().GetStateNBlocks();
1055  }
1056 
1057  /******************************************************/
1058 
1059  std::vector<unsigned int>&
1060  GetControlBlockComponent()
1061  {
1062  return this->GetPDE().GetControlBlockComponent();
1063  }
1064  /******************************************************/
1065 
1066  std::vector<unsigned int>&
1067  GetStateBlockComponent()
1068  {
1069  return this->GetPDE().GetStateBlockComponent();
1070  }
1071 
1072  /******************************************************/
1073  /******************************************************/
1074  /****For the initial values ***************/
1075  template<typename DATACONTAINER>
1076  void
1077  Init_ElementEquation(const DATACONTAINER& edc,
1078  dealii::Vector<double> &local_vector, double scale,
1079  double scale_ico)
1080  {
1081  if (this->GetType() == "adjoint" || this->GetType() == "tangent"
1082  || this->GetType() == "adjoint_hessian")
1083  {
1084  this->GetPDE().Init_ElementEquation(edc, local_vector, scale,
1085  scale_ico);
1086  }
1087  else
1088  {
1089  throw DOpEException("Not implemented",
1090  "OptProblemContainer::Init_ElementEquation");
1091  }
1092  }
1093 
1094  template<typename DATACONTAINER>
1095  void
1096  Init_ElementRhs(const DATACONTAINER& edc,
1097  dealii::Vector<double> &local_vector, double scale)
1098  {
1099  //FIXME: We should take care of the cases of boundary and face functionals...
1100  if (this->GetType() == "adjoint")
1101  {
1102  if (GetFunctional()->NeedTime())
1103  {
1104  if (GetFunctional()->GetType().find("timelocal")
1105  != std::string::npos)
1106  {
1107  if (GetFunctional()->GetType().find("domain")
1108  != std::string::npos)
1109  {
1110  GetFunctional()->ElementValue_U(edc, local_vector, scale);
1111  }
1112  }
1113  }
1114  }
1115  else if (this->GetType() == "tangent")
1116  {
1117  this->GetPDE().Init_ElementRhs_QT(edc, local_vector, scale);
1118  }
1119  else if (this->GetType() == "adjoint_hessian")
1120  {
1121  if (GetFunctional()->NeedTime())
1122  {
1123  if (GetFunctional()->GetType().find("timelocal")
1124  != std::string::npos)
1125  {
1126  if (GetFunctional()->GetType().find("domain")
1127  != std::string::npos)
1128  {
1129  GetFunctional()->ElementValue_UU(edc, local_vector, scale);
1130  }
1131  }
1132  }
1133  }
1134  else
1135  {
1136  throw DOpEException("Not implemented",
1137  "OptProblemContainer::Init_ElementRhs");
1138  }
1139  }
1140 
1141  void
1142  Init_PointRhs(
1143  const std::map<std::string, const dealii::Vector<double>*> &param_values,
1144  const std::map<std::string, const VECTOR*> &domain_values,
1145  VECTOR& rhs_vector, double scale = 1.)
1146  {
1147  if (this->GetType() == "adjoint")
1148  {
1149  if (GetFunctional()->NeedTime())
1150  {
1151  if (GetFunctional()->GetType().find("timelocal")
1152  != std::string::npos)
1153  {
1154  if (GetFunctional()->GetType().find("point")
1155  != std::string::npos)
1156  {
1157  GetFunctional()->PointValue_U(param_values, domain_values,
1158  rhs_vector, scale);
1159  }
1160  }
1161  }
1162  }
1163  else if (this->GetType() == "tangent")
1164  {
1165  //Nothing to do for tangent, since no point sources are allowed in the initial condition.
1166  }
1167  else if (this->GetType() == "adjoint_hessian")
1168  {
1169  if (GetFunctional()->NeedTime())
1170  {
1171  if (GetFunctional()->GetType().find("timelocal")
1172  != std::string::npos)
1173  {
1174  if (GetFunctional()->GetType().find("point")
1175  != std::string::npos)
1176  {
1177  GetFunctional()->PointValue_UU(param_values, domain_values,
1178  rhs_vector, scale);
1179  }
1180  }
1181  }
1182  }
1183  else
1184  {
1185  throw DOpEException("Not implemented",
1186  "OptProblemContainer::Init_PointRhs");
1187  }
1188  }
1189 
1190  template<typename DATACONTAINER>
1191  void
1192  Init_ElementMatrix(const DATACONTAINER& edc,
1193  dealii::FullMatrix<double> &local_entry_matrix, double scale,
1194  double scale_ico)
1195  {
1196  if (this->GetType() == "adjoint" || this->GetType() == "tangent"
1197  || this->GetType() == "adjoint_hessian")
1198  {
1199  this->GetPDE().Init_ElementMatrix(edc, local_entry_matrix, scale,
1200  scale_ico);
1201  }
1202  else
1203  {
1204  throw DOpEException("Not implemented",
1205  "OptProblemContainer::Init_ElementMatrix");
1206  }
1207  }
1208 
1209  /******************************************************/
1210 
1215  bool
1216  EEFunctionalIsCost() const
1217  {
1218  return functional_for_ee_is_cost_;
1219  }
1220 
1221  protected:
1222  FUNCTIONAL*
1223  GetFunctional();
1224  const FUNCTIONAL*
1225  GetFunctional() const;
1226  CONSTRAINTS *
1227  GetConstraints()
1228  {
1229  return constraints_;
1230  }
1231  const CONSTRAINTS *
1232  GetConstraints() const
1233  {
1234  return constraints_;
1235  }
1236 
1237  const VECTOR*
1238  GetBlockVector(const std::map<std::string, const VECTOR*>& values,
1239  std::string name)
1240  {
1241  typename std::map<std::string, const VECTOR*>::const_iterator it =
1242  values.find(name);
1243  if (it == values.end())
1244  {
1245  throw DOpEException("Did not find " + name,
1246  "OptProblemContainer::GetBlockVector");
1247  }
1248  return it->second;
1249  }
1250  const dealii::Vector<double>*
1251  GetVector(const std::map<std::string, const Vector<double>*>& values,
1252  std::string name)
1253  {
1254  typename std::map<std::string, const Vector<double>*>::const_iterator it =
1255  values.find(name);
1256  if (it == values.end())
1257  {
1258  throw DOpEException("Did not find " + name,
1259  "OptProblemContainer::GetVector");
1260  }
1261  return it->second;
1262  }
1263  /******************************************************/
1264  private:
1265  DOpEExceptionHandler<VECTOR>* ExceptionHandler_;
1266  DOpEOutputHandler<VECTOR>* OutputHandler_;
1267  std::string algo_type_;
1268 
1269  bool functional_for_ee_is_cost_;
1270  double c_interval_length_, interval_length_;
1271  unsigned int functional_for_ee_num_;
1272  std::vector<FUNCTIONAL_INTERFACE*> aux_functionals_;
1273  std::map<std::string, unsigned int> functional_position_;
1274  FUNCTIONAL* functional_;
1275  CONSTRAINTS* constraints_;
1276  SpaceTimeHandler<FE, DH, SPARSITYPATTERN, VECTOR, dopedim, dealdim>* STH_;
1277 
1278  std::vector<unsigned int> control_dirichlet_colors_;
1279  std::vector<unsigned int> control_transposed_dirichlet_colors_;
1280  std::vector<std::vector<bool> > control_dirichlet_comps_;
1281  std::vector<const DOpEWrapper::Function<dealdim>*> control_dirichlet_values_;
1282  std::vector<
1283  TransposedGradientDirichletData<DD, VECTOR, dealdim>*> transposed_control_gradient_dirichlet_values_;
1284  std::vector<
1285  TransposedHessianDirichletData<DD, VECTOR, dealdim>*> transposed_control_hessian_dirichlet_values_;
1286 
1287  std::vector<unsigned int> dirichlet_colors_;
1288  std::vector<std::vector<bool> > dirichlet_comps_;
1289  std::vector<PrimalDirichletData<DD, VECTOR, dealdim>*> primal_dirichlet_values_;
1290  std::vector<TangentDirichletData<DD, VECTOR, dealdim>*> tangent_dirichlet_values_;
1291  const dealii::Function<dealdim>* zero_dirichlet_values_;
1292 
1293  const dealii::Function<dealdim>* initial_values_;
1294 
1295  std::vector<unsigned int> control_boundary_equation_colors_;
1296  std::vector<unsigned int> state_boundary_equation_colors_;
1297  std::vector<unsigned int> adjoint_boundary_equation_colors_;
1298 
1299  std::vector<unsigned int> boundary_functional_colors_;
1300 
1301  std::map<std::string, const ControlVector<VECTOR>*> auxiliary_controls_;
1302  std::map<std::string, const StateVector<VECTOR>*> auxiliary_state_;
1303  std::map<std::string, const ConstraintVector<VECTOR>*> auxiliary_constraints_;
1304 
1305  bool initial_; //Do we solve the problem at initial time?
1306 
1307  StateProblem<
1308  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD,
1309  CONSTRAINTS, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>,
1310  PDE, DD, SPARSITYPATTERN, VECTOR, dealdim> * state_problem_;
1311 
1312  friend class StateProblem<
1313  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD,
1314  CONSTRAINTS, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>,
1315  PDE, DD, SPARSITYPATTERN, VECTOR, dealdim> ;
1316  };
1317  /******************************************************/
1318 
1319  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
1320  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
1321  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
1322  template<int, int> class DH>
1323  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
1324  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::OptProblemContainer(
1325  FUNCTIONAL& functional, PDE& pde, CONSTRAINTS& constraints,
1326  SpaceTimeHandler<FE, DH, SPARSITYPATTERN, VECTOR, dopedim, dealdim>& STH) :
1327  ProblemContainerInternal<PDE>(pde), functional_(&functional), constraints_(
1328  &constraints), STH_(&STH), state_problem_(NULL)
1329  {
1330  ExceptionHandler_ = NULL;
1331  OutputHandler_ = NULL;
1332  zero_dirichlet_values_ = new ZeroFunction<dealdim>(
1333  this->GetPDE().GetStateNComponents());
1334  algo_type_ = "";
1335  functional_position_[functional_->GetName()] = 0;
1336  //remember! At functional_values_[0] we store always the cost functional!
1337  functional_for_ee_num_ = dealii::numbers::invalid_unsigned_int;
1338  c_interval_length_ = 1.;
1339  interval_length_ = 1.;
1340  }
1341 
1342  /******************************************************/
1343 
1344  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
1345  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
1346  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
1347  template<int, int> class DH>
1348  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
1349  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::~OptProblemContainer()
1350  {
1351  if (zero_dirichlet_values_ != NULL)
1352  {
1353  delete zero_dirichlet_values_;
1354  }
1355 
1356  for (unsigned int i = 0;
1357  i < transposed_control_gradient_dirichlet_values_.size(); i++)
1358  {
1359  if (transposed_control_gradient_dirichlet_values_[i] != NULL)
1360  delete transposed_control_gradient_dirichlet_values_[i];
1361  }
1362  for (unsigned int i = 0;
1363  i < transposed_control_hessian_dirichlet_values_.size(); i++)
1364  {
1365  if (transposed_control_hessian_dirichlet_values_[i] != NULL)
1366  delete transposed_control_hessian_dirichlet_values_[i];
1367  }
1368  for (unsigned int i = 0; i < primal_dirichlet_values_.size(); i++)
1369  {
1370  if (primal_dirichlet_values_[i] != NULL)
1371  delete primal_dirichlet_values_[i];
1372  }
1373  for (unsigned int i = 0; i < tangent_dirichlet_values_.size(); i++)
1374  {
1375  if (tangent_dirichlet_values_[i] != NULL)
1376  delete tangent_dirichlet_values_[i];
1377  }
1378  if (state_problem_ != NULL)
1379  {
1380  delete state_problem_;
1381  }
1382  }
1383 
1384  /******************************************************/
1385 
1386  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
1387  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
1388  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
1389  template<int, int> class DH>
1390  void
1391  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
1392  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::ReInit(
1393  std::string algo_type)
1394  {
1395  if (state_problem_ != NULL)
1396  {
1397  delete state_problem_;
1398  state_problem_ = NULL;
1399  }
1400 
1401  if (algo_type_ != algo_type && algo_type_ != "")
1402  {
1403  throw DOpEException("Conflicting Algorithms!",
1404  "OptProblemContainer::ReInit");
1405  }
1406  else
1407  {
1408  algo_type_ = algo_type;
1409  this->SetTypeInternal("");
1410 
1411  if (algo_type_ == "reduced")
1412  {
1413  GetSpaceTimeHandler()->ReInit(this->GetPDE().GetControlNBlocks(),
1414  this->GetPDE().GetControlBlockComponent(),
1415  this->GetPDE().GetStateNBlocks(),
1416  this->GetPDE().GetStateBlockComponent());
1417  }
1418  else
1419  {
1420  throw DOpEException("Unknown Algorithm " + algo_type_,
1421  "OptProblemContainer::ReInit");
1422  }
1423  }
1424  }
1425 
1426  /******************************************************/
1427 
1428  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
1429  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
1430  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
1431  template<int, int> class DH>
1432  void
1433  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
1434  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::SetType(
1435  std::string type, unsigned int num)
1436  {
1437  if (this->GetType() != type || this->GetTypeNum() != num)
1438  {
1439  this->SetTypeNumInternal(num);
1440  this->SetTypeInternal(type);
1441  this->GetPDE().SetProblemType(type);
1442  if (functional_for_ee_num_ != dealii::numbers::invalid_unsigned_int)
1443  aux_functionals_[functional_for_ee_num_]->SetProblemType(type);
1444  this->GetConstraints()->SetProblemType(type, num);
1445 
1446 #if dope_dimension > 0
1447  if(dealdim == dopedim)
1448  {
1449  //Prepare DoFHandlerPointer
1450 
1451  {
1452  if(this->GetType() == "state" ||this->GetType() == "adjoint"
1453  || this->GetType() == "adjoint_for_ee" || this->GetType() == "cost_functional"
1454  || this->GetType() == "aux_functional" || this->GetType() == "functional_for_ee"
1455  || this->GetType() == "tangent" || this->GetType() == "adjoint_hessian"
1456  || this->GetType() == "error_evaluation"
1457  || this->GetType().find("constraints") != std::string::npos)
1458  {
1459  GetSpaceTimeHandler()->SetDoFHandlerOrdering(1,0);
1460  }
1461  else if (this->GetType() == "gradient"||this->GetType() == "hessian"||this->GetType() == "hessian_inverse" || this->GetType() == "global_constraint_gradient"|| this->GetType() == "global_constraint_hessian")
1462  {
1463  GetSpaceTimeHandler()->SetDoFHandlerOrdering(0,1);
1464  }
1465  else
1466  {
1467  throw DOpEException("problem_type_ : "+this->GetType()+" not implemented!", "OptProblemContainer::SetType");
1468  }
1469  }
1470  }
1471  else
1472  {
1473  throw DOpEException("dopedim not implemented", "OptProblemContainer::SetType");
1474  }
1475 #else
1476  //dopedim ==0
1477  {
1478  //Prepare DoFHandlerPointer
1479  {
1480 
1481  if (this->GetType() == "state" || this->GetType() == "adjoint"
1482  || this->GetType() == "adjoint_for_ee"
1483  || this->GetType() == "functional_for_ee"
1484  || this->GetType() == "cost_functional"
1485  || this->GetType() == "aux_functional"
1486  || this->GetType() == "tangent"
1487  || this->GetType() == "error_evaluation"
1488  || this->GetType() == "adjoint_hessian")
1489  {
1490  GetSpaceTimeHandler()->SetDoFHandlerOrdering(0, 0);
1491  }
1492  else if (this->GetType() == "gradient"
1493  || this->GetType() == "hessian_inverse"
1494  || this->GetType() == "hessian")
1495  {
1496  GetSpaceTimeHandler()->SetDoFHandlerOrdering(0, 0);
1497  }
1498  else
1499  {
1500  throw DOpEException(
1501  "problem_type_ : " + this->GetType() + " not implemented!",
1502  "OptProblemContainer::SetType");
1503  }
1504  }
1505  }
1506 #endif
1507  }
1508  }
1509 
1510  /******************************************************/
1511 
1512  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
1513  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
1514  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
1515  template<int, int> class DH>
1516  template<typename DATACONTAINER>
1517  double
1518  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD,
1519  CONSTRAINTS, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::ElementFunctional(
1520  const DATACONTAINER& edc)
1521  {
1522 
1523  if (this->GetType() == "cost_functional")
1524  {
1525  // state values in quadrature points
1526  return GetFunctional()->ElementValue(edc);
1527  }
1528  else if (this->GetType() == "aux_functional")
1529  {
1530  // state values in quadrature points
1531  return aux_functionals_[this->GetTypeNum()]->ElementValue(edc);
1532  }
1533  else if (this->GetType() == "functional_for_ee")
1534  { //TODO ist das hier korrekt? Sollten wir eigentlich nicht benoetigen.
1535  return aux_functionals_[functional_for_ee_num_]->ElementValue(edc);
1536  }
1537  else if (this->GetType().find("constraints") != std::string::npos)
1538  {
1539  return GetConstraints()->ElementValue(edc);
1540  }
1541  else
1542  {
1543  throw DOpEException("Not implemented",
1544  "OptProblemContainer::ElementFunctional");
1545  }
1546  }
1547 
1548  /******************************************************/
1549  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
1550  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
1551  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
1552  template<int, int> class DH>
1553  double
1554  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
1555  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::PointFunctional(
1556  const std::map<std::string, const dealii::Vector<double>*> &param_values,
1557  const std::map<std::string, const VECTOR*> &domain_values)
1558  {
1559  if (this->GetType() == "cost_functional")
1560  {
1561  // state values in quadrature points
1562  return GetFunctional()->PointValue(
1563  this->GetSpaceTimeHandler()->GetControlDoFHandler(),
1564  this->GetSpaceTimeHandler()->GetStateDoFHandler(), param_values,
1565  domain_values);
1566 
1567  } //endif cost_functional
1568  else if (this->GetType() == "aux_functional")
1569  {
1570  // state values in quadrature points
1571  return aux_functionals_[this->GetTypeNum()]->PointValue(
1572  this->GetSpaceTimeHandler()->GetControlDoFHandler(),
1573  this->GetSpaceTimeHandler()->GetStateDoFHandler(), param_values,
1574  domain_values);
1575 
1576  } //endif aux_functional
1577  else if (this->GetType() == "functional_for_ee")
1578  {
1579  //TODO ist das hier korrekt? Sollten wir eigentlich nicht benoetigen.
1580  return aux_functionals_[functional_for_ee_num_]->PointValue(
1581  this->GetSpaceTimeHandler()->GetControlDoFHandler(),
1582  this->GetSpaceTimeHandler()->GetStateDoFHandler(), param_values,
1583  domain_values);
1584  } //endif functional_for_ee
1585  else if (this->GetType().find("constraints") != std::string::npos)
1586  {
1587  return GetConstraints()->PointValue(
1588  this->GetSpaceTimeHandler()->GetControlDoFHandler(),
1589  this->GetSpaceTimeHandler()->GetStateDoFHandler(), param_values,
1590  domain_values);
1591 
1592  } //endif constraints
1593  else
1594  {
1595  throw DOpEException("Not implemented",
1596  "OptProblemContainer::PointFunctional");
1597  }
1598  }
1599 
1600  /******************************************************/
1601 
1602  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
1603  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
1604  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
1605  template<int, int> class DH>
1606  template<typename FACEDATACONTAINER>
1607  double
1608  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD,
1609  CONSTRAINTS, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::BoundaryFunctional(
1610  const FACEDATACONTAINER& fdc)
1611  {
1612  if (this->GetType() == "cost_functional")
1613  {
1614  // state values in quadrature points
1615  return GetFunctional()->BoundaryValue(fdc);
1616  }
1617  else if (this->GetType() == "aux_functional")
1618  {
1619  // state values in quadrature points
1620  return aux_functionals_[this->GetTypeNum()]->BoundaryValue(fdc);
1621  }
1622  else if (this->GetType() == "functional_for_ee")
1623  //TODO ist das hier korrekt? Sollten wir eigentlich nicht benoetigen.
1624  {
1625  return aux_functionals_[functional_for_ee_num_]->BoundaryValue(fdc);
1626  }
1627  else if (this->GetType().find("constraints") != std::string::npos)
1628  {
1629  return GetConstraints()->BoundaryValue(fdc);
1630  }
1631  else
1632  {
1633  throw DOpEException("Not implemented",
1634  "OptProblemContainer::BoundaryFunctional");
1635  }
1636  }
1637 
1638  /******************************************************/
1639 
1640  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
1641  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
1642  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
1643  template<int, int> class DH>
1644  template<typename FACEDATACONTAINER>
1645  double
1646  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD,
1647  CONSTRAINTS, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::FaceFunctional(
1648  const FACEDATACONTAINER& fdc)
1649  {
1650  if (this->GetType() == "cost_functional")
1651  {
1652  // state values in quadrature points
1653  return GetFunctional()->FaceValue(fdc);
1654  }
1655  else if (this->GetType() == "aux_functional")
1656  {
1657  // state values in quadrature points
1658  return aux_functionals_[this->GetTypeNum()]->FaceValue(fdc);
1659  }
1660  else if (this->GetType() == "functional_for_ee")
1661  //TODO ist das hier korrekt? Sollten wir eigentlich nicht benoetigen.
1662  {
1663  return aux_functionals_[functional_for_ee_num_]->FaceValue(fdc);
1664  }
1665  else
1666  {
1667  throw DOpEException("Not implemented",
1668  "OptProblemContainer::FaceFunctional");
1669  }
1670  }
1671 
1672  /******************************************************/
1673 
1674  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
1675  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
1676  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
1677  template<int, int> class DH>
1678  double
1679  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
1680  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::AlgebraicFunctional(
1681  const std::map<std::string, const dealii::Vector<double>*> &param_values,
1682  const std::map<std::string, const VECTOR*> &domain_values)
1683  {
1684  if (this->GetType() == "cost_functional")
1685  {
1686  // state values in quadrature points
1687  return GetFunctional()->AlgebraicValue(param_values, domain_values);
1688  }
1689  else if (this->GetType() == "aux_functional")
1690  {
1691  // state values in quadrature points
1692  return aux_functionals_[this->GetTypeNum()]->AlgebraicValue(
1693  param_values, domain_values);
1694  }
1695  else if (this->GetType() == "functional_for_ee")
1696  //TODO ist das hier korrekt? Sollten wir eigentlich nicht benoetigen.
1697  {
1698  return aux_functionals_[functional_for_ee_num_]->AlgebraicValue(
1699  param_values, domain_values);
1700  }
1701  else
1702  {
1703  throw DOpEException("Not implemented",
1704  "OptProblemContainer::AlgebraicFunctional");
1705  }
1706  }
1707 
1708  /******************************************************/
1709 
1710  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
1711  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
1712  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
1713  template<int, int> class DH>
1714  template<typename DATACONTAINER>
1715  void
1716  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD,
1717  CONSTRAINTS, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::ElementEquation(
1718  const DATACONTAINER& edc, dealii::Vector<double> &local_vector,
1719  double scale, double scale_ico)
1720  {
1721 
1722  if (this->GetType() == "state")
1723  {
1724  // state values in quadrature points
1725  this->GetPDE().ElementEquation(edc, local_vector, scale*interval_length_, scale_ico*interval_length_);
1726  }
1727  else if ((this->GetType() == "adjoint")
1728  || (this->GetType() == "adjoint_for_ee"))
1729  {
1730  // state values in quadrature points
1731  this->GetPDE().ElementEquation_U(edc, local_vector, scale*interval_length_,
1732  scale_ico*interval_length_);
1733  }
1734  else if (this->GetType() == "adjoint_hessian")
1735  {
1736  // state values in quadrature points
1737  this->GetPDE().ElementEquation_UTT(edc, local_vector, scale*interval_length_,
1738  scale_ico*interval_length_);
1739  }
1740  else if (this->GetType() == "tangent")
1741  {
1742  // state values in quadrature points
1743  this->GetPDE().ElementEquation_UT(edc, local_vector, scale*interval_length_,
1744  scale_ico*interval_length_);
1745  }
1746  else if ((this->GetType() == "gradient")
1747  || (this->GetType() == "hessian"))
1748  {
1749  // control values in quadrature points
1750  this->GetPDE().ControlElementEquation(edc, local_vector, scale*c_interval_length_);
1751  }
1752  else
1753  {
1754  throw DOpEException("Not implemented",
1755  "OptProblemContainer::ElementEquation");
1756  }
1757  }
1758 
1759  /******************************************************/
1760 
1761  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
1762  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
1763  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
1764  template<int, int> class DH>
1765  void
1766  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
1767  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::AlgebraicResidual(
1768  VECTOR& residual,
1769  const std::map<std::string, const dealii::Vector<double>*> &param_values,
1770  const std::map<std::string, const VECTOR*> &domain_values)
1771  {
1772  if (this->GetType() == "gradient")
1773  {
1774  // state values in quadrature points
1775  return GetFunctional()->AlgebraicGradient_Q(residual, param_values,
1776  domain_values);
1777  }
1778  else
1779  {
1780  throw DOpEException("Not implemented",
1781  "OptProblemContainer::AlgebraicFunctional");
1782  }
1783  }
1784 
1785  /******************************************************/
1786 
1787  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
1788  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
1789  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
1790  template<int, int> class DH>
1791  template<typename DATACONTAINER>
1792  void
1793  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD,
1794  CONSTRAINTS, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::ElementTimeEquation(
1795  const DATACONTAINER& edc, dealii::Vector<double> &local_vector,
1796  double scale)
1797  {
1798 
1799  if (this->GetType() == "state")
1800  {
1801  this->GetPDE().ElementTimeEquation(edc, local_vector, scale);
1802  }
1803  else if (this->GetType() == "adjoint"
1804  || this->GetType() == "adjoint_for_ee")
1805  {
1806  this->GetPDE().ElementTimeEquation_U(edc, local_vector, scale);
1807  }
1808  else if (this->GetType() == "adjoint_hessian")
1809  {
1810  this->GetPDE().ElementTimeEquation_UTT(edc, local_vector, scale);
1811  }
1812  else if (this->GetType() == "tangent")
1813  {
1814  this->GetPDE().ElementTimeEquation_UT(edc, local_vector, scale);
1815  }
1816  else if ((this->GetType() == "gradient")
1817  || (this->GetType() == "hessian"))
1818  {
1819  throw DOpEException("Not implemented",
1820  "OptProblemContainer::ElementTimeEquation");
1821  }
1822  else
1823  {
1824  throw DOpEException("Not implemented",
1825  "OptProblemContainer::ElementTimeEquation");
1826  }
1827  }
1828 
1829  /******************************************************/
1830 
1831  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
1832  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
1833  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
1834  template<int, int> class DH>
1835  template<typename DATACONTAINER>
1836  void
1837  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD,
1838  CONSTRAINTS, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::ElementTimeEquationExplicit(
1839  const DATACONTAINER& edc, dealii::Vector<double> &local_vector,
1840  double scale)
1841  {
1842 
1843  if (this->GetType() == "state")
1844  {
1845  this->GetPDE().ElementTimeEquationExplicit(edc, local_vector,
1846  scale);
1847  }
1848  else if (this->GetType() == "adjoint"
1849  || this->GetType() == "adjoint_for_ee")
1850  {
1851  this->GetPDE().ElementTimeEquationExplicit_U(edc, local_vector,
1852  scale);
1853  }
1854  else if (this->GetType() == "adjoint_hessian")
1855  {
1856  this->GetPDE().ElementTimeEquationExplicit_UTT(edc, local_vector,
1857  scale);
1858  }
1859  else if (this->GetType() == "tangent")
1860  {
1861  this->GetPDE().ElementTimeEquationExplicit_UT(edc, local_vector,
1862  scale);
1863  }
1864  else if ((this->GetType() == "gradient")
1865  || (this->GetType() == "hessian"))
1866  {
1867  throw DOpEException("Not implemented",
1868  "OptProblemContainer::ElementTimeEquationExplicit");
1869  }
1870  else
1871  {
1872  throw DOpEException("Not implemented",
1873  "OptProblemContainer::ElementTimeEquationExplicit");
1874  }
1875  }
1876 
1877  /******************************************************/
1878 
1879  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
1880  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
1881  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
1882  template<int, int> class DH>
1883  template<typename FACEDATACONTAINER>
1884  void
1885  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD,
1886  CONSTRAINTS, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::FaceEquation(
1887  const FACEDATACONTAINER& fdc,
1888  dealii::Vector<double> &local_vector, double scale,
1889  double scale_ico)
1890  {
1891  if (this->GetType() == "state")
1892  {
1893  // state values in quadrature points
1894  this->GetPDE().FaceEquation(fdc, local_vector, scale*interval_length_, scale_ico*interval_length_);
1895  }
1896  else if (this->GetType() == "adjoint"
1897  || this->GetType() == "adjoint_for_ee")
1898  {
1899  // state values in quadrature points
1900  this->GetPDE().FaceEquation_U(fdc, local_vector, scale*interval_length_,
1901  scale_ico*interval_length_);
1902  }
1903  else if (this->GetType() == "adjoint_hessian")
1904  {
1905  // state values in quadrature points
1906  this->GetPDE().FaceEquation_UTT(fdc, local_vector, scale*interval_length_,
1907  scale_ico*interval_length_);
1908  }
1909  else if (this->GetType() == "tangent")
1910  {
1911  // state values in quadrature points
1912  this->GetPDE().FaceEquation_UT(fdc, local_vector, scale*interval_length_,
1913  scale_ico*interval_length_);
1914  }
1915 // else if ((this->GetType() == "gradient") || (this->GetType() == "hessian"))
1916 // {
1917 // // control values in quadrature points
1918 // this->GetPDE().ControlFaceEquation(fdc, local_vector, scale);
1919 // }
1920  else
1921  {
1922  throw DOpEException("Not implemented",
1923  "OptProblemContainer::FaceEquation");
1924  }
1925  }
1926 
1927  /******************************************************/
1928 
1929  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
1930  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
1931  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
1932  template<int, int> class DH>
1933  template<typename FACEDATACONTAINER>
1934  void
1935  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD,
1936  CONSTRAINTS, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::InterfaceEquation(
1937  const FACEDATACONTAINER& fdc,
1938  dealii::Vector<double> &local_vector, double scale,
1939  double scale_ico)
1940  {
1941  if (this->GetType() == "state")
1942  {
1943  this->GetPDE().InterfaceEquation(fdc, local_vector, scale*interval_length_,
1944  scale_ico*interval_length_);
1945  }
1946  else if (this->GetType() == "adjoint"
1947  || this->GetType() == "adjoint_for_ee")
1948  {
1949  // state values in quadrature points
1950  this->GetPDE().InterfaceEquation_U(fdc, local_vector, scale*interval_length_,
1951  scale_ico*interval_length_);
1952  }
1953  else
1954  {
1955  throw DOpEException("Not implemented",
1956  "OptProblemContainer::InterfaceEquation");
1957  }
1958  }
1959  /******************************************************/
1960 
1961  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
1962  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
1963  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
1964  template<int, int> class DH>
1965  template<typename FACEDATACONTAINER>
1966  void
1967  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD,
1968  CONSTRAINTS, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::BoundaryEquation(
1969  const FACEDATACONTAINER& fdc,
1970  dealii::Vector<double> &local_vector, double scale,
1971  double scale_ico)
1972  {
1973  if (this->GetType() == "state")
1974  {
1975  // state values in quadrature points
1976  this->GetPDE().BoundaryEquation(fdc, local_vector, scale*interval_length_,
1977  scale_ico*interval_length_);
1978  }
1979  else if (this->GetType() == "adjoint"
1980  || this->GetType() == "adjoint_for_ee")
1981  {
1982  // state values in quadrature points
1983  this->GetPDE().BoundaryEquation_U(fdc, local_vector, scale*interval_length_,
1984  scale_ico*interval_length_);
1985  }
1986  else if (this->GetType() == "adjoint_hessian")
1987  {
1988  // state values in quadrature points
1989  this->GetPDE().BoundaryEquation_UTT(fdc, local_vector, scale*interval_length_,
1990  scale_ico*interval_length_);
1991  }
1992  else if (this->GetType() == "tangent")
1993  {
1994  // state values in quadrature points
1995  this->GetPDE().BoundaryEquation_UT(fdc, local_vector, scale*interval_length_,
1996  scale_ico*interval_length_);
1997  }
1998 // else if ((this->GetType() == "gradient") || (this->GetType() == "hessian"))
1999 // {
2000 // // control values in quadrature points
2001 // this->GetPDE().ControlBoundaryEquation(fdc, local_vector, scale);
2002 // }
2003  else
2004  {
2005  throw DOpEException("Not implemented",
2006  "OptProblemContainer::ElementBoundaryEquation");
2007  }
2008  }
2009 
2010  /******************************************************/
2011 
2012  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
2013  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
2014  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
2015  template<int, int> class DH>
2016  template<typename DATACONTAINER>
2017  void
2018  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD,
2019  CONSTRAINTS, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::ElementRhs(
2020  const DATACONTAINER& edc, dealii::Vector<double> &local_vector,
2021  double scale)
2022  {
2023  if (this->GetType() == "state")
2024  {
2025  // state values in quadrature points
2026  this->GetPDE().ElementRightHandSide(edc, local_vector, scale*interval_length_);
2027  }
2028  else if (this->GetType() == "adjoint")
2029  {
2030  // state values in quadrature points
2031  if(GetFunctional()->NeedTime())
2032  {
2033  if (GetFunctional()->GetType().find("domain") != std::string::npos)
2034  {
2035  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos)
2036  {
2037  GetFunctional()->ElementValue_U(edc, local_vector, scale*interval_length_);
2038  }
2039  else // Otherwise always local if(GetFunctional()->GetType().find("timelocal") != std::string::npos)
2040  {
2041  GetFunctional()->ElementValue_U(edc, local_vector, scale);
2042  }
2043 
2044  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos && GetFunctional()->GetType().find("timelocal") != std::string::npos)
2045  {
2046  throw DOpEException("Conflicting functional types: "+ GetFunctional()->GetType(),
2047  "OptProblemContainer::ElementRhs");
2048  }
2049  }
2050  }
2051  }
2052  else if (this->GetType() == "adjoint_for_ee")
2053  {
2054  //values of the derivative of the functional for error estimation
2055  if (aux_functionals_[functional_for_ee_num_]->NeedTime())
2056  {
2057  if (aux_functionals_[functional_for_ee_num_]->GetType().find("domain") != std::string::npos)
2058  {
2059  if(aux_functionals_[functional_for_ee_num_]->GetType().find("timedistributed") != std::string::npos)
2060  {
2061  aux_functionals_[functional_for_ee_num_]->ElementValue_U(edc,
2062  local_vector, scale*interval_length_);
2063  }
2064  else
2065  {
2066  aux_functionals_[functional_for_ee_num_]->ElementValue_U(edc, local_vector, scale);
2067  }
2068  if(aux_functionals_[functional_for_ee_num_]->GetType().find("timedistributed") != std::string::npos && aux_functionals_[functional_for_ee_num_]->GetType().find("timelocal") != std::string::npos)
2069  {
2070  throw DOpEException("Conflicting functional types: "+ aux_functionals_[functional_for_ee_num_]->GetType(),
2071  "OptProblemContainer::ElementRhs");
2072  }
2073  }
2074  }
2075  }
2076  else if (this->GetType() == "tangent")
2077  {
2078  // state values in quadrature points
2079  scale *= -1;
2080  this->GetPDE().ElementEquation_QT(edc, local_vector, scale*interval_length_, scale*interval_length_);
2081  }
2082  else if (this->GetType() == "adjoint_hessian")
2083  {
2084  // state values in quadrature points
2085  if(GetFunctional()->NeedTime())
2086  {
2087  if (GetFunctional()->GetType().find("domain") != std::string::npos)
2088  {
2089  if (GetFunctional()->GetType().find("domain") != std::string::npos)
2090  {
2091  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos)
2092  {
2093  GetFunctional()->ElementValue_UU(edc, local_vector, scale*interval_length_);
2094  GetFunctional()->ElementValue_QU(edc, local_vector, scale*interval_length_);
2095  }
2096  else
2097  {
2098  GetFunctional()->ElementValue_UU(edc, local_vector, scale);
2099  GetFunctional()->ElementValue_QU(edc, local_vector, scale);
2100  }
2101  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos && GetFunctional()->GetType().find("timelocal") != std::string::npos)
2102  {
2103  throw DOpEException("Conflicting functional types: "+ GetFunctional()->GetType(),
2104  "OptProblemContainer::ElementRhs");
2105  }
2106  }
2107  }
2108  }
2109  scale *= -1;
2110  this->GetPDE().ElementEquation_UU(edc, local_vector, scale*interval_length_, scale*interval_length_);
2111  this->GetPDE().ElementEquation_QU(edc, local_vector, scale*interval_length_, scale*interval_length_);
2112  //TODO: make some example where this realy matters to check if this is right
2113  this->GetPDE().ElementTimeEquationExplicit_UU(edc, local_vector,
2114  scale);
2115  }
2116  else if (this->GetType() == "gradient")
2117  {
2118  if (GetSpaceTimeHandler()->GetControlType()
2119  == DOpEtypes::ControlType::initial && initial_)
2120  {
2121  this->GetPDE().Init_ElementRhs_Q(edc, local_vector, scale);
2122  }
2123  // state values in quadrature points
2124  if(GetFunctional()->NeedTime())
2125  {
2126  if (GetFunctional()->GetType().find("domain") != std::string::npos)
2127  {
2128  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos)
2129  {
2130  GetFunctional()->ElementValue_Q(edc, local_vector, scale*interval_length_);
2131  }
2132  else
2133  {
2134  GetFunctional()->ElementValue_Q(edc, local_vector, scale);
2135  }
2136  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos && GetFunctional()->GetType().find("timelocal") != std::string::npos)
2137  {
2138  throw DOpEException("Conflicting functional types: "+ GetFunctional()->GetType(),
2139  "OptProblemContainer::ElementRhs");
2140  }
2141  }
2142  }
2143  scale *= -1;
2144  this->GetPDE().ElementEquation_Q(edc, local_vector, scale*interval_length_, scale*interval_length_);
2145  }
2146  else if (this->GetType() == "hessian")
2147  {
2148  if (GetSpaceTimeHandler()->GetControlType()
2149  == DOpEtypes::ControlType::initial && initial_)
2150  {
2151  this->GetPDE().Init_ElementRhs_QTT(edc, local_vector, scale);
2152  this->GetPDE().Init_ElementRhs_QQ(edc, local_vector, scale);
2153  }
2154  if(GetFunctional()->NeedTime())
2155  {
2156  if (GetFunctional()->GetType().find("domain") != std::string::npos)
2157  {
2158  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos)
2159  {
2160  GetFunctional()->ElementValue_QQ(edc, local_vector, scale*interval_length_);
2161  GetFunctional()->ElementValue_UQ(edc, local_vector, scale*interval_length_);
2162  }
2163  else
2164  {
2165  GetFunctional()->ElementValue_QQ(edc, local_vector, scale);
2166  GetFunctional()->ElementValue_UQ(edc, local_vector, scale);
2167  }
2168  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos && GetFunctional()->GetType().find("timelocal") != std::string::npos)
2169  {
2170  throw DOpEException("Conflicting functional types: "+ GetFunctional()->GetType(),
2171  "OptProblemContainer::ElementRhs");
2172  }
2173  }
2174  }
2175 
2176  scale *= -1;
2177  this->GetPDE().ElementEquation_QTT(edc, local_vector, scale*interval_length_, scale*interval_length_);
2178  this->GetPDE().ElementEquation_UQ(edc, local_vector, scale*interval_length_, scale*interval_length_);
2179  this->GetPDE().ElementEquation_QQ(edc, local_vector, scale*interval_length_, scale*interval_length_);
2180  }
2181  else if (this->GetType() == "global_constraint_gradient")
2182  {
2183  assert(interval_length_==1.);
2184  GetConstraints()->ElementValue_Q(edc, local_vector, scale);
2185  }
2186  else if (this->GetType() == "global_constraint_hessian")
2187  {
2188  assert(interval_length_==1.);
2189  GetConstraints()->ElementValue_QQ(edc, local_vector, scale);
2190  }
2191  else
2192  {
2193  throw DOpEException("Not implemented",
2194  "OptProblemContainer::ElementRhs");
2195  }
2196 
2197  }
2198 
2199  /******************************************************/
2200 
2201  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
2202  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
2203  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
2204  template<int, int> class DH>
2205  void
2206  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
2207  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::PointRhs(
2208  const std::map<std::string, const dealii::Vector<double>*> &param_values,
2209  const std::map<std::string, const VECTOR*> &domain_values,
2210  VECTOR& rhs_vector, double scale)
2211  {
2212 
2213  if (this->GetType() == "adjoint")
2214  {
2215  // state values in quadrature points
2216  if(GetFunctional()->NeedTime())
2217  {
2218  if (GetFunctional()->GetType().find("point") != std::string::npos)
2219  {
2220  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos)
2221  {
2222  GetFunctional()->PointValue_U(
2223  this->GetSpaceTimeHandler()->GetControlDoFHandler(),
2224  this->GetSpaceTimeHandler()->GetStateDoFHandler(), param_values,
2225  domain_values, rhs_vector, scale*interval_length_);
2226  }
2227  else
2228  {
2229  GetFunctional()->PointValue_U(
2230  this->GetSpaceTimeHandler()->GetControlDoFHandler(),
2231  this->GetSpaceTimeHandler()->GetStateDoFHandler(), param_values,
2232  domain_values, rhs_vector, scale);
2233  }
2234  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos && GetFunctional()->GetType().find("timelocal") != std::string::npos)
2235  {
2236  throw DOpEException("Conflicting functional types: "+ GetFunctional()->GetType(),
2237  "OptProblemContainer::PointRhs");
2238  }
2239  }
2240  }
2241  }
2242  else if (this->GetType() == "adjoint_for_ee")
2243  {
2244  //values of the derivative of the functional for error estimation
2245  if (aux_functionals_[functional_for_ee_num_]->NeedTime())
2246  {
2247  if (aux_functionals_[functional_for_ee_num_]->GetType().find("point") != std::string::npos)
2248  {
2249  if(aux_functionals_[functional_for_ee_num_]->GetType().find("timedistributed") != std::string::npos)
2250  {
2251  aux_functionals_[functional_for_ee_num_]->PointValue_U(
2252  this->GetSpaceTimeHandler()->GetControlDoFHandler(),
2253  this->GetSpaceTimeHandler()->GetStateDoFHandler(), param_values,
2254  domain_values, rhs_vector, scale*interval_length_);
2255  }
2256  else
2257  {
2258  aux_functionals_[functional_for_ee_num_]->PointValue_U(
2259  this->GetSpaceTimeHandler()->GetControlDoFHandler(),
2260  this->GetSpaceTimeHandler()->GetStateDoFHandler(), param_values,
2261  domain_values, rhs_vector, scale);
2262  }
2263  if(aux_functionals_[functional_for_ee_num_]->GetType().find("timedistributed") != std::string::npos && aux_functionals_[functional_for_ee_num_]->GetType().find("timelocal") != std::string::npos)
2264  {
2265  throw DOpEException("Conflicting functional types: "+ aux_functionals_[functional_for_ee_num_]->GetType(),
2266  "OptProblemContainer::PointRhs");
2267  }
2268  }
2269  }
2270  }
2271  else if (this->GetType() == "adjoint_hessian")
2272  {
2273  // state values in quadrature points
2274  if(GetFunctional()->NeedTime())
2275  {
2276  if (GetFunctional()->GetType().find("point") != std::string::npos)
2277  {
2278  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos)
2279  {
2280  GetFunctional()->PointValue_UU(
2281  this->GetSpaceTimeHandler()->GetControlDoFHandler(),
2282  this->GetSpaceTimeHandler()->GetStateDoFHandler(), param_values,
2283  domain_values, rhs_vector, scale*interval_length_);
2284  GetFunctional()->PointValue_QU(
2285  this->GetSpaceTimeHandler()->GetControlDoFHandler(),
2286  this->GetSpaceTimeHandler()->GetStateDoFHandler(), param_values,
2287  domain_values, rhs_vector, scale*interval_length_);
2288  }
2289  else
2290  {
2291  GetFunctional()->PointValue_UU(
2292  this->GetSpaceTimeHandler()->GetControlDoFHandler(),
2293  this->GetSpaceTimeHandler()->GetStateDoFHandler(), param_values,
2294  domain_values, rhs_vector, scale);
2295  GetFunctional()->PointValue_QU(
2296  this->GetSpaceTimeHandler()->GetControlDoFHandler(),
2297  this->GetSpaceTimeHandler()->GetStateDoFHandler(), param_values,
2298  domain_values, rhs_vector, scale);
2299  }
2300  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos && GetFunctional()->GetType().find("timelocal") != std::string::npos)
2301  {
2302  throw DOpEException("Conflicting functional types: "+ GetFunctional()->GetType(),
2303  "OptProblemContainer::PointRhs");
2304  }
2305  }
2306  }
2307  }
2308  else if (this->GetType() == "gradient")
2309  {
2310  // state values in quadrature points
2311  if(GetFunctional()->NeedTime())
2312  {
2313  if (GetFunctional()->GetType().find("point") != std::string::npos)
2314  {
2315  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos)
2316  {
2317  GetFunctional()->PointValue_Q(
2318  this->GetSpaceTimeHandler()->GetControlDoFHandler(),
2319  this->GetSpaceTimeHandler()->GetStateDoFHandler(), param_values,
2320  domain_values, rhs_vector, scale*interval_length_);
2321  }
2322  else
2323  {
2324  GetFunctional()->PointValue_Q(
2325  this->GetSpaceTimeHandler()->GetControlDoFHandler(),
2326  this->GetSpaceTimeHandler()->GetStateDoFHandler(), param_values,
2327  domain_values, rhs_vector, scale);
2328  }
2329  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos && GetFunctional()->GetType().find("timelocal") != std::string::npos)
2330  {
2331  throw DOpEException("Conflicting functional types: "+ GetFunctional()->GetType(),
2332  "OptProblemContainer::PointRhs");
2333  }
2334  }
2335  }
2336  }
2337  else if (this->GetType() == "hessian")
2338  {
2339  // state values in quadrature points
2340  if(GetFunctional()->NeedTime())
2341  {
2342  if (GetFunctional()->GetType().find("point") != std::string::npos)
2343  {
2344  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos)
2345  {
2346  GetFunctional()->PointValue_QQ(
2347  this->GetSpaceTimeHandler()->GetControlDoFHandler(),
2348  this->GetSpaceTimeHandler()->GetStateDoFHandler(), param_values,
2349  domain_values, rhs_vector, scale*interval_length_);
2350  GetFunctional()->PointValue_UQ(
2351  this->GetSpaceTimeHandler()->GetControlDoFHandler(),
2352  this->GetSpaceTimeHandler()->GetStateDoFHandler(), param_values,
2353  domain_values, rhs_vector, scale*interval_length_);
2354  }
2355  else
2356  {
2357  GetFunctional()->PointValue_QQ(
2358  this->GetSpaceTimeHandler()->GetControlDoFHandler(),
2359  this->GetSpaceTimeHandler()->GetStateDoFHandler(), param_values,
2360  domain_values, rhs_vector, scale);
2361  GetFunctional()->PointValue_UQ(
2362  this->GetSpaceTimeHandler()->GetControlDoFHandler(),
2363  this->GetSpaceTimeHandler()->GetStateDoFHandler(), param_values,
2364  domain_values, rhs_vector, scale);
2365  }
2366  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos && GetFunctional()->GetType().find("timelocal") != std::string::npos)
2367  {
2368  throw DOpEException("Conflicting functional types: "+ GetFunctional()->GetType(),
2369  "OptProblemContainer::PointRhs");
2370  }
2371  }
2372  }
2373  }
2374  else
2375  {
2376  throw DOpEException("Not implemented", "OptProblem::PointRhs");
2377  }
2378 
2379  }
2380 
2381  /******************************************************/
2382 
2383  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
2384  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
2385  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
2386  template<int, int> class DH>
2387  template<typename FACEDATACONTAINER>
2388  void
2389  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD,
2390  CONSTRAINTS, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::FaceRhs(
2391  const FACEDATACONTAINER& fdc,
2392  dealii::Vector<double> &local_vector, double scale)
2393  {
2394 
2395  if (this->GetType() == "state")
2396  {
2397  // state values in face quadrature points
2398  this->GetPDE().FaceRightHandSide(fdc, local_vector, scale*interval_length_);
2399  }
2400  else if (this->GetType() == "adjoint")
2401  {
2402  // state values in quadrature points
2403  if(GetFunctional()->NeedTime())
2404  {
2405  if (GetFunctional()->GetType().find("face") != std::string::npos)
2406  {
2407  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos)
2408  {
2409  GetFunctional()->FaceValue_U(fdc, local_vector, scale*interval_length_);
2410  }
2411  else
2412  {
2413  GetFunctional()->FaceValue_U(fdc, local_vector, scale);
2414  }
2415 
2416  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos && GetFunctional()->GetType().find("timelocal") != std::string::npos)
2417  {
2418  throw DOpEException("Conflicting functional types: "+ GetFunctional()->GetType(),
2419  "OptProblemContainer::FaceRhs");
2420  }
2421  }
2422  }
2423  }
2424  else if (this->GetType() == "adjoint_for_ee")
2425  {
2426  //values of the derivative of the functional for error estimation
2427  if (aux_functionals_[functional_for_ee_num_]->NeedTime())
2428  {
2429  if (aux_functionals_[functional_for_ee_num_]->GetType().find("face")
2430  != std::string::npos)
2431  {
2432  if(aux_functionals_[functional_for_ee_num_]->GetType().find("timedistributed") != std::string::npos)
2433  {
2434  aux_functionals_[functional_for_ee_num_]->FaceValue_U(fdc,
2435  local_vector, scale*interval_length_);
2436  }
2437  else
2438  {
2439  aux_functionals_[functional_for_ee_num_]->FaceValue_U(fdc,
2440  local_vector, scale);
2441  }
2442  if(aux_functionals_[functional_for_ee_num_]->GetType().find("timedistributed") != std::string::npos && aux_functionals_[functional_for_ee_num_]->GetType().find("timelocal") != std::string::npos)
2443  {
2444  throw DOpEException("Conflicting functional types: "+ aux_functionals_[functional_for_ee_num_]->GetType(),
2445  "OptProblemContainer::FaceRhs");
2446  }
2447  }
2448  }
2449 
2450  }
2451  else if (this->GetType() == "tangent")
2452  {
2453  // state values in quadrature points
2454  scale *= -1;
2455  this->GetPDE().FaceEquation_QT(fdc, local_vector, scale*interval_length_, scale*interval_length_);
2456  }
2457  else if (this->GetType() == "adjoint_hessian")
2458  {
2459  // state values in quadrature points
2460  if(GetFunctional()->NeedTime())
2461  {
2462  if (GetFunctional()->GetType().find("face") != std::string::npos)
2463  {
2464  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos)
2465  {
2466  GetFunctional()->FaceValue_UU(fdc, local_vector, scale*interval_length_);
2467  GetFunctional()->FaceValue_QU(fdc, local_vector, scale*interval_length_);
2468  }
2469  else
2470  {
2471  GetFunctional()->FaceValue_UU(fdc, local_vector, scale);
2472  GetFunctional()->FaceValue_QU(fdc, local_vector, scale);
2473  }
2474  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos && GetFunctional()->GetType().find("timelocal") != std::string::npos)
2475  {
2476  throw DOpEException("Conflicting functional types: "+ GetFunctional()->GetType(),
2477  "OptProblemContainer::FaceRhs");
2478  }
2479  }
2480  }
2481 
2482  scale *= -1;
2483  this->GetPDE().FaceEquation_UU(fdc, local_vector, scale*interval_length_, scale*interval_length_);
2484  this->GetPDE().FaceEquation_QU(fdc, local_vector, scale*interval_length_, scale*interval_length_);
2485  }
2486  else if (this->GetType() == "gradient")
2487  {
2488  // state values in quadrature points
2489  if(GetFunctional()->NeedTime())
2490  {
2491  if (GetFunctional()->GetType().find("face") != std::string::npos)
2492  {
2493  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos)
2494  {
2495  GetFunctional()->FaceValue_Q(fdc, local_vector, scale*interval_length_);
2496  }
2497  else
2498  {
2499  GetFunctional()->FaceValue_Q(fdc, local_vector, scale);
2500  }
2501  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos && GetFunctional()->GetType().find("timelocal") != std::string::npos)
2502  {
2503  throw DOpEException("Conflicting functional types: "+ GetFunctional()->GetType(),
2504  "OptProblemContainer::FaceRhs");
2505  }
2506  }
2507  }
2508 
2509  scale *= -1;
2510  this->GetPDE().FaceEquation_Q(fdc, local_vector, scale*interval_length_, scale*interval_length_);
2511  }
2512  else if (this->GetType() == "hessian")
2513  {
2514  // state values in quadrature points
2515  if(GetFunctional()->NeedTime())
2516  {
2517  if (GetFunctional()->GetType().find("face") != std::string::npos)
2518  {
2519  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos)
2520  {
2521  GetFunctional()->FaceValue_QQ(fdc, local_vector, scale*interval_length_);
2522  GetFunctional()->FaceValue_UQ(fdc, local_vector, scale*interval_length_);
2523  }
2524  else
2525  {
2526  GetFunctional()->FaceValue_QQ(fdc, local_vector, scale);
2527  GetFunctional()->FaceValue_UQ(fdc, local_vector, scale);
2528  }
2529  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos && GetFunctional()->GetType().find("timelocal") != std::string::npos)
2530  {
2531  throw DOpEException("Conflicting functional types: "+ GetFunctional()->GetType(),
2532  "OptProblemContainer::FaceRhs");
2533  }
2534  }
2535  }
2536 
2537  scale *= -1;
2538  this->GetPDE().FaceEquation_QTT(fdc, local_vector, scale*interval_length_, scale*interval_length_);
2539  this->GetPDE().FaceEquation_UQ(fdc, local_vector, scale*interval_length_, scale*interval_length_);
2540  this->GetPDE().FaceEquation_QQ(fdc, local_vector, scale*interval_length_, scale*interval_length_);
2541  }
2542  else
2543  {
2544  throw DOpEException("Not implemented",
2545  "OptProblemContainer::FaceRhs");
2546  }
2547 
2548  }
2549 
2550  /******************************************************/
2551 
2552  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
2553  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
2554  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
2555  template<int, int> class DH>
2556  template<typename FACEDATACONTAINER>
2557  void
2558  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD,
2559  CONSTRAINTS, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::BoundaryRhs(
2560  const FACEDATACONTAINER& fdc,
2561  dealii::Vector<double> &local_vector, double scale)
2562  {
2563 
2564  if (this->GetType() == "state")
2565  {
2566  // state values in face quadrature points
2567  this->GetPDE().BoundaryRightHandSide(fdc, local_vector, scale*interval_length_);
2568  }
2569  else if (this->GetType() == "adjoint")
2570  {
2571  // state values in quadrature points
2572  if(GetFunctional()->NeedTime())
2573  {
2574  if (GetFunctional()->GetType().find("boundary") != std::string::npos)
2575  {
2576  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos)
2577  {
2578  GetFunctional()->BoundaryValue_U(fdc, local_vector, scale*interval_length_);
2579  }
2580  else
2581  {
2582  GetFunctional()->BoundaryValue_U(fdc, local_vector, scale);
2583  }
2584  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos && GetFunctional()->GetType().find("timelocal") != std::string::npos)
2585  {
2586  throw DOpEException("Conflicting functional types: "+ GetFunctional()->GetType(),
2587  "OptProblemContainer::BoundaryRhs");
2588  }
2589  }
2590  }
2591  }
2592  else if (this->GetType() == "adjoint_for_ee")
2593  {
2594  //values of the derivative of the functional for error estimation
2595  if (aux_functionals_[functional_for_ee_num_]->NeedTime())
2596  {
2597  if (aux_functionals_[functional_for_ee_num_]->GetType().find(
2598  "boundary") != std::string::npos)
2599  {
2600  if(aux_functionals_[functional_for_ee_num_]->GetType().find("timedistributed") != std::string::npos)
2601  {
2602  aux_functionals_[functional_for_ee_num_]->BoundaryValue_U(fdc,
2603  local_vector, scale*interval_length_);
2604  }
2605  else
2606  {
2607  aux_functionals_[functional_for_ee_num_]->BoundaryValue_U(fdc,
2608  local_vector, scale);
2609  }
2610  if(aux_functionals_[functional_for_ee_num_]->GetType().find("timedistributed") != std::string::npos && aux_functionals_[functional_for_ee_num_]->GetType().find("timelocal") != std::string::npos)
2611  {
2612  throw DOpEException("Conflicting functional types: "+ aux_functionals_[functional_for_ee_num_]->GetType(),
2613  "OptProblemContainer::BoundaryRhs");
2614  }
2615  }
2616  }
2617  }
2618  else if (this->GetType() == "tangent")
2619  {
2620  // state values in quadrature points
2621  scale *= -1;
2622  this->GetPDE().BoundaryEquation_QT(fdc, local_vector, scale*interval_length_,
2623  scale*interval_length_);
2624  }
2625  else if (this->GetType() == "adjoint_hessian")
2626  {
2627  // state values in quadrature points
2628  if(GetFunctional()->NeedTime())
2629  {
2630  if (GetFunctional()->GetType().find("boundary") != std::string::npos)
2631  {
2632  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos)
2633  {
2634  GetFunctional()->BoundaryValue_UU(fdc, local_vector, scale*interval_length_);
2635  GetFunctional()->BoundaryValue_QU(fdc, local_vector, scale*interval_length_);
2636  }
2637  else
2638  {
2639  GetFunctional()->BoundaryValue_UU(fdc, local_vector, scale);
2640  GetFunctional()->BoundaryValue_QU(fdc, local_vector, scale);
2641  }
2642  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos && GetFunctional()->GetType().find("timelocal") != std::string::npos)
2643  {
2644  throw DOpEException("Conflicting functional types: "+ GetFunctional()->GetType(),
2645  "OptProblemContainer::BoundaryRhs");
2646  }
2647  }
2648  }
2649 
2650  scale *= -1;
2651  this->GetPDE().BoundaryEquation_UU(fdc, local_vector, scale*interval_length_,
2652  scale*interval_length_);
2653  this->GetPDE().BoundaryEquation_QU(fdc, local_vector, scale*interval_length_,
2654  scale*interval_length_);
2655  }
2656  else if (this->GetType() == "gradient")
2657  {
2658  // state values in quadrature points
2659  if(GetFunctional()->NeedTime())
2660  {
2661  if (GetFunctional()->GetType().find("boundary") != std::string::npos)
2662  {
2663  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos)
2664  {
2665  GetFunctional()->BoundaryValue_Q(fdc, local_vector, scale*interval_length_);
2666  }
2667  else
2668  {
2669  GetFunctional()->BoundaryValue_Q(fdc, local_vector, scale);
2670  }
2671  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos && GetFunctional()->GetType().find("timelocal") != std::string::npos)
2672  {
2673  throw DOpEException("Conflicting functional types: "+ GetFunctional()->GetType(),
2674  "OptProblemContainer::BoundaryRhs");
2675  }
2676  }
2677  }
2678 
2679  scale *= -1;
2680  this->GetPDE().BoundaryEquation_Q(fdc, local_vector, scale*interval_length_,
2681  scale*interval_length_);
2682  }
2683  else if (this->GetType() == "hessian")
2684  {
2685  // state values in quadrature points
2686  if(GetFunctional()->NeedTime())
2687  {
2688  if (GetFunctional()->GetType().find("boundary") != std::string::npos)
2689  {
2690  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos)
2691  {
2692  GetFunctional()->BoundaryValue_QQ(fdc, local_vector, scale*interval_length_);
2693  GetFunctional()->BoundaryValue_UQ(fdc, local_vector, scale*interval_length_);
2694  }
2695  else
2696  {
2697  GetFunctional()->BoundaryValue_QQ(fdc, local_vector, scale);
2698  GetFunctional()->BoundaryValue_UQ(fdc, local_vector, scale);
2699  }
2700  if(GetFunctional()->GetType().find("timedistributed") != std::string::npos && GetFunctional()->GetType().find("timelocal") != std::string::npos)
2701  {
2702  throw DOpEException("Conflicting functional types: "+ GetFunctional()->GetType(),
2703  "OptProblemContainer::BoundaryRhs");
2704  }
2705  }
2706  }
2707 
2708  scale *= -1;
2709  this->GetPDE().BoundaryEquation_QTT(fdc, local_vector, scale*interval_length_,
2710  scale*interval_length_);
2711  this->GetPDE().BoundaryEquation_UQ(fdc, local_vector, scale*interval_length_,
2712  scale*interval_length_);
2713  this->GetPDE().BoundaryEquation_QQ(fdc, local_vector, scale*interval_length_,
2714  scale*interval_length_);
2715  }
2716  else
2717  {
2718  throw DOpEException("Not implemented",
2719  "OptProblemContainer::BoundaryRhs");
2720  }
2721 
2722  }
2723 
2724  /******************************************************/
2725 
2726  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
2727  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
2728  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
2729  template<int, int> class DH>
2730  template<typename DATACONTAINER>
2731  void
2732  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD,
2733  CONSTRAINTS, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::ElementMatrix(
2734  const DATACONTAINER& edc,
2735  dealii::FullMatrix<double> &local_entry_matrix, double scale,
2736  double scale_ico)
2737  {
2738 
2739  if (this->GetType() == "state" || this->GetType() == "tangent")
2740  {
2741  // state values in quadrature points
2742  this->GetPDE().ElementMatrix(edc, local_entry_matrix, scale*interval_length_, scale_ico*interval_length_);
2743  }
2744  else if (this->GetType() == "adjoint"
2745  || this->GetType() == "adjoint_for_ee"
2746  || this->GetType() == "adjoint_hessian")
2747  {
2748  // state values in quadrature points
2749  this->GetPDE().ElementMatrix_T(edc, local_entry_matrix, scale*interval_length_,
2750  scale_ico*interval_length_);
2751  }
2752  else if ((this->GetType() == "gradient")
2753  || (this->GetType() == "hessian"))
2754  {
2755  // control values in quadrature points
2756  this->GetPDE().ControlElementMatrix(edc, local_entry_matrix, scale*c_interval_length_);
2757  }
2758  else
2759  {
2760  throw DOpEException("Not implemented",
2761  "OptProblemContainer::ElementMatrix");
2762  }
2763 
2764  }
2765 
2766  /******************************************************/
2767 
2768  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
2769  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
2770  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
2771  template<int, int> class DH>
2772  template<typename DATACONTAINER>
2773  void
2774  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD,
2775  CONSTRAINTS, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::ElementTimeMatrix(
2776  const DATACONTAINER& edc, FullMatrix<double> &local_entry_matrix)
2777  {
2778 
2779  if (this->GetType() == "state" || this->GetType() == "tangent")
2780  {
2781  // state values in quadrature points
2782  this->GetPDE().ElementTimeMatrix(edc, local_entry_matrix);
2783  }
2784  else if (this->GetType() == "adjoint"
2785  || this->GetType() == "adjoint_for_ee"
2786  || this->GetType() == "adjoint_hessian")
2787  {
2788  this->GetPDE().ElementTimeMatrix_T(edc, local_entry_matrix);
2789  }
2790  else if ((this->GetType() == "gradient")
2791  || (this->GetType() == "hessian"))
2792  {
2793  throw DOpEException("Not implemented",
2794  "OptProblemContainer::ElementTimeMatrix");
2795  }
2796  else
2797  {
2798  throw DOpEException("Not implemented",
2799  "OptProblemContainer::ElementTimeMatrix");
2800  }
2801 
2802  }
2803 
2804  /******************************************************/
2805 
2806  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
2807  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
2808  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
2809  template<int, int> class DH>
2810  template<typename DATACONTAINER>
2811  void
2812  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD,
2813  CONSTRAINTS, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::ElementTimeMatrixExplicit(
2814  const DATACONTAINER& edc,
2815  dealii::FullMatrix<double> &local_entry_matrix)
2816  {
2817 
2818  if (this->GetType() == "state" || this->GetType() == "tangent")
2819  {
2820  this->GetPDE().ElementTimeMatrixExplicit(edc, local_entry_matrix);
2821  }
2822  else if (this->GetType() == "adjoint"
2823  || this->GetType() == "adjoint_hessian")
2824  {
2825  this->GetPDE().ElementTimeMatrixExplicit_T(edc, local_entry_matrix);
2826  }
2827  else if ((this->GetType() == "gradient")
2828  || (this->GetType() == "hessian"))
2829  {
2830  throw DOpEException("Not implemented",
2831  "OptProblemContainer::ElementTimeMatrixExplicit");
2832  }
2833  else
2834  {
2835  throw DOpEException("Not implemented",
2836  "OptProblemContainer::ElementTimeMatrixExplicit");
2837  }
2838 
2839  }
2840 
2841  /******************************************************/
2842 
2843  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
2844  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
2845  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
2846  template<int, int> class DH>
2847  template<typename FACEDATACONTAINER>
2848  void
2849  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD,
2850  CONSTRAINTS, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::FaceMatrix(
2851  const FACEDATACONTAINER& fdc, FullMatrix<double> &local_entry_matrix,
2852  double scale, double scale_ico)
2853  {
2854  if (this->GetType() == "state" || this->GetType() == "tangent")
2855  {
2856  // state values in face quadrature points
2857  this->GetPDE().FaceMatrix(fdc, local_entry_matrix, scale*interval_length_, scale_ico*interval_length_);
2858  }
2859  else if (this->GetType() == "adjoint"
2860  || this->GetType() == "adjoint_for_ee"
2861  || this->GetType() == "adjoint_hessian")
2862  {
2863  // state values in quadrature points
2864  this->GetPDE().FaceMatrix_T(fdc, local_entry_matrix, scale*interval_length_,
2865  scale_ico*interval_length_);
2866  }
2867 // else if ((this->GetType() == "gradient") || (this->GetType() == "hessian"))
2868 // {
2869 // // control values in quadrature points
2870 // this->GetPDE().ControlFaceMatrix(fdc, local_entry_matrix);
2871 // }
2872  else
2873  {
2874  throw DOpEException("Not implemented",
2875  "OptProblemContainer::NewtonFaceMatrix");
2876  }
2877 
2878  }
2879 
2880  /******************************************************/
2881 
2882  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
2883  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
2884  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
2885  template<int, int> class DH>
2886  template<typename FACEDATACONTAINER>
2887  void
2888  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD,
2889  CONSTRAINTS, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::InterfaceMatrix(
2890  const FACEDATACONTAINER& fdc, FullMatrix<double> &local_entry_matrix,
2891  double scale, double scale_ico)
2892  {
2893  if (this->GetType() == "state")
2894  {
2895  this->GetPDE().InterfaceMatrix(fdc, local_entry_matrix, scale*interval_length_,
2896  scale_ico*interval_length_);
2897  }
2898  else if (this->GetType() == "adjoint"
2899  || this->GetType() == "adjoint_for_ee"
2900  || this->GetType() == "adjoint_hessian")
2901  {
2902  this->GetPDE().InterfaceMatrix_T(fdc, local_entry_matrix, scale*interval_length_,
2903  scale_ico*interval_length_);
2904  }
2905  else
2906  {
2907  throw DOpEException("Not implemented",
2908  "OptProblemContainer::NewtonInterfaceMatrix");
2909  }
2910  }
2911 
2912  /******************************************************/
2913 
2914  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
2915  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
2916  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
2917  template<int, int> class DH>
2918  template<typename FACEDATACONTAINER>
2919  void
2920  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD,
2921  CONSTRAINTS, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::BoundaryMatrix(
2922  const FACEDATACONTAINER& fdc, FullMatrix<double> &local_matrix,
2923  double scale, double scale_ico)
2924  {
2925  if (this->GetType() == "state" || this->GetType() == "tangent")
2926  {
2927  // state values in face quadrature points
2928  this->GetPDE().BoundaryMatrix(fdc, local_matrix, scale*interval_length_,
2929  scale_ico*interval_length_);
2930  }
2931  else if (this->GetType() == "adjoint"
2932  || this->GetType() == "adjoint_for_ee"
2933  || this->GetType() == "adjoint_hessian")
2934  {
2935  // state values in quadrature points
2936  this->GetPDE().BoundaryMatrix_T(fdc, local_matrix, scale*interval_length_,
2937  scale_ico*interval_length_);
2938  }
2939 // else if ((this->GetType() == "gradient") || (this->GetType() == "hessian"))
2940 // {
2941 // // control values in quadrature points
2942 // this->GetPDE().ControlBoundaryMatrix(fdc, local_matrix);
2943 // }
2944  else
2945  {
2946  throw DOpEException("Not implemented",
2947  "OptProblemContainer::ElementBoundaryMatrix");
2948  }
2949 
2950  }
2951 
2952  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
2953  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
2954  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
2955  template<int, int> class DH>
2956  void
2957  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
2958  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::ComputeLocalControlConstraints(
2959  VECTOR& constraints,
2960  const std::map<std::string, const dealii::Vector<double>*> &/*values*/,
2961  const std::map<std::string, const VECTOR*> &block_values)
2962  {
2963  if (this->GetType() == "constraints")
2964  {
2965  if (this->GetSpaceTimeHandler()->GetNLocalConstraints() != 0)
2966  {
2967  const VECTOR& control = *GetBlockVector(block_values, "control");
2968  this->GetConstraints()->EvaluateLocalControlConstraints(control,
2969  constraints);
2970  }
2971  }
2972  else
2973  {
2974  throw DOpEException("Wrong problem type" + this->GetType(),
2975  "OptProblemContainer::ComputeLocalConstraints");
2976  }
2977  }
2978 
2979  /******************************************************/
2980 
2981  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
2982  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
2983  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
2984  template<int, int> class DH>
2985  std::string
2986  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
2987  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::GetDoFType() const
2988  {
2989  if (this->GetType() == "state" || this->GetType() == "adjoint"
2990  || this->GetType() == "adjoint_for_ee" || this->GetType() == "tangent"
2991  || this->GetType() == "adjoint_hessian")
2992  {
2993  return "state";
2994  }
2995  else if ((this->GetType() == "gradient") || (this->GetType() == "hessian")
2996  || (this->GetType() == "hessian_inverse"))
2997  {
2998  return "control";
2999  }
3000  else
3001  {
3002  throw DOpEException("Unknown Type:" + this->GetType(),
3003  "OptProblemContainer::GetDoFType");
3004  }
3005  }
3006 
3007  /******************************************************/
3008 
3009  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3010  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3011  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3012  template<int, int> class DH>
3013  const FE<dealdim, dealdim>&
3014  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3015  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::GetFESystem() const
3016  {
3017  if ((this->GetType() == "state") || (this->GetType() == "adjoint")
3018  || this->GetType() == "adjoint_for_ee" || this->GetType() == "tangent"
3019  || this->GetType() == "adjoint_hessian")
3020  {
3021  return this->GetSpaceTimeHandler()->GetFESystem("state");
3022  }
3023  else if ((this->GetType() == "gradient") || (this->GetType() == "hessian")
3024  || (this->GetType() == "global_constraint_gradient"))
3025  {
3026 #if dope_dimension > 0
3027  if(dopedim == dealdim)
3028  return this->GetSpaceTimeHandler()->GetFESystem("control");
3029  else
3030  throw DOpEException("Non matching dimensions!","OptProblemContainer::GetFESystem");
3031 #else
3032  return this->GetSpaceTimeHandler()->GetFESystem("state");
3033 #endif
3034  }
3035  else
3036  {
3037  throw DOpEException("Unknown Type:" + this->GetType(),
3038  "OptProblemContainer::GetFESystem");
3039  }
3040  }
3041 
3042  /******************************************************/
3043 
3044  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3045  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3046  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3047  template<int, int> class DH>
3048  UpdateFlags
3049  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3050  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::GetUpdateFlags() const
3051  {
3052 
3053  UpdateFlags r;
3054  if (this->GetType().find("aux_functional") != std::string::npos)
3055  {
3056  r = aux_functionals_[this->GetTypeNum()]->GetUpdateFlags();
3057  }
3058  else if (this->GetType().find("functional") != std::string::npos)
3059  {
3060  r = this->GetFunctional()->GetUpdateFlags();
3061  }
3062  else if (this->GetType().find("constraints") != std::string::npos)
3063  {
3064  r = this->GetConstraints()->GetUpdateFlags();
3065  }
3066  else if (this->GetType() == "functional_for_ee")
3067  {
3068  r = aux_functionals_[functional_for_ee_num_]->GetUpdateFlags();
3069  }
3070  else
3071  {
3072  r = this->GetPDE().GetUpdateFlags();
3073  if (this->GetType() == "adjoint_hessian" || this->GetType() == "adjoint"
3074  || (this->GetType() == "hessian")
3075  || (this->GetType() == "gradient"))
3076  {
3077  r = r | this->GetFunctional()->GetUpdateFlags();
3078  }
3079  else if (this->GetType() == "adjoint_for_ee")
3080  {
3081  r = r | aux_functionals_[functional_for_ee_num_]->GetUpdateFlags();
3082  }
3083  }
3084  return r | update_JxW_values;
3085  }
3086 
3087  /******************************************************/
3088 
3089  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3090  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3091  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3092  template<int, int> class DH>
3093  UpdateFlags
3094  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3095  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::GetFaceUpdateFlags() const
3096  {
3097  UpdateFlags r;
3098  if (this->GetType().find("aux_functional") != std::string::npos)
3099  {
3100  r = aux_functionals_[this->GetTypeNum()]->GetFaceUpdateFlags();
3101  }
3102  else if (this->GetType().find("functional") != std::string::npos)
3103  {
3104  r = this->GetFunctional()->GetFaceUpdateFlags();
3105  }
3106  else if (this->GetType().find("constraints") != std::string::npos)
3107  {
3108  r = this->GetConstraints()->GetFaceUpdateFlags();
3109  }
3110  else if (this->GetType() == "functional_for_ee")
3111  {
3112  r = aux_functionals_[functional_for_ee_num_]->GetUpdateFlags();
3113  }
3114  else
3115  {
3116  r = this->GetPDE().GetFaceUpdateFlags();
3117  if (this->GetType() == "adjoint_hessian" || this->GetType() == "adjoint"
3118  || (this->GetType() == "hessian"))
3119  {
3120  r = r | this->GetFunctional()->GetFaceUpdateFlags();
3121  }
3122  else if (this->GetType() == "adjoint_for_ee")
3123  {
3124  r = r
3125  | aux_functionals_[functional_for_ee_num_]->GetFaceUpdateFlags();
3126  }
3127  }
3128  return r | update_JxW_values;
3129  }
3130 
3131  /******************************************************/
3132 
3133  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3134  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3135  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3136  template<int, int> class DH>
3137  std::string
3138  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3139  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::GetFunctionalType() const
3140  {
3141  if (this->GetType() == "aux_functional")
3142  {
3143  return aux_functionals_[this->GetTypeNum()]->GetType();
3144  }
3145  else if (this->GetType() == "functional_for_ee")
3146  {
3147  return aux_functionals_[functional_for_ee_num_]->GetType();
3148  }
3149  return GetFunctional()->GetType();
3150  }
3151 
3152  /******************************************************/
3153 
3154  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3155  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3156  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3157  template<int, int> class DH>
3158  std::string
3159  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3160  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::GetFunctionalName() const
3161  {
3162  if (this->GetType() == "aux_functional")
3163  {
3164  return aux_functionals_[this->GetTypeNum()]->GetName();
3165  }
3166  else if (this->GetType() == "functional_for_ee")
3167  {
3168  return aux_functionals_[functional_for_ee_num_]->GetName();
3169  }
3170  return GetFunctional()->GetName();
3171  }
3172 
3173  /******************************************************/
3174 
3175  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3176  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3177  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3178  template<int, int> class DH>
3179  std::string
3180  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3181  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::GetConstraintType() const
3182  {
3183  return GetConstraints()->GetType();
3184  }
3185 
3186  /******************************************************/
3187 
3188  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3189  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3190  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3191  template<int, int> class DH>
3192  void
3193  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3194  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::SetTime(double time,
3195  const TimeIterator& interval, bool initial)
3196  {
3197  GetSpaceTimeHandler()->SetInterval(interval);
3198  initial_ = initial;
3199  // GetSpaceTimeHandler()->SetTimeDoFNumber(time_point);
3200  interval_length_ = GetSpaceTimeHandler()->GetStepSize();
3201 
3202  if(GetSpaceTimeHandler()->GetControlType() == DOpEtypes::ControlType::initial || GetSpaceTimeHandler()->GetControlType() == DOpEtypes::ControlType::stationary )
3203  {
3204  c_interval_length_ = 1.;
3205  }
3206  else
3207  {
3208  c_interval_length_ = GetSpaceTimeHandler()->GetStepSize();
3209  }
3210 
3211  { //Zeit an Dirichlet Werte uebermitteln
3212  for (unsigned int i = 0;
3213  i < transposed_control_gradient_dirichlet_values_.size(); i++)
3214  transposed_control_gradient_dirichlet_values_[i]->SetTime(time);
3215  for (unsigned int i = 0;
3216  i < transposed_control_hessian_dirichlet_values_.size(); i++)
3217  transposed_control_hessian_dirichlet_values_[i]->SetTime(time);
3218  for (unsigned int i = 0; i < primal_dirichlet_values_.size(); i++)
3219  primal_dirichlet_values_[i]->SetTime(time);
3220  for (unsigned int i = 0; i < tangent_dirichlet_values_.size(); i++)
3221  tangent_dirichlet_values_[i]->SetTime(time);
3222  for (unsigned int i = 0; i < control_dirichlet_values_.size(); i++)
3223  control_dirichlet_values_[i]->SetTime(time);
3224  //Functionals
3225  GetFunctional()->SetTime(time);
3226  for (unsigned int i = 0; i < aux_functionals_.size(); i++)
3227  aux_functionals_[i]->SetTime(time);
3228  //PDE
3229  this->GetPDE().SetTime(time);
3230  }
3231  //Update Auxiliary Control, State and Constraint Vectors
3232  {
3233  typename std::map<std::string, const ControlVector<VECTOR> *>::iterator it =
3234  auxiliary_controls_.begin();
3235  for (; it != auxiliary_controls_.end(); it++)
3236  {
3237  it->second->SetTime(time, interval);
3238  }
3239  }
3240  {
3241  typename std::map<std::string, const StateVector<VECTOR> *>::iterator it =
3242  auxiliary_state_.begin();
3243  for (; it != auxiliary_state_.end(); it++)
3244  {
3245  it->second->SetTime(time, interval);
3246  }
3247  }
3248  {
3249  typename std::map<std::string, const ConstraintVector<VECTOR> *>::iterator it =
3250  auxiliary_constraints_.begin();
3251  for (; it != auxiliary_constraints_.end(); it++)
3252  {
3253  it->second->SetTime(time, interval);
3254  }
3255  }
3256 
3257  }
3258 
3259  /******************************************************/
3260 
3261  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3262  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3263  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3264  template<int, int> class DH>
3265  void
3266  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3267  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::ComputeSparsityPattern(
3268  SPARSITYPATTERN & sparsity) const
3269  {
3270  if (this->GetType() == "state" || this->GetType() == "tangent"
3271  || this->GetType() == "adjoint_for_ee" || this->GetType() == "adjoint"
3272  || this->GetType() == "adjoint_hessian")
3273  {
3274  this->GetSpaceTimeHandler()->ComputeStateSparsityPattern(sparsity);
3275  }
3276  else if ((this->GetType() == "gradient")
3277  || (this->GetType() == "hessian"))
3278  {
3279 #if dope_dimension > 0
3280  this->GetSpaceTimeHandler()->ComputeControlSparsityPattern(sparsity);
3281 #else
3282  throw DOpEException("Wrong dimension",
3283  "OptProblemContainer::ComputeSparsityPattern");
3284 #endif
3285  }
3286  else
3287  {
3288  throw DOpEException("Unknown type " + this->GetType(),
3289  "OptProblemContainer::ComputeSparsityPattern");
3290  }
3291  }
3292 
3293  /******************************************************/
3294 
3295  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3296  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3297  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3298  template<int, int> class DH>
3299  void
3300  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3301  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::PostProcessConstraints(
3302  ConstraintVector<VECTOR>& g) const
3303  {
3304  return this->GetConstraints()->PostProcessConstraints(g);
3305  }
3306  /******************************************************/
3307 
3308  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3309  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3310  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3311  template<int, int> class DH>
3312  void
3313  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3314  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::AddAuxiliaryControl(
3315  const ControlVector<VECTOR>* c, std::string name)
3316  {
3317  if (auxiliary_controls_.find(name) != auxiliary_controls_.end())
3318  {
3319  throw DOpEException(
3320  "Adding multiple Data with name " + name + " is prohibited!",
3321  "OptProblemContainer::AddAuxiliaryControl");
3322  }
3323  auxiliary_controls_.insert(
3324  std::pair<std::string, const ControlVector<VECTOR>*>(name, c));
3325  }
3326 
3327  /******************************************************/
3328 
3329  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3330  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3331  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3332  template<int, int> class DH>
3333  void
3334  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3335  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::AddAuxiliaryState(
3336  const StateVector<VECTOR>* c, std::string name)
3337  {
3338  if (auxiliary_state_.find(name) != auxiliary_state_.end())
3339  {
3340  throw DOpEException(
3341  "Adding multiple Data with name " + name + " is prohibited!",
3342  "OptProblemContainer::AddAuxiliaryState");
3343  }
3344  auxiliary_state_.insert(
3345  std::pair<std::string, const StateVector<VECTOR>*>(name, c));
3346  }
3347  /******************************************************/
3348 
3349  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3350  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3351  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3352  template<int, int> class DH>
3353  void
3354  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3355  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::AddAuxiliaryConstraint(
3356  const ConstraintVector<VECTOR>* c, std::string name)
3357  {
3358  if (auxiliary_constraints_.find(name) != auxiliary_constraints_.end())
3359  {
3360  throw DOpEException(
3361  "Adding multiple Data with name " + name + " is prohibited!",
3362  "OptProblemContainer::AddAuxiliaryConstraint");
3363  }
3364  auxiliary_constraints_.insert(
3365  std::pair<std::string, const ConstraintVector<VECTOR>*>(name, c));
3366  }
3367  /******************************************************/
3368 
3369  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3370  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3371  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3372  template<int, int> class DH>
3373  const ControlVector<VECTOR>*
3374  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3375  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::GetAuxiliaryControl(
3376  std::string name) const
3377  {
3378  typename std::map<std::string, const ControlVector<VECTOR> *>::const_iterator it =
3379  auxiliary_controls_.find(name);
3380  if (it == auxiliary_controls_.end())
3381  {
3382  throw DOpEException("Could not find Data with name " + name,
3383  "OptProblemContainer::GetAuxiliaryControl");
3384  }
3385  return it->second;
3386  }
3387 
3388  /******************************************************/
3389 
3390  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3391  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3392  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3393  template<int, int> class DH>
3394  const StateVector<VECTOR>*
3395  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3396  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::GetAuxiliaryState(
3397  std::string name) const
3398  {
3399  typename std::map<std::string, const StateVector<VECTOR> *>::const_iterator it =
3400  auxiliary_state_.find(name);
3401  if (it == auxiliary_state_.end())
3402  {
3403  throw DOpEException("Could not find Data with name " + name,
3404  "OptProblemContainer::GetAuxiliaryState");
3405  }
3406  return it->second;
3407  }
3408 
3409  /******************************************************/
3410 
3411  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3412  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3413  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3414  template<int, int> class DH>
3415  void
3416  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3417  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::DeleteAuxiliaryControl(
3418  std::string name)
3419  {
3420  typename std::map<std::string, const ControlVector<VECTOR> *>::iterator it =
3421  auxiliary_controls_.find(name);
3422  if (it == auxiliary_controls_.end())
3423  {
3424  throw DOpEException(
3425  "Deleting Data " + name + " is impossible! Data not found",
3426  "OptProblemContainer::DeleteAuxiliaryControl");
3427  }
3428  auxiliary_controls_.erase(it);
3429  }
3430 
3431  /******************************************************/
3432 
3433  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3434  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3435  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3436  template<int, int> class DH>
3437  void
3438  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3439  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::DeleteAuxiliaryState(
3440  std::string name)
3441  {
3442  typename std::map<std::string, const StateVector<VECTOR> *>::iterator it =
3443  auxiliary_state_.find(name);
3444  if (it == auxiliary_state_.end())
3445  {
3446  throw DOpEException(
3447  "Deleting Data " + name + " is impossible! Data not found",
3448  "OptProblemContainer::DeleteAuxiliaryState");
3449  }
3450  auxiliary_state_.erase(it);
3451  }
3452 
3453  /******************************************************/
3454 
3455  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3456  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3457  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3458  template<int, int> class DH>
3459  void
3460  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3461  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::DeleteAuxiliaryConstraint(
3462  std::string name)
3463  {
3464  typename std::map<std::string, const ConstraintVector<VECTOR> *>::iterator it =
3465  auxiliary_constraints_.find(name);
3466  if (it == auxiliary_constraints_.end())
3467  {
3468  throw DOpEException(
3469  "Deleting Data " + name + " is impossible! Data not found",
3470  "OptProblemContainer::DeleteAuxiliaryConstraint");
3471  }
3472  auxiliary_constraints_.erase(it);
3473  }
3474 
3475  /******************************************************/
3476 
3477  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3478  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3479  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3480  template<int, int> class DH>
3481  FUNCTIONAL*
3482  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3483  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::GetFunctional()
3484  {
3485  if (this->GetType() == "aux_functional"
3486  || this->GetType() == "functional_for_ee")
3487  {
3488  //This may no longer happen!
3489  abort();
3490  // return aux_functionals_[this->GetTypeNum()];
3491  }
3492  return functional_;
3493 
3494  }
3495 
3496  /******************************************************/
3497 
3498  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3499  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3500  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3501  template<int, int> class DH>
3502  const FUNCTIONAL*
3503  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3504  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::GetFunctional() const
3505  {
3506  if (this->GetType() == "aux_functional"
3507  || this->GetType() == "functional_for_ee")
3508  {
3509  //This may no longer happen!
3510  abort();
3511  // return aux_functionals_[this->GetTypeNum()];
3512  }
3513  return functional_;
3514  }
3515 
3516  /******************************************************/
3517 
3518  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3519  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3520  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3521  template<int, int> class DH>
3522  bool
3523  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3524  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::HasFaces() const
3525  {
3526  if (this->GetType().find("aux_functional") != std::string::npos)
3527  {
3528  return aux_functionals_[this->GetTypeNum()]->HasFaces();
3529  }
3530  else if (this->GetType().find("functional") != std::string::npos)
3531  {
3532  return this->GetFunctional()->HasFaces();
3533  }
3534  else if (this->GetType().find("constraint") != std::string::npos)
3535  {
3536  return this->GetConstraints()->HasFaces();
3537  }
3538  else
3539  {
3540  if ((this->GetType() == "state") || (this->GetType() == "tangent")
3541  || (this->GetType() == "gradient"))
3542  {
3543  return this->GetPDE().HasFaces();
3544  }
3545  else if ((this->GetType() == "adjoint")
3546  || (this->GetType() == "adjoint_hessian")
3547  || (this->GetType() == "hessian"))
3548  {
3549  return this->GetPDE().HasFaces() || this->GetFunctional()->HasFaces();
3550  }
3551  else if (this->GetType() == "adjoint_for_ee")
3552  {
3553  return this->GetPDE().HasFaces()
3554  || aux_functionals_[functional_for_ee_num_]->HasFaces();
3555  }
3556  else
3557  {
3558  throw DOpEException("Unknown Type: '" + this->GetType() + "'!",
3559  "OptProblemContainer::HasFaces");
3560  }
3561  }
3562  }
3563 
3564  /******************************************************/
3565 
3566  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3567  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3568  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3569  template<int, int> class DH>
3570  bool
3571  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3572  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::HasPoints() const
3573  {
3574  if (this->GetType().find("constraint") != std::string::npos
3575  || (this->GetType() == "functional")
3576  || this->GetType() == "aux_functional" || (this->GetType() == "state")
3577  || (this->GetType() == "tangent"))
3578  {
3579  // We dont need PointRhs in this cases.
3580  return false;
3581  }
3582  else if ((this->GetType() == "adjoint")
3583  || (this->GetType() == "adjoint_hessian")
3584  || (this->GetType() == "hessian"))
3585  {
3586  return this->GetFunctional()->HasPoints();
3587  }
3588  else if (this->GetType() == "adjoint_for_ee")
3589  {
3590  return aux_functionals_[functional_for_ee_num_]->HasPoints();
3591  }
3592  else if (this->GetType() == "gradient")
3593  {
3594  return this->GetFunctional()->HasPoints();
3595  }
3596  else
3597  {
3598  throw DOpEException("Unknown Type: '" + this->GetType() + "'!",
3599  "OptProblem::HasPoints");
3600  }
3601  }
3602 // }
3603 
3604  /******************************************************/
3605 
3606  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3607  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3608  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3609  template<int, int> class DH>
3610  bool
3611  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3612  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::HasInterfaces() const
3613  {
3614  if (this->GetType().find("aux_functional") != std::string::npos)
3615  {
3616  return false;
3617  }
3618  else if (this->GetType().find("functional") != std::string::npos)
3619  {
3620  return false;
3621  }
3622  else if (this->GetType().find("constraint") != std::string::npos)
3623  {
3624  return false;
3625  }
3626  else
3627  {
3628  if ((this->GetType() == "state") || this->GetType() == "adjoint_for_ee"
3629  || (this->GetType() == "tangent") || (this->GetType() == "gradient")
3630  || (this->GetType() == "adjoint")
3631  || (this->GetType() == "adjoint_hessian")
3632  || (this->GetType() == "hessian"))
3633  {
3634  return this->GetPDE().HasInterfaces();
3635  }
3636  else
3637  {
3638  throw DOpEException("Unknown Type: '" + this->GetType() + "'!",
3639  "OptProblemContainer::HasFaces");
3640  }
3641  }
3642  }
3643 
3644  /******************************************************/
3645 
3646  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3647  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3648  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3649  template<int, int> class DH>
3650  void
3651  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3652  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::SetControlDirichletBoundaryColors(
3653  unsigned int color, const std::vector<bool>& comp_mask,
3654  const DOpEWrapper::Function<dealdim>* values)
3655  {
3656  assert(values);
3657 
3658  unsigned int comp = control_dirichlet_colors_.size();
3659  for (unsigned int i = 0; i < control_dirichlet_colors_.size(); ++i)
3660  {
3661  if (control_dirichlet_colors_[i] == color)
3662  {
3663  comp = i;
3664  break;
3665  }
3666  }
3667  if (comp != control_dirichlet_colors_.size())
3668  {
3669  std::stringstream s;
3670  s << "ControlDirichletColor" << color << " has multiple occurences !";
3671  throw DOpEException(s.str(),
3672  "OptProblemContainer::SetControlDirichletBoundary");
3673  }
3674  control_dirichlet_colors_.push_back(color);
3675  control_dirichlet_comps_.push_back(comp_mask);
3676  control_dirichlet_values_.push_back(values);
3677  }
3678 
3679  /******************************************************/
3680 
3681  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3682  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3683  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3684  template<int, int> class DH>
3685  void
3686  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3687  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::SetDirichletBoundaryColors(
3688  unsigned int color, const std::vector<bool>& comp_mask,
3689  const DD* values)
3690  {
3691  assert(values);
3692  assert(values->n_components() == comp_mask.size());
3693 
3694  unsigned int comp = dirichlet_colors_.size();
3695  for (unsigned int i = 0; i < dirichlet_colors_.size(); ++i)
3696  {
3697  if (dirichlet_colors_[i] == color)
3698  {
3699  comp = i;
3700  break;
3701  }
3702  }
3703  if (comp != dirichlet_colors_.size())
3704  {
3705  std::stringstream s;
3706  s << "DirichletColor" << color << " has multiple occurrences !";
3707  throw DOpEException(s.str(),
3708  "OptProblemContainer::SetDirichletBoundary");
3709  }
3710  dirichlet_colors_.push_back(color);
3711  dirichlet_comps_.push_back(comp_mask);
3712  PrimalDirichletData<DD, VECTOR, dealdim>* data =
3713  new PrimalDirichletData<DD, VECTOR, dealdim>(*values);
3714  primal_dirichlet_values_.push_back(data);
3715  TangentDirichletData<DD, VECTOR, dealdim>* tdata =
3716  new TangentDirichletData<DD, VECTOR, dealdim>(*values);
3717  tangent_dirichlet_values_.push_back(tdata);
3718 
3719  if (values->NeedsControl())
3720  {
3721  control_transposed_dirichlet_colors_.push_back(color);
3722  TransposedGradientDirichletData<DD, VECTOR, dealdim> * gdata =
3723  new TransposedGradientDirichletData<DD, VECTOR, dealdim>(
3724  *values);
3725  transposed_control_gradient_dirichlet_values_.push_back(gdata);
3726  TransposedHessianDirichletData<DD, VECTOR, dealdim> * hdata =
3727  new TransposedHessianDirichletData<DD, VECTOR, dealdim>(
3728  *values);
3729  transposed_control_hessian_dirichlet_values_.push_back(hdata);
3730  }
3731  }
3732 
3733  /******************************************************/
3734 
3735  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3736  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3737  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3738  template<int, int> class DH>
3739  const std::vector<unsigned int>&
3740  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3741  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::GetDirichletColors() const
3742  {
3743  if ((this->GetType() == "state") || (this->GetType() == "adjoint")
3744  || this->GetType() == "adjoint_for_ee"
3745  || (this->GetType() == "tangent")
3746  || (this->GetType() == "adjoint_hessian"))
3747  {
3748  return dirichlet_colors_;
3749  }
3750  else if ((this->GetType() == "gradient") || (this->GetType() == "hessian")
3751  || (this->GetType() == "global_constraint_gradient"))
3752  {
3753  return control_dirichlet_colors_;
3754  }
3755  else
3756  {
3757  throw DOpEException("Unknown Type:" + this->GetType(),
3758  "OptProblemContainer::GetDirichletColors");
3759  }
3760  }
3761  /******************************************************/
3762 
3763  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3764  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3765  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3766  template<int, int> class DH>
3767  const std::vector<unsigned int>&
3768  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3769  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::GetTransposedDirichletColors() const
3770  {
3771  if ((this->GetType() == "gradient") || (this->GetType() == "hessian"))
3772  {
3773  return control_transposed_dirichlet_colors_;
3774  }
3775  else
3776  {
3777  throw DOpEException("Unknown Type:" + this->GetType(),
3778  "OptProblemContainer::GetTransposedDirichletColors");
3779  }
3780  }
3781 
3782  /******************************************************/
3783 
3784  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3785  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3786  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3787  template<int, int> class DH>
3788  const std::vector<bool>&
3789  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3790  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::GetDirichletCompMask(
3791  unsigned int color) const
3792  {
3793  if ((this->GetType() == "state") || (this->GetType() == "adjoint")
3794  || this->GetType() == "adjoint_for_ee"
3795  || (this->GetType() == "tangent")
3796  || (this->GetType() == "adjoint_hessian"))
3797  {
3798  unsigned int comp = dirichlet_colors_.size();
3799  for (unsigned int i = 0; i < dirichlet_colors_.size(); ++i)
3800  {
3801  if (dirichlet_colors_[i] == color)
3802  {
3803  comp = i;
3804  break;
3805  }
3806  }
3807  if (comp == dirichlet_colors_.size())
3808  {
3809  std::stringstream s;
3810  s << "DirichletColor" << color << " has not been found !";
3811  throw DOpEException(s.str(),
3812  "OptProblemContainer::GetDirichletCompMask");
3813  }
3814  return dirichlet_comps_[comp];
3815  }
3816  else if ((this->GetType() == "gradient")
3817  || (this->GetType() == "hessian"))
3818  {
3819  unsigned int comp = control_dirichlet_colors_.size();
3820  for (unsigned int i = 0; i < control_dirichlet_colors_.size(); ++i)
3821  {
3822  if (control_dirichlet_colors_[i] == color)
3823  {
3824  comp = i;
3825  break;
3826  }
3827  }
3828  if (comp == control_dirichlet_colors_.size())
3829  {
3830  std::stringstream s;
3831  s << "ControlDirichletColor" << color << " has not been found !";
3832  throw DOpEException(s.str(),
3833  "OptProblemContainer::GetDirichletCompMask");
3834  }
3835  return control_dirichlet_comps_[comp];
3836  }
3837  else
3838  {
3839  throw DOpEException("Unknown Type:" + this->GetType(),
3840  "OptProblemContainer::GetDirichletCompMask");
3841  }
3842  }
3843  /******************************************************/
3844 
3845  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3846  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3847  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3848  template<int, int> class DH>
3849  const std::vector<bool>&
3850  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3851  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::GetTransposedDirichletCompMask(
3852  unsigned int color) const
3853  {
3854  if ((this->GetType() == "gradient") || (this->GetType() == "hessian"))
3855  {
3856  unsigned int comp = dirichlet_colors_.size();
3857  for (unsigned int i = 0; i < dirichlet_colors_.size(); ++i)
3858  {
3859  if (dirichlet_colors_[i] == color)
3860  {
3861  comp = i;
3862  break;
3863  }
3864  }
3865  if (comp == dirichlet_colors_.size())
3866  {
3867  std::stringstream s;
3868  s << "DirichletColor" << color << " has not been found !";
3869  throw DOpEException(s.str(),
3870  "OptProblemContainer::GetDirichletCompMask");
3871  }
3872  return dirichlet_comps_[comp];
3873  }
3874  else
3875  {
3876  throw DOpEException("Unknown Type:" + this->GetType(),
3877  "OptProblemContainer::GetTransposedDirichletCompMask");
3878  }
3879  }
3880 
3881  /******************************************************/
3882 
3883  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3884  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3885  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3886  template<int, int> class DH>
3887  const Function<dealdim>&
3888  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3889  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::GetDirichletValues(
3890  unsigned int color,
3891  const std::map<std::string, const dealii::Vector<double>*> &param_values,
3892  const std::map<std::string, const VECTOR*> &domain_values) const
3893  {
3894 
3895  unsigned int col = dirichlet_colors_.size();
3896  if ((this->GetType() == "state") || (this->GetType() == "adjoint")
3897  || this->GetType() == "adjoint_for_ee"
3898  || (this->GetType() == "tangent")
3899  || (this->GetType() == "adjoint_hessian"))
3900  {
3901  for (unsigned int i = 0; i < dirichlet_colors_.size(); ++i)
3902  {
3903  if (dirichlet_colors_[i] == color)
3904  {
3905  col = i;
3906  break;
3907  }
3908  }
3909  if (col == dirichlet_colors_.size())
3910  {
3911  std::stringstream s;
3912  s << "DirichletColor" << color << " has not been found !";
3913  throw DOpEException(s.str(),
3914  "OptProblemContainer::GetDirichletValues");
3915  }
3916  }
3917  else if (this->GetType() == "gradient" || (this->GetType() == "hessian"))
3918  {
3919  col = control_dirichlet_colors_.size();
3920  for (unsigned int i = 0; i < control_dirichlet_colors_.size(); ++i)
3921  {
3922  if (control_dirichlet_colors_[i] == color)
3923  {
3924  col = i;
3925  break;
3926  }
3927  }
3928  if (col == control_dirichlet_colors_.size())
3929  {
3930  std::stringstream s;
3931  s << "ControlDirichletColor" << color << " has not been found !";
3932  throw DOpEException(s.str(),
3933  "OptProblemContainer::GetDirichletValues");
3934  }
3935  }
3936  else
3937  {
3938  throw DOpEException("Unknown Type:" + this->GetType(),
3939  "OptProblemContainer::GetDirichletValues");
3940  }
3941 
3942  if (this->GetType() == "state")
3943  {
3944  primal_dirichlet_values_[col]->ReInit(param_values, domain_values,
3945  color);
3946  return *(primal_dirichlet_values_[col]);
3947  }
3948  else if (this->GetType() == "tangent")
3949  {
3950  tangent_dirichlet_values_[col]->ReInit(param_values, domain_values,
3951  color);
3952  return *(tangent_dirichlet_values_[col]);
3953  }
3954  else if (this->GetType() == "adjoint"
3955  || this->GetType() == "adjoint_for_ee"
3956  || (this->GetType() == "adjoint_hessian"))
3957  {
3958  return *(zero_dirichlet_values_);
3959  }
3960  else if (this->GetType() == "gradient" || (this->GetType() == "hessian"))
3961  {
3962  return *(control_dirichlet_values_[col]);
3963  }
3964  else
3965  {
3966  throw DOpEException("Unknown Type:" + this->GetType(),
3967  "OptProblemContainer::GetDirichletValues");
3968  }
3969  }
3970  /******************************************************/
3971 
3972  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
3973  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
3974  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
3975  template<int, int> class DH>
3976  const TransposedDirichletDataInterface<dealdim>&
3977  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
3978  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::GetTransposedDirichletValues(
3979  unsigned int color,
3980  const std::map<std::string, const dealii::Vector<double>*> &param_values,
3981  const std::map<std::string, const VECTOR*> &domain_values) const
3982  {
3983  unsigned int col = control_transposed_dirichlet_colors_.size();
3984  if (this->GetType() == "gradient" || (this->GetType() == "hessian"))
3985  {
3986  for (unsigned int i = 0;
3987  i < control_transposed_dirichlet_colors_.size(); ++i)
3988  {
3989  if (control_transposed_dirichlet_colors_[i] == color)
3990  {
3991  col = i;
3992  break;
3993  }
3994  }
3995  if (col == control_transposed_dirichlet_colors_.size())
3996  {
3997  std::stringstream s;
3998  s << "TransposedControlDirichletColor" << color
3999  << " has not been found !";
4000  throw DOpEException(s.str(),
4001  "OptProblemContainer::GetTransposedDirichletValues");
4002  }
4003  }
4004  else
4005  {
4006  throw DOpEException("Unknown Type:" + this->GetType(),
4007  "OptProblemContainer::GetTransposedDirichletValues");
4008  }
4009 
4010  if (this->GetType() == "gradient")
4011  {
4012  transposed_control_gradient_dirichlet_values_[col]->ReInit(param_values,
4013  domain_values, color);
4014  return *(transposed_control_gradient_dirichlet_values_[col]);
4015  }
4016  else if (this->GetType() == "hessian")
4017  {
4018  transposed_control_hessian_dirichlet_values_[col]->ReInit(param_values,
4019  domain_values, color);
4020  return *(transposed_control_hessian_dirichlet_values_[col]);
4021  }
4022  else
4023  {
4024  throw DOpEException("Unknown Type:" + this->GetType(),
4025  "OptProblemContainer::GetTransposedDirichletValues");
4026  }
4027  }
4028 
4029  /******************************************************/
4030 
4031  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
4032  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
4033  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
4034  template<int, int> class DH>
4035  const std::vector<unsigned int>&
4036  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
4037  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::GetBoundaryEquationColors() const
4038  {
4039  if (this->GetType() == "state" || this->GetType() == "tangent")
4040  {
4041  return state_boundary_equation_colors_;
4042  }
4043  else if (this->GetType() == "adjoint"
4044  || this->GetType() == "adjoint_for_ee"
4045  || this->GetType() == "adjoint_hessian")
4046  {
4047  return adjoint_boundary_equation_colors_;
4048  }
4049  else if (this->GetType() == "gradient" || (this->GetType() == "hessian")
4050  || (this->GetType() == "global_constraint_gradient"))
4051  {
4052  return control_boundary_equation_colors_;
4053  }
4054  else
4055  {
4056  throw DOpEException("Unknown Type:" + this->GetType(),
4057  "OptProblemContainer::GetBoundaryEquationColors");
4058  }
4059  }
4060 
4061  /******************************************************/
4062 
4063  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
4064  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
4065  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
4066  template<int, int> class DH>
4067  void
4068  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
4069  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::SetControlBoundaryEquationColors(
4070  unsigned int color)
4071  {
4072  { //Control Boundary Equation colors are simply inserted
4073  unsigned int comp = control_boundary_equation_colors_.size();
4074  for (unsigned int i = 0; i < control_boundary_equation_colors_.size();
4075  ++i)
4076  {
4077  if (control_boundary_equation_colors_[i] == color)
4078  {
4079  comp = i;
4080  break;
4081  }
4082  }
4083  if (comp != control_boundary_equation_colors_.size())
4084  {
4085  std::stringstream s;
4086  s << "Boundary Equation Color" << color
4087  << " has multiple occurences !";
4088  throw DOpEException(s.str(),
4089  "OptProblemContainer::SetControlBoundaryEquationColors");
4090  }
4091  control_boundary_equation_colors_.push_back(color);
4092  }
4093  }
4094  /******************************************************/
4095 
4096  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
4097  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
4098  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
4099  template<int, int> class DH>
4100  void
4101  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
4102  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::SetBoundaryEquationColors(
4103  unsigned int color)
4104  {
4105  { //State Boundary Equation colors are simply inserted
4106  unsigned int comp = state_boundary_equation_colors_.size();
4107  for (unsigned int i = 0; i < state_boundary_equation_colors_.size();
4108  ++i)
4109  {
4110  if (state_boundary_equation_colors_[i] == color)
4111  {
4112  comp = i;
4113  break;
4114  }
4115  }
4116  if (comp != state_boundary_equation_colors_.size())
4117  {
4118  std::stringstream s;
4119  s << "Boundary Equation Color" << color
4120  << " has multiple occurences !";
4121  throw DOpEException(s.str(),
4122  "OptProblemContainer::SetBoundaryEquationColors");
4123  }
4124  state_boundary_equation_colors_.push_back(color);
4125  }
4126  { //For the adjoint they are added with the boundary functional colors
4127  unsigned int comp = adjoint_boundary_equation_colors_.size();
4128  for (unsigned int i = 0; i < adjoint_boundary_equation_colors_.size();
4129  ++i)
4130  {
4131  if (adjoint_boundary_equation_colors_[i] == color)
4132  {
4133  comp = i;
4134  break;
4135  }
4136  }
4137  if (comp != adjoint_boundary_equation_colors_.size())
4138  {
4139  //Seems this color is already added, however it might have been a functional color
4140  //so we don't do anything.
4141  }
4142  else
4143  {
4144  adjoint_boundary_equation_colors_.push_back(color);
4145  }
4146  }
4147  }
4148 
4149  /******************************************************/
4150 
4151  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
4152  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
4153  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
4154  template<int, int> class DH>
4155  const std::vector<unsigned int>&
4156  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
4157  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::GetBoundaryFunctionalColors() const
4158  {
4159  if (this->GetType() == "cost_functional"
4160  || this->GetType() == "aux_functional"
4161  || this->GetType() == "functional_for_ee") //fixme: what about error_evaluation?
4162  {
4163  return boundary_functional_colors_;
4164  }
4165  else
4166  {
4167  throw DOpEException("Unknown Type:" + this->GetType(),
4168  "OptProblemContainer::GetBoundaryFunctionalColors");
4169  }
4170  }
4171 
4172  /******************************************************/
4173 
4174  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
4175  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
4176  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
4177  template<int, int> class DH>
4178  void
4179  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
4180  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::SetBoundaryFunctionalColors(
4181  unsigned int color)
4182  {
4183  { //Boundary Functional colors are simply inserted
4184  unsigned int comp = boundary_functional_colors_.size();
4185  for (unsigned int i = 0; i < boundary_functional_colors_.size(); ++i)
4186  {
4187  if (boundary_functional_colors_[i] == color)
4188  {
4189  comp = i;
4190  break;
4191  }
4192  }
4193  if (comp != boundary_functional_colors_.size())
4194  {
4195  std::stringstream s;
4196  s << "Boundary Functional Color" << color
4197  << " has multiple occurences !";
4198  throw DOpEException(s.str(),
4199  "OptProblemContainer::SetBoundaryFunctionalColors");
4200  }
4201  boundary_functional_colors_.push_back(color);
4202  }
4203  { //For the adjoint they are addeed to the boundary equation colors
4204  unsigned int comp = adjoint_boundary_equation_colors_.size();
4205  for (unsigned int i = 0; i < adjoint_boundary_equation_colors_.size();
4206  ++i)
4207  {
4208  if (adjoint_boundary_equation_colors_[i] == color)
4209  {
4210  comp = i;
4211  break;
4212  }
4213  }
4214  if (comp != adjoint_boundary_equation_colors_.size())
4215  {
4216  //Seems this color is already added, however it might have been a equation color
4217  //so we don't do anything.
4218  }
4219  else
4220  {
4221  adjoint_boundary_equation_colors_.push_back(color);
4222  }
4223  }
4224  }
4225 
4226  /******************************************************/
4227 
4228  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
4229  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
4230  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
4231  template<int, int> class DH>
4232  unsigned int
4233  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
4234  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::GetControlNBlocks() const
4235  {
4236  return this->GetPDE().GetControlNBlocks();
4237  }
4238 
4239  /******************************************************/
4240 
4241  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
4242  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
4243  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
4244  template<int, int> class DH>
4245  unsigned int
4246  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
4247  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::GetStateNBlocks() const
4248  {
4249  return this->GetPDE().GetStateNBlocks();
4250  }
4251 
4252  /******************************************************/
4253 
4254  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
4255  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
4256  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
4257  template<int, int> class DH>
4258  unsigned int
4259  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
4260  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::GetNBlocks() const
4261  {
4262  if ((this->GetType() == "state") || (this->GetType() == "adjoint_for_ee")
4263  || (this->GetType() == "adjoint") || (this->GetType() == "tangent")
4264  || (this->GetType() == "adjoint_hessian"))
4265  {
4266  return this->GetStateNBlocks();
4267  }
4268  else if ((this->GetType() == "gradient")
4269  || (this->GetType() == "hessian"))
4270  {
4271  return this->GetControlNBlocks();
4272  }
4273  else
4274  {
4275  throw DOpEException("Unknown Type:" + this->GetType(),
4276  "OptProblemContainer::GetNBlocks");
4277  }
4278  }
4279 
4280  /******************************************************/
4281 
4282  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
4283  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
4284  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
4285  template<int, int> class DH>
4286  unsigned int
4287  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
4288  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::GetDoFsPerBlock(
4289  unsigned int b) const
4290  {
4291  if ((this->GetType() == "state") || (this->GetType() == "adjoint")
4292  || (this->GetType() == "adjoint_for_ee")
4293  || (this->GetType() == "tangent")
4294  || (this->GetType() == "adjoint_hessian"))
4295  {
4296  return GetSpaceTimeHandler()->GetStateDoFsPerBlock(b);
4297  }
4298  else if ((this->GetType() == "gradient")
4299  || (this->GetType() == "hessian"))
4300  {
4301  return GetSpaceTimeHandler()->GetControlDoFsPerBlock(b);
4302  }
4303  else
4304  {
4305  throw DOpEException("Unknown Type:" + this->GetType(),
4306  "OptProblemContainer::GetDoFsPerBlock");
4307  }
4308  }
4309 
4310  /******************************************************/
4311 
4312  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
4313  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
4314  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
4315  template<int, int> class DH>
4316  const std::vector<unsigned int>&
4317  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
4318  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::GetDoFsPerBlock() const
4319  {
4320  if ((this->GetType() == "state") || (this->GetType() == "adjoint")
4321  || (this->GetType() == "adjoint_for_ee")
4322  || (this->GetType() == "tangent")
4323  || (this->GetType() == "adjoint_hessian"))
4324  {
4325  return GetSpaceTimeHandler()->GetStateDoFsPerBlock();
4326  }
4327  else if ((this->GetType() == "gradient")
4328  || (this->GetType() == "hessian"))
4329  {
4330  return GetSpaceTimeHandler()->GetControlDoFsPerBlock();
4331  }
4332  else
4333  {
4334  throw DOpEException("Unknown Type:" + this->GetType(),
4335  "OptProblemContainer::GetDoFsPerBlock");
4336  }
4337  }
4338 
4339  /******************************************************/
4340 
4341  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
4342  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
4343  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
4344  template<int, int> class DH>
4345  const dealii::ConstraintMatrix&
4346  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
4347  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::GetDoFConstraints() const
4348  {
4349  if ((this->GetType() == "state") || (this->GetType() == "adjoint")
4350  || (this->GetType() == "adjoint_for_ee")
4351  || (this->GetType() == "tangent")
4352  || (this->GetType() == "adjoint_hessian"))
4353  {
4354  return GetSpaceTimeHandler()->GetStateDoFConstraints();
4355  }
4356  else if ((this->GetType() == "gradient") || (this->GetType() == "hessian")
4357  || (this->GetType() == "global_constraint_gradient"))
4358  {
4359  return GetSpaceTimeHandler()->GetControlDoFConstraints();
4360  }
4361  else
4362  {
4363  throw DOpEException("Unknown Type:" + this->GetType(),
4364  "OptProblemContainer::GetDoFConstraints");
4365  }
4366  }
4367 
4368  /******************************************************/
4369 
4370  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
4371  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
4372  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
4373  template<int, int> class DH>
4374  bool
4375  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
4376  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::NeedTimeFunctional() const
4377  {
4378  if (this->GetType() == "cost_functional")
4379  return GetFunctional()->NeedTime();
4380  else if (this->GetType() == "aux_functional")
4381  return aux_functionals_[this->GetTypeNum()]->NeedTime();
4382  else if (this->GetType() == "functional_for_ee")
4383  return aux_functionals_[functional_for_ee_num_]->NeedTime();
4384  else
4385  throw DOpEException("Not implemented",
4386  "OptProblemContainer::NeedTimeFunctional");
4387  }
4388 
4389  /******************************************************/
4390 
4391  template<typename FUNCTIONAL_INTERFACE, typename FUNCTIONAL, typename PDE,
4392  typename DD, typename CONSTRAINTS, typename SPARSITYPATTERN,
4393  typename VECTOR, int dopedim, int dealdim, template<int, int> class FE,
4394  template<int, int> class DH>
4395  bool
4396  OptProblemContainer<FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS,
4397  SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH>::HasControlInDirichletData() const
4398  {
4399  return (!control_transposed_dirichlet_colors_.empty());
4400  }
4401 
4402 }
4403 #endif
DOpE::OptProblemContainer::ElementTimeMatrix
void ElementTimeMatrix(const DATACONTAINER &dc, dealii::FullMatrix< double > &local_entry_matrix)
DOpE::OptProblemContainer::ElementMatrix
void ElementMatrix(const DATACONTAINER &dc, dealii::FullMatrix< double > &local_entry_matrix, double scale=1., double scale_ico=1.)
Definition: optproblemcontainer.h:2733
DOpE::OptProblemContainer::DeleteAuxiliaryState
void DeleteAuxiliaryState(std::string name)
Definition: optproblemcontainer.h:3439
DOpE::OptProblemContainer::GetOutputHandler
DOpEOutputHandler< VECTOR > * GetOutputHandler()
Definition: optproblemcontainer.h:822
DOpE::OptProblemContainer::Init_PointRhs
void Init_PointRhs(const std::map< std::string, const dealii::Vector< double > * > &param_values, const std::map< std::string, const VECTOR * > &domain_values, VECTOR &rhs_vector, double scale=1.)
Definition: optproblemcontainer.h:1142
DOpE::OptProblemContainer::OptProblemContainer
OptProblemContainer(FUNCTIONAL &functional, PDE &pde, CONSTRAINTS &constraints, SpaceTimeHandler< FE, DH, SPARSITYPATTERN, VECTOR, dopedim, dealdim > &STH)
Definition: optproblemcontainer.h:1324
DOpE::OptProblemContainer::BoundaryEquation
void BoundaryEquation(const FACEDATACONTAINER &dc, dealii::Vector< double > &local_vector, double scale, double scale_ico)
Definition: optproblemcontainer.h:1968
DOpE::OptProblemContainer::GetUpdateFlags
dealii::UpdateFlags GetUpdateFlags() const
Definition: optproblemcontainer.h:3050
DOpE::StateProblem
Definition: stateproblem.h:48
DOpE::OptProblemContainer::GetAuxiliaryConstraint
const ConstraintVector< VECTOR > * GetAuxiliaryConstraint(std::string name)
Definition: optproblemcontainer.h:911
DOpE::OptProblemContainer::GetNFunctionals
unsigned int GetNFunctionals() const
Definition: optproblemcontainer.h:754
DOpE::OptProblemContainer::FaceFunctional
double FaceFunctional(const FACEDATACONTAINER &fdc)
Definition: optproblemcontainer.h:1647
DOpE::OptProblemContainer::Init_ElementEquation
void Init_ElementEquation(const DATACONTAINER &edc, dealii::Vector< double > &local_vector, double scale, double scale_ico)
Definition: optproblemcontainer.h:1077
DOpE::OptProblemContainer::RegisterOutputHandler
void RegisterOutputHandler(DOpEOutputHandler< VECTOR > *OH)
Definition: optproblemcontainer.h:164
DOpE::OptProblemContainer::BoundaryRhs
void BoundaryRhs(const FACEDATACONTAINER &dc, dealii::Vector< double > &local_vector, double scale=1.)
Definition: optproblemcontainer.h:2559
DOpE::ConstraintVector
Definition: constraintvector.h:48
DOpE::OptProblemContainer::AlgebraicFunctional
double AlgebraicFunctional(const std::map< std::string, const dealii::Vector< double > * > &values, const std::map< std::string, const VECTOR * > &block_values)
Definition: optproblemcontainer.h:1680
DOpE::OptProblemContainer::ComputeLocalControlConstraints
void ComputeLocalControlConstraints(VECTOR &constraints, const std::map< std::string, const dealii::Vector< double > * > &values, const std::map< std::string, const VECTOR * > &block_values)
Definition: optproblemcontainer.h:2958
DOpE::OptProblemContainer::ElementTimeMatrixExplicit
void ElementTimeMatrixExplicit(const DATACONTAINER &dc, dealii::FullMatrix< double > &local_entry_matrix)
Definition: optproblemcontainer.h:2813
DOpE::OptProblemContainer::ElementTimeEquationExplicit
void ElementTimeEquationExplicit(const DATACONTAINER &dc, dealii::Vector< double > &local_vector, double scale=1.)
Definition: optproblemcontainer.h:1838
DOpE::OptProblemContainer::GetBoundaryFunctionalColors
const std::vector< unsigned int > & GetBoundaryFunctionalColors() const
Definition: optproblemcontainer.h:4157
DOpE::OptProblemContainer::SetTime
void SetTime(double time, const TimeIterator &interval, bool initial=false)
Definition: optproblemcontainer.h:3194
DOpE::OptProblemContainer::ElementRhs
void ElementRhs(const DATACONTAINER &dc, dealii::Vector< double > &local_vector, double scale=1.)
Definition: optproblemcontainer.h:2019
DOpE::OptProblemContainer::ReInit
void ReInit(std::string algo_type)
Definition: optproblemcontainer.h:1392
DOpE::OptProblemContainer::GetBlockVector
const VECTOR * GetBlockVector(const std::map< std::string, const VECTOR * > &values, std::string name)
Definition: optproblemcontainer.h:1238
DOpE::DOpEtypes::initial
Definition: dopetypes.h:106
DOpE::OptProblemContainer::GetInitialValues
const dealii::Function< dealdim > & GetInitialValues() const
Definition: optproblemcontainer.h:662
DOpE::OptProblemContainer::GetFESystem
const FE< dealdim, dealdim > & GetFESystem() const
Definition: optproblemcontainer.h:3015
DOpE::OptProblemContainer::PostProcessConstraints
void PostProcessConstraints(ConstraintVector< VECTOR > &g) const
Definition: optproblemcontainer.h:3301
DOpE::OptProblemContainer::GetFaceUpdateFlags
dealii::UpdateFlags GetFaceUpdateFlags() const
Definition: optproblemcontainer.h:3095
DOpE::OptProblemContainer::GetFunctionalType
std::string GetFunctionalType() const
Definition: optproblemcontainer.h:3139
DOpE::OptProblemContainer::AddFunctional
void AddFunctional(FUNCTIONAL_INTERFACE *F)
Definition: optproblemcontainer.h:692
DOpE::OptProblemContainer::GetFunctional
FUNCTIONAL * GetFunctional()
Definition: optproblemcontainer.h:3483
DOpE::TangentDirichletData
Definition: tangentdirichletdata.h:39
DOpE::OptProblemContainer::HasInterfaces
bool HasInterfaces() const
Definition: optproblemcontainer.h:3612
DOpE::OptProblemContainer::BoundaryFunctional
double BoundaryFunctional(const FACEDATACONTAINER &fdc)
Definition: optproblemcontainer.h:1609
DOpE::OptProblemContainer::ElementFunctional
double ElementFunctional(const DATACONTAINER &edc)
Definition: optproblemcontainer.h:1519
DOpE::ProblemContainerInternal::GetType
std::string GetType() const
Definition: problemcontainer_internal.h:109
DOpE::OptProblemContainer::GetBoundaryEquationColors
const std::vector< unsigned int > & GetBoundaryEquationColors() const
Definition: optproblemcontainer.h:4037
DOpE::OptProblemContainer::SetDirichletBoundaryColors
void SetDirichletBoundaryColors(unsigned int color, const std::vector< bool > &comp_mask, const DD *values)
Definition: optproblemcontainer.h:3687
DOpE::OptProblemContainer::GetDoFsPerBlock
const std::vector< unsigned int > & GetDoFsPerBlock() const
Definition: optproblemcontainer.h:4318
DOpE::OptProblemContainer::SetControlBoundaryEquationColors
void SetControlBoundaryEquationColors(unsigned int color)
Definition: optproblemcontainer.h:4069
DOpE::OptProblemContainer::PointRhs
void PointRhs(const std::map< std::string, const dealii::Vector< double > * > &param_values, const std::map< std::string, const VECTOR * > &domain_values, VECTOR &rhs_vector, double scale=1.)
Definition: optproblemcontainer.h:2207
DOpE::ProblemContainerInternal
Definition: problemcontainer_internal.h:38
DOpE::DOpEOutputHandler
Definition: optproblemcontainer.h:70
DOpE::OptProblemContainer::GetControlBoxConstraints
void GetControlBoxConstraints(VECTOR &lb, VECTOR &ub) const
Definition: optproblemcontainer.h:566
DOpE::OptProblemContainer::AddAuxiliaryControl
void AddAuxiliaryControl(const ControlVector< VECTOR > *c, std::string name)
Definition: optproblemcontainer.h:3314
DOpE::OptProblemContainer::GetDirichletColors
const std::vector< unsigned int > & GetDirichletColors() const
Definition: optproblemcontainer.h:3741
DOpE::OptProblemContainer::HasPoints
bool HasPoints() const
Definition: optproblemcontainer.h:3572
DOpE::OptProblemContainer::GetDoFType
std::string GetDoFType() const
Definition: optproblemcontainer.h:2987
DOpE::OptProblemContainer::GetName
virtual std::string GetName() const
Definition: optproblemcontainer.h:115
DOpE::OptProblemContainer::ElementEquation
void ElementEquation(const DATACONTAINER &edc, dealii::Vector< double > &local_vector, double scale, double scale_ico)
Definition: optproblemcontainer.h:1717
DOpE::OptProblemContainer::ElementTimeEquation
void ElementTimeEquation(const DATACONTAINER &dc, dealii::Vector< double > &local_vector, double scale=1.)
Definition: optproblemcontainer.h:1794
DOpE::OptProblemContainer::GetDirichletCompMask
const std::vector< bool > & GetDirichletCompMask(unsigned int color) const
Definition: optproblemcontainer.h:3790
DOpE::OptProblemContainer::SetType
void SetType(std::string type, unsigned int num=0)
Definition: optproblemcontainer.h:1434
DOpE::OptProblemContainer::PointFunctional
double PointFunctional(const std::map< std::string, const dealii::Vector< double > * > &param_values, const std::map< std::string, const VECTOR * > &domain_values)
Definition: optproblemcontainer.h:1555
DOpE::TimeIterator
Definition: timeiterator.h:63
DOpE::DOpEtypes::stationary
Definition: dopetypes.h:105
DOpE::OptProblemContainer::GetNBlocks
unsigned int GetNBlocks() const
Definition: optproblemcontainer.h:4260
DOpE::OptProblemContainer::GetAuxiliaryState
const StateVector< VECTOR > * GetAuxiliaryState(std::string name) const
Definition: optproblemcontainer.h:3396
DOpE::OptProblemContainer::GetDirichletValues
const dealii::Function< dealdim > & GetDirichletValues(unsigned int color, const std::map< std::string, const dealii::Vector< double > * > &param_values, const std::map< std::string, const VECTOR * > &domain_values) const
Definition: optproblemcontainer.h:3889
DOpE::ControlVector
Definition: controlvector.h:49
DOpE::OptProblemContainer::Init_ElementMatrix
void Init_ElementMatrix(const DATACONTAINER &edc, dealii::FullMatrix< double > &local_entry_matrix, double scale, double scale_ico)
Definition: optproblemcontainer.h:1192
DOpE::OptProblemContainer::GetStateProblem
StateProblem< OptProblemContainer< FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH >, PDE, DD, SPARSITYPATTERN, VECTOR, dealdim > & GetStateProblem()
Definition: optproblemcontainer.h:128
DOpE::ProblemContainerInternal::GetPDE
const PDE & GetPDE() const
Definition: problemcontainer_internal.h:103
DOpE::OptProblemContainer::EEFunctionalIsCost
bool EEFunctionalIsCost() const
Definition: optproblemcontainer.h:1216
DOpE::OptProblemContainer::GetDoFConstraints
const dealii::ConstraintMatrix & GetDoFConstraints() const
Definition: optproblemcontainer.h:4347
DOpE::OptProblemContainer::GetConstraintType
std::string GetConstraintType() const
Definition: optproblemcontainer.h:3181
DOpE::OptProblemContainer
Definition: optproblemcontainer.h:101
problemcontainer_internal.h
DOpE::OptProblemContainer::GetStateNBlocks
unsigned int GetStateNBlocks()
Definition: optproblemcontainer.h:1052
DOpE::OptProblemContainer::FaceEquation
void FaceEquation(const FACEDATACONTAINER &dc, dealii::Vector< double > &local_vector, double scale, double scale_ico)
Definition: optproblemcontainer.h:1886
DOpE::OptProblemContainer::AddAuxiliaryConstraint
void AddAuxiliaryConstraint(const ConstraintVector< VECTOR > *c, std::string name)
Definition: optproblemcontainer.h:3355
DOpE::OptProblemContainer::GetAuxiliaryControl
const ControlVector< VECTOR > * GetAuxiliaryControl(std::string name) const
Definition: optproblemcontainer.h:3375
DOpE::OptProblemContainer::GetExceptionHandler
DOpEExceptionHandler< VECTOR > * GetExceptionHandler()
Definition: optproblemcontainer.h:813
DOpE::OptProblemContainer::FaceRhs
void FaceRhs(const FACEDATACONTAINER &dc, dealii::Vector< double > &local_vector, double scale=1.)
Definition: optproblemcontainer.h:2390
DOpE::OptProblemContainer::GetSpaceTimeHandler
const SpaceTimeHandler< FE, DH, SPARSITYPATTERN, VECTOR, dopedim, dealdim > * GetSpaceTimeHandler() const
Definition: optproblemcontainer.h:844
DOpE::OptProblemContainer::RegisterExceptionHandler
void RegisterExceptionHandler(DOpEExceptionHandler< VECTOR > *OH)
Definition: optproblemcontainer.h:172
DOpE::TransposedHessianDirichletData
Definition: transposedhessiandirichletdata.h:41
DOpE::OptProblemContainer::GetConstraints
CONSTRAINTS * GetConstraints()
Definition: optproblemcontainer.h:1227
DOpE::OptProblemContainer::FaceMatrix
void FaceMatrix(const FACEDATACONTAINER &dc, dealii::FullMatrix< double > &local_entry_matrix, double scale=1., double scale_ico=1.)
DOpE::OptProblemContainer::DeleteAuxiliaryFromIntegrator
void DeleteAuxiliaryFromIntegrator(INTEGRATOR &integrator)
Definition: optproblemcontainer.h:996
DOpE::OptProblemContainer::GetVector
const dealii::Vector< double > * GetVector(const std::map< std::string, const Vector< double > * > &values, std::string name)
Definition: optproblemcontainer.h:1251
DOpE::OptProblemContainer::GetSpaceTimeHandler
SpaceTimeHandler< FE, DH, SPARSITYPATTERN, VECTOR, dopedim, dealdim > * GetSpaceTimeHandler()
Definition: optproblemcontainer.h:852
DOpE::OptProblemContainer::GetConstraints
const CONSTRAINTS * GetConstraints() const
Definition: optproblemcontainer.h:1232
DOpE::OptProblemContainer::HasFaces
bool HasFaces() const
Definition: optproblemcontainer.h:3524
DOpE::PrimalDirichletData
Definition: primaldirichletdata.h:41
DOpE::OptProblemContainer::GetTransposedDirichletColors
const std::vector< unsigned int > & GetTransposedDirichletColors() const
Definition: optproblemcontainer.h:3769
DOpE::StateVector
Definition: statevector.h:50
DOpE::OptProblemContainer::GetTransposedDirichletValues
const TransposedDirichletDataInterface< dealdim > & GetTransposedDirichletValues(unsigned int color, const std::map< std::string, const dealii::Vector< double > * > &param_values, const std::map< std::string, const VECTOR * > &domain_values) const
Definition: optproblemcontainer.h:3978
DOpE::OptProblemContainer::ComputeSparsityPattern
void ComputeSparsityPattern(SPARSITYPATTERN &sparsity) const
Definition: optproblemcontainer.h:3267
DOpE::OptProblemContainer::AlgebraicResidual
void AlgebraicResidual(VECTOR &residual, const std::map< std::string, const dealii::Vector< double > * > &values, const std::map< std::string, const VECTOR * > &block_values)
Definition: optproblemcontainer.h:1767
DOpE::OptProblemContainer::InterfaceEquation
void InterfaceEquation(const FACEDATACONTAINER &dc, dealii::Vector< double > &local_vector, double scale, double scale_ico)
Definition: optproblemcontainer.h:1936
DOpE::SpaceTimeHandler
Definition: spacetimehandler.h:71
fevalues_wrapper.h
DOpE::OptProblemContainer::GetBaseProblem
OptProblemContainer< FUNCTIONAL_INTERFACE, FUNCTIONAL, PDE, DD, CONSTRAINTS, SPARSITYPATTERN, VECTOR, dopedim, dealdim, FE, DH > & GetBaseProblem()
Definition: optproblemcontainer.h:145
DOpE::OptProblemContainer::SetFunctionalForErrorEstimation
void SetFunctionalForErrorEstimation(std::string functional_name)
Definition: optproblemcontainer.h:718
DOpE::OptProblemContainer::GetStateNBlocks
unsigned int GetStateNBlocks() const
Definition: optproblemcontainer.h:4247
DOpE::OptProblemContainer::GetControlBlockComponent
std::vector< unsigned int > & GetControlBlockComponent()
Definition: optproblemcontainer.h:1060
DOpE::OptProblemContainer::AddAuxiliaryState
void AddAuxiliaryState(const StateVector< VECTOR > *c, std::string name)
Definition: optproblemcontainer.h:3335
DOpE::OptProblemContainer::NeedTimeFunctional
bool NeedTimeFunctional() const
Definition: optproblemcontainer.h:4376
dofhandler_wrapper.h
DOpE::OptProblemContainer::GetFunctionalName
std::string GetFunctionalName() const
Definition: optproblemcontainer.h:3160
function_wrapper.h
DOpE::TransposedDirichletDataInterface
Definition: transposeddirichletdatainterface.h:37
DOpE::OptProblemContainer::GetControlNBlocks
unsigned int GetControlNBlocks() const
Definition: optproblemcontainer.h:4234
DOpE::OptProblemContainer::GetTransposedDirichletCompMask
const std::vector< bool > & GetTransposedDirichletCompMask(unsigned int color) const
Definition: optproblemcontainer.h:3851
DOpE::OptProblemContainer::~OptProblemContainer
virtual ~OptProblemContainer()
Definition: optproblemcontainer.h:1349
DOpE::OptProblemContainer::HasControlInDirichletData
bool HasControlInDirichletData() const
Definition: optproblemcontainer.h:4397
DOpE::OptProblemContainer::GetStateBlockComponent
std::vector< unsigned int > & GetStateBlockComponent()
Definition: optproblemcontainer.h:1067
DOpE::OptProblemContainer::SetInitialValues
void SetInitialValues(const dealii::Function< dealdim > *values)
Definition: optproblemcontainer.h:657
DOpE::DOpEException
Definition: dopeexception.h:35
DOpE::TransposedGradientDirichletData
Definition: transposedgradientdirichletdata.h:41
DOpE::OptProblemContainer::AddAuxiliaryToIntegrator
void AddAuxiliaryToIntegrator(INTEGRATOR &integrator)
Definition: optproblemcontainer.h:933
DOpE::OptProblemContainer::SetBoundaryEquationColors
void SetBoundaryEquationColors(unsigned int color)
Definition: optproblemcontainer.h:4102
DOpE::OptProblemContainer::SetControlDirichletBoundaryColors
void SetControlDirichletBoundaryColors(unsigned int color, const std::vector< bool > &comp_mask, const DOpEWrapper::Function< dealdim > *values)
Definition: optproblemcontainer.h:3652
DOpE::OptProblemContainer::DeleteAuxiliaryConstraint
void DeleteAuxiliaryConstraint(std::string name)
Definition: optproblemcontainer.h:3461
DOpE::OptProblemContainer::GetFunctionalPosition
const std::map< std::string, unsigned int > & GetFunctionalPosition() const
Definition: optproblemcontainer.h:1046
DOpE::OptProblemContainer::SetBoundaryFunctionalColors
void SetBoundaryFunctionalColors(unsigned int color)
Definition: optproblemcontainer.h:4180
DOpE::OptProblemContainer::DeleteAuxiliaryControl
void DeleteAuxiliaryControl(std::string name)
Definition: optproblemcontainer.h:3417
DOpE::OptProblemContainer::BoundaryMatrix
void BoundaryMatrix(const FACEDATACONTAINER &dc, dealii::FullMatrix< double > &local_matrix, double scale=1., double scale_ico=1.)
DOpE::OptProblemContainer::Init_ElementRhs
void Init_ElementRhs(const DATACONTAINER &edc, dealii::Vector< double > &local_vector, double scale)
Definition: optproblemcontainer.h:1096
DOpE::OptProblemContainer::InterfaceMatrix
void InterfaceMatrix(const FACEDATACONTAINER &dc, dealii::FullMatrix< double > &local_entry_matrix, double scale=1., double scale_ico=1.)
DOpE::DOpEExceptionHandler
Definition: optproblemcontainer.h:72
Generated on Fri Oct 24 2014 20:06:50 for DOpE by   doxygen 1.8.6