Rheolef  7.1
an efficient C++ finite element environment
basis_visu_gnuplot.icc
Go to the documentation of this file.
1 #ifndef _RHEOLEF_BASIS_VISU_GNUPLOT_ICC
2 #define _RHEOLEF_BASIS_VISU_GNUPLOT_ICC
3 // file automatically generated by: ../../../rheolef/nfem/pbasis/make_basis_list_cxx.sh
24 #include "rheolef/basis_raw.h"
25 #include "rheolef/iorheo.h"
26 #include "rheolef/rheostream.h"
27 namespace rheolef { namespace details {
28 using namespace std;
29 
30 template<class Basis>
31 void
32 put (const Basis& b, ostream& os, reference_element hat_K)
33 {
34  typedef typename Basis::size_type size_type;
35  typedef typename Basis::value_type T;
36  bool verbose = iorheo::getverbose(os);
37  bool clean = iorheo::getclean(os);
38  bool execute = iorheo::getexecute(os);
39  size_type nsub = iorheo::getsubdivide(os);
40  if (nsub <= 1) nsub = 20; // default value
41 
42  size_type loc_ndof = reference_element::n_node(hat_K.variant(), b.degree());
43  string basename = "basis-" + b.name() + "-" + hat_K.name();
44  string filelist;
45  // --------------------
46  // .plot
47  // --------------------
48  string plot_name = basename + ".plot";
49  string gdat_name = basename + ".gdat";
50  ofstream plot (plot_name.c_str());
51  cerr << "! file \"" << plot_name << "\" created" << endl;
52  filelist += " " + plot_name;
53  size_t d = hat_K.dimension();
54  check_macro (d==1 || d==2, "unsupported dimension " << d);
55  plot << "gdat = \"" << gdat_name << "\"" << endl
56  << "set colors classic" << endl
57  << "set size square" << endl
58  ;
59  if (d == 1) {
60  plot << "plot \\" << endl
61  ;
62  for (size_t loc_idof = 0; loc_idof < loc_ndof; ++loc_idof) {
63  plot << " gdat u 1:"<<loc_idof+2 << " t \"L" << loc_idof+1 << "\" w l";
64  if (loc_idof+1 != loc_ndof) {
65  plot << ",\\";
66  }
67  plot << endl;
68  }
69  } else if (d == 2) {
70  plot << "pause_duration = 0.5" << endl
71  << "d = 2" << endl
72  << "i = d+1" << endl
73  << "imax = d+" << loc_ndof << endl
74  << "set hidden3d back # offset 1 trianglepattern 3 undefined 1 altdiagonal bentover" << endl
75  << "set linetype 1 lw 1 lc rgb \"#000000\"" << endl
76  << "set linetype 2 lw 1 lc rgb \"#0000ff\"" << endl
77  << "splot gdat u 1:2:i t sprintf(\"P%d\",i-d) w l lc 0" << endl
78  << "load \"" << basename << ".loop\"" << endl
79  ;
80  string loop_name = basename + ".loop";
81  ofstream loop (loop_name.c_str());
82  cerr << "! file \"" << loop_name << "\" created" << endl;
83  filelist += " " + loop_name;
84  loop << "i = i+1" << endl
85  << "if (i <= imax) \\" << endl
86  << " pause pause_duration; \\" << endl
87  << " replot; \\" << endl
88  << " reread" << endl
89  ;
90  }
91  plot << "pause -1 \"<return>\"" << endl;
92  plot.close();
93 
94  // --------------------
95  // .gdat
96  // --------------------
97  ofstream gdat (gdat_name.c_str());
98  cerr << "! file \"" << gdat_name << "\" created" << endl;
99  filelist += " " + gdat_name;
100  gdat << setprecision(std::numeric_limits<T>::digits10)
101  << "# basis " << b.name() << endl
102  << "# element " << hat_K.name() << endl
103  << "# degree " << b.degree() << endl
104  ;
105  Eigen::Matrix<T,Eigen::Dynamic,1> value;
106  switch (hat_K.variant()) {
107  case reference_element::p: {
108  break;
109  }
110  case reference_element::e: {
111  gdat << "# size " << loc_ndof << endl
112  << "# x";
113  for (size_t loc_idof = 0; loc_idof < loc_ndof; ++loc_idof) {
114  gdat << " L" << loc_idof+1;
115  }
116  gdat << endl;
117  for (size_type i = 0; i <= nsub; i++) {
118  point_basic<T> hat_x (T(int(i))/T(int(nsub)));
119  b.evaluate (hat_K, hat_x, value);
120  gdat << hat_x[0];
121  for (size_t loc_idof = 0, loc_ndof = value.size(); loc_idof < loc_ndof; ++loc_idof) {
122  gdat << " " << value[loc_idof];
123  }
124  gdat << endl;
125  }
126  break;
127  }
128  case reference_element::t: {
129  for (size_type j = 0; j <= nsub; j++) {
130  for (size_type i1 = 0; i1 <= nsub; i1++) {
131  size_type i = std::min(i1, nsub-j);
132  point_basic<T> hat_x (T(int(i))/T(int(nsub)), T(int(j))/T(int(nsub)));
133  b.evaluate (hat_K, hat_x, value);
134  gdat << hat_x[0] << " " << hat_x[1];
135  for (size_t loc_idof = 0, loc_ndof = value.size(); loc_idof < loc_ndof; ++loc_idof) {
136  gdat << " " << value[loc_idof];
137  }
138  gdat << endl;
139  }
140  gdat << endl;
141  }
142  break;
143  }
144  case reference_element::q: {
145  for (size_type j = 0; j <= nsub; j++) {
146  for (size_type i = 0; i <= nsub; i++) {
147  point_basic<T> hat_x (2*T(int(i))/T(int(nsub))-1, 2*T(int(j))/T(int(nsub))-1);
148  b.evaluate (hat_K, hat_x, value);
149  gdat << hat_x[0] << " " << hat_x[1];
150  for (size_t loc_idof = 0, loc_ndof = value.size(); loc_idof < loc_ndof; ++loc_idof) {
151  gdat << " " << value[loc_idof];
152  }
153  gdat << endl;
154  }
155  gdat << endl;
156  }
157  break;
158  }
159 #ifdef TODO // 3D visu ?
160  case reference_element::T: {
161  for (size_type k = 0; k <= nsub; k++) {
162  for (size_type j = 0; j+k <= nsub; j++) {
163  for (size_type i = 0; i+j+k <= nsub; i++) {
164  point_basic<T> hat_x (T(i)/T(nsub), T(j)/T(nsub), T(k)/T(nsub));
165  b.eval_lagrange (hat_x, hat_K, degree, sopt, inv_vdm, value);
166  Lambda = max(Lambda, norm(value,1));
167  }
168  }
169  }
170  return Lambda;
171  }
172  case reference_element::P: {
173  for (size_type k = 0; k <= degree; k++) {
174  for (size_type j = 0; j <= degree; j++) {
175  for (size_type i = 0; i+j <= degree; i++) {
176  point_basic<T> hat_x (T(i)/T(nsub), T(j)/T(nsub), 2*T(k)/T(nsub)-1);
177  b.eval_lagrange (hat_x, hat_K, degree, sopt, inv_vdm, value);
178  Lambda = max(Lambda, norm(value,1));
179  }
180  }
181  }
182  break;
183  }
184  case reference_element::H: {
185  for (size_type k = 0; k <= degree; k++) {
186  for (size_type j = 0; j <= degree; j++) {
187  for (size_type i = 0; i <= degree; i++) {
188  point_basic<T> hat_x (2*T(i)/T(nsub)-1, 2*T(j)/T(nsub), 2*T(k)/T(nsub)-1);
189  b.eval_lagrange (hat_x, hat_K, degree, sopt, inv_vdm, value);
190  Lambda = max(Lambda, norm(value,1));
191  }
192  }
193  }
194  break;
195  }
196 #endif // TODO
197  default:
198  error_macro ("unexpected element type `"<<hat_K.name()<<"'");
199  }
200  // -----------
201  if (execute) {
202  // -----------
203  string command = "gnuplot " + plot_name;
204  if (verbose) clog << "! " << command << endl;
205  int status = system (command.c_str());
206  }
207  // -----------
208  if (clean) {
209  // -----------
210  string command = "/bin/rm -f " + filelist;
211  if (verbose) clog << "! " << command << endl;
212  int status = system (command.c_str());
213  }
214 }
215 
216 }} // namespace rheolef::details
217 #endif // _RHEOLEF_BASIS_VISU_GNUPLOT_ICC
rheolef::reference_element::e
static const variant_type e
Definition: reference_element.h:76
rheolef::reference_element::H
static const variant_type H
Definition: reference_element.h:81
rheolef::reference_element::n_node
static size_type n_node(variant_type variant, size_type order)
Definition: reference_element.cc:201
mkgeo_couette.basename
string basename
Definition: mkgeo_couette.sh:73
rheolef::point_basic
Definition: point.h:87
check_macro
check_macro(expr1.have_homogeneous_space(Xh1), "dual(expr1,expr2); expr1 should have homogeneous space. HINT: use dual(interpolate(Xh, expr1),expr2)")
mkgeo_ball.b
int b
Definition: mkgeo_ball.sh:152
rheolef::value
rheolef::std value
rheolef::reference_element::T
static const variant_type T
Definition: reference_element.h:79
mkgeo_ball.d
int d
Definition: mkgeo_ball.sh:154
rheolef::size_type
size_t size_type
Definition: basis_get.cc:76
rheolef::norm
T norm(const vec< T, M > &x)
norm(x): see the expression page for the full documentation
Definition: vec.h:387
rheolef::details::put
void put(const Basis &b, ostream &os, reference_element hat_K)
Definition: basis_visu_gnuplot.icc:32
rheolef::reference_element
see the reference_element page for the full documentation
Definition: reference_element.h:66
mkgeo_ball.command
string command
Definition: mkgeo_ball.sh:136
rheolef::reference_element::name
char name() const
Definition: reference_element.h:100
mkgeo_ball.clean
clean
Definition: mkgeo_ball.sh:335
rheolef::reference_element::variant
variant_type variant() const
Definition: reference_element.h:99
rheolef
This file is part of Rheolef.
Definition: compiler_eigen.h:37
error_macro
#define error_macro(message)
Definition: dis_macros.h:49
rheolef::reference_element::dimension
size_type dimension() const
Definition: reference_element.h:101
mkgeo_ball.verbose
verbose
Definition: mkgeo_ball.sh:133
rheolef::reference_element::p
static const variant_type p
Definition: reference_element.h:75
rheolef::reference_element::q
static const variant_type q
Definition: reference_element.h:78
rheolef::reference_element::P
static const variant_type P
Definition: reference_element.h:80
rheolef::reference_element::t
static const variant_type t
Definition: reference_element.h:77
size_type
field::size_type size_type
Definition: branch.cc:425
mkgeo_contraction.status
status
Definition: mkgeo_contraction.sh:290
T
Expr1::float_type T
Definition: field_expr.h:261