← Back to branch summary

~ubuntu-branches/ubuntu/saucy/rheolef/saucy

« back to all changes in this revision

Viewing changes to nfem/plib/field.cc

  • Committer: Bazaar Package Importer
  • Author(s): Pierre Saramito
  • Date: 2011-03-23 11:14:26 UTC
  • mfrom: (1.1.4 upstream)
  • Revision ID: james.westby@ubuntu.com-20110323111426-cjvhey7lxt6077ty
Tags: 5.93-1
* New upstream release (minor changes):
  - some extra warning message deleted in heap_allocator
  - graphic output with mayavi2 fixed
  - add doc refman in .info and .pdf format

Show diffs side-by-side

added added

removed removed

Lines of Context:
nfem/plib/field.cc
19
19
/// 
20
20
/// =========================================================================
21
21
#include "rheolef/field.h"
 
22
#include <complex>
22
23
 
23
24
namespace rheolef {
24
25
 
25
 
template <class T>
26
 
field_basic<T>::field_basic (
27
 
        const space& V,
28
 
        const Float& init_value)
 
26
template <class T, class M>
 
27
field_basic<T,M>::field_basic (
 
28
        const space_type& V,
 
29
        const T& init_value)
29
30
 : _V (V), 
30
 
   _u (),
31
 
