VlachosGroup / pMuTT

Python Multiscale Thermochemistry Toolbox (pMuTT)
https://vlachosgroup.github.io/pMuTT/
40 stars 23 forks source link

help wanted set_vib_wavenumbers_from_outcar #169

Closed MCVifpen closed 4 years ago

MCVifpen commented 4 years ago

Hello,

I need some help with the set_vib_wavenumbers_from_outcar method to parse VASP frequencies from an OUTCAR file to a dictionary. I need a way to get automatically the right frequencies for a thermochemistry calculation.

VASP yields 3N frequencies, though for thermochemistry calculations you only need 3N-6 or 3N-5 eigenvalues, i.e.: it computes redundant values from translations and rotations. For example, if I apply the set_vib_wavenumbers_from_outcar method to a H2 frequency calculation (3N-5 = 1 expected frequency) I get

{'vib_wavenumbers': [4300.073749, 64.352654, 63.163151, 4300.073749, 64.352654, 63.163151]}

If I apply what seems to me a sensible cutoff to remove unwanted frequencies, such as set_vib_wavenumbers_from_outcar("OUTCAR",,min_frequency_cutoff=100) I still get {'vib_wavenumbers': [4300.073749, 4300.073749]} that is , two frequencies instead of one.

What would you to get one single frequency and in a robust way that you can apply to any molecule ?

Thank you for your help.

Kind regards, Manuel

jonlym commented 4 years ago

This is an interesting question. I would advise to always be careful when eliminating frequencies as symmetry may result in duplicate frequencies that are physical (e.g. methane). That being said, the example you provided is legitimate.

You can modify your code to check if the molecular is linear/non-linear (perhaps by using pmutt.statmech.rot.get_geometry_from_atoms) and if the number of frequencies does not match your array, filter them further using a function like the following:

def get_unique_float(in_data, tol=1.e-6):
    """Returns unique values. Solution from https://stackoverflow.com/a/5429049/6865095

    Parameters
    ----------
        in_data : numpy.ndarray
            Array with duplicate values
        tol : float, optional
            Tolerance to use to determine if values are unique. Default value
            is 1.e-6
    Returns
    -------
        out_data : numpy.ndarray
            Array without duplicate values
    """
    b = in_data.copy()
    b.sort()
    d = np.append(True, np.diff(b))
    out_data = b[d>tol]
    return out_data

When I fed the function the frequencies you posted, it returned the following:

[  63.163151   64.352654 4300.073749]

Let me know if you think that solves your problem.

MCVifpen commented 4 years ago

Hello,

Thank you very much for the hint. Based on your answer, I’ve being working on the below piece of code. I’m afraid that something is wrong in the get_unique_float function (see error message below), but I can’t figure out how to fix it.

I would appreciate if you could give me a hand, all the more as it has to do more with my python skills than with pmutt this time.

Manuel

Entrée [13]:

N=len(H2_molecule.get_atomic_numbers()) if N>1 and get_geometry_from_atoms(H2_molecule)=='linear': int_dof=3N-5 elif N>1 and get_geometry_from_atoms(H2_molecule)=='nonlinear': int_dof=3N-6 else: int_dof=0

‘here I compute the number of expected vibrations’ ‘then I check against the number of VASP vibrations and try to clean the array’

if len(H2_vib['vib_wavenumbers'])>int_dof: get_unique_float(H2_vib['vib_wavenumbers'], tol=1.e-6)

TypeError Traceback (most recent call last)

