line_visualiser.h
Go to the documentation of this file.
1 //LIC// ====================================================================
2 //LIC// This file forms part of oomph-lib, the object-oriented,
3 //LIC// multi-physics finite-element library, available
4 //LIC// at http://www.oomph-lib.org.
5 //LIC//
6 //LIC// Version 1.0; svn revision $LastChangedRevision: 1182 $
7 //LIC//
8 //LIC// $LastChangedDate: 2016-05-20 16:50:20 +0100 (Fri, 20 May 2016) $
9 //LIC//
10 //LIC// Copyright (C) 2006-2016 Matthias Heil and Andrew Hazel
11 //LIC//
12 //LIC// This library is free software; you can redistribute it and/or
13 //LIC// modify it under the terms of the GNU Lesser General Public
14 //LIC// License as published by the Free Software Foundation; either
15 //LIC// version 2.1 of the License, or (at your option) any later version.
16 //LIC//
17 //LIC// This library is distributed in the hope that it will be useful,
18 //LIC// but WITHOUT ANY WARRANTY; without even the implied warranty of
19 //LIC// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
20 //LIC// Lesser General Public License for more details.
21 //LIC//
22 //LIC// You should have received a copy of the GNU Lesser General Public
23 //LIC// License along with this library; if not, write to the Free Software
24 //LIC// Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
25 //LIC// 02110-1301 USA.
26 //LIC//
27 //LIC// The authors may be contacted at oomph-lib@maths.man.ac.uk.
28 //LIC//
29 //LIC//====================================================================
30 //Header file for line visualiser
31 
32 //Include guard to prevent multiple inclusions of the header
33 #ifndef OOMPH_LINE_VISUALISER_HEADER
34 #define OOMPH_LINE_VISUALISER_HEADER
35 
36 // Config header generated by autoconfig
37 #ifdef HAVE_CONFIG_H
38 #include <oomph-lib-config.h>
39 #endif
40 
41 #include<fstream>
42 #include<iostream>
43 
44 
45 namespace oomph
46 {
47 
48 //====================================================================
49 /// \short Class to aid visualisation of the values on a set
50 /// of points. NOTE: in a distributed problem, output is only done
51 /// on processor 0.
52 //====================================================================
54  {
55 
56  public :
57 
58 
59  /// \short Constructor: Pass pointer to mesh and coordinates of
60  /// desired plot points: coord_vec[j][i] is the i-th coordinate of
61  /// the j-th plot point. Optional final parameter specifies the
62  /// maximum search radius in bin when locating the plot points.
63  /// Defaults to DBL_MAX and will therefore keep searching through
64  /// the entire bin structure if a point cannot be found. It's worth
65  /// limiting this to the order of the size of a typical element
66  /// in the mesh in parallel computations with fine meshes, otherwise
67  /// the setup can take forever.
68  LineVisualiser(Mesh* mesh_pt, const Vector<Vector<double> >& coord_vec,
69  const double& max_search_radius=DBL_MAX) :
70  Max_search_radius(max_search_radius),
71  Comm_pt(mesh_pt->communicator_pt())
72  {
73  // Do the actual work
74  setup(mesh_pt,coord_vec);
75  }
76 
77 
78  /// \short Constructor reading centerline file
79  /// - Open "file_name" and extract 3 first doubles of each line
80  /// - Skip lines which does not begin with a number. Scaling
81  /// factor allows points defined in input file to be scaled.
82  LineVisualiser(Mesh* mesh_pt, const std::string file_name,
83  const double &scale = 1.0) :
84  Max_search_radius(DBL_MAX),
85  Comm_pt(mesh_pt->communicator_pt())
86  {
87  setup_from_file(mesh_pt,file_name,scale);
88  }
89 
90 
91  /// \short Constructor reading centerline file
92  /// - Open "file_name" and extract 3 first doubles of each line
93  /// - Skip lines which does not begin with a number. Scaling
94  /// factor allows points defined in input file to be scaled.
95  /// Second parameter specifies the
96  /// maximum search radius in bin when locating the plot points. It's worth
97  /// setting this to the order of the size of a typical element
98  /// in the mesh in parallel computations with fine meshes, otherwise
99  /// the setup can take forever.
100  LineVisualiser(Mesh* mesh_pt, const double& max_search_radius,
101  const std::string file_name, const double &scale = 1.0) :
102  Max_search_radius(max_search_radius),
103  Comm_pt(mesh_pt->communicator_pt())
104  {
105  setup_from_file(mesh_pt,file_name,scale);
106  }
107 
108 
109 
110 
111  /// \short Output function: output each plot point.
112  /// NOTE: in a distributed problem, output is only done
113  /// on processor 0.
114  void output(std::ostream &outfile)
115  {
116  // Get data in array
118  get_output_data(data);
119 
120 
121  // Loop over the points
122  for (unsigned i=0; i<Nplot_points; i++)
123  {
124  // Get the size of the line
125  unsigned n=data[i].size();
126 
127  // Loop over the values on the line
128  for (unsigned j=0;j<n;j++)
129  {
130  outfile << data[i][j] << " ";
131  }
132  if (n > 0)
133  {
134  outfile << std::endl;
135  }
136  }
137  }
138 
139  /// \short Output data function: store data associated with each
140  /// plot point in data array
142  {
143  // Resize output data array
144  data.resize(Nplot_points);
145 
146  // Check if mesh is distributed and a communication pointer
147  // exists.
148  if(Comm_pt!=0)
149  {
150  int nproc=Comm_pt->nproc();
151 
152  if (nproc>1)
153  {
154 
155 #ifdef OOMPH_HAS_MPI
156 
157  // Declaration of MPI variables
158  MPI_Status stat;
159  int tag=0;
160  int my_rank=Comm_pt->my_rank();
161 
162 
163  // Buffer
164  unsigned buff_size;
165 
166  // Create array which contains data found in every process
168 
169  // Loop over the points to fill in vec
170  for (unsigned i=0; i<Nplot_points; i++)
171  {
172  // Check if the point was found in the mesh
173  if (Plot_point[i].first != NULL) // success
174  {
175  // Check if the point is halo
176  if (!((*Plot_point[i].first).is_halo()))
177  {
178  // Get the line of output data from the element
179  // (specified by .first), at its local coordinate
180  // (specified by .second)
181  Plot_point[i].first->point_output_data(Plot_point[i].second,
182  vec[i]);
183  }
184  }
185  }
186 
187 
188  // Analyse which plot points have been found
189  // locally and concatenate the data:
190 
191  // This contains the flat-packed doubles to be sent
192  // for all located plot points
193  Vector<double> local_values;
194 
195  // Number of values to be sent for each plot point
196  // (almost certainly the same for all plot points, but...)
197  // size_values[i] gives the number of doubles to be
198  // sent for plot point i.
199  Vector<unsigned> size_values;
200 
201  // Each processor indicates if it has found a given plot point.
202  // Once this is gathered on the root processor we know
203  // exactly which data we'll receive from where.
204  Vector<unsigned> tmp_proc_point_found_plus_one(Nplot_points,0);
205 
206  // Loop over the plot points
207  for (unsigned i=0; i<Nplot_points; i++)
208  {
209  unsigned ndata=vec[i].size();
210  if (ndata!=0)
211  {
212  // Store the number of fields
213  size_values.push_back(ndata);
214 
215  // Update found vector
216  tmp_proc_point_found_plus_one[i]=my_rank+1;
217 
218  // Store values
219  for (unsigned j=0;j<ndata;j++)
220  {
221  local_values.push_back(vec[i][j]);
222  }
223  }
224  }
225 
226  // hierher: Get rid of mpi Helpers
227 
228  // Gather information on root
229 
230  // Find out who's found the points
231  Vector<unsigned> proc_point_found_plus_one(Nplot_points,0);
232  MPI_Reduce(&tmp_proc_point_found_plus_one[0],
233  &proc_point_found_plus_one[0],
234  Nplot_points, MPI_UNSIGNED, MPI_MAX, 0,
235  Comm_pt->mpi_comm());
236 
237 
238  // Main process write data
239  if (my_rank == 0)
240  {
241  // Collect all the data
242  Vector<Vector<double> > received_data(nproc-1);
243  Vector<Vector<unsigned> > received_size(nproc-1);
244  Vector<unsigned> counter_d(nproc-1,0);
245  Vector<unsigned> counter_s(nproc-1,0);
246 
247  // Loop over processors that send their points
248  for (int i=1; i<nproc; i++)
249  {
250  // Receive sizes of data
251  MPI_Recv(&buff_size, 1, MPI_UNSIGNED, i,
252  tag, Comm_pt->mpi_comm(), &stat);
253  received_size[i-1].resize(std::max(unsigned(1),buff_size));
254  MPI_Recv(&received_size[i-1][0], buff_size, MPI_UNSIGNED, i,
255  tag, Comm_pt->mpi_comm(), &stat);
256 
257  // Receive actual data
258  MPI_Recv(&buff_size, 1, MPI_UNSIGNED, i,
259  tag, Comm_pt->mpi_comm(), &stat);
260  received_data[i-1].resize(std::max(unsigned(1),buff_size));
261  MPI_Recv(&received_data[i-1][0], buff_size, MPI_DOUBLE, i,
262  tag, Comm_pt->mpi_comm(), &stat);
263  }
264 
265  // Analyse data for each point
266  for (unsigned i=0; i<Nplot_points; i++)
267  {
268  // Somebody has found it
269  if (proc_point_found_plus_one[i] != 0)
270  {
271  // Root processor has found it
272  if (proc_point_found_plus_one[i] == 1)
273  {
274  // Copy directly from vec vector
275  data[i]=vec[i];
276  }
277  // Another (non-root) processor has found it
278  else
279  {
280  unsigned line_i=proc_point_found_plus_one[i]-2;
281 
282  // Resize data line
283  data[i].resize(received_size[line_i][counter_s[line_i] ]);
284 
285  // Copy values
286  for (unsigned j=0;j<received_size[line_i][counter_s[line_i] ];
287  j++)
288  {
289  data[i][j]=received_data[line_i][counter_d[line_i]+j];
290  }
291 
292  //Increase counter
293  counter_d[line_i]+=received_size[line_i][counter_s[line_i] ];
294  counter_s[line_i]++;
295  }
296  } // end somebody has found it -- no output at all if nobody
297  // has found the point (e.g. outside mesh)
298  }
299  }
300  // Send data to root
301  else
302  {
303  //Send the number of fields to the main process
304  buff_size = size_values.size();
305  MPI_Send(&buff_size, 1, MPI_UNSIGNED, 0, tag,
306  Comm_pt->mpi_comm());
307 
308  //Send the sizes of fields to the main process
309  if (buff_size==0) size_values.resize(1);
310  MPI_Send(&size_values[0], buff_size, MPI_UNSIGNED, 0, tag,
311  Comm_pt->mpi_comm());
312 
313  //Send the number of data fields to the main process
314  buff_size = local_values.size();
315  MPI_Send(&buff_size, 1, MPI_UNSIGNED, 0, tag,
316  Comm_pt->mpi_comm());
317 
318  //Send the data to the main process
319  if (buff_size==0) local_values.resize(1);
320  MPI_Send(&local_values[0], buff_size, MPI_DOUBLE, 0, tag,
321  Comm_pt->mpi_comm());
322  }
323 
324 #endif // Serial version
325  }
326  }
327  else
328  {
329 
330  // Loop over the points
331  for (unsigned i=0; i<Nplot_points; i++)
332  {
333  // Check if the point was found in the mesh
334  if (Plot_point[i].first != NULL) // success
335  {
336  // Copy line into data array
337  Plot_point[i].first->point_output_data(Plot_point[i].second,
338  data[i]);
339  }
340  else // not found -- keep empty block there for debugging
341  {
342  //oomph_info << "Point " << i << " not found\n";
343  }
344  }
345  }
346  }
347 
348 
349  /// \short Update plot points coordinates (in preparation of remesh,
350  /// say).
352  {
353  // Resize coord_vec
354  coord_vec.resize(Nplot_points);
355 
356  // Check that the communication pointer is initialised and the
357  // problem is distributed.
358  if(Comm_pt!=0)
359  {
360  int nproc=Comm_pt->nproc();
361  if (nproc>1)
362  {
363 
364 #ifdef OOMPH_HAS_MPI
365 
366  // Declaration of MPI variables
367  MPI_Status stat;
368  int tag; // cgj: tag should be initialised before use
369  int my_rank=Comm_pt->my_rank();
370 
371 
372  // Buffer
373  unsigned buff_size;
374 
375  // Create array which contains data found in every process
377 
378  for (unsigned i=0; i<Nplot_points; i++)
379  {
380  if (Plot_point[i].first != NULL)
381  {
382  if (!((*Plot_point[i].first).is_halo()))
383  {
384  unsigned dim = Plot_point[i].second.size();
385 
386  vec[i].resize(dim);
387 
388  for (unsigned j=0; j<dim; j++)
389  {
390  vec[i][j]=Plot_point[i].first->
391  interpolated_x(Plot_point[i].second,j);
392  }
393  }
394  }
395  }
396 
397 
398  // Analyse which plot points have been found
399  // locally and concatenate the data:
400 
401  // This contains the flat-packed doubles to be sent
402  // for all located plot points
403  Vector<double> local_values;
404 
405  // Number of values to be sent for each plot point
406  // (almost certainly the same for all plot points, but...)
407  // size_values[i] gives the number of doubles to be
408  // sent for plot point i.
409  Vector<unsigned> size_values;
410 
411  // Each processor indicates if it has found a given plot point.
412  // Once this is gathered on the root processor we know
413  // exactly which data we'll receive from where.
414  Vector<unsigned> tmp_proc_point_found_plus_one(Nplot_points,0);
415 
416  // Loop over the plot points
417  for (unsigned i=0; i<Nplot_points; i++)
418  {
419  unsigned ndata=vec[i].size();
420  if (ndata!=0)
421  {
422  // Store the number of fields
423  size_values.push_back(ndata);
424 
425  // Update found vector
426  tmp_proc_point_found_plus_one[i]=my_rank+1;
427 
428 
429  // Store values
430  for (unsigned j=0;j<ndata;j++)
431  {
432  local_values.push_back(vec[i][j]);
433  }
434  }
435  }
436 
437  // Gather information on root
438 
439  // Find out who's found the points
440  Vector<unsigned> proc_point_found_plus_one(Nplot_points,0);
441  MPI_Reduce(&tmp_proc_point_found_plus_one[0],
442  &proc_point_found_plus_one[0],
443  Nplot_points, MPI_UNSIGNED, MPI_MAX, 0,
444  Comm_pt->mpi_comm());
445 
446  // Main process write data
447  if (my_rank == 0)
448  {
449  // Collect all the data
450  Vector<Vector<double> > received_data(nproc-1);
451  Vector<Vector<unsigned> > received_size(nproc-1);
452  Vector<unsigned> counter_d(nproc-1,0);
453  Vector<unsigned> counter_s(nproc-1,0);
454 
455  // Loop over processors that send their points
456  for (int i=1; i<nproc; i++)
457  {
458  // Receive sizes of data
459  MPI_Recv(&buff_size, 1, MPI_UNSIGNED, i,
460  tag, Comm_pt->mpi_comm(), &stat);
461  received_size[i-1].resize(std::max(unsigned(1),buff_size));
462  MPI_Recv(&received_size[i-1][0], buff_size, MPI_UNSIGNED, i,
463  tag, Comm_pt->mpi_comm(), &stat);
464 
465  // Receive actual data
466  MPI_Recv(&buff_size, 1, MPI_UNSIGNED, i,
467  tag, Comm_pt->mpi_comm(), &stat);
468  received_data[i-1].resize(std::max(unsigned(1),buff_size));
469  MPI_Recv(&received_data[i-1][0], buff_size, MPI_DOUBLE, i,
470  tag, Comm_pt->mpi_comm(), &stat);
471  }
472 
473  // Analyse data for each point
474  for (unsigned i=0; i<Nplot_points; i++)
475  {
476  // Somebody has found it
477  if (proc_point_found_plus_one[i] != 0)
478  {
479  // Root processor has found it
480  if (proc_point_found_plus_one[i] == 1)
481  {
482  // Copy directly from vec vector
483  coord_vec[i]=vec[i];
484  }
485  // Another (non-root) processor has found it
486  else
487  {
488  unsigned line_i=proc_point_found_plus_one[i]-2;
489 
490  // Resize data line
491  coord_vec[i].resize(received_size[line_i][counter_s[line_i] ]);
492 
493  // Copy values
494  for (unsigned j=0;j<received_size[line_i][counter_s[line_i] ];
495  j++)
496  {
497  coord_vec[i][j]=received_data[line_i][counter_d[line_i]+j];
498  }
499 
500  //Increase counter
501  counter_d[line_i]+=received_size[line_i][counter_s[line_i] ];
502  counter_s[line_i]++;
503  }
504  } // end somebody has found it -- no output at all if nobody
505  // has found the point (e.g. outside mesh)
506  }
507  }
508  // Send data to root
509  else
510  {
511  //Send the number of fields to the main process
512  buff_size = size_values.size();
513  MPI_Send(&buff_size, 1, MPI_UNSIGNED, 0, tag,
514  Comm_pt->mpi_comm());
515 
516  //Send the sizes of fields to the main process
517  if (buff_size==0) size_values.resize(1);
518  MPI_Send(&size_values[0], buff_size, MPI_UNSIGNED, 0, tag,
519  Comm_pt->mpi_comm());
520 
521  //Send the number of data fields to the main process
522  buff_size = local_values.size();
523  MPI_Send(&buff_size, 1, MPI_UNSIGNED, 0, tag,
524  Comm_pt->mpi_comm());
525 
526  //Send the data to the main process
527  if (buff_size==0) local_values.resize(1);
528  MPI_Send(&local_values[0], buff_size, MPI_DOUBLE, 0, tag,
529  Comm_pt->mpi_comm());
530  }
531 
532 #endif // Serial version
533  }
534  }
535  else
536  {
538  }
539  }
540 
541  private:
542 
543  /// \short Max radius beyond which we stop searching the bin. Initialised
544  /// to DBL_MAX so keep going until the point is found or until
545  /// we've searched every single bin. Overwriting this means we won't search
546  /// in bins whose closest vertex is at a distance greater than
547  /// Max_search_radius from the point to be located.
549 
550  /// \short Pointer to communicator -- allows us to collect data on
551  /// processor 0 if the mesh is distributed.
553 
554  /// Helper function to setup from file
555  void setup_from_file(Mesh* mesh_pt, const std::string file_name,
556  const double &scale)
557  {
558  // Open file use ifstream
559  std::ifstream file_input(file_name.c_str(),std::ios_base::in);
560  if (!file_input)
561  {
562  std::ostringstream error_message;
563  error_message << "Cannot open file " << file_name << "\n";
564  throw OomphLibError(error_message.str(),
565  OOMPH_CURRENT_FUNCTION,
566  OOMPH_EXCEPTION_LOCATION);
567  }
568  if (!file_input.is_open())
569  {
570  std::ostringstream error_message;
571  error_message << "Cannot open file " << file_name << "\n";
572  throw OomphLibError(error_message.str(),
573  OOMPH_CURRENT_FUNCTION,
574  OOMPH_EXCEPTION_LOCATION);
575  }
576 
577  // Declaration of variables
578  std::string line;
579  Vector<Vector<double> > coord_vec_tmp; // Coord array
580 
581  // Loop over the lines of the input file
582  while(getline(file_input,line)/* != 0*/)
583  {
584  // Test if the first char of the line is a number
585  // using ascii enumeration of chars
586  if (isdigit(line[0]))
587  {
588  Vector<double> tmp(3);
589 
590  // Read the 3 first doubles of the line
591  // Return 3 if success and less if error
592  int n=sscanf(line.c_str(),"%lf %lf %lf",&tmp[0],&tmp[1],&tmp[2]);
593 
594  if (n == 3) // success
595  {
596  // Rescaling
597  for (unsigned i=0; i<3; i++)
598  {
599  tmp[i]*=scale;
600  }
601 
602  // Add the new point to the list
603  coord_vec_tmp.push_back(tmp);
604  }
605  else // error
606  {
607  oomph_info << "Line ignored \n";
608  }
609  }
610  }
611 
612  // Call to the helper function
613  setup(mesh_pt,coord_vec_tmp);
614  }
615 
616 
617 
618  /// \short Helper function to setup the output structures
619  void setup(Mesh* mesh_pt, const Vector<Vector<double> > &coord_vec)
620  {
621  // Read out number of plot points
622  Nplot_points=coord_vec.size();
623 
624  if (Nplot_points==0) return;
625 
626  // Keep track of unlocated plot points
627  unsigned count_not_found_local=0;
628 
629  // Dimension
630  unsigned dim=coord_vec[0].size();
631 
632  // Make space
633  Plot_point.resize(Nplot_points);
634 
635  // Transform mesh into a geometric object
636  // hierher change this to the mesh's own commmunicator
637  unsigned suppress_synchronisation_of_bins_flag=1;
638  MeshAsGeomObject mesh_geom_tmp(mesh_pt,Comm_pt,
639  suppress_synchronisation_of_bins_flag);
640 
641  // Limit the search radius
642  mesh_geom_tmp.max_search_radius()=Max_search_radius;
643 
644  // Loop over input points
645  double tt_start=TimingHelpers::timer();
646  //oomph_info << "Looking for " << Nplot_points << " plot points\n";
647 
648  for (unsigned i=0; i<Nplot_points; i++)
649  {
650  // Local coordinate of the plot point with its element
651  Vector<double> s(dim,0.0);
652 
653  // Pointer to GeomObject that contains the plot point
654  GeomObject* geom_pt=0;
655 
656  // Locate zeta
657  mesh_geom_tmp.locate_zeta(coord_vec[i],geom_pt,s);
658 
659  // Upcast GeomElement as a FiniteElement
660  FiniteElement* fe_pt=dynamic_cast<FiniteElement*>(geom_pt);
661 
662  // Another one not found locally...
663  if (fe_pt==0)
664  {
665  count_not_found_local++;
666  /* oomph_info << "NOT Found the one at " */
667  /* << coord_vec[i][0] << " " */
668  /* << coord_vec[i][1] << "\n"; */
669  }
670  else
671  {
672  /* oomph_info << "Found the one at " */
673  /* << coord_vec[i][0] << " " */
674  /* << coord_vec[i][1] << "\n"; */
675  }
676 
677  // Save result in a pair
678  Plot_point[i]=std::pair<FiniteElement*,Vector<double> >(fe_pt,s);
679  }
680 
681 
682  oomph_info << "Number of points not found locally: "
683  << count_not_found_local << std::endl;
684 
685  // Global equivalent (is overwritten below if mpi)
686  unsigned count_not_found=count_not_found_local;
687 
688  // Check communication pointer exists and problem is
689  // distributed.
690  if(Comm_pt!=0)
691  {
692  int nproc=Comm_pt->nproc();
693  if (nproc>1)
694  {
695 
696 #ifdef OOMPH_HAS_MPI
697 
698  // Declaration of MPI variables
699  int my_rank=Comm_pt->my_rank();
700 
701  // Each processor indicates if it has found a given plot point.
702  // Once this is gathered on the root processor we know
703  // exactly which data we'll receive from where.
704  Vector<unsigned> tmp_proc_point_found_plus_one(Nplot_points,0);
705 
706  // Loop over the plot points
707  for (unsigned i=0; i<Nplot_points; i++)
708  {
709  // Found locally?
710  if (Plot_point[i].first!=0)
711  {
712  tmp_proc_point_found_plus_one[i]=my_rank+1;
713  }
714  }
715 
716  // hierher: Get rid of mpi Helpers
717 
718  // Gather information on root
719 
720  // Find out who's found the points
721  Vector<unsigned> proc_point_found_plus_one(Nplot_points,0);
722  MPI_Reduce(&tmp_proc_point_found_plus_one[0],
723  &proc_point_found_plus_one[0],
724  Nplot_points, MPI_UNSIGNED, MPI_MAX, 0,
725  Comm_pt->mpi_comm());
726 
727 
728  // Main process analyses data
729  if (my_rank == 0)
730  {
731  // Analyse data for each point
732  count_not_found=0;
733  for (unsigned i=0; i<Nplot_points; i++)
734  {
735  // Nobody has found it
736  if (proc_point_found_plus_one[i] == 0)
737  {
738  count_not_found++;
739  }
740  }
741  }
742 
743  // Now tell everybody about it
744  MPI_Bcast(&count_not_found,1,MPI_UNSIGNED,0,
745  Comm_pt->mpi_comm());
746 
747 #endif
748  }
749  }
750  double tt_end=TimingHelpers::timer();
751  oomph_info
752  << "Total time for location of plot points in LineVisualiser: "
753  << tt_end-tt_start << " [" << count_not_found
754  << " plot points were not found (with max search radius="
755  << Max_search_radius << ")]\n";
756  }
757 
758  // Get coordinates of found points
760  {
761  data.resize(Nplot_points);
762  for (unsigned i=0; i<Nplot_points; i++)
763  {
764  if (Plot_point[i].first != NULL)
765  {
766  unsigned dim = Plot_point[i].second.size();
767 
768  data[i].resize(dim);
769 
770  for (unsigned j=0; j<dim; j++)
771  {
772  data[i][j]=Plot_point[i].first->interpolated_x(
773  Plot_point[i].second,j);
774  }
775  }
776  }
777 
778  }
779 
780 
781  /// \short Vector of pairs containing points to finite elements and
782  /// local coordinates
784 
785  /// Number of plot points
786  unsigned Nplot_points;
787 
788  }; //end of class
789 
790 }
791 
792 #endif
793 
void setup_from_file(Mesh *mesh_pt, const std::string file_name, const double &scale)
Helper function to setup from file.
unsigned Nplot_points
Number of plot points.
cstr elem_len * i
Definition: cfortran.h:607
LineVisualiser(Mesh *mesh_pt, const double &max_search_radius, const std::string file_name, const double &scale=1.0)
Constructor reading centerline file.
A general Finite Element class.
Definition: elements.h:1271
double Max_search_radius
Max radius beyond which we stop searching the bin. Initialised to DBL_MAX so keep going until the poi...
void locate_zeta(const Vector< double > &zeta, GeomObject *&sub_geom_object_pt, Vector< double > &s, const bool &use_coordinate_as_initial_guess=false)
Find the sub geometric object and local coordinate therein that corresponds to the intrinsic coordina...
void update_plot_points_coordinates(Vector< Vector< double > > &coord_vec)
Update plot points coordinates (in preparation of remesh, say).
OomphInfo oomph_info
Class to aid visualisation of the values on a set of points. NOTE: in a distributed problem...
void get_output_data(Vector< Vector< double > > &data)
Output data function: store data associated with each plot point in data array.
void get_local_plot_points_coordinates(Vector< Vector< double > > &data)
static char t char * s
Definition: cfortran.h:572
LineVisualiser(Mesh *mesh_pt, const std::string file_name, const double &scale=1.0)
Constructor reading centerline file.
double timer()
returns the time in seconds after some point in past
void setup(Mesh *mesh_pt, const Vector< Vector< double > > &coord_vec)
Helper function to setup the output structures.
double & max_search_radius()
Set maximum search radius for locate zeta. This is initialised do DBL_MAX so we brutally search throu...
std::string string(const unsigned &i)
Return the i-th string or "" if the relevant string hasn't been defined.
OomphCommunicator * Comm_pt
Pointer to communicator – allows us to collect data on processor 0 if the mesh is distributed...
Vector< std::pair< FiniteElement *, Vector< double > > > Plot_point
Vector of pairs containing points to finite elements and local coordinates.
void output(std::ostream &outfile)
Output function: output each plot point. NOTE: in a distributed problem, output is only done on proce...
A general mesh class.
Definition: mesh.h:74
LineVisualiser(Mesh *mesh_pt, const Vector< Vector< double > > &coord_vec, const double &max_search_radius=DBL_MAX)
Constructor: Pass pointer to mesh and coordinates of desired plot points: coord_vec[j][i] is the i-th...
An oomph-lib wrapper to the MPI_Comm communicator object. Just contains an MPI_Comm object (which is ...
Definition: communicator.h:57