28# include "rheolef/field.h"
29# include "rheolef/field_wdof_indirect.h"
30# include "rheolef/field_evaluate.h"
31# include "rheolef/piola_util.h"
32# include "rheolef/rheostream.h"
34# include "rheolef/rounder.h"
35# include "rheolef/interpolate.h"
40odiststream&
visu_gnuplot (odiststream& ops,
const geo_basic<T,sequential>& omega);
50 const geo_basic<T,sequential>& omega,
52 const field_basic<T,sequential>& uh,
53 const fem_on_pointset<T>& fops,
58 typedef typename geo_basic<T,sequential>::size_type
size_type;
59 typedef point_basic<size_type> ilat;
61 const Eigen::Matrix<piola<T>,Eigen::Dynamic,1>& piola = fops.get_piola_on_pointset().get_piola (omega, K);
62 Eigen::Matrix<T,Eigen::Dynamic,1>
value;
64 for (
size_type loc_idof = 0, loc_ndof = piola.size(); loc_idof < loc_ndof; loc_idof++) {
65 bbox.update (piola[loc_idof].F, value[loc_idof]);
67 for (
size_type i = 0; i <= my_order; i++) {
69 piola[loc_idof].F.put (gdat, dim); gdat <<
" " << value[loc_idof] << endl;
79 const geo_basic<T,sequential>& omega,
81 const field_basic<T,sequential>& uh,
82 const fem_on_pointset<T>& fops,
87 typedef typename geo_basic<T,sequential>::size_type
size_type;
88 typedef point_basic<size_type> ilat;
90 const Eigen::Matrix<piola<T>,Eigen::Dynamic,1>& piola = fops.get_piola_on_pointset().get_piola (omega, K);
91 Eigen::Matrix<T,Eigen::Dynamic,1> value;
93 for (
size_type loc_idof = 0, loc_ndof = piola.size(); loc_idof < loc_ndof; loc_idof++) {
94 bbox.update (piola[loc_idof].F, value[loc_idof]);
96 for (
size_type j = 0; j <= my_order; j++) {
97 for (
size_type i1 = 0; i1 <= my_order; i1++) {
100 piola[iloc].F.put (gdat, dim); gdat <<
" " << value[iloc] << endl;
104 gdat << endl << endl;
111 const geo_basic<T,sequential>& omega,
112 const geo_element& K,
113 const field_basic<T,sequential>& uh,
114 const fem_on_pointset<T>& fops,
119 typedef typename geo_basic<T,sequential>::size_type
size_type;
120 typedef point_basic<size_type> ilat;
122 const Eigen::Matrix<piola<T>,Eigen::Dynamic,1>& piola = fops.get_piola_on_pointset().get_piola (omega, K);
123 Eigen::Matrix<T,Eigen::Dynamic,1> value;
125 for (
size_type loc_idof = 0, loc_ndof = piola.size(); loc_idof < loc_ndof; loc_idof++) {
126 bbox.update (piola[loc_idof].F, value[loc_idof]);
128 for (
size_type j = 0; j < my_order+1; j++) {
129 for (
size_type i = 0; i < my_order+1; i++) {
131 piola[loc_idof00].F.put (gdat, dim); gdat <<
" " << value[loc_idof00] << endl;
135 gdat << endl << endl;
166 ostream& os = ods.
os();
167 bool verbose = iorheo::getverbose(os);
168 bool clean = iorheo::getclean(os);
169 bool execute = iorheo::getexecute(os);
170 bool fill = iorheo::getfill(os);
171 bool elevation = iorheo::getelevation(os);
172 bool color = iorheo::getcolor(os);
173 bool gray = iorheo::getgray(os);
174 bool black_and_white = iorheo::getblack_and_white(os);
175 bool reader_on_stdin = iorheo::getreader_on_stdin(os);
176 string format = iorheo::getimage_format(os);
177 string basename = iorheo::getbasename(os);
178 size_type subdivide = iorheo::getsubdivide(os);
179 size_type n_isovalue = iorheo::getn_isovalue(os);
180 size_type n_isovalue_negative = iorheo::getn_isovalue_negative(os);
181 string outfile_fmt =
"";
183 if (!clean) tmp =
"";
187 size_type map_dim = omega.map_dimension();
188 size_type nv = omega.sizes().ownership_by_dimension[0].size();
189 size_type nedg = omega.sizes().ownership_by_dimension[1].size();
190 size_type nfac = omega.sizes().ownership_by_dimension[2].size();
191 size_type nvol = omega.sizes().ownership_by_dimension[3].size();
192 size_type ne = omega.sizes().ownership_by_dimension[map_dim].size();
195 if (subdivide == 0) {
196 subdivide = std::max(omega.order(), subdivide);
197 subdivide = std::max(b.degree (), subdivide);
199 basis_basic<T> subdivide_pointset (
"P"+std::to_string(subdivide));
204 bbox.
xmin = omega.xmin();
205 bbox.
xmax = omega.xmax();
209 std::vector<T> values;
210 if (n_isovalue_negative != 0) {
211 for (
size_t i = 0; i <= n_isovalue_negative; i++) {
212 values.push_back (bbox.
umin*(n_isovalue_negative - i)/n_isovalue_negative);
214 for (
size_t i = 1; i <= n_isovalue - n_isovalue_negative; i++) {
215 values.push_back (bbox.
umax*i/(n_isovalue - n_isovalue_negative));
222 string filename = tmp+basename +
".gdat";
223 string gdatname = filename;
224 filelist = filelist +
" " + filename;
225 ofstream gdat (filename.c_str());
226 if (verbose) clog <<
"! file \"" << filename <<
"\" created.\n";
227 gdat << setprecision(numeric_limits<float_type>::digits10);
228 size_type used_dim = (fill ? map_dim : 1);
229 for (
size_type ie = 0, ne = omega.size(used_dim); ie < ne; ie++) {
230 const geo_element& K = omega.get_geo_element(used_dim,ie);
231 put (gdat, omega, K, uh, fops, subdivide, bbox);
240 if (fabs(bbox.
umax - bbox.
umin) < eps) {
247 filename = tmp+basename +
".plot";
248 filelist = filelist +
" " + filename;
249 ofstream plot (filename.c_str());
250 if (verbose) clog <<
"! file \"" << filename <<
"\" created.\n";
252 plot <<
"#!gnuplot" << endl
253 << setprecision(numeric_limits<float_type>::digits10);
255 outfile_fmt = basename +
"." + format;
256 string terminal = format;
257 if (terminal ==
"ps") {
258 terminal =
"postscript eps";
259 if (color) terminal +=
" color";
261 if (terminal ==
"jpg") terminal =
"jpeg";
262 if (terminal ==
"jpeg" || terminal ==
"png" || terminal ==
"gif") {
265 plot <<
"set terminal " << terminal << endl
266 <<
"set output \"" << outfile_fmt <<
"\"" << endl;
268 plot <<
"umin = " << bbox.
umin << endl
269 <<
"umax = " << bbox.
umax << endl
270 <<
"n_isovalue = " << n_isovalue << endl
271 <<
"uincr = 1.0*(umax-umin)/n_isovalue" << endl;
273 if (bbox.
xmin[0] >= bbox.
xmax[0]) plot <<
"#";
274 plot <<
"set xrange [" << bbox.
xmin[0] <<
":" << bbox.
xmax[0] <<
"]" << endl;
275 if (bbox.
xmin[1] >= bbox.
xmax[1]) plot <<
"#";
276 plot <<
"set yrange [" << bbox.
xmin[1] <<
":" << bbox.
xmax[1] <<
"]" << endl;
277 if (bbox.
umin >= bbox.
umax) plot <<
"#";
278 plot <<
"set zrange [umin:umax]" << endl;
280 plot <<
"set size square" << endl;
282 plot <<
"set size ratio -1 # equal scales" << endl
283 <<
"set noxtics" << endl
284 <<
"set noytics" << endl;
287 plot <<
"set xyplane at umin-0.1*(umax-umin)" << endl;
289 plot <<
"set view map" << endl;
292 plot <<
"unset surface" << endl
293 <<
"set contour base" << endl;
294 if (values.size() != 0) {
295 plot <<
"set cntrparam levels discrete ";
296 for (
size_t i = 0, n = values.size(); i < n; i++) {
298 if (i+1 != n) plot <<
", ";
302 plot <<
"eps = (umax-umin)*1e-4"<< endl;
303 plot <<
"uincr_eps = 1.0*(umax-umin-2*eps)/n_isovalue"<< endl;
304 plot <<
"set cntrparam levels incremental umin+eps,uincr_eps,umax-eps"<< endl;
307 if (black_and_white) {
308 plot <<
"set cbtics out scale 0.5" << endl;
310 plot <<
"set cbtics uincr" << endl
311 <<
"set cbrange [umin:umax]" << endl;
313 plot <<
"set palette gray" << endl;
315 plot <<
"set palette rgbformulae 33,13,-4" << endl;
317 plot <<
"set palette rgbformulae 0,0,0" << endl;
319 plot <<
"set palette maxcolors n_isovalue+1" << endl
320 <<
"set nokey" << endl;
321 }
else if (dim == 3 && map_dim == 2) {
324 T dx_max = max(dx[0],max(dx[1],dx[2]));
325 if (dx_max == 0) dx_max = 0.1;
326 dx[0] = max(dx[0],dx_max);
327 if (omega.dimension() >= 2) dx[1] = max(dx[1],dx_max);
328 if (omega.dimension() == 3) dx[2] = max(dx[2],dx_max);
332 plot <<
"set xrange [" << xmin[0] <<
":" << xmax[0] <<
"]" << endl
333 <<
"set yrange [" << xmin[1] <<
":" << xmax[1] <<
"]" << endl
334 <<
"set zrange [" << xmin[2] <<
":" << xmax[2] <<
"]" << endl
335 <<
"set xyplane at " << xmin[2] << endl
336 <<
"set view equal xyz # equal scales" << endl
337 <<
"set view 70,120" << endl
338 <<
"n_isovalue = " << n_isovalue << endl
339 <<
"n_subdivide = 40 " << endl
340 <<
"uincr = 1.0*(umax-umin)/n_isovalue" << endl
341 <<
"set cbtics uincr" << endl
342 <<
"set pm3d interpolate n_subdivide,n_subdivide corners2color mean" << endl
343 <<
"set palette rgbformulae 33,13,-4 maxcolors n_isovalue" << endl;
345 plot <<
"set noxlabel" << endl
346 <<
"set noylabel" << endl
347 <<
"set nozlabel" << endl;
349 plot <<
"set xlabel \"x\"" << endl
350 <<
"set ylabel \"y\"" << endl
351 <<
"set zlabel \"z\"" << endl;
356 plot <<
"set colors classic" << endl;
358 plot <<
"plot \"" << gdatname <<
"\" notitle with lines lw 2";
359 if (
gray || black_and_white) {
360 plot <<
" lc 0" << endl;
364 if (!fill && dim == 2) {
369 plot <<
" \"" << gdatname <<
"\" notitle";
371 if (black_and_white) {
372 plot <<
" with lines palette lw 2" << endl;
374 plot <<
" with lines palette" << endl;
377 plot <<
" with lines palette lw 2" << endl;
383 if (format ==
"" && !reader_on_stdin) {
384 plot <<
"pause -1 \"<return>\"\n";
393 command =
"gnuplot ";
394 if (reader_on_stdin) command +=
"-persist ";
395 command += tmp + basename +
".plot";
396 if (verbose) clog <<
"! " << command << endl;
398 status = system (command.c_str());
401 if (verbose) clog <<
"! file \"" << outfile_fmt <<
"\" created" << endl;
408 command =
"/bin/rm -f " + filelist;
409 if (verbose) clog <<
"! " << command << endl;
410 status = system (command.c_str());
422 Float vscale = iorheo::getvectorscale(ods.
os());
426 "gnuplot vector: unsupported non-isoparametric approx " << uh.
get_space().get_basis().name());
429 T diam_omega =
norm (omega.xmax() - omega.xmin());
430 T h_moy = diam_omega/
pow(n_elt,1./
d);
434 if (norm_max_uh + 1 == 1) norm_max_uh = 1;
436 T scale = vscale*(h_moy/norm_max_uh);
440 for (
size_type inod = 0, nnod = x.size(); inod < nnod; inod++) {
442 for (
size_type i_comp = 0; i_comp < n_comp; i_comp++) {
444 vi[i_comp] = uh.
dof(idof);
446 x[inod] += vscale*vi;
449 deformed_omega.set_nodes(x);
466 default:
error_macro (
"do not known how to print " << uh.
valued() <<
"-valued field");
473#define _RHEOLEF_instanciate(T) \
474template odiststream& visu_gnuplot<T> (odiststream&, const field_basic<T,sequential>&);
477#undef _RHEOLEF_instanciate
field::size_type size_type
see the Float page for the full documentation
see the disarray page for the full documentation
void initialize(const basis_basic< T > &fem_basis, const piola_on_pointset< T > &pops)
const geo_type & get_geo() const
const space_type & get_space() const
const std::string & valued() const
typename float_traits< T >::type float_type
generic mesh with rerefence counting
see the geo_element page for the full documentation
variant_type variant() const
see the integrate_option page for the full documentation
odiststream: see the diststream page for the full documentation
void initialize(const basis_basic< T > &piola_basis, const quadrature< T > &quad, const integrate_option &iopt)
static size_type ilat2loc_inod(size_type order, const point_basic< size_type > &ilat)
static size_type ilat2loc_inod(size_type order, const point_basic< size_type > &ilat)
static size_type ilat2loc_inod(size_type order, const point_basic< size_type > &ilat)
static const variant_type q
static const variant_type e
static const variant_type t
#define _RHEOLEF_instanciate(T)
#define error_macro(message)
check_macro(expr1.have_homogeneous_space(Xh1), "dual(expr1,expr2); expr1 should have homogeneous space. HINT: use dual(interpolate(Xh, expr1),expr2)")
verbose clean transpose logscale grid shrink ball stereo iso volume skipvtk deformation fastfieldload lattice reader_on_stdin gray
This file is part of Rheolef.
void 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, Eigen::Matrix< T, Eigen::Dynamic, 1 > &value)
odiststream & visu_gnuplot_vector(odiststream &ods, const field_basic< T, sequential > &uh)
std::string get_tmpdir()
get_tmpdir: see the rheostream page for the full documentation
space_mult_list< T, M > pow(const space_basic< T, M > &X, size_t n)
void put(std::ostream &out, std::string name, const tiny_matrix< T > &a)
odiststream & visu_gnuplot_scalar(odiststream &ods, const field_basic< T, sequential > &uh)
floorer_type< T > floorer(const T &prec)
ceiler_type< T > ceiler(const T &prec)
field_basic< T, M > interpolate(const space_basic< T, M > &V2h, const field_basic< T, M > &u1h)
see the interpolate page for the full documentation
bool file_exists(const std::string &filename)
file_exists: see the rheostream page for the full documentation
odiststream & visu_gnuplot(odiststream &, const field_basic< T, sequential > &)
T norm(const vec< T, M > &x)
norm(x): see the expression page for the full documentation