00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013 #ifndef __deal2__mesh_worker_local_results_h
00014 #define __deal2__mesh_worker_local_results_h
00015
00016 #include <deal.II/base/config.h>
00017 #include <deal.II/base/std_cxx1x/function.h>
00018 #include <deal.II/base/geometry_info.h>
00019 #include <deal.II/lac/matrix_block.h>
00020 #include <deal.II/lac/block_vector.h>
00021 #include <deal.II/meshworker/vector_selector.h>
00022
00023 DEAL_II_NAMESPACE_OPEN
00024
00025 class BlockIndices;
00026 template<int,int> class DoFHandler;
00027 template<int,int> class MGDoFHandler;
00028
00174 namespace MeshWorker
00175 {
00208 template <typename number>
00209 class LocalResults
00210 {
00211 public:
00219 unsigned int n_values () const;
00220
00229 unsigned int n_vectors () const;
00230
00234 unsigned int n_matrices () const;
00235
00240 unsigned int n_quadrature_points() const;
00241
00247 unsigned int n_quadrature_values() const;
00248
00253 number& value(unsigned int i);
00254
00259 number value(unsigned int i) const;
00260
00264 BlockVector<number>& vector(unsigned int i);
00265
00269 const BlockVector<number>& vector(unsigned int i) const;
00270
00280 MatrixBlock<FullMatrix<number> >& matrix(unsigned int i, bool external = false);
00281
00291 const MatrixBlock<FullMatrix<number> >& matrix(unsigned int i, bool external = false) const;
00292
00302 Table<2, number>& quadrature_values();
00303
00308 number& quadrature_value(unsigned int k, unsigned int i);
00309
00314 number quadrature_value(unsigned int k, unsigned int i) const;
00315
00324 void initialize_numbers(const unsigned int n);
00333 void initialize_vectors(const unsigned int n);
00347 void initialize_matrices(unsigned int n, bool both);
00348
00363 template <class MATRIX>
00364 void initialize_matrices(const MatrixBlockVector<MATRIX>& matrices,
00365 bool both);
00366
00381 template <class MATRIX>
00382 void initialize_matrices(const MGMatrixBlockVector<MATRIX>& matrices,
00383 bool both);
00384
00390 void initialize_quadrature(unsigned int np, unsigned int nv);
00391
00402 void reinit(const BlockIndices& local_sizes);
00403
00404 template <class STREAM>
00405 void print_debug(STREAM& os) const;
00406
00410 std::size_t memory_consumption () const;
00411
00412 private:
00418 void initialize_local(MatrixBlock<FullMatrix<number> >& M,
00419 const unsigned int row,
00420 const unsigned int col);
00421
00427 std::vector<number> J;
00428
00435 std::vector<BlockVector<number> > R;
00436
00444 std::vector<MatrixBlock<FullMatrix<number> > > M1;
00445
00456 std::vector<MatrixBlock<FullMatrix<number> > > M2;
00457
00464 Table<2, number> quadrature_data;
00465 };
00466
00467
00468
00469 template <typename number>
00470 inline void
00471 LocalResults<number>::initialize_numbers(unsigned int n)
00472 {
00473 J.resize(n);
00474 }
00475
00476
00477 template <typename number>
00478 inline void
00479 LocalResults<number>::initialize_vectors(const unsigned int n)
00480 {
00481 R.resize(n);
00482 }
00483
00484
00485 template <typename number>
00486 template <class MATRIX>
00487 inline void
00488 LocalResults<number>::initialize_matrices(
00489 const MatrixBlockVector<MATRIX>& matrices,
00490 bool both)
00491 {
00492 M1.resize(matrices.size());
00493 if (both)
00494 M2.resize(matrices.size());
00495 for (unsigned int i=0;i<matrices.size();++i)
00496 {
00497 const unsigned int row = matrices.block(i).row;
00498 const unsigned int col = matrices.block(i).column;
00499
00500 M1[i].row = row;
00501 M1[i].column = col;
00502 if (both)
00503 {
00504 M2[i].row = row;
00505 M2[i].column = col;
00506 }
00507 }
00508 }
00509
00510
00511 template <typename number>
00512 template <class MATRIX>
00513 inline void
00514 LocalResults<number>::initialize_matrices(
00515 const MGMatrixBlockVector<MATRIX>& matrices,
00516 bool both)
00517 {
00518 M1.resize(matrices.size());
00519 if (both)
00520 M2.resize(matrices.size());
00521 for (unsigned int i=0;i<matrices.size();++i)
00522 {
00523 const MGLevelObject<MatrixBlock<MATRIX> >& o = matrices.block(i);
00524 const unsigned int row = o[o.min_level()].row;
00525 const unsigned int col = o[o.min_level()].column;
00526
00527 M1[i].row = row;
00528 M1[i].column = col;
00529 if (both)
00530 {
00531 M2[i].row = row;
00532 M2[i].column = col;
00533 }
00534 }
00535 }
00536
00537
00538 template <typename number>
00539 inline void
00540 LocalResults<number>::initialize_matrices(const unsigned int n,
00541 const bool both)
00542 {
00543 M1.resize(n);
00544 if (both)
00545 M2.resize(n);
00546 for (unsigned int i=0;i<n;++i)
00547 {
00548 M1[i].row = 0;
00549 M1[i].column = 0;
00550 if (both)
00551 {
00552 M2[i].row = 0;
00553 M2[i].column = 0;
00554 }
00555 }
00556 }
00557
00558
00559 template <typename number>
00560 inline void
00561 LocalResults<number>::initialize_quadrature(unsigned int np, unsigned int nv)
00562 {
00563 quadrature_data.reinit(np, nv);
00564 }
00565
00566
00567 template <typename number>
00568 inline
00569 unsigned int
00570 LocalResults<number>::n_values() const
00571 {
00572 return J.size();
00573 }
00574
00575
00576 template <typename number>
00577 inline
00578 unsigned int
00579 LocalResults<number>::n_vectors() const
00580 {
00581 return R.size();
00582 }
00583
00584
00585 template <typename number>
00586 inline
00587 unsigned int
00588 LocalResults<number>::n_matrices() const
00589 {
00590 return M1.size();
00591 }
00592
00593
00594 template <typename number>
00595 inline
00596 unsigned int
00597 LocalResults<number>::n_quadrature_points() const
00598 {
00599 return quadrature_data.n_rows();
00600 }
00601
00602
00603 template <typename number>
00604 inline
00605 unsigned int
00606 LocalResults<number>::n_quadrature_values() const
00607 {
00608 return quadrature_data.n_cols();
00609 }
00610
00611
00612 template <typename number>
00613 inline
00614 number&
00615 LocalResults<number>::value(unsigned int i)
00616 {
00617 AssertIndexRange(i,J.size());
00618 return J[i];
00619 }
00620
00621
00622 template <typename number>
00623 inline
00624 BlockVector<number>&
00625 LocalResults<number>::vector(unsigned int i)
00626 {
00627 AssertIndexRange(i,R.size());
00628 return R[i];
00629 }
00630
00631
00632 template <typename number>
00633 inline
00634 MatrixBlock<FullMatrix<number> >&
00635 LocalResults<number>::matrix(unsigned int i, bool external)
00636 {
00637 if (external)
00638 {
00639 AssertIndexRange(i,M2.size());
00640 return M2[i];
00641 }
00642 AssertIndexRange(i,M1.size());
00643 return M1[i];
00644 }
00645
00646
00647 template <typename number>
00648 inline
00649 number&
00650 LocalResults<number>::quadrature_value(unsigned int k, unsigned int i)
00651 {
00652 return quadrature_data(k,i);
00653 }
00654
00655
00656 template <typename number>
00657 inline
00658 Table<2, number>&
00659 LocalResults<number>::quadrature_values()
00660 {
00661 return quadrature_data;
00662 }
00663
00664
00665 template <typename number>
00666 inline
00667 number
00668 LocalResults<number>::value(unsigned int i) const
00669 {
00670 AssertIndexRange(i,J.size());
00671 return J[i];
00672 }
00673
00674
00675 template <typename number>
00676 inline
00677 const BlockVector<number>&
00678 LocalResults<number>::vector(unsigned int i) const
00679 {
00680 AssertIndexRange(i,R.size());
00681 return R[i];
00682 }
00683
00684
00685 template <typename number>
00686 inline
00687 const MatrixBlock<FullMatrix<number> >&
00688 LocalResults<number>::matrix(unsigned int i, bool external) const
00689 {
00690 if (external)
00691 {
00692 AssertIndexRange(i,M2.size());
00693 return M2[i];
00694 }
00695 AssertIndexRange(i,M1.size());
00696 return M1[i];
00697 }
00698
00699
00700 template <typename number>
00701 inline
00702 number
00703 LocalResults<number>::quadrature_value(unsigned int k, unsigned int i) const
00704 {
00705 return quadrature_data(k,i);
00706 }
00707
00708
00709 template <typename number>
00710 template <class STREAM>
00711 void
00712 LocalResults<number>::print_debug(STREAM& os) const
00713 {
00714 os << "J: " << J.size() << std::endl;
00715 os << "R: " << R.size() << std::endl;
00716 for (unsigned int i=0;i<R.size();++i)
00717 {
00718 os << " " << R[i].n_blocks() << " -";
00719 for (unsigned int j=0;j<R[i].n_blocks();++j)
00720 os << ' ' << R[i].block(j).size();
00721 os << std::endl;
00722 }
00723 os << "M: " << M1.size() << " face " << M2.size() << std::endl;
00724 for (unsigned int i=0;i<M1.size();++i)
00725 {
00726 os << " " << M1[i].row << "," << M1[i].column
00727 << " " << M1[i].matrix.m() << 'x' << M1[i].matrix.n();
00728 if (i < M2.size())
00729 os << " face " << M2[i].row << "," << M2[i].column
00730 << " " << M2[i].matrix.m() << 'x' << M2[i].matrix.n();
00731 os << std::endl;
00732 }
00733 }
00734
00735 }
00736
00737
00738 DEAL_II_NAMESPACE_CLOSE
00739
00740 #endif
00741