IRA interface

class ira_mod.IRA(shlib='')[source]

The Iterative Rotatitions and Assignments (IRA) algorithm is a shape matching algorithm designed for matching generic atomic structures. It solves the problem:

\[P_B B = \mathbf{R} A + \mathbf{t}\]

where \(A\) and \(B\) are two atomic structures, \(\mathbf{R}\) is a rotation matrix, \(\mathbf{t}\) is a translation vector, and \(P_B\) is a permutation matrix.

For futher details please see the publication:

M Gunde, N Salles, A Hemeryck, L Martin Samos: J. Chem. Inf. Model. 2021, 61, 11, 5446–5457 DOI: https://doi.org/10.1021/acs.jcim.1c00567 arXiv: https://export.arxiv.org/abs/2111.00939

or the thesis: M Gunde: “Development of IRA : a shape matching algorithm, its implementation, and utility in a general off-lattice kMC kernel”, 2021 url: https://theses.hal.science/tel-03635139/

In order to use the module, the path to this file should be added to the environment variable PYTHONPATH:

export PYTHONPATH=/your/path/to/IRA/interface

Example import:

>>> import ira_mod
>>> sofi = ira_mod.IRA( )
cshda(nat1, typ1, coords1, nat2, typ2, coords2, lat=None, thr=None)[source]

Constrained Shortest Distance Assignment (CShDA) algorithm, is a LAP (Linear Assignment Problem) solver, which constrains the atoms to be assigned such that distances are minimised locally. This is sometimes referred to as the “inverse bottleneck assignment”.

CShDA can be used to compute the assignment of atoms and the corresponding distances between atoms in permuted structures. It is used inside the IRA algorithm.

This is a wrapper to the libira_cshda and libira_cshda_pbc routines from library_ira.f90

Note

Requirement: nat1 \(\le\) nat2

== input: ==

Parameters

  • nat1 (int) – number of atoms of structure 1

  • typ1 (np.ndarray( nat1, dtype = int or string )) – atomic types of structure 1

  • coords1 (np.ndarray( (nat1, 3), dtype = float )) – atomic positions of structure 1

  • nat2 (int) – number of atoms of structure 2

  • typ2 (atomic types of structure 2) – np.ndarray( nat2, dtype = int or string )

  • coords2 (np.ndarray( (nat2, 3), dtype = float)) – atomic positions of structure 2

== optional ==

Parameters

  • lat (np.ndarray( (3, 3), dtype = float )) – matrix of lattice vectors in rows. If not None, it triggers cshda_pbc routine. Use with care.

  • thr (float) – threshold for early exit of cshda. Use with care.

== output: ==

Parameters

  • found (np.ndarray( nat2, dtype = int )) – array of assignments of atoms of structure1 to structure2, e.g.: found[3]=4 means the atom index 3 from structure1 has been assigned to atom index 4 from structure2. When nat1 != nat2, the first nat1 indices give the permutation, the rest stay fixed.

  • dists (np.ndarray( nat1, dtype = float )) – array of distances from atom index i in structure1 to atom found[i] in structure2. The Hausdorff distance can be obtained as dH = np.max( dists )

Raises

ValueError – when the requirement nat1 =< nat2 is not met.

match(nat1, typ1, coords1, nat2, typ2, coords2, kmax_factor, candidate1=None, candidate2=None)[source]

The Iterative Rotatitions and Assignments (IRA) procedure to match two structures, including the SVD correction at the end.

This function is a wrapper to the libira_match routine in library_ira.f90

It returns the solution to:

\[P_B B = \mathbf{R} A + \mathbf{t}\]

where \(A\) and \(B\) are two atomic structures, \(\mathbf{R}\) is a rotation matrix, \(\mathbf{t}\) is a translation vector, and \(P_B\) is a permutation matrix.

Solution is to be applied to struc2, as:

>>> for i in range(nat2):
>>>     coords2[i] = np.matmul( rotation, coords2[i] ) + translation

or alternatively to struc1 as:

>>> for i in range(nat1):
>>>     coords1[i] = np.matmul( rotation.T, coords1[i] ) - np.matmul( rotation.T, translation )

Apply permutation:

>>> coords2[:] = coords2[permutation]
>>> typ2[:] = typ2[permutation]

Note

Requirement: nat1 \(\le\) nat2

== input: ==

Parameters

  • nat1 (int) – number of atoms of structure 1

  • typ1 (np.ndarray( nat1, dtype = int or string )) – atomic types of structure 1

  • coords1 (np.ndarray( (nat1, 3), dtype = float )) – atomic positions of structure 1

  • nat2 (int) – number of atoms of structure 2

  • typ2 (atomic types of structure 2) – np.ndarray( nat2, dtype = int or string )

  • coords2 (np.ndarray( (nat2, 3), dtype = float)) – atomic positions of structure 2

  • kmax_factor (float) – factor for multiplication of search radius, the value should be > 1.0, for very non-congruent structures a higher value is needed. Lower value speeds up the algorithm, but can give non-optimal results.

== optional ==

Parameters

  • candidate1 (int, or np.ndarray( AnySize, dtype = int )) – list of candidate central atoms of structure1. This is useful when some additional information is known about the structures, e.g. we want to impose some specific atoms as possible central atoms. Particularly useful when matching structures with different number of atoms. If None, geometric center is used as central point. NOTE: the indices should be starting from 0. Use with care.

  • candidate2 (int, or np.ndarray( AnySize, dtype = int )) – similar as candidate1, except for structure2

Raises

  • ValueError – when kmax_factor is =< 1.0

  • ValueError – when the requirement nat1 =< nat2 is not met.

== output: ==

Parameters

  • rotation (np.ndarray( (3,3), dtype = float)) – rotation matrix

  • translation (np.ndarray( 3, dtype = float)) – translation vector

  • permutation (np.ndarray( (nat2), dtype = int )) – permutation of atoms

  • hd (float) – Hausdorff distance value of the match