2D Placement of Streamlines

This chapter describes the CGAL's 2D streamline placement package. Basic definitions and notions are given in Section 41.1. Section 41.2 gives a description of the integration process. Section 41.3 provides a brief description of the algorithm. Section 41.4 presents the implementation of the package, and Section 41.5 details two example placements.

In physics, a *field* is an assignment of a quantity to every
point in space. For example, a gravitational field assigns a
gravitational potential to each point in space.

Vector and direction fields are commonly used for modeling physical
phenomena, where a direction and magnitude, namely a vector is assigned to
each point inside a domain (such as the magnitude and direction of the
force at each point in a magnetic field).

Streamlines are important tools for visualizing flow fields. A
*streamline* is a curve everywhere tangent to the field. In practice,
a streamline is often represented as a polyline (series of points)
iteratively elongated by bidirectional numerical integration started
from a *seed point*, until it comes close to another streamline
(according to a specified distance called *the separating distance*), hits the domain boundary, reaches a critical point or
generates a closed path.

A *valid* placement of streamlines consists of saturating the
domain with a set of tangential streamlines in accordance with a
specified *density*, determined by the *separating distance*
between the streamlines.

A streamline can be considered as the path traced by an imaginary
massless particle dropped into a steady fluid flow described by the
field. The construction of this path consists in the solving an
ordinary differential equation for successive time intervals. In this
way, we obtain a series of points $$*p _{k}, 0<k<n* which allow visualizing
the streamline. The differential equation is defined as follows :

$$*(dp)/(dt) = v(p(t)), p(0) = p _{0}*

where *p(t)* is the
position of the particle at time *t*, *v* is a function which
assigns a vector value at each point in the domain (possibly by
interpolation), and $$*p _{0}* is the initial position. The position after
a given interval $$

$$*p(t + t) = p(t) +
_{t}^{t+t} v(p(t)) dt*

. Several numeric methods have been proposed to solve this equation. In this package, the Euler, and the Second Order Runge Kutta algorithm are implemented.

This algorithm approximates the point computation by this formula

$$*p _{k+1} = p_{k} + hv(p_{k})*

where *h* specifies the
*integration step* (see Figure 41.2.1). The integration
can be done forward (resp. backward) by specifying a positive
(resp. negative) integration step *h*. The streamline is then
constructed by successive integration from a seed point both forward
and backward.

This method introduces an intermediate point *p'_k* between $$*p _{k}* and $$

$$*
*

See [PTVF02] for further details about numerical integration.

The algorithm implemented in this package [MAD05]
consists of placing one streamline at a time by numerical integration
starting farthest away from all previously placed
streamlines.

The input of our algorithm is given by (i) a flow field, (ii) a
*density* specified either globally, by the inverse of the
ideal spacing distance, or locally by a density field, and (iii) a
*saturation* ratio over the desired spacing required to trigger
the seeding of a new streamline.

The input flow field is given by a discrete set of vectors or
directions sampled within a domain, associated with an interpolation
scheme (*e.g.* bilinear interpolation over a regular grid, or
natural neighbor interpolation over an irregular point set to allow
for an evaluation at each point coordinate within the domain).

The *output* is a streamline placement, represented as a list
of streamlines. The core idea of our algorithm consists of placing
one streamline at a time by numerical integration seeded at the
farthest point from all previously placed streamlines.

The streamlines are approximated by polylines, whose points are
inserted to a 2D Delaunay
triangulation (see figure 41). The empty circumscribed
circles of the Delaunay triangles provide us with a good approximation
of the cavities in the domain.

After each streamline integration, all incident triangles whose
circumcircle diameter is larger (within the saturation ratio) than the
desired spacing distance are pushed to a priority queue sorted by the
triangle circumcircle diameter. To start each new streamline
integration, the triangle with largest circumcircle diameter (and
hence the biggest cavity) is popped out of the queue. We first test if
it is still a valid triangle of the triangulation, since it could have
been destroyed by a streamline previously added to the
triangulation. If it is not, we pop another triangle out of the
queue. If it is, we use the center of its circumcircle as seed point
to integrate a new streamline.

Our algorithm terminates when the priority queue is empty. The size of the biggest cavity being monotonically decreasing, our algorithm guarantees the domain saturation.

Streamlines are represented as polylines, and are obtained by
iterative integration from the seed point. A polyline is represented as a range of points. The computation is
processed via a list of Delaunay triangulation vertices.

