updated_top.py

import foyer
import mbuild as mb
import unyt as u
import gmso
from gmso.lib.potential_templates import PotentialTemplateLibrary
from gmso.external import from_parmed
from gmso.formats.top import write_top
from gmso.formats.gro import write_gro

def gmso_top():
    ethane = mb.load('CC', smiles=True)
    ff = foyer.forcefields.load_OPLSAA()
    pmd_ethane = ff.apply(ethane)
    typed_ethane = from_parmed(pmd_ethane)
    typed_ethane.name = "ethane"
    
    per_torsion = PotentialTemplateLibrary()["RyckaertBellemansTorsionPotential"]
    params = {"c0": 0.627 * u.Unit("kJ / mol"),
              "c1": 1.8828 * u.Unit("kJ / mol"),
              "c2": 0.0 * u.Unit("kJ / mol"),
              "c3": -2.5104 * u.Unit("kJ / mol"),
              "c4": 0.0 * u.Unit("kJ / mol"),
              "c5": 0.0 * u.Unit("kJ / mol"),
              }

    periodic_dihedral_type = gmso.core.potential.Potential.from_template(
            per_torsion, params)
    for dihedral in typed_ethane.dihedrals:
        dihedral.connection_type = periodic_dihedral_type
    
    typed_ethane.update_connections()
    
    write_top(typed_ethane, 'gmso.top')
    write_gro(typed_ethane, 'gmso.gro')

def parmed_top():
    ethane = mb.load('CC', smiles=True)
    ff = foyer.forcefields.load_OPLSAA()
    pmd_ethane = ff.apply(ethane)
    pmd_ethane.save('pmd.top', overwrite=True)
    pmd_ethane.save('pmd.gro', overwrite=True)

gmso_top()
parmed_top()