Hermes2D  3.0
vector_view.cpp
1 // This file is part of Hermes2D.
2 //
3 // Hermes2D is free software: you can redistribute it and/or modify
4 // it under the terms of the GNU General Public License as published by
5 // the Free Software Foundation, either version 2 of the License, or
6 // (at your option) any later version.
7 //
8 // Hermes2D is distributed in the hope that it will be useful,
9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
11 // GNU General Public License for more details.
12 //
13 // You should have received a copy of the GNU General Public License
14 // along with Hermes2D. If not, see <http://www.gnu.org/licenses/>.
15 
16 #ifndef NOGLUT
17 
18 #include <GL/freeglut.h>
19 #include "global.h"
20 #include "vector_view.h"
21 
22 namespace Hermes
23 {
24  namespace Hermes2D
25  {
26  namespace Views
27  {
28  VectorView::VectorView(const char* title, WinGeom* wg)
29  : View(title, wg)
30  {
31  vec = new Vectorizer(OpenGL);
32  gx = gy = 0.0;
33  gs = 20.0;
34  hexa = true;
35  mode = 0;
36  lines = false;
37  pmode = false;
38  length_coef = 1.0;
39  }
40 
41  VectorView::VectorView(char* title, WinGeom* wg)
42  : View(title, wg)
43  {
44  vec = new Vectorizer(OpenGL);
45  gx = gy = 0.0;
46  gs = 20.0;
47  hexa = true;
48  mode = 0;
49  lines = false;
50  pmode = false;
51  length_coef = 1.0;
52  }
53 
54  VectorView::~VectorView()
55  {
56  delete vec;
57  }
58 
59  void VectorView::show(MeshFunctionSharedPtr<double> vsln)
60  {
61  if (vsln->get_num_components() < 2)
62  throw Hermes::Exceptions::Exception("The single-argument version of show() is only for vector-valued solutions.");
63  show(vsln, vsln, H2D_FN_VAL_0, H2D_FN_VAL_1);
64  }
65 
66  void VectorView::show(MeshFunctionSharedPtr<double> xsln, MeshFunctionSharedPtr<double> ysln, int xitem, int yitem, MeshFunctionSharedPtr<double> xdisp, MeshFunctionSharedPtr<double> ydisp, double dmult)
67  {
68  if (xsln != nullptr && ysln != nullptr && xsln == ysln)
69  this->warn("Identical solutions passed to the two-argument version of show(). Most likely this is a mistake.");
70 
71  vec->set_displacement(xdisp, ydisp, dmult);
72 
73  vec->lock_data();
74  MeshFunctionSharedPtr<double> slns[2] = { xsln, ysln };
75  int items[2] = { xitem, yitem };
76  vec->process_solution(slns, items);
77 
78  if (range_auto)
79  {
80  range_min = vec->get_min_value();
81  range_max = vec->get_max_value();
82  }
83  vec->calc_vertices_aabb(&vertices_min_x, &vertices_max_x, &vertices_min_y, &vertices_max_y);
84  vec->unlock_data();
85 
86  create();
87  update_layout();
88  reset_view(false);
89 
90  refresh();
91 
92  wait_for_draw();
93  }
94 
95  static int n_vert(int i) { return (i + 1) % 3; }
96  static int p_vert(int i) { return (i + 2) % 3; }
97 
98  void VectorView::set_mode(int mode)
99  {
100  this->mode = mode % 3;
101  refresh();
102  }
103 
104  Vectorizer* VectorView::get_vectorizer()
105  {
106  return this->vec;
107  }
108 
109  void VectorView::set_vectorizer_criterion(LinearizerCriterion criterion)
110  {
111  this->vec->set_criterion(criterion);
112  }
113 
114  void VectorView::plot_arrow(double x, double y, double xval, double yval, double max, double min, double gs)
115  {
116  if (mode == 1)
117  glColor3f(0.0f, 0.0f, 0.0f);
118  else
119  glColor3f(0.5f, 0.5f, 0.5f);
120 
121  // magnitude
122  double Real_mag = sqrt(sqr(xval) + sqr(yval));
123  double mag = Real_mag;
124  if (Real_mag > max) mag = max;
125  double length = mag / max * gs * length_coef;
126  double width = 0.1 * gs;
127  if (mode == 1) width *= 1.2;
128  double xnew = x + gs * xval * mag / (max * Real_mag) * length_coef;
129  double ynew = y - gs * yval * mag / (max * Real_mag) * length_coef;
130 
131  if ((mag) / (max - min) < 1e-5)
132  {
133  glTranslated(x, y, 0.0);
134 
135  glBegin(GL_QUADS);
136  glVertex2d(width, width);
137  glVertex2d(width, -width);
138  glVertex2d(-width, -width);
139  glVertex2d(-width, width);
140  glEnd();
141  }
142  else
143  {
144  glBegin(GL_LINES);
145  glVertex2d(x, y);
146  glVertex2d(xnew, ynew);
147  glEnd();
148 
149  glTranslated(x, y, 0.0);
150  glRotated(atan2(-yval, xval) * 180.0 / M_PI, 0.0, 0.0, 1.0);
151 
152  glBegin(GL_TRIANGLES);
153  glVertex2d(length + 3 * width, 0.0);
154  glVertex2d(length - 2 * width, width);
155  glVertex2d(length - 2 * width, -width);
156  glEnd();
157  }
158 
159  glLoadIdentity();
160 
161  if (mode == 1)
162  {
163  float color[3];
164  // 0.0 -- 1.0
165  get_palette_color((mag - min) / (max - min), color);
166  glColor3f(color[0], color[1], color[2]);
167 
168  if (mag / (max - min) < 1e-5)
169  {
170  glBegin(GL_QUADS);
171  glVertex2d(width, width);
172  glVertex2d(width, -width);
173  glVertex2d(-width, -width);
174  glVertex2d(-width, width);
175  glEnd();
176  }
177  else
178  {
179  glBegin(GL_LINES);
180  glVertex2d(x, y);
181  glVertex2d(xnew, ynew);
182  glEnd();
183 
184  glTranslated(x - 1, y, 0.0);
185  glRotated(atan2(-yval, xval) * 180.0 / M_PI, 0.0, 0.0, 1.0);
186 
187  glBegin(GL_TRIANGLES);
188  glVertex2d(length + 3 * width, 0.0);
189  glVertex2d(length - 2 * width, width);
190  glVertex2d(length - 2 * width, -width);
191  glEnd();
192 
193  glLoadIdentity();
194  }
195  }
196  }
197 
198  void VectorView::draw_edges_2d()
199  {
200  glColor3f(0.5, 0.5, 0.5);
201  glBegin(GL_LINES);
202  for (Vectorizer::Iterator<VectorLinearizerDataDimensions<LINEARIZER_DATA_TYPE>::edge_t> it = vec->edges_begin(); !it.end; ++it)
203  {
205  int& edge_marker = it.get_marker();
206 
207  if (lines || edge_marker)
208  {
209  glVertex2d(transform_x(edge[0][0]), transform_y(edge[0][1]));
210  glVertex2d(transform_x(edge[1][0]), transform_y(edge[1][1]));
211  }
212  }
213  glEnd();
214  }
215 
216  void VectorView::on_display()
217  {
218  set_ortho_projection();
219  glDisable(GL_LIGHTING);
220  glDisable(GL_DEPTH_TEST);
221  glDisable(GL_TEXTURE_1D);
222  glPolygonMode(GL_FRONT_AND_BACK, pmode ? GL_LINE : GL_FILL);
223 
224  // initial grid point and grid step
225  double gt = gs;
226  if (hexa) gt *= sqrt(3.0) / 2.0;
227 
228  double max_length = 0.0;
229 
230  // transform all vertices
231  vec->lock_data();
232 
233  // value range
234  double min = range_min, max = range_max;
235  if (range_auto) { min = vec->get_min_value(); max = vec->get_max_value(); }
236  double irange = 1.0 / (max - min);
237  // special case: constant solution
238  if (fabs(min - max) < Hermes::HermesSqrtEpsilon) { irange = 1.0; min -= 0.5; }
239 
240  // draw all triangles
241  if (mode != 1)
242  glEnable(GL_TEXTURE_1D);
243  glBindTexture(GL_TEXTURE_1D, gl_pallete_tex_id);
244  glBegin(GL_TRIANGLES);
245  glColor3f(0.95f, 0.95f, 0.95f);
246 
247  for (Vectorizer::Iterator<VectorLinearizerDataDimensions<LINEARIZER_DATA_TYPE>::triangle_t> it = vec->triangles_begin(); !it.end; ++it)
248  {
249  VectorLinearizerDataDimensions<LINEARIZER_DATA_TYPE>::triangle_t& triangle = it.get();
250 
251  double mag = sqrt(sqr(triangle[0][2]) + sqr(triangle[0][3]));
252  glTexCoord2d((mag - min) * irange * tex_scale + tex_shift, 0.0);
253  glVertex2d(transform_x(triangle[0][0]), transform_y(triangle[0][1]));
254 
255  double length = sqr(triangle[0][2]) + sqr(triangle[0][3]);
256  if (length > max_length)
257  max_length = length;
258 
259  mag = sqrt(sqr(triangle[1][2]) + sqr(triangle[1][3]));
260  glTexCoord2d((mag - min) * irange * tex_scale + tex_shift, 0.0);
261  glVertex2d(transform_x(triangle[1][0]), transform_y(triangle[1][1]));
262 
263  length = sqr(triangle[1][2]) + sqr(triangle[1][3]);
264  if (length > max_length)
265  max_length = length;
266 
267  mag = sqrt(sqr(triangle[2][2]) + sqr(triangle[2][3]));
268  glTexCoord2d((mag - min) * irange * tex_scale + tex_shift, 0.0);
269  glVertex2d(transform_x(triangle[2][0]), transform_y(triangle[2][1]));
270 
271  length = sqr(triangle[2][2]) + sqr(triangle[2][3]);
272  if (length > max_length)
273  max_length = length;
274  }
275 
276  max_length = sqrt(max_length);
277  glEnd();
278  glDisable(GL_TEXTURE_1D);
279 
280  draw_edges_2d();
281 
282  // draw arrows
283  if (mode != 2)
284  {
285  for (Vectorizer::Iterator<VectorLinearizerDataDimensions<LINEARIZER_DATA_TYPE>::triangle_t> it = vec->triangles_begin(); !it.end; ++it)
286  {
287  VectorLinearizerDataDimensions<LINEARIZER_DATA_TYPE>::triangle_t& triangle = it.get();
288  {
289  double miny = 1e100;
290  int idx = -1, k, l1, l2, r2, r1, s;
291  double lry, x;
292  double mr, ml, lx, rx, xval, yval;
293 
294  double wh = output_height + gt, ww = output_width + gs;
295 
296  double vert[3][4] = {
297  {
298  triangle[0][0],
299  triangle[0][1],
300  triangle[0][2],
301  triangle[0][3]
302  },
303  {
304  triangle[1][0],
305  triangle[1][1],
306  triangle[1][2],
307  triangle[1][3]
308  },
309  {
310  triangle[2][0],
311  triangle[2][1],
312  triangle[2][2],
313  triangle[2][3]
314  }
315  };
316 
317  double tvert[3][2] = {
318  {
319  transform_x(triangle[0][0]),
320  transform_y(triangle[0][1])
321  },
322  {
323  transform_x(triangle[1][0]),
324  transform_y(triangle[1][1])
325  },
326  {
327  transform_x(triangle[2][0]),
328  transform_y(triangle[2][1])
329  }
330  };
331 
332  if ((tvert[0][0] < -gs) && (tvert[1][0] < -gs) && (tvert[2][0] < -gs)) continue;
333  if ((tvert[0][0] > ww) && (tvert[1][0] > ww) && (tvert[2][0] > ww)) continue;
334  if ((tvert[0][1] < -gt) && (tvert[1][1] < -gt) && (tvert[2][1] < -gt)) continue;
335  if ((tvert[0][1] > wh) && (tvert[1][1] > wh) && (tvert[2][1] > wh)) continue;
336 
337  // find vertex with min y-coordinate
338  for (k = 0; k < 3; k++)
339  {
340  if (tvert[k][1] < miny)
341  miny = tvert[idx = k][1];
342  }
343  l1 = r1 = idx;
344  l2 = n_vert(idx);
345  r2 = p_vert(idx);
346 
347  // plane of x and y values on triangle
348  double a[2], b[2], c[2], d[2];
349  for (int n = 0; n < 2; n++)
350  {
351  a[n] = (tvert[l1][1] - tvert[l2][1])*(vert[r1][2 + n] - vert[r2][2 + n]) - (vert[l1][2 + n] - vert[l2][2 + n])*(tvert[r1][1] - tvert[r2][1]);
352  b[n] = (vert[l1][2 + n] - vert[l2][2 + n])*(tvert[r1][0] - tvert[r2][0]) - (tvert[l1][0] - tvert[l2][0])*(vert[r1][2 + n] - vert[r2][2 + n]);
353  c[n] = (tvert[l1][0] - tvert[l2][0])*(tvert[r1][1] - tvert[r2][1]) - (tvert[l1][1] - tvert[l2][1])*(tvert[r1][0] - tvert[r2][0]);
354  d[n] = -a[n] * tvert[l1][0] - b[n] * tvert[l1][1] - c[n] * vert[l1][2 + n];
355  a[n] /= c[n]; b[n] /= c[n]; d[n] /= c[n];
356  }
357 
358  // first step
359  s = (int)ceil((tvert[l1][1] - gy) / gt);
360  lry = gy + s*gt;
361  bool shift = hexa && (s & 1);
362 
363  // if there are two points with min y-coordinate, switch to the next segment
364  if ((tvert[l1][1] == tvert[l2][1]) || (tvert[r1][1] == tvert[r2][1]))
365  {
366  if (tvert[l1][1] == tvert[l2][1])
367  {
368  l1 = l2;
369  l2 = r2;
370  }
371  else if (tvert[r1][1] == tvert[r2][1])
372  {
373  r1 = r2;
374  r2 = l2;
375  }
376  }
377 
378  // slope of the left and right segment
379  ml = (tvert[l1][0] - tvert[l2][0]) / (tvert[l1][1] - tvert[l2][1]);
380  mr = (tvert[r1][0] - tvert[r2][0]) / (tvert[r1][1] - tvert[r2][1]);
381  // x-coordinates of the endpoints of the first line
382  lx = tvert[l1][0] + ml * (lry - (tvert[l1][1]));
383  rx = tvert[r1][0] + mr * (lry - (tvert[r1][1]));
384 
385  if (lry < -gt)
386  {
387  k = (int)floor(-lry / gt);
388  lry += gt * k;
389  lx += k * ml * gt;
390  rx += k * mr * gt;
391  }
392 
393  // while we are in triangle
394  while (((lry < tvert[l2][1]) || (lry < tvert[r2][1])) && (lry < wh))
395  {
396  // while we are in the segment
397  while (((lry <= tvert[l2][1]) && (lry <= tvert[r2][1])) && (lry < wh))
398  {
399  double gz = gx;
400  if (shift) gz -= 0.5*gs;
401  s = (int)ceil((lx - gz) / gs);
402  x = gz + s*gs;
403  if (hexa) shift = !shift;
404 
405  if (x < -gs)
406  {
407  k = (int)floor(-x / gs);
408  x += gs * k;
409  }
410  // go along the line
411  while ((x < rx) && (x < ww))
412  {
413  // plot the arrow
414  xval = -a[0] * x - b[0] * lry - d[0];
415  yval = -a[1] * x - b[1] * lry - d[1];
416  plot_arrow(x, lry, xval, yval, max, min, gs);
417  x += gs;
418  }
419  // move to the next line
420  lx += ml*gt;
421  rx += mr*gt;
422  lry += gt;
423  }
424  // change segment
425  if (lry >= tvert[l2][1])
426  {
427  l1 = l2; l2 = r2;
428  ml = (tvert[l1][0] - tvert[l2][0]) / (tvert[l1][1] - tvert[l2][1]);
429  lx = tvert[l1][0] + ml * (lry - (tvert[l1][1]));
430  }
431  else
432  {
433  r1 = r2; r2 = l2;
434  mr = (tvert[r1][0] - tvert[r2][0]) / (tvert[r1][1] - tvert[r2][1]);
435  rx = tvert[r1][0] + mr * (lry - (tvert[r1][1]));
436  }
437  }
438  }
439  }
440 
441  vec->unlock_data();
442  }
443  }
444 
445  void VectorView::on_mouse_move(int x, int y)
446  {
447  if (dragging)
448  {
449  gx += (x - mouse_x);
450  gy += (y - mouse_y);
451  }
452  View::on_mouse_move(x, y);
453  }
454 
455  void VectorView::on_key_down(unsigned char key, int x, int y)
456  {
457  switch (key)
458  {
459  case 'm':
460  lines = !lines;
461  refresh();
462  break;
463 
464  case 'l':
465  pmode = !pmode;
466  refresh();
467  break;
468 
469  case 'c':
470  reset_view(true);
471  refresh();
472  break;
473 
474  case 'f':
475  set_palette_filter(pal_filter != GL_LINEAR);
476  break;
477 
478  case 'x':
479  set_grid_type(!hexa);
480  break;
481 
482  case 'b':
483  mode++;
484  if (mode > 2) mode = 0;
485  refresh();
486  break;
487 
488  case '*':
489  case '/':
490  if (key == '*') length_coef *= 1.1; else length_coef /= 1.1;
491  refresh();
492  break;
493 
494  default:
495  View::on_key_down(key, x, y);
496  break;
497  }
498  }
499 
500  const char* VectorView::get_help_text() const
501  {
502  return
503  "VectorView\n\n"
504  "Controls:\n"
505  " Left mouse - pan\n"
506  " Right mouse - zoom\n"
507  " B - toggle view mode (type of arrows x no arrows)\n"
508  " * - extend arrows\n"
509  " / - shorten arrows\n"
510  " C - center image\n"
511  " F - toggle smooth palette\n"
512  " X - toggle hexagonal grid\n"
513  " H - render high-quality frame\n"
514  " M - toggle mesh\n"
515  " P - cycle palettes\n"
516  " S - save screenshot\n"
517  " F1 - this help\n"
518  " Esc, Q - quit";
519  }
520  }
521  }
522 }
523 #endif
Definition: adapt.h:24
Common definitions for Hermes2D.
LinearizerMultidimensional< VectorLinearizerDataDimensions< LINEARIZER_DATA_TYPE > > Vectorizer
Linearizer for vector cases - historically called Vectorizer.
Definition: linearizer.h:203
Typedefs used throughout the Linearizer functionality.
void set_criterion(LinearizerCriterion criterion)
Definition: linearizer.cpp:68
Abstract class for criterion according to which the linearizer stops dividing elements at some point ...