Line data Source code
1 : /*
2 : * Copyright (c) 2022: G-CSC, Goethe University Frankfurt
3 : * Author: Lukas Larisch
4 : *
5 : * This file is part of UG4.
6 : *
7 : * UG4 is free software: you can redistribute it and/or modify it under the
8 : * terms of the GNU Lesser General Public License version 3 (as published by the
9 : * Free Software Foundation) with the following additional attribution
10 : * requirements (according to LGPL/GPL v3 §7):
11 : *
12 : * (1) The following notice must be displayed in the Appropriate Legal Notices
13 : * of covered and combined works: "Based on UG4 (www.ug4.org/license)".
14 : *
15 : * (2) The following notice must be displayed at a prominent place in the
16 : * terminal output of covered works: "Based on UG4 (www.ug4.org/license)".
17 : *
18 : * (3) The following bibliography is recommended for citation and must be
19 : * preserved in all covered files:
20 : * "Reiter, S., Vogel, A., Heppner, I., Rupp, M., and Wittum, G. A massively
21 : * parallel geometric multigrid solver on hierarchically distributed grids.
22 : * Computing and visualization in science 16, 4 (2013), 151-164"
23 : * "Vogel, A., Reiter, S., Rupp, M., Nägel, A., and Wittum, G. UG4 -- a novel
24 : * flexible software system for simulating pde based models on high performance
25 : * computers. Computing and visualization in science 16, 4 (2013), 165-179"
26 : *
27 : * This program is distributed in the hope that it will be useful,
28 : * but WITHOUT ANY WARRANTY; without even the implied warranty of
29 : * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
30 : * GNU Lesser General Public License for more details.
31 : */
32 :
33 : #ifndef __UG__LIB_DISC__ORDERING_STRATEGIES_ALGORITHMS_DIRECTIONAL_ORDERING__
34 : #define __UG__LIB_DISC__ORDERING_STRATEGIES_ALGORITHMS_DIRECTIONAL_ORDERING__
35 :
36 : #include <boost/graph/adjacency_list.hpp>
37 : #include <boost/graph/graph_traits.hpp>
38 : #include <boost/graph/properties.hpp>
39 :
40 : #include <boost/graph/cuthill_mckee_ordering.hpp>
41 :
42 : #include <set>
43 : #include <algorithm> //reverse
44 : #include <utility> //pair
45 :
46 : #include "lib_disc/domain.h"
47 : #include "lib_disc/function_spaces/grid_function.h"
48 : #include "lib_disc/spatial_disc/user_data/user_data.h"
49 :
50 : #include "lib_algebra/ordering_strategies/algorithms/IOrderingAlgorithm.h"
51 : #include "lib_algebra/ordering_strategies/algorithms/util.h"
52 :
53 : #include <assert.h>
54 : #include "common/error.h"
55 :
56 :
57 : namespace ug{
58 :
59 0 : bool CompareScalar(const std::pair<number, size_t> &p1,
60 : const std::pair<number, size_t> &p2)
61 : {
62 0 : return p1.first<p2.first;
63 : }
64 :
65 :
66 : template <typename TAlgebra, typename TDomain, typename O_t>
67 : class DirectionalOrdering : public IOrderingAlgorithm<TAlgebra, O_t>
68 : {
69 : public:
70 : typedef typename TAlgebra::matrix_type M_t;
71 : typedef typename TAlgebra::vector_type V_t;
72 : typedef IOrderingAlgorithm<TAlgebra, O_t> baseclass;
73 :
74 : typedef typename std::pair<MathVector<TDomain::dim>, size_t> Position_t;
75 : typedef typename std::pair<number, size_t> Scalar_t;
76 : typedef MathVector<TDomain::dim> small_vec_t;
77 : typedef GridFunction<TDomain, TAlgebra> GridFunc_t;
78 :
79 : typedef SmartPtr<UserData<MathVector<TDomain::dim>, TDomain::dim> > TSpUserData;
80 :
81 0 : DirectionalOrdering(){}
82 :
83 : /// clone constructor
84 0 : DirectionalOrdering( const DirectionalOrdering<TAlgebra, TDomain, O_t> &parent )
85 0 : : baseclass(){}
86 :
87 0 : SmartPtr<IOrderingAlgorithm<TAlgebra, O_t> > clone()
88 : {
89 0 : return make_sp(new DirectionalOrdering<TAlgebra, TDomain, O_t>(*this));
90 : }
91 :
92 0 : void compute(){
93 0 : m_vScalars.resize(m_vPositions.size());
94 : small_vec_t pos;
95 0 : for(size_t i = 0; i < m_vPositions.size(); ++i){
96 : pos = m_vPositions[i].first;
97 0 : m_vScalars[i].first = pos*(*m_dir); //scalar product
98 0 : m_vScalars[i].second = m_vPositions[i].second;
99 : }
100 :
101 0 : std::sort(m_vScalars.begin(), m_vScalars.end(), CompareScalar);
102 :
103 0 : for(size_t i = 0; i < m_vScalars.size(); ++i){
104 0 : o[m_vScalars[i].second] = i;
105 : }
106 :
107 : #ifdef UG_DEBUG
108 : check();
109 : #endif
110 0 : }
111 :
112 0 : void init(M_t* A, const V_t& V){
113 : const GridFunc_t* pGridF;
114 : unsigned n;
115 :
116 : try{
117 0 : if((pGridF = dynamic_cast<const GridFunc_t*>(&V)) == 0){
118 0 : UG_THROW(name() << "::init: No DoFDistribution specified.");
119 : }
120 :
121 : SmartPtr<DoFDistribution> dd = ((GridFunc_t*) pGridF)->dof_distribution();
122 :
123 : n = dd->num_indices();
124 :
125 0 : if(n != A->num_rows ()){
126 0 : UG_THROW(name() << "::init: #indices != #rows");
127 : }
128 :
129 0 : o.resize(n);
130 0 : ExtractPositions(pGridF->domain(), dd, m_vPositions);
131 : }
132 0 : catch(...){
133 0 : throw;
134 : }
135 :
136 : #ifdef UG_ENABLE_DEBUG_LOGS
137 : UG_LOG("Using " << name() << " (direction " << *m_dir << ")\n");
138 : #endif
139 0 : }
140 :
141 0 : void init(M_t* A){
142 0 : UG_THROW(name() << "::init: Cannot initialize smoother without a geometry. Specify the 2nd argument for init!");
143 : }
144 :
145 0 : void init(M_t* A, const V_t&, const O_t& inv_map){
146 0 : UG_THROW(name() << "::init: induced subgraph version not implemented yet!");
147 : }
148 :
149 0 : void init(M_t* A, const O_t& inv_map){
150 0 : UG_THROW(name() << "::init: induced subgraph version not implemented yet!");
151 : }
152 :
153 0 : void check(){
154 0 : if(!is_permutation(o)){
155 0 : UG_THROW(name() << "::check: Not a permutation!");
156 : }
157 0 : }
158 :
159 0 : O_t& ordering(){
160 0 : return o;
161 : }
162 :
163 0 : void set_direction(small_vec_t *dir){
164 0 : m_dir = dir;
165 0 : }
166 :
167 0 : virtual const char* name() const {return "DirectionalOrdering";}
168 :
169 : private:
170 : O_t o;
171 :
172 : small_vec_t *m_dir;
173 :
174 : std::vector<Position_t> m_vPositions;
175 : std::vector<Scalar_t> m_vScalars;
176 : };
177 :
178 : } //namespace
179 :
180 :
181 : #endif //guard
|
