← Back to branch summary

~ubuntu-branches/ubuntu/trusty/rheolef/trusty

« back to all changes in this revision

Viewing changes to nfem/plib/form_concat.h

  • Committer: Package Import Robot
  • Author(s): Pierre Saramito
  • Date: 2012-04-06 09:12:21 UTC
  • mfrom: (1.1.5)
  • Revision ID: package-import@ubuntu.com-20120406091221-m58me99p1nxqui49
Tags: 6.0-1
http://bugs.debian.org/667356
http://bugs.debian.org/661689
* New upstream release 6.0 (major changes):
  - massively distributed and parallel support
  - full FEM characteristic method (Lagrange-Gakerkin method) support
  - enhanced users documentation 
  - source code supports g++-4.7 (closes: #667356)
* debian/control: dependencies for MPI distributed solvers added
* debian/rules: build commands simplified
* debian/librheolef-dev.install: man1/* to man9/* added
* debian/changelog: package description rewritted (closes: #661689)

Show diffs side-by-side

added added

removed removed

Lines of Context:
nfem/plib/form_concat.h
 
1
#ifndef _RHEOLEF_FORM_CONCAT_H
 
2
#define _RHEOLEF_FORM_CONCAT_H
 
3
///
 
4
/// This file is part of Rheolef.
 
5
///
 
6
/// Copyright (C) 2000-2009 Pierre Saramito <Pierre.Saramito@imag.fr>
 
7
///
 
8
/// Rheolef is free software; you can redistribute it and/or modify
 
9
/// it under the terms of the GNU General Public License as published by
 
10
/// the Free Software Foundation; either version 2 of the License, or
 
11
/// (at your option) any later version.
 
12
///
 
13
/// Rheolef is distributed in the hope that it will be useful,
 
14
/// but WITHOUT ANY WARRANTY; without even the implied warranty of
 
15
/// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 
16
/// GNU General Public License for more details.
 
17
///
 
18
/// You should have received a copy of the GNU General Public License
 
19
/// along with Rheolef; if not, write to the Free Software
 
20
/// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
 
21
///
 
22
/// =========================================================================
 
23
// build form from initializer list (c++ 2011)
 
24
//
 
25
#include "rheolef/form.h"
 
26
#include "rheolef/csr_concat.h"
 
27
 
 
28
namespace rheolef {
 
29
 
 
30
// =========================================================================
 
31
// 1rst case : one-line matrix initializer
 
32
//  A = {a, b};                 // matrix & vector
 
33
// =========================================================================
 
34
template <class T, class M>
 
35
struct form_concat_value {
 
36
// typedef:
 
37
 typedef enum { scalar, field, field_transpose, form} variant_type;
 
38
// allocators:
 
39
 form_concat_value (const T& x)               : s(x), v(), m(),  variant(scalar) {}
 
40
 form_concat_value (const form_basic<T,M>& x) : s(),  v(), m(x), variant(form) {}
 
41
// io/debug:
 
42
 friend std::ostream& operator<< (std::ostream& o, const form_concat_value<T,M>& x) {
 
43
       if (x.variant == scalar) return o << "s";
 
44
  else if (x.variant == field)  return o << "f";
 
45
  else if (x.variant == field_transpose) return o << "ft";
 
46
  else return o << "m";
 
47
 }
 
48
// data:
 
49
public:
 
50
 T                 s;
 
51
 field_basic<T,M>  v;
 
52
 form_basic<T,M>   m;
 
53
 variant_type      variant;
 
54
};
 
55
template <class T, class M>
 
56
struct form_concat_line {
 
57
// typedef:
 
58
 typedef typename form_basic<T,M>::size_type   size_type;
 
59
 typedef form_concat_value<T,M>                value_type;
 
60
 typedef typename std::list<value_type>::const_iterator const_iterator;
 
61
// allocators:
 
62
 form_concat_line () : _l() {}
 
63
 form_concat_line (const std::initializer_list<value_type>& il) : _l() {
 
64
#ifdef _RHEOLEF_HAVE_STD_INITIALIZER_ITERATOR
 
65
    typedef typename std::initializer_list<value_type>::const_iterator const_iterator;
 
66
#else // _RHEOLEF_HAVE_STD_INITIALIZER_ITERATOR
 
67
    typedef const value_type* const_iterator;
 
68
#endif // _RHEOLEF_HAVE_STD_INITIALIZER_ITERATOR
 
69
    for(const_iterator iter = il.begin(); iter != il.end(); ++iter) {
 
70
        _l.push_back(*iter);
 
71
    }
 
72
 }
 
73
// accessors:
 
74
 const_iterator begin() const { return _l.begin(); }
 
75
 const_iterator end()   const { return _l.end(); }
 
76
 
 
77
 friend std::ostream& operator<< (std::ostream& o, const form_concat_line<T,M>& x) {
 
78
    std::cout << "{";
 
79
    for(typename std::list<value_type>::const_iterator iter = x._l.begin(); iter != x._l.end(); ++iter) {
 
80
        std::cout << *iter << " ";
 
81
    }
 
82
    return std::cout << "}";
 
83
 }
 
84
// internals:
 
85
 
 
86
 void build_form_pass0 (std::vector<std::pair<bool,space_basic<T,M> > >& l_Xh, space_basic<T,M>& Yh, size_t i_comp = 0) const;
 
87
 static void build_first_space (const std::vector<std::pair<bool,space_basic<T,M> > >& l_Xh, space_basic<T,M>& Xh);
 
88
 void build_form_pass1 (space_basic<T,M>& Xh, space_basic<T,M>& Yh) const;
 
89
 form_basic<T,M> build_form_pass2 (const space_basic<T,M>& Xh, const space_basic<T,M>& Yh) const;
 
90
 form_basic<T,M> build_form () const;
 
91
 
 
92
// data:
 
93
protected:
 
94
 std::list<value_type> _l;
 
95
};
 
96
template <class T, class M>
 
97
void
 
98
form_concat_line<T,M>::build_form_pass0 (std::vector<std::pair<bool,space_basic<T,M> > >& l_Xh, space_basic<T,M>& Yh, size_t i_comp) const
 
99
{
 
100
  // ------------------------------------------------------------
 
101
  // pass 0 : lazy first space computation, compute second space
 
102
  // ------------------------------------------------------------
 
103
  size_t j_comp = 0;
 
104
  bool have_Yh = false;
 
105
  typename std::vector<std::pair<bool,space_basic<T,M> > >::iterator xh_iter = l_Xh.begin();
 
106
  for (typename std::list<value_type>::const_iterator iter = _l.begin(); iter != _l.end(); ++iter, xh_iter++, j_comp++) {
 
107
    const value_type& x = *iter;
 
108
    switch (x.variant) {
 
109
      case form_concat_value<T,M>::scalar: {
 
110
        check_macro (x.s == 0, "unsupported non-nul scalar `"<<x.s<<"' in form concatenation"
 
111
                   << " at ("<<i_comp<<","<<j_comp<<")");
 
112
        break;
 
113
      }
 
114
      case form_concat_value<T,M>::form: {
 
115
        if (!(*xh_iter).first) {
 
116
          (*xh_iter).first  = true;
 
117
          (*xh_iter).second = x.m.get_first_space();
 
118
        } else {
 
119
          check_macro (x.m.get_first_space() == (*xh_iter).second, "form initializer: invalid second space `"
 
120
                << x.m.get_first_space().stamp() << "': expect `" << (*xh_iter).second.stamp() << "'"
 
121
                << " at ("<<i_comp<<","<<j_comp<<")");
 
122
        }
 
123
        if (!have_Yh) {
 
124
          have_Yh = true;
 
125
          Yh      = x.m.get_second_space();
 
126
        } else {
 
127
          check_macro (x.m.get_second_space() == Yh, "form initializer: invalid second space `"
 
128
                << x.m.get_second_space().stamp() << "': expect `" << Yh.stamp() << "'"
 
129
                << " at ("<<i_comp<<","<<j_comp<<")");
 
130
        }
 
131
        break;
 
132
      }
 
133
      default: error_macro ("non-form or scalar concatenation not yet supported"
 
134
                << " at ("<<i_comp<<","<<j_comp<<")");
 
135
    }
 
136
  }
 
137
  check_macro (have_Yh, "form concatenation: "<<i_comp<<"th row space remains undefined");
 
138
}
 
139
template <class T, class M>
 
140
void
 
141
form_concat_line<T,M>::build_first_space (const std::vector<std::pair<bool,space_basic<T,M> > >& l_Xh, space_basic<T,M>& Xh)
 
142
{
 
143
  // ------------------------------------------------------------
 
144
  // pass 0b : first space computation
 
145
  // ------------------------------------------------------------
 
146
  space_mult_list<T,M> sml_X;
 
147
  size_t j_comp = 0;
 
148
  for (typename std::vector<std::pair<bool,space_basic<T,M> > >::const_iterator
 
149
        xh_iter = l_Xh.begin(),
 
150
        xh_last = l_Xh.end(); xh_iter != xh_last; xh_iter++, j_comp++) {
 
151
    check_macro ((*xh_iter).first, "form concatenation: "<<j_comp<<"th column space remains undefined");
 
152
    sml_X *= (*xh_iter).second;
 
153
  }
 
154
  Xh = space_basic<T,M>(sml_X);
 
155
}
 
156
template <class T, class M>
 
157
void
 
158
form_concat_line<T,M>::build_form_pass1 (space_basic<T,M>& Xh, space_basic<T,M>& Yh) const
 
159
{
 
160
  // --------------------------------
 
161
  // pass 1 : both spaces computation
 
162
  // --------------------------------
 
163
  std::vector<std::pair<bool,space_basic<T,M> > > l_Xh (_l.size(), std::pair<bool,space_basic<T,M> >(false, space_basic<T,M>()));
 
164
  build_form_pass0 (l_Xh, Yh);
 
165
  build_first_space (l_Xh, Xh);
 
166
}
 
167
template <class T, class M>
 
168
form_basic<T,M>
 
169
form_concat_line<T,M>::build_form_pass2 (const space_basic<T,M>& Xh, const space_basic<T,M>& Yh) const
 
170
{
 
171
  // -----------------------
 
172
  // pass 2 : compute values
 
173
  // -----------------------
 
174
  form_basic<T,M> a (Xh, Yh);
 
175
  csr_concat_line<T,M> uu, ub, bu, bb;
 
176
  for(typename std::list<value_type>::const_iterator iter = _l.begin(); iter != _l.end(); ++iter) {
 
177
    const value_type& x = *iter;
 
178
    switch (x.variant) {
 
179
      case form_concat_value<T,M>::form: {
 
180
        uu.push_back (x.m.uu());
 
181
        ub.push_back (x.m.ub());
 
182
        bu.push_back (x.m.bu());
 
183
        bb.push_back (x.m.bb());
 
184
        break;
 
185
      }
 
186
      default: error_macro ("non-form concatenation not yet supported");
 
187
    }
 
188
  }
 
189
  a.set_uu() = uu.build_csr();
 
190
  a.set_ub() = ub.build_csr();
 
191
  a.set_bu() = bu.build_csr();
 
192
  a.set_bb() = bb.build_csr();
 
193
  return a;
 
194
}
 
195
template <class T, class M>
 
196
form_basic<T,M>
 
197
form_concat_line<T,M>::build_form() const
 
198
{
 
199
  space_basic<T,M> Xh, Yh;
 
200
  build_form_pass1 (Xh, Yh);
 
201
  return build_form_pass2 (Xh, Yh);
 
202
}
 
203
// -------------------------------
 
204
// form cstor from std::initializer
 
205
// -------------------------------
 
206
#ifdef _RHEOLEF_HAVE_STD_INITIALIZER_LIST
 
207
template <class T, class M>
 
208
inline
 
209
form_basic<T,M>::form_basic (const std::initializer_list<form_concat_value<T,M> >& init_list)
 
210
{
 
211
  form_concat_line<T,M> cc (init_list);
 
212
  form_basic<T,M>::operator= (cc.build_form());
 
213
}
 
214
#endif // _RHEOLEF_HAVE_STD_INITIALIZER_LIST
 
215
// =========================================================================
 
216
// 2nd case : multi-line form initializer
 
217
//  A = { {a, b  },
 
218
//        {c, d} };
 
219
// =========================================================================
 
220
template <class T, class M>
 
221
struct form_concat {
 
222
// typedef:
 
223
 typedef typename form_basic<T,M>::size_type   size_type;
 
224
 typedef form_concat_line<T,M>           line_type;
 
225
 typedef form_concat_value<T,M>          value_type;
 
226
// allocators:
 
227
 form_concat () : _l() {}
 
228
 form_concat (const std::initializer_list<line_type>& il) : _l() {
 
229
#ifdef _RHEOLEF_HAVE_STD_INITIALIZER_ITERATOR
 
230
    typedef typename std::initializer_list<line_type>::const_iterator const_iterator;
 
231
#else // _RHEOLEF_HAVE_STD_INITIALIZER_ITERATOR
 
232
    typedef const line_type* const_iterator;
 
233
#endif // _RHEOLEF_HAVE_STD_INITIALIZER_ITERATOR
 
234
    for(const_iterator iter = il.begin(); iter != il.end(); ++iter) {
 
235
        _l.push_back(*iter);
 
236
    }
 
237
 }
 
238
 friend std::ostream& operator<< (std::ostream& o, const form_concat<T,M>& x) {
 
239
    std::cout << "{";
 
240
    for(typename std::list<line_type>::const_iterator iter = x._l.begin(); iter != x._l.end(); ++iter) {
 
241
        std::cout << *iter << " ";
 
242
    }
 
243
    return std::cout << "}";
 
244
 }
 
245
// internals:
 
246
 form_basic<T,M> build_form () const;
 
247
 
 
248
// data:
 
249
protected:
 
250
 std::list<line_type> _l;
 
251
};
 
252
template <class T, class M>
 
253
form_basic<T,M>
 
254
form_concat<T,M>::build_form() const
 
255
{
 
256
  // ---------------------------
 
257
  // pass 1 : compute spaces
 
258
  // ---------------------------
 
259
  size_t i_comp = 0;
 
260
  space_mult_list<T,M> sml_Y;
 
261
  std::vector<std::pair<bool,space_basic<T,M> > > l_Xh (_l.size(), std::pair<bool,space_basic<T,M> >(false, space_basic<T,M>()));
 
262
  for (typename std::list<line_type>::const_iterator iter = _l.begin(); iter != _l.end(); ++iter, i_comp++) {
 
263
    const line_type& line = *iter;
 
264
    space_basic<T,M> Yih;
 
265
    line.build_form_pass0 (l_Xh, Yih, i_comp);
 
266
    sml_Y *= Yih;
 
267
  }
 
268
  space_basic<T,M> Xh;
 
269
  form_concat_line<T,M>::build_first_space (l_Xh, Xh);
 
270
  space_basic<T,M> Yh (sml_Y);
 
271
  // ------------------------
 
272
  // pass 2 : copy
 
273
  // ------------------------
 
274
  csr_concat<T,M> uu, ub, bu, bb;
 
275
  for (typename std::list<line_type>::const_iterator iter = _l.begin(); iter != _l.end(); ++iter) {
 
276
    const line_type& line = *iter;
 
277
    csr_concat_line<T,M> uu_i, ub_i, bu_i, bb_i;
 
278
    for (typename std::list<value_type>::const_iterator jter = line.begin(); jter != line.end(); ++jter) {
 
279
      const value_type& x = *jter;
 
280
      switch (x.variant) {
 
281
        case form_concat_value<T,M>::scalar: {
 
282
          check_macro (x.s == 0, "unsupported non-nul scalar `"<<x.s<<"' in form concatenation");
 
283
          // zero: no values to insert in the sparse structure => nothing more to do
 
284
          uu_i.push_back (x.s);
 
285
          ub_i.push_back (x.s);
 
286
          bu_i.push_back (x.s);
 
287
          bb_i.push_back (x.s);
 
288
          break;
 
289
        }
 
290
        case form_concat_value<T,M>::form: {
 
291
          uu_i.push_back (x.m.uu());
 
292
          ub_i.push_back (x.m.ub());
 
293
          bu_i.push_back (x.m.bu());
 
294
          bb_i.push_back (x.m.bb());
 
295
          break;
 
296
        }
 
297
        default: error_macro ("non-form or scalar concatenation not yet supported");
 
298
      }
 
299
    }
 
300
    uu.push_back (uu_i);
 
301
    ub.push_back (ub_i);
 
302
    bu.push_back (bu_i);
 
303
    bb.push_back (bb_i);
 
304
  }
 
305
  form_basic<T,M> a(Xh, Yh);
 
306
  a.set_uu() = uu.build_csr();
 
307
  a.set_ub() = ub.build_csr();
 
308
  a.set_bu() = bu.build_csr();
 
309
  a.set_bb() = bb.build_csr();
 
310
  return a;
 
311
}
 
312
 
 
313
#ifdef _RHEOLEF_HAVE_STD_INITIALIZER_LIST
 
314
template <class T, class M>
 
315
inline
 
316
form_basic<T,M>::form_basic (const std::initializer_list<form_concat_line<T,M> >& init_list)
 
317
{
 
318
  form_concat<T,M> cc (init_list);
 
319
  form_basic<T,M>::operator= (cc.build_form());
 
320
}
 
321
#endif // _RHEOLEF_HAVE_STD_INITIALIZER_LIST
 
322
 
 
323
} // namespace rheolef
 
324
#endif // _RHEOLEF_FORM_CONCAT_H