   _b ()
32
 
{
33
 
   _u.resize (_V.iu_ownership(), init_value);
34
 
   _b.resize (_V.ib_ownership(), init_value);
35
 
}
36
 
template <class T>
37
 
oparstream& 
38
 
field_basic<T>::put (oparstream& ops) const
 
31
   u (),
 
32
   b ()
 
33
{
 
34
   u.resize (_V.iu_ownership(), init_value);
 
35
   b.resize (_V.ib_ownership(), init_value);
 
36
}
 
37
template <class T, class M>
 
38
void 
 
39
field_basic<T,M>::resize (
 
40
        const space_type& V,
 
41
        const T& init_value)
 
42
{
 
43
   if (_V == V) return;
 
44
   _V = V;
 
45
   u.resize (_V.iu_ownership(), init_value);
 
46
   b.resize (_V.ib_ownership(), init_value);
 
47
}
 
48
template <class T, class M>
 
49
idiststream& 
 
50
field_basic<T,M>::get (idiststream& ips)
 
51
{
 
52
warning_macro ("field::get...");
 
53
  communicator comm = ips.comm();
 
54
  if ( ! dis_scatch (ips, comm, "\nfield")) {
 
55
    error_macro ("read field failed");
 
56
    return ips;
 
57
  }
 
58
  size_type      version,  dis_size1;
 
59
  std::string    geo_name, approx;
 
60
  ips >> version >> dis_size1
 
61
      >> geo_name
 
62
      >> approx;
 
63
  geo_type omega (geo_name);
 
64
 
 
65
  if (!(_V.get_geo().name() == geo_name && _V.get_element().name() == approx)) {
 
66
    // cannot re-use _V: build a new space and resize the field
 
67
warning_macro ("field::get: resize...");
 
68
    resize (space_type (omega,approx));
 
69
warning_macro ("field::get: resize done");
 
70
  }
 
71
  // TODO: read and permut in one pass ? avoid 2 passes of global exchanges
 
72
  vec<T,M> u_io  (_V.ios_ownership(), std::numeric_limits<T>::max());
 
73
  u_io.get_values (ips);
 
74
  vec<T,M> u_new (_V.ownership(),     std::numeric_limits<T>::max());
 
75
  for (size_type ios_idof = 0, ios_ndof = _V.ios_ownership().size(); ios_idof < ios_ndof; ios_idof++) {
 
76
    T value = u_io [ios_idof];
 
77
    size_type dis_idof = _V.ios_idof2dis_idof (ios_idof);
 
78
    u_new.dis_entry (dis_idof) = value;
 
79
  }
 
80
  u_new.dis_entry_assembly();
 
81
 
 
82
  for (size_type idof = 0, ndof = _V.ownership().size(); idof < ndof; idof++) {
 
83
    dof (idof) = u_new [idof];
 
84
  }
 
85
warning_macro ("field::get done");
 
86
  return ips;
 
87
}
 
88
template <class T, class M>
 
89
odiststream& 
 
90
field_basic<T,M>::put (odiststream& ops) const
39
91
{
40
92
    using namespace std;
41
 
    // merge distributed blocked and non-blocked onto proc=0
42
 
    communicator_type comm = _V.comm();
43
 
    size_type io_proc = oparstream::io_proc();
44
 
    size_type my_proc = comm.rank();
45
 
    distributor io_ownership (_V.par_size(), comm, (my_proc == io_proc ? _V.par_size() : 0)); 
46
 
    vec<T> u_io (io_ownership, std::numeric_limits<T>::max());
47
 
    for (size_type idof = 0, n_dof = _V.size(); idof < n_dof; idof++) {
48
 
        T value = operator[] (idof);
49
 
        size_type par_old_idof = _V.idof2par_old_idof (idof);
50
 
        u_io.set (par_old_idof, value);
51
 
        trace_macro ("idof="<<idof<<", par_old_idof="<<par_old_idof<<", value="<<value);
 
93
    // merge distributed blocked and non-blocked and apply iso_dof permutation:
 
94
    vec<T,M> u_io (_V.ios_ownership(), std::numeric_limits<T>::max());
 
95
    warning_macro ("V.ios_ownership.size=" << _V.ios_ownership().size());
 
96
    warning_macro ("V.ownership.size=" << _V.ownership().size());
 
97
    for (size_type idof = 0, ndof = _V.size(); idof < ndof; idof++) {
 
98
        T value = dof (idof);
 
99
        size_type ios_dis_idof = _V.idof2ios_dis_idof (idof);
 
100
        u_io.dis_entry (ios_dis_idof) = value;
52
101
    }
53
 
    u_io.assembly();  
 
102
    u_io.dis_entry_assembly();  
54
103
 
55
104
    // then output
56
 
    ops << "field" << endl
57
 
        << "1 " << u_io.par_size() << endl 
 
105
    ops << setprecision(std::numeric_limits<Float>::digits10)
 
106
        << "field" << endl
 
107
        << "1 " << u_io.dis_size() << endl 
58
108
        << _V.get_geo().name() << endl
59
 
        << "P1" << endl
 
109
        << _V.get_element().name() << endl
60
110
        << endl
61
111
        << u_io;
 
112
    
62
113
    return ops;
63
114
}
64
115
// ----------------------------------------------------------------------------
65
116
// instanciation in library
66
117
// ----------------------------------------------------------------------------
67
 
template class field_basic<Float>;
68
 
template oparstream& operator << (oparstream& s, const field_basic<Float>& x);
 
118
#define _RHEOLEF_instanciation_base(T,M)                                \
 
119
template class field_basic<T,M>;                                        \
 
120
template odiststream& operator << (odiststream&, const field_basic<T,M>&);
 
121
 
 
122
#ifdef TODO
 
123
#define _RHEOLEF_instanciation(T,M)                                     \
 
124
_RHEOLEF_instanciation_base(T,M)                                        \
 
125
_RHEOLEF_instanciation_base(std::complex<T>,M)
 
126
#else // TODO
 
127
#define _RHEOLEF_instanciation(T,M)                                     \
 
128
_RHEOLEF_instanciation_base(T,M)
 
129
#endif // TODO
 
130
 
 
131
_RHEOLEF_instanciation(Float,sequential)
 
132
#ifdef _RHEOLEF_HAVE_MPI
 
133
_RHEOLEF_instanciation(Float,distributed)
 
134
#endif // _RHEOLEF_HAVE_MPI
69
135
 
70
136
} // namespace rheolef