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__
45 template<
typename TDomain>
53 template<
typename TDomain>
58 m_srcCoord = srcCoord;
62 template<
typename TDomain>
72 template<
typename TDomain>
75 { m_snkTransportData = snkTransportData; }
77 template<
typename TDomain>
84 template<
typename TDomain>
92 template<
typename TDomain>
97 for(
size_t i = 0; i < vLfeID.size(); ++i)
99 UG_THROW(
"DiracSourceDisc: Lagrange polynomials expected.");
102 m_bNonRegularGrid = bNonRegularGrid;
103 register_all_funcs();
106 template<
typename TDomain>
114 template<
typename TDomain>
115 template<
typename TElem,
typename TFVGeom>
119 if (m_snkTransportData.invalid())
return;
120 UG_LOG(
"'DiracSourceDisc::add_jac_A_elem' called for " << vCornerCoords[0] << std::endl);
124 (*m_snkTransportData)(snkSVal, m_srcCoord, this->time(), -1);
125 J(_C_, co, _C_, co) -= snkSVal;
130 template<
typename TDomain>
131 template<
typename TElem,
typename TFVGeom>
135 if (m_snkTransportData.invalid())
return;
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;
144 template<
typename TDomain>
145 template<
typename TElem,
typename TFVGeom>
149 UG_LOG(
"'DiracSourceDisc::add_rhs_elem' called for " << vCornerCoords[0]<< std::endl);
150 UG_LOG(
"elem:"<< elem << std::endl);
152 if (m_srcData.invalid())
return;
154 UG_ASSERT(vCornerCoords[0] == m_srcCoord,
"Source must be located at element corner!");
159 (*m_srcData)(srcVal, m_srcCoord, this->time (), -1);
160 rhs(0, co) += srcVal;
173 template<
typename TDomain>
181 register_func<RegularVertex, FV1Geometry<RegularVertex, dim> >();
184 template<
typename TDomain>
185 template<
typename TElem,
typename TFVGeom>
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>);
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>);
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
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