\(\newcommand{\AA}{\unicode{x212B}}\)

\(\Delta\)SCF Post-Processing (1D)#

Now that the LumiWork is complete, it is time to analyze the results. This first section shows how to do that following the so-called “single effective phonon mode model” or “one-dimensional configuration coordinate diagram” (1D-CCM). We recall that, in this model, it is assumed that there exists a fictitious effective phonon mode whose eigenvectors follow exactly the ground-state to excited-state atomic relaxation, with an eigenfrequency computed with equation (48). For this analysis, we use the \(\Delta\)SCF computations shown in the previous tutorial to instantiate a DeltaSCF object:

import warnings
warnings.filterwarnings('ignore')

from abipy.lumi.deltaSCF import DeltaSCF
scf_files=[
    "../workflows_data/flow_deltaSCF/w0/t2/outdata/out_GSR.nc",
    "../workflows_data/flow_deltaSCF/w0/t3/outdata/out_GSR.nc",
    "../workflows_data/flow_deltaSCF/w0/t4/outdata/out_GSR.nc",
    "../workflows_data/flow_deltaSCF/w0/t5/outdata/out_GSR.nc",
]

results = DeltaSCF.from_four_points_file(scf_files)

# or
#
# results = DeltaSCF.from_json_file("../workflows_data/flow_deltaSCF/w0/outdata/lumi.json")
#
# or only two relaxations
#
# results = DeltaSCF.from_relax_file(["../workflows_data/flow_deltaSCF/w0/t0/outdata/out_GSR.nc",
#                                     "../workflows_data/flow_deltaSCF/w0/t1/outdata/out_GSR.nc"])

Note

Energy units are given in eV, distances are in \(\AA\).

Total energies and electronic levels#

The energies of the four relevant states are accessible with:

print(results.ag_energy,results.ag_star_energy, results.ae_energy, results.ag_star_energy)
-12934.332201652453 eV -12932.701536166687 eV -12934.277595698035 eV -12932.701536166687 eV

The electronic eigenenergies, if computed at a single k-point, typically \(\Gamma\), can be plotted with (spin up in black, spin down in red):

results.plot_eigen_energies(scf_files);
../../_images/ec28860a43c259d5fd725bb8a26a0286d42265e8fb6efccc0dd0ac142a206dc3.png

In the ground state, notice the 7 Eu\(_{4f}\) states located in the gap. In the excited state (simulated with constrained occupation, as shown with the (un)filled markers), the created 4f hole lowers the energy of an occupied 5d state, which is now located at the top of the gap. The 6 remaining occupied 4f states are pushed down in the VB. If you have computed the four band structures associated to each point, you can use results.plot_four_BandStructures(nscf_files) where nscf_files is the list of four band structures .nc files.

nscf_files = [
    "../workflows_data/flow_deltaSCF/w0/t6/outdata/out_GSR.nc",
    "../workflows_data/flow_deltaSCF/w0/t7/outdata/out_GSR.nc",
    "../workflows_data/flow_deltaSCF/w0/t8/outdata/out_GSR.nc",
    "../workflows_data/flow_deltaSCF/w0/t9/outdata/out_GSR.nc",
]

results.plot_four_BandStructures(nscf_files);
../../_images/27f2b5089c8ca3c93f4eac62d093880c6bb76cda847ad08a7e0fb49f783c0ae2.png

Notice the strong dispersion of the bands close the CB bottom due to an interaction between Eu\(_{5d}\) replica. Increasing the supercell size reduces this dispersion.

Atomic relaxation#

The ground/excited states structures are accessible with:

results.structure_gs()
#results.structures_ex()
Structure Summary
Lattice
    abc : 8.004850964704389 8.004850964704389 6.399859775548332
 angles : 90.0 90.0 90.0
 volume : 410.08790413783584
      A : np.float64(8.004850964704389) np.float64(0.0) np.float64(0.0)
      B : np.float64(0.0) np.float64(8.004850964704389) np.float64(0.0)
      C : np.float64(0.0) np.float64(0.0) np.float64(6.399859775548332)
    pbc : True True True
