pyforce.tools package

Submodules

pyforce.tools.backends module

class pyforce.tools.backends.IntegralCalculator(grid: pyvista.UnstructuredGrid, gdim=3)[source]

Bases: object

Class to compute integrals and norms on a PyVista UnstructuredGrid. This class supports both 2D and 3D grids and computes cell sizes accordingly. It provides methods to compute integrals, averages, \(L^1\) norms, \(L^2\) inner products, and \(L^2\) norms.

Parameters:
  • grid (pv.UnstructuredGrid) – The PyVista UnstructuredGrid on which the operations will be performed.

  • gdim (int, optional (Default = 3)) – The geometric dimension of the grid. It can be either 2 or 3. If set to 2, the area of the cells will be computed; if set to 3, the volume of the cells will be computed.

grid

The PyVista UnstructuredGrid on which the operations will be performed.

Type:

pv.UnstructuredGrid

gdim

The geometric dimension of the grid (2 or 3).

Type:

int

cell_sizes

The sizes of the cells in the grid, computed as areas for 2D grids or volumes for 3D grids.

Type:

np.ndarray

n_points

The number of points in the grid.

Type:

int

n_cells

The number of cells in the grid.

Type:

int

L1_norm(u)[source]

Computes the \(L^1\) norm of a function u over the domain

\[\|u\|_{L^1}=\int_\Omega |u| \,d\Omega\]
Parameters:

u (np.ndarray) – Function belonging to the grid

Returns:

value\(L^1\) norm of the function

Return type:

float

L2_inner_product(u, v)[source]

Compute the L2 inner product of two functions \(\left( u, v \right)\) over the domain

\[\left( u, v \right)_{L^2} = \int_\Omega u \cdot v \,d\Omega\]
Parameters:
  • u (np.ndarray) – First function belonging to the grid.

  • v (np.ndarray) – Second function belonging to the grid.

Returns:

inner_product – The L2 inner product of the two functions.

Return type:

float

L2_norm(u)[source]

Computes the \(L^2\) norm of a function u over the domain

\[\|u\|_{L^2}=\sqrt{\left( u, u \right)_{L^2}}\]
Parameters:

u (np.ndarray) – Function belonging to the grid

Returns:

norm\(L^2\) norm of the function

Return type:

float

average(u)[source]

Computes the integral average of a given scalar function u over the domain

\[\langle u \rangle = \frac{1}{|\Omega|}\int_\Omega u \,d\Omega\]
Parameters:

u (np.ndarray) – Function belonging to the grid

Returns:

ave_value – Average over the domain

Return type:

float

check_input(u)[source]

Check if the input array u is valid for the grid. It should have the same number of elements as the number of points or cells in the grid.

Parameters:

u (np.ndarray) – The input array to be checked.

integral(u)[source]

Computes the integral of a given scalar function u over the domain

\[\int_\Omega u \,d\Omega\]
Parameters:

u (np.ndarray) – Function belonging to the grid

Returns:

value – Integral over the domain

Return type:

float

class pyforce.tools.backends.LoopProgress(msg: str, final: float = 100)[source]

Bases: object

A class to make progress bar.

Parameters:
  • msg (str) – Message to be displayed

  • final (float, optional (Default = 100)) – Maximum value for the iterations

update(step: float, percentage: bool = False)[source]

Update message to display and clears the previous one.

Parameters:
  • step (float) – Interger or float value to add at the counter.

  • percentage (boolean, optional (Default = False)) – Indicates if the bar should be displayed in %.

class pyforce.tools.backends.Timer[source]

Bases: object

start()[source]

Start a new timer

stop()[source]

Stop the timer, and report the elapsed time

exception pyforce.tools.backends.TimerError[source]

Bases: Exception

A custom exception used to report errors in use of Timer class

pyforce.tools.functions_list module

class pyforce.tools.functions_list.FunctionsList(dofs: int | None = None, snap_matrix: ndarray | None = None)[source]

Bases: object

A class wrapping a list of functions. They are stored as a list of np.ndarray.

Parameters:

dofs (int) – Degrees of freedom of the functions \(\mathcal{N}_h\).

_list

List containing the functions as np.ndarray.

Type:

list

fun_shape

Number of degrees of freedom of the functions \(\mathcal{N}_h\).

Type:

int

append(dofs_fun: ndarray) None[source]

Extend the current list by adding a new function. The snapshot is stored in a list as a numpy array, to avoid problems when the number of elements becomes large.