in 7 int_dof=0 8 if len(H2_vib['vib_wavenumbers'])>int_dof: ----> 9 get_unique_float(H2_vib['vib_wavenumbers'], tol=1.e-6) in get_unique_float(in_data, tol) 18 b.sort() 19 d = np.append(True, np.diff(b)) ---> 20 out_data = b[d>tol] 21 return out_data TypeError: only integer scalar arrays can be converted to a scalar index De : Jonathan Lym [mailto:notifications@github.com] Envoyé : mardi 7 juillet 2020 15:37 À : VlachosGroup/pMuTT Cc : CORRAL VALERO Manuel; Author Objet : Re: [VlachosGroup/pMuTT] help wanted set_vib_wavenumbers_from_outcar (#169) This is an interesting question. I would advise to always be careful when eliminating frequencies as symmetry may result in duplicate frequencies that are physical (e.g. methane). That being said, the example you provided is legitimate. You can modify your code to check if the molecular is linear/non-linear (perhaps by using pmutt.statmech.rot.get_geometry_from_atoms) and if the number of frequencies does not match your array, filter them further using a function like the following: def get_unique_float(in_data, tol=1.e-6): """Returns unique values. Solution from https://stackoverflow.com/a/5429049/6865095 using-a-d Parameters ---------- in_data : numpy.ndarray Array with duplicate values tol : float, optional Tolerance to use to determine if values are unique. Default value is 1.e-6 Returns ------- out_data : numpy.ndarray Array without duplicate values """ b = in_data.copy() b.sort() d = np.append(True, np.diff(b)) out_data = b[d>tol] return out_data When I fed the function the frequencies you posted, it returned the following: [ 63.163151 64.352654 4300.073749] Let me know if you think that solves your problem. — You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub, or unsubscribe. __________________________ Avant d'imprimer, pensez à l'environnement ! Please consider the environment before printing ! Ce message et toutes ses pièces jointes sont confidentiels et établis à l'intention exclusive de ses destinataires. Toute utilisation non conforme à sa destination, toute diffusion ou toute publication, totale ou partielle, est interdite, sauf autorisation expresse. IFP Energies nouvelles décline toute responsabilité au titre de ce message. This message and any attachments are confidential and intended solely for the addressees. Any unauthorised use or dissemination is prohibited. IFP Energies nouvelles should not be liable for this message. __________________________
jonlym commented 4 years ago

My guess is H2_vib['vib_wavenumbers'] is not a NumPy array. Please correct me if I'm wrong. However, the input to the function needs to be a NumPy array in order to work.

jonlym commented 4 years ago

Hi Manuel, Did this solve your problem? Let me know so I can help further or close the issue. Thank you again for your interest in pMuTT!

MCVifpen commented 4 years ago

Dear Jonathan,

I tried this morning and it solved my problem, indeed.

I’ve parsed the H2_vib['vib_wavenumbers'] dictionnary into the get_unique_float function as get_unique_float(np.array(H2_vib['vib_wavenumbers']), tol=1.e-6) and it worked fine.

Thank you again ! Manuel

De : Jonathan Lym [mailto:notifications@github.com] Envoyé : mercredi 8 juillet 2020 16:53 À : VlachosGroup/pMuTT Cc : CORRAL VALERO Manuel; Author Objet : Re: [VlachosGroup/pMuTT] help wanted set_vib_wavenumbers_from_outcar (#169)

Hi Manuel, Did this solve your problem? Let me know so I can help further or close the issue. Thank you again for your interest in pMuTT!

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHubhttps://github.com/VlachosGroup/pMuTT/issues/169#issuecomment-655569167, or unsubscribehttps://github.com/notifications/unsubscribe-auth/AQDTCESFB4UYU76VBT772YTR2SB3TANCNFSM4OSSK3UA.


Avant d'imprimer, pensez à l'environnement ! Please consider the environment before printing ! Ce message et toutes ses pièces jointes sont confidentiels et établis à l'intention exclusive de ses destinataires. Toute utilisation non conforme à sa destination, toute diffusion ou toute publication, totale ou partielle, est interdite, sauf autorisation expresse. IFP Energies nouvelles décline toute responsabilité au titre de ce message. This message and any attachments are confidential and intended solely for the addressees. Any unauthorised use or dissemination is prohibited. IFP Energies nouvelles should not be liable for this message.


jonlym commented 4 years ago

Great! I'll close the issue then.