Class: FlowField

class pykitPIV.flowfield.FlowFieldSpecs(n_images=1, size=(512, 512), size_buffer=10, time_separation=1, random_seed=None, flowfield_type='random smooth', gaussian_filters=(10, 10), n_gaussian_filter_iter=10, displacement=(2, 2), constant_u_magnitude=(1, 4), constant_v_magnitude=(1, 4), radial_source=True, radial_imposed_source_location=None, radial_sigma=20, radial_epsilon=1e-06, directionally_divergent_imposed_origin=(0, 0), directionally_divergent_directions=['x', 'y', 'both'], directionally_divergent_sources=[True, False], taylor_green_k=1, imposed_origin=(0, 0), apply_SLM=False, integral_time_scale=1, sigma=1, n_stochastic_particles=1000000, n_iterations=100, dtype=<class 'numpy.float64'>)

Configuration object for the FlowField class.

Example:

from pykitPIV import FlowFieldSpecs

# Instantiate an object of FlowFieldSpecs class:
flowfield_spec = FlowFieldSpecs()

# Change one field of flowfield_spec:
flowfield_spec.gaussian_filters = (4, 6)

# You can print the current values of all attributes:
print(flowfield_spec)

class pykitPIV.flowfield.FlowField(n_images, size=(512, 512), size_buffer=10, time_separation=1, dtype=<class 'numpy.float64'>, random_seed=None)

Generates or uploads velocity field(s) for a set of n_images number of PIV image pairs. The velocity field(s) can later be used to advect particles between two consecutive PIV images.

The velocity vector, \(\vec{V} = [u, v]\), is represented as a four-dimensional tensor array of size \((N, 2, H, W)\), where the second index corresponds to \(u\) and \(v\) velocity component, respectively.

The velocity magnitude is computed as \(|\vec{V}| = \sqrt{u^2 + v^2}\) and is also represented as a four-dimensional tensor array of size \((N, 1, H, W)\).

The user also specifies the time separation, \(\Delta t\), between two PIV image pairs with the parameter time_separation. With this, the user can populate the corresponding displacement field, \(d\vec{\mathbf{s}} = [dx, dy] = [u \Delta t, v \Delta t]\), and the displacement field magnitude, \(|d\vec{\mathbf{s}}|\).

Note

Only two-dimensional velocity fields are supported at the moment.

Example:

from pykitPIV import FlowField
import numpy as np

# We are going to generate 10 flow fields for 10 PIV image pairs:
n_images = 10

# Specify size in pixels for each image:
image_size = (128, 512)

# Initialize a flow field object:
flowfield = FlowField(n_images=n_images,
                      size=image_size,
                      size_buffer=10,
                      time_separation=1,
                      dtype=np.float32,
                      random_seed=100)

Parameters:

• n_images – int specifying the number of PIV image pairs, \(N\), to create.

• size – (optional) tuple of two int elements specifying the size of each image in pixels. The first number is image height, \(h\), the second number is image width, \(w\).

• size_buffer – (optional) int specifying the buffer, \(b\), in pixels \([\text{px}]\) to add to the image size in the width and height direction. This number should be approximately equal to the maximum displacement that particles are subject to in order to allow for new particles to arrive into the image area.

• time_separation – (optional) float or int specifying the time separation, \(\Delta t\), between two consecutive PIV images.

• dtype – (optional) numpy.dtype specifying the data type for produced vector fields. To reduce memory, you can switch from the default numpy.float64 to numpy.float32.

• random_seed – (optional) int specifying the random seed for random number generation in numpy. If specified, all operations are reproducible.

Attributes:

• n_images - (read-only) as per user input.

• size - (read-only) as per user input.

• size_buffer - (read-only) as per user input.

• time_separation - (can be re-set) as per user input.

• dtype - (read-only) as per user input.

• random_seed - (read-only) as per user input.

• size_with_buffer - (read-only) tuple specifying the size of each image in pixels with buffer added.

• displacement - (read-only) tuple of two numerical elements specifying the minimum (first element) and maximum (second element) maximum displacement on each image. This attribute becomes available after calling one of the velocity field generators. The assumption is that this displacement value corresponds to time being unity.

• displacement_per_image - (read-only) numpy.ndarray specifying the template for the maximum displacements per each image. Template displacements are random numbers between minimum and maximum displacement. This attribute becomes available after calling one of the velocity field generators. The assumption is that this displacement value corresponds to time being unity.

• velocity_field - (read-only) numpy.ndarray specifying the velocity vector, \(\vec{V} = [u, v]\), per each image. It has size \((N, 2, H+2b, W+2b)\). The second index corresponds to \(u\) and \(v\) velocity component, respectively. This attribute becomes available after calling one of the velocity field generators.

• velocity_field_magnitude - (read-only) numpy.ndarray specifying the velocity field magnitude, \(|\vec{V}| = \sqrt{u^2 + v^2}\), per each image. It has size \((N, 1, H+2b, W+2b)\). This attribute becomes available after calling one of the velocity field generators.

• displacement_field - (read-only) numpy.ndarray specifying the displacement field, \(d\vec{\mathbf{s}} = [dx, dy]\), in the \(x\) and \(y\) direction. It is computed as the velocity component multiplied by time separation and has a unit of \(\text{px}\). It has size \((N, 2, H+2b, W+2b)\). The second index corresponds to \(dx\) and \(dy\) displacement, respectively.

• displacement_field_magnitude - (read-only) numpy.ndarray specifying the displacement field magnitude, \(|d\vec{\mathbf{s}}| = \sqrt{dx^2 + dy^2}\). It has a unit of \(\text{px}\). It has size \((N, 1, H+2b, W+2b)\).

