Skip to content

Commit

Permalink
improvements
Browse files Browse the repository at this point in the history 
* added progress option to track CENSO progress using printout to stderr
* updated NMR kt2 reference shielding constants
* fixed wb97x-V in TM ($disp3)
* updated README
* added information functional availability in -tutorial
* added load_balancing automatical O and P balance
* corrected NMR SD of sigma
* created censo_nmr_ref.json for user editable shift references
  • Loading branch information
@fabothch
fabothch authored Mar 1, 2021
1 parent cc57c95 commit 80bb98f
Show file tree
Hide file tree
Showing 15 changed files with 308 additions and 57 deletions.
10 changes: 8 additions & 2 deletions README.rst
View file
Original file line number Diff line number Diff line change
Expand Up @@ -58,19 +58,24 @@ settings and provide paths to the external programs e.g. `xtb`, `crest`, `orca`
vi /home/$USER/.ensorc
Interactive Documentation can be accessed:
**Interactive Documentation can be accessed:**

.. code::
$ censo -tutorial
Explainations on the commandline arguments can be printed by:

.. code::
$ censo --help
Online Documentation:
---------------------

Can be found following: https://fbohle.gitbook.io/censo/


The molecule numbering from the input structure ensemble is kept throughout the
entire program. There are several program parts which can be used to filter a structure
ensemble:
Expand All @@ -96,6 +101,7 @@ ensemble:


Usage:
------

.. code::
Expand Down
80 changes: 78 additions & 2 deletions censo_qm/cfg.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -5,7 +5,8 @@
"""
import os

__version__ = "1.0.4"
__version__ = "1.0.5"

DESCR = f"""
______________________________________________________________
| |
Expand Down Expand Up @@ -325,6 +326,16 @@

class NmrRef():
"""nmrreference data in the format: [reference-molecule][func-geometry][funcS][solvent]
h_tm_shieldings
c_tm_shieldings
f_tm_shieldings
si_tm_shieldings
p_tm_shieldings
h_orca_shieldings
c_orca_shieldings
f_orca_shieldings
si_orca_shieldings
p_orca_shieldings
"""
h_tm_shieldings = {
"TMS": {
Expand Down Expand Up @@ -388,6 +399,18 @@ class NmrRef():
"methanol": 31.644428750000003,
"thf": 31.633562058333336,
"toluene": 31.6153103
},
"kt2": {
"gas": 31.667160058333337,
"acetone": 31.629431858333334,
"chcl3": 31.634026916666674,
"acetonitrile": 31.627110300000002,
"ch2cl2": 31.645328925,
"dmso": 31.629352433333338,
"h2o": 31.635154825,
"methanol": 31.642592158333333,
"thf": 31.631994574999997,
"toluene": 31.612610516666663
}
},
"pbeh-3c": {
Expand Down Expand Up @@ -794,6 +817,18 @@ class NmrRef():
"methanol": 188.355007975,
"thf": 188.4689729,
"toluene": 188.594034275
},
"kt2": {
"gas": 189.78494644999998,
"acetone": 190.329502875,
"chcl3": 190.204013175,
"acetonitrile": 190.397052075,
"ch2cl2": 190.06665505,
"dmso": 190.4107424,
"h2o": 190.40970589999998,
"methanol": 190.188391875,
"thf": 190.2872299,
"toluene": 190.41299607500002
}
},
"pbeh-3c": {
Expand Down Expand Up @@ -2228,6 +2263,18 @@ class NmrRef():
"methanol": 331.3449225,
"thf": 331.3305993,
"toluene": 331.3160006
},
"kt2": {
"gas": 340.923509,
"acetone": 340.3208483,
"chcl3": 340.3430966,
"acetonitrile": 340.3155874,
"ch2cl2": 340.3599861,
"dmso": 340.3553012,
"h2o": 340.4180652,
"methanol": 340.3808603,
"thf": 340.348664,
"toluene": 340.3406705
}
},
"pbeh-3c": {
Expand Down Expand Up @@ -2572,4 +2619,33 @@ class NmrRef():
},
}
}

def NMRRef_to_dict(self):
"""Convert NMRRef data to a dict object"""
dict_ret = dict(
h_tm_shieldings=self.h_tm_shieldings,
c_tm_shieldings=self.c_tm_shieldings,
f_tm_shieldings=self.f_tm_shieldings,
si_tm_shieldings=self.si_tm_shieldings,
p_tm_shieldings=self.p_tm_shieldings,
h_orca_shieldings=self.h_orca_shieldings,
c_orca_shieldings=self.c_orca_shieldings,
f_orca_shieldings=self.f_orca_shieldings,
si_orca_shieldings=self.si_orca_shieldings,
p_orca_shieldings=self.p_orca_shieldings,
)
return dict_ret