PeriodicSite: Eu (4.002, 0.0, 1.6) [0.5, 0.0, 0.25]
PeriodicSite: Sr (4.002, 0.0, 4.8) [0.5, 0.0, 0.75]
PeriodicSite: Sr (0.0, 4.002, 0.0006355) [0.0, 0.5, 9.93e-05]
PeriodicSite: Sr (0.0, 4.002, 3.199) [0.0, 0.5, 0.4999]
PeriodicSite: Li (1.006, 6.616, 1.6) [0.1257, 0.8265, 0.25]
PeriodicSite: Li (1.003, 6.614, 4.8) [0.1253, 0.8263, 0.75]
PeriodicSite: Li (6.999, 1.389, 1.6) [0.8743, 0.1735, 0.25]
PeriodicSite: Li (7.001, 1.391, 4.8) [0.8747, 0.1737, 0.75]
PeriodicSite: Li (1.389, 1.006, 0.0007597) [0.1735, 0.1257, 0.0001187]
PeriodicSite: Li (1.389, 1.006, 3.199) [0.1735, 0.1257, 0.4999]
PeriodicSite: Li (6.616, 6.999, 0.0007597) [0.8265, 0.8743, 0.0001187]
PeriodicSite: Li (6.616, 6.999, 3.199) [0.8265, 0.8743, 0.4999]
PeriodicSite: Al (2.808, 5.463, 0.0008531) [0.3508, 0.6824, 0.0001333]
PeriodicSite: Al (2.808, 5.463, 3.199) [0.3508, 0.6824, 0.4999]
PeriodicSite: Al (5.196, 2.542, 0.0008531) [0.6492, 0.3176, 0.0001333]
PeriodicSite: Al (5.196, 2.542, 3.199) [0.6492, 0.3176, 0.4999]
PeriodicSite: Al (2.544, 2.806, 1.6) [0.3178, 0.3505, 0.25]
PeriodicSite: Al (2.544, 2.808, 4.8) [0.3178, 0.3508, 0.75]
PeriodicSite: Al (5.461, 5.199, 1.6) [0.6822, 0.6495, 0.25]
PeriodicSite: Al (5.461, 5.197, 4.8) [0.6822, 0.6492, 0.75]
PeriodicSite: N (4.669, 5.856, -0.001125) [0.5833, 0.7316, -0.0001758]
PeriodicSite: N (4.669, 5.856, 3.201) [0.5833, 0.7316, 0.5002]
PeriodicSite: N (3.336, 2.149, -0.001125) [0.4167, 0.2684, -0.0001758]
PeriodicSite: N (3.336, 2.149, 3.201) [0.4167, 0.2684, 0.5002]
PeriodicSite: N (2.155, 4.669, 1.6) [0.2692, 0.5832, 0.25]
PeriodicSite: N (2.156, 4.67, 4.8) [0.2693, 0.5835, 0.75]
PeriodicSite: N (5.85, 3.336, 1.6) [0.7308, 0.4168, 0.25]
PeriodicSite: N (5.849, 3.334, 4.8) [0.7307, 0.4165, 0.75]
PeriodicSite: O (1.993, 7.104, -0.001526) [0.249, 0.8874, -0.0002385]
PeriodicSite: O (1.993, 7.104, 3.201) [0.249, 0.8874, 0.5002]
PeriodicSite: O (6.011, 0.9013, -0.001526) [0.751, 0.1126, -0.0002385]
PeriodicSite: O (6.011, 0.9013, 3.201) [0.751, 0.1126, 0.5002]
PeriodicSite: O (0.8989, 2.0, 1.6) [0.1123, 0.2499, 0.25]
PeriodicSite: O (0.8984, 2.003, 4.8) [0.1122, 0.2503, 0.75]
PeriodicSite: O (7.106, 6.004, 1.6) [0.8877, 0.7501, 0.25]
PeriodicSite: O (7.106, 6.002, 4.8) [0.8878, 0.7497, 0.75]

It is sometimes interesting to decompose the gs-ex displacements per specie:

results.get_dataframe_species()
symbol mass $\Delta R^2$ $\Delta Q^2$
0 Sr 87.620000 2.492076e-07 0.000022
1 Li 6.941000 1.707388e-03 0.011851
2 Eu 151.964000 0.000000e+00 0.000000
3 Al 26.981539 2.010499e-03 0.054246
4 N 14.006700 8.944237e-03 0.125279
5 O 15.999400 5.898415e-03 0.094371

