chromo.util package

chromo.util.linalg

Utility functions for linear algebra calculations.

chromo.util.linalg.angle_between_vectors()

Angle between two vectors, v and w, returned in radians.

chromo.util.linalg.arbitrary_axis_rotation()

Rotate about an axis defined by two points.

Notes

Generate a transformation matrix for counterclockwise (right handed convention) rotation of angle rot_angle about an arbitrary axis from points r0 to r1.

TODO: Store rot_vec as attribute of poly (to avoid recreating vector)

Parameters:

• axis (array_like (3,) of double) – Rotation axis

• point (array_like (3,) of double) – Rotation fulcrum

• rot_angle (double) – Angle of rotation. Positive rotation is counterclockwise when rotation axis is pointing directly out of the screen

• poly (poly.TransformedObject) – Container holding the transformation matrix

chromo.util.linalg.generate_translation_mat()

Generate translation matrix.

Notes

Generate the homogeneous transformation matrix for a translation of distance delta_x, delta_y, delta_z in the x, y, and z directions, respectively.

Parameters:

• delta (array_like (3,) of double) – Array of translation magnitudes in x, y, z directions, in form [dx, dy, dz]

• polymer (poly.TransformedObject) – Container holding the transformation matrix

chromo.util.linalg.get_arbitrary_axis_rotation_matrix()

Rotate a point around an arbitrary axis in 3D space.

chromo.util.linalg.get_prism_verticies()

Get the cartesian coordinates of verticies for specified prism geometry.

Parameters:

• num_sides (int) – Number of sides on the face of the prism used to represent the nucleosome’s geometry; this determines the locations of verticies of the Prism

• width (float) – Determines the shape of the prism defining the location of the nucleosome’s verticies. The width gives the diameter of the circle circumscribing the base of the prism in the simulation’s distance units.

• height (float) – Determines the shape of the prism defining the location of the nucleosome’s verticies. The height gives the height of the prism in the simulation’s distance units.

Returns:

Cartesian coordinates of the M verticies of the prism geometry.

Return type:

np.ndarray (M, 3) of double

chromo.util.linalg.get_rotation_matrix()

Find the rotation matrix that aligns two pairs of vectors.

Notes

This function provides a transformation matrix that rotates t3_local to align with t3 while simultaneously rotating t2_local to align with t2. To do this, the function first finds the rotation matrix that aligns t3_local to t3. Once that first rotation is applied, the function then finds the rotation matrix that aligns the intermediately-rotated t2_local to t2. The final rotation matrix is the product of these two rotations.

This function does not work with all possible input vectors. It is possible that the resulting rotation matrix will generate a rotated t3_local vector that is not aligned with t3, but rather is aligned with -t3. This is because the rotation matrix is not unique. This function checks for this possibility and corrects the rotation matrix accordingly.

chromo.util.linalg.return_arbitrary_axis_rotation()

Generate & return homogeneous transform matrix of arbitrary axis rot.

Notes

Use this function when you want to return a transformation matrix, not just update the transformation_mat attribute of a polymer.

This function is not yet tested.

Parameters:

• axis (array_like (3,) of double) – Rotation axis

• point (array_like (3,) of double) – Rotation fulcrum

• rot_angle (double) – Angle of rotation. Positive rotation is counterclockwise when rotation axis is pointing directly out of the screen

Returns:

double[ – Homogeneous transformation matrix encoding arbitrary axis rotation.

Return type:

, ::1]

chromo.util.linalg.rotation_matrix_from_vectors()

Find the rotation matrix that aligns vec1 to vec2

Notes

This method is based on the following Stack Overflow answer: https://stackoverflow.com/a/67767180 from the forum: “Calculate rotation matrix to align two vectors in 3D space.”

Parameters:

• vec1 (array_like (3,) of double) – A 3D “source” vector

• vec2 (array_like (3,) of double) – A 3D “destination” vector

Returns:

A 3x3 rotation matrix that rotates vec1 to vec2

Return type:

array_like (3, 3) of double

chromo.util.reproducibility

A module with routines designed to make code “automatically” reproducible.

Notes

Reproducible here should be taken in the scientific sense. We want to, given the output of our simulation, always be able to

1. Reconstruct the inputs to the simulation exactly.