def dict_to_NMRRef(self, dictionary):
"""Convert dict object to NMRRef data """
NmrRef_object = NmrRef()
NmrRef_object.h_tm_shieldings = dictionary.get('h_tm_shieldings', NmrRef_object.h_tm_shieldings)
NmrRef_object.c_tm_shieldings = dictionary.get('c_tm_shieldings', NmrRef_object.c_tm_shieldings)
NmrRef_object.f_tm_shieldings = dictionary.get('f_tm_shieldings', NmrRef_object.f_tm_shieldings)
NmrRef_object.si_tm_shieldings = dictionary.get('si_tm_shieldings', NmrRef_object.si_tm_shieldings)
NmrRef_object.p_tm_shieldings = dictionary.get('p_tm_shieldings', NmrRef_object.p_tm_shieldings)
NmrRef_object.h_orca_shieldings = dictionary.get('h_orca_shieldings', NmrRef_object.h_orca_shieldings)
NmrRef_object.c_orca_shieldings = dictionary.get('c_orca_shieldings', NmrRef_object.c_orca_shieldings)
NmrRef_object.f_orca_shieldings = dictionary.get('f_orca_shieldings', NmrRef_object.f_orca_shieldings)
NmrRef_object.si_orca_shieldings = dictionary.get('si_orca_shieldings', NmrRef_object.si_orca_shieldings)
NmrRef_object.p_orca_shieldings = dictionary.get('p_orca_shieldings', NmrRef_object.p_orca_shieldings)
return NmrRef_object
14 changes: 10 additions & 4 deletions censo_qm/cheapscreening.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -148,7 +148,6 @@ def part0(config, conformers, ensembledata):
"unpaired": config.unpaired,
"solvent": config.solvent,
"sm": config.sm_rrho,
"omp": config.omp,
"gfn_version": config.part0_gfnv,
"energy": 0.0,
"energy2": 0.0,
Expand Down Expand Up @@ -190,7 +189,6 @@ def part0(config, conformers, ensembledata):
"unpaired": config.unpaired,
"solvent": "gas",
"sm": "gas-phase",
"omp": config.omp,
"energy": 0.0,
"energy2": 0.0,
"success": False,
Expand Down Expand Up @@ -231,10 +229,18 @@ def part0(config, conformers, ensembledata):
calculate, store_confs, save_errors = ensemble2coord(
config, folder, calculate, store_confs, save_errors
)

# parallel calculation:
calculate = run_in_parallel(
config, q, resultq, job, config.maxthreads, calculate, instruction, folder
config,
q,
resultq,
job,
config.maxthreads,
config.omp,
calculate,
instruction,
config.balance,
folder
)

for conf in list(calculate):
Expand Down
14 changes: 14 additions & 0 deletions censo_qm/inputhandling.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -845,6 +845,16 @@ def cml(startup_description, options, argv=None):
" threads with each (omp) 4 cores --> 20 cores need to be available on "
"the machine.",
)
group7.add_argument(
"-balance",
"--balance",
dest="balance",
choices=["on", "off"],
action="store",
metavar="",
help="Automatically balance the number of threads and cores when a low number"
"of conformers is left. (never exceed O*P cores).",
)
group11 = parser.add_argument_group("Concerning overall mRRHO calculations")
group11.add_argument(
"-imagthr",
Expand Down Expand Up @@ -1066,6 +1076,7 @@ class internal_settings:
"scale": "scale",
"cosmorsparam": "cosmorsparam",
"progress": "progress",
"balance": "balance",
}
knownbasissets3 = [
"SVP",
Expand Down Expand Up @@ -1504,6 +1515,7 @@ def __init__(self):
("basis", {"default": "automatic", "type": str}),
("maxthreads", {"default": 1, "type": int}),
("omp", {"default": 1, "type": int}),
("balance", {"default": False, "type": bool}),
("cosmorsparam", {"default": "automatic", "type": str}),

]
Expand Down Expand Up @@ -1782,6 +1794,7 @@ def __init__(self, path=os.getcwd(), *args, **kwargs):
self.basis = "automatic"
self.maxthreads = 1
self.omp = 1
self.balance=False
self.cosmorsparam = "automatic"
# part0
self.part0 = False
Expand Down Expand Up @@ -2725,6 +2738,7 @@ def print_parameters(self):
info.append(["crestcheck", "checking the DFT-ensemble using CREST"])
info.append(["maxthreads", "maxthreads"])
info.append(["omp", "omp"])
info.append(["balance", "automatically balance maxthreads and omp"])

