OpenFCST: The open-source Fuel Cell Simulation Toolbox
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
mesh_worker_info.h
Go to the documentation of this file.
1 // ------------------------------------------------------------------
2 //
3 // $Id: mesh_worker_info.h 1348 2013-08-16 00:45:07Z secanell $
4 //
5 // Copyright (C) 2006, 2007, 2008, 2009 by Guido Kanschat.
6 // Copyright (C) 2012 by Valentin N. Zingan, University of Alberta.
7 //
8 // This file is subject to QPL and may not be distributed
9 // without copyright and license information. Please refer
10 // to the file deal.II/doc/license.html for the text and
11 // further information on this license.
12 //
13 // ------------------------------------------------------------------
14 
15 #ifndef _DEALII_APPFRAME_MESH_WORKER_INFO_H_
16 #define _DEALII_APPFRAME_MESH_WORKER_INFO_H_
17 
18 #include <base/geometry_info.h>
19 #include <base/quadrature_lib.h>
20 #include <base/mg_level_object.h>
21 
22 #include <lac/block_indices.h>
23 
24 #include <fe/fe.h>
25 #include <fe/fe_tools.h>
26 
27 #include <multigrid/mg_tools.h>
28 
29 #include <appframe/matrix_base.h>
30 
31 #include <boost/shared_ptr.hpp>
32 
33 namespace dealii
34 {
35  template<int,int> class DoFHandler;
36  template<int,int> class MGDoFHandler;
37 }
38 
39 using namespace dealii;
40 using namespace AppFrame;
41 
42 // ################################################################################################
43 
58 // ##### CATEGORY-1 ################################################# INFO_OBJECTS ################
59 //
60 // ##### MeshWorker::BlockInfo #####
61 // ##### MeshWorker::VectorSelector #####
62 // ##### MeshWorker::InfoObjects::DoFInfo #####
63 // ##### MeshWorker::InfoObjects::IntegrationInfo #####
64 // ##### MeshWorker::InfoObjects::IntegrationInfoBox #####
65 
77 namespace MeshWorker
78 {
79 
86 struct BlockInfo : public Subscriptor
87 {
91  BlockIndices global;
92 
96  BlockIndices local;
97 
101  std::vector<BlockIndices> levels;
102 
108  std::vector<unsigned int> base_element;
109 
118  std::vector<unsigned int> local_renumbering;
119 
123  template<int dim, int spacedim>
124  void initialize(const DoFHandler<dim, spacedim>& dof_handler)
125  {
126  const FiniteElement<dim, spacedim>& fe = dof_handler.get_fe();
127  std::vector<unsigned int> sizes(fe.n_blocks());
128  DoFTools::count_dofs_per_block(dof_handler, sizes);
129  global.reinit(sizes);
130  }
131 
135  template<int dim, int spacedim>
136  void initialize_local(const DoFHandler<dim, spacedim>& dof_handler)
137  {
138  const FiniteElement<dim, spacedim>& fe = dof_handler.get_fe();
139  std::vector<unsigned int> sizes(fe.n_blocks());
140 
141  base_element.resize(fe.n_blocks());
142 
143  for(unsigned int i = 0; i < base_element.size(); ++i)
144  base_element[i] = fe.block_to_base_index(i).first;
145 
146  local_renumbering.resize(fe.n_dofs_per_cell());
147  FETools::compute_block_renumbering(fe,
148  local_renumbering,
149  sizes,
150  false);
151  local.reinit(sizes);
152  }
153 
157  template<int dim, int spacedim>
158  void initialize(const MGDoFHandler<dim, spacedim>& mg_dof_handler)
159  {
160  initialize((DoFHandler<dim, spacedim>&) mg_dof_handler);
161 
162  std::vector< std::vector<unsigned int> > sizes(mg_dof_handler.get_tria().n_levels(),
163  std::vector<unsigned int>(mg_dof_handler.get_fe().n_blocks())
164  );
165 
166  MGTools::count_dofs_per_block(mg_dof_handler, sizes);
167 
168  levels.resize(sizes.size());
169  for(unsigned int i = 0; i < sizes.size(); ++i)
170  levels[i].reinit(sizes[i]);
171  }
172 };
173 
182 class VectorSelector
183 {
184 public:
185 
190  VectorSelector() { }
191 
198  VectorSelector(const VectorSelector& other);
199 
209  void add_vector(const std::string& name,
210  bool values = true,
211  bool gradients = false,
212  bool hessians = false);
213 
222  void cache(const FEVectors& src);
223 
227  bool empty() const
228  {
229  return (selection.size() == 0);
230  }
231 
235  bool has_values() const
236  {
237  return selected_values;
238  }
239 
243  bool has_gradients() const
244  {
245  return selected_gradients;
246  }
247 
251  bool has_hessians() const
252  {
253  return selected_hessians;
254  }
255 
259  unsigned int n_values() const
260  {
261  return cache_values.size();
262  }
263 
267  unsigned int n_gradients() const
268  {
269  return cache_gradients.size();
270  }
271 
275  unsigned int n_hessians() const
276  {
277  return cache_hessians.size();
278  }
279 
283  unsigned int value_index(unsigned int i) const
284  {
285  AssertIndexRange(i, cache_values.size());
286  return cache_values[i];
287  }
288 
292  unsigned int gradient_index(unsigned int i) const
293  {
294  AssertIndexRange(i, cache_gradients.size());
295  return cache_gradients[i];
296  }
297 
301  unsigned int hessian_index(unsigned int i) const
302  {
303  AssertIndexRange(i, cache_hessians.size());
304  return cache_hessians[i];
305  }
306 
312  template<typename STREAM>
313  void print(STREAM& s,
314  const FEVectors& v) const;
315 
316 private:
317 
328  struct ListEntry
329  {
331  std::string name;
332 
334  bool values;
335 
337  bool gradients;
338 
340  bool hessians;
341  };
342 
346  std::vector<ListEntry> selection;
347 
351  bool selected_values;
352 
356  bool selected_gradients;
357 
361  bool selected_hessians;
362 
366  std::vector<unsigned int> cache_values;
367 
371  std::vector<unsigned int> cache_gradients;
372 
376  std::vector<unsigned int> cache_hessians;
377 };
378 
379 // --- implementation ---
380 
381 inline
382 VectorSelector::VectorSelector(const VectorSelector& other)
383 :
384 selection(other.selection),
385 selected_values(other.selected_values),
386 selected_gradients(other.selected_gradients),
387 selected_hessians(other.selected_hessians)
388 { }
389 
390 inline
391 void
392 VectorSelector::add_vector(const std::string& name,
393  bool values,
394  bool gradients,
395  bool hessians)
396 {
397  ListEntry e;
398  e.name = name;
399  e.values = values;
400  e.gradients = gradients;
401  e.hessians = hessians;
402 
403  selected_values |= values;
404  selected_gradients |= gradients;
405  selected_hessians |= hessians;
406 
407  selection.push_back(e);
408 }
409 
410 inline
411 void
412 VectorSelector::cache(const FEVectors& src)
413 {
414  cache_values.resize(0);
415  cache_gradients.resize(0);
416  cache_hessians.resize(0);
417 
418  for(std::vector<ListEntry>::const_iterator iter = selection.begin();
419  iter != selection.end();
420  ++iter)
421  {
422  const unsigned int k = src.find_vector(iter->name);
423 
424  Assert( k != deal_II_numbers::invalid_unsigned_int,
425  ApplicationData::ExcNotFound("vector", iter->name) );
426 
427  if(iter->values) cache_values.push_back(k);
428  if(iter->gradients) cache_gradients.push_back(k);
429  if(iter->hessians) cache_hessians.push_back(k);
430  }
431 }
432 
433 template<typename STREAM>
434 inline
435 void
436 VectorSelector::print(STREAM& s,
437  const FEVectors& v) const
438 {
439  s << "values: ";
440  for(unsigned int i = 0; i < cache_values.size(); ++i)
441  s << " '" << v.vector_name(cache_values[i]) << '\'';
442  s << std::endl << "gradients: ";
443  for(unsigned int i = 0; i < cache_gradients.size(); ++i)
444  s << " '" << v.vector_name(cache_gradients[i]) << '\'';
445  s << std::endl << "hessians: ";
446  for(unsigned int i = 0; i < cache_hessians.size(); ++i)
447  s << " '" << v.vector_name(cache_hessians[i]) << '\'';
448  s << std::endl;
449 }
450 
464 inline
465 unsigned int
466 selector_size(const VectorSelector& s,
467  const std::vector<double>&)
468 {
469  return s.n_values();
470 }
471 
472 template<int dim>
473 inline
474 unsigned int
475 selector_size(const VectorSelector& s,
476  const std::vector< Tensor<1, dim> >&)
477 {
478  return s.n_gradients();
479 }
480 
481 template<int dim>
482 inline
483 unsigned int
484 selector_size(const VectorSelector& s,
485  const std::vector< Tensor<2, dim> >&)
486 {
487  return s.n_hessians();
488 }
489 
490 inline
491 unsigned int
492 selector_index(const VectorSelector& s,
493  unsigned int i,
494  const std::vector<double>&)
495 {
496  return s.value_index(i);
497 }
498 
499 template<int dim>
500 inline
501 unsigned int
502 selector_index(const VectorSelector& s,
503  unsigned int i,
504  const std::vector< Tensor<1, dim> >&)
505 {
506  return s.gradient_index(i);
507 }
508 
509 template<int dim>
510 inline
511 unsigned int
512 selector_index(const VectorSelector& s,
513  unsigned int i,
514  const std::vector< Tensor<2, dim> >&)
515 {
516  return s.hessian_index(i);
517 }
518 
519 namespace InfoObjects
520 {
521 
553 template<int dim, typename TYPE>
554 inline
555 void
556 fill_data(const FEValuesBase<dim>& fe_values,
557  const FEVector& fe_vector,
558  const std::vector<unsigned int>& local_dof_indices,
559  unsigned int first_index,
560  unsigned int n_indices,
561  std::vector< std::vector<TYPE> >& result)
562 {
563  Assert(false, ExcNotImplemented());
564 }
565 
566 template<int dim>
567 inline
568 void
569 fill_data(const FEValuesBase<dim>& fe_values,
570  const FEVector& fe_vector,
571  const std::vector<unsigned int>& local_dof_indices,
572  unsigned int first_index,
573  unsigned int n_indices,
574  std::vector< std::vector<double> >& result)
575 {
576  fe_values.get_function_values(fe_vector,
577  make_slice(local_dof_indices, first_index, n_indices),
578  result,
579  true);
580 }
581 
582 template<int dim>
583 inline
584 void
585 fill_data(const FEValuesBase<dim>& fe_values,
586  const FEVector& fe_vector,
587  const std::vector<unsigned int>& local_dof_indices,
588  unsigned int first_index,
589  unsigned int n_indices,
590  std::vector< std::vector< Tensor<1,dim> > >& result)
591 {
592  fe_values.get_function_grads(fe_vector,
593  make_slice(local_dof_indices, first_index, n_indices),
594  result,
595  true);
596 }
597 
598 template<int dim>
599 inline
600 void
601 fill_data(const FEValuesBase<dim>& fe_values,
602  const FEVector& fe_vector,
603  const std::vector<unsigned int>& local_dof_indices,
604  unsigned int first_index,
605  unsigned int n_indices,
606  std::vector< std::vector< Tensor<2,dim> > >& result)
607 {
608  fe_values.get_function_hessians(fe_vector,
609  make_slice(local_dof_indices, first_index, n_indices),
610  result,
611  true);
612 }
613 
640 template<int dim, int spacedim = dim>
641 class DoFInfo
642 {
643 public:
644 
648  DoFInfo(const BlockInfo& block_info);
649 
653  DoFInfo(const FEVectors&,
654  const BlockInfo& block_info);
655 
660  template<typename DHCellIterator>
661  void reinit(const DHCellIterator& c)
662  {
663  get_indices(c);
664 
665  dof_cell = c;
666  dof_active_cell = c;
667  cell = static_cast<typename Triangulation<dim>::cell_iterator> (c);
668 
669  face_number = deal_II_numbers::invalid_unsigned_int;
670  sub_number = deal_II_numbers::invalid_unsigned_int;
671  }
672 
678  template<typename DHCellIterator, typename DHFaceIterator>
679  void reinit(const DHCellIterator& c,
680  const DHFaceIterator& f,
681  const unsigned int n)
682  {
683  if( cell.state() != IteratorState::valid || cell != static_cast<typename Triangulation<dim>::cell_iterator> (c) )
684  {
685  get_indices(c);
686  }
687 
688  dof_cell = c;
689  dof_active_cell = c;
690  cell = static_cast<typename Triangulation<dim>::cell_iterator> (c);
691 
692  dof_face = f;
693  face = static_cast<typename Triangulation<dim>::face_iterator> (f);
694 
695  face_number = n;
696  sub_number = deal_II_numbers::invalid_unsigned_int;
697  }
698 
704  template<typename DHCellIterator, typename DHFaceIterator>
705  void reinit(const DHCellIterator& c,
706  const DHFaceIterator& f,
707  const unsigned int n,
708  const unsigned int s)
709  {
710  if( cell.state() != IteratorState::valid || cell != static_cast<typename Triangulation<dim>::cell_iterator> (c) )
711  {
712  get_indices(c);
713  }
714 
715  dof_cell = c;
716  dof_active_cell = c;
717  cell = static_cast<typename Triangulation<dim>::cell_iterator> (c);
718 
719  dof_face = f;
720  face = static_cast<typename Triangulation<dim>::face_iterator> (f);
721 
722  face_number = n;
723  sub_number = s;
724  }
725 
729  typename DoFHandler<dim>::cell_iterator dof_cell;
730 
734  typename DoFHandler<dim>::active_cell_iterator dof_active_cell;
735 
739  typename DoFHandler<dim>::face_iterator dof_face;
740 
744  typename Triangulation<dim>::cell_iterator cell;
745 
749  typename Triangulation<dim>::face_iterator face;
750 
760  unsigned int face_number;
761 
772  unsigned int sub_number;
773 
778  std::vector<unsigned int> indices;
779 
783  SmartPointer<const BlockInfo> block_info;
784 
785 private:
786 
788  void get_indices(const typename DoFHandler<dim, spacedim>::cell_iterator c);
789 
791  void get_indices(const typename MGDoFHandler<dim, spacedim>::cell_iterator c);
792 
796  std::vector<unsigned int> indices_org;
797 };
798 
799 // --- implementation ---
800 
801 template<int dim, int spacedim>
802 inline
803 DoFInfo<dim, spacedim>::DoFInfo(const BlockInfo& block_info)
804 :
805 block_info(&block_info,
806  typeid(*this).name())
807 { }
808 
809 template<int dim, int spacedim>
810 inline
811 DoFInfo<dim, spacedim>::DoFInfo(const FEVectors&,
812  const BlockInfo& block_info)
813 :
814 block_info(&block_info,
815  typeid(*this).name())
816 { }
817 
818 template<int dim, int spacedim>
819 inline
820 void
821 DoFInfo<dim, spacedim>::get_indices(const typename DoFHandler<dim, spacedim>::cell_iterator c)
822 {
823  indices.resize(c->get_fe().dofs_per_cell);
824 
825  if(block_info->local_renumbering.size() == 0)
826  {
827  c->get_dof_indices(indices);
828  }
829  else
830  {
831  indices_org.resize(c->get_fe().dofs_per_cell);
832  c->get_dof_indices(indices_org);
833 
834  for(unsigned int i = 0; i < indices.size(); ++i)
835  indices[this->block_info->local_renumbering[i]] = indices_org[i];
836  }
837 }
838 
839 template<int dim, int spacedim>
840 inline
841 void
842 DoFInfo<dim, spacedim>::get_indices(const typename MGDoFHandler<dim, spacedim>::cell_iterator c)
843 {
844  indices.resize(c->get_fe().dofs_per_cell);
845 
846  if(block_info->local_renumbering.size() == 0)
847  {
848  c->get_mg_dof_indices(indices);
849  }
850  else
851  {
852  indices_org.resize(c->get_fe().dofs_per_cell);
853  c->get_mg_dof_indices(indices_org);
854 
855  for(unsigned int i = 0; i < indices.size(); ++i)
856  indices[this->block_info->local_renumbering[i]] = indices_org[i];
857  }
858 }
859 
906 template<int dim, typename FEVALUESBASE>
907 class IntegrationInfo : public DoFInfo<dim>
908 {
909 public:
910 
914  IntegrationInfo(const BlockInfo& block_info);
915 
919  IntegrationInfo(const FEVectors& data,
920  const BlockInfo& block_info);
921 
964  template<typename FEVALUES>
965  void initialize(const FEVALUES* FE_VALUES,
966  const FiniteElement<dim>& fe,
967  const Mapping<dim>& mapping,
968  const Quadrature<FEVALUES::integral_dimension>& quadrature,
969  const UpdateFlags flags);
970 
977  void initialize_selector(const VectorSelector& s)
978  {
979  selector = s;
980 
981  if(global_data)
982  selector.cache(*global_data);
983  }
984 
991  void initialize_data(const FEVectors& data)
992  {
993  global_data = &data;
994  selector.cache(*global_data);
995  }
996 
1001  void clear()
1002  {
1003  fevalv.resize(0);
1004  }
1005 
1007  const FEVALUESBASE& fe() const
1008  {
1009  AssertDimension(fevalv.size(), 1);
1010  return *fevalv[0];
1011  }
1012 
1014  const FEVALUESBASE& fe(unsigned int i) const
1015  {
1016  Assert( i < fevalv.size(), ExcIndexRange(i, 0, fevalv.size()) );
1017  return *fevalv[i];
1018  }
1019 
1021  const FEVALUESBASE& get_fe_val_unsplit() const
1022  {
1023  return *fe_val_unsplit;
1024  }
1025 
1030  template<typename DHCellIterator>
1031  void reinit(const DHCellIterator& c);
1032 
1037  template<typename DHCellIterator, typename DHFaceIterator>
1038  void reinit(const DHCellIterator& c,
1039  const DHFaceIterator& f,
1040  const unsigned int fn);
1041 
1046  template<typename DHCellIterator, typename DHFaceIterator>
1047  void reinit(const DHCellIterator& c,
1048  const DHFaceIterator& f,
1049  const unsigned int fn,
1050  const unsigned int sn);
1051 
1064  template<typename TYPE>
1065  void fill_local_data(std::vector< std::vector< std::vector<TYPE> > >& data,
1066  bool split_fevalues) const;
1067 
1069  bool multigrid;
1070 
1080  SmartPointer<const FEVectors> global_data;
1081 
1090  VectorSelector selector;
1091 
1100  std::vector< std::vector< std::vector<double> > > values;
1101 
1110  std::vector< std::vector< std::vector< Tensor<1,dim> > > > gradients;
1111 
1113  std::vector< std::vector< std::vector< Tensor<1,dim> > > >& derivatives;
1114 
1123  std::vector< std::vector< std::vector< Tensor<2,dim> > > > hessians;
1124 
1125 private:
1126 
1135  template<typename TYPE>
1136  void fill_local_data(std::vector< std::vector< std::vector<TYPE> > >& data,
1137  VectorSelector& selector,
1138  bool split_fevalues) const;
1139 
1141  std::vector< boost::shared_ptr<FEVALUESBASE> > fevalv;
1142 
1149  boost::shared_ptr<FEVALUESBASE> fe_val_unsplit;
1150 };
1151 
1152 // --- implementation ---
1153 
1154 template<int dim, typename FEVALUESBASE>
1155 inline
1156 IntegrationInfo<dim, FEVALUESBASE>::IntegrationInfo(const BlockInfo& block_info)
1157 :
1158 DoFInfo<dim>(block_info),
1159 multigrid(false),
1160 global_data(0, typeid(*this).name()),
1161 derivatives(gradients),
1162 fevalv(0)
1163 { }
1164 
1165 template<int dim, typename FEVALUESBASE>
1166 inline
1167 IntegrationInfo<dim, FEVALUESBASE>::IntegrationInfo(const FEVectors& data,
1168  const BlockInfo& block_info)
1169 :
1170 DoFInfo<dim>(block_info),
1171 multigrid(false),
1172 global_data(&data),
1173 derivatives(gradients),
1174 fevalv(0)
1175 { }
1176 
1177 template<int dim, typename FEVALUESBASE>
1178 template<typename FEVALUES>
1179 inline
1180 void
1181 IntegrationInfo<dim, FEVALUESBASE>::initialize(const FEVALUES*,
1182  const FiniteElement<dim>& fe,
1183  const Mapping<dim>& mapping,
1184  const Quadrature<FEVALUES::integral_dimension>& quadrature,
1185  const UpdateFlags flags)
1186 {
1187  if( this->block_info->local_renumbering.size() != 0 ) // mode # 1 - the main mode
1188  {
1189  fevalv.resize(fe.n_base_elements());
1190  for(unsigned int i = 0; i < fevalv.size(); ++i)
1191  {
1192  fevalv[i] = boost::shared_ptr<FEVALUESBASE>( new FEVALUES(mapping,
1193  fe.base_element(i),
1194  quadrature,
1195  flags)
1196  );
1197  }
1198 
1199  // --- (ZVN) fe_val_unsplit initialization ---
1200 
1201  fe_val_unsplit = boost::shared_ptr<FEVALUESBASE>( new FEVALUES(mapping,
1202  fe,
1203  quadrature,
1204  flags)
1205  );
1206  }
1207  else // mode # 2
1208  {
1209  fevalv.resize(1);
1210  fevalv[0] = boost::shared_ptr<FEVALUESBASE>( new FEVALUES(mapping,
1211  fe,
1212  quadrature,
1213  flags)
1214  );
1215  }
1216 
1217  values.resize(global_data->n_vectors(),
1218  std::vector< std::vector<double> >(fe.n_components(),
1219  std::vector<double>(quadrature.size())
1220  )
1221  );
1222  gradients.resize(global_data->n_vectors(),
1223  std::vector< std::vector< Tensor<1,dim> > >(fe.n_components(),
1224  std::vector< Tensor<1,dim> >(quadrature.size())
1225  )
1226  );
1227  hessians.resize(global_data->n_vectors(),
1228  std::vector< std::vector< Tensor<2,dim> > >(fe.n_components(),
1229  std::vector< Tensor<2,dim> >(quadrature.size())
1230  )
1231  );
1232 }
1233 
1234 template<int dim, typename FEVALUESBASE>
1235 template<typename DHCellIterator>
1236 inline
1237 void
1238 IntegrationInfo<dim, FEVALUESBASE>::reinit(const DHCellIterator& c)
1239 {
1240  DoFInfo<dim>::reinit(c);
1241 
1242  for(unsigned int i = 0; i < fevalv.size(); ++i)
1243  {
1244  FEVALUESBASE& fe_values_base = *fevalv[i];
1245  FEValues<dim>& fe_values = dynamic_cast<FEValues<dim>&>(fe_values_base);
1246  fe_values.reinit(this->cell);
1247  }
1248 
1249  // --- (ZVN) ---
1250 
1251  FEVALUESBASE& fe_values_base_unsplit = *fe_val_unsplit;
1252  FEValues<dim>& fe_values_unsplit = dynamic_cast<FEValues<dim>&>(fe_values_base_unsplit);
1253  fe_values_unsplit.reinit(this->dof_active_cell);
1254 
1255  // --- (ZVN) ---
1256 
1257  const bool split_fevalues = this->block_info->local_renumbering.size() != 0;
1258  fill_local_data(values,
1259  selector,
1260  split_fevalues);
1261  fill_local_data(gradients,
1262  selector,
1263  split_fevalues);
1264  fill_local_data(hessians,
1265  selector,
1266  split_fevalues);
1267 }
1268 
1269 template<int dim, typename FEVALUESBASE>
1270 template<typename DHCellIterator, typename DHFaceIterator>
1271 inline
1272 void
1273 IntegrationInfo<dim, FEVALUESBASE>::reinit(const DHCellIterator& c,
1274  const DHFaceIterator& f,
1275  const unsigned int fn)
1276 {
1277  DoFInfo<dim>::reinit(c,
1278  f,
1279  fn);
1280 
1281  for(unsigned int i = 0; i < fevalv.size(); ++i)
1282  {
1283  FEVALUESBASE& fe_values_base = *fevalv[i];
1284  FEFaceValues<dim>& fe_face_values = dynamic_cast<FEFaceValues<dim>&>(fe_values_base);
1285  fe_face_values.reinit(this->cell,
1286  fn);
1287  }
1288 
1289  // --- (ZVN) ---
1290 
1291  FEVALUESBASE& fe_values_base_unsplit = *fe_val_unsplit;
1292  FEFaceValues<dim>& fe_face_values_unsplit = dynamic_cast<FEFaceValues<dim>&>(fe_values_base_unsplit);
1293  fe_face_values_unsplit.reinit(this->dof_active_cell,
1294  fn);
1295 
1296  // --- (ZVN) ---
1297 
1298  const bool split_fevalues = this->block_info->local_renumbering.size() != 0;
1299  fill_local_data(values,
1300  selector,
1301  split_fevalues);
1302  fill_local_data(gradients,
1303  selector,
1304  split_fevalues);
1305  fill_local_data(hessians,
1306  selector,
1307  split_fevalues);
1308 }
1309 
1310 template<int dim, typename FEVALUESBASE>
1311 template<typename DHCellIterator, typename DHFaceIterator>
1312 inline
1313 void
1314 IntegrationInfo<dim, FEVALUESBASE>::reinit(const DHCellIterator& c,
1315  const DHFaceIterator& f,
1316  const unsigned int fn,
1317  const unsigned int sn)
1318 {
1319  DoFInfo<dim>::reinit(c,
1320  f,
1321  fn,
1322  sn);
1323 
1324  for(unsigned int i = 0; i < fevalv.size(); ++i)
1325  {
1326  FEVALUESBASE& fe_values_base = *fevalv[i];
1327  FESubfaceValues<dim>& fe_subface_values = dynamic_cast<FESubfaceValues<dim>&>(fe_values_base);
1328  fe_subface_values.reinit(this->cell,
1329  fn,
1330  sn);
1331  }
1332 
1333  // --- (ZVN) ---
1334 
1335  FEVALUESBASE& fe_values_base_unsplit = *fe_val_unsplit;
1336  FESubfaceValues<dim>& fe_subface_values_unsplit = dynamic_cast<FESubfaceValues<dim>&>(fe_values_base_unsplit);
1337  fe_subface_values_unsplit.reinit(this->dof_active_cell,
1338  fn,
1339  sn);
1340 
1341  // --- (ZVN) ---
1342 
1343  const bool split_fevalues = this->block_info->local_renumbering.size() != 0;
1344  fill_local_data(values,
1345  selector,
1346  split_fevalues);
1347  fill_local_data(gradients,
1348  selector,
1349  split_fevalues);
1350  fill_local_data(hessians,
1351  selector,
1352  split_fevalues);
1353 }
1354 
1355 template<int dim, typename FEVALUESBASE>
1356 template<typename TYPE>
1357 inline
1358 void
1359 IntegrationInfo<dim, FEVALUESBASE>::fill_local_data(std::vector< std::vector< std::vector<TYPE> > >& data,
1360  bool split_fevalues) const
1361 {
1362  if(split_fevalues)
1363  {
1364  std::vector< std::vector<TYPE> > local_data;
1365 
1366  unsigned int comp = 0;
1367  for(unsigned int b = 0; b < this->block_info->base_element.size(); ++b)
1368  {
1369  const unsigned int no_fe = this->block_info->base_element[b];
1370  const FEValuesBase<dim>& fe_values_base = this->fe(no_fe);
1371 
1372  const unsigned int n_comp = fe_values_base.get_fe().n_components();
1373  local_data.resize(n_comp,
1374  std::vector<TYPE>(fe_values_base.n_quadrature_points));
1375 
1376  for(unsigned int i = 0; i < this->global_data->n_vectors(); ++i)
1377  {
1378  const FEVector& src = this->global_data->vector(i);
1379  MeshWorker::InfoObjects::fill_data(fe_values_base,
1380  src,
1381  this->indices,
1382  this->block_info->local.block_start(b),
1383  fe_values_base.dofs_per_cell,
1384  local_data);
1385 
1386  for(unsigned int c = 0; c < local_data.size(); ++c)
1387  for(unsigned int k = 0; k < local_data[c].size(); ++k)
1388  data[i][comp+c][k] = local_data[c][k];
1389  }
1390 
1391  comp += n_comp;
1392  }
1393  }
1394  else
1395  {
1396  for(unsigned int i = 0; i < this->global_data->n_vectors(); ++i)
1397  {
1398  const FEVector& src = this->global_data->vector(i);
1399  const FEValuesBase<dim>& fe_values_base = this->fe();
1400 
1401  MeshWorker::InfoObjects::fill_data(fe_values_base,
1402  src,
1403  this->indices,
1404  0,
1405  fe_values_base.get_fe().dofs_per_cell,
1406  data[i]);
1407  }
1408  }
1409 }
1410 
1411 template<int dim, typename FEVALUESBASE>
1412 template<typename TYPE>
1413 inline
1414 void
1415 IntegrationInfo<dim, FEVALUESBASE>::fill_local_data(std::vector< std::vector< std::vector<TYPE> > >& data,
1416  VectorSelector& selector,
1417  bool split_fevalues) const
1418 {
1419  if(split_fevalues)
1420  {
1421  std::vector< std::vector<TYPE> > local_data;
1422 
1423  unsigned int comp = 0;
1424  for(unsigned int b = 0; b < this->block_info->base_element.size(); ++b)
1425  {
1426  const unsigned int no_fe = this->block_info->base_element[b];
1427  const FEValuesBase<dim>& fe_values_base = this->fe(no_fe);
1428 
1429  const unsigned int n_comp = fe_values_base.get_fe().n_components();
1430  local_data.resize(n_comp,
1431  std::vector<TYPE>(fe_values_base.n_quadrature_points));
1432 
1433  const unsigned int n = selector_size(selector,
1434  local_data[0]);
1435 
1436  for(unsigned int i = 0; i < n; ++i)
1437  {
1438  const FEVector& src = this->global_data->vector( selector_index(selector, i, local_data[0]) );
1439  MeshWorker::InfoObjects::fill_data(fe_values_base,
1440  src,
1441  this->indices,
1442  this->block_info->local.block_start(b),
1443  fe_values_base.dofs_per_cell,
1444  local_data);
1445 
1446  for(unsigned int c = 0; c < local_data.size(); ++c)
1447  for(unsigned int k = 0; k < local_data[c].size(); ++k)
1448  data[i][comp+c][k] = local_data[c][k];
1449  }
1450 
1451  comp += n_comp;
1452  }
1453  }
1454  else
1455  {
1456  for(unsigned int i = 0; i < this->global_data->n_vectors(); ++i)
1457  {
1458  const FEVector& src = this->global_data->vector(i);
1459  const FEValuesBase<dim>& fe_values_base = this->fe();
1460 
1461  MeshWorker::InfoObjects::fill_data(fe_values_base,
1462  src,
1463  this->indices,
1464  0,
1465  fe_values_base.get_fe().dofs_per_cell,
1466  data[i]);
1467  }
1468  }
1469 }
1470 
1478 template<int dim>
1479 class IntegrationInfoBox
1480 {
1481 public:
1482 
1486  IntegrationInfoBox(const BlockInfo& block_info);
1487 
1493  template<typename WORKER>
1494  void initialize(const WORKER& integrator,
1495  const FiniteElement<dim>& fe,
1496  const Mapping<dim>& mapping,
1497  const FEVectors& global_data);
1498 
1503  IntegrationInfo<dim, FEValuesBase<dim> > cell_info;
1504  IntegrationInfo<dim, FEFaceValuesBase<dim> > bdry_info;
1505  IntegrationInfo<dim, FEFaceValuesBase<dim> > face_info;
1506  IntegrationInfo<dim, FEFaceValuesBase<dim> > subface_info;
1507  IntegrationInfo<dim, FEFaceValuesBase<dim> > neighbor_info;
1508 };
1509 
1510 // --- implementation ---
1511 
1512 template<int dim>
1513 inline
1514 IntegrationInfoBox<dim>::IntegrationInfoBox(const BlockInfo& block_info)
1515 :
1516 cell_info(block_info),
1517 bdry_info(block_info),
1518 face_info(block_info),
1519 subface_info(block_info),
1520 neighbor_info(block_info)
1521 { }
1522 
1523 template<int dim>
1524 template<typename WORKER>
1525 inline
1526 void
1527 IntegrationInfoBox<dim>::initialize(const WORKER& integrator,
1528  const FiniteElement<dim>& fe,
1529  const Mapping<dim>& mapping,
1530  const FEVectors& global_data)
1531 {
1532  FEValues<dim>* fe_values = 0;
1533  FEFaceValues<dim>* fe_face_values = 0;
1534  FESubfaceValues<dim>* fe_subface_values = 0;
1535 
1536  cell_info.initialize_selector(integrator.cell_selector);
1537  bdry_info.initialize_selector(integrator.bdry_selector);
1538  face_info.initialize_selector(integrator.face_selector);
1539  subface_info.initialize_selector(integrator.face_selector);
1540  neighbor_info.initialize_selector(integrator.face_selector);
1541 
1542  cell_info.initialize_data(global_data);
1543  bdry_info.initialize_data(global_data);
1544  face_info.initialize_data(global_data);
1545  subface_info.initialize_data(global_data);
1546  neighbor_info.initialize_data(global_data);
1547 
1548  cell_info.initialize(fe_values,
1549  fe,
1550  mapping,
1551  integrator.cell_quadrature,
1552  integrator.cell_flags);
1553 
1554  bdry_info.initialize(fe_face_values,
1555  fe,
1556  mapping,
1557  integrator.bdry_quadrature,
1558  integrator.face_flags);
1559 
1560  face_info.initialize(fe_face_values,
1561  fe,
1562  mapping,
1563  integrator.face_quadrature,
1564  integrator.face_flags);
1565 
1566  subface_info.initialize(fe_subface_values,
1567  fe,
1568  mapping,
1569  integrator.face_quadrature,
1570  integrator.face_flags);
1571 
1572  neighbor_info.initialize(fe_face_values,
1573  fe,
1574  mapping,
1575  integrator.face_quadrature,
1576  integrator.ngbr_flags);
1577 }
1578 
1579 } // namespace InfoObjects
1580 
1581 } // namespace MeshWorker
1582 
1583 #endif
Generated on Fri Aug 30 2013 17:15:54 for OpenFCST: The open-source Fuel Cell Simulation Toolbox by   doxygen 1.8.2