# coding: utf-8 from msspec.utils import * from ase.build import bulk from ase.visualize import view import numpy as np from msspec.calculator import MSSPEC, XRaySource from msspec.iodata import Data from itertools import product DATA = None def AlN_cluster(emitter_tag, emitter_plane, diameter=0, planes=0, term_anion=True, tetra_down=True): """ This function is a kind of overload of the hemispherical_cluster function with specific attributes to the AlN structure :param emitter_tag: An integer that allows to identify the kind of atom to use as the emitter :param emitter_plane: The plane where the emitter is. 0 means the surface. :param diameter: The diameter of the cluster (in Angstroms). :param planes: The total number of planes. :param term_anion: True if the surface plane is anion terminated, False otherwise :param tetra_down: The orientation of the tetrahedral :return: """ # create the initial cluster of AlN cluster = bulk('AlN', crystalstructure='wurtzite', a=3.11, c=4.975) # tag each atom in the unit cell because they are all in a different chemical/geometrical environment # (0 and 2 for Aluminum's atoms and 1 and 3 for Nitride's atoms) [atom.set('tag', i) for i, atom in enumerate(cluster)] # change the orientation of the tetrahedron if not tetra_down: cluster.rotate(180, 'y') # From this base pattern, creat an hemispherically shaped cluster cluster = hemispherical_cluster(cluster, emitter_tag=emitter_tag, emitter_plane=emitter_plane, diameter=diameter, planes=planes) # Depending on the number of planes above the emitter, the termination may not be the one you wish, # we test this and raise an exception in such a case # Get the termination by finding the kind of one atom located at the topmost z coordinate termination = cluster[np.argsort(cluster.get_positions()[:, 2])[-1]].symbol # test if the combination of parameters is possible assert (termination == 'N' and term_anion) or (termination == 'Al' and not term_anion), ( "This termination isn't compatible with your others parameters, you must change the tag of " "your emitter, the plane of your emitter or your termination") return cluster def create_clusters(side='Al', emitter='Al', diameter=15, planes=6): clusters = [] if side == 'Al': term_anion = tetra_down = False elif side == 'N': term_anion = tetra_down = True if emitter == 'Al': tags = (0, 2) level = '2p' ke = 1407. elif emitter == 'N': tags = (1, 3) level = '1s' ke = 1083. for emitter_tag, emitter_plane in product(tags, range(0, planes)): nplanes = emitter_plane + 2 try: cluster = AlN_cluster(emitter_tag, emitter_plane, diameter=diameter, planes=nplanes, term_anion=term_anion, tetra_down=tetra_down) except AssertionError: continue cluster.absorber = get_atom_index(cluster, 0, 0, 0) cluster.info.update({'emitter_plane': emitter_plane, 'emitter_tag': emitter_tag, 'emitter': emitter, 'side': side, 'level': level, 'ke': ke}) clusters.append(cluster) return clusters def compute_polar_scans(clusters, theta=np.arange(-20., 80., 1.), phi=0., data=DATA): for ic, cluster in enumerate(clusters): # select the spectroscopy of the calculation and create a new folder for each calculation side, emitter, tag, plane, level, ke = [cluster.info[k] for k in ('side', 'emitter', 'emitter_tag', 'emitter_plane', 'level', 'ke')] calc = MSSPEC(spectroscopy='PED', folder='calc{:0>3d}_S{}_E{}_T{:d}_P{:d}'.format(ic, side, emitter, tag, plane)) calc.calculation_parameters.scattering_order = max(1, min(4, plane)) calc.source_parameters.energy = XRaySource.AL_KALPHA calc.source_parameters.theta = -35 calc.detector_parameters.angular_acceptance = 4. calc.detector_parameters.average_sampling = 'medium' calc.calculation_parameters.path_filtering = 'forward_scattering' calc.calculation_parameters.off_cone_events = 1 [a.set('forward_angle', 30.) for a in cluster] calc.calculation_parameters.vibrational_damping = 'averaged_tl' [a.set('mean_square_vibration', 0.006) for a in cluster] calc.set_atoms(cluster) data = calc.get_theta_scan(level=level, theta=theta, phi=phi, kinetic_energy=ke, data=data) dset = data[-1] dset.title = 'Polar scan {emitter:s}({level:s} tag {tag:d}) in plane #{plane:d}'.format(emitter=emitter, level=level, tag=tag, plane=plane) for name, value in cluster.info.items(): dset.add_parameter(group='AlnTuto', name=name, value=str(value), unit="") #data.save('all_polar_scans.hdf5', append=True) data.save('all_polar_scans.hdf5') def analysis(filename='all_polar_scans.hdf5'): data=Data.load(filename) # create datasets to store the sum of all emitter tmp_data = {} for dset in data: # retrieve some info side = dset.get_parameter(group='AlnTuto', name='side')['value'] emitter = dset.get_parameter(group='AlnTuto', name='emitter')['value'] try: key = '{}_{}'.format(side, emitter) tmp_data[key] += dset.cross_section except KeyError: tmp_data[key] = dset.cross_section.copy() # get the index of theta = 0; it0 = np.where(dset.theta == 0)[0][0] for key, cs in tmp_data.items(): tmp_data[key + '_norm'] = cs / cs[it0] tmp_data['Al_side'] = tmp_data['Al_Al_norm'] / tmp_data['Al_N_norm'] tmp_data['N_side'] = tmp_data['N_Al_norm'] / tmp_data['N_N_norm'] # now add all columns substrate_dset = data.add_dset('Total substrate signal') substrate_dset.add_columns(theta=dset.theta.copy()) substrate_dset.add_columns(**tmp_data) view = substrate_dset.add_view('Al side', title=r'AlN Polar scan in the (10$\bar{1}$0) azimuthal plane - Al side polarity', xlabel=r'$\Theta (\degree$)', ylabel='Signal (a.u.)') view.select('theta', 'Al_Al_norm', legend='Al 2p') view.select('theta', 'Al_N_norm', legend='N 1s') view.set_plot_options(autoscale=True) view = substrate_dset.add_view('N side', title=r'AlN Polar scan in the (10$\bar{1}$0) azimuthal plane - N side polarity', xlabel=r'$\Theta (\degree$)', ylabel='Signal (a.u.)') view.select('theta', 'N_Al_norm', legend='Al 2p') view.select('theta', 'N_N_norm', legend='N 1s') view.set_plot_options(autoscale=True) view = substrate_dset.add_view('Ratios', title=r'Al(2p)/N(1s) ratios on both polar sides of AlN in the (10$\bar{1}$0) ' r'azimuthal plane', xlabel=r'$\Theta (\degree$)', ylabel='Intenisty ratio') view.select('theta', 'Al_side', legend='Al side', where="theta >= 0 and theta <=70") view.select('theta', 'N_side', legend='N side', where="theta >= 0 and theta <=70") view.set_plot_options(autoscale=True) data.save('analysis.hdf5') data.view() DIAMETER = 10 PLANES = 4 clusters = create_clusters(side='Al', emitter='Al', diameter=DIAMETER, planes=PLANES) + \ create_clusters(side='Al', emitter='N', diameter=DIAMETER, planes=PLANES) + \ create_clusters(side='N', emitter='Al', diameter=DIAMETER, planes=PLANES) + \ create_clusters(side='N', emitter='N', diameter=DIAMETER, planes=PLANES) compute_polar_scans(clusters, phi=0.) analysis()