or per atom:

results.get_dataframe_atoms(defect_symbol="Eu")
symbol mass $\Delta R$ $\Delta Q^2$ $\Delta F$ dist. from defect
0 Eu 151.964000 0.000000 0.000000 3.021136e-18 0.000000
1 Sr 87.620000 0.000000 0.000000 3.021136e-18 3.199930
2 Sr 87.620000 0.000353 0.000011 2.404364e-02 5.881894
3 Sr 87.620000 0.000353 0.000011 2.404364e-02 5.881894
4 Li 6.941000 0.022438 0.003495 3.219426e-02 3.302428
5 Li 6.941000 0.006908 0.000331 1.195830e-02 4.600811
6 Li 6.941000 0.022438 0.003495 3.219426e-02 3.302428
7 Li 6.941000 0.006908 0.000331 1.195830e-02 4.600811
8 Li 6.941000 0.012298 0.001050 3.689746e-02 3.224950
9 Li 6.941000 0.012298 0.001050 3.689746e-02 3.224950
10 Li 6.941000 0.012298 0.001050 3.689746e-02 3.224950
11 Li 6.941000 0.012298 0.001050 3.689746e-02 3.224950
12 Al 26.981539 0.011001 0.003266 6.143308e-02 3.231813
13 Al 26.981539 0.011001 0.003266 6.143308e-02 3.231813
14 Al 26.981539 0.011001 0.003266 6.143308e-02 3.231813
15 Al 26.981539 0.011001 0.003266 6.143308e-02 3.231813
16 Al 26.981539 0.010330 0.002879 7.510210e-02 3.162003
17 Al 26.981539 0.025622 0.017713 1.495512e-01 4.500308
18 Al 26.981539 0.010330 0.002879 7.510210e-02 3.162003
19 Al 26.981539 0.025622 0.017713 1.495512e-01 4.500308
20 N 14.006700 0.038935 0.021233 4.005216e-01 2.761134
21 N 14.006700 0.038935 0.021233 4.005216e-01 2.761134
22 N 14.006700 0.038935 0.021233 4.005216e-01 2.761134
23 N 14.006700 0.038935 0.021233 4.005216e-01 2.761134
24 N 14.006700 0.028770 0.011594 1.530554e-01 3.813568
25 N 14.006700 0.024751 0.008580 1.248930e-02 4.976787
26 N 14.006700 0.028770 0.011594 1.530554e-01 3.813568
27 N 14.006700 0.024751 0.008580 1.248930e-02 4.976787
28 O 15.999400 0.037741 0.022789 3.681187e-01 2.722718
29 O 15.999400 0.037741 0.022789 3.681187e-01 2.722718
30 O 15.999400 0.037741 0.022789 3.681187e-01 2.722718
31 O 15.999400 0.037741 0.022789 3.681187e-01 2.722718
32 O 15.999400 0.002985 0.000143 1.467203e-01 3.692334
33 O 15.999400 0.009567 0.001464 9.676330e-02 4.887502
34 O 15.999400 0.002985 0.000143 1.467203e-01 3.692334
35 O 15.999400 0.009567 0.001464 9.676330e-02 4.887502

You can plot these displacements (or ground state forces at the excited state atomic positions) as a function of the distance with respect to the defect. This allows to check the convergence of your calculation with respect to the supercell size. We see in our toy example that the results are not converged: the displacements are underestimated because o f cancelation error due to periodic defect replica. Note that the forces decay faster with respect to the distance (this will be important in later tutorials).

results.plot_delta_R_distance(defect_symbol="Eu");
../../_images/8154c14532e6523e436abb6c400850ed6bf7636990045cf5992d459f5a90c1a4.png
results.plot_delta_F_distance(defect_symbol="Eu");
../../_images/4fade43e7aa93653dcf269fb19fec8bd6e40a25e605d044ea2950ee89a745c62.png

In order to visualize these displacements on a VESTA structure, you can follow these steps.

