00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012 #ifndef __deal2__hp_dof_objects_h
00013 #define __deal2__hp_dof_objects_h
00014
00015 #include <deal.II/base/config.h>
00016 #include <deal.II/hp/fe_collection.h>
00017 #include <deal.II/hp/dof_handler.h>
00018
00019 #include <vector>
00020
00021 DEAL_II_NAMESPACE_OPEN
00022
00023 namespace internal
00024 {
00025 namespace hp
00026 {
00027
00095 template <int dim>
00096 class DoFObjects
00097 {
00098 public:
00104 std::vector<unsigned int> dof_offsets;
00105
00112 std::vector<unsigned int> dofs;
00113
00139 template <int dimm, int spacedim>
00140 void
00141 set_dof_index (const ::hp::DoFHandler<dimm,spacedim> &dof_handler,
00142 const unsigned int obj_index,
00143 const unsigned int fe_index,
00144 const unsigned int local_index,
00145 const unsigned int global_index,
00146 const unsigned int obj_level);
00147
00171 template <int dimm, int spacedim>
00172 unsigned int
00173 get_dof_index (const ::hp::DoFHandler<dimm,spacedim> &dof_handler,
00174 const unsigned int obj_index,
00175 const unsigned int fe_index,
00176 const unsigned int local_index,
00177 const unsigned int obj_level) const;
00178
00202 template <int dimm, int spacedim>
00203 unsigned int
00204 n_active_fe_indices (const ::hp::DoFHandler<dimm,spacedim> &dof_handler,
00205 const unsigned int obj_index) const;
00206
00212 template <int dimm, int spacedim>
00213 unsigned int
00214 nth_active_fe_index (const ::hp::DoFHandler<dimm,spacedim> &dof_handler,
00215 const unsigned int obj_level,
00216 const unsigned int obj_index,
00217 const unsigned int n) const;
00218
00225 template <int dimm, int spacedim>
00226 bool
00227 fe_index_is_active (const ::hp::DoFHandler<dimm,spacedim> &dof_handler,
00228 const unsigned int obj_index,
00229 const unsigned int fe_index,
00230 const unsigned int obj_level) const;
00231
00237 std::size_t memory_consumption () const;
00238 };
00239
00240
00241
00242
00243 template <int dim>
00244 template <int dimm, int spacedim>
00245 inline
00246 unsigned int
00247 DoFObjects<dim>::
00248 get_dof_index (const ::hp::DoFHandler<dimm,spacedim> &dof_handler,
00249 const unsigned int obj_index,
00250 const unsigned int fe_index,
00251 const unsigned int local_index,
00252 const unsigned int obj_level) const
00253 {
00254 Assert ((fe_index != ::hp::DoFHandler<dimm,spacedim>::default_fe_index),
00255 ExcMessage ("You need to specify a FE index when working "
00256 "with hp DoFHandlers"));
00257 Assert (&dof_handler != 0,
00258 ExcMessage ("No DoFHandler is specified for this iterator"));
00259 Assert (&dof_handler.get_fe() != 0,
00260 ExcMessage ("No finite element collection is associated with "
00261 "this DoFHandler"));
00262 Assert (fe_index < dof_handler.get_fe().size(),
00263 ExcIndexRange (fe_index, 0, dof_handler.get_fe().size()));
00264 Assert (local_index <
00265 dof_handler.get_fe()[fe_index].template n_dofs_per_object<dim>(),
00266 ExcIndexRange(local_index, 0,
00267 dof_handler.get_fe()[fe_index]
00268 .template n_dofs_per_object<dim>()));
00269 Assert (obj_index < dof_offsets.size(),
00270 ExcIndexRange (obj_index, 0, dof_offsets.size()));
00271
00272
00273
00274
00275 Assert (dof_offsets[obj_index] != numbers::invalid_unsigned_int,
00276 ExcMessage ("You are trying to access degree of freedom "
00277 "information for an object on which no such "
00278 "information is available"));
00279
00280 if (dim == dimm)
00281 {
00282
00283
00284
00285
00286
00287
00288 Assert (fe_index == dof_handler.levels[obj_level]->active_fe_indices[obj_index],
00289 ExcMessage ("FE index does not match that of the present cell"));
00290 return dofs[dof_offsets[obj_index]+local_index];
00291 }
00292 else
00293 {
00294
00295
00296
00297
00298
00299
00300
00301
00302
00303
00304
00305
00306
00307 const unsigned int starting_offset = dof_offsets[obj_index];
00308 const unsigned int *pointer = &dofs[starting_offset];
00309 while (true)
00310 {
00311 Assert (*pointer != numbers::invalid_unsigned_int,
00312 ExcInternalError());
00313 if (*pointer == fe_index)
00314 return *(pointer + 1 + local_index);
00315 else
00316 pointer += dof_handler.get_fe()[*pointer]
00317 .template n_dofs_per_object<dim>() + 1;
00318 }
00319 }
00320 }
00321
00322
00323
00324 template <int dim>
00325 template <int dimm, int spacedim>
00326 inline
00327 void
00328 DoFObjects<dim>::
00329 set_dof_index (const ::hp::DoFHandler<dimm,spacedim> &dof_handler,
00330 const unsigned int obj_index,
00331 const unsigned int fe_index,
00332 const unsigned int local_index,
00333 const unsigned int global_index,
00334 const unsigned int obj_level)
00335 {
00336 Assert ((fe_index != ::hp::DoFHandler<dimm,spacedim>::default_fe_index),
00337 ExcMessage ("You need to specify a FE index when working "
00338 "with hp DoFHandlers"));
00339 Assert (&dof_handler != 0,
00340 ExcMessage ("No DoFHandler is specified for this iterator"));
00341 Assert (&dof_handler.get_fe() != 0,
00342 ExcMessage ("No finite element collection is associated with "
00343 "this DoFHandler"));
00344 Assert (fe_index < dof_handler.get_fe().size(),
00345 ExcIndexRange (fe_index, 0, dof_handler.get_fe().size()));
00346 Assert (local_index <
00347 dof_handler.get_fe()[fe_index].template n_dofs_per_object<dim>(),
00348 ExcIndexRange(local_index, 0,
00349 dof_handler.get_fe()[fe_index]
00350 .template n_dofs_per_object<dim>()));
00351 Assert (obj_index < dof_offsets.size(),
00352 ExcIndexRange (obj_index, 0, dof_offsets.size()));
00353
00354
00355
00356
00357 Assert (dof_offsets[obj_index] != numbers::invalid_unsigned_int,
00358 ExcMessage ("You are trying to access degree of freedom "
00359 "information for an object on which no such "
00360 "information is available"));
00361
00362 if (dim == dimm)
00363 {
00364
00365
00366
00367
00368
00369
00370 Assert (fe_index == dof_handler.levels[obj_level]->active_fe_indices[obj_index],
00371 ExcMessage ("FE index does not match that of the present cell"));
00372 dofs[dof_offsets[obj_index]+local_index] = global_index;
00373 }
00374 else
00375 {
00376
00377
00378
00379
00380
00381
00382
00383
00384
00385
00386
00387
00388
00389 const unsigned int starting_offset = dof_offsets[obj_index];
00390 unsigned int *pointer = &dofs[starting_offset];
00391 while (true)
00392 {
00393 Assert (*pointer != numbers::invalid_unsigned_int,
00394 ExcInternalError());
00395 if (*pointer == fe_index)
00396 {
00397 *(pointer + 1 + local_index) = global_index;
00398 return;
00399 }
00400 else
00401 pointer += dof_handler.get_fe()[*pointer]
00402 .template n_dofs_per_object<dim>() + 1;
00403 }
00404 }
00405 }
00406
00407
00408
00409 template <int dim>
00410 template <int dimm, int spacedim>
00411 inline
00412 unsigned int
00413 DoFObjects<dim>::
00414 n_active_fe_indices (const ::hp::DoFHandler<dimm,spacedim> &dof_handler,
00415 const unsigned int obj_index) const
00416 {
00417 Assert (dim <= dimm, ExcInternalError());
00418 Assert (&dof_handler != 0,
00419 ExcMessage ("No DoFHandler is specified for this iterator"));
00420 Assert (&dof_handler.get_fe() != 0,
00421 ExcMessage ("No finite element collection is associated with "
00422 "this DoFHandler"));
00423 Assert (obj_index < dof_offsets.size(),
00424 ExcIndexRange (obj_index, 0, dof_offsets.size()));
00425
00426
00427
00428
00429 if (dof_offsets[obj_index] == numbers::invalid_unsigned_int)
00430 return 0;
00431
00432
00433
00434
00435
00436 if (dim == dimm)
00437 return 1;
00438 else
00439 {
00440
00441
00442
00443
00444
00445
00446
00447
00448
00449
00450
00451
00452 const unsigned int starting_offset = dof_offsets[obj_index];
00453 const unsigned int *pointer = &dofs[starting_offset];
00454 unsigned int counter = 0;
00455 while (true)
00456 {
00457 if (*pointer == numbers::invalid_unsigned_int)
00458
00459 return counter;
00460 else
00461 {
00462 ++counter;
00463 pointer += dof_handler.get_fe()[*pointer]
00464 .template n_dofs_per_object<dim>() + 1;
00465 }
00466 }
00467 }
00468 }
00469
00470
00471
00472 template <int dim>
00473 template <int dimm, int spacedim>
00474 inline
00475 unsigned int
00476 DoFObjects<dim>::
00477 nth_active_fe_index (const ::hp::DoFHandler<dimm,spacedim> &dof_handler,
00478 const unsigned int obj_level,
00479 const unsigned int obj_index,
00480 const unsigned int n) const
00481 {
00482 Assert (dim <= dimm, ExcInternalError());
00483 Assert (&dof_handler != 0,
00484 ExcMessage ("No DoFHandler is specified for this iterator"));
00485 Assert (&dof_handler.get_fe() != 0,
00486 ExcMessage ("No finite element collection is associated with "
00487 "this DoFHandler"));
00488 Assert (obj_index < dof_offsets.size(),
00489 ExcIndexRange (obj_index, 0, dof_offsets.size()));
00490
00491
00492
00493
00494 Assert (dof_offsets[obj_index] != numbers::invalid_unsigned_int,
00495 ExcMessage ("You are trying to access degree of freedom "
00496 "information for an object on which no such "
00497 "information is available"));
00498
00499 if (dim == dimm)
00500 {
00501
00502
00503
00504 Assert (n == 0, ExcIndexRange (n, 0, 1));
00505
00506 return dof_handler.levels[obj_level]->active_fe_indices[obj_index];
00507 }
00508 else
00509 {
00510 Assert (n < n_active_fe_indices(dof_handler, obj_index),
00511 ExcIndexRange (n, 0,
00512 n_active_fe_indices(dof_handler, obj_index)));
00513
00514
00515
00516
00517
00518
00519
00520
00521
00522
00523
00524
00525
00526
00527 const unsigned int starting_offset = dof_offsets[obj_index];
00528 const unsigned int *pointer = &dofs[starting_offset];
00529 unsigned int counter = 0;
00530 while (true)
00531 {
00532 Assert (*pointer != numbers::invalid_unsigned_int,
00533 ExcInternalError());
00534
00535 const unsigned int fe_index = *pointer;
00536
00537 Assert (fe_index < dof_handler.get_fe().size(),
00538 ExcInternalError());
00539
00540 if (counter == n)
00541 return fe_index;
00542
00543 ++counter;
00544 pointer += dof_handler.get_fe()[fe_index]
00545 .template n_dofs_per_object<dim>() + 1;
00546 }
00547 }
00548 }
00549
00550
00551
00552 template <int dim>
00553 template <int dimm, int spacedim>
00554 inline
00555 bool
00556 DoFObjects<dim>::
00557 fe_index_is_active (const ::hp::DoFHandler<dimm,spacedim> &dof_handler,
00558 const unsigned int obj_index,
00559 const unsigned int fe_index,
00560 const unsigned int obj_level) const
00561 {
00562 Assert (&dof_handler != 0,
00563 ExcMessage ("No DoFHandler is specified for this iterator"));
00564 Assert (&dof_handler.get_fe() != 0,
00565 ExcMessage ("No finite element collection is associated with "
00566 "this DoFHandler"));
00567 Assert (obj_index < dof_offsets.size(),
00568 ExcIndexRange (obj_index, 0, dof_offsets.size()));
00569 Assert ((fe_index != ::hp::DoFHandler<dimm,spacedim>::default_fe_index),
00570 ExcMessage ("You need to specify a FE index when working "
00571 "with hp DoFHandlers"));
00572 Assert (fe_index < dof_handler.get_fe().size(),
00573 ExcIndexRange (fe_index, 0, dof_handler.get_fe().size()));
00574
00575
00576
00577
00578 Assert (dof_offsets[obj_index] != numbers::invalid_unsigned_int,
00579 ExcMessage ("You are trying to access degree of freedom "
00580 "information for an object on which no such "
00581 "information is available"));
00582
00583 if (dim == dimm)
00584 {
00585
00586
00587
00588
00589
00590 Assert (obj_index < dof_handler.levels[obj_level]->active_fe_indices.size(),
00591 ExcInternalError());
00592 return (fe_index == dof_handler.levels[obj_level]->active_fe_indices[obj_index]);
00593 }
00594 else
00595 {
00596
00597
00598
00599
00600
00601
00602
00603
00604
00605
00606
00607
00608
00609 const unsigned int starting_offset = dof_offsets[obj_index];
00610 const unsigned int *pointer = &dofs[starting_offset];
00611 while (true)
00612 {
00613 if (*pointer == numbers::invalid_unsigned_int)
00614
00615 return false;
00616 else
00617 if (*pointer == fe_index)
00618 return true;
00619 else
00620 pointer += dof_handler.get_fe()[*pointer]
00621 .template n_dofs_per_object<dim>() + 1;
00622 }
00623 }
00624 }
00625
00626 }
00627 }
00628
00629 DEAL_II_NAMESPACE_CLOSE
00630
00631 #endif
00632