Interpolation

Description

This treatment interpolates a multiblock source CFD domain on another target CFD domain.

The interpolation method is the Inverse Distance Weighting.

The complexity is N*logN thanks to the use of a fast spatial search structure (kd-tree).

The set of N known data points: \([ \ldots, (x_i, y_i), \ldots]\)

Interpolated value u at a given point x based on samples \(\displaystyle u_{i}=u(x_{i})\):

\(\displaystyle u(x) = \frac{\sum_1^N \frac{1}{d(x,x_i)^p} u_i} {\sum_1^N \frac{1}{d(x,x_i)^p}}\) if \(\displaystyle d(x,x_i) \neq 0 \forall i\)

\(\displaystyle u(x) = u_i\) if \(\displaystyle \exists i, d(x,x_i) = 0\)

x is the arbitrary interpolated point. \(\displaystyle x_i\) is the known interpolating point. \(\displaystyle d(x,x_i)\) is the distance from the known point \(\displaystyle x_i\) to the unknown point x. N is the total number of known points used in interpolation. p is the power parameter.

Warning

dependency on scipy.spatial

Parameters

  • source: Base

    The source mesh with the variables you want to interpolate.

  • target: Base

    The mesh to interpolate on. The base is modified in-place at the output.

  • coordinates: list(str), default= None

    The variable names that define the set of coordinates used for interpolation. The KDTree extracts closest points based on these variables. If None, use the coordinates from the source or target base.

  • variables: list(str), default= None

    List of variables to apply for the interpolation. If None, all variables will be interpolated, except the coordinates. Using this parameter allow to call the interpolation treatment on bases with different variable locations.

  • tolerance: float, default= 1e-10

    The value that tells if the distance from the closest point is smaller than this tolerance, the interpolated value is equal to the value of that closest point.

  • duplicated_tolerance: int, default= None

    Number of decimal places to round to for duplicated points detection. If decimals is negative, it specifies the number of positions to the left of the decimal point. If None, the detection will be exact.

  • invdist_power: int, default= 1

    The power of the inverse distance, this acts like a norm (1, 2, or infinity etc.) and can be used to smooth/refined the interpolation.

  • nb_points: int, default= None

    The number of closest points used for the interpolation. If not provided by the user, it will automatically be set depending on the type of cells in the mesh.

  • with_boundaries: bool, default= False

    Whether or not to use data from the boundaries.

Preconditions

Postconditions

The target base is modified in-place.

Example

import antares
myt = antares.Treatment('interpolation')
myt['source'] = source_base
myt['target'] = target_base
myt.execute()

Main functions

class antares.treatment.TreatmentInterpolation.TreatmentInterpolation

Class for Interpolation Treatment.

execute()

Interpolate data from a base onto data from another base.

Returns:

the target base containing the results

Return type:

Base

Example

"""
This example illustrates the interpolation of a
cloud of points against another cloud of points
"""

from copy import copy
import os
if not os.path.isdir('OUTPUT'):
    os.makedirs('OUTPUT')


from antares import Reader, Treatment, Writer

# ------------------
# Reading the files
# ------------------
reader = Reader('bin_tp')
reader['filename'] = os.path.join('..', 'data', 'ROTOR37', 'GENERIC', 'flow_<zone>_<instant>.dat')
source_base = reader.read()

print(source_base)

# --------------------
# Build a target base
# --------------------
# (same mesh but at cell centers)
tmp_base = copy(source_base[:, :, ('x', 'y', 'z')])
tmp_base.node_to_cell()
target_base = tmp_base.get_location('cell')

# --------------
# Interpolation
# --------------
treatment = Treatment('interpolation')
treatment['source'] = source_base
treatment['target'] = target_base
result = treatment.execute()

# Note that the result contains the same instants as the source and
# each instant contains the variable interpolated and the coordinates
# used for interpolation
print(result[0])
# >>>  Zone object
#      - instants : ['ite0', 'ite1', 'ite2']
print(result[0][0])
# >>>  Instant object
#      - variables : [('x', 'cell'), ('y', 'cell'), ('z', 'cell'), ('ro', 'cell'), ('rovx', 'cell'), ('rovy', 'cell'), ('rovz', 'cell'), ('roE', 'cell')]
#      - shape : None

# -------------------
# Writing the result
# -------------------
writer = Writer('bin_tp')
writer['base'] = result
writer['filename'] = os.path.join('OUTPUT', 'ex_interpolation.plt')
writer.dump()