ug4
lagrange_dirac_source_impl.h
Go to the documentation of this file.
1 /*
2  * Copyright (c) 2019-2020: G-CSC, Goethe University Frankfurt
3  * Author: Arne Naegel
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 #ifndef __H__UG__LIB_DISC__SPACIAL_DISCRETIZATION__ELEM_DISC__DIRAC_SOURCE__LAGRANGE_DIRAC_SOURCE_IMPL__
34 #define __H__UG__LIB_DISC__SPACIAL_DISCRETIZATION__ELEM_DISC__DIRAC_SOURCE__LAGRANGE_DIRAC_SOURCE_IMPL__
35 
36 #include "lagrange_dirac_source.h"
39 
40 
41 
42 namespace ug
43 {
44 
45 template<typename TDomain>
47 add_source(number scale, MathVector<dim> &coord)
48 {
49  add_source(make_sp(new ConstUserNumber<dim>(scale)), coord);
50 }
51 
52 
53 template<typename TDomain>
56 {
57  m_srcData = srcData;
58  m_srcCoord = srcCoord;
59 }
60 
61 #ifdef UG_FOR_LUA
62 template<typename TDomain>
64 add_source(const char* luaFctName, MathVector<dim> &coord)
65 {
66  add_source(LuaUserDataFactory<number,dim>::create(luaFctName), coord);
67 }
68 #endif
69 
70 
71 
72 template<typename TDomain>
75 { m_snkTransportData = snkTransportData; }
76 
77 template<typename TDomain>
80 { add_transport_sink(make_sp(new ConstUserNumber<dim>(scale))); }
81 
82 
83 #ifdef UG_FOR_LUA
84 template<typename TDomain>
86 add_transport_sink(const char* luaFctName)
87 { add_transport_sink(LuaUserDataFactory<number,dim>::create(luaFctName)); }
88 #endif
89 
90 
91 
92 template<typename TDomain>
94 prepare_setting(const std::vector<LFEID>& vLfeID, bool bNonRegularGrid)
95 {
96  // check that Lagrange 1st order
97  for(size_t i = 0; i < vLfeID.size(); ++i)
98  if(vLfeID[i].type() != LFEID::LAGRANGE )
99  UG_THROW("DiracSourceDisc: Lagrange polynomials expected.");
100 
101  // update assemble functions
102  m_bNonRegularGrid = bNonRegularGrid;
103  register_all_funcs();
104 }
105 
106 template<typename TDomain>
108 use_hanging() const
109 {
110  return true;
111 }
112 
113 
114 template<typename TDomain>
115 template<typename TElem, typename TFVGeom>
117 add_jac_A_elem(LocalMatrix& J, const LocalVector& u, GridObject* elem, const MathVector<dim> vCornerCoords[])
118 {
119  if (m_snkTransportData.invalid()) return;
120  UG_LOG("'DiracSourceDisc::add_jac_A_elem' called for " << vCornerCoords[0] << std::endl);
121  // Request source strength form user data.
122  const int co = 0;
123  double snkSVal;
124  (*m_snkTransportData)(snkSVal, m_srcCoord, this->time(), -1);
125  J(_C_, co, _C_, co) -= snkSVal;
126 }
127 
128 
130 template<typename TDomain>
131 template<typename TElem, typename TFVGeom>
133 add_def_A_elem(LocalVector& d, const LocalVector& u, GridObject* elem, const MathVector<dim> vCornerCoords[])
134 {
135  if (m_snkTransportData.invalid()) return;
136  // Request source strength form user data.
137  const int co = 0;
138  double snkSVal;
139  UG_LOG("'DiracSourceDisc::add_def_A_elem' called for " << vCornerCoords[0] << std::endl);
140  (*m_snkTransportData)(snkSVal, m_srcCoord, this->time (), -1);
141  d(_C_, co) -= u(_C_, co) * snkSVal;
142 }
143 
144 template<typename TDomain>
145 template<typename TElem, typename TFVGeom>
147 add_rhs_elem(LocalVector& rhs, GridObject* elem, const MathVector<dim> vCornerCoords[])
148 {
149  UG_LOG("'DiracSourceDisc::add_rhs_elem' called for " << vCornerCoords[0]<< std::endl);
150  UG_LOG("elem:"<< elem << std::endl);
151 
152  if (m_srcData.invalid()) return;
153 
154  UG_ASSERT(vCornerCoords[0] == m_srcCoord, "Source must be located at element corner!");
155 
156  // Request source strength form user data.
157  const int co = 0;
158  double srcVal;
159  (*m_srcData)(srcVal, m_srcCoord, this->time (), -1);
160  rhs(0, co) += srcVal;
161 
162 }
163 
164 
165 
166 
167 
169 // register assemble functions
171 
172 // register for 1D
173 template<typename TDomain>
176 {
177 
178 
179 // switch assemble functions
180  // boost::mpl::for_each<ElemList>(RegisterFV1(this));
181  register_func<RegularVertex, FV1Geometry<RegularVertex, dim> >();
182 }
183 
184 template<typename TDomain>
185 template<typename TElem, typename TFVGeom>
188 {
190  typedef this_type T;
191 
192  this->clear_add_fct(id);
193  this->set_prep_elem_loop_fct( id, &T::template prep_elem_loop<TElem, TFVGeom>);
194  this->set_prep_elem_fct( id, &T::template prep_elem<TElem, TFVGeom>);
195  this->set_fsh_elem_loop_fct( id, &T::template fsh_elem_loop<TElem, TFVGeom>);
196  this->set_add_jac_A_elem_fct( id, &T::template add_jac_A_elem<TElem, TFVGeom>);
197  this->set_add_jac_M_elem_fct( id, &T::template add_jac_M_elem<TElem, TFVGeom>);
198  this->set_add_def_A_elem_fct( id, &T::template add_def_A_elem<TElem, TFVGeom>);
199  this->set_add_def_M_elem_fct( id, &T::template add_def_M_elem<TElem, TFVGeom>);
200  this->set_add_rhs_elem_fct( id, &T::template add_rhs_elem<TElem, TFVGeom>);
201 
202  // error estimator parts
203  this->set_prep_err_est_elem_loop( id, &T::template prep_err_est_elem_loop<TElem, TFVGeom>);
204  this->set_prep_err_est_elem( id, &T::template prep_err_est_elem<TElem, TFVGeom>);
205  this->set_compute_err_est_A_elem( id, &T::template compute_err_est_A_elem<TElem, TFVGeom>);
206  this->set_fsh_err_est_elem_loop( id, &T::template fsh_err_est_elem_loop<TElem, TFVGeom>);
207 }
208 
210 // explicit template instantiations
212 
213 #ifdef UG_DIM_1
214 template class DiracSourceDisc<Domain1d>;
215 #endif
216 #ifdef UG_DIM_2
217 template class DiracSourceDisc<Domain2d>;
218 #endif
219 #ifdef UG_DIM_3
220 template class DiracSourceDisc<Domain3d>;
221 #endif
222 
223 
224 } // namespace ug
225 
226 
227 #endif /*__H__UG__LIB_DISC__SPACIAL_DISCRETIZATION__ELEM_DISC__NEUMANN_BOUNDARY__FV1__INNER_BOUNDARY_IMPL__*/
Definition: smart_pointer.h:108
constant scalar user data
Definition: const_user_data.h:153
Definition: lagrange_dirac_source.h:73
virtual void prepare_setting(const std::vector< LFEID > &vLfeID, bool bNonRegularGrid)
type of trial space for each function used
Definition: lagrange_dirac_source_impl.h:94
void add_jac_A_elem(LocalMatrix &J, const LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[])
assembles the local stiffness matrix using a finite volume scheme
Definition: lagrange_dirac_source_impl.h:117
void add_source(number scale, MathVector< dim > &srcCoord)
Setting a scaling factor for the flux.
Definition: lagrange_dirac_source_impl.h:47
virtual bool use_hanging() const
returns if hanging nodes are used
Definition: lagrange_dirac_source_impl.h:108
void register_func()
Definition: lagrange_dirac_source_impl.h:187
void add_transport_sink(SmartPtr< UserData< number, dim > > snkData)
Setting a scaling factor for the flux.
Definition: lagrange_dirac_source_impl.h:74
void register_all_funcs()
Definition: lagrange_dirac_source_impl.h:175
void add_def_A_elem(LocalVector &d, const LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[])
assembles the stiffness part of the local defect
Definition: lagrange_dirac_source_impl.h:133
void add_rhs_elem(LocalVector &rhs, GridObject *elem, const MathVector< dim > vCornerCoords[])
assembles the local right hand side
Definition: lagrange_dirac_source_impl.h:147
The base class for all geometric objects, such as vertices, edges, faces, volumes,...
Definition: grid_base_objects.h:157
@ LAGRANGE
Definition: local_finite_element_id.h:104
Definition: local_algebra.h:422
Definition: local_algebra.h:198
Factory providing LuaUserData.
Definition: lua_user_data.h:180
Type based UserData.
Definition: user_data.h:143
Definition: grid_base_object_traits.h:68
#define UG_ASSERT(expr, msg)
Definition: assert.h:70
#define UG_THROW(msg)
Definition: error.h:57
#define UG_LOG(msg)
Definition: log.h:367
double number
Definition: types.h:124
the ug namespace
ReferenceObjectID
these ids are used to identify the shape of a geometric object.
Definition: grid_base_objects.h:74
SmartPtr< T, FreePolicy > make_sp(T *inst)
returns a SmartPtr for the passed raw pointer
Definition: smart_pointer.h:836