Parameters:

dofs_fun (np.ndarray) – Functions to be appended.

build_from_matrix(matrix: ndarray) None[source]

Build the list from a matrix \(S\in\mathbb{R}^{\mathcal{N}_h\times N_s}\). The matrix is transposed to match the shape of the functions.

Parameters:

matrix (np.ndarray) – Matrix containing the functions as columns.

clear() None[source]

Clear the storage.

copy()[source]

Defining the copy of the _list of elements

delete(idx: int) None[source]

Delete a single element in position idx.

Parameters:

idx (int) – Integers indicating the position inside the _list to delete.

lin_combine(vec: ndarray) ndarray[source]

Linearly combine functions (a_i = vec[i]) in the list \(\{\phi_i\}_{i=0}^N\):

\[\sum_{i=0}^N a_i\cdot \phi_i\]

given N = len(vec) and \(\mathbf{a}\in\mathbb{R}^N\).

Parameters:

vec (np.ndarray) – Iterable containing the coefficients of the linear combination.

Returns:

combination – Function object storing the result of the linear combination

Return type:

np.ndarray or Function

max(axis: int | None = None) ndarray[source]

Returns the maximum of the functions in the list along the specified axis.

Parameters:

axis (int, optional (Default=None)) – Axis along which to compute the maximum. If None, computes the maximum over the entire array.

Returns:

max_values – Maximum values of the functions in the list.

Return type:

np.ndarray

mean(axis: int | None = None) ndarray[source]

Returns the mean of the functions in the list along the specified axis.

Parameters:

axis (int, optional (Default=None)) – Axis along which to compute the mean. If None, computes the mean over the entire array.

Returns:

mean_values – Mean values of the functions in the list.

Return type:

np.ndarray

min(axis: int | None = None) ndarray[source]

Returns the minimum of the functions in the list along the specified axis.

Parameters:

axis (int, optional (Default=None)) – Axis along which to compute the minimum. If None, computes the minimum over the entire array.

Returns:

min_values – Minimum values of the functions in the list.

Return type:

np.ndarray

plot(grid: pyvista.UnstructuredGrid, idx_to_plot: int, varname: str = 'u', clim: tuple | None = None, cmap: str = 'jet', resolution: int = [1080, 720], title: str | None = None, **kwargs) pyvista.Plotter[source]

Plot a function from the list.

Parameters:
  • grid (pyvista.UnstructuredGrid) – The PyVista grid on which the function will be plotted.

  • idx_to_plot (int) – Index of the function to plot.

  • varname (str) – The name to assign to the data in the PyVista grid (default is ‘u’).

  • clim (tuple, optional) – The color limits for the plot. If None, the limits will be automatically determined from the data.

  • cmap (str, optional) – The colormap to use for the plot (default is ‘jet’). Other options include ‘viridis’, ‘plasma’, ‘inferno’, etc.

  • resolution (list, optional) – The resolution of the plot in pixels (default is [1080, 720]).

  • zoom (float, optional) – The zoom level for the plot (default is 1.0).

  • title (str, optional) – The title of the plot. If None, no title will be displayed.

  • **kwargs (dict, optional) – Additional keyword arguments passed to the PyVista plotting functions.

plot_sequence(grid: pyvista.UnstructuredGrid, sampling: int = 1, varname: str = 'u', clim: tuple | None = None, cmap: str = 'jet', resolution: int = [1080, 720], title: str | None = None, title_size: int = 20, view: str = 'xy', **kwargs) pyvista.Plotter[source]

Plot a sequence of functions from the list with a specified sampling rate.

Parameters:
  • grid (pyvista.UnstructuredGrid) – The PyVista grid on which the functions will be plotted.

  • sampling (int, optional (Default = 1)) – The sampling rate for the functions to be plotted. If sampling=1, all functions are plotted; if sampling=2, every second function is plotted, and so on.

  • varname (str, optional (Default = 'u')) – The name to assign to the data in the PyVista grid.

  • clim (tuple, optional) – The color limits for the plot. If None, the limits will be automatically determined from the data.

  • cmap (str, optional (Default = 'jet')) – The colormap to use for the plot. Other options include ‘viridis’, ‘plasma’, ‘inferno’, etc.

  • resolution (list, optional (Default = [1080, 720])) – The resolution of the plot in pixels.

  • title (str, optional) – The title of the plot. If None, no title will be displayed.

  • title_size (int, optional (Default = 20)) – The font size of the title.

  • view (str, optional (Default = 'xy')) – The view direction for the plot. Options include ‘xy’, ‘xz’, ‘yz’.

  • **kwargs (dict, optional) – Additional keyword arguments passed to the PyVista plotting functions.