2. Run our code again with those inputs to generate that exact same output.

For our purposes, this means saving

1. The Git commit hash of the simulation used.

2. All random seeds used.

3. Serializing all inputs used.

If in doubt, you probably want to simply wrap any simulation function you write with make_reproducible.

chromo.util.reproducibility.builtin_types = [<class 'bool'>, <class 'bytes'>, <class 'bytearray'>, <class 'complex'>, <class 'float'>, <class 'int'>, <class 'str'>]

Types to be saved as CSV entries by default.

Notes

These are basically just the built-ins that permit simple save/load cycling via CSV entries with Pandas. Some extras are included in the split_builtin_params function itself. Please see that code for a full description.

The following strategy would also work, but would lead to us having to pickle some builtin things, like “map”, “Error”, etc.. builtin_types = tuple(getattr(builtins, t) for t in dir(builtins) if isinstance(getattr(builtins, t), type))

These additional types are not allowed by default because they make loading in the CSV a little trickier: …, frozenset, tuple]

chromo.util.reproducibility.get_unique_subfolder(root)

Create and return a thread-safe unique folder.

Notes

Uses the fact that mkdir is required to be atomic on POSIX-compliant system to make sure that two threads aren’t given the same folder.

Parameters:

root (str or Path) – The base of the filename which we should try appending numbers to until we find a number which doesn’t already exist.

Returns:

The folder which was created.

Return type:

Path

chromo.util.reproducibility.get_unique_subfolder_name(root)

Return name of a thread-safe unique folder.

Parameters:

root (str or Path) – The base of the filename which we should try appending numbers to until we find a number which doesn’t already exist.

Returns:

The path which would provide a unique subfolder

Return type:

Path

chromo.util.reproducibility.parse_type(arg_name: str, value: Any) → Tuple[List[str], List[Any]]

Parse argument data types (for handling kwargs)

Parameters:

• arg_name (str) – Name of the argument

• value (Any) – Value associated with arg_name

Returns:

• List[str] – List of argument names after parsing

• List[Any] – List of argument values after parsing

chromo.util.reproducibility.replicate_chromo_essentials(**kwargs)

Create a compressed copy of essential modules for running a simulation.

Notes

Modules are not compiled; only the essential modules are included in the copy; if a path to the modification patterns and simulation run script are included, they will be stored in the copy; if the run command is specified, a Bash file will be automatically generated for re-creating the chromo Conda environment, reinstalling dependencies, recompiling the essential modules, and re-running the simulation.

The function begins by loading the essential modules from essentials.txt. Then, each essential module is copied to a folder in the output directory. If provided, the run script, modification patterns, and run command will be copied to the simulator copy. The copy of the simulator will then be compressed. A Bash script for installing modules, compiling the simulator code, and re-running the simulation will then be automatically generate.

With this level of reproducibility and acceptable runtimes of 1-3 days per simulation, it will not be necessary to share exact results; we can simply share code for reproducing results when we go to publish this material.

chromo.util.reproducibility.split_builtin_params(sim: Callable[[Any], Any], *args, **kwargs) → Tuple[Dict[str, Any], Dict[str, Any]]

Split parameters between built in and non-built in data types.

Parameters:

sim (Callable[[Any], Any]) – Any function on which make_reproducible is being called.

Returns:

• Dict[str, Any] – Dictionary of parameters for built in data types

• Dict[str, Any] – Dictionary of parameters for non-built in data types

chromo.util.reproducibility.to_file_params(non_builtins_kwargs, folder, suffix='')

Call .to_file for each parameter, handling sequences.

Parameters:

• non_builtins_kwargs (Dict[str, RoundTrippable]) – A RoundTrippable must implement the .name/.to_file/.from_file interface, or be a standard type, such as a numpy array or a DataFrame.

• folder (Path) – The folder to save the output to.

chromo.util.__init__

chromo.util.combine_repeat(a, idx)

Combine the repeat entries in an array of bin densities.

Flatten one dimension of the value array a so each column represents a different reader protein, but densities in all bins are combined into the same column.

Flatten the index array idx so that indices for each of the eight bins representing the density distribution of a bead are in the same column.

Sort the index array idx and the value array a by the values in idx

Identify the unique bins and the indices where they are first encountered in the sourted index array.