pykitPIV.flowfield.FlowField.print_available_fields(self)

Prints the available velocity fields and points to the function from FlowField that generates them.

Example:

from pykitPIV import FlowField

# Initialize a flow field object:
flowfield = FlowField(n_images=10)

# Print the available velocity fields:
flowfield.print_available_fields()

The above will print:

Velocity fields available in pykitPIV:

- constant
    Use function: generate_constant_velocity_field

- shearing
    Use function: generate_shearing_velocity_field

- random smooth
    Use function: generate_random_velocity_field

- sinusoidal
    Use function: generate_sinusoidal_velocity_field

- checkered
    Use function: generate_checkered_velocity_field

- Chebyshev
    Use function: generate_chebyshev_velocity_field

- spherical harmonics
    Use function: generate_spherical_harmonics_velocity_field

- radial
    Use function: generate_radial_velocity_field

- potential
    Use function: generate_potential_velocity_field

- taylor green vortex
    Use function: generate_taylor_green_vortex_velocity_field

- Langevin
    Use function: generate_langevin_velocity_field

pykitPIV.flowfield.FlowField.generate_constant_velocity_field(self, u_magnitude=(1, 4), v_magnitude=(1, 4))

Generates a constant velocity field.

Example:

from pykitPIV import FlowField

# We are going to generate 10 flow fields for 10 PIV image pairs:
n_images = 10

# Specify size in pixels for each image:
image_size = (128, 128)

# Initialize a flow field object:
flowfield = FlowField(n_images=n_images,
                      size=image_size,
                      size_buffer=10,
                      time_separation=1,
                      random_seed=100)

# Generate a constant velocity field:
flowfield.generate_constant_velocity_field(u_magnitude=(1, 4),
                                           v_magnitude=(1, 4))

# Access the velocity components tensor:
flowfield.velocity_field

# Access the velocity field magnitude:
flowfield.velocity_field_magnitude

Parameters:

• u_magnitude – (optional) tuple of two numerical elements specifying the minimum (first element) and maximum (second element) magnitude of the \(u\) component of velocity.

• v_magnitude – (optional) tuple of two numerical elements specifying the minimum (first element) and maximum (second element) magnitude of the \(v\) component of velocity.

pykitPIV.flowfield.FlowField.generate_shearing_velocity_field(self, u_magnitude=(1, 4), v_magnitude=(1, 4), shear_type='linear')

Generates a shearing flow velocity field.

Example:

from pykitPIV import FlowField

# We are going to generate 10 flow fields for 10 PIV image pairs:
n_images = 10

# Specify size in pixels for each image:
image_size = (128, 128)

# Initialize a flow field object:
flowfield = FlowField(n_images=n_images,
                      size=image_size,
                      size_buffer=10,
                      time_separation=1,
                      random_seed=100)

# Generate a constant velocity field:
flowfield.generate_shearing_velocity_field(u_magnitude=(1, 4),
                                        v_magnitude=(1, 4),
                                        shear_type='linear')

# Access the velocity components tensor:
flowfield.velocity_field

# Access the velocity field magnitude:
flowfield.velocity_field_magnitude

Parameters:

• u_magnitude – (optional) tuple of two numerical elements specifying the minimum (first element) and maximum (second element) magnitude of the \(u\) component of velocity.

• v_magnitude – (optional) tuple of two numerical elements specifying the minimum (first element) and maximum (second element) magnitude of the \(v\) component of velocity.

• shear_type – (optional) str specifying the type of shearing. It can be ‘linear’ or ‘quadratic’.

pykitPIV.flowfield.FlowField.generate_random_velocity_field(self, displacement=(0, 10), gaussian_filters=(10, 30), n_gaussian_filter_iter=6)

Generates random velocity field by smoothing a random initialization of pixels with a series of Gaussian filters.

Example:

from pykitPIV import FlowField

# We are going to generate 10 flow fields for 10 PIV image pairs:
n_images = 10

# Specify size in pixels for each image:
image_size = (128, 512)

# Initialize a flow field object:
flowfield = FlowField(n_images=n_images,
                      size=image_size,
                      size_buffer=10,
                      time_separation=1,
                      random_seed=100)

# Generate random velocity field:
flowfield.generate_random_velocity_field(displacement=(0, 10),
                                         gaussian_filters=(10, 30),
                                         n_gaussian_filter_iter=6)

# Access the velocity components tensor:
flowfield.velocity_field

# Access the velocity field magnitude:
flowfield.velocity_field_magnitude

Parameters:

• displacement – (optional) tuple of two numerical elements specifying the minimum (first element) and maximum (second element) displacement.

• gaussian_filters – (optional) tuple of two numerical elements specifying the minimum (first element) and maximum (second element) Gaussian filter size (bandwidth) for smoothing out the random velocity fields to randomly sample from.

• n_gaussian_filter_iter – (optional) int specifying the number of iterations applying a Gaussian filter to the random velocity field to eventually arrive at a smoothed velocity map. With no iterations, each pixel attains a random velocity component value.

pykitPIV.flowfield.FlowField.generate_sinusoidal_velocity_field(self, amplitudes=(2, 4), wavelengths=(8, 128), components='both')

Generates a sinusoidal velocity field as per Scarano and Riethmuller (2000). Each velocity component is computed as:

\[u(w) = A \sin(2 \pi w / \Lambda)\]

\[v(h) = A \sin(2 \pi h / \Lambda)\]

where \(A\) is the amplitude in pixels \([\text{px}]\) and \(\Lambda\) is the wavelength in pixels \([\text{px}]\).