Returns:

The PyVista plotter object used for the visualization.

Return type:

pv.Plotter

return_matrix() ndarray[source]

Returns the list of arrays as a matrix \(S\in\mathbb{R}^{\mathcal{N}_h\times N_s}\).

shape() tuple[source]

Returns the shape of the list as a tuple (dofs, number of functions).

sort(order: list) None[source]

Sort the list according to the given order (list of indices).

Parameters:

order (list) – List of indices for the sorting.

std(axis: int | None = None) ndarray[source]

Returns the standard deviation of the functions in the list along the specified axis.

Parameters:

axis (int, optional (Default=None)) – Axis along which to compute the standard deviation. If None, computes the standard deviation over the entire array.

Returns:

std_values – Standard deviation values of the functions in the list.

Return type:

np.ndarray

store(var_name: str, filename: str, order: list | None = None, format: str = 'h5', compression: bool = True) None[source]

Store the functions in the list to a file.

Parameters:
  • var_name (str) – Name of the variable.

  • filename (str) – Name of the file to save.

  • order (list, optional (Default=None)) – List of integers containing the ordered indices.

  • format (str, optional (Default='h5')) – Format of the file to save. It can be either ‘h5’ or ‘npz’.

  • compression (bool, optional (Default=True)) – If True, the data is compressed when saved in h5/npz format.

pyforce.tools.functions_list.train_test_split(params: list, fun_list: FunctionsList, test_size: float = 0.33, random_state: int = 42)[source]

This function can be used to create a train and test using sklearn capabilities.

Parameters:
  • params (list) – List of parameters to be split.

  • fun_list (FunctionsList) – Object containing the functions as a list of arrays to convert.

  • test_size (float) – DimensionFunctional space of the functions (the dofs should be compliant).

  • random_state (int, optional (Default = 42)) – Random seed for the splitting algorithm.

Returns:

  • train_params (list) – List of the train parameters.

  • test_params (list) – List of the test parameters.

  • train_fun (list) – List of the train functions.

  • test_fun (list) – List of the test functions.

pyforce.tools.plotting module

pyforce.tools.plotting.plot_function(grid: pyvista.UnstructuredGrid, snap: ndarray, varname: str = 'u', filename: str | None = None, clim: tuple | None = None, cmap: str = 'jet', resolution: int = [1080, 720], title: str | None = None, **kwargs) pyvista.Plotter[source]

Plots a function on a PyVista grid and saves the plot to a file.

Parameters:
  • grid (pyvista.UnstructuredGrid) – The PyVista grid on which the function will be plotted.

  • snap (np.ndarray) – The function values to be plotted, should match the number of points in the grid.

  • varname (str) – The name to assign to the data in the PyVista grid (default is ‘u’).

  • filename (str, optional) – The name of the file to save the plot. If None, the plot will not be saved.

  • clim (tuple, optional) – The color limits for the plot. If None, the limits will be automatically determined from the data.

  • cmap (str, optional) – The colormap to use for the plot (default is ‘jet’). Other options include ‘viridis’, ‘plasma’, ‘inferno’, etc.

  • resolution (list, optional) – The resolution of the plot in pixels (default is [1080, 720]).

  • zoom (float, optional) – The zoom level for the plot (default is 1.0).

  • title (str, optional) – The title of the plot. If None, no title will be displayed.

  • **kwargs (dict, optional) – Additional keyword arguments passed to the PyVista plotting functions.

Returns:

The PyVista plotter object used for the visualization.

Return type:

pv.Plotter

pyforce.tools.scalers module

class pyforce.tools.scalers.MinMaxScaler(wrt_ic: bool = False, global_scale: bool = False)[source]

Bases: ScalerWrapper

Scaler that scales the data to a specified range (default is [0, 1]).

Parameters:
  • wrt_ic (bool, optional) – If True, the scaling is done with respect to the initial condition - element 0 of the input data (default is False).

  • global_scale (bool, optional) – If True, the scaling is done globally across all data (default is False).

fit(train_X: FunctionsList)[source]

Fit the MinMaxScaler to the training data.

Parameters:

train_X (FunctionsList) – The training data to fit the scaler to.

inverse_transform(X: FunctionsList)[source]

Inverse transform the data using the fitted MinMaxScaler.

