PDOS for molecules

[Mohammad Maymoun]

Greetings,

Hope are you wells,

I am trying to construct a PDOS for a molecule systeme, using the in-built functionalities in ASE, using KOOPMANS code, and I am generated a script.py as you see below. However, I did not find the PDOS when I ran it.

Can you please appreciate your assistances in resolving this issue?

script.py:
import pickle
import numpy as np
import matplotlib.pyplot as plt
from koopmans import io
from ase_koopmans.dft.dos import DOS
from scipy.ndimage import gaussian_filter1d

# === Load Koopmans workflow ===
wf = io.read("aa.pkl")

# === Access the final KI calculation ===
final_calc = wf.calculations[-1]
atoms = final_calc.atoms

# === Construct DOS object ===
dos = DOS(final_calc)

# === Extract projected DOS ===
energies, pdos = dos.get_projected_dos()  # pdos shape: (n_atoms, n_energies)

# === Group PDOS by element symbol ===
atom_indices_by_symbol = {}
for i, atom in enumerate(atoms):
    atom_indices_by_symbol.setdefault(atom.symbol, []).append(i)

# === Plot PDOS with optional Gaussian smoothing ===
plt.figure(figsize=(8, 5))

sigma = 1.0  # Gaussian broadening in energy points (adjust if needed)

for symbol, indices in atom_indices_by_symbol.items():
    pdos_sum = np.sum([pdos[i] for i in indices], axis=0)
    pdos_smoothed = gaussian_filter1d(pdos_sum, sigma=sigma)
    plt.plot(energies, pdos_smoothed, label=f"{symbol} atoms")

# === Customize plot ===
plt.xlabel("Energy (eV)")
plt.ylabel("Projected DOS (a.u.)")
plt.title("Projected DOS (Koopmans KI) - aa")
plt.axhline(0, color="black", linewidth=0.5)
plt.grid(True)
plt.legend()
plt.tight_layout()
plt.show()

[Edward Linscott]

Dear Mohammad Maymoun,

plt.show() just shows the figure on your screen (but does not save it); plt.savefig(...) is the function you want to use to save it to file

Best wishes,

Edward

[Mohammad Maymoun]

Dear Edward,

Thank you for your response. I replaced plt.show() with plt.savefig() and ran plot.py. Unfortunately, I received this error: AttributeError: 'DOS' object has no attribute 'get_projected_dos'.

Some editing must be done to calculate the PDOS for a molecular system in the input file (ozone.json)?.
Could you please share an example of how to construct PDOS using ASE with Koopmans?

Here is the plot.py script I use:

import pickle
import numpy as np
import matplotlib.pyplot as plt
from koopmans import io
from ase_koopmans.dft.dos import DOS
from scipy.ndimage import gaussian_filter1d

# === Load Koopmans workflow ===

wf = io.read("ozone.pkl")

# === Access the final KI calculation ===
final_calc = wf.calculations[-1]
atoms = final_calc.atoms

# === Construct DOS object ===
dos = DOS(final_calc)

# === Extract projected DOS ===

energies = dos.energies
pdos = dos.get_projected_dos()  

# === Group PDOS by element symbol ===
atom_indices_by_symbol = {}
for i, atom in enumerate(atoms):
    atom_indices_by_symbol.setdefault(atom.symbol, []).append(i)

# === Plot PDOS with optional Gaussian smoothing ===
plt.figure(figsize=(8, 5))

sigma = 1.0  # Gaussian broadening in energy points (adjust if needed)

for symbol, indices in atom_indices_by_symbol.items():
    pdos_sum = np.sum([pdos[i] for i in indices], axis=0)
    pdos_smoothed = gaussian_filter1d(pdos_sum, sigma=sigma)
    plt.plot(energies, pdos_smoothed, label=f"{symbol} atoms")

# === Customize plot ===
plt.xlabel("Energy (eV)")
plt.ylabel("Projected DOS (a.u.)")

plt.title("Projected DOS (Koopmans KI) - ozone")

plt.axhline(0, color="black", linewidth=0.5)
plt.grid(True)
plt.legend()
plt.tight_layout()

plt.savefig("ozone_koopmans_pdos.png", dpi=300)

--
You received this message because you are subscribed to the Google Groups "koopmans-users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to koopmans-user...@googlegroups.com.
To view this discussion, visit https://groups.google.com/d/msgid/koopmans-users/0009892e-bfc8-4a5d-83f3-2f5ecc2a32d2n%40googlegroups.com.

--

Best regards,

Mohammad Maymoun

******************************

PhD in Materials Physics and Energy

[Edward Linscott]

Dear Mohammad,

I am not sure where you got this script — as per the error, get_projected_dos is not a method of the DOS class.

As per an earlier message of mine, you can plot the total DOS with

# Plot and save the DOS
plt.plot(dos.get_energies(), dos.get_dos())
plt.savefig('dos.png')

Plotting a pDOS is tricker. Unlike the total DOS, it is not possible to reconstruct the pDOS from the orbital energies alone, so generating a pDOS is not as straightforward. In fact, I don't think this is currently implemented in `kcp.x`. I have opened an issue on GitHub and will look into implementing it.

Best wishes,

Edward

[Mohammad Maymoun]

Dear Edward,

Thank you very much for your feedback. In fact, I generated this script myself. My task is to compute  the partial density of states (PDOS) from Koopmans' code for molecular systems (for instance, the Kohn-Sham functional or KI ).

Best regards 

Mohammad