33 #ifndef __H__UG__LIB_ALGEBRA__OPERATOR__PRECONDITIONER__PROJECTED_GAUSS_SEIDEL__OBSTACLE_IN_NORMAL_DIR_IMPL__
34 #define __H__UG__LIB_ALGEBRA__OPERATOR__PRECONDITIONER__PROJECTED_GAUSS_SEIDEL__OBSTACLE_IN_NORMAL_DIR_IMPL__
44 template <
typename TDomain,
typename TAlgebra>
51 if (
dim > 1) normal[1] = 1.5;
52 if (
dim > 2) normal[2] = -3.0;
55 transform_eulerian_coord_sys(transformedONB, normal);
59 template <
typename TDomain,
typename TAlgebra>
67 if (
dim > 1) unityX[1] = 0.0;
68 if (
dim > 2) unityX[2] = 0.0;
73 VecScale(normalizedFirstBaseVec, firstTransformedBaseVec, 1.0/normOfFirstBaseVec);
77 VecScaleAdd(orthoVec, 0.5, unityX, -0.5, normalizedFirstBaseVec);
85 for(
size_t i = 0; i < (size_t)
dim; ++i)
87 for(
size_t j = 0; j < (size_t)
dim; ++j){
88 transBaseVec[j] = hMat(i,j);
90 transformedONB[i] = transBaseVec;
122 template <
typename TDomain,
typename TAlgebra>
123 template <
typename TElem,
typename TIterator>
134 vector<MathVector<dim> > vCorner;
135 vector<MathVector<dim> > vCoPos;
140 for(TIterator iter = iterBegin; iter != iterEnd; ++iter)
155 int nCorner = (int)vCorner.size();
156 for (
int i = 0; i < nCorner; ++i)
157 vCoPos.push_back(vCorner[rRefElem.id(
dim-1, 0, 0, i)]);
159 if ((
int)vCorner.size() ==
dim)
160 ElementNormal<face_type0, dim>(normal, &vCoPos[0]);
162 ElementNormal<face_type1, dim>(normal, &vCoPos[0]);
174 const size_t comp = dof[1];
180 const number obsVal = m_mObstacleValues[dof];
184 if (!(tmpSol < obsVal))
187 m_vActiveDofs.push_back(dof);
196 template <
typename TDomain,
typename TAlgebra>
203 for (vector<int>::iterator obsSI = m_vObsSubsets.begin();
204 obsSI != m_vObsSubsets.end(); ++obsSI)
206 const int si = *obsSI;
207 UG_LOG(
"si: " <<si<<
"\n");
215 adjust_sol_and_cor_elem<RegularEdge>
216 (m_spDD->template begin<RegularEdge>(si), m_spDD->template end<RegularEdge>(si),
217 sol_i, c_i, dofIsActive, dof);
220 adjust_sol_and_cor_elem<Triangle>
221 (m_spDD->template begin<Triangle>(si), m_spDD->template end<Triangle>(si),
222 sol_i, c_i, dofIsActive, dof);
223 adjust_sol_and_cor_elem<Quadrilateral>
224 (m_spDD->template begin<Quadrilateral>(si), m_spDD->template end<Quadrilateral>(si),
225 sol_i, c_i, dofIsActive, dof);
228 UG_THROW(
"ObstacleInNormalDir::adjust_sol_and_cor:"
229 "SubsetDimension "<< subsetDim <<
" (subset="<< si <<
") not supported.");
234 template <
typename TDomain,
typename TAlgebra>
239 for (vector<
MultiIndex<2> >::iterator itActiveInd = m_vActiveDofs.begin();
240 itActiveInd < m_vActiveDofs.end(); ++itActiveInd)
244 number defect =
BlockRef(d[(*itActiveInd)[0]], (*itActiveInd)[1]);
246 BlockRef(d[(*itActiveInd)[0]], (*itActiveInd)[1]) = 0.0;
250 template <
typename TDomain,
typename TAlgebra>
Definition: fe_geom.h:183
dimension dependent base class for reference elements
Definition: reference_element.h:183
vector_type::value_type value_type
Value type.
Definition: obstacle_constraint_interface.h:88
algebra_type::vector_type vector_type
Vector type.
Definition: obstacle_constraint_interface.h:85
Definition: multi_index.h:50
void adjust_sol_and_cor(value_type &sol_i, value_type &c_i, bool &dofIsActive, const DoFIndex &dof)
projects the i-th index of the solution onto the admissible set and adjusts the correction
Definition: obstacle_in_normal_dir_impl.h:199
void transform_eulerian_coord_sys(MathVector< dim > transformedONB[], const MathVector< dim > &firstTransformedBaseVec)
Definition: obstacle_in_normal_dir_impl.h:61
void adjust_defect_to_constraint(vector_type &d)
the defect needs to be adjusted for the active indices (those indices, which are in contact)
Definition: obstacle_in_normal_dir_impl.h:237
void preprocess()
preprocess is useful to attach the normals for every obstacle DoF
Definition: obstacle_in_normal_dir_impl.h:46
void restrict_obs_values()
restricts the obstacle values to a coarser grid in a multigrid hierarchy
Definition: obstacle_in_normal_dir_impl.h:253
void adjust_sol_and_cor_elem(TIterator iterBegin, TIterator iterEnd, value_type &sol_i, value_type &c_i, bool &dofIsAdmissible, const DoFIndex &dof)
Definition: obstacle_in_normal_dir_impl.h:126
Definition: reference_element.h:345
static const DimReferenceElement< dim > & get(ReferenceObjectID roid)
returns a dimension dependent Reference Element
Definition: reference_element.h:280
Definition: reference_element.h:445
Definition: reference_element.h:394
Definition: reference_element.h:303
void CollectCornerCoordinates(std::vector< typename TAAPos::ValueType > &vCornerCoordsOut, const TElem &elem, const TAAPos &aaPos, bool clearContainer=true)
returns the corner coordinates of a geometric object
Definition: domain_util_impl.h:75
void MatHouseholder(matrix_t &mOut, const vector_t &orthoVec)
Definition: math_matrix_functions_common_impl.hpp:942
#define UG_THROW(msg)
Definition: error.h:57
#define UG_LOG(msg)
Definition: log.h:367
double number
Definition: types.h:124
vector_t::value_type VecLength(const vector_t &v)
returns the length of v. Slower than VecLengthSq.
Definition: math_vector_functions_common_impl.hpp:341
void VecScaleAdd(vector_t &vOut, typename vector_t::value_type s1, const vector_t &v1, typename vector_t::value_type s2, const vector_t &v2)
Scales two Vectors, adds them and returns the sum in a third vector.
Definition: math_vector_functions_common_impl.hpp:265
void VecScale(vector_t &vOut, const vector_t &v, typename vector_t::value_type s)
scales a MathVector<N>
Definition: math_vector_functions_common_impl.hpp:252
ReferenceObjectID
these ids are used to identify the shape of a geometric object.
Definition: grid_base_objects.h:74
double & BlockRef(T &vec, size_t i)
Definition: blocks.h:66
int DimensionOfSubset(const ISubsetHandler &sh, int si)
returns the current dimension of the subset
Definition: subset_dim_util.cpp:89
ReferenceVertex face_type1
Definition: obstacle_in_normal_dir_impl.h:104
ReferenceVertex face_type0
Definition: obstacle_in_normal_dir_impl.h:103
DimFEGeometry< 1, 1 > DimFEGeo
Definition: obstacle_in_normal_dir_impl.h:105
ReferenceEdge face_type0
Definition: obstacle_in_normal_dir_impl.h:110
ReferenceEdge face_type1
Definition: obstacle_in_normal_dir_impl.h:111
DimFEGeometry< 2, 1 > DimFEGeo
Definition: obstacle_in_normal_dir_impl.h:112
ReferenceTriangle face_type0
Definition: obstacle_in_normal_dir_impl.h:117
ReferenceQuadrilateral face_type1
Definition: obstacle_in_normal_dir_impl.h:118
DimFEGeometry< 3, 2 > DimFEGeo
Definition: obstacle_in_normal_dir_impl.h:119
Definition: active_set_impl.h:45
void face_type0
Definition: obstacle_in_normal_dir_impl.h:96
void face_type1
Definition: obstacle_in_normal_dir_impl.h:97
void DimFEGeo
Definition: obstacle_in_normal_dir_impl.h:98