rheolef  6.6
field_seq_put_vtk.cc
Go to the documentation of this file.
1 #include "rheolef/field.h"
2 #include "rheolef/field_expr_ops.h"
3 #include "rheolef/field_nonlinear_expr_ops.h"
4 #include "rheolef/interpolate.h"
5 #include "rheolef/piola.h"
6 #include "rheolef/rheostream.h"
7 #include "rheolef/iorheo.h"
8 #include "rheolef/field_evaluate.h"
9 #include "rheolef/space_component.h"
10 namespace rheolef {
11 using namespace std;
12 
13 template <class T> odiststream& geo_put_vtk (odiststream& ods, const geo_basic<T,sequential>& omega,
14  const numbering<T,sequential>& my_numb, const array<point_basic<T>,sequential>& my_node, bool append_data);
15 template <class T> odiststream& geo_put_vtk (odiststream& ods, const geo_basic<T,sequential>& omega,
16  size_t my_order, const array<point_basic<T>,sequential>& my_node, bool append_data);
17 
18 template <class T>
19 static
20 field_basic<T,sequential>
22 {
24  check_macro (Vh.get_geo() == uh.get_geo(), "interpolate: incompatible geos");
26  const space_basic<T,sequential>& Uh = uh.get_space();
27  const geo_basic<T,sequential>& omega = uh.get_geo();
28  std::vector<bool> dof_done (Vh.ndof(), false);
29  std::vector<size_type> u_dis_idof;
30  std::vector<size_type> v_dis_idof;
31  std::vector<point_basic<T> > hat_x;
32  basis_basic<T> fem_basis = Uh.get_numbering().get_basis();
33  basis_basic<T> pointset = Vh.get_numbering().get_basis();
34  basis_on_pointset<T> fem_basis_on_geo_nodes (pointset, fem_basis);
35  for (size_type ie = 0, ne = omega.size(); ie < ne; ie++) {
36  const geo_element& K = omega[ie];
37  Vh.get_numbering().get_basis().hat_node (K.variant(), hat_x);
38  Uh.dis_idof (K, u_dis_idof);
39  Vh.dis_idof (K, v_dis_idof);
40  for (size_type v_loc_idof = 0, v_loc_ndof = v_dis_idof.size(); v_loc_idof < v_loc_ndof; v_loc_idof++) {
41  if (dof_done[v_dis_idof[v_loc_idof]]) continue;
42  dof_done[v_dis_idof[v_loc_idof]] = true;
43  vh.dof (v_dis_idof[v_loc_idof]) = field_evaluate (uh, fem_basis_on_geo_nodes, K, u_dis_idof, v_loc_idof);
44  }
45  }
46  return vh;
47 }
48 template <class T>
49 static
50 field_basic<T,sequential>
52 {
53  if (uh.size() == 0) {
54  return interpolate_scalar (Vh, uh);
55  }
58  for (size_type i_comp = 0, n_comp = uh.size(); i_comp < n_comp; i_comp++) {
59  vh [i_comp] = interpolate_scalar (space_basic<T,sequential>(Vh[i_comp]), field_basic<T,sequential>(uh[i_comp]));
60  }
61  return vh;
62 }
63 template <class T>
64 odiststream&
65 put_vtk_scalar_values (odiststream& ods, const field_basic<T,sequential>& uh, std::string name, bool put_header)
66 {
68  ostream& vtk = ods.os();
69  size_type degree = uh.get_space().get_numbering().get_basis().degree();
70  vtk << setprecision(numeric_limits<T>::digits10);
71  if (put_header) {
72  std::string data_type = (degree == 0) ? "CELL_DATA" : "POINT_DATA";
73  vtk << data_type << " " << uh.ndof() << endl;
74  }
75  vtk << "SCALARS " << name << " float" << endl
76  << "LOOKUP_TABLE default" << endl;
77  for (size_type idof = 0, ndof = uh.ndof(); idof < ndof; idof++)
78  vtk << uh.dof(idof) << endl;
79  vtk << endl;
80  return ods;
81 }
82 template <class T>
83 odiststream&
84 put_vtk_vector_values (odiststream& ods, const field_basic<T,sequential>& uh, std::string name, bool put_header)
85 {
87  ostream& vtk = ods.os();
88  size_type n_comp = uh.size();
89  std::vector<field_component_const<T,sequential> > uh_comp (n_comp);
90  for (size_type i_comp = 0; i_comp < n_comp; i_comp++) {
91  uh_comp[i_comp] = uh[i_comp];
92  }
94  field_basic<T,sequential> norm_uh = interpolate (Xh, norm(uh));
95  put_vtk_scalar_values (ods, norm_uh, name+"_norm", put_header);
96  vtk << setprecision(numeric_limits<T>::digits10)
97  << "VECTORS " << name << " float" << endl;
98  std::vector<T> u_dof (n_comp);
99  for (size_type idof = 0, ndof = uh_comp[0].ndof(); idof < ndof; idof++) {
100  for (size_type i_comp = 0; i_comp < n_comp; i_comp++) {
101  vtk << uh_comp[i_comp].dof (idof);
102  if (i_comp != 2) vtk << " ";
103  }
104  for (size_type i_comp = n_comp; i_comp < 3; i_comp++) {
105  vtk << "0";
106  if (i_comp != 2) vtk << " ";
107  }
108  vtk << endl;
109  }
110  vtk << endl;
111  return ods;
112 }
113 template <class T>
114 odiststream&
115 put_vtk_tensor_values (odiststream& ods, const field_basic<T,sequential>& tau_h, std::string name, bool put_header)
116 {
118  ostream& vtk = ods.os();
119  space_basic<T,sequential> Xh (tau_h.get_geo(), tau_h.get_approx());
120  field_basic<T,sequential> norm_tau_h = interpolate (Xh, norm(tau_h));
121  put_vtk_scalar_values (ods, norm_tau_h, name+"_norm", put_header);
122  vtk << setprecision(numeric_limits<T>::digits10)
123  << "TENSORS " << name << " float" << endl;
124  size_type d = tau_h.get_geo().dimension();
125  switch (d) {
126  case 1: {
127  field_basic<T,sequential> t00 = tau_h(0,0);
128  for (size_type idof = 0, ndof = t00.ndof(); idof < ndof; idof++) {
129  vtk << t00.dof(idof) << " 0 0" << endl
130  << "0 0 0" << endl
131  << "0 0 0" << endl;
132  }
133  break;
134  }
135  case 2: {
136  field_basic<T,sequential> t00 = tau_h(0,0);
137  field_basic<T,sequential> t01 = tau_h(0,1);
138  field_basic<T,sequential> t11 = tau_h(1,1);
139  for (size_type idof = 0, ndof = t00.ndof(); idof < ndof; idof++) {
140  vtk << t00.dof(idof) << " " << t01.dof(idof) << " 0" << endl
141  << t01.dof(idof) << " " << t11.dof(idof) << " 0" << endl
142  << "0 0 0" << endl;
143  }
144  break;
145  }
146  default: {
147  field_basic<T,sequential> t00 = tau_h(0,0);
148  field_basic<T,sequential> t01 = tau_h(0,1);
149  field_basic<T,sequential> t11 = tau_h(1,1);
150  field_basic<T,sequential> t02 = tau_h(0,2);
151  field_basic<T,sequential> t12 = tau_h(1,2);
152  field_basic<T,sequential> t22 = tau_h(2,2);
153  for (size_type idof = 0, ndof = t00.ndof(); idof < ndof; idof++) {
154  vtk << t00.dof(idof) << " " << t01.dof(idof) << " " << t02.dof(idof) << endl
155  << t01.dof(idof) << " " << t11.dof(idof) << " " << t12.dof(idof) << endl
156  << t02.dof(idof) << " " << t12.dof(idof) << " " << t22.dof(idof) << endl;
157  }
158  }
159  }
160  return ods;
161 }
162 template <class T>
163 odiststream&
164 field_put_vtk (odiststream& ods, const field_basic<T,sequential>& uh, std::string name, bool put_geo)
165 {
167  size_type degree = uh.get_space().get_numbering().get_basis().degree();
168  size_type order = uh.get_space().get_geo().order();
170  if (degree == 0) {
171  vh = uh;
172  if (put_geo) {
173  geo_put_vtk (ods, vh.get_geo(), order, vh.get_geo().get_nodes(), false);
174  }
175  } else if (order <= degree && ! uh.get_space().get_numbering().has_compact_support_inside_element()) {
176  vh = uh;
177  if (put_geo) {
178  geo_put_vtk (ods, vh.get_geo(), degree, vh.get_space().get_xdofs(), false);
179  }
180  } else { // order > degree or discontinuous
181  std::string approx = "P" + itos(order);
182  if (uh.get_space().get_numbering().has_compact_support_inside_element()) approx += "d";
183  space_basic<T,sequential> Vh (uh.get_geo(), approx, uh.get_space().valued());
184  vh = interpolate (Vh, uh);
185  if (put_geo) {
186  geo_put_vtk (ods, vh.get_geo(), vh.get_space().get_numbering(), vh.get_space().get_xdofs(), false);
187  }
188  }
189  bool put_header = put_geo;
190  if (name == "") { name = vh.get_space().valued(); }
191  switch (vh.get_space().valued_tag()) {
192  case space_constant::scalar: put_vtk_scalar_values (ods, vh, name, put_header); break;
193  case space_constant::vector: put_vtk_vector_values (ods, vh, name, put_header); break;
194  case space_constant::tensor: put_vtk_tensor_values (ods, vh, name, put_header); break;
195  default: error_macro ("put_vtk: do not known how to print " << vh.valued() << "-valued field");
196  }
197  return ods;
198 }
199 template <class T>
200 odiststream&
202 {
203  return field_put_vtk (ods, uh, "", true);
204 }
205 template odiststream& field_put_vtk<Float> (odiststream&, const field_basic<Float,sequential>&, std::string, bool);
206 template odiststream& field_put_vtk<Float> (odiststream&, const field_basic<Float,sequential>&);
207 
208 }// namespace rheolef
const basis_basic< T > & get_basis() const
Definition: numbering.h:76
#define error_macro(message)
Definition: compiler.h:100
std::size_t size_type
Definition: field.h:228
field - piecewise polynomial finite element field
Definition: field.h:225
const geo_type & get_geo() const
Definition: field.h:272
const std::string & valued() const
Definition: field.h:276
size_type size(size_type dim) const
Definition: geo.h:1025
idiststream, odiststream - distributed interface for large data streams
Definition: diststream.h:68
size_type size() const
Definition: field.h:297
odiststream & put_vtk_tensor_values(odiststream &ods, const field_basic< T, sequential > &tau_h, std::string name, bool put_header)
std::string get_approx() const
Definition: field.h:274
STL namespace.
irheostream, orheostream - large data streams
Definition: compiler.h:7
T & dof(size_type idof)
Definition: field.h:527
static field_basic< T, sequential > interpolate(const space_basic< T, sequential > &Vh, const field_basic< T, sequential > &uh)
string itos(string::size_type i)
Definition: rheostream.cc:38
variant_type variant() const
odiststream & put_vtk_scalar_values(odiststream &ods, const field_basic< T, sequential > &uh, std::string name, bool put_header)
T norm(const vec< T, M > &x)
Definition: vec_expr.h:326
size_type ndof() const
Definition: field.h:306
#define check_macro(ok_condition, message)
Definition: compiler.h:104
size_t d
const numbering< T, sequential > & get_numbering() const
Definition: space.h:381
void put_header(odiststream &out, const branch_basic< T, sequential > &b)
Definition: branch.cc:30
odiststream & field_put_vtk(odiststream &, const field_basic< T, sequential > &, std::string, bool)
basis - polynomial basis
Definition: basis.h:77
odiststream & geo_put_vtk(odiststream &ods, const geo_basic< T, sequential > &omega, const numbering< T, sequential > &my_numb, const array< point_basic< T >, sequential > &my_node, bool append_data)
geo_element - element of a mesh
const space_type & get_space() const
Definition: field.h:271
odiststream & put_vtk_vector_values(odiststream &ods, const field_basic< T, sequential > &uh, std::string name, bool put_header)
static field_basic< T, sequential > interpolate_scalar(const space_basic< T, sequential > &Vh, const field_basic< T, sequential > &uh)
void dis_idof(const geo_element &K, std::vector< size_type > &dis_idof) const
Definition: space.h:451
reference_element_H::size_type size_type
template odiststream & field_put_vtk< Float >(odiststream &, const field_basic< Float, sequential > &, std::string, bool)
std::ostream & os()
Definition: diststream.h:166
T field_evaluate(const field_basic< T, M > &uh, const basis_on_pointset< T > &bops, reference_element hat_K, const std::vector< size_t > &dis_idof, size_t q)
verbose clean transpose logscale grid shrink ball stereo iso volume velocity elevation bezieradapt label color color format format format format format format format format format format format format format format format format vtk
const geo_basic< T, sequential > & get_geo() const
Definition: space.h:374