To implement the triangular grid, the CGAL *Delaunay_triangulation_2* class is used. The priority queue used
to store candidate seed points is taken from the Standard Template
Library [Sil97].

The first example illustrates the generation of a 2D streamline placement from a vector field defined on a regular grid.

#include <iostream> #include <fstream> #include <CGAL/Cartesian.h> #include <CGAL/Filtered_kernel.h> #include <CGAL/Stream_lines_2.h> #include <CGAL/Runge_kutta_integrator_2.h> #include <CGAL/Regular_grid_2.h> #include <CGAL/Triangular_field_2.h> typedef double coord_type; typedef CGAL::Cartesian<coord_type> K1; typedef CGAL::Filtered_kernel<K1> K; typedef CGAL::Regular_grid_2<K> Field; typedef CGAL::Runge_kutta_integrator_2<Field> Runge_kutta_integrator; typedef CGAL::Stream_lines_2<Field, Runge_kutta_integrator> Stl; typedef Stl::Point_iterator_2 Point_iterator; typedef Stl::Stream_line_iterator_2 Stl_iterator; typedef Stl::Point_2 Point_2; typedef Stl::Vector_2 Vector_2; int main() { Runge_kutta_integrator runge_kutta_integrator; /*data.vec.cin is an ASCII file containing the vector field values*/ std::ifstream infile("data/vnoise.vec.cin", std::ios::in); double iXSize, iYSize; unsigned int x_samples, y_samples; iXSize = iYSize = 512; infile >> x_samples; infile >> y_samples; Field regular_grid_2(x_samples, y_samples, iXSize, iYSize); /*fill the grid with the appropreate values*/ for (unsigned int i=0;i<x_samples;i++) for (unsigned int j=0;j<y_samples;j++) { double xval, yval; infile >> xval; infile >> yval; regular_grid_2.set_field(i, j, Vector_2(xval, yval)); } infile.close(); /* the placement of streamlines */ std::cout << "processing...\n"; double dSep = 3.5; double dRat = 1.6; Stl Stream_lines(regular_grid_2, runge_kutta_integrator,dSep,dRat); std::cout << "placement generated\n"; /*writing streamlines to streamlines_on_regular_grid_1.stl */ std::ofstream fw("streamlines_on_regular_grid_1.stl",std::ios::out); fw << Stream_lines.number_of_lines() << "\n"; for(Stl_iterator sit = Stream_lines.begin(); sit != Stream_lines.end(); sit++) { fw << "\n"; for(Point_iterator pit = sit->first; pit != sit->second; pit++){ Point_2 p = *pit; fw << p.x() << " " << p.y() << "\n"; } } fw.close(); }The second example depicts the generation of a streamline placement from a vector field defined on a triangular grid.

#include <iostream> #include <fstream> #include <CGAL/Cartesian.h> #include <CGAL/Filtered_kernel.h> #include <CGAL/Stream_lines_2.h> #include <CGAL/Runge_kutta_integrator_2.h> #include <CGAL/Triangular_field_2.h> typedef double coord_type; typedef CGAL::Cartesian<coord_type> K1; typedef CGAL::Filtered_kernel<K1> K; typedef K::Point_2 Point; typedef K::Vector_2 Vector; typedef CGAL::Triangular_field_2<K> Field; typedef CGAL::Runge_kutta_integrator_2<Field> Runge_kutta_integrator; typedef CGAL::Stream_lines_2<Field, Runge_kutta_integrator> Stl; typedef Stl::Stream_line_iterator_2 stl_iterator; int main() { Runge_kutta_integrator runge_kutta_integrator(1); /*datap.tri.cin and datav.tri.cin are ascii files where are stored the vector velues*/ std::ifstream inp("data/datap.tri.cin"); std::ifstream inv("data/datav.tri.cin"); std::istream_iterator<Point> beginp(inp); std::istream_iterator<Vector> beginv(inv); std::istream_iterator<Point> endp; Field triangular_field(beginp, endp, beginv); /* the placement of streamlines */ std::cout << "processing...\n"; double dSep = 30.0; double dRat = 1.6; Stl Stream_lines(triangular_field, runge_kutta_integrator,dSep,dRat); std::cout << "placement generated\n"; /*writing streamlines to streamlines.stl */ std::cout << "streamlines.stl\n"; std::ofstream fw("streamlines.stl",std::ios::out); Stream_lines.print_stream_lines(fw); }

Next: Reference Manual

CGAL Open Source Project.
Release 3.2.1.
13 July 2006.