pyCRS : Conformer usage for Database and CRSManager¶

from scm.plams import Settings, init, finish, from_smiles
from scm.plams.recipes.adfcosmorsconformers import ADFCOSMORSConfJob, ADFCOSMORSConfFilter
from scm.conformers import ConformersJob

init()

mol = from_smiles("CC(=O)O")

conf_sett = Settings()
conf_sett.input.AMS.Generator.RDKit
conf_sett.input.AMS.Generator.RDKit.InitialNConformers = 50
conf_job = ConformersJob(name="conformers_uff", molecule=mol, settings=conf_sett)

dftb_sett = Settings()
dftb_sett.input.AMS.Task = "Optimize"
dftb_sett.input.DFTB

# ADFCOSMORSConfFilter(max number of conformers, max energy range)
fil1 = ADFCOSMORSConfFilter(20, 22)  # applied to UFF
fil2 = ADFCOSMORSConfFilter(10, 12)  # applied to DFTB
fil3 = ADFCOSMORSConfFilter(5, 7)  # applied to ADF gas phase

mol_info = {"CAS":"64-19-7", "IUPAC":"Acetic acid", "SMILES":"CC(=O)O"}

job = ADFCOSMORSConfJob(
    mol,
    conf_gen=conf_job,
    first_filter=fil1,
    additional=[(dftb_sett, fil2)],
    final_filter=fil3,
    coskf_name="acetic_acid",
    coskf_dir="coskf_acetic_acid",
    mol_info=mol_info,
)
job.run()

finish()

from pyCRS.Database import COSKFDatabase
import os

##################  Note: Ensure to download the coskf_acetic_acid or generating conformers for acetic acid before running the script  ##################

conformer_path = os.path.join(os.getcwd(), "coskf_acetic_acid")

db = COSKFDatabase("my_coskf_db.db")

db.add_compound("acetic_acid_0.coskf", coskf_path=conformer_path)
db.add_compound("acetic_acid_1.coskf", coskf_path=conformer_path)

print("Displaying all ConformerRow of acetic acid")
for key, rows in db.get_conformers("Acetic acid").items():
    for row in rows:
        print(f" {row}")

print("Displaying the CompoundRow of acetic acid with the lowest energy structure (default)")
for key, rows in db.get_compounds("Acetic acid").items():
    for row in rows:
        print(f" {row}")

db.visualize_conformers(identifier="acetic acid")

from scm.plams import Settings, init, finish
from pyCRS.Database import COSKFDatabase
from pyCRS.CRSManager import CRSSystem
import matplotlib.pyplot as plt
import numpy as np
import os

########  Note: Please ensure to install the ADFCRS-2018 database via amspackages (SCM->Packages) before running the below script  ########

db = COSKFDatabase("my_coskf_db.db")

db.add_compound("Water.coskf")

cal = CRSSystem()

mixture = {}
x_range = np.linspace(0, 1, 11)
for x in x_range:
    mixture["Acetic acid"] = x
    mixture["Water"] = 1 - x
    cal.add_Mixture(mixture, database="my_coskf_db.db", temperature=298.15, problem_type="activitycoef", conformer=True)

init()
cal.runCRSJob()
finish()

gamma_water = []
gamma_acid = []
comp_acetic = {}
for out in cal.outputs:
    res = out.get_results()
    gamma_acid.append(res["gamma"][0])
    gamma_water.append(res["gamma"][1])
    if comp_acetic == {}:
        comp_acetic = out.get_multispecies_dist()[0]
    else:
        tmp_comp_acetic = out.get_multispecies_dist()[0]
        for key in comp_acetic:
            comp_acetic[key].extend(tmp_comp_acetic[key])

fig, axs = plt.subplots(2)
axs[0].plot(x_range, gamma_acid, label="$\gamma_1$ (acetic)")
axs[0].plot(x_range, gamma_water, label="$\gamma_2$ (water)")

for struct, val in comp_acetic.items():
    axs[1].plot(x_range, val, label=os.path.basename(struct))

plt.setp(axs[0], ylabel="Activity coefficients")
plt.setp(axs[1], ylabel="Distribution")
for i in range(2):
    axs[i].legend()
    axs[i].grid()

plt.xlabel("$x_1$ (acetic acid)")
# plt.savefig('./pyCRS_activitycoef_conformer',dpi=300)
plt.show()