Example:

from pykitPIV import FlowField

# We are going to generate 10 flow fields for 10 PIV image pairs:
n_images = 10

# Specify size in pixels for each image:
image_size = (128, 512)

# Initialize a flow field object:
flowfield = FlowField(n_images=n_images,
                      size=image_size,
                      size_buffer=10,
                      time_separation=1,
                      random_seed=100)

# Generate sinusoidal velocity field:
flowfield.generate_sinusoidal_velocity_field(amplitudes=(2, 4),
                                             wavelengths=(20, 40),
                                             components='u')

# Access the velocity components tensor:
flowfield.velocity_field

# Access the velocity field magnitude:
flowfield.velocity_field_magnitude

Parameters:

• amplitudes – (optional) tuple of two numerical elements specifying the minimum (first element) and maximum (second element) amplitude, \(A\), in pixels \([\text{px}]\).

• wavelengths – (optional) tuple of two numerical elements specifying the minimum (first element) and maximum (second element) wavelength, \(\Lambda\), in pixels \([\text{px}]\).

• components – (optional) str specifying which velocity components should be generated. It can be one of the following: 'both', 'u', or 'v'. When only one component is generated, the other one is zero.

pykitPIV.flowfield.FlowField.generate_checkered_velocity_field(self, displacement=(1, 10), m=10, n=10, rotation=None)

Generates a checkered velocity field. Each velocity component is computed as:

\[u(h, w), v(h, w) = \sin(m \cdot h) \cdot \cos(n \cdot w)\]

where \(m\) and \(n\) free parameters.

Optionally, the checkered pattern can be rotated with the additional rotation term:

\[u(h, w), v(h, w) = \sin(m \cdot h) \cdot \cos(n \cdot w) \cdot \sin(r \cdot (h + w))\]

were \(r\) is a free parameter.

Example:

from pykitPIV import FlowField

# We are going to generate 10 flow fields for 10 PIV image pairs:
n_images = 10

# Specify size in pixels for each image:
image_size = (128,512)

# Initialize a flow field object:
flowfield = FlowField(n_images=n_images,
                      size=image_size,
                      size_buffer=10,
                      time_separation=1,
                      random_seed=100)

# Generate checkered velocity field:
flowfield.generate_checkered_velocity_field(displacement=(1, 10),
                                            m=10,
                                            n=10,
                                            rotation=10)

# Access the velocity components tensor:
flowfield.velocity_field

# Access the velocity field magnitude:
flowfield.velocity_field_magnitude

Parameters:

• displacement – (optional) tuple of two numerical elements specifying the minimum (first element) and maximum (second element) displacement

• m – (optional) int or float specifying the scaling for the sine function, \(m\).

• n – (optional) int or float specifying the scaling for the cosine function, \(n\).

• rotation – (optional) int or float specifying the scaling for the rotation with the sine function, \(r\). If set to None, the checkered pattern will not be rotated.

pykitPIV.flowfield.FlowField.generate_chebyshev_velocity_field(self, displacement=(1, 10), start=0.3, stop=0.8, order=10)

Generates a velocity field using Chebyshev polynomials of the second kind. Each velocity component is computed as:

\[u(h, w), v(h, w) = U_n(h) U_n(w) \sin(n(h+w))\]

where \(U_n\) is the Chebyshev polynomial of the second kind and of order \(n\).

Example:

from pykitPIV import FlowField

# We are going to generate 10 flow fields for 10 PIV image pairs:
n_images = 10

# Specify size in pixels for each image:
image_size = (128, 512)

# Initialize a flow field object:
flowfield = FlowField(n_images=n_images,
                      size=image_size,
                      size_buffer=10,
                      time_separation=1,
                      random_seed=100)

# Generate Chebyshev velocity field:
flowfield.generate_chebyshev_velocity_field(displacement=(0, 10),
                                            start=0.3,
                                            stop=0.8,
                                            order=10)

# Access the velocity components tensor:
flowfield.velocity_field

# Access the velocity field magnitude:
flowfield.velocity_field_magnitude

Parameters:

• displacement – (optional) tuple of two numerical elements specifying the minimum (first element) and maximum (second element) displacement.

• start – (optional) int or float specifying the start value for generating the coordinate grid.

• stop – (optional) int or float specifying the end value for generating the coordinate grid.

• order – (optional) int specifying the order of the Chebyshev polynomial.

pykitPIV.flowfield.FlowField.generate_spherical_harmonics_velocity_field(self, displacement=(1, 10), start=0.3, stop=0.8, order=1, degree=1)

Generates a velocity field using spherical harmonics. Each velocity component is computed as:

\[u(h, w), v(h, w) = Y_n^k \sin(n(h+w))\]

where \(Y_n^k\) is the spherical harmonics function of order \(n\) and degree \(k\).

Note

This function uses scipy.special.sph_harm_y.

Example:

from pykitPIV import FlowField

# We are going to generate 10 flow fields for 10 PIV image pairs:
n_images = 10

# Specify size in pixels for each image:
image_size = (128, 512)

# Initialize a flow field object:
flowfield = FlowField(n_images=n_images,
                      size=image_size,
                      size_buffer=10,
                      time_separation=1,
                      random_seed=100)

# Generate spherical harmonics velocity field:
flowfield.generate_spherical_harmonics_velocity_field(displacement=(0, 10),
                                                      start=0.3,
                                                      stop=0.8,
                                                      order=1,
                                                      degree=1)

# Access the velocity components tensor:
flowfield.velocity_field

# Access the velocity field magnitude:
flowfield.velocity_field_magnitude