Add the density within each bin (separately for each reader protein).

Parameters:

• a (double[:, :, ::1]) – 3D memoryview representing the distribution of densities for each bead (dimension 0) in the nearest eight neighboring bins (dimension 2) for each reader protein (dimension 1)

• idx (long[:, ::1]) – 2D memoryview of bin indices; the first dimension represents each bead affected, and the second dimension contains the indices of the eight bins containing the bead

Returns:

a_combine – 2D memoryview of total density in each affected bin (dimension 0) of the polymer and each reader protein (dimension 1)

Return type:

double[:, ::1]

chromo.util.poly_paths

Functions for specification of initial polymer paths.

chromo.util.poly_paths.confined_gaussian_walk(num_points: int, step_sizes: ndarray, confine_type: str, confine_length: float) → ndarray

Generate coordinates for Gaussian random walk w/ fixed path length.

Parameters:

• num_points (int) – Number of points in the Gaussian random walk

• step_sizes (np.ndarray (N-1,) of float) – Distance between each point in the random walk

• confine_type (str) – Name of the confining boundary. To indicate model w/o confinement, enter a blank string for this argument

• confine_length (double) – The lengthscale associated with the confining boundary. Length representation specified in function associated w/ confine_type

Returns:

Coordinates of each point in the Gaussian random walk, where rows represent individual points and columns give x, y, z coordinates

Return type:

np.ndarray (N, 3)

chromo.util.poly_paths.confined_uniform_random_walk(num_points: int, step_size: float, confine_type: str, confine_length: float) → ndarray

Generate coordinates for uniform random walk w/ fixed path length.

Parameters:

• num_points (int) – Number of points in the uniform random walk

• step_size (float) – Distance between each point in the random walk

• confine_type (str) – Name of the confining boundary. To indicate model w/o confinement, enter a blank string for this argument

• confine_length (double) – The lengthscale associated with the confining boundary. Length representation specified in function associated w/ confine_type

Returns:

Coordinates of each point in the uniform random walk, where rows represent individual points and columns give x, y, z coordinates

Return type:

np.ndarray (N, 3)

chromo.util.poly_paths.coordinates_in_x_y(num_beads: int, bead_length: float, shape_func: Callable[[float], float], step_size: float) → ndarray

Find bead coordiantes about path in x, y plane w/ fixed contour length.

Notes

Get coordinates in the x, y plane which splits the path of the polymer into segments of equal path length.

TODO: If we want to have variable linker lengths, the spacing of monomeric units must be different, and this must be accounted for when selecting x-positions of beads.

Parameters:

• num_beads (int) – Number of monomer units on the polymer

• bead_length (float or (N,) array_like of float) – The amount of polymer path length between this bead and the next bead. For now, a constant value is assumed (the first value if an array is passed).

• shape_func (Callable[[float], float]) – Shape of the polymer where z = 0 and y = f(x)

• step_size (float) – Step size for numerical evaluation of contour length when domain of the shape function used to initialize polymer

Returns:

Coordinates of each point with path defined by shape_func

Return type:

np.ndarray

chromo.util.poly_paths.coordinates_in_x_y_z(num_beads: int, bead_length: float, shape_func_x: Callable[[float], float], shape_func_y: Callable[[float], float], shape_func_z: Callable[[float], float], step_size: float) → Tuple[ndarray, List[float]]

Generate coordinates for 3D initialization of a polymer.

Parameters:

• num_beads (int) – Number of monomer units on the polymer

• bead_length (float or (N,) array_like of float) – The amount of polymer path length between this bead and the next bead. For now, a constant value is assumed (the first value if an array is passed).

• shape_func_x (Callable[[float], float]) – Parametric functions to obtain the x coordinates of the path

• shape_func_y (Callable[[float], float]) – Parametric functions to obtain the y coordinates of the path

• shape_func_z (Callable[[float], float]) – Parametric functions to obtain the z coordinates of the path

• step_size (float) – Step size for numerical evaluation of contour length when domain of the shape function used to initialize polymer

Returns:

• np.ndarray – Coordinates of points with path defined by parametric shape functions

• List[float] – Parameter values corresponding to (x, y, z) points obtained from shape functions

chromo.util.poly_paths.estimate_tangents_from_coordinates(coordinates: ndarray) → Tuple[ndarray, ndarray]