if self.part0:
# PART0:
Expand Down
58 changes: 41 additions & 17 deletions censo_qm/nmrproperties.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -5,13 +5,15 @@
import os
import shutil
import sys
import json
from random import normalvariate
from multiprocessing import JoinableQueue as Queue
from .cfg import (
PLENGTH,
DIGILEN,
AU2KCAL,
WARNLEN,
CODING,
NmrRef
)
from .parallel import run_in_parallel
Expand Down Expand Up @@ -90,7 +92,7 @@ def average_shieldings(config, calculate, element_ref_shield, energy, solv, rrho
sigma_std_dev[i] = []

for conf in calculate:
# get shielding constants
# get averaged shielding constants
if not element:
element = get_atom(
os.path.normpath(os.path.join(config.cwd, "CONF" + str(conf.id), "NMR"))
Expand All @@ -100,7 +102,7 @@ def average_shieldings(config, calculate, element_ref_shield, energy, solv, rrho
[conf.shieldings.get(eq_atom, 0.0) for eq_atom in chemeq[atom]]
) / len(chemeq[atom])
averaged[atom] = conf.bm_weight * sigma + averaged.get(atom, 0.0)

# get SD on shieldings based on SD Gsolv
if solv is not None:
get_std_dev = (
lambda conf: conf.lowlevel_gsolv_compare_info["std_dev"]
Expand All @@ -114,21 +116,18 @@ def average_shieldings(config, calculate, element_ref_shield, energy, solv, rrho
0.0, get_std_dev(conf)
)
calculate = calc_boltzmannweights(calculate, "free_energy", config.temperature)
tmp_sigma = {}
for i in range(1, config.nat + 1):
tmp_sigma[i] = 0
for conf in calculate:
# get shielding constants
if not element:
element = get_atom(
os.path.normpath(
os.path.join(config.cwd, "CONF" + str(conf.id), "NMR")
)
)
for atom in conf.shieldings.keys():
sigma = sum(
[conf.shieldings.get(eq_atom, 0.0) for eq_atom in chemeq[atom]]
) / len(chemeq[atom])
sigma_std_dev[atom].append(
conf.bm_weight * sigma + averaged.get(atom, 0.0)
)
tmp_sigma[atom] += conf.bm_weight * sigma
for atom in conf.shieldings.keys():
sigma_std_dev[atom].append(tmp_sigma[atom])


print("\nAveraged shielding constants:")
print("# in coord element σ(sigma) SD(σ based on SD Gsolv) shift σ_ref")
Expand All @@ -151,7 +150,6 @@ def average_shieldings(config, calculate, element_ref_shield, energy, solv, rrho
f"{std_dev:^ 24.6f} {make_shift(atom):>5} {element_ref_shield.get(element[atom], 0.0)}"
)
except Exception as e:
#print(e)
print(f"{atom:< {10}} {element[atom]:^{7}} {sigma:> {maxsigma}.2f}")
print("".ljust(int(80), "-"))

Expand Down Expand Up @@ -247,7 +245,22 @@ def write_anmrrc(config):
):
for element in ref_decision.keys():
ref_decision[element]["active"] = True
nmrref = NmrRef()

# read NMR_references
nmr_ref_user_path = os.path.expanduser(
os.path.join("~/.censo_assets/", "censo_nmr_ref.json")
)
if os.path.isfile(nmr_ref_user_path):
try:
with open(nmr_ref_user_path, "r", encoding=CODING, newline=None) as inp:
tmp_ref = json.load(inp)
except (ValueError, TypeError, FileNotFoundError):
print(
f"{'ERROR:':{WARNLEN}}Your censo_nmr_ref.json file in {nmr_ref_user_path} is corrupted!\n"
)
tmp_ref = {}
nmrref = NmrRef().dict_to_NMRRef(tmp_ref)
# end reading NMR_references
for element, value in ref_decision.items():
if value.get("active", False):
try:
Expand Down Expand Up @@ -746,7 +759,6 @@ def part4(config, conformers, store_confs, ensembledata):
"solvent": config.solvent,
"sm": config.sm4_j,
"success": False,
"omp": config.omp,
# nmractive nuclei
"h_active": config.h_active,
"c_active": config.c_active,
Expand Down Expand Up @@ -796,8 +808,10 @@ def part4(config, conformers, store_confs, ensembledata):
resultq,
job,
config.maxthreads,
config.omp,
calculate,
instruction_j,
config.balance,
folder,
)
for conf in list(calculate):
Expand Down Expand Up @@ -880,7 +894,6 @@ def part4(config, conformers, store_confs, ensembledata):
"solvent": config.solvent,
"sm": config.sm4_s,
"success": False,
"omp": config.omp,
# nmractive nuclei
"h_active": config.h_active,
"c_active": config.c_active,
Expand Down Expand Up @@ -950,8 +963,10 @@ def part4(config, conformers, store_confs, ensembledata):
resultq,
job,
config.maxthreads,
config.omp,
calculate,
instruction_s,
config.balance,
folder,
)
for conf in list(calculate):
Expand Down Expand Up @@ -1053,7 +1068,16 @@ def part4(config, conformers, store_confs, ensembledata):
# write generic:
instructgeneric = {"jobtype": "genericout", "nat": int(config.nat)}
calculate = run_in_parallel(
config, q, resultq, job, config.maxthreads, calculate, instructgeneric, folder
config,
q,
resultq,
job,
config.maxthreads,
config.omp,
calculate,
instructgeneric,
False,
folder
)

# printout the averaged shielding constants
Expand Down
Loading

0 comments on commit 80bb98f

Please sign in to comment.