(1) create a cif file with the ground state structure. (2) Open the structure with VESTA software and save it in a .vesta format (this should be done manually). (3) Use draw_displacements_vesta() method.

results.structure_gs().to(filename="gs_stru.cif")
# then open this cif with vesta, save it as .vesta file format
"# generated using pymatgen\ndata_Sr3Li8EuAl8(NO)8\n_symmetry_space_group_name_H-M   'P 1'\n_cell_length_a   8.00485096\n_cell_length_b   8.00485096\n_cell_length_c   6.39985978\n_cell_angle_alpha   90.00000000\n_cell_angle_beta   90.00000000\n_cell_angle_gamma   90.00000000\n_symmetry_Int_Tables_number   1\n_chemical_formula_structural   Sr3Li8EuAl8(NO)8\n_chemical_formula_sum   'Sr3 Li8 Eu1 Al8 N8 O8'\n_cell_volume   410.08790414\n_cell_formula_units_Z   1\nloop_\n _symmetry_equiv_pos_site_id\n _symmetry_equiv_pos_as_xyz\n  1  'x, y, z'\nloop_\n _atom_site_type_symbol\n _atom_site_label\n _atom_site_symmetry_multiplicity\n _atom_site_fract_x\n _atom_site_fract_y\n _atom_site_fract_z\n _atom_site_occupancy\n  Eu  Eu0  1  0.50000000  0.00000000  0.25000000  1\n  Sr  Sr1  1  0.50000000  0.00000000  0.75000000  1\n  Sr  Sr2  1  0.00000000  0.50000000  0.00009930  1\n  Sr  Sr3  1  0.00000000  0.50000000  0.49990070  1\n  Li  Li4  1  0.12568471  0.82654005  0.25000000  1\n  Li  Li5  1  0.12534406  0.82629134  0.75000000  1\n  Li  Li6  1  0.87431529  0.17345995  0.25000000  1\n  Li  Li7  1  0.87465594  0.17370866  0.75000000  1\n  Li  Li8  1  0.17351474  0.12571140  0.00011871  1\n  Li  Li9  1  0.17351474  0.12571140  0.49988129  1\n  Li  Li10  1  0.82648526  0.87428860  0.00011871  1\n  Li  Li11  1  0.82648526  0.87428860  0.49988129  1\n  Al  Al12  1  0.35084176  0.68244058  0.00013330  1\n  Al  Al13  1  0.35084176  0.68244058  0.49986670  1\n  Al  Al14  1  0.64915824  0.31755942  0.00013330  1\n  Al  Al15  1  0.64915824  0.31755942  0.49986670  1\n  Al  Al16  1  0.31784161  0.35050234  0.25000000  1\n  Al  Al17  1  0.31781068  0.35081909  0.75000000  1\n  Al  Al18  1  0.68215839  0.64949766  0.25000000  1\n  Al  Al19  1  0.68218932  0.64918091  0.75000000  1\n  N  N20  1  0.58325401  0.73159491  -0.00017577  1\n  N  N21  1  0.58325401  0.73159491  0.50017577  1\n  N  N22  1  0.41674599  0.26840509  -0.00017577  1\n  N  N23  1  0.41674599  0.26840509  0.50017577  1\n  N  N24  1  0.26919588  0.58323483  0.25000000  1\n  N  N25  1  0.26928648  0.58345522  0.75000000  1\n  N  N26  1  0.73080412  0.41676517  0.25000000  1\n  N  N27  1  0.73071352  0.41654478  0.75000000  1\n  O  O28  1  0.24902582  0.88740903  -0.00023850  1\n  O  O29  1  0.24902582  0.88740903  0.50023850  1\n  O  O30  1  0.75097418  0.11259097  -0.00023850  1\n  O  O31  1  0.75097418  0.11259097  0.50023850  1\n  O  O32  1  0.11229476  0.24989457  0.25000000  1\n  O  O33  1  0.11223134  0.25025797  0.75000000  1\n  O  O34  1  0.88770524  0.75010543  0.25000000  1\n  O  O35  1  0.88776866  0.74974203  0.75000000  1\n"
results.draw_displacements_vesta(in_path="gs_stru.vesta",color_vector=[0, 0, 0])
Vesta files created and stored in : 
 /home/runner/work/abipy_book/abipy_book/abipy_book/lumabi/post_process_1D/VESTA_FILES