Parameters:

• displacement – (optional) tuple of two numerical elements specifying the minimum (first element) and maximum (second element) displacement.

• start – (optional) int or float specifying the start value for generating the azimuthal/polar coordinate grid.

• stop – (optional) int or float specifying the end value for generating the azimuthal/polar coordinate grid.

• order – (optional) int specifying the order, \(n\), of the spherical harmonics.

• degree – (optional) int specifying the degree, \(k\), of the spherical harmonics.

pykitPIV.flowfield.FlowField.generate_directionally_divergent_velocity_field(self, direction='x', source=True, imposed_origin=(0, 0), displacement=(1, 4))

Generates a velocity field that has non-zero derivative in only the selected direction(s). The non-zero components of the velocity field have the form:

\[u(h, w) = (w - w_o)^2 (h - h_o)\]

\[v(h, w) = (h - h_o)^2 (w - w_o)\]

where \((h_o, w_o)\) is the origin location on the PIV image.

This flow can control the direction that contributes to the divergence.

Example:

from pykitPIV import FlowField

# We are going to generate 10 flow fields for 10 PIV image pairs:
n_images = 10

# Specify size in pixels for each image:
image_size = (128, 512)

# Initialize a flow field object:
flowfield = FlowField(n_images=n_images,
                      size=image_size,
                      size_buffer=10,
                      time_separation=1,
                      random_seed=100)

# Generate a potential velocity field:
flowfield.generate_directionally_divergent_velocity_field(direction='x',
                                                          source=True,
                                                          imposed_origin=None,
                                                          displacement=(2, 2))

# Access the velocity components tensor:
flowfield.velocity_field

# Access the velocity field magnitude:
flowfield.velocity_field_magnitude

Parameters:

• direction – (optional) str specifying the direction where there is a non-zero velocity field. It can be 'x' or 'y' or 'both'.

• source – (optional) bool specifying generation of the source (True) or sink (False).

• imposed_origin – (optional) tuple specifying the user-imposed origin location, \((h_o, w_o)\). Note that you need to account for the buffer size when specifying the values \((h_o, w_o)\). If set to None, the origin location will be randomized in each image.

• displacement – (optional) tuple of two numerical elements specifying the minimum (first element) and maximum (second element) displacement.

pykitPIV.flowfield.FlowField.generate_radial_velocity_field(self, source=True, displacement=(1, 4), imposed_source_location=None, sigma=20.0, epsilon=1e-06)

Generates a radial velocity field that has a single source or sink. The velocity components are computed as:

\[u(h, w) = \Big( 1 - \exp{ \big( \frac{-r^2}{2 \sigma_1^2} \big)} \Big) \frac{(w - w_s)}{r + \varepsilon}\]

\[v(h, w) = \Big( 1 - \exp{ \big( \frac{-r^2}{2 \sigma_2^2} \big)} \Big) \frac{(h - h_s)}{r + \varepsilon}\]

where \(\sigma\) is the size of source/sink, \(r\) is the radial distance from the source/sink computed as:

\[r = \sqrt{((w - w_s)^2 + (h - h_s)^2)}\]

and where \(w_s\) and \(h_s\) are the location of the source/sink.

\(\varepsilon\) is a small value added to avoid dividing by zero.

Example:

from pykitPIV import FlowField

# We are going to generate 10 flow fields for 10 PIV image pairs:
n_images = 10

# Specify size in pixels for each image:
image_size = (200, 200)

# Initialize a flow field object:
flowfield = FlowField(n_images=n_images,
                      size=image_size,
                      size_buffer=10,
                      time_separation=1,
                      random_seed=100)

# Generate radial velocity field:
flowfield.generate_radial_velocity_field(source=True,
                                         displacement=(1, 4),
                                         imposed_source_location=(50, 50),
                                         sigma=5.0,
                                         epsilon=1e-6)

# Access the velocity components tensor:
flowfield.velocity_field

# Access the velocity field magnitude:
flowfield.velocity_field_magnitude

Parameters:

• source – (optional) bool specifying generation of the source (True) or sink (False).

• displacement – (optional) tuple of two numerical elements specifying the minimum (first element) and maximum (second element) displacement.

• imposed_source_location – (optional) tuple specifying the user-imposed source location, \((h_s, w_s)\). Note that you need to account for the buffer size when specifying the values \((h_s, w_s)\).

• sigma – (optional) `int or float or tuple of two `int or float specifying the size of the source/sink in each direction, \(\sigma_1\) and \(\sigma_2\).

• epsilon – (optional) float specifying the \(\varepsilon\).

pykitPIV.flowfield.FlowField.generate_potential_velocity_field(self, imposed_origin=(0, 0), displacement=(1, 4))

Generates a potential velocity field of the form:

\[u(h, w) = 2 (w - w_o)\]

\[v(h, w) = - 2 (h - h_o)\]

where \((h_o, w_o)\) is the origin location on the PIV image.

This potential flow is irrotational and divergence-free and can therefore be used for generating BOS images.

Example:

from pykitPIV import FlowField

# We are going to generate 10 flow fields for 10 PIV image pairs:
n_images = 10

# Specify size in pixels for each image:
image_size = (128, 512)

# Initialize a flow field object:
flowfield = FlowField(n_images=n_images,
                      size=image_size,
                      size_buffer=10,
                      time_separation=1,
                      random_seed=100)

# Generate a potential velocity field:
flowfield.generate_potential_velocity_field(imposed_origin=None,
                                            displacement=(2, 2))

# Access the velocity components tensor:
flowfield.velocity_field

