Loading [MathJax]/extensions/tex2jax.js
ug4
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
inner_boundary_impl.h
Go to the documentation of this file.
1/*
2 * Copyright (c) 2011-2015: G-CSC, Goethe University Frankfurt
3 * Author: Markus Breit
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__NEUMANN_BOUNDARY__FV1__INNER_BOUNDARY_IMPL__
34#define __H__UG__LIB_DISC__SPACIAL_DISCRETIZATION__ELEM_DISC__NEUMANN_BOUNDARY__FV1__INNER_BOUNDARY_IMPL__
35
36#include "inner_boundary.h"
40
41
42namespace ug
43{
44
45template <typename TImpl, typename TDomain>
51
52
53template <typename TImpl, typename TDomain>
56{
57 m_fluxScale.set_data(scaleFct);
58 this->register_import(m_fluxScale);
59}
60
61#ifdef UG_FOR_LUA
62template <typename TImpl, typename TDomain>
64set_flux_scale(const char* luaScaleFctName)
65{
66 set_flux_scale(LuaUserDataFactory<number,dim>::create(luaScaleFctName));
67}
68#endif
69
70
71
72template <typename TImpl, typename TDomain>
74prepare_setting(const std::vector<LFEID>& vLfeID, bool bNonRegularGrid)
75{
76 // check that Lagrange 1st order
77 for(size_t i = 0; i < vLfeID.size(); ++i)
78 if(vLfeID[i].type() != LFEID::LAGRANGE || vLfeID[i].order() != 1)
79 UG_THROW("FV1InnerBoundary: 1st order lagrange expected.");
80
81 // update assemble functions
82 m_bNonRegularGrid = bNonRegularGrid;
83 register_all_fv1_funcs();
84}
85
86template <typename TImpl, typename TDomain>
88use_hanging() const
89{
90 return true;
91}
92
93
94template <typename TImpl, typename TDomain>
96{
97 m_bPrevSolRequired = b;
98}
99
100
101template <typename TImpl, typename TDomain>
102template <typename TAlgebra>
104(
105 number future_time,
106 number time,
107 VectorProxyBase* upb
108)
109{
110 if (!m_bPrevSolRequired)
111 return;
112
113 // the proxy provided here will not survive, but the vector it wraps will,
114 // so we re-wrap it in a new proxy for ourselves
116 tProxy* proxy = static_cast<tProxy*>(upb);
117 m_spOldSolutionProxy = make_sp(new tProxy(proxy->m_v));
118}
119
120
121template <typename TImpl, typename TDomain>
122template<typename TElem, typename TFVGeom>
124prep_elem_loop(const ReferenceObjectID roid, const int si)
125{
126 m_si = si;
127
128 // set local positions
129 if (!TFVGeom::usesHangingNodes)
130 {
131 static const int refDim = TElem::dim;
132 static const TFVGeom& geo = GeomProvider<TFVGeom>::get();
133 const MathVector<refDim>* vBFip = geo.bf_local_ips();
134 const size_t numBFip = geo.num_bf_local_ips();
135
136 m_fluxScale.template set_local_ips<refDim>(vBFip, numBFip, false);
137 }
138}
139
140
141template <typename TImpl, typename TDomain>
142template<typename TElem, typename TFVGeom>
146
147
148template <typename TImpl, typename TDomain>
149template<typename TElem, typename TFVGeom>
151prep_elem(const LocalVector& u, GridObject* elem, const ReferenceObjectID roid, const MathVector<dim> vCornerCoords[])
152{
153#ifdef UG_PARALLEL
154 DistributedGridManager& dgm = *this->approx_space()->domain()->grid()->distributed_grid_manager();
155 m_bCurrElemIsHSlave = dgm.get_status(elem) & ES_H_SLAVE;
156#endif
157
158 // on horizontal interfaces: only treat hmasters
159 if (m_bCurrElemIsHSlave) return;
160
161 // update Geometry for this element
162 static TFVGeom& geo = GeomProvider<TFVGeom>::get();
163 try {geo.update(elem, vCornerCoords, &(this->subset_handler()));}
164 UG_CATCH_THROW("FV1InnerBoundaryElemDisc::prep_elem: "
165 "Cannot update Finite Volume Geometry.");
166
167 // set local positions
168 if (TFVGeom::usesHangingNodes)
169 {
170 const int refDim = TElem::dim;
171 const MathVector<refDim>* vBFip = geo.bf_local_ips();
172 const size_t numBFip = geo.num_bf_local_ips();
173
174 m_fluxScale.template set_local_ips<refDim>(vBFip, numBFip);
175 }
176
177 // set global positions
178 const MathVector<dim>* vBFip = geo.bf_global_ips();
179 const size_t numBFip = geo.num_bf_global_ips();
180
181 m_fluxScale.set_global_ips(vBFip, numBFip);
182
183 // create a local vector with values of the previous solution
184 if (!m_bPrevSolRequired)
185 return;
186
187 // copy local vector of current solution to get indices and map access right,
188 // then copy correct values
189 m_locUOld = u;
190 const LocalIndices& ind = m_locUOld.get_indices();
191 for (size_t fct = 0; fct < m_locUOld.num_all_fct(); ++fct)
192 {
193 for(size_t dof = 0; dof < m_locUOld.num_all_dof(fct); ++dof)
194 {
195 const DoFIndex di = ind.multi_index(fct, dof);
196 m_locUOld.value(fct, dof) = m_spOldSolutionProxy->evaluate(di);
197 }
198 }
199}
200
201// assemble stiffness part of Jacobian
202template <typename TImpl, typename TDomain>
203template<typename TElem, typename TFVGeom>
205add_jac_A_elem(LocalMatrix& J, const LocalVector& u, GridObject* elem, const MathVector<dim> vCornerCoords[])
206{
207 // on horizontal interfaces: only treat hmasters
208 if (m_bCurrElemIsHSlave) return;
209
210 // get finite volume geometry
211 const static TFVGeom& fvgeom = GeomProvider<TFVGeom>::get();
212
213 FluxDerivCond fdc;
214 size_t nFct = u.num_fct();
215 m_uAtCorner.resize(nFct);
216 m_uOldAtCorner.resize(nFct);
217 for (size_t i = 0; i < fvgeom.num_bf(); ++i)
218 {
219 // get current BF
220 const typename TFVGeom::BF& bf = fvgeom.bf(i);
221
222 // get associated node and subset index
223 const int co = bf.node_id();
224
225 // get corner coordinates
226 const MathVector<dim>& cc = bf.global_corner(0);
227
228 // get solution at the corner of the bf
229 for (size_t fct = 0; fct < nFct; fct++)
230 m_uAtCorner[fct] = u(fct,co);
231
232 if (!m_bPrevSolRequired)
233 {
234 if (!static_cast<TImpl*>(this)->fluxDensityDerivFct(m_uAtCorner, elem, cc, m_si, fdc))
235 UG_THROW("FV1InnerBoundaryElemDisc::add_jac_A_elem:"
236 " Call to fluxDensityDerivFct did not succeed.");
237 }
238 else
239 {
240 for (size_t fct = 0; fct < nFct; fct++)
241 m_uOldAtCorner[fct] = m_locUOld(fct,co);
242
243 if (!static_cast<TImpl*>(this)->fluxDensityDerivFct(m_uAtCorner, m_uOldAtCorner, elem, cc, m_si, fdc))
244 UG_THROW("FV1InnerBoundaryElemDisc::add_jac_A_elem:"
245 " Call to fluxDensityDerivFct did not succeed.");
246 }
247
248 // scale with volume of BF and flux scale (if applicable)
249 number scale = bf.volume();
250 if (m_fluxScale.data_given())
251 scale *= m_fluxScale[i];
252
253 for (size_t j=0; j<fdc.fluxDeriv.size(); j++)
254 for (size_t k=0; k<fdc.fluxDeriv[j].size(); k++)
255 fdc.fluxDeriv[j][k].second *= scale;
256
257 // add to Jacobian
258 for (size_t j=0; j<fdc.fluxDeriv.size(); j++)
259 {
260 for (size_t k=0; k<fdc.fluxDeriv[j].size(); k++)
261 {
263 J(fdc.from[j], co, fdc.fluxDeriv[j][k].first, co) += fdc.fluxDeriv[j][k].second;
265 J(fdc.to[j], co, fdc.fluxDeriv[j][k].first, co) -= fdc.fluxDeriv[j][k].second;
266 }
267 }
268 }
269}
270
271template <typename TImpl, typename TDomain>
272template<typename TElem, typename TFVGeom>
276
277
278template <typename TImpl, typename TDomain>
279template<typename TElem, typename TFVGeom>
281add_def_A_elem(LocalVector& d, const LocalVector& u, GridObject* elem, const MathVector<dim> vCornerCoords[])
282{
283 // on horizontal interfaces: only treat hmasters
284 if (m_bCurrElemIsHSlave) return;
285
286 // get finite volume geometry
287 static TFVGeom& fvgeom = GeomProvider<TFVGeom>::get();
288
289 FluxCond fc;
290 size_t nFct = u.num_fct();
291 m_uAtCorner.resize(nFct);
292 m_uOldAtCorner.resize(nFct);
293
294 // loop Boundary Faces
295 for (size_t i = 0; i < fvgeom.num_bf(); ++i)
296 {
297 // get current BF
298 const typename TFVGeom::BF& bf = fvgeom.bf(i);
299
300 // get associated node and subset index
301 const int co = bf.node_id();
302
303 // get solution at the corner of the bf
304 for (size_t fct = 0; fct < nFct; fct++)
305 m_uAtCorner[fct] = u(fct,co);
306
307 // get corner coordinates
308 const MathVector<dim>& cc = bf.global_corner(0);
309
310 // get flux densities in that node
311 if (!m_bPrevSolRequired)
312 {
313 if (!static_cast<TImpl*>(this)->fluxDensityFct(m_uAtCorner, elem, cc, m_si, fc))
314 {
315 UG_THROW("FV1InnerBoundaryElemDisc::add_def_A_elem:"
316 " Call to fluxDensityFct did not succeed.");
317 }
318 }
319 else
320 {
321 // get solution at the corner of the bf
322 for (size_t fct = 0; fct < nFct; fct++)
323 m_uOldAtCorner[fct] = m_locUOld(fct,co);
324
325 if (!static_cast<TImpl*>(this)->fluxDensityFct(m_uAtCorner, m_uOldAtCorner, elem, cc, m_si, fc))
326 UG_THROW("FV1InnerBoundaryElemDisc::add_def_A_elem:"
327 " Call to fluxDensityFct did not succeed.");
328 }
329
330 // scale with volume of BF and flux scale (if applicable)
331 number scale = bf.volume();
332 if (m_fluxScale.data_given())
333 scale *= m_fluxScale[i];
334
335 for (size_t j=0; j<fc.flux.size(); j++)
336 fc.flux[j] *= scale;
337
338 // add to defect
339 for (size_t j=0; j<fc.flux.size(); j++)
340 {
341 if (fc.from[j] != InnerBoundaryConstants::_IGNORE_) d(fc.from[j], co) += fc.flux[j];
342 if (fc.to[j] != InnerBoundaryConstants::_IGNORE_) d(fc.to[j], co) -= fc.flux[j];
343 }
344 }
345}
346
347template <typename TImpl, typename TDomain>
348template<typename TElem, typename TFVGeom>
352
353template <typename TImpl, typename TDomain>
354template<typename TElem, typename TFVGeom>
358
359
360// ////////////////////////////////////////////////////////////////////////////////////////////////
361// error estimation (begin) //
362
363// prepares the loop over all elements of one type for the computation of the error estimator
364template <typename TImpl, typename TDomain>
365template <typename TElem, typename TFVGeom>
367prep_err_est_elem_loop(const ReferenceObjectID roid, const int si)
368{
369 m_si = si;
370
371 // get the error estimator data object and check that it is of the right type
372 // we check this at this point in order to be able to dispense with this check later on
373 // (i.e. in prep_err_est_elem and compute_err_est_A_elem())
374 if (this->m_spErrEstData.get() == NULL)
375 {
376 UG_THROW("No ErrEstData object has been given to this ElemDisc!");
377 }
378
379 err_est_type* err_est_data = dynamic_cast<err_est_type*>(this->m_spErrEstData.get());
380
381 if (!err_est_data)
382 {
383 UG_THROW("Dynamic cast to MultipleSideAndElemErrEstData failed."
384 << std::endl << "Make sure you handed the correct type of ErrEstData to this discretization.");
385 }
386
387 if (!err_est_data->equal_side_order())
388 {
389 UG_THROW("The underlying error estimator data objects of this discretization's "
390 "error estimator do not all have the same integration orders. This case "
391 "is not supported by the implementation. If you need it, implement!");
392 }
393
394 if (err_est_data->num() < 1)
395 {
396 UG_THROW("No underlying error estimator data objects present. No IPs can be determined.");
397 }
398
399 // set local positions
400 if (!TFVGeom::usesHangingNodes)
401 {
402 static const int refDim = TElem::dim;
403
404 // get local IPs
405 size_t numSideIPs;
406 const MathVector<refDim>* sideIPs;
407 try
408 {
409 numSideIPs = err_est_data->get(0)->num_side_ips(roid);
410 sideIPs = err_est_data->get(0)->template side_local_ips<refDim>(roid);
411 }
412 UG_CATCH_THROW("Integration points for error estimator cannot be set.");
413
414 // set local IPs in import
415 m_fluxScale.template set_local_ips<refDim>(sideIPs, numSideIPs, false);
416
417 // store values of shape functions in local IPs
420
421 m_shapeValues.resize(numSideIPs, trialSpace.num_sh());
422 for (size_t ip = 0; ip < numSideIPs; ip++)
423 trialSpace.shapes(m_shapeValues.shapesAtSideIP(ip), sideIPs[ip]);
424 }
425}
426
427template <typename TImpl, typename TDomain>
428template<typename TElem, typename TFVGeom>
430prep_err_est_elem(const LocalVector& u, GridObject* elem, const MathVector<dim> vCornerCoords[])
431{
432#ifdef UG_PARALLEL
433 DistributedGridManager& dgm = *this->approx_space()->domain()->grid()->distributed_grid_manager();
434 m_bCurrElemIsHSlave = dgm.get_status(elem) & ES_H_SLAVE;
435#endif
436
437 // on horizontal interfaces: only treat hmasters
438 if (m_bCurrElemIsHSlave) return;
439
440 // get error estimator
441 err_est_type* err_est_data = dynamic_cast<err_est_type*>(this->m_spErrEstData.get());
442
443 // roid
445
446 // set local positions
447 if (TFVGeom::usesHangingNodes)
448 {
449 static const int refDim = TElem::dim;
450
451 size_t numSideIPs;
452 const MathVector<refDim>* sideIPs;
453 try
454 {
455 numSideIPs = err_est_data->get(0)->num_side_ips(roid);
456 sideIPs = err_est_data->get(0)->template side_local_ips<refDim>(roid);
457 }
458 UG_CATCH_THROW("Integration points for error estimator cannot be set.");
459
460 // set local IPs in import
461 m_fluxScale.template set_local_ips<refDim>(sideIPs, numSideIPs);
462
463 // store values of shape functions in local IPs
466
467 m_shapeValues.resize(numSideIPs, trialSpace.num_sh());
468 for (size_t ip = 0; ip < numSideIPs; ip++)
469 trialSpace.shapes(m_shapeValues.shapesAtSideIP(ip), sideIPs[ip]);
470 }
471
472 // set global positions
473 size_t numSideIPs;
474 MathVector<dim>* sideIPs;
475
476 try
477 {
478 numSideIPs = err_est_data->get(0)->num_side_ips(roid);
479 sideIPs = err_est_data->get(0)->side_global_ips(elem, vCornerCoords);
480 }
481 UG_CATCH_THROW("Global integration points for error estimator cannot be set.");
482
483 // set local IPs in import
484 m_fluxScale.set_global_ips(&sideIPs[0], numSideIPs);
485}
486
487// computes the error estimator contribution for one element
488template <typename TImpl, typename TDomain>
489template <typename TElem, typename TFVGeom>
491compute_err_est_A_elem(const LocalVector& u, GridObject* elem, const MathVector<dim> vCornerCoords[], const number& scale)
492{
493 // on horizontal interfaces: only treat hmasters
494 if (m_bCurrElemIsHSlave) return;
495
496 // get error estimator
497 err_est_type* err_est_data = dynamic_cast<err_est_type*>(this->m_spErrEstData.get());
498
499 // cast this elem to side_type of error estimator
501 dynamic_cast<typename SideAndElemErrEstData<domain_type>::side_type*>(elem);
502 if (!side)
503 {
504 UG_THROW("Error in DependentNeumannBoundaryFV1::compute_err_est_A_elem():\n"
505 "Element that error assembling routine is called for has the wrong type.");
506 }
507
508 // global IPs
509 ReferenceObjectID roid = side->reference_object_id();
510 size_t numSideIPs = err_est_data->get(0)->num_side_ips(roid);
511 MathVector<dim>* globIPs = err_est_data->get(0)->side_global_ips(side, vCornerCoords);
512
513 // loop IPs
514 try
515 {
516 for (size_t sip = 0; sip < numSideIPs; sip++)
517 {
518 // get values of u at ip (interpolate)
519 size_t nFct = u.num_fct();
520 std::vector<LocalVector::value_type> uAtIP = std::vector<LocalVector::value_type>(nFct);
521
522 for (size_t fct = 0; fct < nFct; fct++)
523 {
524 uAtIP[fct] = 0.0;
525 for (size_t sh = 0; sh < m_shapeValues.num_sh(); sh++)
526 uAtIP[fct] += u(fct,sh) * m_shapeValues.shapeAtSideIP(sh,sip);
527 }
528
529 // ip coordinates
530 const MathVector<dim>& ipCoords = globIPs[sip];
531
532 FluxCond fc;
533 if (!static_cast<TImpl*>(this)->fluxDensityFct(uAtIP, elem, ipCoords, m_si, fc))
534 {
535 UG_THROW("FV1InnerBoundaryElemDisc::compute_err_est_A_elem:"
536 " Call to fluxDensityFct did not succeed.");
537 }
538
539 // scale fluxes if applicable
540 if (m_fluxScale.data_given())
541 for (size_t j = 0; j < fc.flux.size(); ++j)
542 fc.flux[j] *= m_fluxScale[sip];
543
544
545 // subtract from estimator values
546 // sign must be opposite of that in add_def_A_elem
547 // as the difference between this and the actual flux of the unknown is calculated
548 for (size_t j=0; j<fc.flux.size(); j++)
549 {
551 (*err_est_data->get(this->m_fctGrp[fc.from[j]])) (side,sip) -= scale * fc.flux[j];
553 (*err_est_data->get(this->m_fctGrp[fc.to[j]])) (side,sip) += scale * fc.flux[j];
554 }
555 }
556 }
557 UG_CATCH_THROW("Values for the error estimator could not be assembled at every IP." << std::endl
558 << "Maybe wrong type of ErrEstData object? This implementation needs: MultipleSideAndElemErrEstData.");
559}
560
561// postprocesses the loop over all elements of one type in the computation of the error estimator
562template <typename TImpl, typename TDomain>
563template <typename TElem, typename TFVGeom>
569
570// error estimation (end) //
571// ////////////////////////////////////////////////////////////////////////////////////////////////
572
573
575// register assemble functions
577
578// register for 1D
579template <typename TImpl, typename TDomain>
582{
583// register prep_timestep function for all known algebra types
585
586// register assembling functions for all element types
587 typedef typename domain_traits<dim>::ManifoldElemList ElemList;
588 boost::mpl::for_each<ElemList>(RegisterAssemblingFuncs(this));
589}
590
591template <typename TImpl, typename TDomain>
592template<typename TElem, typename TFVGeom>
595{
597 typedef this_type T;
598
599 this->clear_add_fct(id);
600 this->set_prep_elem_loop_fct( id, &T::template prep_elem_loop<TElem, TFVGeom>);
601 this->set_prep_elem_fct( id, &T::template prep_elem<TElem, TFVGeom>);
602 this->set_fsh_elem_loop_fct( id, &T::template fsh_elem_loop<TElem, TFVGeom>);
603 this->set_add_jac_A_elem_fct( id, &T::template add_jac_A_elem<TElem, TFVGeom>);
604 this->set_add_jac_M_elem_fct( id, &T::template add_jac_M_elem<TElem, TFVGeom>);
605 this->set_add_def_A_elem_fct( id, &T::template add_def_A_elem<TElem, TFVGeom>);
606 this->set_add_def_M_elem_fct( id, &T::template add_def_M_elem<TElem, TFVGeom>);
607 this->set_add_rhs_elem_fct( id, &T::template add_rhs_elem<TElem, TFVGeom>);
608
609 // error estimator parts
610 this->set_prep_err_est_elem_loop( id, &T::template prep_err_est_elem_loop<TElem, TFVGeom>);
611 this->set_prep_err_est_elem( id, &T::template prep_err_est_elem<TElem, TFVGeom>);
612 this->set_compute_err_est_A_elem( id, &T::template compute_err_est_A_elem<TElem, TFVGeom>);
613 this->set_fsh_err_est_elem_loop( id, &T::template fsh_err_est_elem_loop<TElem, TFVGeom>);
614}
615
616} // namespace ug
617
618
619#endif /*__H__UG__LIB_DISC__SPACIAL_DISCRETIZATION__ELEM_DISC__NEUMANN_BOUNDARY__FV1__INNER_BOUNDARY_IMPL__*/
Definition smart_pointer.h:108
void shapes(shape_type *vShape, const MathVector< dim > &x) const
returns all shape functions evaluated at a point
Definition local_shape_function_set.h:220
constant scalar user data
Definition const_user_data.h:153
Type based UserData.
Definition user_data.h:501
manages the layouts and interfaces which are associated with a distributed grid.
Definition distributed_grid.h:88
byte get_status(GridObject *go) const
returns the status of the given object.
Definition distributed_grid.cpp:379
MultipleSideAndElemErrEstData< domain_type > err_est_type
error estimator type
Definition inner_boundary.h:122
void prep_elem_loop(const ReferenceObjectID roid, const int si)
prepares the loop over all elements
Definition inner_boundary_impl.h:124
void prep_err_est_elem_loop(const ReferenceObjectID roid, const int si)
prepares the loop over all elements of one type for the computation of the error estimator
Definition inner_boundary_impl.h:367
void register_all_fv1_funcs()
Definition inner_boundary_impl.h:581
virtual void prepare_setting(const std::vector< LFEID > &vLfeID, bool bNonRegularGrid) override
type of trial space for each function used
Definition inner_boundary_impl.h:74
void register_fv1_func()
Definition inner_boundary_impl.h:594
void prep_timestep(number future_time, number time, VectorProxyBase *upb)
prepare the time step
Definition inner_boundary_impl.h:104
void add_def_M_elem(LocalVector &d, const LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[])
assembles the mass part of the local defect
Definition inner_boundary_impl.h:350
FV1InnerBoundaryElemDisc< TImpl, TDomain > this_type
own type
Definition inner_boundary.h:112
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 inner_boundary_impl.h:281
void add_jac_M_elem(LocalMatrix &J, const LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[])
assembles the local mass matrix using a finite volume scheme
Definition inner_boundary_impl.h:274
void previous_solution_required(bool b)
Definition inner_boundary_impl.h:95
virtual bool use_hanging() const override
returns if hanging nodes are used
Definition inner_boundary_impl.h:88
void add_rhs_elem(LocalVector &rhs, GridObject *elem, const MathVector< dim > vCornerCoords[])
assembles the local right hand side
Definition inner_boundary_impl.h:356
void fsh_elem_loop()
finishes the loop over all elements
Definition inner_boundary_impl.h:144
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 inner_boundary_impl.h:205
void set_flux_scale(number scale)
Setting the flux function.
Definition inner_boundary_impl.h:47
void compute_err_est_A_elem(const LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[], const number &scale)
computes the error estimator contribution for one element
Definition inner_boundary_impl.h:491
void prep_elem(const LocalVector &u, GridObject *elem, const ReferenceObjectID roid, const MathVector< dim > vCornerCoords[])
prepares the element for assembling
Definition inner_boundary_impl.h:151
void prep_err_est_elem(const LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[])
prepares the element for assembling the error estimator
Definition inner_boundary_impl.h:430
void fsh_err_est_elem_loop()
postprocesses the loop over all elements of one type in the computation of the error estimator
Definition inner_boundary_impl.h:565
static TGeom & get()
returns a singleton based on the identifier
Definition geom_provider.h:141
The base class for all geometric objects, such as vertices, edges, faces, volumes,...
Definition grid_base_objects.h:157
virtual ReferenceObjectID reference_object_id() const =0
@ LAGRANGE
Definition local_finite_element_id.h:104
Lagrange shape functions of first order.
Definition lagrangep1.h:55
size_t num_sh() const
Definition lagrangep1.h:87
Definition local_algebra.h:50
const DoFIndex & multi_index(size_t fct, size_t dof) const
global algebra multi-index for (fct, dof)
Definition local_algebra.h:132
Definition local_algebra.h:422
Definition local_algebra.h:198
size_t num_fct() const
returns the number of currently accessible functions
Definition local_algebra.h:328
Factory providing LuaUserData.
Definition lua_user_data.h:180
a mathematical Vector with N entries.
Definition math_vector.h:97
Definition multi_index.h:50
Provider, holding a single instance of an object.
Definition provider.h:48
domain_traits< dim >::side_type side_type
type of the sides (face, edge) and the elems (volume, face)
Definition err_est_data.h:230
size_t num_side_ips(const side_type *pSide)
get number of side IPs of a specific side
Definition err_est_data_impl.h:419
Definition grid_base_object_traits.h:68
@ ES_H_SLAVE
Definition distributed_grid.h:59
#define UG_CATCH_THROW(msg)
Definition error.h:64
#define UG_THROW(msg)
Definition error.h:57
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
static const size_t _IGNORE_
index value for which a flux is ignored
Definition inner_boundary.h:69
struct that holds information about the flux densities and from where to where the flux occurs
Definition inner_boundary.h:75
std::vector< size_t > from
Definition inner_boundary.h:78
std::vector< number > flux
Definition inner_boundary.h:77
std::vector< size_t > to
Definition inner_boundary.h:79
Definition inner_boundary.h:85
std::vector< size_t > from
Definition inner_boundary.h:90
std::vector< std::vector< std::pair< size_t, number > > > fluxDeriv
Definition inner_boundary.h:89
std::vector< size_t > to
Definition inner_boundary.h:91
Proxy struct for generic passing of any vector type.
Definition elem_disc_interface.h:70
Definition elem_disc_interface.h:77
Definition domain_traits.h:53