Return the gyration tensor of a molecule.

    The gyration tensor should be invariant to the molecule's position.
    The known formulas for the gyration tensor have the correction for the
    centre of mass of the molecule, therefore, the coordinates are first
    corrected for the centre of mass and essentially shifted to the origin.

    Parameters
    ----------
    elements : numpy.ndarray
        The array containing the molecule's elemental data.

    coordinates :  numpy.ndarray
        The array containing the Cartesian coordinates of the molecule.

    Returns
    -------
    numpy.ndarray
        The gyration tensor of a molecule invariant to the molecule's position.

def get_gyration_tensor(elements, coordinates):
    """
    Return the gyration tensor of a molecule.

    The gyration tensor should be invariant to the molecule's position.
    The known formulas for the gyration tensor have the correction for the
    centre of mass of the molecule, therefore, the coordinates are first
    corrected for the centre of mass and essentially shifted to the origin.

    Parameters
    ----------
    elements : numpy.ndarray
        The array containing the molecule's elemental data.

    coordinates :  numpy.ndarray
        The array containing the Cartesian coordinates of the molecule.

    Returns
    -------
    numpy.ndarray
        The gyration tensor of a molecule invariant to the molecule's position.

    """
    # First calculate COM for correction.
    com = centre_of_mass(elements, coordinates)
    # Correct the coordinates for the COM.
    coordinates = coordinates - com
    # Calculate diagonal and then other values of the matrix.
    diag = np.sum(coordinates**2, axis=0)
    xy = np.sum(coordinates[:, 0] * coordinates[:, 1])
    xz = np.sum(coordinates[:, 0] * coordinates[:, 2])
    yz = np.sum(coordinates[:, 1] * coordinates[:, 2])
    S = np.array([[diag[0], xy, xz], [xy, diag[1], yz],
                  [xz, yz, diag[2]]]) / coordinates.shape[0]
    return (S)