Estimate t3 and t2 tangent vectors from array of coordinates.

Parameters:

coordinates (np.ndarray) – Ordered coordinates representing path of a polymer

Returns:

• np.ndarray (N, 3) – Matrix representing t3 tangents, where rows represent x, y, z components of tangent at each point.

• np.ndarray (N, 3) – Matrix representing t2 tangents, where rows represent x, y, z coordinates of tangent vector at each point

chromo.util.poly_paths.gaussian_walk(num_steps: int, step_sizes: ndarray) → ndarray

Generate coordinates for Gaussian random walk w/ fixed path length.

Parameters:

• num_steps (int) – Number of steps in the Gaussian random walk

• step_sizes (np.ndarray) – Distances between subsequent points in the random walk

Returns:

Coordinates of each point in the Gaussian random walk, where rows represent individual points and columns give x, y, z coordinates

Return type:

np.ndarray (N, 3)

chromo.util.poly_paths.get_tangent_vals_x_y(x: ndarray, shape_func: Callable[[float], float], step_size: float) → Tuple[ndarray, ndarray]

Get the t3 and t2 vectors at a fixed position along shape_func.

Parameters:

• x (np.ndarray) – Vector of independent variable positions at which to evaluate tangent

• shape_func (Callable[[float], float]) – Function defining the shape on which to obtain tangent

• step_size (float) – Step size to use when numerically evaluating tangent vectors

Returns:

• np.ndarray (N, 3) – Matrix representing t3 tangents, where rows represent x, y, z components of tangent at each point.

• np.ndarray (N, 3) – Matrix representing t2 tangents, where rows represent x, y, z coordinates of tangent vector at each point

chromo.util.poly_paths.get_tangent_vals_x_y_z(t: ndarray, shape_func_x: Callable[[float], float], shape_func_y: Callable[[float], float], shape_func_z: Callable[[float], float], step_size: float, r: Optional[ndarray] = None) → Tuple[ndarray, ndarray]

Get the t3 and t2 vectors at a fixed position along shape_func.

Parameters:

• t (np.ndarray) – Vector of parameter values at which to evaluate the tangents

• r (Optional[np.ndarray]) – Coordinates of points corresponding to parameter values

• shape_func_x (Callable[[float], float]) – Parametric functions to obtain the x coordinates of the path

• shape_func_y (Callable[[float], float]) – Parametric functions to obtain the y coordinates of the path

• shape_func_z (Callable[[float], float]) – Parametric functions to obtain the z coordinates of the path

• step_size (float) – Step size to use when numerically evaluating tangent vectors

Returns:

• array_like (N, 3) of double – Matrix representing t3 tangents, where rows represent x, y, z components of tangent at each point.

• array_like (N, 3) of double – Matrix representing t2 tangents, where rows represent x, y, z coordinates of tangent vector at each point

chromo.util.poly_paths.in_confinement(point, confine_type, confine_length)

Check if a proposed point lies inside a specified confinement.

Parameters:

• point (array_like (1, 3)) – Point for which to evaluate position in confinement

• confine_type (str) – Name of the confining boundary. To indicate model without confinement, enter a blank string for this argument.

• confine_length (double) – The lengthscale associated with the confining boundary. What the length represents is specified in the function associated with confine_type.

Returns:

Indicator for whether the point lies in the confinement (True) or outside the confinement (False)

Return type:

bool

chromo.util.poly_paths.in_cubical_confinement(point, confine_length)

Check if a proposed point lies inside a spherical confinement.

Parameters:

• point (array_like (1, 3)) – Point for which to evaluate position in confinement.

• confine_length (double) – Edge length of cubical confimement.

Returns:

Indicator for whether the point lies in the confinement (True) or outside the confinement (False)

Return type:

bool

chromo.util.poly_paths.in_spherical_confinement(point, boundary_radius)

Check if a proposed point lies inside a spherical confinement.

Parameters:

• point (array_like (1, 3)) – Point for which to evaluate position in confinement.

• boundary_radius (double) – Radial distance of the confining boundary from the origin.

Returns:

Indicator for whether the point lies in the confinement (True) or outside the confinement (False)

Return type:

bool

