rheolef  6.5
field_seq_put_gmsh_pos.cc
Go to the documentation of this file.
1 #include "rheolef/field.h"
2 #include "rheolef/field_expr_ops.h"
3 #include "rheolef/piola.h"
4 #include "rheolef/rheostream.h"
5 #include "rheolef/iorheo.h"
6 #include "rheolef/field_evaluate.h"
7 #include "rheolef/space_component.h"
8 namespace rheolef {
9 using namespace std;
10 
11 template <class T>
12 odiststream&
13 field_put_gmsh_pos (odiststream& ods, const field_basic<T,sequential>& uh, std::string name)
14 {
16  if (name == "") { name = uh.get_space().valued(); }
17  ostream& gmsh = ods.os();
18  gmsh << setprecision(numeric_limits<T>::digits10);
19  check_macro (uh.get_space().get_numbering().name() == "P1",
20  "gmsh: unsupported " << uh.get_space().get_numbering().name() << " approximation");
21  check_macro (uh.get_geo().order() == 1,
22  "gmsh: unsupported geo order > 1");
23  /*
24  type #list-of-coords #list-of-values
25  --------------------------------------------------------------------
26  Scalar point SP 3 1 * nb-time-steps
27  Vector point VP 3 3 * nb-time-steps
28  Tensor point TP 3 9 * nb-time-steps
29  Scalar line SL 6 2 * nb-time-steps
30  Vector line VL 6 6 * nb-time-steps
31  Tensor line TL 6 18 * nb-time-steps
32  Scalar triangle ST 9 3 * nb-time-steps
33  Vector triangle VT 9 9 * nb-time-steps
34  Tensor triangle TT 9 27 * nb-time-steps
35  Scalar quadrangle SQ 12 4 * nb-time-steps
36  Vector quadrangle VQ 12 12 * nb-time-steps
37  Tensor quadrangle TQ 12 36 * nb-time-steps
38  Scalar tetrahedron SS 12 4 * nb-time-steps
39  Vector tetrahedron VS 12 12 * nb-time-steps
40  Tensor tetrahedron TS 12 36 * nb-time-steps
41  Scalar hexahedron SH 24 8 * nb-time-steps
42  Vector hexahedron VH 24 24 * nb-time-steps
43  Tensor hexahedron TH 24 72 * nb-time-steps
44  Scalar prism SI 18 6 * nb-time-steps
45  Vector prism VI 18 18 * nb-time-steps
46  Tensor prism TI 18 54 * nb-time-steps
47  Scalar pyramid SY 15 5 * nb-time-steps
48  Vector pyramid VY 15 15 * nb-time-steps
49  Tensor pyramid TY 15 45 * nb-time-steps
50  */
51  gmsh << "View \"" << name << "\" {" << endl;
52  static char gmsh_pos_type [reference_element::max_variant];
53  gmsh_pos_type [reference_element::p] = 'P';
54  gmsh_pos_type [reference_element::e] = 'L';
55  gmsh_pos_type [reference_element::t] = 'T';
56  gmsh_pos_type [reference_element::q] = 'Q';
57  gmsh_pos_type [reference_element::T] = 'S';
58  gmsh_pos_type [reference_element::H] = 'H';
59  gmsh_pos_type [reference_element::P] = 'I';
60  const geo_basic<T,sequential>& omega = uh.get_geo();
61  switch (uh.get_space().valued_tag()) {
63  for (size_type ie = 0, ne = omega.size(); ie < ne; ie++) {
64  const geo_element& K = omega[ie];
65  gmsh << "S" << gmsh_pos_type[K.variant()] << "(";
66  for (size_type iloc = 0, nloc = K.size(); iloc < nloc; iloc++) {
67  const point_basic<T>& x = omega.node(K[iloc]);
68  for (size_type i_comp = 0; i_comp < 3; i_comp++) {
69  gmsh << x[i_comp];
70  if (i_comp != 2) gmsh << ",";
71  }
72  if (iloc != nloc-1) gmsh << ",";
73  }
74  gmsh << "){";
75  for (size_type iloc = 0, nloc = K.size(); iloc < nloc; iloc++) {
76  gmsh << uh.dof(K[iloc]);
77  if (iloc != nloc-1) gmsh << ",";
78  }
79  gmsh << "};" << endl;
80  }
81  break;
82  }
83 #ifdef TODO
85  break;
86  }
87 #endif // TODO
90  size_type d = uh.get_geo().dimension();
91  for (size_type i_comp = 0; i_comp < d; i_comp++) {
92  for (size_type j_comp = 0; j_comp < d; j_comp++) {
93  uh_comp[i_comp][j_comp] = uh(i_comp,j_comp);
94  }}
95  tensor_basic<T> value;
96 #define HAVE_ID_DEFAULT_TENSOR
97 #ifdef HAVE_ID_DEFAULT_TENSOR
98  value (0,0) = value (1,1) = value (2,2) = 1;
99 #endif
100  for (size_type ie = 0, ne = omega.size(); ie < ne; ie++) {
101  const geo_element& K = omega[ie];
102  gmsh << "T" << gmsh_pos_type[K.variant()] << "(";
103  for (size_type iloc = 0, nloc = K.size(); iloc < nloc; iloc++) {
104  const point_basic<T>& x = omega.node(K[iloc]);
105  for (size_type i_comp = 0; i_comp < 3; i_comp++) {
106  gmsh << x[i_comp];
107  if (i_comp != 2) gmsh << ",";
108  }
109  if (iloc != nloc-1) gmsh << ",";
110  }
111  gmsh << "){";
112  for (size_type iloc = 0, nloc = K.size(); iloc < nloc; iloc++) {
113  size_type inod = K[iloc];
114  const point_basic<T>& x = omega.node(inod);
115  for (size_type i_comp = 0; i_comp < d; i_comp++) {
116  for (size_type j_comp = 0; j_comp < d; j_comp++) {
117  value(i_comp,j_comp) = uh_comp [i_comp][j_comp].dof (inod);
118  }}
119  for (size_type i_comp = 0; i_comp < 3; i_comp++) {
120  for (size_type j_comp = 0; j_comp < 3; j_comp++) {
121  gmsh << value(i_comp,j_comp);
122  if (i_comp != 2 || j_comp != 2 || iloc != nloc-1) gmsh << ",";
123  }}
124  }
125  gmsh << "};" << endl;
126  }
127  break;
128  }
129  default: error_macro ("put_gmsh: do not known how to print " << uh.valued() << "-valued field");
130  }
131  gmsh << "};" << endl;
132  return ods;
133 }
134 template <class T>
135 odiststream&
137 {
138  return field_put_gmsh_pos (ods, uh, "");
139 }
140 template odiststream& field_put_gmsh_pos<Float> (odiststream&, const field_basic<Float,sequential>&, std::string);
141 template odiststream& field_put_gmsh_pos<Float> (odiststream&, const field_basic<Float,sequential>&);
142 
143 }// namespace rheolef
144