# Access the velocity field magnitude:
flowfield.velocity_field_magnitude

Parameters:

• imposed_origin – (optional) tuple specifying the user-imposed origin location, \((h_o, w_o)\). Note that you need to account for the buffer size when specifying the values \((h_o, w_o)\). If set to None, the origin location will be randomized in each image.

• displacement – (optional) tuple of two numerical elements specifying the minimum (first element) and maximum (second element) displacement.

pykitPIV.flowfield.FlowField.generate_taylor_green_vortex_velocity_field(self, k=(0.01, 0.05), imposed_origin=(0, 0), displacement=(1, 4))

Generates a Taylor-Green vortex of the form:

\[u(h, w) = \sin(k (w - w_o)) \cos(k (h - h_o))\]

\[v(h, w) = - \cos(k (w - w_o)) \sin(k (h - h_o))\]

where \((h_o, w_o)\) is the origin location on the PIV image and \(k\) is the wavelength that can be where \((h_o, w_o)\) is the origin location on the PIV image and \(k\) is the wavelength that can be used to adjust the sizes of vortices.

This flow field divergence-free.

Example:

from pykitPIV import FlowField

# We are going to generate 10 flow fields for 10 PIV image pairs:
n_images = 10

# Specify size in pixels for each image:
image_size = (256, 256)

# Initialize a flow field object:
flowfield = FlowField(n_images=n_images,
                      size=image_size,
                      size_buffer=10,
                      time_separation=1,
                      random_seed=100)

# Generate a Taylor-Green vortex velocity field:
flowfield.generate_taylor_green_vortex_velocity_field(k=(0.01, 0.05),
                                                      imposed_origin=None,
                                                      displacement=(2, 2))

# Access the velocity components tensor:
flowfield.velocity_field

# Access the velocity field magnitude:
flowfield.velocity_field_magnitude

Parameters:

• k – tuple of two numerical elements specifying the minimum (first element) and maximum (second element) wavelength, \(k\).

• imposed_origin – (optional) tuple specifying the user-imposed origin location, \((h_o, w_o)\). Note that you need to account for the buffer size when specifying the values \((h_o, w_o)\). If set to None, the origin location will be randomized in each image.

• displacement – (optional) tuple of two numerical elements specifying the minimum (first element) and maximum (second element) displacement.

pykitPIV.flowfield.FlowField.generate_langevin_velocity_field(self, mean_field, integral_time_scale=1, sigma=1, n_stochastic_particles=10000, n_iterations=100, verbose=False)

Generates a velocity field using the simplified Langevin model (SLM) by solving the following stochastic partial differential equation:

\[dU^*(t) = - \frac{1}{T_L} U^*(t) dt + \sqrt{\frac{1}{T_L} 2 \sigma^2 } d W(t)\]

where \(U^*(t)\) is one component of particle’s velocity, \(T_L\) is the Lagrangian integral time scale, \(\sigma\) is the turbulent kinetic energy, and \(W(t)\) is the stochastic variable.

In essence, the SLM model simulates stationary isotropic turbulence in the Lagrangian sense, where stochastic particles are allowed to drift and carry fluid velocity to a different location.

In practice, we solve the following finite-difference set of equations for both components of velocity:

\[ \begin{align}\begin{aligned}u^*(t + \Delta t) = u^*(t) - (u^*(t) - u_m^*) \frac{\Delta t}{T_L} + \sqrt{\frac{1}{T_L} 2 \sigma^2 \Delta t} \xi(t)\\x(t + \Delta t) = x(t) + u^*(t) \Delta t\end{aligned}\end{align} \]

\[ \begin{align}\begin{aligned}v^*(t + \Delta t) = v^*(t) - (v^*(t) - v_m^*) \frac{\Delta t}{T_L} + \sqrt{\frac{1}{T_L} 2 \sigma^2 \Delta t} \xi(t)\\y(t + \Delta t) = y(t) + v^*(t) \Delta t\end{aligned}\end{align} \]

where \(u_m^*\) and \(v_m^*\) are the components of the mean velocity, \(\xi(t)\) is a random variable \(\xi(t) \in \mathcal{N}(0,1)\), and \(x(t)\) and \(y(t)\) are coordinates of the stochastic particles.

The above equations define the update to velocity components and particle positions over \(n\) iterations (set by the parameter n_iterations). \(\Delta t\) is then computed as:

\[\Delta t = \frac{T_L}{n}\]

Next, a per-pixel velocity is found as an average over the ensemble of stochastic particles present within that pixel. In the case where particle drift makes any pixel empty, velocity component is interpolated.

For a thorough discussion of SLM see section 12.3 of S.B. Pope - Turbulent Flows (1995).

Note

Note that the stochastic particles for solving the SLM are in practice generated using the Particle class but they are independent of the seeded particles used for PIV.

Note

If the Lagrangian integral time scale is large, you may want to create a large image buffer to account for the stochastic particles moving a considerable distance compared to image size.

Example:

from pykitPIV import FlowField

# We are going to generate 10 flow fields for 10 PIV image pairs:
n_images = 10

# Specify size in pixels for each image:
image_size = (128, 128)

# Initialize a mean flow field object:
mean_flowfield = FlowField(n_images=n_images,
                           size=image_size,
                           size_buffer=10,
                           time_separation=1,
                           random_seed=100)

# Generate sinusoidal velocity field that will serve as the mean velocity field:
mean_flowfield.generate_sinusoidal_velocity_field(amplitudes=(2, 4),
                                                  wavelengths=(20, 40),
                                                  components='both')