chromo.util.poly_paths.looped_confined_gaussian_walk(num_points: int, step_size: float, confine_type: str, confine_length: float) → ndarray

Generate looped, confined Gaussian random walk w/ fixed path length.

Parameters:

• num_points (int) – Number of points in the Gaussian random walk

• step_size (float) – Distance between each point in the random walk

• confine_type (str) – Name of the confining boundary. To indicate model w/o confinement, enter a blank string for this argument

• confine_length (double) – The lengthscale associated with the confining boundary. Length representation specified in function associated w/ confine_type

Returns:

Coordinates of each point in the Gaussian random walk, where rows represent individual points and columns give x, y, z coordinates

Return type:

np.ndarray (N, 3)

chromo.util.poly_paths.numerical_derivative(shape_func: Callable[[float], float], point: float, step_size: float) → float

Numerically evaluate the derivative of shape_func at a point.

Parameters:

• shape_funct (Callable[[float], float]) – Function from which to evaluate derivative

• point (float) – Point at which to evaluate derivative

• step_size (float) – Step-size to apply when numerically evaluating derivative

Returns:

Numerical approximation of derivative at point

Return type:

float

chromo.util.bead_selection

Utility functions for random bead selection.

chromo.util.bead_selection.from_left()

Randomly select index exponentially decaying from left.

Parameters:

• window (int) – Bead window size for selection (must be less than polymer length N_beads)

• N_beads (int) – Number of beads in the polymer chain

Returns:

Index of a bead selected with exponentially decaying probability from first position in the chain.

Return type:

int

chromo.util.bead_selection.from_point()

Randomly select index w/ exponentially decaying probability from point.

Parameters:

• window (int) – Bead window size for selection (must be less than polymer length N_beads)

• N_beads (int) – Number of beads in the polymer chain

• ind0 (int) – Index of first point from which to exponentially decaying probability with distance

Returns:

Index of new point selected based on distance from ind0

Return type:

int

chromo.util.bead_selection.from_right()

Randomly select index exponentially decaying from right.

Parameters:

• window (int) – Bead window size for selection (must be less than polymer length N_beads)

• N_beads (int) – Number of beads in the polymer chain

Returns:

Index of a bead selected with exponentially decaying probability from last position in the chain.

Return type:

int

chromo.util.poly_stat

Average Polymer Statistics

Generate average polymer statistics from Monte Carlo simulation output. This module defines a PolyStat object, which loads polymer configurations from a Polymer object. PolyStat can be used to sample beads from the polymer and calculate basic polymer statistics such as mean squared end-to-end distance and mean 4th power end-to-end distance.

Bead sampling methods include overlapping sliding windows and non-overlapping sliding windows. Overlapping sliding windows offer the benefit of increased data, though the results are biased by central beads which exist in multiple bins of the average. Non-overlapping sliding windows reduce the bias in the results, but include fewer samples in the average.

Joseph Wakim Spakowitz Lab Modified: June 17, 2021

class chromo.util.poly_stat.PolyStats(r, lp, sampling_method: str)

Bases: object

Class representation of the polymer statistics analysis toolkit.

calc_r2(windows: ndarray) → float

Calculate the average squared end-to-end distance of the polymer.

Mean squared end-to-end distance is non-dimensionalized by the persistence length of the polymer, dividing the dimensional quantity by ((2 * self.lp) ** 2).

Parameters:

windows (np.ndarray (N, 2)) – Windows of bead indices for which the average squared end-to-end distance will be calculated

Returns:

Average squared end-to-end distance for specified bead windows

Return type:

float

calc_r4(windows: ndarray) → float

Calculate the average 4th power end-to-end distance of the polymer.

Mean 4th power end-to-end distance is non-dimensionalized by the persistence length of the polymer, dividing the dimensional quantity by ((2 * self.lp) ** 4).

Parameters:

windows (np.ndarray (N, 2)) – Windows of bead indices for which the average 4th moment end-to-end distance will be calculated

Returns:

Average 4th power end-to-end distance for specified bead windows

Return type:

float

load_indices(bead_separation: int) → ndarray

Load bead indices for windows in the average.

Parameters:

bead_separation (int) – Separation of beads in windows for which average is calcualated.

Returns:

Pairs of bead indicies for windows of statistics

Return type:

np.ndarray (N, 2)