The resulting vesta file should look like:

../../_images/draw_displacements_vesta.png

You can modify the vectors drawing by changing the default arguments

help(results.draw_displacements_vesta)
Help on method draw_displacements_vesta in module abipy.lumi.deltaSCF:

draw_displacements_vesta(in_path, mass_weighted=False, scale_vector=20, width_vector=0.3, color_vector=(255, 0, 0), centered=True, factor_keep_vectors=0.1, out_path='VESTA_FILES', out_filename='gs_ex_relaxation') method of abipy.lumi.deltaSCF.DeltaSCF instance
    Draw the ground state to excited state atomic relaxation on a vesta structure.

    Args:
        in_path : path where the initial .vesta structure in stored, should correspond to the ground state relaxed structure.
        mass_weighted : If True, weight the displacements by the atomic masses. Draw the \Delta Q in that case.
        scale_vector : scaling factor of the vector modulus
        width_vector : vector width
        color_vector : color in rgb format
        centered : center the vesta structure around [0,0,0]
        factor_keep_vectors : draw only the eigenvectors with magnitude > factor_keep_vectors * max(magnitude)
        out_path : path where .vesta files with vector are stored

Luminescent properties following the 1D-CCM#

One can visualize the 1D-CCM and associated displaced parabolas with

results.draw_displaced_parabolas();
../../_images/ffec15046780f321efb50b5f0a8c687fb849d1c5b83efd7a93d75400e3a9a32b.png

or get a dataframe (or dictionary) with the main 1D-CCM parameters.

results.get_dataframe()
#results.get_dict_results()
E_em E_abs E_zpl E_fc_gs E_fc_ex Delta_S Delta_R Delta_Q Eff_freq_gs Eff_freq_ex S_em S_abs
None 1.507326 1.630665 1.561932 0.054606 0.068734 0.12334 0.136238 0.534574 0.039969 0.044842 1.366207 1.532787

Finally, one can plot the luminescence lineshape at 0 K:

results.plot_lineshape_1D_zero_temp(energy_range=[1,2]);
help(results.plot_lineshape_1D_zero_temp)
../../_images/c79d23c2f35daa6df3633e223005dbd0ff08a429eb0078f88a01d03915598f13.png
Help on method plot_lineshape_1D_zero_temp in module abipy.lumi.deltaSCF:

plot_lineshape_1D_zero_temp(energy_range=[0.5, 5], max_m=25, phonon_width=0.01, with_omega_cube='True', normalized='Area', ax=None, **kwargs) -> 'Figure' method of abipy.lumi.deltaSCF.DeltaSCF instance
    Plot the the emission lineshape following the effective phonon 1D-CCM at T=0K.
    NOT based on the generating function.

    Args:
        ax: |matplotlib-Axes| or None if a new figure should be created.
        energy_range:  Energy range at which the intensities are computed, ex : [0.5,5]
        max_m: Maximal vibrational state m considered
        phonon_width: fwhm of each phonon peak, in eV
        with_omega_cube: Considered or not the omega^3 dependence of the intensity
        normlized: Normalisation procedure. 'Area' if Area under the curve = 1. 'Sum' if maximum of the curve = 1.

    Returns: |matplotlib-Figure|




    Keyword arguments controlling the display of the figure:

    ================  ====================================================
    kwargs            Meaning
    ================  ====================================================
    title             Title of the plot. Default: None.
    show              True to show the figure. Default: True.
    savefig           "abc.png" or "abc.svg" to save the figure to a file.
    size_kwargs       Dictionary with options passed to fig.set_size_inches
                      e.g. size_kwargs=dict(w=3, h=4).
    tight_layout      True to call fig.tight_layout. Default: False.
    ax_grid           True (False) to add (remove) grid from all axes in fig.
                      Default: None i.e. fig is left unchanged.
    ax_annotate       Add labels to subplots e.g. (a), (b). Default: False
    fig_close         Close figure. Default: False.
    plotly            Try to convert mpl figure to plotly: Default: False
    ================  ====================================================