# Initialize a flow field object:
flowfield = FlowField(n_images=n_images,
                      size=image_size,
                      size_buffer=10,
                      time_separation=1,
                      random_seed=100)

# Solve the SLM for the mean velocity fields:
flowfield.generate_langevin_velocity_field(mean_field=mean_flowfield.velocity_field,
                                           integral_time_scale=1,
                                           sigma=1,
                                           n_stochastic_particles=10000,
                                           n_iterations=100,
                                           verbose=True)

# Access the velocity components tensor:
flowfield.velocity_field

# Access the velocity field magnitude:
flowfield.velocity_field_magnitude

Parameters:

• mean_field – numpy.ndarray specifying the mean velocity field with components \(u_m^*\) and \(v_m^*\), respectively. It should be of size \((N, 2, H+2 \cdot b, W+2 \cdot b)\). The SLM model essentially provides a stochastic drift over that user-specified mean velocity, preserving the Lagrangian integral time scale and the turbulent kinetic energy. The mean velocity can be synthetically generated using one of the velocity field generators provided in this class.

• integral_time_scale – (optional) int or float specifying the Lagrangian integral time scale, \(T_L\).

• sigma – (optional) int or float specifying the turbulent kinetic energy.

• n_stochastic_particles – (optional) int specifying the number of stochastic particles for which the SLM will be solved.

• n_iterations – (optional) int specifying the number of iterations, \(n\), during which the finite-difference equations are updated.

• verbose – (optional) bool for printing verbose details.

pykitPIV.flowfield.FlowField.compute_displacement_field(self)

Computes the displacement field, \(d \vec{\mathbf{s}}\), and the displacement field magnitude, \(|d \vec{\mathbf{s}}|\), based on the current value of the time separation, \(\Delta t\). Calling this function populates the class attributes displacement_field and displacement_field_magnitude.

Example:

from pykitPIV import FlowField

# We are going to generate 10 flow fields for 10 PIV image pairs:
n_images = 10

# Specify size in pixels for each image:
image_size = (128, 512)

# Initialize a flow field object:
flowfield = FlowField(n_images=n_images,
                      size=image_size,
                      size_buffer=10,
                      time_separation=1,
                      random_seed=100)

# Generate random velocity field:
flowfield.generate_random_velocity_field(displacement=(0, 10),
                                         gaussian_filters=(10, 30),
                                         n_gaussian_filter_iter=6)

# Compute the displacement field:
flowfield.compute_displacement_field()

pykitPIV.flowfield.FlowField.upload_velocity_field(self, velocity_field)

Uploads a custom velocity field, e.g., generated from synthetic turbulence or from turbulence databases.

An example synthetic turbulence generator implemented in Python can be found here: Saad et al. - TurboGenPY. You can also access its web-based version.

Another popular source of turbulent velocity fields is the Johns Hopkins Turbulence Database (JHTD).

Example:

from pykitPIV import FlowField

# We are going to use 10 flow fields for 10 PIV image pairs:
n_images = 10

# Specify size in pixels for each image:
image_size = (128, 512)

# Initialize a flow field object:
flowfield = FlowField(n_images=n_images,
                      size=image_size,
                      size_buffer=10,
                      time_separation=1,
                      random_seed=100)

# Importing an external velocity field file:
# ...
# ...

# Prepare the velocity_field variable that has shape (10, 2, 148, 532):
# velocity_field = ...

# Upload velocity field:
flowfield.upload_velocity_field(velocity_field)

Parameters:

velocity_field – numpy.ndarray specifying the velocity components. It should be of size \((N, 2, H, W)\), where \(N\) is the number of PIV image pairs, \(2\) refers to each velocity component \(u\) and \(v\) respectively, \(H\) is the height and \(W\) is the width of each PIV image. It can also be of size \((1, 2, H, W)\), in which case the same velocity field will be applied to all PIV image pairs.

Flow field utilities

pykitPIV.flowfield.compute_divergence(vector_field, edge_order=1)

Computes the divergence of the specified vector field.

If the vector field is the velocity field, \(\vec{V} = [u, v]\):

\[\nabla \cdot \vec{V} = \frac{\partial u}{\partial x} + \frac{\partial v}{\partial y}\]

Note

Derivatives are computed using numpy.gradient.

Example:

from pykitPIV import FlowField, compute_divergence

# Initialize a flow field object:
flowfield = FlowField(10,
                      size=(200, 200),
                      size_buffer=0,
                      time_separation=1,
                      random_seed=100)

# Generate random velocity field:
flowfield.generate_random_velocity_field(gaussian_filters=(10, 11),
                                         n_gaussian_filter_iter=20,
                                         displacement=(1, 2))

# Extract the velocity field components:
velocity_field = flowfield.velocity_field

# Compute the divergence of the velocity field:
divergence = compute_divergence(vector_field=velocity_field,
                                edge_order=1)

Parameters:

• vector_field – numpy.ndarray specifying the vector field components. It should be of size \((N, 2, H, W)\), where \(N\) is the number of PIV image pairs, \(2\) refers to each vector field component and \(H\) is the height and \(W\) is the width of each PIV image. For example, it can be the velocity field with components \(u\) and \(v\), or the displacement field with components \(dx\) and \(dy\).

• edge_order –

  (optional) int specifying the order for the gradient computation at image boundaries as per numpy.gradient.

Returns:

• divergence - divergence of the velocity field, \(\nabla \cdot \vec{V}\). It has size \((N, H, W)\).

pykitPIV.flowfield.compute_vorticity(vector_field, edge_order=1)

Computes the vorticity of the specified vector field.

If the vector field is the velocity field, \(\vec{V} = [u, v]\):