chromo.util.poly_stat.calc_mean_r2(output_dir: str, output_files: List[str], window_size: int) → float

Calculate the mean squared end-to-end distance from configuration files.

Parameters:

• output_dir (str) – Path to the configuration output directory for a particular simulation being analyzed

• output_files (List[str]) – List of output files from which to generate polymer statistics

• window_size (int) – Spacing between beads for which to calculate average squared end-to-end distance

chromo.util.poly_stat.find_polymers_in_output_dir(directory: Optional[str] = '.') → List[str]

Obtain the number of polymers in directory.

Parameters:

directory (Optional[str]) – Path to directory into which the number of polymers will be retrieved

Returns:

List of polymer names

Return type:

List[str]

chromo.util.poly_stat.get_latest_configuration(polymer_prefix: Optional[str] = 'Chr-1', directory: Optional[str] = '.') → str

Identify file path to latest polymer configuration in output directory.

Parameters:

• polymer_prefix (str) – Polymer identifier in a configuration file, preceeding snapshot name

• directory (str) – Path to the output directory containing polymer configuration files

Returns:

Path to the file containing the latest polymer configuration file

Return type:

str

chromo.util.poly_stat.get_latest_simulation(directory: Optional[str] = '.')

Get the latest simulation file from directory.

Parameters:

directory (Optional[str]) – Path to directory into which the latest simulation will be retrieved

Returns:

Most recent simulation output directory.

Return type:

str

chromo.util.poly_stat.jump_sample(bead_separation, num_beads)

Generate list of bead index pairs for non-overlaping sampling scheme.

If the end of the polymer does not complete a bin, it is excluded from the average.

Parameters:

• bead_separation (int) – Number of beads in window for average calculation

• num_beads (int) – Number of beads in the polymer chain

Returns:

windows – Pairs of bead indicies for windows of statistics

Return type:

array_like (N, 2)

chromo.util.poly_stat.list_output_files(output_dir: str, equilibration_steps: int) → List[str]

List configuration files in output directory following equilibration.

Begins by listing all files in the output directory, then filters list to configuration files, and finally selects configuration snapshots following equilibration.

Parameters:

• output_dir (str) – Output directory containing simulated configuration files

• equilibration_steps (int) – Number of Monte Carlo steps to ignore due to equilibration

Returns:

List of files in the output directory corresponding to configurations following equilibration.

Return type:

List[str]

chromo.util.poly_stat.overlap_sample(bead_separation, num_beads)

Generate list of bead index pairs for sliding window sampling scheme.

Parameters:

• bead_separation (int) – Number of beads in window for average calculation

• num_beads (int) – Number of beads in the polymer chain

Returns:

windows – Pairs of bead indicies for windows of statistics

Return type:

array_like (N, 2)

chromo.util.poly_stat.save_summary_statistics(window_sizes: List[int], polymer_stats: List[float], save_path: str)

Save summary polymer statistics to an output file.

Parameters:

• window_sizes (List[int]) – Sequence of window sizes for which summary statistics are reported.

• polymer_stats (List[float]) – Polymer statistics being saved

• save_path (str) – Path to file in which to save summary statistics

chromo.util.mc_stat

Utilities for tracking MC move acceptance rates.

class chromo.util.mc_stat.AcceptanceTracker(log_dir: str, log_file_prefix: str, moves_in_average: float)

Bases: Tracker

Class representation of MC acceptance rate tracker.

Tracks acceptance rate of attempted MC moves of a particular type. Downweights historical acceptance rates with each step to more heavily consider recent acceptance rates. Logs the acceptance rate, move amplitude, and bead selection amplitude with each step.

This class will provide metrics to the MC adapter, which will dynamically adjust move and bead amplitudes.

create_log_file(ind: int)

Create the log file tracking amplitudes and acceptance rates.

The log file will have five columns: snapshot count, iteration count bead selection amplitude, move amplitude, and move acceptance rate.

Snapshot count and will be recorded as constant values for each snapshot. The iteration count, bead selection amplitude, move amplitude, and move acceptance rate will be recorded for each iteration in the snapshot.

This method simply creates the log file and adds the column labels. If additional properties are to be logged, add them as column labels in this method, and output them using the save_move_log method.

Parameters:

ind (int) – Index with which to append log file name