Parameters:

X (FunctionsList) – The data to inverse transform.

Returns:

The inverse transformed data.

Return type:

FunctionsList

transform(X: FunctionsList)[source]

Transform the data using the fitted MinMaxScaler.

Parameters:

X (FunctionsList) – The data to transform.

Returns:

The transformed data.

Return type:

FunctionsList

class pyforce.tools.scalers.ScalerWrapper[source]

Bases: object

Base class for scaling functions in the pyforce library.

fit(train_X: FunctionsList, **kwargs)[source]

Fit the scaler to the data.

Parameters:

train_X (FunctionsList) – The data to fit the scaler to.

inverse_transform(X: FunctionsList)[source]

Inverse transform the data using the fitted scaler.

Parameters:

X (FunctionsList) – The data to inverse transform.

Returns:

The inverse transformed data.

Return type:

FunctionsList

transform(X: FunctionsList)[source]

Transform the data using the fitted scaler.

Parameters:

X (FunctionsList) – The data to transform.

Returns:

The transformed data.

Return type:

FunctionsList

class pyforce.tools.scalers.StandardScaler(wrt_ic: bool = False, global_scale: bool = False, with_std: bool = True)[source]

Bases: ScalerWrapper

Scaler that centers the data to have zero mean and unit variance (if with_std is True).

Parameters:
  • wrt_ic (bool, optional) – If True, the scaling is done with respect to the initial condition - element 0 of the input data (default is False).

  • global_scale (bool, optional) – If True, the scaling is done globally across all data (default is False).

  • with_std (bool, optional) – If True, the data is scaled to have unit variance (default is True).

fit(train_X: FunctionsList)[source]

Fit the StandardScaler to the training data.

Parameters:

train_X (FunctionsList) – The training data to fit the scaler to.

inverse_transform(X: FunctionsList)[source]

Inverse transform the data using the fitted StandardScaler.

Parameters:

X (FunctionsList) – The data to inverse transform.

Returns:

The inverse transformed data.

Return type:

FunctionsList

transform(X: FunctionsList)[source]

Transform the data using the fitted StandardScaler.

Parameters:

X (FunctionsList) – The data to transform.

Returns:

The transformed data.

Return type:

FunctionsList

pyforce.tools.write_read module

pyforce.tools.write_read.ImportFunctionsList(filename: str, format: str = 'h5', return_var_name: bool = False)[source]

This function can be used to load from xdmf/h5 files scalar or vector list of functions.

Parameters:
  • filename (str) – Name of the file to read as xdmf/h5 file.

  • format (str, optional (Default = 'h5')) – Format of the file to read. It can be either ‘h5’, ‘npy’, or ‘npz’.

  • return_var_name (bool, optional (Default = False)) – If True, the variable name is returned along with the FunctionsList.

Returns:

snap – Imported list of functions.

Return type:

FunctionsList

class pyforce.tools.write_read.ReadFromOF(path: str, skip_zero_time: bool = False, decomposed_case: bool = False)[source]

Bases: object

A class used to import data from OpenFOAM - either decomposed or not, single region. Note: “Decomposed” here refers exclusively to parallel MPI domain decomposition (i.e., data split into processor* directories), not multi-region setups.