\[\nabla \times \vec{V} = \frac{\partial v }{\partial x} - \frac{\partial u}{\partial y}\]

Note

Derivatives are computed using numpy.gradient.

Example:

from pykitPIV import FlowField, compute_vorticity

# Initialize a flow field object:
flowfield = FlowField(10,
                      size=(200, 200),
                      size_buffer=0,
                      time_separation=1,
                      random_seed=100)

# Generate random velocity field:
flowfield.generate_random_velocity_field(gaussian_filters=(10, 11),
                                         n_gaussian_filter_iter=20,
                                         displacement=(1, 2))

# Extract the velocity field components:
velocity_field = flowfield.velocity_field

# Compute the vorticity of the velocity field:
vorticity = compute_vorticity(vector_field=velocity_field,
                              edge_order=1)

Parameters:

• vector_field – numpy.ndarray specifying the vector field components. It should be of size \((N, 2, H, W)\), where \(N\) is the number of PIV image pairs, \(2\) refers to each vector field component and \(H\) is the height and \(W\) is the width of each PIV image. For example, it can be the velocity field with components \(u\) and \(v\), or the displacement field with components \(dx\) and \(dy\).

• edge_order –

  (optional) int specifying the order for the gradient computation at image boundaries as per numpy.gradient.

Returns:

• vorticity - vorticity of the velocity field, \(\nabla \times \vec{V}\). It has size \((N, H, W)\).

pykitPIV.flowfield.compute_rate_of_strain_tensor(vector_field, edge_order=1)

Computes the components of the rate-of-strain tensor of the specified vector field.

If the vector field is the velocity field, \(\vec{V} = [u, v]\), the velocity gradient can be decomposed into the rate-of-strain tensor and the rate-of-rotation (vorticity) tensor:

\[\nabla \vec{V} = \mathbf{S} + \pmb{\Omega}\]

where \(\mathbf{S} = \frac{1}{2} [ (\nabla \vec{V}) + (\nabla \vec{V})^{\top}]\) is the rate-of-strain tensor and \(\pmb{\Omega} = \frac{1}{2} [ (\nabla \vec{V}) - (\nabla \vec{V})^{\top}]\) is the rate-of-rotation (vorticity) tensor.

The components of the rate-of-strain tensor are computed as

\[\begin{split}\mathbf{S} = \begin{bmatrix} S_{11} & S_{12} \\ S_{21} & S_{22} \end{bmatrix} = \begin{bmatrix} \frac{\partial u}{\partial x} & \frac{1}{2} \left( \frac{\partial v}{\partial x} + \frac{\partial u}{\partial y} \right) \\ \frac{1}{2} \left( \frac{\partial u}{\partial y} + \frac{\partial v}{\partial x} \right) & \frac{\partial v}{\partial y} \end{bmatrix}\end{split}\]

Note

Derivatives are computed using numpy.gradient.

Example:

from pykitPIV import FlowField, compute_rate_of_strain_tensor

# Initialize a flow field object:
flowfield = FlowField(10,
                      size=(200, 200),
                      size_buffer=0,
                      time_separation=1,
                      random_seed=100)

# Generate random velocity field:
flowfield.generate_random_velocity_field(gaussian_filters=(10, 11),
                                         n_gaussian_filter_iter=20,
                                         displacement=(1, 2))

# Extract the velocity field components:
velocity_field = flowfield.velocity_field

# Compute the components of the rate-of-strain tensor:
S_11, S_22, S_12 = compute_rate_of_strain_tensor(vector_field=velocity_field,
                                                 edge_order=1)

Parameters:

• vector_field – numpy.ndarray specifying the vector field components. It should be of size \((N, 2, H, W)\), where \(N\) is the number of PIV image pairs, \(2\) refers to each vector field component and \(H\) is the height and \(W\) is the width of each PIV image. For example, it can be the velocity field with components \(u\) and \(v\), or the displacement field with components \(dx\) and \(dy\).

• edge_order –

  (optional) int specifying the order for the gradient computation at image boundaries as per numpy.gradient.

Returns:

• S_11 - the \(S_{11}\) component of the rate-of-strain tensor. It has size \((N, H, W)\).

• S_22 - the \(S_{22}\) component of the rate-of-strain tensor. It has size \((N, H, W)\).

• S_12 - the \(S_{12}\) (equal to \(S_{21}\)) component of the rate-of-strain tensor. It has size \((N, H, W)\).

pykitPIV.flowfield.compute_rate_of_rotation_tensor(vector_field, edge_order=1)

Computes the components of the rate-of-rotation (vorticity) tensor of the specified vector field.

If the vector field is the velocity field, \(\vec{V} = [u, v]\), the velocity gradient can be decomposed into the rate-of-strain tensor and the rate-of-rotation (vorticity) tensor:

\[\nabla \vec{V} = \mathbf{S} + \pmb{\Omega}\]

where \(\mathbf{S} = \frac{1}{2} [ (\nabla \vec{V}) + (\nabla \vec{V})^{\top}]\) is the rate-of-strain tensor and \(\pmb{\Omega} = \frac{1}{2} [ (\nabla \vec{V}) - (\nabla \vec{V})^{\top}]\) is the rate-of-rotation (vorticity) tensor.

The components of the rate-of-rotation (vorticity) tensor are computed as

\[\begin{split}\pmb{\Omega} = \begin{bmatrix} \omega_{11} & \omega_{12} \\ \omega_{21} & \omega_{22} \end{bmatrix} = \begin{bmatrix} 0 & \frac{1}{2} \left( \frac{\partial v}{\partial x} - \frac{\partial u}{\partial y} \right) \\ \frac{1}{2} \left( \frac{\partial u}{\partial y} - \frac{\partial v}{\partial x} \right) & 0 \end{bmatrix}\end{split}\]