log_move(amp_move_limit: float, amp_bead_limit: int, amp_move: float, amp_bead: int, dE: float)

Add a proposed move and current acceptance rate to the log.

Parameters:

• amp_move_limit (float) – Maximum move amplitude allowed at current MC step

• amp_bead_limit (int) – Maximum bead selection amplitude allowed at current MC step

• amp_move (float) – Amplitude of the proposed move

• amp_bead (int) – Selection amplitude of the proposed move

• dE (float) – Change in energy associated with move

save_move_log(snapshot: int)

Save the move and bead amplitude log to a file and clear lists.

Parameters:

snapshot (int) – Current snapshot number, with which to label

update_acceptance_rate(accept: float, log_update: int)

Incrementally update acceptance rate based on move acceptance.

Apply an exponentially weighted moving average to maintain a running measure of move acceptance weight.

Parameters:

• accept (int) – Binary indicator of move acceptance (1.) or rejection (0.)

• log_update (bint) – Indicator for whether (1) or not (2) to record the updated acceptance rate after the MC move

class chromo.util.mc_stat.ConfigurationTracker(log_path: str, initial_config: ndarray)

Bases: Tracker

Track progressive reconfiguraion of the polymer.

create_log_file()

Create an output file for the tracker.

log_move(config: ndarray)

Log the change in polymer configuration.

Parameters:

config (np.ndarray (N, 3)) – Current polymer configuration; rows represent individual beads and columns represent (x, y, z) coordinates

save_move_log(snapshot: int)

Save to the ConfigurationTracker output file and reset log.

Parameters:

snapshot (int) – MC snapshot number

class chromo.util.mc_stat.Tracker

Bases: ABC

Abstract class representation of a MC performance tracker.

abstract create_log_file()

Create an output file for the tracker.

abstract log_move()

Log performance after an MC move.

abstract save_move_log()

Save the log to the output file, then reset log.

chromo.util.mc_stat.initialize_table(path: str, columns: List[str])

Initialize a log file with empty columns.

Parameters:

• path (str) – Path to the log file being initiailzed

• columns (List[str]) – Column names for log file

chromo.util.gjk

Gilbert-Johnson-Keerthi (GJK) Collision Detection

Notes

Adapted from matlab function provided by Matthew Sheen (2016) at: https://github.com/mws262/MATLAB-GJK-Collision-Detection/blob/master/GJK.m

In our implementation of the GJK algorithm, point c will be the oldest support point selected, followed by point b, and lastly point a. This way, between iterations, point a will become point b and point b will become point c.

chromo.util.gjk.get_furthest_in_dir(vertices: ndarray, dir_: ndarray) → ndarray

Get the furthest point of a shape from its center in a given direction.

Parameters:

• vertices (np.ndarray (M, 3)) – Verticies of shape

• dir (np.ndarray (3, )) – Direction in which to get the furthest point

Returns:

Cartesian (x, y, z) coordinates of the furthest point in the specified direction.

Return type:

np.ndarray

chromo.util.gjk.gjk_collision(vertices_1: ndarray, vertices_2: ndarray, max_iters: int) → bool

Test for collision by GJK algorithm given arrays of vertex coordinates.

Parameters:

• vertices_1 (array_like (M, 3) of double) – Array of vertices for first object, where the rows represent individual vertices and columns represent cartesian (x, y, z) coordinates for the corresponding vertices

• vertices_2 (array_like (M, 3) of double) – Array of vertices for second object

• max_iters (int) – Maximum iterations of the GJK algorithm to try

Returns:

Flag for collision between objects (True = collision, False = no collision)

Return type:

bool

chromo.util.gjk.pick_line(vertices_1: ndarray, vertices_2: ndarray, dir_: ndarray) → Tuple[ndarray, ndarray]

Pick an initial line along the Minkowski difference between the two shapes.

Parameters:

• vertices_1 (array_like (M, 3) of double) – Array of vertices for first object, where the rows represent individual vertices and columns represent cartesian (x, y, z) coordinates for the corresponding vertices

• vertices_2 (array_like (M, 3) of double) – Array of vertices for second object

• dir (np.ndarray (3, )) – Initial direction with which to select support points

Returns:

Two support points as Cartesian (x, y, z) coordinates

Return type:

np.ndarray (2, 3)

chromo.util.gjk.pick_tetrahedron(a: ndarray, b: ndarray, c: ndarray, vertices_1: ndarray, vertices_2: ndarray, max_iters: int) → bool

Check if 3D simplex (tetrahedron) contains origin, indicating collision.

If two 3D, convex shapes collide, then a 3D simplex can be drawn between the Minkowski difference of the two shapes which contains the origin. This function only needs to be evaluated if a triangle case contains the origin.

Parameters:

• a (array_like (3,) of double) – Coordinates of the newest support point at current iteration

• b (array_like (3,) of double) – Coordinates of the second-newest support point at current iteration

• c (array_like (3,) of double) – Coordinates of the oldest support point at current iteration

• vertices_1 (array_like (M, 3) of double) – Array of vertices for first object, where the rows represent individual vertices and columns represent cartesian (x, y, z) coordinates for the corresponding vertices

• vertices_2 (array_like (M, 3) of double) – Array of vertices for second object

• max_iters (int) – Maximum iterations of the GJK algorithm to try

Returns:

Flag for tetrahedron containing origin (True = contains origin, False = does not contain origin) -> indicator of whether 3D collision occured

Return type:

bool

chromo.util.gjk.pick_triangle(a: ndarray, b: ndarray, vertices_1: ndarray, vertices_2: ndarray, max_iters: int) → Tuple[ndarray, ndarray, ndarray, bool]

Check if 2D simplex (triangle) contains origin.

Notes

If two convex shapes collide, then there exists a triangle formed from vertices of the Minkowski difference which contains the origin.

Begin by initializing a flag indicating whether the generated triangle contains the origin – initialize as false.

Generate the first support point using the first two points specified as arguments. This support point must be selected in the direction towards the origin, which is obtained using a triple product.

Any triangle can be used to partition 2D space into seven Voronoi regions: a, b, c, ab, ac, bc, abc. If our triangle contains the origin, then the orgin is located in region abc. Because of the way we selected our three points of our triangle, we can eliminate certain regions as candidates for containing the origin – we eliminate regions a, b, c, bc. If our triangle contains the origin, then regions ab and ac do NOT contain the origin; to check that our triangle contains the origin, verify that this is the case.

If region ab contain the origin, then the dot product of abp (perpendicular to side ab towards the origin) and ao (the vector from the origin to a) will be positive. Likewise, if region ac contains the origin, then the dot product between acp (perpendicular to side ac towards the origin) and ao will be positive.

Parameters:

• a (array_like (3,) of double) – Coordinates of the newer support point

• b (array_like (3,) of double) – Coordinates of the newer older point

• vertices_1 (array_like (M, 3) of double) – Array of vertices for first object, where the rows represent individual vertices and columns represent cartesian (x, y, z) coordinates for the corresponding vertices

• vertices_2 (array_like (M, 3) of double) – Array of vertices for second object

• max_iters (int) – Maximum iterations of the GJK algorithm to try

Returns:

• np.ndarray (3, ) x3 – Verticies of new support points on triangle

• bool – Flag for triangle containing origin (True = contains origin, False = does not contain origin)

chromo.util.gjk.shape_overlap(s1: ndarray, s2: ndarray) → bool

Check if any two vertices of the two shapes overlap

Parameters:

• s1 (array_like (M, 3) of double) – Array of vertices for first shape, where the rows represent individual vertices and columns represent cartesian (x, y, z) coordinates for the corresponding vertices

• s2 (array_like (M, 3) of double) – Array of vertices for second shaoe

Returns:

Flag indicating overlapping vertices (True = overlapping vertices -> collision, False = non-overlapping vertices -> more analysis needed)

Return type:

bool

chromo.util.gjk.support_function(vertices_1, vertices_2, dir_)

Support function for obtaining support points on Minkowski difference.

Parameters:

• vertices_1 (array_like (M, 3) of double) – Array of vertices for first object, where the rows represent individual vertices and columns represent cartesian (x, y, z) coordinates for the corresponding vertices

• vertices_2 (array_like (M, 3) of double) – Array of vertices for second object

• dir (np.ndarray (3, )) – Initial direction with which to select support points

Returns:

Support point in given direction given vertex

Return type:

np.ndarray (3, )