Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
dof_accessor.templates.h
1 // ---------------------------------------------------------------------
2 // @f$Id: dof_accessor.templates.h 31932 2013-12-08 02:15:54Z heister @f$
3 //
4 // Copyright (C) 1999 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef __deal2__dof_accessor_templates_h
18 #define __deal2__dof_accessor_templates_h
19 
20 
21 #include <deal.II/base/config.h>
22 #include <deal.II/base/types.h>
23 #include <deal.II/lac/constraint_matrix.h>
24 #include <deal.II/dofs/dof_accessor.h>
25 #include <deal.II/dofs/dof_levels.h>
26 #include <deal.II/dofs/dof_faces.h>
27 #include <deal.II/hp/dof_level.h>
28 #include <deal.II/hp/dof_faces.h>
29 #include <deal.II/grid/tria_iterator.h>
30 #include <deal.II/grid/tria_iterator.templates.h>
31 
32 #include <vector>
33 #include <limits>
34 
35 DEAL_II_NAMESPACE_OPEN
36 
37 
38 /*------------------------- Functions: DoFAccessor ---------------------------*/
39 
40 
41 template <int structdim, class DH, bool lda>
42 inline
44 {
45  Assert (false, ExcInvalidObject());
46 }
47 
48 
49 
50 template <int structdim, class DH, bool lda>
51 inline
54  const int level,
55  const int index,
56  const DH *dof_handler)
57  :
58  ::internal::DoFAccessor::Inheritance<structdim,DH::dimension,
59  DH::space_dimension>::BaseClass (tria,
60  level,
61  index),
62  dof_handler(const_cast<DH *>(dof_handler))
63 {}
64 
65 
66 
67 template <int structdim, class DH, bool lda>
68 template <int structdim2, int dim2, int spacedim2>
70 {
71  Assert (false, ExcInvalidObject());
72 }
73 
74 
75 
76 template <int structdim, class DH, bool lda>
77 template <int dim2, class DH2, bool lda2>
78 inline
80  : BaseClass(other), dof_handler(const_cast<DH *>(other.dof_handler))
81 {
82  Assert (false, ExcInvalidObject());
83 }
84 
85 
86 
87 template <int structdim, class DH, bool lda>
88 template <bool lda2>
89 inline
91  : BaseClass(other), dof_handler(const_cast<DH *>(other.dof_handler))
92 {
93 }
94 
95 
96 
97 template <int structdim, class DH, bool lda>
98 inline
99 void
101 {
102  Assert (dh != 0, ExcInvalidObject());
103  this->dof_handler = dh;
104 }
105 
106 
107 
108 template <int structdim, class DH, bool lda>
109 inline
110 const DH &
112 {
113  return *this->dof_handler;
114 }
115 
116 
117 
118 template <int structdim, class DH, bool lda>
119 inline
120 void
123 {
124  Assert (this->dof_handler != 0, ExcInvalidObject());
125  BaseClass::copy_from(da);
126 }
127 
128 
129 
130 template <int structdim, class DH, bool lda>
131 template <bool lda2>
132 inline
133 void
135 {
136  BaseClass::copy_from (a);
137  set_dof_handler (a.dof_handler);
138 }
139 
140 
141 
142 template <int structdim, class DH, bool lda>
143 template <int dim2, class DH2, bool lda2>
144 inline
145 bool
147 {
148  Assert (structdim == dim2, ExcCantCompareIterators());
149  Assert (this->dof_handler == a.dof_handler, ExcCantCompareIterators());
150  return (BaseClass::operator == (a));
151 }
152 
153 
154 
155 template <int structdim, class DH, bool lda>
156 template <int dim2, class DH2, bool lda2>
157 inline
158 bool
160 {
161  Assert (structdim == dim2, ExcCantCompareIterators());
162  Assert (this->dof_handler == a.dof_handler, ExcCantCompareIterators());
163  return (BaseClass::operator != (a));
164 }
165 
166 
167 
168 template <int structdim, class DH, bool lda>
169 inline
171 DoFAccessor<structdim,DH,lda>::child (const unsigned int i) const
172 {
173  Assert (static_cast<unsigned int>(this->level()) < this->dof_handler->levels.size(),
174  ExcMessage ("DoFHandler not initialized"));
175 
178 
179  TriaIterator<DoFAccessor<structdim,DH,lda> > q (*t, this->dof_handler);
180  return q;
181 }
182 
183 
184 template <int structdim, class DH, bool lda>
185 inline
188 {
189  Assert (static_cast<unsigned int>(this->level()) < this->dof_handler->levels.size(),
190  ExcMessage ("DoFHandler not initialized"));
191  Assert (this->level () > 0,
192  ExcMessage ("Cell is at coarsest level."));
193 
194  int previous_level;
195 
196  if (DH::dimension==structdim)
197  previous_level = this->level () - 1;
198 
199  else
200  previous_level = 0;
201 
203  previous_level,
204  this->parent_index (),
205  this->dof_handler);
206 
207  return q;
208 }
209 
210 
211 
212 namespace internal
213 {
214  namespace DoFAccessor
215  {
222  {
227  template <int spacedim>
228  static
230  get_dof_index (const ::DoFHandler<1,spacedim> &dof_handler,
231  const unsigned int obj_level,
232  const unsigned int obj_index,
233  const unsigned int fe_index,
234  const unsigned int local_index,
236  {
237  return dof_handler.levels[obj_level]->dof_object.
238  get_dof_index (dof_handler,
239  obj_index,
240  fe_index,
241  local_index);
242  }
243 
244 
245  template <int spacedim>
246  static
247  void
248  set_dof_index (const ::DoFHandler<1,spacedim> &dof_handler,
249  const unsigned int obj_level,
250  const unsigned int obj_index,
251  const unsigned int fe_index,
252  const unsigned int local_index,
254  const types::global_dof_index global_index)
255  {
256  dof_handler.levels[obj_level]->dof_object.
257  set_dof_index (dof_handler,
258  obj_index,
259  fe_index,
260  local_index,
261  global_index);
262  }
263 
264 
265  template <int spacedim>
266  static
268  get_dof_index (const ::DoFHandler<2,spacedim> &dof_handler,
269  const unsigned int obj_level,
270  const unsigned int obj_index,
271  const unsigned int fe_index,
272  const unsigned int local_index,
274  {
275  // faces have no levels
276  Assert (obj_level == 0, ExcInternalError());
277  return dof_handler.faces->lines.
278  get_dof_index (dof_handler,
279  obj_index,
280  fe_index,
281  local_index);
282  }
283 
284 
285  template <int spacedim>
286  static
287  void
288  set_dof_index (const ::DoFHandler<2,spacedim> &dof_handler,
289  const unsigned int obj_level,
290  const unsigned int obj_index,
291  const unsigned int fe_index,
292  const unsigned int local_index,
294  const types::global_dof_index global_index)
295  {
296  // faces have no levels
297  Assert (obj_level == 0, ExcInternalError());
298  dof_handler.faces->lines.
299  set_dof_index (dof_handler,
300  obj_index,
301  fe_index,
302  local_index,
303  global_index);
304  }
305 
306 
307  template <int spacedim>
308  static
310  get_dof_index (const ::DoFHandler<2,spacedim> &dof_handler,
311  const unsigned int obj_level,
312  const unsigned int obj_index,
313  const unsigned int fe_index,
314  const unsigned int local_index,
316  {
317  return dof_handler.levels[obj_level]->dof_object.
318  get_dof_index (dof_handler,
319  obj_index,
320  fe_index,
321  local_index);
322  }
323 
324 
325  template <int spacedim>
326  static
327  void
328  set_dof_index (const ::DoFHandler<2,spacedim> &dof_handler,
329  const unsigned int obj_level,
330  const unsigned int obj_index,
331  const unsigned int fe_index,
332  const unsigned int local_index,
334  const types::global_dof_index global_index)
335  {
336  dof_handler.levels[obj_level]->dof_object.
337  set_dof_index (dof_handler,
338  obj_index,
339  fe_index,
340  local_index,
341  global_index);
342  }
343 
344 
345  template <int spacedim>
346  static
348  get_dof_index (const ::DoFHandler<3,spacedim> &dof_handler,
349  const unsigned int obj_level,
350  const unsigned int obj_index,
351  const unsigned int fe_index,
352  const unsigned int local_index,
354  {
355  // faces have no levels
356  Assert (obj_level == 0, ExcInternalError());
357  return dof_handler.faces->lines.
358  get_dof_index (dof_handler,
359  obj_index,
360  fe_index,
361  local_index);
362  }
363 
364 
365  template <int spacedim>
366  static
367  void
368  set_dof_index (const ::DoFHandler<3,spacedim> &dof_handler,
369  const unsigned int obj_level,
370  const unsigned int obj_index,
371  const unsigned int fe_index,
372  const unsigned int local_index,
374  const types::global_dof_index global_index)
375  {
376  // faces have no levels
377  Assert (obj_level == 0, ExcInternalError());
378  dof_handler.faces->lines.
379  set_dof_index (dof_handler,
380  obj_index,
381  fe_index,
382  local_index,
383  global_index);
384  }
385 
386 
387 
388  template <int spacedim>
389  static
391  get_dof_index (const ::DoFHandler<3,spacedim> &dof_handler,
392  const unsigned int obj_level,
393  const unsigned int obj_index,
394  const unsigned int fe_index,
395  const unsigned int local_index,
397  {
398  // faces have no levels
399  Assert (obj_level == 0, ExcInternalError());
400  return dof_handler.faces->quads.
401  get_dof_index (dof_handler,
402  obj_index,
403  fe_index,
404  local_index);
405  }
406 
407 
408  template <int spacedim>
409  static
410  void
411  set_dof_index (const ::DoFHandler<3,spacedim> &dof_handler,
412  const unsigned int obj_level,
413  const unsigned int obj_index,
414  const unsigned int fe_index,
415  const unsigned int local_index,
417  const types::global_dof_index global_index)
418  {
419  // faces have no levels
420  Assert (obj_level == 0, ExcInternalError());
421  dof_handler.faces->quads.
422  set_dof_index (dof_handler,
423  obj_index,
424  fe_index,
425  local_index,
426  global_index);
427  }
428 
429 
430 
431  template <int spacedim>
432  static
434  get_dof_index (const ::DoFHandler<3,spacedim> &dof_handler,
435  const unsigned int obj_level,
436  const unsigned int obj_index,
437  const unsigned int fe_index,
438  const unsigned int local_index,
440  {
441  return dof_handler.levels[obj_level]->dof_object.
442  get_dof_index (dof_handler,
443  obj_index,
444  fe_index,
445  local_index);
446  }
447 
448 
449  template <int spacedim>
450  static
451  void
452  set_dof_index (const ::DoFHandler<3,spacedim> &dof_handler,
453  const unsigned int obj_level,
454  const unsigned int obj_index,
455  const unsigned int fe_index,
456  const unsigned int local_index,
458  const types::global_dof_index global_index)
459  {
460  dof_handler.levels[obj_level]->dof_object.
461  set_dof_index (dof_handler,
462  obj_index,
463  fe_index,
464  local_index,
465  global_index);
466  }
467 
468 
469  template <int spacedim>
470  static
472  get_dof_index (const ::hp::DoFHandler<1,spacedim> &dof_handler,
473  const unsigned int obj_level,
474  const unsigned int obj_index,
475  const unsigned int fe_index,
476  const unsigned int local_index,
477  const ::internal::int2type<1> &)
478  {
479  return dof_handler.levels[obj_level]->
480  get_dof_index (obj_index,
481  fe_index,
482  local_index);
483  }
484 
485 
486  template <int spacedim>
487  static
488  void
489  set_dof_index (const ::hp::DoFHandler<1,spacedim> &dof_handler,
490  const unsigned int obj_level,
491  const unsigned int obj_index,
492  const unsigned int fe_index,
493  const unsigned int local_index,
494  const ::internal::int2type<1> &,
495  const types::global_dof_index global_index)
496  {
497  dof_handler.levels[obj_level]->
498  set_dof_index (obj_index,
499  fe_index,
500  local_index,
501  global_index);
502  }
503 
504 
505  template <int spacedim>
506  static
508  get_dof_index (const ::hp::DoFHandler<2,spacedim> &dof_handler,
509  const unsigned int obj_level,
510  const unsigned int obj_index,
511  const unsigned int fe_index,
512  const unsigned int local_index,
513  const ::internal::int2type<1> &)
514  {
515  return dof_handler.faces->lines.
516  get_dof_index (dof_handler,
517  obj_index,
518  fe_index,
519  local_index,
520  obj_level);
521  }
522 
523 
524  template <int spacedim>
525  static
526  void
527  set_dof_index (const ::hp::DoFHandler<2,spacedim> &dof_handler,
528  const unsigned int obj_level,
529  const unsigned int obj_index,
530  const unsigned int fe_index,
531  const unsigned int local_index,
532  const ::internal::int2type<1> &,
533  const types::global_dof_index global_index)
534  {
535  dof_handler.faces->lines.
536  set_dof_index (dof_handler,
537  obj_index,
538  fe_index,
539  local_index,
540  global_index,
541  obj_level);
542  }
543 
544 
545  template <int spacedim>
546  static
548  get_dof_index (const ::hp::DoFHandler<2,spacedim> &dof_handler,
549  const unsigned int obj_level,
550  const unsigned int obj_index,
551  const unsigned int fe_index,
552  const unsigned int local_index,
553  const ::internal::int2type<2> &)
554  {
555  return dof_handler.levels[obj_level]->
556  get_dof_index (obj_index,
557  fe_index,
558  local_index);
559  }
560 
561 
562  template <int spacedim>
563  static
564  void
565  set_dof_index (const ::hp::DoFHandler<2,spacedim> &dof_handler,
566  const unsigned int obj_level,
567  const unsigned int obj_index,
568  const unsigned int fe_index,
569  const unsigned int local_index,
570  const ::internal::int2type<2> &,
571  const types::global_dof_index global_index)
572  {
573  dof_handler.levels[obj_level]->
574  set_dof_index (obj_index,
575  fe_index,
576  local_index,
577  global_index);
578  }
579 
580 
581  template <int spacedim>
582  static
584  get_dof_index (const ::hp::DoFHandler<3,spacedim> &dof_handler,
585  const unsigned int obj_level,
586  const unsigned int obj_index,
587  const unsigned int fe_index,
588  const unsigned int local_index,
589  const ::internal::int2type<1> &)
590  {
591  return dof_handler.faces->lines.
592  get_dof_index (dof_handler,
593  obj_index,
594  fe_index,
595  local_index,
596  obj_level);
597  }
598 
599 
600  template <int spacedim>
601  static
602  void
603  set_dof_index (const ::hp::DoFHandler<3,spacedim> &dof_handler,
604  const unsigned int obj_level,
605  const unsigned int obj_index,
606  const unsigned int fe_index,
607  const unsigned int local_index,
608  const ::internal::int2type<1> &,
609  const types::global_dof_index global_index)
610  {
611  dof_handler.faces->lines.
612  set_dof_index (dof_handler,
613  obj_index,
614  fe_index,
615  local_index,
616  global_index,
617  obj_level);
618  }
619 
620 
621  template <int spacedim>
622  static
624  get_dof_index (const ::hp::DoFHandler<3,spacedim> &dof_handler,
625  const unsigned int obj_level,
626  const unsigned int obj_index,
627  const unsigned int fe_index,
628  const unsigned int local_index,
629  const ::internal::int2type<2> &)
630  {
631  return dof_handler.faces->quads.
632  get_dof_index (dof_handler,
633  obj_index,
634  fe_index,
635  local_index,
636  obj_level);
637  }
638 
639 
640  template <int spacedim>
641  static
642  void
643  set_dof_index (const ::hp::DoFHandler<3,spacedim> &dof_handler,
644  const unsigned int obj_level,
645  const unsigned int obj_index,
646  const unsigned int fe_index,
647  const unsigned int local_index,
648  const ::internal::int2type<2> &,
649  const types::global_dof_index global_index)
650  {
651  dof_handler.faces->quads.
652  set_dof_index (dof_handler,
653  obj_index,
654  fe_index,
655  local_index,
656  global_index,
657  obj_level);
658  }
659 
660 
661  template <int spacedim>
662  static
664  get_dof_index (const ::hp::DoFHandler<3,spacedim> &dof_handler,
665  const unsigned int obj_level,
666  const unsigned int obj_index,
667  const unsigned int fe_index,
668  const unsigned int local_index,
669  const ::internal::int2type<3> &)
670  {
671  return dof_handler.levels[obj_level]->
672  get_dof_index (obj_index,
673  fe_index,
674  local_index);
675  }
676 
677 
678  template <int spacedim>
679  static
680  void
681  set_dof_index (const ::hp::DoFHandler<3,spacedim> &dof_handler,
682  const unsigned int obj_level,
683  const unsigned int obj_index,
684  const unsigned int fe_index,
685  const unsigned int local_index,
686  const ::internal::int2type<3> &,
687  const types::global_dof_index global_index)
688  {
689  dof_handler.levels[obj_level]->
690  set_dof_index (obj_index,
691  fe_index,
692  local_index,
693  global_index);
694  }
695 
696 
697  template <int structdim, int dim, int spacedim>
698  static
699  bool
700  fe_index_is_active (const ::DoFHandler<dim,spacedim> &,
701  const unsigned int,
702  const unsigned int,
703  const unsigned int fe_index,
704  const ::internal::int2type<structdim> &)
705  {
706  return (fe_index == 0);
707  }
708 
709 
710 
711  template <int structdim, int dim, int spacedim>
712  static
713  unsigned int
714  n_active_fe_indices (const ::DoFHandler<dim,spacedim> &dof_handler,
715  const unsigned int obj_level,
716  const unsigned int obj_index,
717  const ::internal::int2type<structdim> &)
718  {
719  // check that the object we look
720  // at is in fact active. the
721  // problem is that we have
722  // templatized on the
723  // dimensionality of the object,
724  // so it may be a cell, a face,
725  // or a line. we have a bit of
726  // trouble doing this all in the
727  // generic case, so only check if
728  // it is either a cell or a
729  // line. the only case this
730  // leaves out is faces in 3d --
731  // let's hope that this never is
732  // a problem
733  Assert ((dim==structdim
734  ?
735  typename
737  raw_cell_iterator (&dof_handler.get_tria(),
738  obj_level,
739  obj_index)->used()
740  :
741  (structdim==1
742  ?
743  typename
744  internal::Triangulation::Iterators<dim,spacedim>::
745  raw_line_iterator (&dof_handler.get_tria(),
746  obj_level,
747  obj_index)->used()
748  :
749  true))
750  == true,
751  ExcMessage ("This cell is not active and therefore can't be "
752  "queried for its active FE indices"));
753  return 1;
754  }
755 
756 
757 
758  template <int structdim, int dim, int spacedim>
759  static
760  unsigned int
761  nth_active_fe_index (const ::DoFHandler<dim,spacedim> &dof_handler,
762  const unsigned int obj_level,
763  const unsigned int obj_index,
764  const unsigned int n,
765  const ::internal::int2type<structdim> &)
766  {
767  // check that the object we look
768  // at is in fact active. the
769  // problem is that we have
770  // templatized on the
771  // dimensionality of the object,
772  // so it may be a cell, a face,
773  // or a line. we have a bit of
774  // trouble doing this all in the
775  // generic case, so only check if
776  // it is either a cell or a
777  // line. the only case this
778  // leaves out is faces in 3d --
779  // let's hope that this never is
780  // a problem
781  Assert ((dim==structdim
782  ?
783  typename
785  raw_cell_iterator (&dof_handler.get_tria(),
786  obj_level,
787  obj_index)->used()
788  :
789  (structdim==1
790  ?
791  typename
792  internal::Triangulation::Iterators<dim,spacedim>::
793  raw_line_iterator (&dof_handler.get_tria(),
794  obj_level,
795  obj_index)->used()
796  :
797  true))
798  == true,
799  ExcMessage ("This cell is not active and therefore can't be "
800  "queried for its active FE indices"));
801  Assert (n == 0, ExcIndexRange (n, 0, 1));
802 
803  return ::DoFHandler<dim,spacedim>::default_fe_index;
804  }
805 
806 
807  template <int spacedim>
808  static
809  bool
810  fe_index_is_active (const ::hp::DoFHandler<1,spacedim> &dof_handler,
811  const unsigned int obj_level,
812  const unsigned int obj_index,
813  const unsigned int fe_index,
814  const ::internal::int2type<1> &)
815  {
816  return dof_handler.levels[obj_level]->fe_index_is_active(obj_index,
817  fe_index);
818  }
819 
820 
821  template <int spacedim>
822  static
823  unsigned int
824  n_active_fe_indices (const ::hp::DoFHandler<1,spacedim> &,
825  const unsigned int obj_level,
826  const unsigned int obj_index,
827  const ::internal::int2type<1> &)
828  {
829  // on a cell, the number of active elements is one
830  return 1;
831  }
832 
833 
834 
835  template <int spacedim>
836  static
837  unsigned int
838  nth_active_fe_index (const ::hp::DoFHandler<1,spacedim> &dof_handler,
839  const unsigned int obj_level,
840  const unsigned int obj_index,
841  const unsigned int n,
842  const ::internal::int2type<1> &)
843  {
844  Assert (n==0, ExcMessage("On cells, there can only be one active FE index"));
845  return dof_handler.levels[obj_level]->active_fe_index (obj_index);
846  }
847 
848 
849  template <int spacedim>
850  static
851  bool
852  fe_index_is_active (const ::hp::DoFHandler<2,spacedim> &dof_handler,
853  const unsigned int obj_level,
854  const unsigned int obj_index,
855  const unsigned int fe_index,
856  const ::internal::int2type<1> &)
857  {
858  return dof_handler.faces->lines.fe_index_is_active(dof_handler,
859  obj_index,
860  fe_index,
861  obj_level);
862  }
863 
864 
865  template <int spacedim>
866  static
867  unsigned int
868  n_active_fe_indices (const ::hp::DoFHandler<2,spacedim> &dof_handler,
869  const unsigned int ,
870  const unsigned int obj_index,
871  const ::internal::int2type<1> &)
872  {
873  return dof_handler.faces->lines.n_active_fe_indices (dof_handler,
874  obj_index);
875  }
876 
877 
878  template <int spacedim>
879  static
880  unsigned int
881  nth_active_fe_index (const ::hp::DoFHandler<2,spacedim> &dof_handler,
882  const unsigned int obj_level,
883  const unsigned int obj_index,
884  const unsigned int n,
885  const ::internal::int2type<1> &)
886  {
887  return dof_handler.faces->lines.nth_active_fe_index (dof_handler,
888  obj_level,
889  obj_index,
890  n);
891  }
892 
893 
894 
895  template <int spacedim>
896  static
897  bool
898  fe_index_is_active (const ::hp::DoFHandler<2,spacedim> &dof_handler,
899  const unsigned int obj_level,
900  const unsigned int obj_index,
901  const unsigned int fe_index,
902  const ::internal::int2type<2> &)
903  {
904  return dof_handler.levels[obj_level]->fe_index_is_active(obj_index,
905  fe_index);
906  }
907 
908 
909  template <int spacedim>
910  static
911  unsigned int
912  n_active_fe_indices (const ::hp::DoFHandler<2,spacedim> &,
913  const unsigned int obj_level,
914  const unsigned int obj_index,
915  const ::internal::int2type<2> &)
916  {
917  // on a cell, the number of active elements is one
918  return 1;
919  }
920 
921 
922 
923  template <int spacedim>
924  static
925  unsigned int
926  nth_active_fe_index (const ::hp::DoFHandler<2,spacedim> &dof_handler,
927  const unsigned int obj_level,
928  const unsigned int obj_index,
929  const unsigned int n,
930  const ::internal::int2type<2> &)
931  {
932  Assert (n==0, ExcMessage("On cells, there can only be one active FE index"));
933  return dof_handler.levels[obj_level]->active_fe_index (obj_index);
934  }
935 
936 
937 
938  template <int spacedim>
939  static
940  bool
941  fe_index_is_active (const ::hp::DoFHandler<3,spacedim> &dof_handler,
942  const unsigned int obj_level,
943  const unsigned int obj_index,
944  const unsigned int fe_index,
945  const ::internal::int2type<1> &)
946  {
947  return dof_handler.faces->lines.fe_index_is_active(dof_handler,
948  obj_index,
949  fe_index,
950  obj_level);
951  }
952 
953 
954  template <int spacedim>
955  static
956  unsigned int
957  n_active_fe_indices (const ::hp::DoFHandler<3,spacedim> &dof_handler,
958  const unsigned int,
959  const unsigned int obj_index,
960  const ::internal::int2type<1> &)
961  {
962  return dof_handler.faces->lines.n_active_fe_indices (dof_handler,
963  obj_index);
964  }
965 
966 
967 
968  template <int spacedim>
969  static
970  unsigned int
971  nth_active_fe_index (const ::hp::DoFHandler<3,spacedim> &dof_handler,
972  const unsigned int obj_level,
973  const unsigned int obj_index,
974  const unsigned int n,
975  const ::internal::int2type<1> &)
976  {
977  return dof_handler.faces->lines.nth_active_fe_index (dof_handler,
978  obj_level,
979  obj_index,
980  n);
981  }
982 
983 
984 
985  template <int spacedim>
986  static
987  bool
988  fe_index_is_active (const ::hp::DoFHandler<3,spacedim> &dof_handler,
989  const unsigned int obj_level,
990  const unsigned int obj_index,
991  const unsigned int fe_index,
992  const ::internal::int2type<2> &)
993  {
994  return dof_handler.faces->quads.fe_index_is_active(dof_handler,
995  obj_index,
996  fe_index,
997  obj_level);
998  }
999 
1000  template <int spacedim>
1001  static
1002  bool
1003  fe_index_is_active (const ::hp::DoFHandler<3,spacedim> &dof_handler,
1004  const unsigned int obj_level,
1005  const unsigned int obj_index,
1006  const unsigned int fe_index,
1007  const ::internal::int2type<3> &)
1008  {
1009  return dof_handler.levels[obj_level]->fe_index_is_active(obj_index,
1010  fe_index);
1011  }
1012 
1013 
1014  template <int spacedim>
1015  static
1016  unsigned int
1017  n_active_fe_indices (const ::hp::DoFHandler<3,spacedim> &dof_handler,
1018  const unsigned int ,
1019  const unsigned int obj_index,
1020  const ::internal::int2type<2> &)
1021  {
1022  return dof_handler.faces->quads.n_active_fe_indices (dof_handler,
1023  obj_index);
1024  }
1025 
1026 
1027 
1028  template <int spacedim>
1029  static
1030  unsigned int
1031  nth_active_fe_index (const ::hp::DoFHandler<3,spacedim> &dof_handler,
1032  const unsigned int obj_level,
1033  const unsigned int obj_index,
1034  const unsigned int n,
1035  const ::internal::int2type<2> &)
1036  {
1037  return dof_handler.faces->quads.nth_active_fe_index (dof_handler,
1038  obj_level,
1039  obj_index,
1040  n);
1041  }
1042 
1043 
1044 
1045  template <int spacedim>
1046  static
1047  unsigned int
1048  n_active_fe_indices (const ::hp::DoFHandler<3,spacedim> &,
1049  const unsigned int obj_level,
1050  const unsigned int obj_index,
1051  const ::internal::int2type<3> &)
1052  {
1053  // on a cell, the number of active elements is one
1054  return 1;
1055  }
1056 
1057 
1058 
1059  template <int spacedim>
1060  static
1061  unsigned int
1062  nth_active_fe_index (const ::hp::DoFHandler<3,spacedim> &dof_handler,
1063  const unsigned int obj_level,
1064  const unsigned int obj_index,
1065  const unsigned int n,
1066  const ::internal::int2type<3> &)
1067  {
1068  Assert (n==0, ExcMessage("On cells, there can only be one active FE index"));
1069  return dof_handler.levels[obj_level]->active_fe_index (obj_index);
1070  }
1071 
1082  template <int dim, int spacedim>
1083  static
1084  void
1086  const unsigned int vertex_index,
1087  const unsigned int fe_index,
1088  const unsigned int local_index,
1089  const types::global_dof_index global_index)
1090  {
1092  ExcMessage ("Only the default FE index is allowed for non-hp DoFHandler objects"));
1093  Assert (dof_handler.selected_fe != 0,
1094  ExcMessage ("No finite element collection is associated with "
1095  "this DoFHandler"));
1096  Assert (local_index < dof_handler.selected_fe->dofs_per_vertex,
1097  ExcIndexRange(local_index, 0,
1098  dof_handler.selected_fe->dofs_per_vertex));
1099 
1100  dof_handler.vertex_dofs[vertex_index *
1101  dof_handler.selected_fe->dofs_per_vertex
1102  + local_index]
1103  = global_index;
1104  }
1105 
1106 
1107  template <int dim, int spacedim>
1108  static
1109  void
1111  const unsigned int vertex_index,
1112  const unsigned int fe_index,
1113  const unsigned int local_index,
1114  const types::global_dof_index global_index)
1115  {
1117  ExcMessage ("You need to specify a FE index when working "
1118  "with hp DoFHandlers"));
1119  Assert (dof_handler.finite_elements != 0,
1120  ExcMessage ("No finite element collection is associated with "
1121  "this DoFHandler"));
1122  Assert (local_index < (*dof_handler.finite_elements)[fe_index].dofs_per_vertex,
1123  ExcIndexRange(local_index, 0,
1124  (*dof_handler.finite_elements)[fe_index].dofs_per_vertex));
1125  Assert (fe_index < dof_handler.finite_elements->size(),
1126  ExcInternalError());
1127  Assert (dof_handler.vertex_dofs_offsets[vertex_index] !=
1129  ExcMessage ("This vertex is unused and has no DoFs associated with it"));
1130 
1131  // hop along the list of index
1132  // sets until we find the one
1133  // with the correct fe_index, and
1134  // then poke into that
1135  // part. trigger an exception if
1136  // we can't find a set for this
1137  // particular fe_index
1138  const types::global_dof_index starting_offset = dof_handler.vertex_dofs_offsets[vertex_index];
1139  types::global_dof_index *pointer = &dof_handler.vertex_dofs[starting_offset];
1140  while (true)
1141  {
1142  Assert (pointer <= &dof_handler.vertex_dofs.back(), ExcInternalError());
1143 
1144  // a fe index is always small
1145  Assert((*pointer)<std::numeric_limits<unsigned int>::max(), ExcInternalError());
1146  const types::global_dof_index this_fe_index = *pointer;
1147 
1148  Assert (this_fe_index != numbers::invalid_dof_index,
1149  ExcInternalError());
1150  Assert (this_fe_index < dof_handler.finite_elements->size(),
1151  ExcInternalError());
1152 
1153  if (this_fe_index == fe_index)
1154  {
1155  *(pointer + 1 + local_index) = global_index;
1156  return;
1157  }
1158  else
1159  pointer += static_cast<types::global_dof_index>(
1160  (*dof_handler.finite_elements)[this_fe_index].dofs_per_vertex + 1);
1161  }
1162  }
1163 
1164 
1176  template <int dim, int spacedim>
1177  static
1179  get_vertex_dof_index (const ::DoFHandler<dim,spacedim> &dof_handler,
1180  const unsigned int vertex_index,
1181  const unsigned int fe_index,
1182  const unsigned int local_index)
1183  {
1185  ExcMessage ("Only the default FE index is allowed for non-hp DoFHandler objects"));
1186  Assert (dof_handler.selected_fe != 0,
1187  ExcMessage ("No finite element collection is associated with "
1188  "this DoFHandler"));
1189  Assert (local_index < dof_handler.selected_fe->dofs_per_vertex,
1190  ExcIndexRange(local_index, 0,
1191  dof_handler.selected_fe->dofs_per_vertex));
1192 
1193  return
1194  dof_handler.vertex_dofs[vertex_index *
1195  dof_handler.selected_fe->dofs_per_vertex
1196  + local_index];
1197  }
1198 
1199 
1200  template<int dim, int spacedim>
1201  static
1203  get_vertex_dof_index (const ::hp::DoFHandler<dim,spacedim> &dof_handler,
1204  const unsigned int vertex_index,
1205  const unsigned int fe_index,
1206  const unsigned int local_index)
1207  {
1209  ExcMessage ("You need to specify a FE index when working "
1210  "with hp DoFHandlers"));
1211  Assert (dof_handler.finite_elements != 0,
1212  ExcMessage ("No finite element collection is associated with "
1213  "this DoFHandler"));
1214  Assert (local_index < (*dof_handler.finite_elements)[fe_index].dofs_per_vertex,
1215  ExcIndexRange(local_index, 0,
1216  (*dof_handler.finite_elements)[fe_index].dofs_per_vertex));
1217  Assert (vertex_index < dof_handler.vertex_dofs_offsets.size(),
1218  ExcIndexRange (vertex_index, 0,
1219  dof_handler.vertex_dofs_offsets.size()));
1220  Assert (dof_handler.vertex_dofs_offsets[vertex_index] !=
1222  ExcMessage ("This vertex is unused and has no DoFs associated with it"));
1223 
1224  // hop along the list of index
1225  // sets until we find the one
1226  // with the correct fe_index, and
1227  // then poke into that
1228  // part. trigger an exception if
1229  // we can't find a set for this
1230  // particular fe_index
1231  const types::global_dof_index starting_offset = dof_handler.vertex_dofs_offsets[vertex_index];
1232  const types::global_dof_index *pointer = &dof_handler.vertex_dofs[starting_offset];
1233  while (true)
1234  {
1235  Assert (pointer <= &dof_handler.vertex_dofs.back(), ExcInternalError());
1236 
1237  Assert((*pointer)<std::numeric_limits<types::global_dof_index>::max(), ExcInternalError());
1238  const types::global_dof_index this_fe_index = *pointer;
1239 
1240  Assert (this_fe_index != numbers::invalid_dof_index,
1241  ExcInternalError());
1242  Assert (this_fe_index < dof_handler.finite_elements->size(),
1243  ExcInternalError());
1244 
1245  if (this_fe_index == fe_index)
1246  return *(pointer + 1 + local_index);
1247  else
1248  pointer += static_cast<types::global_dof_index>(
1249  (*dof_handler.finite_elements)[this_fe_index].dofs_per_vertex + 1);
1250  }
1251  }
1252 
1253 
1260  template<int dim, int spacedim>
1261  static
1262  unsigned int
1263  n_active_vertex_fe_indices (const ::hp::DoFHandler<dim,spacedim> &dof_handler,
1264  const unsigned int vertex_index)
1265  {
1266  Assert (dof_handler.finite_elements != 0,
1267  ExcMessage ("No finite element collection is associated with "
1268  "this DoFHandler"));
1269 
1270  // if this vertex is unused, return 0
1271  if (dof_handler.vertex_dofs_offsets[vertex_index] == numbers::invalid_dof_index)
1272  return 0;
1273 
1274  // hop along the list of index
1275  // sets and count the number of
1276  // hops
1277  const types::global_dof_index starting_offset = dof_handler.vertex_dofs_offsets[vertex_index];
1278  const types::global_dof_index *pointer = &dof_handler.vertex_dofs[starting_offset];
1279 
1280  Assert (*pointer != numbers::invalid_dof_index,
1281  ExcInternalError());
1282 
1283  unsigned int counter = 0;
1284  while (true)
1285  {
1286  Assert (pointer <= &dof_handler.vertex_dofs.back(), ExcInternalError());
1287 
1288  const types::global_dof_index this_fe_index = *pointer;
1289 
1290  if (this_fe_index == numbers::invalid_dof_index)
1291  return counter;
1292  else
1293  {
1294  pointer += static_cast<types::global_dof_index>(
1295  (*dof_handler.finite_elements)[this_fe_index].dofs_per_vertex + 1);
1296  ++counter;
1297  }
1298  }
1299  }
1300 
1301 
1302 
1308  template<int dim, int spacedim>
1309  static
1310  unsigned int
1311  nth_active_vertex_fe_index (const ::hp::DoFHandler<dim,spacedim> &dof_handler,
1312  const unsigned int vertex_index,
1313  const unsigned int n)
1314  {
1315  Assert (dof_handler.finite_elements != 0,
1316  ExcMessage ("No finite element collection is associated with "
1317  "this DoFHandler"));
1318  Assert (n < n_active_vertex_fe_indices(dof_handler, vertex_index),
1319  ExcIndexRange (n, 0, n_active_vertex_fe_indices(dof_handler,
1320  vertex_index)));
1321  // make sure we don't ask on
1322  // unused vertices
1323  Assert (dof_handler.vertex_dofs_offsets[vertex_index] !=
1325  ExcInternalError());
1326 
1327  // hop along the list of index
1328  // sets and count the number of
1329  // hops
1330  const types::global_dof_index starting_offset = dof_handler.vertex_dofs_offsets[vertex_index];
1331  const types::global_dof_index *pointer = &dof_handler.vertex_dofs[starting_offset];
1332 
1333  Assert (*pointer != numbers::invalid_dof_index,
1334  ExcInternalError());
1335 
1336  unsigned int counter = 0;
1337  while (true)
1338  {
1339  Assert (pointer <= &dof_handler.vertex_dofs.back(), ExcInternalError());
1340 
1341  Assert((*pointer)<std::numeric_limits<unsigned int>::max(), ExcInternalError());
1342  const types::global_dof_index this_fe_index = *pointer;
1343 
1344  Assert (this_fe_index < dof_handler.finite_elements->size(),
1345  ExcInternalError());
1346 
1347  if (counter == n)
1348  return this_fe_index;
1349 
1350  Assert (this_fe_index != numbers::invalid_dof_index,
1351  ExcInternalError());
1352 
1353  pointer += static_cast<types::global_dof_index>(
1354  (*dof_handler.finite_elements)[this_fe_index].dofs_per_vertex + 1);
1355  ++counter;
1356  }
1357  }
1358 
1359 
1360 
1367  template<int dim, int spacedim>
1368  static
1369  bool
1370  fe_is_active_on_vertex (const ::hp::DoFHandler<dim,spacedim> &dof_handler,
1371  const unsigned int vertex_index,
1372  const unsigned int fe_index)
1373  {
1375  ExcMessage ("You need to specify a FE index when working "
1376  "with hp DoFHandlers"));
1377  Assert (dof_handler.finite_elements != 0,
1378  ExcMessage ("No finite element collection is associated with "
1379  "this DoFHandler"));
1380  Assert (fe_index < dof_handler.finite_elements->size(),
1381  ExcInternalError());
1382 
1383  // make sure we don't ask on
1384  // unused vertices
1385  Assert (dof_handler.vertex_dofs_offsets[vertex_index] !=
1387  ExcInternalError());
1388 
1389  // hop along the list of index
1390  // sets and see whether we find
1391  // the given index
1392  const types::global_dof_index starting_offset = dof_handler.vertex_dofs_offsets[vertex_index];
1393  const types::global_dof_index *pointer = &dof_handler.vertex_dofs[starting_offset];
1394 
1395  Assert (*pointer != numbers::invalid_dof_index,
1396  ExcInternalError());
1397 
1398  while (true)
1399  {
1400  Assert (pointer <= &dof_handler.vertex_dofs.back(), ExcInternalError());
1401 
1402  Assert((*pointer)<std::numeric_limits<types::global_dof_index>::max(), ExcInternalError());
1403  const types::global_dof_index this_fe_index = *pointer;
1404 
1405  Assert (this_fe_index < dof_handler.finite_elements->size(),
1406  ExcInternalError());
1407 
1408  if (this_fe_index == numbers::invalid_dof_index)
1409  return false;
1410  else if (this_fe_index == fe_index)
1411  return true;
1412  else
1413  pointer += (*dof_handler.finite_elements)[this_fe_index].dofs_per_vertex + 1;
1414  }
1415  }
1416 
1417  template<class DH, bool lda>
1418  static
1419  void set_mg_dof_indices (const ::DoFAccessor<1,DH,lda> &,
1420  const int,
1421  const std::vector<types::global_dof_index> &,
1422  const unsigned int)
1423  {
1424  AssertThrow (false, ExcNotImplemented ()); //TODO[TH]: implement
1425  }
1426 
1427 
1428 
1429  template<class DH, bool lda>
1430  static
1431  void set_mg_dof_indices (::DoFAccessor<2, DH,lda> &accessor,
1432  const int level,
1433  const std::vector<types::global_dof_index> &dof_indices,
1434  const unsigned int fe_index)
1435  {
1436  const FiniteElement<DH::dimension, DH::space_dimension> &fe = accessor.get_dof_handler ().get_fe ()[fe_index];
1437  std::vector<types::global_dof_index>::const_iterator next = dof_indices.begin ();
1438 
1439  for (unsigned int vertex = 0; vertex < GeometryInfo<2>::vertices_per_cell; ++vertex)
1440  for (unsigned int dof = 0; dof < fe.dofs_per_vertex; ++dof)
1441  accessor.set_mg_vertex_dof_index(level, vertex, dof, *next++, fe_index);
1442 
1443  for (unsigned int line = 0; line < GeometryInfo<2>::lines_per_cell; ++line)
1444  for (unsigned int dof = 0; dof < fe.dofs_per_line; ++dof)
1445  accessor.line(line)->set_mg_dof_index(level, dof, *next++);
1446 
1447  for (unsigned int dof = 0; dof < fe.dofs_per_quad; ++dof)
1448  accessor.set_mg_dof_index(level, dof, *next++);
1449 
1450  Assert (next == dof_indices.end (), ExcInternalError ());
1451  }
1452 
1453 
1454 
1455  template<class DH, bool lda>
1456  static
1457  void set_mg_dof_indices (const ::DoFAccessor<3, DH,lda> &,
1458  const int,
1459  const std::vector<types::global_dof_index> &,
1460  const unsigned int)
1461  {
1462  AssertThrow (false, ExcNotImplemented ()); //TODO[TH]: implement
1463  }
1464 
1465  };
1466  }
1467 }
1468 
1469 
1470 
1471 template <int dim, class DH, bool lda>
1472 inline
1475  const unsigned int fe_index) const
1476 {
1477  // access the respective DoF
1478  return ::internal::DoFAccessor::Implementation::get_dof_index (*this->dof_handler,
1479  this->level(),
1480  this->present_index,
1481  fe_index,
1482  i,
1484 }
1485 
1486 
1487 template<int structdim, class DH, bool lda>
1488 inline
1491  const unsigned int i) const
1492 {
1493  return this->dof_handler->template get_dof_index<structdim> (level, this->present_index, 0, i);
1494 }
1495 
1496 
1497 template <int dim, class DH, bool lda>
1498 inline
1499 void
1501  const types::global_dof_index index,
1502  const unsigned int fe_index) const
1503 {
1504  // access the respective DoF
1505  ::internal::DoFAccessor::Implementation::set_dof_index (*this->dof_handler,
1506  this->level(),
1507  this->present_index,
1508  fe_index,
1509  i,
1511  index);
1512 }
1513 
1514 
1515 
1516 template <int dim, class DH, bool lda>
1517 inline
1518 unsigned int
1520 {
1521  // access the respective DoF
1522  return
1523  ::internal::DoFAccessor::Implementation::
1524  n_active_fe_indices (*this->dof_handler,
1525  this->level(),
1526  this->present_index,
1528 }
1529 
1530 
1531 
1532 template <int dim, class DH, bool lda>
1533 inline
1534 unsigned int
1536 {
1537  // access the respective DoF
1538  return
1539  ::internal::DoFAccessor::Implementation::
1540  nth_active_fe_index (*this->dof_handler,
1541  this->level(),
1542  this->present_index,
1543  n,
1545 }
1546 
1547 
1548 
1549 template <int dim, class DH, bool lda>
1550 inline
1551 bool
1552 DoFAccessor<dim,DH,lda>::fe_index_is_active (const unsigned int fe_index) const
1553 {
1554  // access the respective DoF
1555  return
1556  ::internal::DoFAccessor::Implementation::
1557  fe_index_is_active (*this->dof_handler,
1558  this->level(),
1559  this->present_index,
1560  fe_index,
1562 }
1563 
1564 
1565 
1566 template <int structdim, class DH, bool lda>
1567 inline
1570  const unsigned int i,
1571  const unsigned int fe_index) const
1572 {
1573  return
1574  ::internal::DoFAccessor::Implementation::get_vertex_dof_index
1575  (*this->dof_handler,
1576  this->vertex_index(vertex),
1577  fe_index,
1578  i);
1579 }
1580 
1581 
1582 template<int structdim, class DH, bool lda>
1583 inline
1586  const unsigned int vertex,
1587  const unsigned int i,
1588  const unsigned int fe_index) const
1589 {
1590  Assert (this->dof_handler != 0, ExcInvalidObject ());
1591  Assert (&this->dof_handler->get_fe () != 0, ExcInvalidObject ());
1593  Assert (i < this->dof_handler->get_fe ()[fe_index].dofs_per_vertex, ExcIndexRange (i, 0, this->dof_handler->get_fe ()[fe_index].dofs_per_vertex));
1594  return this->dof_handler->mg_vertex_dofs[this->vertex_index (vertex)].get_index (level, i);
1595 }
1596 
1597 
1598 template <int structdim, class DH, bool lda>
1599 inline
1600 void
1602  const unsigned int i,
1603  const types::global_dof_index index,
1604  const unsigned int fe_index) const
1605 {
1607  (*this->dof_handler,
1608  this->vertex_index(vertex),
1609  fe_index,
1610  i,
1611  index);
1612 }
1613 
1614 
1615 template<int structdim, class DH, bool lda>
1616 inline
1617 void
1619  const unsigned int vertex,
1620  const unsigned int i,
1621  const types::global_dof_index index,
1622  const unsigned int fe_index) const
1623 {
1624  Assert (this->dof_handler != 0, ExcInvalidObject ());
1625  Assert (&this->dof_handler->get_fe () != 0, ExcInvalidObject ());
1627  Assert (i < this->dof_handler->get_fe ()[fe_index].dofs_per_vertex, ExcIndexRange (i, 0, this->dof_handler->get_fe ()[fe_index].dofs_per_vertex));
1628  this->dof_handler->mg_vertex_dofs[this->vertex_index (vertex)].set_index (level, i, index);
1629 }
1630 
1631 
1632 template<int structdim, class DH, bool lda>
1633 inline
1634 void
1636  const unsigned int i,
1637  const types::global_dof_index index) const
1638 {
1639  this->dof_handler->template set_dof_index<structdim> (level, this->present_index, 0, i, index);
1640 }
1641 
1642 
1643 namespace internal
1644 {
1645  namespace DoFAccessor
1646  {
1647  template <int dim, int spacedim>
1648  inline
1650  get_fe (const FiniteElement<dim,spacedim> &fe,
1651  const unsigned int)
1652  {
1653  return fe;
1654  }
1655 
1656 
1657 
1658  template <int dim, int spacedim>
1659  inline
1661  get_fe (const ::hp::FECollection<dim,spacedim> &fe,
1662  const unsigned int index)
1663  {
1664  return fe[index];
1665  }
1666  }
1667 }
1668 
1669 
1670 template <int dim, class DH, bool lda>
1671 inline
1673 DoFAccessor<dim,DH,lda>::get_fe (const unsigned int fe_index) const
1674 {
1675  Assert (fe_index_is_active (fe_index) == true,
1676  ExcMessage ("This function can only be called for active fe indices"));
1677 
1678  return ::internal::DoFAccessor::get_fe (this->dof_handler->get_fe(), fe_index);
1679 }
1680 
1681 
1682 
1683 namespace internal
1684 {
1685  namespace DoFAccessor
1686  {
1687  template <class DH, bool lda>
1688  void get_dof_indices (const ::DoFAccessor<1,DH,lda> &accessor,
1689  std::vector<types::global_dof_index> &dof_indices,
1690  const unsigned int fe_index)
1691  {
1692  const unsigned int dofs_per_vertex = accessor.get_fe(fe_index).dofs_per_vertex,
1693  dofs_per_line = accessor.get_fe(fe_index).dofs_per_line;
1694  std::vector<types::global_dof_index>::iterator next = dof_indices.begin();
1695  for (unsigned int vertex=0; vertex<2; ++vertex)
1696  for (unsigned int d=0; d<dofs_per_vertex; ++d)
1697  *next++ = accessor.vertex_dof_index(vertex,d,fe_index);
1698  for (unsigned int d=0; d<dofs_per_line; ++d)
1699  *next++ = accessor.dof_index(d,fe_index);
1700  }
1701 
1702 
1703 
1704  template <class DH, bool lda>
1705  void get_dof_indices (const ::DoFAccessor<2,DH,lda> &accessor,
1706  std::vector<types::global_dof_index> &dof_indices,
1707  const unsigned int fe_index)
1708  {
1709  const unsigned int dofs_per_vertex = accessor.get_fe(fe_index).dofs_per_vertex,
1710  dofs_per_line = accessor.get_fe(fe_index).dofs_per_line,
1711  dofs_per_quad = accessor.get_fe(fe_index).dofs_per_quad;
1712  std::vector<types::global_dof_index>::iterator next = dof_indices.begin();
1713  for (unsigned int vertex=0; vertex<4; ++vertex)
1714  for (unsigned int d=0; d<dofs_per_vertex; ++d)
1715  *next++ = accessor.vertex_dof_index(vertex,d,fe_index);
1716  // now copy dof numbers from the line. for
1717  // lines with the wrong orientation (which
1718  // might occur in 3d), we have already made
1719  // sure that we're ok by picking the correct
1720  // vertices (this happens automatically in
1721  // the vertex() function). however, if the
1722  // line is in wrong orientation, we look at
1723  // it in flipped orientation and we will have
1724  // to adjust the shape function indices that
1725  // we see to correspond to the correct
1726  // (face-local) ordering.
1727  for (unsigned int line=0; line<4; ++line)
1728  for (unsigned int d=0; d<dofs_per_line; ++d)
1729  *next++ = accessor.line(line)->dof_index(accessor.get_fe(fe_index).
1730  adjust_line_dof_index_for_line_orientation(d,
1731  accessor.line_orientation(line)),
1732  fe_index);
1733  for (unsigned int d=0; d<dofs_per_quad; ++d)
1734  *next++ = accessor.dof_index(d,fe_index);
1735  }
1736 
1737 
1738 
1739  template <class DH, bool lda>
1740  void get_dof_indices (const ::DoFAccessor<3,DH,lda> &accessor,
1741  std::vector<types::global_dof_index> &dof_indices,
1742  const unsigned int fe_index)
1743  {
1744  const unsigned int dofs_per_vertex = accessor.get_fe(fe_index).dofs_per_vertex,
1745  dofs_per_line = accessor.get_fe(fe_index).dofs_per_line,
1746  dofs_per_quad = accessor.get_fe(fe_index).dofs_per_quad,
1747  dofs_per_hex = accessor.get_fe(fe_index).dofs_per_hex;
1748  std::vector<types::global_dof_index>::iterator next = dof_indices.begin();
1749  for (unsigned int vertex=0; vertex<8; ++vertex)
1750  for (unsigned int d=0; d<dofs_per_vertex; ++d)
1751  *next++ = accessor.vertex_dof_index(vertex,d,fe_index);
1752  // now copy dof numbers from the line. for
1753  // lines with the wrong orientation, we have
1754  // already made sure that we're ok by picking
1755  // the correct vertices (this happens
1756  // automatically in the vertex()
1757  // function). however, if the line is in
1758  // wrong orientation, we look at it in
1759  // flipped orientation and we will have to
1760  // adjust the shape function indices that we
1761  // see to correspond to the correct
1762  // (cell-local) ordering.
1763  for (unsigned int line=0; line<12; ++line)
1764  for (unsigned int d=0; d<dofs_per_line; ++d)
1765  *next++ = accessor.line(line)->dof_index(accessor.get_fe(fe_index).
1766  adjust_line_dof_index_for_line_orientation(d,
1767  accessor.line_orientation(line)),fe_index);
1768  // now copy dof numbers from the face. for
1769  // faces with the wrong orientation, we
1770  // have already made sure that we're ok by
1771  // picking the correct lines and vertices
1772  // (this happens automatically in the
1773  // line() and vertex() functions). however,
1774  // if the face is in wrong orientation, we
1775  // look at it in flipped orientation and we
1776  // will have to adjust the shape function
1777  // indices that we see to correspond to the
1778  // correct (cell-local) ordering. The same
1779  // applies, if the face_rotation or
1780  // face_orientation is non-standard
1781  for (unsigned int quad=0; quad<6; ++quad)
1782  for (unsigned int d=0; d<dofs_per_quad; ++d)
1783  *next++ = accessor.quad(quad)->dof_index(accessor.get_fe(fe_index).
1784  adjust_quad_dof_index_for_face_orientation(d,
1785  accessor.face_orientation(quad),
1786  accessor.face_flip(quad),
1787  accessor.face_rotation(quad)),
1788  fe_index);
1789  for (unsigned int d=0; d<dofs_per_hex; ++d)
1790  *next++ = accessor.dof_index(d,fe_index);
1791  }
1792 
1793 
1794 
1795  template<class DH, bool lda>
1796  void get_mg_dof_indices (const ::DoFAccessor<1, DH,lda> &accessor,
1797  const int level,
1798  std::vector<types::global_dof_index> &dof_indices,
1799  const unsigned int fe_index)
1800  {
1801  const DH &handler = accessor.get_dof_handler();
1802  Assert(handler.n_dofs(level) != numbers::invalid_dof_index,
1803  ExcNotInitialized());
1804 
1806  = handler.get_fe ()[fe_index];
1807  std::vector<types::global_dof_index>::iterator next = dof_indices.begin ();
1808 
1809  for (unsigned int vertex = 0; vertex < GeometryInfo<1>::vertices_per_cell; ++vertex)
1810  for (unsigned int dof = 0; dof < fe.dofs_per_vertex; ++dof)
1811  *next++ = accessor.mg_vertex_dof_index (level, vertex, dof);
1812 
1813  for (unsigned int dof = 0; dof < fe.dofs_per_line; ++dof)
1814  *next++ = accessor.mg_dof_index (level, dof);
1815 
1816  Assert (next == dof_indices.end (), ExcInternalError ());
1817  }
1818 
1819 
1820 
1821  template<class DH, bool lda>
1822  void get_mg_dof_indices (const ::DoFAccessor<2, DH,lda> &accessor,
1823  const int level,
1824  std::vector<types::global_dof_index> &dof_indices,
1825  const unsigned int fe_index)
1826  {
1827  const DH &handler = accessor.get_dof_handler();
1828  Assert(handler.n_dofs(level) != numbers::invalid_dof_index,
1829  ExcNotInitialized());
1830 
1831  const FiniteElement<DH::dimension, DH::space_dimension> &fe = handler.get_fe ()[fe_index];
1832  std::vector<types::global_dof_index>::iterator next = dof_indices.begin ();
1833 
1834  for (unsigned int vertex = 0; vertex < GeometryInfo<2>::vertices_per_cell; ++vertex)
1835  for (unsigned int dof = 0; dof < fe.dofs_per_vertex; ++dof)
1836  *next++ = accessor.mg_vertex_dof_index (level, vertex, dof);
1837 
1838  for (unsigned int line = 0; line < GeometryInfo<2>::lines_per_cell; ++line)
1839  for (unsigned int dof = 0; dof < fe.dofs_per_line; ++dof)
1840  *next++ = accessor.line (line)->mg_dof_index (level, dof);
1841 
1842  for (unsigned int dof = 0; dof < fe.dofs_per_quad; ++dof)
1843  *next++ = accessor.mg_dof_index (level, dof);
1844 
1845  Assert (next == dof_indices.end (), ExcInternalError ());
1846  }
1847 
1848 
1849 
1850  template<class DH, bool lda>
1851  void get_mg_dof_indices (const ::DoFAccessor<3, DH,lda> &accessor,
1852  const int level,
1853  std::vector<types::global_dof_index> &dof_indices,
1854  const unsigned int fe_index)
1855  {
1856  const DH &handler = accessor.get_dof_handler();
1857  Assert(handler.n_dofs(level) != numbers::invalid_dof_index,
1858  ExcNotInitialized());
1859 
1860  const FiniteElement<DH::dimension, DH::space_dimension> &fe = handler.get_fe ()[fe_index];
1861  std::vector<types::global_dof_index>::iterator next = dof_indices.begin ();
1862 
1863  for (unsigned int vertex = 0; vertex < GeometryInfo<3>::vertices_per_cell; ++vertex)
1864  for (unsigned int dof = 0; dof < fe.dofs_per_vertex; ++dof)
1865  *next++ = accessor.mg_vertex_dof_index (level, vertex, dof);
1866 
1867  for (unsigned int line = 0; line < GeometryInfo<3>::lines_per_cell; ++line)
1868  for (unsigned int dof = 0; dof < fe.dofs_per_line; ++dof)
1869  *next++ = accessor.line (line)->mg_dof_index (level, dof);
1870 
1871  for (unsigned int quad = 0; quad < GeometryInfo<3>::quads_per_cell; ++quad)
1872  for (unsigned int dof = 0; dof < fe.dofs_per_quad; ++dof)
1873  *next++ = accessor.quad (quad)->mg_dof_index (level, dof);
1874 
1875  for (unsigned int dof = 0; dof < fe.dofs_per_hex; ++dof)
1876  *next++ = accessor.mg_dof_index (level, dof);
1877 
1878  Assert (next == dof_indices.end (), ExcInternalError ());
1879  }
1880 
1881 
1882  }
1883 }
1884 
1885 
1886 template <int structdim, class DH, bool lda>
1887 inline
1888 void
1889 DoFAccessor<structdim,DH,lda>::get_dof_indices (std::vector<types::global_dof_index> &dof_indices,
1890  const unsigned int fe_index) const
1891 {
1892  Assert (this->dof_handler != 0, ExcNotInitialized());
1893  Assert (&this->dof_handler->get_fe() != 0, ExcMessage ("DoFHandler not initialized"));
1894  Assert (static_cast<unsigned int>(this->level()) < this->dof_handler->levels.size(),
1895  ExcMessage ("DoFHandler not initialized"));
1896 
1897  switch (structdim)
1898  {
1899  case 1:
1900  Assert (dof_indices.size() ==
1901  (2*this->dof_handler->get_fe()[fe_index].dofs_per_vertex +
1902  this->dof_handler->get_fe()[fe_index].dofs_per_line),
1903  ExcVectorDoesNotMatch());
1904  break;
1905  case 2:
1906  Assert (dof_indices.size() ==
1907  (4*this->dof_handler->get_fe()[fe_index].dofs_per_vertex +
1908  4*this->dof_handler->get_fe()[fe_index].dofs_per_line +
1909  this->dof_handler->get_fe()[fe_index].dofs_per_quad),
1910  ExcVectorDoesNotMatch());
1911  break;
1912  case 3:
1913  Assert (dof_indices.size() ==
1914  (8*this->dof_handler->get_fe()[fe_index].dofs_per_vertex +
1915  12*this->dof_handler->get_fe()[fe_index].dofs_per_line +
1916  6*this->dof_handler->get_fe()[fe_index].dofs_per_quad +
1917  this->dof_handler->get_fe()[fe_index].dofs_per_hex),
1918  ExcVectorDoesNotMatch());
1919  break;
1920  default:
1921  Assert (false, ExcNotImplemented());
1922  }
1923 
1924 
1925  // this function really only makes
1926  // sense if either a) there are
1927  // degrees of freedom defined on
1928  // the present object, or b) the
1929  // object is non-active objects but
1930  // all degrees of freedom are
1931  // located on vertices, since
1932  // otherwise there are degrees of
1933  // freedom on sub-objects which are
1934  // not allocated for this
1935  // non-active thing
1936  Assert (this->fe_index_is_active (fe_index)
1937  ||
1938  (this->dof_handler->get_fe()[fe_index].dofs_per_cell ==
1940  this->dof_handler->get_fe()[fe_index].dofs_per_vertex),
1941  ExcInternalError());
1942 
1943  // now do the actual work
1944  ::internal::DoFAccessor::get_dof_indices (*this, dof_indices, fe_index);
1945 }
1946 
1947 
1948 
1949 template<int structdim, class DH, bool lda>
1950 inline
1952  std::vector<types::global_dof_index> &dof_indices,
1953  const unsigned int fe_index) const
1954 {
1955  Assert (this->dof_handler != 0, ExcInvalidObject ());
1956  Assert (&this->dof_handler->get_fe () != 0, ExcInvalidObject ());
1957 
1958  switch (structdim)
1959  {
1960  case 1:
1961  {
1962  Assert (dof_indices.size () ==
1963  2 * this->dof_handler->get_fe ()[fe_index].dofs_per_vertex +
1964  this->dof_handler->get_fe ()[fe_index].dofs_per_line,
1965  ExcVectorDoesNotMatch ());
1966  break;
1967  }
1968 
1969  case 2:
1970  {
1971  Assert (dof_indices.size () ==
1972  4 * (this->dof_handler->get_fe ()[fe_index].dofs_per_vertex +
1973  this->dof_handler->get_fe ()[fe_index].dofs_per_line) +
1974  this->dof_handler->get_fe ()[fe_index].dofs_per_quad,
1975  ExcVectorDoesNotMatch ());
1976  break;
1977  }
1978 
1979  case 3:
1980  {
1981  Assert (dof_indices.size () ==
1982  8 * this->dof_handler->get_fe ()[fe_index].dofs_per_vertex +
1983  12 * this->dof_handler->get_fe ()[fe_index].dofs_per_line +
1984  6 * this->dof_handler->get_fe ()[fe_index].dofs_per_quad +
1985  this->dof_handler->get_fe ()[fe_index].dofs_per_hex,
1986  ExcVectorDoesNotMatch ());
1987  break;
1988  }
1989 
1990  default:
1991  Assert (false, ExcNotImplemented ());
1992  }
1993 
1994  internal::DoFAccessor::get_mg_dof_indices (*this,
1995  level,
1996  dof_indices,
1997  fe_index);
1998 }
1999 
2000 
2001 template<int structdim, class DH, bool lda>
2002 inline
2004  const std::vector<types::global_dof_index> &dof_indices,
2005  const unsigned int fe_index)
2006 {
2007  Assert (this->dof_handler != 0, ExcInvalidObject ());
2008  Assert (&this->dof_handler->get_fe () != 0, ExcInvalidObject ());
2009 
2010  switch (structdim)
2011  {
2012  case 1:
2013  {
2014  Assert (dof_indices.size () ==
2015  2 * this->dof_handler->get_fe ()[fe_index].dofs_per_vertex +
2016  this->dof_handler->get_fe ()[fe_index].dofs_per_line,
2017  ExcVectorDoesNotMatch ());
2018  break;
2019  }
2020 
2021  case 2:
2022  {
2023  Assert (dof_indices.size () ==
2024  4 * (this->dof_handler->get_fe ()[fe_index].dofs_per_vertex +
2025  this->dof_handler->get_fe ()[fe_index].dofs_per_line) +
2026  this->dof_handler->get_fe ()[fe_index].dofs_per_quad,
2027  ExcVectorDoesNotMatch ());
2028  break;
2029  }
2030 
2031  case 3:
2032  {
2033  Assert (dof_indices.size () ==
2034  8 * this->dof_handler->get_fe ()[fe_index].dofs_per_vertex +
2035  12 * this->dof_handler->get_fe ()[fe_index].dofs_per_line +
2036  6 * this->dof_handler->get_fe ()[fe_index].dofs_per_quad +
2037  this->dof_handler->get_fe ()[fe_index].dofs_per_hex,
2038  ExcVectorDoesNotMatch ());
2039  break;
2040  }
2041 
2042  default:
2043  Assert (false, ExcNotImplemented ());
2044  }
2045 
2046  internal::DoFAccessor::Implementation::set_mg_dof_indices (*this,
2047  level,
2048  dof_indices,
2049  fe_index);
2050 }
2051 
2052 
2053 namespace internal
2054 {
2055  namespace DoFAccessor
2056  {
2057  template <bool lda, class DH>
2058  inline
2059  typename ::internal::DoFHandler::Iterators<DH, lda>::quad_iterator
2060  get_quad(const ::Triangulation<DH::dimension, DH::space_dimension> *tria,
2061  unsigned int index,
2062  DH *dof_handler)
2063  {
2064  }
2065 
2066 
2067  template<bool lda>
2068  inline
2069  typename ::internal::DoFHandler::Iterators<::DoFHandler<2,2>, lda>::quad_iterator
2070  get_quad(const ::Triangulation<2,2> *,
2071  unsigned int,
2072  ::DoFHandler<2,2> *)
2073  {
2074  Assert(false, ExcNotImplemented());
2075  return typename ::internal::DoFHandler::Iterators<::DoFHandler<2,2>, lda>::line_iterator();
2076  }
2077 
2078  template<bool lda>
2079  inline
2080  typename ::internal::DoFHandler::Iterators<::DoFHandler<2,3>, lda>::quad_iterator
2081  get_quad(const ::Triangulation<2,3> *,
2082  unsigned int,
2083  ::DoFHandler<2,3> *)
2084  {
2085  Assert(false, ExcNotImplemented());
2086  return typename ::internal::DoFHandler::Iterators<::DoFHandler<2,3>, lda>::line_iterator();
2087  }
2088 
2089  template<bool lda>
2090  inline
2091  typename ::internal::DoFHandler::Iterators<::hp::DoFHandler<2,2>, lda>::quad_iterator
2092  get_quad(const ::Triangulation<2,2> *,
2093  unsigned int,
2094  ::hp::DoFHandler<2,2> *)
2095  {
2096  Assert(false, ExcNotImplemented());
2097  return typename ::internal::DoFHandler::Iterators<::hp::DoFHandler<2,2>, lda>::line_iterator();
2098  }
2099 
2100  template<bool lda>
2101  inline
2102  typename ::internal::DoFHandler::Iterators<::hp::DoFHandler<2,3>, lda>::quad_iterator
2103  get_quad(const ::Triangulation<2,3> *,
2104  unsigned int,
2105  ::hp::DoFHandler<2,3> *)
2106  {
2107  Assert(false, ExcNotImplemented());
2108  return typename ::internal::DoFHandler::Iterators<::hp::DoFHandler<2,3>, lda>::line_iterator();
2109  }
2110  }
2111 }
2112 
2113 
2114 template <int structdim, class DH, bool lda>
2115 inline
2116 typename ::internal::DoFHandler::Iterators<DH,lda>::line_iterator
2117 DoFAccessor<structdim,DH,lda>::line (const unsigned int i) const
2118 {
2119  // if we are asking for a particular line and this object refers to
2120  // a line, then the only valid index is i==0 and we should return
2121  // *this
2122  if (structdim == 1)
2123  {
2124  Assert (i==0, ExcMessage ("You can only ask for line zero if the "
2125  "current object is a line itself."));
2126  return
2127  typename ::internal::DoFHandler::Iterators<DH,lda>::cell_iterator
2128  (&this->get_triangulation(),
2129  this->level(),
2130  this->index(),
2131  &this->get_dof_handler());
2132  }
2133 
2134  // otherwise we need to be in structdim>=2
2135  Assert (structdim > 1, ExcImpossibleInDim(structdim));
2136  Assert (DH::dimension > 1, ExcImpossibleInDim(DH::dimension));
2137 
2138  // checking of 'i' happens in line_index(i)
2139  return typename ::internal::DoFHandler::Iterators<DH,lda>::line_iterator
2140  (this->tria,
2141  0, // only sub-objects are allowed, which have no level
2142  this->line_index(i),
2143  this->dof_handler);
2144 }
2145 
2146 
2147 template <int structdim, class DH, bool lda>
2148 inline
2149 typename ::internal::DoFHandler::Iterators<DH,lda>::quad_iterator
2150 DoFAccessor<structdim,DH,lda>::quad (const unsigned int i) const
2151 {
2152  // if we are asking for a
2153  // particular quad and this object
2154  // refers to a quad, then the only
2155  // valid index is i==0 and we
2156  // should return *this
2157  if (structdim == 2)
2158  {
2159  Assert (i==0, ExcMessage ("You can only ask for quad zero if the "
2160  "current object is a quad itself."));
2161  return
2162  typename ::internal::DoFHandler::Iterators<DH>::cell_iterator
2163  (&this->get_triangulation(),
2164  this->level(),
2165  this->index(),
2166  &this->get_dof_handler());
2167  }
2168 
2169  // otherwise we need to be in structdim>=3
2170  Assert (structdim > 2, ExcImpossibleInDim(structdim));
2171  Assert (DH::dimension > 2, ExcImpossibleInDim(DH::dimension));
2172 
2173  // checking of 'i' happens in quad_index(i)
2174  return typename ::internal::DoFHandler::Iterators<DH,lda>::quad_iterator
2175  (this->tria,
2176  0, // only sub-objects are allowed, which have no level
2177  this->quad_index(i),
2178  this->dof_handler);
2179 }
2180 
2181 
2182 /*------------------------- Functions: DoFAccessor<0,1,spacedim> ---------------------------*/
2183 
2184 
2185 template <template <int, int> class DH, int spacedim, bool lda>
2186 inline
2188 {
2189  Assert (false, ExcInvalidObject());
2190 }
2191 
2192 
2193 
2194 template <template <int, int> class DH, int spacedim, bool lda>
2195 inline
2196 DoFAccessor<0,DH<1,spacedim>, lda>::
2198  const typename TriaAccessor<0,1,spacedim>::VertexKind vertex_kind,
2199  const unsigned int vertex_index,
2200  const DH<1,spacedim> *dof_handler)
2201  :
2202  BaseClass (tria,
2203  vertex_kind,
2204  vertex_index),
2205  dof_handler(const_cast<DH<1,spacedim>*>(dof_handler))
2206 {}
2207 
2208 
2209 
2210 template <template <int, int> class DH, int spacedim, bool lda>
2211 inline
2212 DoFAccessor<0,DH<1,spacedim>, lda>::
2214  const int ,
2215  const int ,
2216  const DH<1,spacedim> *)
2217  :
2218  dof_handler(0)
2219 {
2220  Assert (false,
2221  ExcMessage ("This constructor can not be called for face iterators in 1d."));
2222 }
2223 
2224 
2225 
2226 template <template <int, int> class DH, int spacedim, bool lda>
2227 template <int structdim2, int dim2, int spacedim2>
2228 DoFAccessor<0,DH<1,spacedim>, lda>::DoFAccessor (const InvalidAccessor<structdim2,dim2,spacedim2> &)
2229 {
2230  Assert (false, ExcInvalidObject());
2231 }
2232 
2233 
2234 
2235 template <template <int, int> class DH, int spacedim, bool lda>
2236 template <int dim2, class DH2, bool lda2>
2237 inline
2238 DoFAccessor<0,DH<1,spacedim>, lda>::DoFAccessor (const DoFAccessor<dim2, DH2, lda2> &)
2239 {
2240  Assert (false, ExcInvalidObject());
2241 }
2242 
2243 
2244 
2245 template <template <int, int> class DH, int spacedim, bool lda>
2246 inline
2247 void
2248 DoFAccessor<0,DH<1,spacedim>, lda>::set_dof_handler (DH<1,spacedim> *dh)
2249 {
2250  Assert (dh != 0, ExcInvalidObject());
2251  this->dof_handler = dh;
2252 }
2253 
2254 
2255 
2256 template <template <int, int> class DH, int spacedim, bool lda>
2257 inline
2258 const DH<1,spacedim> &
2259 DoFAccessor<0,DH<1,spacedim>, lda>::get_dof_handler () const
2260 {
2261  return *this->dof_handler;
2262 }
2263 
2264 
2265 
2266 template <template <int, int> class DH, int spacedim, bool lda>
2267 inline
2268 void
2269 DoFAccessor<0,DH<1,spacedim>, lda>::get_dof_indices (
2270  std::vector<types::global_dof_index> &dof_indices,
2271  const unsigned int fe_index) const
2272 {
2273  for (unsigned int i=0; i<dof_indices.size(); ++i)
2274  dof_indices[i]
2276  *dof_handler,
2277  this->global_vertex_index,
2278  fe_index,
2279  i);
2280 }
2281 
2282 
2283 
2284 template <template <int, int> class DH, int spacedim, bool lda>
2285 inline
2287 DoFAccessor<0,DH<1,spacedim>, lda>::
2288 vertex_dof_index (const unsigned int vertex,
2289  const unsigned int i,
2290  const unsigned int fe_index) const
2291 {
2292  Assert (vertex == 0, ExcIndexRange (vertex, 0, 1));
2293  return ::internal::DoFAccessor::Implementation::get_vertex_dof_index (
2294  *dof_handler,
2295  this->global_vertex_index,
2296  fe_index,
2297  i);
2298 }
2299 
2300 
2301 
2302 template <template <int, int> class DH, int spacedim, bool lda>
2303 inline
2304 void
2305 DoFAccessor<0,DH<1,spacedim>, lda>::copy_from (const TriaAccessorBase<0,1,spacedim> &da)
2306 {
2307  Assert (this->dof_handler != 0, ExcInvalidObject());
2308  BaseClass::copy_from(da);
2309 }
2310 
2311 
2312 
2313 template <template <int, int> class DH, int spacedim, bool lda>
2314 template <bool lda2>
2315 inline
2316 void
2317 DoFAccessor<0,DH<1,spacedim>, lda>::copy_from (const DoFAccessor<0,DH<1,spacedim>, lda2> &a)
2318 {
2319  BaseClass::copy_from (a);
2320  set_dof_handler (a.dof_handler);
2321 }
2322 
2323 
2324 
2325 template <template <int, int> class DH, int spacedim, bool lda>
2326 template <int dim2, class DH2, bool lda2>
2327 inline
2328 bool
2329 DoFAccessor<0,DH<1,spacedim>, lda>::operator == (const DoFAccessor<dim2,DH2,lda2> &a) const
2330 {
2331  Assert (dim2 == 0, ExcCantCompareIterators());
2332  Assert (this->dof_handler == a.dof_handler, ExcCantCompareIterators());
2333  return (BaseClass::operator == (a));
2334 }
2335 
2336 
2337 
2338 template <template <int, int> class DH, int spacedim, bool lda>
2339 template <int dim2, class DH2, bool lda2>
2340 inline
2341 bool
2342 DoFAccessor<0,DH<1,spacedim>, lda>::operator != (const DoFAccessor<dim2,DH2,lda2> &a) const
2343 {
2344  Assert (dim2 == 0, ExcCantCompareIterators());
2345  Assert (this->dof_handler == a.dof_handler, ExcCantCompareIterators());
2346  return (BaseClass::operator != (a));
2347 }
2348 
2349 
2350 
2351 /*------------------------- Functions: DoFCellAccessor -----------------------*/
2352 
2353 
2354 namespace internal
2355 {
2356  namespace DoFCellAccessor
2357  {
2358  // make sure we refer to class
2359  // ::DoFCellAccessor, not
2360  // namespace
2361  // ::internal::DoFCellAccessor
2362  using ::DoFCellAccessor;
2363  using ::DoFHandler;
2364 
2370  {
2377  template <int spacedim, bool lda>
2378  static
2379  void
2381  {
2382  // check as in documentation that
2383  // cell is either active, or dofs
2384  // are only in vertices. otherwise
2385  // simply don't update the cache at
2386  // all. the get_dof_indices
2387  // function will then make sure we
2388  // don't access the invalid data
2389  if (accessor.has_children()
2390  &&
2391  (accessor.get_fe().dofs_per_cell !=
2392  accessor.get_fe().dofs_per_vertex * GeometryInfo<1>::vertices_per_cell))
2393  return;
2394 
2395  const unsigned int dofs_per_vertex = accessor.get_fe().dofs_per_vertex,
2396  dofs_per_line = accessor.get_fe().dofs_per_line,
2397  dofs_per_cell = accessor.get_fe().dofs_per_cell;
2398 
2399  // make sure the cache is at least
2400  // as big as we need it when
2401  // writing to the last element of
2402  // this cell
2403  Assert (accessor.present_index * dofs_per_cell + dofs_per_cell
2404  <=
2405  accessor.dof_handler->levels[accessor.present_level]
2406  ->cell_dof_indices_cache.size(),
2407  ExcInternalError());
2408 
2409  std::vector<types::global_dof_index>::iterator next
2410  = (accessor.dof_handler->levels[accessor.present_level]
2411  ->cell_dof_indices_cache.begin() + accessor.present_index * dofs_per_cell);
2412 
2413  for (unsigned int vertex=0; vertex<2; ++vertex)
2414  for (unsigned int d=0; d<dofs_per_vertex; ++d)
2415  *next++ = accessor.vertex_dof_index(vertex,d);
2416  for (unsigned int d=0; d<dofs_per_line; ++d)
2417  *next++ = accessor.dof_index(d);
2418  }
2419 
2420 
2421 
2422  template <int spacedim, bool lda>
2423  static
2424  void
2426  {
2427  // check as in documentation that
2428  // cell is either active, or dofs
2429  // are only in vertices. otherwise
2430  // simply don't update the cache at
2431  // all. the get_dof_indices
2432  // function will then make sure we
2433  // don't access the invalid data
2434  if (accessor.has_children()
2435  &&
2436  (accessor.get_fe().dofs_per_cell !=
2437  accessor.get_fe().dofs_per_vertex * GeometryInfo<2>::vertices_per_cell))
2438  return;
2439 
2440  const unsigned int dofs_per_vertex = accessor.get_fe().dofs_per_vertex,
2441  dofs_per_line = accessor.get_fe().dofs_per_line,
2442  dofs_per_quad = accessor.get_fe().dofs_per_quad,
2443  dofs_per_cell = accessor.get_fe().dofs_per_cell;
2444 
2445  // make sure the cache is at least
2446  // as big as we need it when
2447  // writing to the last element of
2448  // this cell
2449  Assert (accessor.present_index * dofs_per_cell + dofs_per_cell
2450  <=
2451  accessor.dof_handler->levels[accessor.present_level]
2452  ->cell_dof_indices_cache.size(),
2453  ExcInternalError());
2454 
2455  std::vector<types::global_dof_index>::iterator next
2456  = (accessor.dof_handler->levels[accessor.present_level]
2457  ->cell_dof_indices_cache.begin() + accessor.present_index * dofs_per_cell);
2458 
2459  for (unsigned int vertex=0; vertex<4; ++vertex)
2460  for (unsigned int d=0; d<dofs_per_vertex; ++d)
2461  *next++ = accessor.vertex_dof_index(vertex,d);
2462  for (unsigned int line=0; line<4; ++line)
2463  for (unsigned int d=0; d<dofs_per_line; ++d)
2464  *next++ = accessor.line(line)->dof_index(d);
2465  for (unsigned int d=0; d<dofs_per_quad; ++d)
2466  *next++ = accessor.dof_index(d);
2467  }
2468 
2469 
2470  template <int spacedim, bool lda>
2471  static
2472  void
2474  {
2475  // check as in documentation that
2476  // cell is either active, or dofs
2477  // are only in vertices. otherwise
2478  // simply don't update the cache at
2479  // all. the get_dof_indices
2480  // function will then make sure we
2481  // don't access the invalid data
2482  if (accessor.has_children()
2483  &&
2484  (accessor.get_fe().dofs_per_cell !=
2485  accessor.get_fe().dofs_per_vertex * GeometryInfo<3>::vertices_per_cell))
2486  return;
2487 
2488  const unsigned int dofs_per_vertex = accessor.get_fe().dofs_per_vertex,
2489  dofs_per_line = accessor.get_fe().dofs_per_line,
2490  dofs_per_quad = accessor.get_fe().dofs_per_quad,
2491  dofs_per_hex = accessor.get_fe().dofs_per_hex,
2492  dofs_per_cell = accessor.get_fe().dofs_per_cell;
2493 
2494  // make sure the cache is at least
2495  // as big as we need it when
2496  // writing to the last element of
2497  // this cell
2498  Assert (accessor.present_index * dofs_per_cell + dofs_per_cell
2499  <=
2500  accessor.dof_handler->levels[accessor.present_level]
2501  ->cell_dof_indices_cache.size(),
2502  ExcInternalError());
2503 
2504  std::vector<types::global_dof_index>::iterator next
2505  = (accessor.dof_handler->levels[accessor.present_level]
2506  ->cell_dof_indices_cache.begin() + accessor.present_index * dofs_per_cell);
2507 
2508  for (unsigned int vertex=0; vertex<8; ++vertex)
2509  for (unsigned int d=0; d<dofs_per_vertex; ++d)
2510  *next++ = accessor.vertex_dof_index(vertex,d);
2511  // now copy dof numbers from the line. for
2512  // lines with the wrong orientation, we have
2513  // already made sure that we're ok by picking
2514  // the correct vertices (this happens
2515  // automatically in the vertex()
2516  // function). however, if the line is in
2517  // wrong orientation, we look at it in
2518  // flipped orientation and we will have to
2519  // adjust the shape function indices that we
2520  // see to correspond to the correct
2521  // (cell-local) ordering.
2522  for (unsigned int line=0; line<12; ++line)
2523  for (unsigned int d=0; d<dofs_per_line; ++d)
2524  *next++ = accessor.line(line)->dof_index(accessor.dof_handler->get_fe().
2525  adjust_line_dof_index_for_line_orientation(d,
2526  accessor.line_orientation(line)));
2527  // now copy dof numbers from the face. for
2528  // faces with the wrong orientation, we
2529  // have already made sure that we're ok by
2530  // picking the correct lines and vertices
2531  // (this happens automatically in the
2532  // line() and vertex() functions). however,
2533  // if the face is in wrong orientation, we
2534  // look at it in flipped orientation and we
2535  // will have to adjust the shape function
2536  // indices that we see to correspond to the
2537  // correct (cell-local) ordering. The same
2538  // applies, if the face_rotation or
2539  // face_orientation is non-standard
2540  for (unsigned int quad=0; quad<6; ++quad)
2541  for (unsigned int d=0; d<dofs_per_quad; ++d)
2542  *next++ = accessor.quad(quad)->dof_index(accessor.dof_handler->get_fe().
2543  adjust_quad_dof_index_for_face_orientation(d,
2544  accessor.face_orientation(quad),
2545  accessor.face_flip(quad),
2546  accessor.face_rotation(quad)));
2547  for (unsigned int d=0; d<dofs_per_hex; ++d)
2548  *next++ = accessor.dof_index(d);
2549  }
2550 
2551 
2552  // implementation for the case of
2553  // hp::DoFHandler objects. it's
2554  // not implemented there, for no
2555  // space dimension
2556  template <int dim, int spacedim, bool lda>
2557  static
2558  void
2560  {
2561  // caches are only for cells with DoFs, i.e., for active ones
2562  if (accessor.has_children())
2563  return;
2564 
2565  const unsigned int dofs_per_cell = accessor.get_fe().dofs_per_cell;
2566 
2567  // make sure the cache is at least
2568  // as big as we need it when
2569  // writing to the last element of
2570  // this cell
2571  Assert (static_cast<unsigned int>(accessor.present_index)
2572  <
2573  accessor.dof_handler->levels[accessor.present_level]
2574  ->cell_cache_offsets.size(),
2575  ExcInternalError());
2576  Assert (accessor.dof_handler->levels[accessor.present_level]
2577  ->cell_cache_offsets[accessor.present_index]
2578  <=
2579  accessor.dof_handler->levels[accessor.present_level]
2580  ->cell_dof_indices_cache.size(),
2581  ExcInternalError());
2582 
2583  std::vector<types::global_dof_index> dof_indices (dofs_per_cell);
2584  static_cast<const ::DoFAccessor<dim,::hp::DoFHandler<dim,spacedim>,lda> &>
2585  (accessor).get_dof_indices (dof_indices, accessor.active_fe_index());
2586 
2587  types::global_dof_index *next_dof_index
2588  = &accessor.dof_handler->levels[accessor.present_level]
2589  ->cell_dof_indices_cache[accessor.dof_handler->levels[accessor.present_level]
2590  ->cell_cache_offsets[accessor.present_index]];
2591  for (unsigned int i=0; i<dofs_per_cell; ++i, ++next_dof_index)
2592  *next_dof_index = dof_indices[i];
2593  }
2594 
2595 
2596 
2604  template <int spacedim, bool lda>
2605  static
2606  void
2608  const std::vector<types::global_dof_index> &local_dof_indices)
2609  {
2610  Assert (accessor.has_children() == false,
2611  ExcInternalError());
2612 
2613  const unsigned int dofs_per_vertex = accessor.get_fe().dofs_per_vertex,
2614  dofs_per_line = accessor.get_fe().dofs_per_line,
2615  dofs_per_cell = accessor.get_fe().dofs_per_cell;
2616 
2617  Assert (local_dof_indices.size() == dofs_per_cell,
2618  ExcInternalError());
2619 
2620  unsigned int index = 0;
2621 
2622  for (unsigned int vertex=0; vertex<2; ++vertex)
2623  for (unsigned int d=0; d<dofs_per_vertex; ++d, ++index)
2624  accessor.set_vertex_dof_index(vertex,d,
2625  local_dof_indices[index]);
2626 
2627  for (unsigned int d=0; d<dofs_per_line; ++d, ++index)
2628  accessor.set_dof_index(d, local_dof_indices[index]);
2629 
2630  Assert (index == dofs_per_cell,
2631  ExcInternalError());
2632  }
2633 
2634 
2635 
2636  template <int spacedim, bool lda>
2637  static
2638  void
2640  const std::vector<types::global_dof_index> &local_dof_indices)
2641  {
2642  Assert (accessor.has_children() == false,
2643  ExcInternalError());
2644 
2645  const unsigned int dofs_per_vertex = accessor.get_fe().dofs_per_vertex,
2646  dofs_per_line = accessor.get_fe().dofs_per_line,
2647  dofs_per_quad = accessor.get_fe().dofs_per_quad,
2648  dofs_per_cell = accessor.get_fe().dofs_per_cell;
2649 
2650  Assert (local_dof_indices.size() == dofs_per_cell,
2651  ExcInternalError());
2652 
2653  unsigned int index = 0;
2654 
2655  for (unsigned int vertex=0; vertex<4; ++vertex)
2656  for (unsigned int d=0; d<dofs_per_vertex; ++d, ++index)
2657  accessor.set_vertex_dof_index(vertex,d,
2658  local_dof_indices[index]);
2659  for (unsigned int line=0; line<4; ++line)
2660  for (unsigned int d=0; d<dofs_per_line; ++d, ++index)
2661  accessor.line(line)->set_dof_index(d, local_dof_indices[index]);
2662 
2663  for (unsigned int d=0; d<dofs_per_quad; ++d, ++index)
2664  accessor.set_dof_index(d, local_dof_indices[index]);
2665 
2666  Assert (index == dofs_per_cell,
2667  ExcInternalError());
2668  }
2669 
2670 
2671 
2672  template <int spacedim, bool lda>
2673  static
2674  void
2676  const std::vector<types::global_dof_index> &local_dof_indices)
2677  {
2678  Assert (accessor.has_children() == false,
2679  ExcInternalError());
2680 
2681  const unsigned int dofs_per_vertex = accessor.get_fe().dofs_per_vertex,
2682  dofs_per_line = accessor.get_fe().dofs_per_line,
2683  dofs_per_quad = accessor.get_fe().dofs_per_quad,
2684  dofs_per_hex = accessor.get_fe().dofs_per_hex,
2685  dofs_per_cell = accessor.get_fe().dofs_per_cell;
2686 
2687  Assert (local_dof_indices.size() == dofs_per_cell,
2688  ExcInternalError());
2689 
2690  unsigned int index = 0;
2691 
2692  for (unsigned int vertex=0; vertex<8; ++vertex)
2693  for (unsigned int d=0; d<dofs_per_vertex; ++d, ++index)
2694  accessor.set_vertex_dof_index(vertex,d,
2695  local_dof_indices[index]);
2696  // now copy dof numbers into the line. for
2697  // lines with the wrong orientation, we have
2698  // already made sure that we're ok by picking
2699  // the correct vertices (this happens
2700  // automatically in the vertex()
2701  // function). however, if the line is in
2702  // wrong orientation, we look at it in
2703  // flipped orientation and we will have to
2704  // adjust the shape function indices that we
2705  // see to correspond to the correct
2706  // (cell-local) ordering.
2707  for (unsigned int line=0; line<12; ++line)
2708  for (unsigned int d=0; d<dofs_per_line; ++d, ++index)
2709  accessor.line(line)->set_dof_index(accessor.dof_handler->get_fe().
2710  adjust_line_dof_index_for_line_orientation(d,
2711  accessor.line_orientation(line)),
2712  local_dof_indices[index]);
2713  // now copy dof numbers into the face. for
2714  // faces with the wrong orientation, we
2715  // have already made sure that we're ok by
2716  // picking the correct lines and vertices
2717  // (this happens automatically in the
2718  // line() and vertex() functions). however,
2719  // if the face is in wrong orientation, we
2720  // look at it in flipped orientation and we
2721  // will have to adjust the shape function
2722  // indices that we see to correspond to the
2723  // correct (cell-local) ordering. The same
2724  // applies, if the face_rotation or
2725  // face_orientation is non-standard
2726  for (unsigned int quad=0; quad<6; ++quad)
2727  for (unsigned int d=0; d<dofs_per_quad; ++d, ++index)
2728  accessor.quad(quad)->set_dof_index(accessor.dof_handler->get_fe().
2729  adjust_quad_dof_index_for_face_orientation(d,
2730  accessor.face_orientation(quad),
2731  accessor.face_flip(quad),
2732  accessor.face_rotation(quad)),
2733  local_dof_indices[index]);
2734  for (unsigned int d=0; d<dofs_per_hex; ++d, ++index)
2735  accessor.set_dof_index(d, local_dof_indices[index]);
2736 
2737  Assert (index == dofs_per_cell,
2738  ExcInternalError());
2739  }
2740 
2741 
2742  // implementation for the case of
2743  // hp::DoFHandler objects. it's
2744  // not implemented there, for no
2745  // space dimension
2746  template <int dim, int spacedim, bool lda>
2747  static
2748  void
2750  const std::vector<types::global_dof_index> &)
2751  {
2752  Assert (false, ExcNotImplemented());
2753  }
2754 
2755 
2756 
2762  template <int dim, int spacedim, bool lda>
2763  static
2764  unsigned int
2766  {
2767  // ::DoFHandler only supports a
2768  // single active fe with index
2769  // zero
2770  return 0;
2771  }
2772 
2773 
2774 
2775  template <int dim, int spacedim, bool lda>
2776  static
2777  unsigned int
2779  {
2780  Assert (static_cast<unsigned int>(accessor.level()) < accessor.dof_handler->levels.size(),
2781  ExcMessage ("DoFHandler not initialized"));
2782 
2783  return accessor.dof_handler->levels[accessor.level()]
2784  ->active_fe_index(accessor.present_index);
2785  }
2786 
2787 
2788 
2795  template <int dim, int spacedim, bool lda>
2796  static
2797  void
2799  const unsigned int i)
2800  {
2801  // ::DoFHandler only supports a
2802  // single active fe with index
2803  // zero
2804  typedef ::DoFAccessor<dim,DoFHandler<dim,spacedim>, lda> BaseClass;
2805  Assert (i == 0, typename BaseClass::ExcInvalidObject());
2806  }
2807 
2808 
2809 
2810  template <int dim, int spacedim, bool lda>
2811  static
2812  void
2814  const unsigned int i)
2815  {
2816  typedef ::DoFAccessor<dim,DoFHandler<dim,spacedim>, lda> BaseClass;
2817  Assert (accessor.dof_handler != 0,
2818  typename BaseClass::ExcInvalidObject());
2819  Assert (static_cast<unsigned int>(accessor.level()) <
2820  accessor.dof_handler->levels.size(),
2821  ExcMessage ("DoFHandler not initialized"));
2822 
2823  accessor.dof_handler->levels[accessor.level()]
2824  ->set_active_fe_index (accessor.present_index, i);
2825  }
2826 
2827 
2828 
2829  template <int dim, int spacedim, bool lda, typename ForwardIterator, class OutputVector>
2830  static
2831  void
2832  distribute_local_to_global (const DoFCellAccessor<::DoFHandler<dim,spacedim>, lda> &accessor,
2833  ForwardIterator local_source_begin,
2834  ForwardIterator local_source_end,
2835  OutputVector &global_destination)
2836  {
2837  typedef ::DoFAccessor<dim,DoFHandler<dim,spacedim>, lda> BaseClass;
2838  Assert (accessor.dof_handler != 0,
2839  typename BaseClass::ExcInvalidObject());
2840  Assert (&accessor.get_fe() != 0,
2841  typename BaseClass::ExcInvalidObject());
2842  Assert (static_cast<unsigned int>(local_source_end-local_source_begin)
2843  ==
2844  accessor.get_fe().dofs_per_cell,
2845  typename BaseClass::ExcVectorDoesNotMatch());
2846  Assert (accessor.dof_handler->n_dofs() == global_destination.size(),
2847  typename BaseClass::ExcVectorDoesNotMatch());
2848 
2849  Assert (!accessor.has_children(),
2850  ExcMessage ("Cell must be active"));
2851 
2852  const unsigned int n_dofs = local_source_end - local_source_begin;
2853 
2854  types::global_dof_index *dofs = &accessor.dof_handler->levels[accessor.level()]
2855  ->cell_dof_indices_cache[accessor.present_index * n_dofs];
2856 
2857  // distribute cell vector
2858  global_destination.add(n_dofs, dofs, local_source_begin);
2859  }
2860 
2861 
2862 
2863  template <int dim, int spacedim, bool lda, typename ForwardIterator, class OutputVector>
2864  static
2865  void
2866  distribute_local_to_global (const DoFCellAccessor<::hp::DoFHandler<dim,spacedim>, lda> &accessor,
2867  ForwardIterator local_source_begin,
2868  ForwardIterator local_source_end,
2869  OutputVector &global_destination)
2870  {
2871  typedef ::DoFAccessor<dim,DoFHandler<dim,spacedim>, lda> BaseClass;
2872  Assert (accessor.dof_handler != 0,
2873  typename BaseClass::ExcInvalidObject());
2874  Assert (&accessor.get_fe() != 0,
2875  typename BaseClass::ExcInvalidObject());
2876  Assert (local_source_end-local_source_begin == accessor.get_fe().dofs_per_cell,
2877  typename BaseClass::ExcVectorDoesNotMatch());
2878  Assert (accessor.dof_handler->n_dofs() == global_destination.size(),
2879  typename BaseClass::ExcVectorDoesNotMatch());
2880 
2881  const unsigned int n_dofs = local_source_end - local_source_begin;
2882 
2883 //TODO[WB/MK]: This function could me made more efficient because it allocates memory, which could be avoided by passing in another argument as a scratch array. This should be fixed eventually
2884 
2885  // get indices of dofs
2886  std::vector<types::global_dof_index> dofs (n_dofs);
2887  accessor.get_dof_indices (dofs);
2888 
2889  // distribute cell vector
2890  global_destination.add (n_dofs, dofs.begin(), local_source_begin);
2891  }
2892 
2893 
2894 
2895  template <int dim, int spacedim, bool lda, typename ForwardIterator, class OutputVector>
2896  static
2897  void
2898  distribute_local_to_global (const DoFCellAccessor<::DoFHandler<dim,spacedim>, lda> &accessor,
2899  const ConstraintMatrix &constraints,
2900  ForwardIterator local_source_begin,
2901  ForwardIterator local_source_end,
2902  OutputVector &global_destination)
2903  {
2904  typedef ::DoFAccessor<dim,DoFHandler<dim,spacedim>, lda> BaseClass;
2905  Assert (accessor.dof_handler != 0,
2906  typename BaseClass::ExcInvalidObject());
2907  Assert (&accessor.get_fe() != 0,
2908  typename BaseClass::ExcInvalidObject());
2909  Assert (local_source_end-local_source_begin == accessor.get_fe().dofs_per_cell,
2910  typename BaseClass::ExcVectorDoesNotMatch());
2911  Assert (accessor.dof_handler->n_dofs() == global_destination.size(),
2912  typename BaseClass::ExcVectorDoesNotMatch());
2913 
2914  Assert (!accessor.has_children(),
2915  ExcMessage ("Cell must be active."));
2916 
2917  const unsigned int n_dofs = local_source_end - local_source_begin;
2918 
2919  types::global_dof_index *dofs = &accessor.dof_handler->levels[accessor.level()]
2920  ->cell_dof_indices_cache[accessor.present_index * n_dofs];
2921 
2922  // distribute cell vector
2923  constraints.distribute_local_to_global (local_source_begin, local_source_end,
2924  dofs, global_destination);
2925  }
2926 
2927 
2928 
2929  template <int dim, int spacedim, bool lda, typename ForwardIterator, class OutputVector>
2930  static
2931  void
2932  distribute_local_to_global (const DoFCellAccessor<::hp::DoFHandler<dim,spacedim>, lda> &accessor,
2933  const ConstraintMatrix &constraints,
2934  ForwardIterator local_source_begin,
2935  ForwardIterator local_source_end,
2936  OutputVector &global_destination)
2937  {
2938  typedef ::DoFAccessor<dim,DoFHandler<dim,spacedim>, lda> BaseClass;
2939  Assert (accessor.dof_handler != 0,
2940  typename BaseClass::ExcInvalidObject());
2941  Assert (&accessor.get_fe() != 0,
2942  typename BaseClass::ExcInvalidObject());
2943  Assert (local_source_end-local_source_begin == accessor.get_fe().dofs_per_cell,
2944  typename BaseClass::ExcVectorDoesNotMatch());
2945  Assert (accessor.dof_handler->n_dofs() == global_destination.size(),
2946  typename BaseClass::ExcVectorDoesNotMatch());
2947 
2948  const unsigned int n_dofs = local_source_end - local_source_begin;
2949 
2950 //TODO[WB/MK]: This function could me made more efficient because it allocates memory, which could be avoided by passing in another argument as a scratch array. This should be fixed eventually
2951 
2952  // get indices of dofs
2953  std::vector<types::global_dof_index> dofs (n_dofs);
2954  accessor.get_dof_indices (dofs);
2955 
2956  // distribute cell vector
2957  constraints.distribute_local_to_global (local_source_begin, local_source_end,
2958  dofs.begin(), global_destination);
2959  }
2960 
2961 
2962 
2963  template <int dim, int spacedim, bool lda, typename number, class OutputMatrix>
2964  static
2965  void
2966  distribute_local_to_global (const DoFCellAccessor<::DoFHandler<dim,spacedim>, lda> &accessor,
2967  const ::FullMatrix<number> &local_source,
2968  OutputMatrix &global_destination)
2969  {
2970  typedef ::DoFAccessor<dim,DoFHandler<dim,spacedim>, lda> BaseClass;
2971  Assert (accessor.dof_handler != 0,
2972  typename BaseClass::ExcInvalidObject());
2973  Assert (&accessor.get_fe() != 0,
2974  typename BaseClass::ExcInvalidObject());
2975  Assert (local_source.m() == accessor.get_fe().dofs_per_cell,
2976  typename BaseClass::ExcMatrixDoesNotMatch());
2977  Assert (local_source.n() == accessor.get_fe().dofs_per_cell,
2978  typename BaseClass::ExcMatrixDoesNotMatch());
2979  Assert (accessor.dof_handler->n_dofs() == global_destination.m(),
2980  typename BaseClass::ExcMatrixDoesNotMatch());
2981  Assert (accessor.dof_handler->n_dofs() == global_destination.n(),
2982  typename BaseClass::ExcMatrixDoesNotMatch());
2983 
2984  Assert (!accessor.has_children(),
2985  ExcMessage ("Cell must be active."));
2986 
2987  const unsigned int n_dofs = local_source.m();
2988 
2989  types::global_dof_index *dofs = &accessor.dof_handler->levels[accessor.level()]
2990  ->cell_dof_indices_cache[accessor.present_index * n_dofs];
2991 
2992  // distribute cell matrix
2993  for (unsigned int i=0; i<n_dofs; ++i)
2994  global_destination.add(dofs[i], n_dofs, dofs,
2995  &local_source(i,0));
2996  }
2997 
2998 
2999 
3000  template <int dim, int spacedim, bool lda, typename number, class OutputMatrix>
3001  static
3002  void
3003  distribute_local_to_global (const DoFCellAccessor<::hp::DoFHandler<dim,spacedim>, lda> &accessor,
3004  const ::FullMatrix<number> &local_source,
3005  OutputMatrix &global_destination)
3006  {
3007  typedef ::DoFAccessor<dim,DoFHandler<dim,spacedim>, lda> BaseClass;
3008  Assert (accessor.dof_handler != 0,
3009  typename BaseClass::ExcInvalidObject());
3010  Assert (&accessor.get_fe() != 0,
3011  typename BaseClass::ExcInvalidObject());
3012  Assert (local_source.m() == accessor.get_fe().dofs_per_cell,
3013  typename BaseClass::ExcMatrixDoesNotMatch());
3014  Assert (local_source.n() == accessor.get_fe().dofs_per_cell,
3015  typename BaseClass::ExcVectorDoesNotMatch());
3016  Assert (accessor.dof_handler->n_dofs() == global_destination.m(),
3017  typename BaseClass::ExcMatrixDoesNotMatch());
3018  Assert (accessor.dof_handler->n_dofs() == global_destination.n(),
3019  typename BaseClass::ExcMatrixDoesNotMatch());
3020 
3021  const unsigned int n_dofs = local_source.size();
3022 
3023 //TODO[WB/MK]: This function could me made more efficient because it allocates memory, which could be avoided by passing in another argument as a scratch array.
3024 
3025  // get indices of dofs
3026  std::vector<types::global_dof_index> dofs (n_dofs);
3027  accessor.get_dof_indices (dofs);
3028 
3029  // distribute cell matrix
3030  global_destination.add(dofs,local_source);
3031  }
3032 
3033 
3034 
3035  template <int dim, int spacedim, bool lda, typename number,
3036  class OutputMatrix, typename OutputVector>
3037  static
3038  void
3039  distribute_local_to_global (const DoFCellAccessor<::DoFHandler<dim,spacedim>, lda> &accessor,
3040  const ::FullMatrix<number> &local_matrix,
3041  const ::Vector<number> &local_vector,
3042  OutputMatrix &global_matrix,
3043  OutputVector &global_vector)
3044  {
3045  typedef ::DoFAccessor<dim,DoFHandler<dim,spacedim>, lda> BaseClass;
3046  Assert (accessor.dof_handler != 0,
3047  typename BaseClass::ExcInvalidObject());
3048  Assert (&accessor.get_fe() != 0,
3049  typename BaseClass::ExcInvalidObject());
3050  Assert (local_matrix.m() == accessor.get_fe().dofs_per_cell,
3051  typename BaseClass::ExcMatrixDoesNotMatch());
3052  Assert (local_matrix.n() == accessor.get_fe().dofs_per_cell,
3053  typename BaseClass::ExcVectorDoesNotMatch());
3054  Assert (accessor.dof_handler->n_dofs() == global_matrix.m(),
3055  typename BaseClass::ExcMatrixDoesNotMatch());
3056  Assert (accessor.dof_handler->n_dofs() == global_matrix.n(),
3057  typename BaseClass::ExcMatrixDoesNotMatch());
3058  Assert (local_vector.size() == accessor.get_fe().dofs_per_cell,
3059  typename BaseClass::ExcVectorDoesNotMatch());
3060  Assert (accessor.dof_handler->n_dofs() == global_vector.size(),
3061  typename BaseClass::ExcVectorDoesNotMatch());
3062 
3063  Assert (!accessor.has_children(),
3064  ExcMessage ("Cell must be active."));
3065 
3066  const unsigned int n_dofs = accessor.get_fe().dofs_per_cell;
3067  types::global_dof_index *dofs = &accessor.dof_handler->levels[accessor.level()]
3068  ->cell_dof_indices_cache[accessor.present_index *n_dofs];
3069 
3070  // distribute cell matrices
3071  for (unsigned int i=0; i<n_dofs; ++i)
3072  {
3073  global_matrix.add(dofs[i], n_dofs, dofs, &local_matrix(i,0));
3074  global_vector(dofs[i]) += local_vector(i);
3075  }
3076  }
3077 
3078 
3079 
3080  template <int dim, int spacedim, bool lda, typename number,
3081  class OutputMatrix, typename OutputVector>
3082  static
3083  void
3084  distribute_local_to_global (const DoFCellAccessor<::hp::DoFHandler<dim,spacedim>, lda> &accessor,
3085  const ::FullMatrix<number> &local_matrix,
3086  const ::Vector<number> &local_vector,
3087  OutputMatrix &global_matrix,
3088  OutputVector &global_vector)
3089  {
3090  typedef ::DoFAccessor<dim,DoFHandler<dim,spacedim>, lda> BaseClass;
3091  Assert (accessor.dof_handler != 0,
3092  typename BaseClass::ExcInvalidObject());
3093  Assert (&accessor.get_fe() != 0,
3094  typename BaseClass::ExcInvalidObject());
3095  Assert (local_matrix.m() == accessor.get_fe().dofs_per_cell,
3096  typename BaseClass::ExcMatrixDoesNotMatch());
3097  Assert (local_matrix.n() == accessor.get_fe().dofs_per_cell,
3098  typename BaseClass::ExcVectorDoesNotMatch());
3099  Assert (accessor.dof_handler->n_dofs() == global_matrix.m(),
3100  typename BaseClass::ExcMatrixDoesNotMatch());
3101  Assert (accessor.dof_handler->n_dofs() == global_matrix.n(),
3102  typename BaseClass::ExcMatrixDoesNotMatch());
3103  Assert (local_vector.size() == accessor.get_fe().dofs_per_cell,
3104  typename BaseClass::ExcVectorDoesNotMatch());
3105  Assert (accessor.dof_handler->n_dofs() == global_vector.size(),
3106  typename BaseClass::ExcVectorDoesNotMatch());
3107 
3108  const unsigned int n_dofs = local_matrix.size();
3109 
3110 //TODO[WB/MK]: This function could me made more efficient because it
3111 //allocates memory, which could be avoided by passing in another
3112 //argument as a scratch array. Comment(GK) Do not bother and leva this
3113 //to ConstraintMatrix or MeshWorker::Assembler
3114 
3115  // get indices of dofs
3116  std::vector<types::global_dof_index> dofs (n_dofs);
3117  accessor.get_dof_indices (dofs);
3118 
3119  // distribute cell matrix and vector
3120  global_matrix.add(dofs,local_matrix);
3121  global_vector.add(dofs,local_vector);
3122  }
3123  };
3124  }
3125 }
3126 
3127 
3128 template <class DH, bool lda>
3129 inline
3132  const int level,
3133  const int index,
3134  const AccessorData *local_data)
3135  :
3136  DoFAccessor<DH::dimension,DH,lda> (tria,level,index, local_data)
3137 {}
3138 
3139 
3140 template <class DH, bool lda>
3141 template <int structdim2, int dim2, int spacedim2>
3142 inline
3144 {
3145  Assert (false, typename BaseClass::ExcInvalidObject());
3146 }
3147 
3148 
3149 
3150 template <class DH, bool lda>
3151 template <int dim2, class DH2, bool lda2>
3152 inline
3154  :
3155  BaseClass(other)
3156 {}
3157 
3158 
3159 template <class DH, bool lda>
3160 inline
3162 DoFCellAccessor<DH,lda>::neighbor (const unsigned int i) const
3163 {
3165  q (this->tria,
3166  this->neighbor_level (i),
3167  this->neighbor_index (i),
3168  this->dof_handler);
3169 
3170 #ifdef DEBUG
3172  Assert (q->used(), TriaAccessorExceptions::ExcUnusedCellAsNeighbor());
3173 #endif
3174  return q;
3175 }
3176 
3177 
3178 template <class DH, bool lda>
3179 inline
3181 DoFCellAccessor<DH,lda>::child (const unsigned int i) const
3182 {
3184  q (this->tria,
3185  this->level()+1,
3186  this->child_index (i),
3187  this->dof_handler);
3188 
3189 #ifdef DEBUG
3191  Assert (q->used(), TriaAccessorExceptions::ExcUnusedCellAsChild());
3192 #endif
3193  return q;
3194 }
3195 
3196 
3197 template <class DH, bool lda>
3198 inline
3201 {
3203  q (this->tria,
3204  this->level() - 1,
3205  this->parent_index (),
3206  this->dof_handler);
3207 
3208  return q;
3209 }
3210 
3211 
3212 namespace internal
3213 {
3214  namespace DoFCellAccessor
3215  {
3216  template <class DH, bool lda>
3217  inline
3218  TriaIterator<::DoFAccessor<DH::dimension-1,DH,lda> >
3219  get_face (const ::DoFCellAccessor<DH,lda> &cell,
3220  const unsigned int i,
3221  const ::internal::int2type<1>)
3222  {
3224  a (&cell.get_triangulation(),
3225  ((i == 0) && cell.at_boundary(0)
3226  ?
3228  :
3229  ((i == 1) && cell.at_boundary(1)
3230  ?
3232  :
3234  cell.vertex_index(i),
3235  &cell.get_dof_handler());
3236  return ::TriaIterator<::DoFAccessor<0,DH,lda> > (a);
3237  }
3238 
3239 
3240  template <class DH, bool lda>
3241  inline
3242  TriaIterator<::DoFAccessor<DH::dimension-1,DH,lda> >
3243  get_face (const ::DoFCellAccessor<DH,lda> &cell,
3244  const unsigned int i,
3245  const ::internal::int2type<2>)
3246  {
3247  return cell.line(i);
3248  }
3249 
3250 
3251  template <class DH, bool lda>
3252  inline
3253  TriaIterator<::DoFAccessor<DH::dimension-1,DH,lda> >
3254  get_face (const ::DoFCellAccessor<DH,lda> &cell,
3255  const unsigned int i,
3256  const ::internal::int2type<3>)
3257  {
3258  return cell.quad(i);
3259  }
3260  }
3261 }
3262 
3263 
3264 template <class DH, bool lda>
3265 inline
3266 TriaIterator<DoFAccessor<DH::dimension-1,DH,lda> >
3267 DoFCellAccessor<DH,lda>::face (const unsigned int i) const
3268 {
3270  Assert (static_cast<unsigned int>(this->level()) < this->dof_handler->levels.size(),
3271  ExcMessage ("DoFHandler not initialized"));
3272 
3273  const unsigned int dim = DH::dimension;
3274  return ::internal::DoFCellAccessor::get_face (*this, i, ::internal::int2type<dim>());
3275 }
3276 
3277 
3278 
3279 template <class DH, bool lda>
3280 inline
3281 void
3282 DoFCellAccessor<DH,lda>::get_dof_indices (std::vector<types::global_dof_index> &dof_indices) const
3283 {
3284  Assert (this->active(), ExcMessage ("get_dof_indices() only works on active cells."));
3285  Assert (this->is_artificial() == false,
3286  ExcMessage ("Can't ask for DoF indices on artificial cells."));
3287  AssertDimension (dof_indices.size(), this->get_fe().dofs_per_cell);
3288 
3289  const types::global_dof_index *cache
3290  = this->dof_handler->levels[this->present_level]
3291  ->get_cell_cache_start (this->present_index, this->get_fe().dofs_per_cell);
3292  for (unsigned int i=0; i<this->get_fe().dofs_per_cell; ++i, ++cache)
3293  dof_indices[i] = *cache;
3294 }
3295 
3296 
3297 
3298 template<class DH, bool lda>
3299 inline
3300 void DoFCellAccessor<DH,lda>::get_mg_dof_indices (std::vector<types::global_dof_index> &dof_indices) const
3301 {
3302  DoFAccessor<dim, DH,lda>::get_mg_dof_indices (this->level (), dof_indices);
3303 }
3304 
3305 
3306 
3307 template<class DH, bool lda>
3308 inline
3309 void DoFCellAccessor<DH,lda>::set_mg_dof_indices (const std::vector<types::global_dof_index> &dof_indices)
3310 {
3311  DoFAccessor<dim, DH,lda>::set_mg_dof_indices (this->level (), dof_indices);
3312 }
3313 
3314 
3315 
3316 template<class DH, bool lda>
3317 inline
3318 void DoFCellAccessor<DH,lda>::get_active_or_mg_dof_indices (std::vector<types::global_dof_index> &dof_indices) const
3319 {
3320  if (lda)
3321  get_mg_dof_indices (dof_indices);
3322  else
3323  get_dof_indices (dof_indices);
3324 }
3325 
3326 
3327 
3328 template <class DH, bool lda>
3329 template <class InputVector, typename number>
3330 inline
3331 void
3332 DoFCellAccessor<DH,lda>::get_dof_values (const InputVector &values,
3333  Vector<number> &local_values) const
3334 {
3335  get_dof_values (values, local_values.begin(), local_values.end());
3336 }
3337 
3338 
3339 
3340 template <class DH, bool lda>
3341 template <class InputVector, typename ForwardIterator>
3342 inline
3343 void
3344 DoFCellAccessor<DH,lda>::get_dof_values (const InputVector &values,
3345  ForwardIterator local_values_begin,
3346  ForwardIterator local_values_end) const
3347 {
3348  Assert (this->is_artificial() == false,
3349  ExcMessage ("Can't ask for DoF indices on artificial cells."));
3350  Assert (!this->has_children(),
3351  ExcMessage ("Cell must be active."));
3352 
3353  Assert (static_cast<unsigned int>(local_values_end-local_values_begin)
3354  == this->get_fe().dofs_per_cell,
3355  typename DoFCellAccessor::ExcVectorDoesNotMatch());
3356  Assert (values.size() == this->get_dof_handler().n_dofs(),
3357  typename DoFCellAccessor::ExcVectorDoesNotMatch());
3358 
3359  const types::global_dof_index *cache
3360  = this->dof_handler->levels[this->present_level]
3361  ->get_cell_cache_start (this->present_index, this->get_fe().dofs_per_cell);
3362 
3363  values.extract_subvector_to (cache,
3364  cache + this->get_fe().dofs_per_cell,
3365  local_values_begin);
3366 }
3367 
3368 
3369 
3370 template <class DH, bool lda>
3371 template <class InputVector, typename ForwardIterator>
3372 inline
3373 void
3375  const InputVector &values,
3376  ForwardIterator local_values_begin,
3377  ForwardIterator local_values_end) const
3378 {
3379  Assert (this->is_artificial() == false,
3380  ExcMessage ("Can't ask for DoF indices on artificial cells."));
3381  Assert (!this->has_children(),
3382  ExcMessage ("Cell must be active."));
3383 
3384  Assert (static_cast<unsigned int>(local_values_end-local_values_begin)
3385  == this->get_fe().dofs_per_cell,
3386  typename DoFCellAccessor::ExcVectorDoesNotMatch());
3387  Assert (values.size() == this->get_dof_handler().n_dofs(),
3388  typename DoFCellAccessor::ExcVectorDoesNotMatch());
3389 
3390 
3391  const types::global_dof_index *cache
3392  = this->dof_handler->levels[this->present_level]
3393  ->get_cell_cache_start (this->present_index, this->get_fe().dofs_per_cell);
3394 
3395  constraints.get_dof_values(values, *cache, local_values_begin,
3396  local_values_end);
3397 }
3398 
3399 
3400 
3401 template <class DH, bool lda>
3402 template <class OutputVector, typename number>
3403 inline
3404 void
3406  OutputVector &values) const
3407 {
3408  Assert (this->is_artificial() == false,
3409  ExcMessage ("Can't ask for DoF indices on artificial cells."));
3410  Assert (!this->has_children(),
3411  ExcMessage ("Cell must be active."));
3412 
3413  Assert (static_cast<unsigned int>(local_values.size())
3414  == this->get_fe().dofs_per_cell,
3415  typename DoFCellAccessor::ExcVectorDoesNotMatch());
3416  Assert (values.size() == this->get_dof_handler().n_dofs(),
3417  typename DoFCellAccessor::ExcVectorDoesNotMatch());
3418 
3419 
3420  const types::global_dof_index *cache
3421  = this->dof_handler->levels[this->present_level]
3422  ->get_cell_cache_start (this->present_index, this->get_fe().dofs_per_cell);
3423 
3424  for (unsigned int i=0; i<this->get_fe().dofs_per_cell; ++i, ++cache)
3425  values(*cache) = local_values(i);
3426 }
3427 
3428 
3429 
3430 template <class DH, bool lda>
3431 inline
3434 {
3435  return ::internal::DoFAccessor::get_fe (this->dof_handler->get_fe(), active_fe_index());
3436 }
3437 
3438 
3439 
3440 template <class DH, bool lda>
3441 inline
3442 unsigned int
3444 {
3445  return ::internal::DoFCellAccessor::Implementation::active_fe_index (*this);
3446 }
3447 
3448 
3449 
3450 template <class DH, bool lda>
3451 inline
3452 void
3454 {
3456 }
3457 
3458 
3459 
3460 
3461 template <class DH, bool lda>
3462 template <typename number, typename OutputVector>
3463 inline
3464 void
3466  const Vector<number> &local_source,
3467  OutputVector &global_destination) const
3468 {
3469  ::internal::DoFCellAccessor::Implementation::
3470  distribute_local_to_global (*this, local_source.begin(),
3471  local_source.end(), global_destination);
3472 }
3473 
3474 
3475 
3476 template <class DH, bool lda>
3477 template <typename ForwardIterator, typename OutputVector>
3478 inline
3479 void
3481  ForwardIterator local_source_begin,
3482  ForwardIterator local_source_end,
3483  OutputVector &global_destination) const
3484 {
3485  ::internal::DoFCellAccessor::Implementation::
3486  distribute_local_to_global (*this, local_source_begin,
3487  local_source_end, global_destination);
3488 }
3489 
3490 
3491 
3492 template <class DH, bool lda>
3493 template <typename ForwardIterator, typename OutputVector>
3494 inline
3495 void
3497  const ConstraintMatrix &constraints,
3498  ForwardIterator local_source_begin,
3499  ForwardIterator local_source_end,
3500  OutputVector &global_destination) const
3501 {
3502  ::internal::DoFCellAccessor::Implementation::
3503  distribute_local_to_global (*this, constraints, local_source_begin,
3504  local_source_end, global_destination);
3505 }
3506 
3507 
3508 
3509 template <class DH, bool lda>
3510 template <typename number, typename OutputMatrix>
3511 inline
3512 void
3514  const FullMatrix<number> &local_source,
3515  OutputMatrix &global_destination) const
3516 {
3517  ::internal::DoFCellAccessor::Implementation::
3518  distribute_local_to_global (*this,local_source,global_destination);
3519 }
3520 
3521 
3522 
3523 template <class DH, bool lda>
3524 template <typename number, typename OutputMatrix, typename OutputVector>
3525 inline
3526 void
3528  const FullMatrix<number> &local_matrix,
3529  const Vector<number> &local_vector,
3530  OutputMatrix &global_matrix,
3531  OutputVector &global_vector) const
3532 {
3533  ::internal::DoFCellAccessor::Implementation::
3534  distribute_local_to_global (*this,local_matrix,local_vector,
3535  global_matrix,global_vector);
3536 }
3537 
3538 
3539 DEAL_II_NAMESPACE_CLOSE
3540 
3541 #endif
static void update_cell_dof_indices_cache(const DoFCellAccessor< DoFHandler< 1, spacedim >, lda > &accessor)
types::global_dof_index mg_dof_index(const int level, const unsigned int i) const
TriaIterator< DoFCellAccessor< DH, level_dof_access > > child(const unsigned int) const
TriaIterator< DoFCellAccessor< DH, level_dof_access > > neighbor(const unsigned int) const
void distribute_local_to_global(const Vector< number > &local_source, OutputVector &global_destination) const
const DH & get_dof_handler() const
void set_active_fe_index(const unsigned int i)
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:858
SmartPointer< const hp::FECollection< dim, spacedim >, hp::DoFHandler< dim, spacedim > > finite_elements
Definition: dof_handler.h:645
TriaIterator< DoFAccessor< structdim, DH, level_dof_access > > parent() const
unsigned int active_fe_index() const
static const unsigned int dimension
Definition: dof_handler.h:223
bool operator==(const DoFAccessor< dim2, DH2, level_dof_access2 > &) const
void set_dof_handler(DH *dh)
void get_active_or_mg_dof_indices(std::vector< types::global_dof_index > &dof_indices) const
::ExceptionBase & ExcMessage(std::string arg1)
typename::internal::DoFAccessor::Inheritance< structdim, dimension, space_dimension >::BaseClass BaseClass
Definition: dof_accessor.h:215
static void set_active_fe_index(const DoFCellAccessor< DoFHandler< dim, spacedim >, lda > &, const unsigned int i)
const unsigned int dofs_per_quad
Definition: fe_base.h:226
SmartPointer< const FiniteElement< dim, spacedim >, DoFHandler< dim, spacedim > > selected_fe
Definition: dof_handler.h:935
SmartPointer< const Triangulation< dim, spacedim >, DoFHandler< dim, spacedim > > tria
Definition: dof_handler.h:920
static types::global_dof_index get_vertex_dof_index(const ::DoFHandler< dim, spacedim > &dof_handler, const unsigned int vertex_index, const unsigned int fe_index, const unsigned int local_index)
typename::internal::DoFHandler::Iterators< DH, level_dof_access >::line_iterator line(const unsigned int i) const
#define AssertThrow(cond, exc)
Definition: exceptions.h:362
iterator end()
const unsigned int dofs_per_line
Definition: fe_base.h:220
static void set_dof_indices(DoFCellAccessor< DoFHandler< 1, spacedim >, lda > &accessor, const std::vector< types::global_dof_index > &local_dof_indices)
static bool fe_is_active_on_vertex(const ::hp::DoFHandler< dim, spacedim > &dof_handler, const unsigned int vertex_index, const unsigned int fe_index)
DH * dof_handler
Definition: dof_accessor.h:675
static types::global_dof_index get_dof_index(const ::DoFHandler< 1, spacedim > &dof_handler, const unsigned int obj_level, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index,::internal::int2type< 1 >)
types::global_dof_index mg_vertex_dof_index(const int level, const unsigned int vertex, const unsigned int i, const unsigned int fe_index=DH::default_fe_index) const
TriaIterator< DoFAccessor< structdim, DH, level_dof_access > > child(const unsigned int c) const
void set_mg_dof_indices(const int level, const std::vector< types::global_dof_index > &dof_indices, const unsigned int fe_index=DH::default_fe_index)
types::global_dof_index dof_index(const unsigned int i, const unsigned int fe_index=DH::default_fe_index) const
static void set_vertex_dof_index(::DoFHandler< dim, spacedim > &dof_handler, const unsigned int vertex_index, const unsigned int fe_index, const unsigned int local_index, const types::global_dof_index global_index)
static unsigned int active_fe_index(const DoFCellAccessor< DoFHandler< dim, spacedim >, lda > &)
void distribute_local_to_global(const InVector &local_vector, const std::vector< size_type > &local_dof_indices, OutVector &global_vector) const
static unsigned int nth_active_vertex_fe_index(const ::hp::DoFHandler< dim, spacedim > &dof_handler, const unsigned int vertex_index, const unsigned int n)
unsigned int global_dof_index
Definition: types.h:100
void get_mg_dof_indices(const int level, std::vector< types::global_dof_index > &dof_indices, const unsigned int fe_index=DH::default_fe_index) const
#define Assert(cond, exc)
Definition: exceptions.h:299
void get_dof_values(const InputVector &values, Vector< number > &local_values) const
void set_mg_dof_indices(const std::vector< types::global_dof_index > &dof_indices)
void set_dof_values(const Vector< number > &local_values, OutputVector &values) const
iterator begin()
::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
bool fe_index_is_active(const unsigned int fe_index) const
types::global_dof_index vertex_dof_index(const unsigned int vertex, const unsigned int i, const unsigned int fe_index=DH::default_fe_index) const
::ExceptionBase & ExcImpossibleInDim(int arg1)
IteratorState::IteratorStates state() const
static unsigned int n_active_vertex_fe_indices(const ::hp::DoFHandler< dim, spacedim > &dof_handler, const unsigned int vertex_index)
const unsigned int dofs_per_cell
Definition: fe_base.h:271
const FiniteElement< DH::dimension, DH::space_dimension > & get_fe(const unsigned int fe_index) const
TriaIterator< TriaAccessor< structdim, dim, spacedim > > child(const unsigned int i) const
unsigned int nth_active_fe_index(const unsigned int n) const
void set_dof_index(const unsigned int i, const types::global_dof_index index, const unsigned int fe_index=DH::default_fe_index) const
TriaIterator< DoFCellAccessor< DH, level_dof_access > > parent() const
TriaIterator< DoFAccessor< DH::dimension-1, DH, level_dof_access > > face(const unsigned int i) const
std::size_t size() const
bool operator!=(const DoFAccessor< dim2, DH2, level_dof_access2 > &) const
::ExceptionBase & ExcNotImplemented()
::ExceptionBase & ExcNotInitialized()
Iterator reached end of container.
const unsigned int dofs_per_vertex
Definition: fe_base.h:214
const FiniteElement< DH::dimension, DH::space_dimension > & get_fe() const
void get_dof_indices(std::vector< types::global_dof_index > &dof_indices) const
static const unsigned int space_dimension
Definition: dof_handler.h:229
std::vector< types::global_dof_index > vertex_dofs_offsets
Definition: dof_handler.h:863
void get_dof_indices(std::vector< types::global_dof_index > &dof_indices, const unsigned int fe_index=DH::default_fe_index) const
const types::global_dof_index invalid_dof_index
Definition: types.h:217
::ExceptionBase & ExcInternalError()
std::vector< types::global_dof_index > vertex_dofs
Definition: dof_handler.h:840
void get_mg_dof_indices(std::vector< types::global_dof_index > &dof_indices) const
DoFCellAccessor(const Triangulation< DH::dimension, DH::space_dimension > *tria, const int level, const int index, const AccessorData *local_data)
void copy_from(const DoFAccessor< structdim, DH, level_dof_access2 > &a)
void get_dof_values(const VectorType &global_vector, ForwardIteratorInd local_indices_begin, ForwardIteratorVec local_vector_begin, ForwardIteratorVec local_vector_end) const
std::vector< types::global_dof_index > vertex_dofs
Definition: dof_handler.h:911
void set_vertex_dof_index(const unsigned int vertex, const unsigned int i, const types::global_dof_index index, const unsigned int fe_index=DH::default_fe_index) const
typename::internal::DoFHandler::Iterators< DH, level_dof_access >::quad_iterator quad(const unsigned int i) const
unsigned int n_active_fe_indices() const