Either the fluidfoam library (see https://github.com/fluiddyn/fluidfoam) or the foamlib (see https://github.com/gerlero/foamlib) or pyvista are exploited for the import process

Both org and com version of OpenFOAM are supported.

Parameters:
  • path (str) – Path of the OpenFOAM case directory.

  • skip_zero_time (bool, optional) – If True, the zero time folder is skipped. Default is False.

  • decomposed_case (bool, optional) – If True, the case is considered as decomposed. Default is False.

import_field(var_name: str, import_mode: str = 'pyvista', target_times: list[str] | None = None, extract_cell_data: bool = True, verbose: bool = True) tuple[FunctionsList, ndarray][source]

Importing time instances from OpenFOAM directory, if not specified, for all time folders.

Three methods are available for the import process: pyvista, fluidfoam, and foamlib. The latter two are typically faster than pyvista, especially for large cases.

If you want to import specific time instances, provide their folder names in the target_times list.

If you want to import point data instead of cell data, set extract_cell_data to False, only for pyvista method.

Parameters:
  • var_name (str) – Name of the field to import.

  • import_mode (str, (Default: 'pyvista')) – The mode of import, either ‘pyvista’ or ‘fluidfoam’ or ‘foamlib’.

  • target_times (list[str], optional) – List of time folders to import. If None, all time instances are imported.

  • extract_cell_data (boolean, (Default: True)) – If True, the cell data from centroids is extracted; otherwise, point data is extracted.

  • verbose (boolean, (Default = True)) – If True, printing is enabled.

Returns:

  • field (FunctionsList) – Imported list of functions (each element is a numpy.ndarray), sorted in time.

  • time_instants (np.ndarray) – Sorted list of time.

mesh()[source]

Returns the mesh of the OpenFOAM case using pyvista capabilities.

Returns:

mesh – The mesh of the OpenFOAM case.

Return type:

pyvista.UnstructuredGrid

save_mesh(filename: str)[source]

Saves the mesh of the OpenFOAM case to a vtk-file.

Parameters:

filename (str) – Name of the file to save the mesh.

class pyforce.tools.write_read.ReadFromOFMultiRegion(path: str, skip_zero_time: bool = False, decomposed_case: bool = False)[source]

Bases: object

A class used to import data from OpenFOAM - either decomposed or not, multi region. Note: “Decomposed” here refers exclusively to parallel MPI domain decomposition (i.e., data split into processor* directories), not multi-region setups.

Either the fluidfoam library (see https://github.com/fluiddyn/fluidfoam) or the foamlib (see https://github.com/gerlero/foamlib) or pyvista are exploited for the import process.

Both org and com version of OpenFOAM are supported.

Parameters:
  • path (str) – Path of the OpenFOAM case directory.

  • skip_zero_time (bool, optional) – If True, the zero time folder is skipped. Default is False.

  • decomposed_case (bool, optional) – If True, the case is considered as decomposed. Default is False (pyvista only).

import_field(var_name: str, import_mode: str = 'pyvista', target_times: list[str] | None = None, regions_to_import: str | list[str] | None = None, verbose: bool = True, inner_verbose: bool = False) tuple[FunctionsList, ndarray][source]

Importing all time instances (skipping zero folder) from OpenFOAM directory for all available regions, if not skip.

Three methods are available for the import process: pyvista, fluidfoam, and foamlib. The latter two are typically faster than pyvista, especially for large cases.

If you want to import specific time instances, provide their folder names in the target_times list.

Only cell data extraction is supported for multi-region cases: to prevent data misalignment issues.

Note on parallel-decomposed cases: For multi-region cases that are also parallel-decomposed, the import process will automatically handle the reconstruction of the field across all processors for each region, there might be visualisation issues for ‘foamlib’ and ‘fluidfoam’ if the decomposition is not uniform and more advanced techniques are used. If needed, exploit ‘reconstructPar’ from OpenFOAM. Decomposition Methods tested: * hierarchical * simple * … (to be extended if needed)

Parameters:
  • var_name (str) – Name of the field to import.

  • import_mode (str, optional) – Method to use for importing the data. It can be either ‘pyvista’ or ‘fluidfoam’. Default is ‘pyvista’.

  • target_times (list[str], optional) – List of time folders to import. If None, all time instances are imported.

  • regions_to_import (str | list[str], optional) – If specified, only the given region(s) are imported. If None, all valid regions for the specified field are imported.

  • verbose (boolean, (Default = True)) – If True, printing is enabled

  • inner_verbose (boolean, (Default = False)) – If True, verbose mode is enabled for the inner import functions (e.g., _import_region_field_pyvista), otherwise only the outer loop progress is printed.

Returns:

  • snaps (FunctionsList) – Imported list of functions (each element is a numpy.ndarray), sorted in time.

  • time_instants (np.ndarray) – Sorted list of time.

mesh(regions: list[str] | None = None, decomposed_mode: bool = False)[source]

Returns the combines mesh of regions of the OpenFOAM case.

Moreover, the cumulative number of cells for each region is stored in self.ncells_cumulative to facilitate the data extraction process.

Parameters:
  • regions (str, optional) – If specified, only the mesh of the given region is returned. If None, the combined mesh of all regions is returned.

  • decomposed_mode (bool) – If True, the mesh is read from processor* folders. Default is False.

Returns:

mesh – The mesh of the OpenFOAM case.

Return type:

pyvista.UnstructuredGrid

save_mesh(filename: str, region: str | list[str] | None = None)[source]

Saves the mesh (either all regions or a specific list) of the OpenFOAM case to a vtk-file.

Parameters:
  • filename (str) – Name of the file to save the mesh.

  • regions (str | list[str], optional) – If specified, only the mesh of the given region is saved. If None, the combined mesh of all regions is saved.

pyforce.tools.write_read.convert_cell_data_to_point_data(grid: pyvista.UnstructuredGrid, snaps: FunctionsList)[source]

Convert cell data to point data for a given pyvista UnstructuredGrid and corresponding FunctionsList.

Parameters:
  • grid (pyvista.UnstructuredGrid) – The input mesh with cell data.

  • snaps (FunctionsList) – The FunctionsList containing cell data.

Returns:

new_snaps – The FunctionsList containing point data.

Return type:

FunctionsList

pyforce.tools.write_read.convert_point_data_to_cell_data(grid: pyvista.UnstructuredGrid, snaps: FunctionsList)[source]

Convert point data to cell data for a given pyvista UnstructuredGrid and corresponding FunctionsList.

Parameters:
  • grid (pyvista.UnstructuredGrid) – The input mesh with point data.

  • snaps (FunctionsList) – The FunctionsList containing point data.

Returns:

new_snaps – The FunctionsList containing cell data.

Return type:

FunctionsList

pyforce.tools.write_read.get_candidate_channel_all_points(all_grids: pyvista.UnstructuredGrid, candidate_grid: pyvista.UnstructuredGrid | List[pyvista.UnstructuredGrid], candidate_xy: ndarray | list[ndarray], tol=1e-06)[source]

Get all points along a channel (fixed XY, any Z) and corresponding indices in all_grids. This returns all points with XY matching the candidate_xy within a tolerance.

Parameters:
  • all_grids (pv.UnstructuredGrid) – The global mesh.

  • candidate_grid (pv.UnstructuredGrid | list of pv.UnstructuredGrid) – Target region grid(s). Used to ensure the candidate point(s) are belonged to this(these) region(s).

  • candidate_xy (np.ndarray or list of np.ndarray, shape (2,) or (N,2)) – Target XY coordinates of points. Can be a single point or multiple points.

  • tol (float) – Tolerance for matching XY coordinates. Notice: it should be tuned according to different geometric discretisation.

Returns:

  • candidate_points (np.ndarray, shape (M,3)) – All points in the global mesh with matching XY.

  • candidate_idx (np.ndarray, shape (M,)) – Indices of candidate_points in all_grids.

pyforce.tools.write_read.get_candidate_probes(all_grids: pyvista.UnstructuredGrid, candidate_grid: pyvista.UnstructuredGrid | List[pyvista.UnstructuredGrid], candidate_points: ndarray | list[ndarray])[source]

Get target points and corresponding indices within all regions. Step 1: Map the input candidate point(s) to the nearest points in the specified candidate_grid. Step 2: Map the points obtained in Step 1 to the corresponding points in all_grids. This two-step mapping avoids assigning candidate points to incorrect regions due to differences in mesh discretisation.

Parameters:
  • all_grids (pv.UnstructuredGrid) – The grid containing all valid regions.

  • candidate_grid (pv.UnstructuredGrid | list of pv.UnstructuredGrid) – Target region grid(s). Used to ensure the candidate point(s) are belonged to this(these) region(s).

  • candidate_points (np.ndarray | list of np.ndarray, shape (N, 3,) | (3,) each) – Target point(s). Can be a single point or a list of points.

Returns:

  • candidate_points (np.ndarray of shape (N, 3)) – Coordinates of the target points, used for training or post-processing.

  • candidate_idx (np.ndarray shape (M,)) – Indices of target points within all_grids, used for training or post-processing.

pyforce.tools.write_read.get_candidate_regions_points(all_grids: pyvista.UnstructuredGrid, candidate: pyvista.UnstructuredGrid | List[pyvista.UnstructuredGrid], tolerance=1e-06)[source]

Get target region points and corresponding indices within all regions.

Parameters:
  • all_grids (pv.UnstructuredGrid) – The grid containing all valid regions.

  • candidate (pv.UnstructuredGrid | list of pv.UnstructuredGrid) – Target region grid(s). Can be a single grid or a list of grids for multiple regions.

  • tolerance (float) – Distance tolerance used during KDTree matching. Default is 1e-6.

Returns:

  • candidate_points (np.ndarray of shape (N, 3)) – Coordinates of the target points, used for training or post-processing.

  • candidate_idx (np.ndarray shape (M,)) – Indices of target points within all_grids, used for training or post-processing.

Module contents

pyforce/tools.

Basic tools for the pyforce library.