78 lines
3.1 KiB
Python
78 lines
3.1 KiB
Python
|
# coding: utf-8
|
||
|
|
|
||
|
|
# import all we need and start by msspec
|
||
|
|
from msspec.calculator import MSSPEC
|
||
|
|
|
||
|
|
# we will build a simple atomic chain
|
||
|
|
from ase import Atom, Atoms
|
||
|
|
|
||
|
|
# we need some numpy functions
|
||
|
|
import numpy as np
|
||
|
|
|
||
|
|
|
||
|
|
symbol = 'Ni' # The kind of atom for the chain
|
||
|
|
orders = (1, 5) # We will run the calculation for single scattering
|
||
|
|
# and multiple scattering (5th diffusion order)
|
||
|
|
chain_lengths = (2,3,5) # We will run the calculation for differnt lengths
|
||
|
|
# of the atomic chain
|
||
|
|
a = 3.5 # The distance bewteen 2 atoms
|
||
|
|
|
||
|
|
# Define an empty variable to store all the results
|
||
|
|
all_data = None
|
||
|
|
|
||
|
|
# 2 for nested loops over the chain length and the order of diffusion
|
||
|
|
for chain_length in chain_lengths:
|
||
|
|
for order in orders:
|
||
|
|
# We build the atomic chain by
|
||
|
|
# 1- stacking each atom one by one along the z axis
|
||
|
|
chain = Atoms([Atom(symbol, position = (0., 0., i*a)) for i in
|
||
|
|
range(chain_length)])
|
||
|
|
# 2- rotating the chain by 45 degrees with respect to the y axis
|
||
|
|
chain.rotate('y', np.radians(45.))
|
||
|
|
# 3- setting a custom Muffin-tin radius of 1.5 angstroms for all
|
||
|
|
# atoms (needed if you want to enlarge the distance between
|
||
|
|
# the atoms while keeping the radius constant)
|
||
|
|
[atom.set('mt_radius', 1.5) for atom in chain]
|
||
|
|
# 4- defining the absorber to be the first atom in the chain at
|
||
|
|
# x = y = z = 0
|
||
|
|
chain.absorber = 0
|
||
|
|
|
||
|
|
# We define a new PED calculator
|
||
|
|
calc = MSSPEC(spectroscopy = 'PED')
|
||
|
|
calc.set_atoms(chain)
|
||
|
|
# Here is how to tweak the scattering order
|
||
|
|
calc.calculation_parameters.scattering_order = order
|
||
|
|
# This line below is where we actually run the calculation
|
||
|
|
all_data = calc.get_theta_scan(level='3s', kinetic_energy=1000.,
|
||
|
|
theta=np.arange(0., 80.), data=all_data)
|
||
|
|
|
||
|
|
# OPTIONAL, to improve the display of the data we will change the dataset
|
||
|
|
# default title as well as the plot title
|
||
|
|
t = "order {:d}, n = {:d}".format(order, chain_length) # A useful title
|
||
|
|
dset = all_data[-1] # get the last dataset
|
||
|
|
dset.title = t # change its title
|
||
|
|
# get its last view (there is only one defined for each dataset)
|
||
|
|
v = dset.views()[-1]
|
||
|
|
v.set_plot_options(title=t) # change the title of the figure
|
||
|
|
|
||
|
|
|
||
|
|
|
||
|
|
# OPTIONAL, set the same scale for all plots
|
||
|
|
# 1. iterate over all computed cross_sections to find the absolute minimum and
|
||
|
|
# maximum of the data
|
||
|
|
min_cs = max_cs = 0
|
||
|
|
for dset in all_data:
|
||
|
|
min_cs = min(min_cs, np.min(dset.cross_section))
|
||
|
|
max_cs = max(max_cs, np.max(dset.cross_section))
|
||
|
|
|
||
|
|
# 2. for each view in each dataset, change the scale accordingly
|
||
|
|
for dset in all_data:
|
||
|
|
v = dset.views()[-1]
|
||
|
|
v.set_plot_options(ylim=[min_cs, max_cs])
|
||
|
|
|
||
|
|
# Pop up the graphical window
|
||
|
|
all_data.view()
|
||
|
|
# You can end your script with the line below to remove the temporary
|
||
|
|
# folder needed for the calculation
|
||
|
|
calc.shutdown()
|