Note

Derivatives are computed using numpy.gradient.

Example:

from pykitPIV import FlowField, compute_rate_of_rotation_tensor

# Initialize a flow field object:
flowfield = FlowField(10,
                      size=(200, 200),
                      size_buffer=0,
                      time_separation=1,
                      random_seed=100)

# Generate random velocity field:
flowfield.generate_random_velocity_field(gaussian_filters=(10, 11),
                                         n_gaussian_filter_iter=20,
                                         displacement=(1, 2))

# Extract the velocity field components:
velocity_field = flowfield.velocity_field

# Compute the components of the rate-of-rotation (vorticity) tensor:
omega_12, omega_21 = compute_rate_of_rotation_tensor(vector_field=velocity_field,
                                                     edge_order=1)

Parameters:

• vector_field – numpy.ndarray specifying the vector field components. It should be of size \((N, 2, H, W)\), where \(N\) is the number of PIV image pairs, \(2\) refers to each vector field component and \(H\) is the height and \(W\) is the width of each PIV image. For example, it can be the velocity field with components \(u\) and \(v\), or the displacement field with components \(dx\) and \(dy\).

• edge_order –

  (optional) int specifying the order for the gradient computation at image boundaries as per numpy.gradient.

Returns:

• omega_12 - the \(\omega_{12}\) component of the rate-of-rotation (vorticity) tensor. It has size \((N, H, W)\).

• omega_21 - the \(\omega_{21}\) component of the rate-of-rotation (vorticity) tensor. It has size \((N, H, W)\).

pykitPIV.flowfield.compute_q_criterion(vector_field, edge_order=1)

Computes the Q-criterion for the specified vector field, \(\vec{V} = [u, v]\).

The Q-criterion comes from decomposing the velocity gradient into the rate-of-strain tensor and the rate-of-rotation (vorticity) tensor:

\[\nabla \vec{V} = \mathbf{S} + \pmb{\Omega}\]

where \(\mathbf{S} = \frac{1}{2} [ (\nabla \vec{V}) + (\nabla \vec{V})^{\top}]\) is the rate-of-strain tensor and \(\pmb{\Omega} = \frac{1}{2} [ (\nabla \vec{V}) - (\nabla \vec{V})^{\top}]\) is the rate-of-rotation (vorticity) tensor.

The Q-criterion is defined as regions where the vorticity tensor magnitude dominates over the strain rate tensor magnitude:

\[Q = \frac{1}{2} \big( \| \pmb{\Omega} \|^2 - \| \mathbf{S} \|^2 \big) > 0\]

where \(\| \bullet \|\) is the Euclidean (or Frobenius) norm.

In 2D, the rate-of-rotation (vorticity) tensor is

\[\begin{split}\pmb{\Omega} = \begin{bmatrix} 0 & \frac{1}{2} \left( \frac{\partial v}{\partial x} - \frac{\partial u}{\partial y} \right) \\ \frac{1}{2} \left( \frac{\partial u}{\partial y} - \frac{\partial v}{\partial x} \right) & 0 \end{bmatrix}\end{split}\]

and the rate-of-strain tensor is

\[\begin{split}\mathbf{S} = \begin{bmatrix} \frac{\partial u}{\partial x} & \frac{1}{2} \left( \frac{\partial v}{\partial x} + \frac{\partial u}{\partial y} \right) \\ \frac{1}{2} \left( \frac{\partial u}{\partial y} + \frac{\partial v}{\partial x} \right) & \frac{\partial v}{\partial y} \end{bmatrix}\end{split}\]

Hence, in 2D, we can write the Q-criterion as:

\[Q = \frac{1}{2} \left( \frac{1}{2} \left( \frac{\partial v}{\partial x} - \frac{\partial u}{\partial y} \right)^2 - \left( \frac{\partial u}{\partial x} \right)^2 - \left( \frac{\partial v}{\partial y} \right)^2 - \frac{1}{2} \left( \frac{\partial u}{\partial y} + \frac{\partial v}{\partial x} \right)^2 \right)\]

See Jeong & Hussain and Haller for more info on the Q-criterion.

Note

Derivatives are computed using numpy.gradient.

Example:

from pykitPIV import FlowField, compute_q_criterion

# Initialize a flow field object:
flowfield = FlowField(10,
                      size=(200, 200),
                      size_buffer=0,
                      time_separation=1,
                      random_seed=100)

# Generate random velocity field:
flowfield.generate_random_velocity_field(gaussian_filters=(10, 11),
                                         n_gaussian_filter_iter=20,
                                         displacement=(1, 2))

# Extract the velocity field components:
velocity_field = flowfield.velocity_field

# Compute the Q-criterion for the velocity field:
q_criterion = compute_q_criterion(vector_field=velocity_field,
                                  edge_order=1)

Parameters:

• vector_field – numpy.ndarray specifying the vector field components. It should be of size \((N, 2, H, W)\), where \(N\) is the number of PIV image pairs, \(2\) refers to each vector field component and \(H\) is the height and \(W\) is the width of each PIV image. For example, it can be the velocity field with components \(u\) and \(v\), or the displacement field with components \(dx\) and \(dy\).

• edge_order –

  (optional) int specifying the order for the gradient computation at image boundaries as per numpy.gradient.

Returns:

• q_criterion - Q-criterion of the velocity field, \(Q\). It has size \((N, H, W)\).