Cluster space information

The purpose of this example is to demonstrate how to obtain basic information about a cluster space.

Import modules

First, one needs to import the ClusterSpace class, which is used to store information regarding a given cluster space. In addition, one needs the methods get_singlet_info() and view_singlets() for extracting specific details regarding the singlets. Moreover, the ASE function bulk() will be needed for constructing structures.

from ase.build import bulk
from icet import ClusterSpace, get_singlet_info, view_singlets

Generate prototype structure

The next step is to build a prototype structure, here a bulk unit cell of rhenium. It is furthermore decided that the cluster space will be created through substitution of the Re atoms with several other elements. Also, the cutoffs for pairs, triplets and quadruplets are set to 10 Å, 7 Å, and 5 Å, respectively.

prototype = bulk('Re')
subelements = ['Re', 'Ti', 'W', 'Mo']
cutoffs = [10.0, 7.0, 5.0]

Create cluster space

The cluster space is created by initiating a ClusterSpace object and providing the prototype structure, cutoffs and list of elements defined previously as arguments. Next, the print() method is used to print all relevant information regarding the cluster space in tabular format.

cluster_space = ClusterSpace(prototype, cutoffs, subelements)
print(cluster_space)

The final call should produce the following (partial) output:

------------------------- Cluster Space -------------------------
 subelements: Ti Mo W Re
 cutoffs: 10.0 7.0 5.0
 number of orbits: 3768
-----------------------------------------------------------------
order |  radius  | multiplicity | index | orbit |    MC vector
-----------------------------------------------------------------
  1   |   0.0000 |        2     |    0  |    0  |    [0]
  1   |   0.0000 |        2     |    1  |    0  |    [1]
  1   |   0.0000 |        2     |    2  |    0  |    [2]
  2   |   1.3699 |        6     |    3  |    1  |  [0, 0]
...
  4   |   2.7094 |       12     | 3767  |  135  | [2, 2, 2, 2]
-----------------------------------------------------------------

Information regarding singlets

Additonal information regarding the singlets is extracted with help of the get_singlet_info function. Afterwards, the corresponding clusters are printed by calling view_singlets. One should note that both functions take the prototype ASE Atoms object created earlier as input argument.

print('\nSinglets:')
cluster_data = get_singlet_info(prototype)
for singlet in cluster_data:
    for key in singlet.keys():
        print(' {:22} : {}'.format(key, singlet[key]))
view_singlets(prototype)

These lines ought to yield the following result:

Singlets:
 orbit_index            : 0
 sites                  : [[0 : [ 0.  0.  0.]], [1 : [ 0.  0.  0.]]]
 multiplicity           : 2
 representative_site    : [0 : [ 0.  0.  0.]]

Source code

The complete source code is available in examples/get_cluster_space_info.py
'''
This example demonstrates how to obtain basic information about a cluster
space.
'''

# Import modules
from ase.build import bulk
from icet import ClusterSpace, get_singlet_info, view_singlets

# Create a prototype structure, decide which additional elements to populate
# it with (Re, Ti, W and Mo) and set the cutoffs for pairs (10.0 A),
# triplets (7.0 A) and quadruplets (5.0 A).
prototype = bulk('Re')
subelements = ['Re', 'Ti', 'W', 'Mo']
cutoffs = [10.0, 7.0, 5.0]

# Generate and print the cluster space.
cluster_space = ClusterSpace(prototype, cutoffs, subelements)
print(cluster_space)

# Extract and print additional information regarding the singlets.
print('\nSinglets:')
cluster_data = get_singlet_info(prototype)
for singlet in cluster_data:
    for key in singlet.keys():
        print(' {:22} : {}'.format(key, singlet[key]))
view_singlets(prototype)