Line data Source code
1 : /*
2 : * Copyright (c) 2013-2015: G-CSC, Goethe University Frankfurt
3 : * Author: Dmitry Logashenko
4 : *
5 : * This file is part of UG4.
6 : *
7 : * UG4 is free software: you can redistribute it and/or modify it under the
8 : * terms of the GNU Lesser General Public License version 3 (as published by the
9 : * Free Software Foundation) with the following additional attribution
10 : * requirements (according to LGPL/GPL v3 §7):
11 : *
12 : * (1) The following notice must be displayed in the Appropriate Legal Notices
13 : * of covered and combined works: "Based on UG4 (www.ug4.org/license)".
14 : *
15 : * (2) The following notice must be displayed at a prominent place in the
16 : * terminal output of covered works: "Based on UG4 (www.ug4.org/license)".
17 : *
18 : * (3) The following bibliography is recommended for citation and must be
19 : * preserved in all covered files:
20 : * "Reiter, S., Vogel, A., Heppner, I., Rupp, M., and Wittum, G. A massively
21 : * parallel geometric multigrid solver on hierarchically distributed grids.
22 : * Computing and visualization in science 16, 4 (2013), 151-164"
23 : * "Vogel, A., Reiter, S., Rupp, M., Nägel, A., and Wittum, G. UG4 -- a novel
24 : * flexible software system for simulating pde based models on high performance
25 : * computers. Computing and visualization in science 16, 4 (2013), 165-179"
26 : *
27 : * This program is distributed in the hope that it will be useful,
28 : * but WITHOUT ANY WARRANTY; without even the implied warranty of
29 : * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
30 : * GNU Lesser General Public License for more details.
31 : */
32 :
33 : /*
34 : * A linker that projects a given vector to the manifold of given elements.
35 : */
36 : #ifndef __H__UG__LIB_DISC__SPATIAL_DISC__PROJECT_LINKER__
37 : #define __H__UG__LIB_DISC__SPATIAL_DISC__PROJECT_LINKER__
38 :
39 : /* ug4 headers */
40 : #include "common/common.h"
41 :
42 : #include "lib_disc/reference_element/reference_mapping_provider.h"
43 : #include "linker.h"
44 :
45 : namespace ug {
46 :
47 : /**
48 : * User data projecting a given vector filed to the plane of the low-dim. elements.
49 : */
50 : template <int dim>
51 : class ProjectionLinker
52 : : public StdDataLinker<ProjectionLinker<dim>, MathVector<dim>, dim>
53 : {
54 : /// Base class type
55 : typedef StdDataLinker<ProjectionLinker<dim>, MathVector<dim>, dim> base_type;
56 :
57 : public:
58 :
59 : /// Constructor
60 0 : ProjectionLinker
61 : (
62 : SmartPtr<CplUserData<MathVector<dim>, dim> > spVector ///< vector to project
63 : )
64 0 : {
65 : this->set_num_input (1);
66 0 : m_spVector = spVector;
67 0 : m_spDVector = spVector.template cast_dynamic<DependentUserData<MathVector<dim>, dim> > ();
68 0 : this->set_input (0, spVector, spVector);
69 0 : }
70 :
71 : /// Constructor
72 : ProjectionLinker
73 : (
74 : MathVector<dim> vector ///< vector to project
75 : )
76 : {
77 : this->set_num_input (1);
78 : m_spVector = make_sp (new ConstUserVector<dim> (vector));
79 : m_spDVector = m_spVector.template cast_dynamic<DependentUserData<MathVector<dim>, dim> > ();
80 : this->set_input (0, m_spVector, m_spVector);
81 : }
82 :
83 : /// Returns true because without a grid function, we do not get the element to project to!
84 0 : virtual bool requires_grid_fct() const {return true;}
85 :
86 : /// Evaluation with no element is impossible
87 0 : inline void evaluate
88 : (
89 : MathVector<dim>& value,
90 : const MathVector<dim>& glob_ip,
91 : number time,
92 : int si
93 : ) const
94 : {
95 0 : UG_THROW ("ProjectionLinker: Cannot evaluate without any specification of the element!");
96 : }
97 :
98 : /// Computation only of the projection
99 : template <int refDim>
100 0 : inline void evaluate
101 : (
102 : MathVector<dim> vValue[],
103 : const MathVector<dim> vGlobIP[],
104 : number time,
105 : int si,
106 : GridObject* elem,
107 : const MathVector<dim> vCornerCoords[],
108 : const MathVector<refDim> vLocIP[],
109 : const size_t nip,
110 : LocalVector* u,
111 : const MathMatrix<refDim, dim>* vJT = NULL
112 : ) const
113 : {
114 : // 1. Get the vectors themselves
115 0 : (*m_spVector) (vValue, vGlobIP, time, si, elem, vCornerCoords, vLocIP, nip, u, vJT);
116 :
117 : // 2. Prepare the data of the element
118 0 : const ReferenceObjectID roid = elem->reference_object_id ();
119 0 : DimReferenceMapping<refDim, dim>& rMapping = ReferenceMappingProvider::get<refDim, dim> (roid);
120 0 : rMapping.update (vCornerCoords);
121 :
122 : // 3. At every integration point
123 0 : for (size_t i = 0; i < nip; i++)
124 : {
125 : // 3a. Get the Jacobian
126 : MathMatrix<dim, refDim> J;
127 0 : rMapping.jacobian (J, vLocIP[i]);
128 : // 3b. Project the vector to the subspace spanned by the columns of the Jacobian
129 0 : OrthogProjectVec (vValue[i], J);
130 : }
131 0 : }
132 :
133 : /// Computation of the projection and its derivatives
134 : template <int refDim>
135 0 : void eval_and_deriv
136 : (
137 : MathVector<dim> vValue[],
138 : const MathVector<dim> vGlobIP[],
139 : number time,
140 : int si,
141 : GridObject* elem,
142 : const MathVector<dim> vCornerCoords[],
143 : const MathVector<refDim> vLocIP[],
144 : const size_t nip,
145 : LocalVector* u,
146 : bool bDeriv,
147 : int s,
148 : std::vector<std::vector<MathVector<dim> > > vvvDeriv[],
149 : const MathMatrix<refDim, dim>* vJT = NULL
150 : ) const
151 : {
152 : // 1. Get the vectors themselves
153 0 : const MathVector<dim>* vVector = m_spVector->values (s);
154 :
155 : // 2a. Prepare the data of the element
156 0 : const ReferenceObjectID roid = elem->reference_object_id ();
157 0 : DimReferenceMapping<refDim, dim>& rMapping = ReferenceMappingProvider::get<refDim, dim> (roid);
158 0 : rMapping.update (vCornerCoords);
159 :
160 : // 2b. Check if we should compute the derivatives
161 0 : if (this->zero_derivative ())
162 : bDeriv = false;
163 : else
164 0 : this->set_zero (vvvDeriv, nip);
165 :
166 : // 3. At every integration point
167 0 : for (size_t i = 0; i < nip; i++)
168 : {
169 : // 3a. Get the Jacobian
170 : MathMatrix<dim, refDim> J;
171 0 : rMapping.jacobian (J, vLocIP[i]);
172 : // 3b. Project the vector to the subspace spanned by the columns of the Jacobian
173 0 : vValue[i] = vVector[i];
174 0 : OrthogProjectVec (vValue[i], J);
175 : // 3c. Project the derivatives in the same way
176 0 : if (! bDeriv) continue;
177 0 : for (size_t fct = 0; fct < m_spDVector->num_fct(); fct++)
178 : {
179 : const MathVector<dim>* vDVector = m_spDVector->deriv (s, i, fct);
180 : const size_t c_fct = this->input_common_fct (0, fct);
181 0 : for (size_t sh = 0; sh < this->num_sh (c_fct); sh++)
182 : {
183 0 : vvvDeriv[i][c_fct][sh] = vDVector [sh];
184 0 : OrthogProjectVec (vvvDeriv[i][c_fct][sh], J);
185 : }
186 : }
187 : }
188 0 : }
189 :
190 : private:
191 :
192 : /// data to project
193 : SmartPtr<CplUserData<MathVector<dim>, dim> > m_spVector;
194 : SmartPtr<DependentUserData<MathVector<dim>, dim> > m_spDVector;
195 : };
196 :
197 : } // end namespace ug
198 :
199 : #endif // __H__UG__LIB_DISC__SPATIAL_DISC__PROJECT_LINKER__
200 :
201 : /* End of File */
|
