Adds Formatter to Project

This commit is contained in:
2023-06-11 16:04:25 -03:00
parent 82f3092f3e
commit c4dae5e8d1
29 changed files with 1151 additions and 721 deletions

View File

@@ -1,4 +1,3 @@
from diceplayer.shared.utils.logger import Logger
logger = Logger(__name__)

View File

@@ -1,12 +1,12 @@
from diceplayer.shared.interface.dice_interface import DiceInterface
from diceplayer.player import Player
from diceplayer import logger
from diceplayer.player import Player
from diceplayer.shared.interface.dice_interface import DiceInterface
from pathlib import Path
import argparse
import logging
import pickle
import sys
from pathlib import Path
__VERSION = "v0.0.1"
@@ -25,18 +25,20 @@ def main():
"-v", "--version", action="version", version="diceplayer-" + __VERSION
)
parser.add_argument(
"-i", "--input",
"-i",
"--input",
dest="infile",
default="control.yml",
metavar="INFILE",
help="input file of diceplayer [default = control.in]"
help="input file of diceplayer [default = control.in]",
)
parser.add_argument(
"-o", "--output",
"-o",
"--output",
dest="outfile",
default="run.log",
metavar="OUTFILE",
help="output file of diceplayer [default = run.log]"
help="output file of diceplayer [default = run.log]",
)
args = parser.parse_args()

View File

@@ -1,24 +1,24 @@
from diceplayer.shared.interface.gaussian_interface import GaussianInterface
from diceplayer.shared.interface.dice_interface import DiceInterface
from diceplayer.shared.utils.dataclass_protocol import Dataclass
from diceplayer import logger
from diceplayer.shared.config.dice_config import DiceConfig
from diceplayer.shared.config.gaussian_config import GaussianDTO
from diceplayer.shared.config.player_config import PlayerConfig
from diceplayer.shared.config.dice_config import DiceConfig
from diceplayer.shared.utils.misc import weekday_date_time
from diceplayer.shared.environment.atom import Atom
from diceplayer.shared.environment.molecule import Molecule
from diceplayer.shared.environment.system import System
from diceplayer.shared.environment.atom import Atom
from diceplayer.shared.interface.dice_interface import DiceInterface
from diceplayer.shared.interface.gaussian_interface import GaussianInterface
from diceplayer.shared.utils.dataclass_protocol import Dataclass
from diceplayer.shared.utils.misc import weekday_date_time
from diceplayer.shared.utils.ptable import atomsymb
from diceplayer import logger
from dataclasses import fields
from typing import Type, Tuple
from pathlib import Path
import pickle
import yaml
import sys
import os
import os
import pickle
import sys
from dataclasses import fields
from pathlib import Path
from typing import Tuple, Type
ENV = ["OMP_STACKSIZE"]
@@ -29,9 +29,7 @@ class Player:
raise ValueError("Must specify either infile or optimization")
elif infile is not None:
self.config = self.set_config(
self.read_keywords(infile)
)
self.config = self.set_config(self.read_keywords(infile))
self.system = System()
@@ -60,7 +58,6 @@ class Player:
)
for cycle in range(self.initial_cycle, self.initial_cycle + self.config.maxcyc):
logger.info(
f"------------------------------------------------------------------------------------------\n"
f" Step # {cycle}\n"
@@ -78,9 +75,7 @@ class Player:
def prepare_system(self):
for i, mol in enumerate(self.system.molecule):
logger.info(
f"Molecule {i + 1} - {mol.molname}"
)
logger.info(f"Molecule {i + 1} - {mol.molname}")
mol.print_mol_info()
logger.info(
@@ -113,7 +108,6 @@ class Player:
geoms_file_path.touch()
def print_keywords(self) -> None:
def log_keywords(config: Dataclass, dto: Type[Dataclass]):
for key in sorted(list(map(lambda f: f.name, fields(dto)))):
if getattr(config, key) is not None:
@@ -162,9 +156,7 @@ class Player:
with open(self.config.dice.ljname) as file:
ljc_data = file.readlines()
else:
raise RuntimeError(
f"Potential file {self.config.dice.ljname} not found."
)
raise RuntimeError(f"Potential file {self.config.dice.ljname} not found.")
combrule = ljc_data.pop(0).split()[0]
if combrule not in ("*", "+"):
@@ -191,7 +183,6 @@ class Player:
)
for i in range(ntypes):
try:
nsites, molname = ljc_data.pop(0).split()[:2]
except ValueError:
@@ -209,10 +200,7 @@ class Player:
atom_fields = ["lbl", "na", "rx", "ry", "rz", "chg", "eps", "sig"]
for j in range(nsites):
new_atom = dict(zip(
atom_fields,
ljc_data.pop(0).split()
))
new_atom = dict(zip(atom_fields, ljc_data.pop(0).split()))
self.system.molecule[i].add_atom(
Atom(**self.validate_atom_dict(i, j, new_atom))
)
@@ -227,16 +215,12 @@ class Player:
)
logger.info(f"Combination rule: {self.config.dice.combrule}")
logger.info(
f"Types of molecules: {len(self.system.molecule)}\n"
)
logger.info(f"Types of molecules: {len(self.system.molecule)}\n")
i = 0
for mol in self.system.molecule:
i += 1
logger.info(
"{} atoms in molecule type {}:".format(len(mol.atom), i)
)
logger.info("{} atoms in molecule type {}:".format(len(mol.atom), i))
logger.info(
"---------------------------------------------------------------------------------"
)
@@ -285,42 +269,37 @@ class Player:
self.gaussian_interface.reset()
if self.config.opt:
if 'position' not in result:
raise RuntimeError(
'Optimization failed. No position found in result.'
)
if "position" not in result:
raise RuntimeError("Optimization failed. No position found in result.")
self.system.update_molecule(result['position'])
self.system.update_molecule(result["position"])
else:
if 'charges' not in result:
if "charges" not in result:
raise RuntimeError(
'Charges optimization failed. No charges found in result.'
"Charges optimization failed. No charges found in result."
)
diff = self.system.molecule[0] \
.update_charges(result['charges'])
diff = self.system.molecule[0].update_charges(result["charges"])
self.system.print_charges_and_dipole(cycle)
self.print_geoms(cycle)
if diff < self.config.gaussian.chg_tol:
logger.info(
f'Charges converged after {cycle} cycles.'
)
logger.info(f"Charges converged after {cycle} cycles.")
raise StopIteration()
def print_geoms(self, cycle: int):
with open(self.config.geoms_file, 'a') as file:
file.write(f'Cycle # {cycle}\n')
with open(self.config.geoms_file, "a") as file:
file.write(f"Cycle # {cycle}\n")
for atom in self.system.molecule[0].atom:
symbol = atomsymb[atom.na]
file.write(
f'{symbol:<2s} {atom.rx:>10.6f} {atom.ry:>10.6f} {atom.rz:>10.6f}\n'
f"{symbol:<2s} {atom.rx:>10.6f} {atom.ry:>10.6f} {atom.rz:>10.6f}\n"
)
file.write('\n')
file.write("\n")
@staticmethod
def validate_atom_dict(molecule_type, molecule_site, atom_dict: dict) -> dict:
@@ -329,69 +308,69 @@ class Player:
if len(atom_dict) < 8:
raise ValueError(
f'Invalid number of fields for site {molecule_site} for molecule type {molecule_type}.'
f"Invalid number of fields for site {molecule_site} for molecule type {molecule_type}."
)
try:
atom_dict['lbl'] = int(atom_dict['lbl'])
atom_dict["lbl"] = int(atom_dict["lbl"])
except Exception:
raise ValueError(
f'Invalid lbl fields for site {molecule_site} for molecule type {molecule_type}.'
f"Invalid lbl fields for site {molecule_site} for molecule type {molecule_type}."
)
try:
atom_dict['na'] = int(atom_dict['na'])
atom_dict["na"] = int(atom_dict["na"])
except Exception:
raise ValueError(
f'Invalid na fields for site {molecule_site} for molecule type {molecule_type}.'
f"Invalid na fields for site {molecule_site} for molecule type {molecule_type}."
)
try:
atom_dict['rx'] = float(atom_dict['rx'])
atom_dict["rx"] = float(atom_dict["rx"])
except Exception:
raise ValueError(
f'Invalid rx fields for site {molecule_site} for molecule type {molecule_type}. '
f'Value must be a float.'
f"Invalid rx fields for site {molecule_site} for molecule type {molecule_type}. "
f"Value must be a float."
)
try:
atom_dict['ry'] = float(atom_dict['ry'])
atom_dict["ry"] = float(atom_dict["ry"])
except Exception:
raise ValueError(
f'Invalid ry fields for site {molecule_site} for molecule type {molecule_type}. '
f'Value must be a float.'
f"Invalid ry fields for site {molecule_site} for molecule type {molecule_type}. "
f"Value must be a float."
)
try:
atom_dict['rz'] = float(atom_dict['rz'])
atom_dict["rz"] = float(atom_dict["rz"])
except Exception:
raise ValueError(
f'Invalid rz fields for site {molecule_site} for molecule type {molecule_type}. '
f'Value must be a float.'
f"Invalid rz fields for site {molecule_site} for molecule type {molecule_type}. "
f"Value must be a float."
)
try:
atom_dict['chg'] = float(atom_dict['chg'])
atom_dict["chg"] = float(atom_dict["chg"])
except Exception:
raise ValueError(
f'Invalid chg fields for site {molecule_site} for molecule type {molecule_type}. '
f'Value must be a float.'
f"Invalid chg fields for site {molecule_site} for molecule type {molecule_type}. "
f"Value must be a float."
)
try:
atom_dict['eps'] = float(atom_dict['eps'])
atom_dict["eps"] = float(atom_dict["eps"])
except Exception:
raise ValueError(
f'Invalid eps fields for site {molecule_site} for molecule type {molecule_type}. '
f'Value must be a float.'
f"Invalid eps fields for site {molecule_site} for molecule type {molecule_type}. "
f"Value must be a float."
)
try:
atom_dict['sig'] = float(atom_dict['sig'])
atom_dict["sig"] = float(atom_dict["sig"])
except Exception:
raise ValueError(
f'Invalid sig fields for site {molecule_site} for molecule type {molecule_type}. '
f'Value must be a float.'
f"Invalid sig fields for site {molecule_site} for molecule type {molecule_type}. "
f"Value must be a float."
)
return atom_dict
@@ -400,9 +379,7 @@ class Player:
formatstr = "{:<3d} {:>3d} {:>10.5f} {:>10.5f} {:>10.5f} {:>10.6f} {:>9.5f} {:>7.4f} {:>9.4f}"
mol = self.system.molecule[0]
logger.info(
"{} atoms in molecule type {}:".format(len(mol.atom), 1)
)
logger.info("{} atoms in molecule type {}:".format(len(mol.atom), 1))
logger.info(
"---------------------------------------------------------------------------------"
)
@@ -432,28 +409,21 @@ class Player:
def save_run_in_pickle(self, cycle):
try:
with open('latest-step.pkl', 'wb') as pickle_file:
pickle.dump(
(self.config, self.system, cycle),
pickle_file
)
with open("latest-step.pkl", "wb") as pickle_file:
pickle.dump((self.config, self.system, cycle), pickle_file)
except Exception:
raise RuntimeError(
f'Could not save pickle file latest-step.pkl.'
)
raise RuntimeError(f"Could not save pickle file latest-step.pkl.")
@staticmethod
def load_run_from_pickle() -> Tuple[PlayerConfig, System, int]:
pickle_path = Path("latest-step.pkl")
try:
with open(pickle_path, 'rb') as pickle_file:
with open(pickle_path, "rb") as pickle_file:
save = pickle.load(pickle_file)
return save[0], save[1], save[2] + 1
except Exception:
raise RuntimeError(
f'Could not load pickle file {pickle_path}.'
)
raise RuntimeError(f"Could not load pickle file {pickle_path}.")
@staticmethod
def set_config(data: dict) -> PlayerConfig:
@@ -461,18 +431,16 @@ class Player:
@staticmethod
def read_keywords(infile) -> dict:
with open(infile, 'r') as yml_file:
with open(infile, "r") as yml_file:
config = yaml.load(yml_file, Loader=yaml.SafeLoader)
if "diceplayer" in config:
return config.get("diceplayer")
raise RuntimeError(
f'Could not find diceplayer section in {infile}.'
)
raise RuntimeError(f"Could not find diceplayer section in {infile}.")
@classmethod
def from_file(cls, infile: str) -> 'Player':
def from_file(cls, infile: str) -> "Player":
return cls(infile=infile)
@classmethod

View File

@@ -1,7 +1,8 @@
from diceplayer.shared.utils.dataclass_protocol import Dataclass
from dataclasses import dataclass
from dacite import from_dict
from dataclasses import dataclass
from typing import List
@@ -10,6 +11,7 @@ class DiceConfig(Dataclass):
"""
Data Transfer Object for the Dice configuration.
"""
ljname: str
outname: str
dens: float
@@ -22,24 +24,17 @@ class DiceConfig(Dataclass):
press: float = 1.0
temp: float = 300.0
progname: str = "dice"
randominit: str = 'first'
randominit: str = "first"
def __post_init__(self):
if not isinstance(self.ljname, str):
raise ValueError(
"Error: 'ljname' keyword not specified in config file"
)
raise ValueError("Error: 'ljname' keyword not specified in config file")
if not isinstance(self.outname, str):
raise ValueError(
"Error: 'outname' keyword not specified in config file"
)
raise ValueError("Error: 'outname' keyword not specified in config file")
if not isinstance(self.dens, float):
raise ValueError(
"Error: 'dens' keyword not specified in config file"
)
raise ValueError("Error: 'dens' keyword not specified in config file")
if not isinstance(self.nmol, list):
raise ValueError(

View File

@@ -1,31 +1,29 @@
from diceplayer.shared.utils.dataclass_protocol import Dataclass
from dataclasses import dataclass
from dacite import from_dict
from dataclasses import dataclass
@dataclass
class GaussianDTO(Dataclass):
"""
Data Transfer Object for the Gaussian configuration.
"""
level: str
qmprog: str
chgmult = [0, 1]
pop: str = 'chelpg'
pop: str = "chelpg"
chg_tol: float = 0.01
keywords: str = None
def __post_init__(self):
if self.qmprog not in ("g03", "g09", "g16"):
raise ValueError(
"Error: invalid qmprog value."
)
raise ValueError("Error: invalid qmprog value.")
if self.level is None:
raise ValueError(
"Error: 'level' keyword not specified in config file."
)
raise ValueError("Error: 'level' keyword not specified in config file.")
@classmethod
def from_dict(cls, param: dict):

View File

@@ -1,9 +1,10 @@
from diceplayer.shared.utils.dataclass_protocol import Dataclass
from diceplayer.shared.config.gaussian_config import GaussianDTO
from diceplayer.shared.config.dice_config import DiceConfig
from diceplayer.shared.config.gaussian_config import GaussianDTO
from diceplayer.shared.utils.dataclass_protocol import Dataclass
from dacite import from_dict
from dataclasses import dataclass
from dacite import from_dict
@dataclass
@@ -11,6 +12,7 @@ class PlayerConfig(Dataclass):
"""
Data Transfer Object for the player configuration.
"""
opt: bool
maxcyc: int
nprocs: int
@@ -21,10 +23,10 @@ class PlayerConfig(Dataclass):
mem: int = None
switchcyc: int = 3
qmprog: str = 'g16'
qmprog: str = "g16"
altsteps: int = 20000
geoms_file = 'geoms.xyz'
simulation_dir = 'simfiles'
geoms_file = "geoms.xyz"
simulation_dir = "simfiles"
def __post_init__(self):
MIN_STEP = 20000
@@ -33,16 +35,12 @@ class PlayerConfig(Dataclass):
@classmethod
def from_dict(cls, param: dict):
if param['dice'] is None:
raise ValueError(
"Error: 'dice' keyword not specified in config file."
)
param['dice'] = DiceConfig.from_dict(param['dice'])
if param["dice"] is None:
raise ValueError("Error: 'dice' keyword not specified in config file.")
param["dice"] = DiceConfig.from_dict(param["dice"])
if param['gaussian'] is None:
raise ValueError(
"Error: 'gaussian' keyword not specified in config file."
)
param['gaussian'] = GaussianDTO.from_dict(param['gaussian'])
if param["gaussian"] is None:
raise ValueError("Error: 'gaussian' keyword not specified in config file.")
param["gaussian"] = GaussianDTO.from_dict(param["gaussian"])
return from_dict(PlayerConfig, param)

View File

@@ -17,15 +17,15 @@ class Atom:
"""
def __init__(
self,
lbl: int,
na: int,
rx: float,
ry: float,
rz: float,
chg: float,
eps: float,
sig: float,
self,
lbl: int,
na: int,
rx: float,
ry: float,
rz: float,
chg: float,
eps: float,
sig: float,
) -> None:
"""
The constructor function __init__ is used to create new instances of the Atom class.

View File

@@ -1,18 +1,18 @@
from __future__ import annotations
from diceplayer.shared.utils.ptable import ghost_number
from diceplayer import logger
from diceplayer.shared.environment.atom import Atom
from diceplayer.shared.utils.misc import BOHR2ANG
from diceplayer import logger
from diceplayer.shared.utils.ptable import ghost_number
from nptyping import NDArray, Shape, Float
from numpy.linalg import linalg
import numpy as np
from nptyping import Float, NDArray, Shape
from numpy.linalg import linalg
from typing import List, Any, Tuple, Union
from copy import deepcopy
import logging
import math
from copy import deepcopy
from typing import Any, List, Tuple, Union
class Molecule:
@@ -127,9 +127,9 @@ class Molecule:
dx = atom1.rx - atom2.rx
dy = atom1.ry - atom2.ry
dz = atom1.rz - atom2.rz
distances.append(math.sqrt(dx ** 2 + dy ** 2 + dz ** 2))
distances.append(math.sqrt(dx**2 + dy**2 + dz**2))
return np.array(distances).reshape(dim, dim-1)
return np.array(distances).reshape(dim, dim - 1)
def inertia_tensor(self) -> NDArray[Shape["3, 3"], Float]:
"""
@@ -147,9 +147,9 @@ class Molecule:
dy = atom.ry - self.com[1]
dz = atom.rz - self.com[2]
Ixx += atom.mass * (dy ** 2 + dz ** 2)
Iyy += atom.mass * (dz ** 2 + dx ** 2)
Izz += atom.mass * (dx ** 2 + dy ** 2)
Ixx += atom.mass * (dy**2 + dz**2)
Iyy += atom.mass * (dz**2 + dx**2)
Izz += atom.mass * (dx**2 + dy**2)
Ixy += atom.mass * dx * dy * -1
Ixz += atom.mass * dx * dz * -1
@@ -169,7 +169,9 @@ class Molecule:
try:
evals, evecs = linalg.eigh(self.inertia_tensor())
except ValueError:
raise RuntimeError("Error: diagonalization of inertia tensor did not converge")
raise RuntimeError(
"Error: diagonalization of inertia tensor did not converge"
)
return evals, evecs
@@ -355,7 +357,7 @@ class Molecule:
)
)
def minimum_distance(self, molec: 'Molecule') -> float:
def minimum_distance(self, molec: "Molecule") -> float:
"""
Return the minimum distance between two molecules
@@ -374,6 +376,6 @@ class Molecule:
dx = atom1.rx - atom2.rx
dy = atom1.ry - atom2.ry
dz = atom1.rz - atom2.rz
distances.append(math.sqrt(dx ** 2 + dy ** 2 + dz ** 2))
distances.append(math.sqrt(dx**2 + dy**2 + dz**2))
return min(distances)
return min(distances)

View File

@@ -1,13 +1,14 @@
from diceplayer.shared.environment.molecule import Molecule
from diceplayer.shared.utils.ptable import atomsymb
from diceplayer.shared.utils.misc import BOHR2ANG
from diceplayer import logger
from diceplayer.shared.environment.molecule import Molecule
from diceplayer.shared.utils.misc import BOHR2ANG
from diceplayer.shared.utils.ptable import atomsymb
from typing import List, Tuple, TextIO
from copy import deepcopy
from numpy import linalg
import numpy as np
from numpy import linalg
import math
from copy import deepcopy
from typing import List, TextIO, Tuple
class System:
@@ -59,7 +60,6 @@ class System:
logger.info(f"RMSD = {rmsd:>8.5f} Angstrom")
def rmsd_fit(self, p_index: int, r_index: int) -> Tuple[float, Molecule]:
projecting_mol = self.molecule[p_index]
reference_mol = self.molecule[r_index]
@@ -117,9 +117,9 @@ class System:
rmsd = 0
for i in range(dim):
rmsd += (
(x[i, 0] - y[i, 0]) ** 2
+ (x[i, 1] - y[i, 1]) ** 2
+ (x[i, 2] - y[i, 2]) ** 2
(x[i, 0] - y[i, 0]) ** 2
+ (x[i, 1] - y[i, 1]) ** 2
+ (x[i, 2] - y[i, 2]) ** 2
)
rmsd = math.sqrt(rmsd / dim)
@@ -214,12 +214,12 @@ class System:
#
def print_charges_and_dipole(self, cycle: int) -> None:
"""
Print the charges and dipole of the molecule in the Output file
Print the charges and dipole of the molecule in the Output file
Args:
cycle (int): Number of the cycle
fh (TextIO): Output file
"""
Args:
cycle (int): Number of the cycle
fh (TextIO): Output file
"""
logger.info("Cycle # {}\n".format(cycle))
logger.info("Number of site: {}\n".format(len(self.molecule[0].atom)))
@@ -228,6 +228,10 @@ class System:
logger.info(
"{:>10.6f} {:>10.6f} {:>10.6f} {:>10.6f} {:>10.6f}\n".format(
chargesAndDipole[0], chargesAndDipole[1], chargesAndDipole[2], chargesAndDipole[3], chargesAndDipole[4]
chargesAndDipole[0],
chargesAndDipole[1],
chargesAndDipole[2],
chargesAndDipole[3],
chargesAndDipole[4],
)
)

View File

@@ -7,10 +7,7 @@ from abc import ABC, abstractmethod
class Interface(ABC):
__slots__ = [
'step',
'system'
]
__slots__ = ["step", "system"]
def __init__(self):
self.system: System | None = None

View File

@@ -1,20 +1,21 @@
from __future__ import annotations
from diceplayer import logger
from diceplayer.shared.config.player_config import PlayerConfig
from diceplayer.shared.environment.system import System
from diceplayer.shared.interface import Interface
from diceplayer import logger
from multiprocessing import Process, connection
from setproctitle import setproctitle
from typing import Final, TextIO
from pathlib import Path
import subprocess
import shutil
import random
import time
import sys
import os
import random
import shutil
import subprocess
import sys
import time
from multiprocessing import Process, connection
from pathlib import Path
from typing import Final, TextIO
DICE_END_FLAG: Final[str] = "End of simulation"
DICE_FLAG_LINE: Final[int] = -2
@@ -74,19 +75,11 @@ class DiceInterface(Interface):
if not simulation_dir.exists():
simulation_dir.mkdir(parents=True)
proc_dir = Path(
simulation_dir,
f"step{cycle:02d}",
f"p{proc:02d}"
)
proc_dir = Path(simulation_dir, f"step{cycle:02d}", f"p{proc:02d}")
proc_dir.mkdir(parents=True, exist_ok=True)
def _make_dice_inputs(self, cycle, proc):
proc_dir = Path(
self.step.simulation_dir,
f"step{cycle:02d}",
f"p{proc:02d}"
)
proc_dir = Path(self.step.simulation_dir, f"step{cycle:02d}", f"p{proc:02d}")
self._make_potentials(proc_dir)
@@ -94,17 +87,17 @@ class DiceInterface(Interface):
# This is logic is used to make the initial configuration file
# for the next cycle using the last.xyz file from the previous cycle.
if self.step.dice.randominit == 'first' and cycle > 1:
if self.step.dice.randominit == "first" and cycle > 1:
last_xyz = Path(
self.step.simulation_dir,
f"step{(cycle - 1):02d}",
f"p{proc:02d}",
"last.xyz"
"last.xyz",
)
if not last_xyz.exists():
raise FileNotFoundError(f"File {last_xyz} not found.")
with open(last_xyz, 'r') as last_xyz_file:
with open(last_xyz, "r") as last_xyz_file:
self._make_init_file(proc_dir, last_xyz_file)
last_xyz_file.seek(0)
self.step.dice.dens = self._new_density(last_xyz_file)
@@ -122,11 +115,7 @@ class DiceInterface(Interface):
def _run_dice(self, cycle: int, proc: int):
working_dir = os.getcwd()
proc_dir = Path(
self.step.simulation_dir,
f"step{cycle:02d}",
f"p{proc:02d}"
)
proc_dir = Path(self.step.simulation_dir, f"step{cycle:02d}", f"p{proc:02d}")
logger.info(
f"Simulation process {str(proc_dir)} initiated with pid {os.getpid()}"
@@ -134,7 +123,7 @@ class DiceInterface(Interface):
os.chdir(proc_dir)
if not (self.step.dice.randominit == 'first' and cycle > 1):
if not (self.step.dice.randominit == "first" and cycle > 1):
self.run_dice_file(cycle, proc, "NVT.ter")
if len(self.step.dice.nstep) == 2:
@@ -168,17 +157,16 @@ class DiceInterface(Interface):
SECONDARY_MOLECULE_LENGTH = 0
for i in range(1, len(self.step.dice.nmol)):
SECONDARY_MOLECULE_LENGTH += self.step.dice.nmol[i] * len(self.system.molecule[i].atom)
SECONDARY_MOLECULE_LENGTH += self.step.dice.nmol[i] * len(
self.system.molecule[i].atom
)
xyz_lines = xyz_lines[-SECONDARY_MOLECULE_LENGTH:]
input_file = Path(proc_dir, self.step.dice.outname + ".xy")
with open(input_file, 'w') as f:
with open(input_file, "w") as f:
for atom in self.system.molecule[0].atom:
f.write(
f"{atom.rx:>10.6f} {atom.ry:>10.6f} {atom.rz:>10.6f}\n"
)
f.write(f"{atom.rx:>10.6f} {atom.ry:>10.6f} {atom.rz:>10.6f}\n")
for line in xyz_lines:
atom = line.split()
@@ -205,7 +193,7 @@ class DiceInterface(Interface):
def _make_nvt_ter(self, cycle, proc_dir):
file = Path(proc_dir, "NVT.ter")
with open(file, 'w') as f:
with open(file, "w") as f:
f.write(f"title = {self.title} - NVT Thermalization\n")
f.write(f"ncores = {self.step.ncores}\n")
f.write(f"ljname = {self.step.dice.ljname}\n")
@@ -236,9 +224,8 @@ class DiceInterface(Interface):
f.write(f"upbuf = {self.step.dice.upbuf}")
def _make_nvt_eq(self, cycle, proc_dir):
file = Path(proc_dir, "NVT.eq")
with open(file, 'w') as f:
with open(file, "w") as f:
f.write(f"title = {self.title} - NVT Production\n")
f.write(f"ncores = {self.step.ncores}\n")
f.write(f"ljname = {self.step.dice.ljname}\n")
@@ -269,9 +256,8 @@ class DiceInterface(Interface):
f.write("seed = {}\n".format(seed))
def _make_npt_ter(self, cycle, proc_dir):
file = Path(proc_dir, "NPT.ter")
with open(file, 'w') as f:
with open(file, "w") as f:
f.write(f"title = {self.title} - NPT Thermalization\n")
f.write(f"ncores = {self.step.ncores}\n")
f.write(f"ljname = {self.step.dice.ljname}\n")
@@ -303,7 +289,7 @@ class DiceInterface(Interface):
def _make_npt_eq(self, proc_dir):
file = Path(proc_dir, "NPT.eq")
with open(file, 'w') as f:
with open(file, "w") as f:
f.write(f"title = {self.title} - NPT Production\n")
f.write(f"ncores = {self.step.ncores}\n")
f.write(f"ljname = {self.step.dice.ljname}\n")
@@ -332,7 +318,7 @@ class DiceInterface(Interface):
fstr = "{:<3d} {:>3d} {:>10.5f} {:>10.5f} {:>10.5f} {:>10.6f} {:>9.5f} {:>7.4f}\n"
file = Path(proc_dir, self.step.dice.ljname)
with open(file, 'w') as f:
with open(file, "w") as f:
f.write(f"{self.step.dice.combrule}\n")
f.write(f"{len(self.step.dice.nmol)}\n")
@@ -370,7 +356,9 @@ class DiceInterface(Interface):
)
def run_dice_file(self, cycle: int, proc: int, file_name: str):
with open(Path(file_name), 'r') as infile, open(Path(file_name + ".out"), 'w') as outfile:
with open(Path(file_name), "r") as infile, open(
Path(file_name + ".out"), "w"
) as outfile:
if shutil.which("bash") is not None:
exit_status = subprocess.call(
[
@@ -385,11 +373,15 @@ class DiceInterface(Interface):
)
if exit_status != 0:
raise RuntimeError(f"Dice process step{cycle:02d}-p{proc:02d} did not exit properly")
raise RuntimeError(
f"Dice process step{cycle:02d}-p{proc:02d} did not exit properly"
)
with open(Path(file_name + ".out"), 'r') as outfile:
with open(Path(file_name + ".out"), "r") as outfile:
flag = outfile.readlines()[DICE_FLAG_LINE].strip()
if flag != DICE_END_FLAG:
raise RuntimeError(f"Dice process step{cycle:02d}-p{proc:02d} did not exit properly")
raise RuntimeError(
f"Dice process step{cycle:02d}-p{proc:02d} did not exit properly"
)
logger.info(f"Dice {file_name} - step{cycle:02d}-p{proc:02d} exited properly")

View File

@@ -1,24 +1,23 @@
from __future__ import annotations
from diceplayer import logger
from diceplayer.shared.config.player_config import PlayerConfig
from diceplayer.shared.environment.atom import Atom
from diceplayer.shared.environment.molecule import Molecule
from diceplayer.shared.environment.system import System
from diceplayer.shared.environment.atom import Atom
from diceplayer.shared.interface import Interface
from diceplayer.shared.utils.misc import date_time
from diceplayer.shared.utils.ptable import atomsymb
from diceplayer.shared.interface import Interface
from diceplayer import logger
from typing import Tuple, List, Dict, Any
from nptyping import NDArray
import numpy as np
from nptyping import NDArray
from pathlib import Path
import os
import shutil
import subprocess
import textwrap
import shutil
import os
from pathlib import Path
from typing import Any, Dict, List, Tuple
class GaussianInterface(Interface):
@@ -33,10 +32,7 @@ class GaussianInterface(Interface):
self._copy_chk_file_from_previous_step(cycle)
asec_charges = self.populate_asec_vdw(cycle)
self._make_gaussian_input_file(
cycle,
asec_charges
)
self._make_gaussian_input_file(cycle, asec_charges)
self._run_gaussian(cycle)
self._run_formchk(cycle)
@@ -49,9 +45,7 @@ class GaussianInterface(Interface):
raise NotImplementedError("Optimization not implemented yet.")
else:
return_value['charges'] = np.array(
self._read_charges_from_fchk(cycle)
)
return_value["charges"] = np.array(self._read_charges_from_fchk(cycle))
return return_value
@@ -60,36 +54,22 @@ class GaussianInterface(Interface):
del self.system
def _make_qm_dir(self, cycle: int):
qm_dir_path = Path(
self.step.simulation_dir,
f"step{cycle:02d}",
"qm"
)
qm_dir_path = Path(self.step.simulation_dir, f"step{cycle:02d}", "qm")
if not qm_dir_path.exists():
qm_dir_path.mkdir()
def _copy_chk_file_from_previous_step(self, cycle: int):
current_chk_file_path = Path(
self.step.simulation_dir,
f"step{cycle:02d}",
"qm",
f"asec.chk"
self.step.simulation_dir, f"step{cycle:02d}", "qm", f"asec.chk"
)
if current_chk_file_path.exists():
raise FileExistsError(
f"File {current_chk_file_path} already exists."
)
raise FileExistsError(f"File {current_chk_file_path} already exists.")
previous_chk_file_path = Path(
self.step.simulation_dir,
f"step{(cycle - 1):02d}",
"qm",
f"asec.chk"
self.step.simulation_dir, f"step{(cycle - 1):02d}", "qm", f"asec.chk"
)
if not previous_chk_file_path.exists():
raise FileNotFoundError(
f"File {previous_chk_file_path} does not exist."
)
raise FileNotFoundError(f"File {previous_chk_file_path} does not exist.")
shutil.copy(previous_chk_file_path, current_chk_file_path)
@@ -104,17 +84,19 @@ class GaussianInterface(Interface):
for proc in range(1, self.step.nprocs + 1):
proc_charges.append(self._read_charges_from_last_step(cycle, proc))
asec_charges, thickness, picked_mols = \
self._evaluate_proc_charges(nsites_total, proc_charges)
asec_charges, thickness, picked_mols = self._evaluate_proc_charges(
nsites_total, proc_charges
)
logger.info(f"In average, {(sum(picked_mols) / norm_factor):^7.2f} molecules\n"
f"were selected from each of the {len(picked_mols)} configurations\n"
f"of the production simulations to form the ASEC, comprising a shell with\n"
f"minimum thickness of {(sum(thickness) / norm_factor):>6.2f} Angstrom\n"
)
logger.info(
f"In average, {(sum(picked_mols) / norm_factor):^7.2f} molecules\n"
f"were selected from each of the {len(picked_mols)} configurations\n"
f"of the production simulations to form the ASEC, comprising a shell with\n"
f"minimum thickness of {(sum(thickness) / norm_factor):>6.2f} Angstrom\n"
)
for charge in asec_charges:
charge['chg'] = charge['chg'] / norm_factor
charge["chg"] = charge["chg"] / norm_factor
return asec_charges
@@ -135,23 +117,19 @@ class GaussianInterface(Interface):
def _read_charges_from_last_step(self, cycle: int, proc: int) -> list[str]:
last_xyz_file_path = Path(
self.step.simulation_dir,
f"step{cycle:02d}",
f"p{proc:02d}",
"last.xyz"
self.step.simulation_dir, f"step{cycle:02d}", f"p{proc:02d}", "last.xyz"
)
if not last_xyz_file_path.exists():
raise FileNotFoundError(
f"File {last_xyz_file_path} does not exist."
)
raise FileNotFoundError(f"File {last_xyz_file_path} does not exist.")
with open(last_xyz_file_path, 'r') as last_xyz_file:
with open(last_xyz_file_path, "r") as last_xyz_file:
lines = last_xyz_file.readlines()
return lines
def _evaluate_proc_charges(self, total_nsites: int, proc_charges: list[list[str]]) -> Tuple[
List[Dict[str, float | Any]], List[float], List[int]]:
def _evaluate_proc_charges(
self, total_nsites: int, proc_charges: list[list[str]]
) -> Tuple[List[Dict[str, float | Any]], List[float], List[int]]:
asec_charges = []
thickness = []
@@ -164,9 +142,7 @@ class GaussianInterface(Interface):
f"Number of sites does not match total number of sites."
)
thickness.append(
self._calculate_proc_thickness(charges)
)
thickness.append(self._calculate_proc_thickness(charges))
nsites_ref_mol = len(self.system.molecule[0].atom)
charges = charges[nsites_ref_mol:]
@@ -183,7 +159,12 @@ class GaussianInterface(Interface):
for site in range(len(self.system.molecule[type].atom)):
line = charges.pop(0).split()
if line[0].title() != atomsymb[self.system.molecule[type].atom[site].na].strip():
if (
line[0].title()
!= atomsymb[
self.system.molecule[type].atom[site].na
].strip()
):
raise SyntaxError(
f"Error: Invalid Dice Output. Atom type does not match."
)
@@ -201,13 +182,18 @@ class GaussianInterface(Interface):
)
)
distance = self.system.molecule[0] \
.minimum_distance(new_molecule)
distance = self.system.molecule[0].minimum_distance(new_molecule)
if distance < thickness[-1]:
for atom in new_molecule.atom:
asec_charges.append(
{"lbl": atomsymb[atom.na], "rx": atom.rx, "ry": atom.ry, "rz": atom.rz, "chg": atom.chg}
{
"lbl": atomsymb[atom.na],
"rx": atom.rx,
"ry": atom.ry,
"rz": atom.rz,
"chg": atom.chg,
}
)
mol_count += 1
@@ -230,18 +216,17 @@ class GaussianInterface(Interface):
def _make_gaussian_input_file(self, cycle: int, asec_charges: list[dict]) -> None:
gaussian_input_file_path = Path(
self.step.simulation_dir,
f"step{cycle:02d}",
"qm",
f"asec.gjf"
self.step.simulation_dir, f"step{cycle:02d}", "qm", f"asec.gjf"
)
with open(gaussian_input_file_path, 'w') as gaussian_input_file:
with open(gaussian_input_file_path, "w") as gaussian_input_file:
gaussian_input_file.writelines(
self._generate_gaussian_input(cycle, asec_charges)
)
def _generate_gaussian_input(self, cycle: int, asec_charges: list[dict]) -> list[str]:
def _generate_gaussian_input(
self, cycle: int, asec_charges: list[dict]
) -> list[str]:
gaussian_input = ["%Chk=asec.chk\n"]
if self.step.mem is not None:
@@ -268,7 +253,9 @@ class GaussianInterface(Interface):
gaussian_input.append("\nForce calculation - Cycle number {}\n".format(cycle))
gaussian_input.append("\n")
gaussian_input.append(f"{self.step.gaussian.chgmult[0]},{self.step.gaussian.chgmult[1]}\n")
gaussian_input.append(
f"{self.step.gaussian.chgmult[0]},{self.step.gaussian.chgmult[1]}\n"
)
for atom in self.system.molecule[0].atom:
symbol = atomsymb[atom.na]
@@ -283,7 +270,7 @@ class GaussianInterface(Interface):
for charge in asec_charges:
gaussian_input.append(
"{:>10.5f} {:>10.5f} {:>10.5f} {:>11.8f}\n".format(
charge['rx'], charge['ry'], charge['rz'], charge['chg']
charge["rx"], charge["ry"], charge["rz"], charge["chg"]
)
)
@@ -292,11 +279,7 @@ class GaussianInterface(Interface):
return gaussian_input
def _run_gaussian(self, cycle: int) -> None:
qm_dir = Path(
self.step.simulation_dir,
f"step{(cycle):02d}",
"qm"
)
qm_dir = Path(self.step.simulation_dir, f"step{(cycle):02d}", "qm")
working_dir = os.getcwd()
os.chdir(qm_dir)
@@ -319,7 +302,10 @@ class GaussianInterface(Interface):
"bash",
"-c",
"exec -a {}-step{} {} {}".format(
self.step.gaussian.qmprog, cycle, self.step.gaussian.qmprog, infile
self.step.gaussian.qmprog,
cycle,
self.step.gaussian.qmprog,
infile,
),
]
)
@@ -334,18 +320,16 @@ class GaussianInterface(Interface):
os.chdir(working_dir)
def _run_formchk(self, cycle: int):
qm_dir = Path(
self.step.simulation_dir,
f"step{(cycle):02d}",
"qm"
)
qm_dir = Path(self.step.simulation_dir, f"step{(cycle):02d}", "qm")
work_dir = os.getcwd()
os.chdir(qm_dir)
logger.info("Formatting the checkpoint file... \n")
exit_status = subprocess.call(["formchk", "asec.chk"], stdout=subprocess.DEVNULL)
exit_status = subprocess.call(
["formchk", "asec.chk"], stdout=subprocess.DEVNULL
)
if exit_status != 0:
raise SystemError("Formchk process did not exit properly")
@@ -355,12 +339,7 @@ class GaussianInterface(Interface):
os.chdir(work_dir)
def _read_charges_from_fchk(self, cycle: int):
fchk_file_path = Path(
"simfiles",
f"step{cycle:02d}",
"qm",
"asec.fchk"
)
fchk_file_path = Path("simfiles", f"step{cycle:02d}", "qm", "asec.fchk")
with open(fchk_file_path) as fchk:
fchkfile = fchk.readlines()

View File

@@ -1,4 +1,4 @@
from typing import runtime_checkable, Protocol
from typing import Protocol, runtime_checkable
@runtime_checkable

View File

@@ -1,12 +1,11 @@
from pathlib import Path
import logging
from pathlib import Path
def valid_logger(func):
def wrapper(*args, **kwargs):
logger = args[0]
assert logger._was_set, \
"Logger is not set. Please call set_logger() first."
assert logger._was_set, "Logger is not set. Please call set_logger() first."
return func(*args, **kwargs)
@@ -24,7 +23,7 @@ class Logger:
if self._logger is None:
self._logger = logging.getLogger(logger_name)
def set_logger(self, outfile='run.log', level=logging.INFO, stream=None):
def set_logger(self, outfile="run.log", level=logging.INFO, stream=None):
outfile_path = None
if outfile is not None and stream is None:
outfile_path = Path(outfile)
@@ -57,14 +56,14 @@ class Logger:
def _create_handlers(self, outfile_path: Path, stream):
handlers = []
if outfile_path is not None:
handlers.append(logging.FileHandler(outfile_path, mode='a+'))
handlers.append(logging.FileHandler(outfile_path, mode="a+"))
elif stream is not None:
handlers.append(logging.StreamHandler(stream))
else:
handlers.append(logging.StreamHandler())
for handler in handlers:
handler.setFormatter(logging.Formatter('%(message)s'))
handler.setFormatter(logging.Formatter("%(message)s"))
self._logger.addHandler(handler)
def close(self):

View File

@@ -14,6 +14,7 @@ ANG2BOHR: Final[float] = 1 / BOHR2ANG
####################################### functions ######################################
def weekday_date_time():
return time.strftime("%A, %d %b %Y at %H:%M:%S")
@@ -31,8 +32,8 @@ def compress_files_1mb(path):
if os.path.getsize(file) > 1024 * 1024: ## If bigger than 1MB
filegz = file + ".gz"
try:
with open(file, 'rb') as f_in:
with gzip.open(filegz, 'wb') as f_out:
with open(file, "rb") as f_in:
with gzip.open(filegz, "wb") as f_out:
shutil.copyfileobj(f_in, f_out)
except:
sys.exit("Error: cannot compress file {}".format(file))

View File

@@ -2,34 +2,222 @@
dice_ghost_label = "Xx"
#### Tuple of atom symbols
atomsymb = ( "00",
atomsymb = (
"00",
"H ",
"He",
"Li",
"Be",
"B ",
"C ",
"N ",
"O ",
"F ",
"Ne",
"Na",
"Mg",
"Al",
"Si",
"P ",
"S ",
"Cl",
"Ar",
"K ",
"Ca",
"Sc",
"Ti",
"V ",
"Cr",
"Mn",
"Fe",
"Co",
"Ni",
"Cu",
"Zn",
"Ga",
"Ge",
"As",
"Se",
"Br",
"Kr",
"Rb",
"Sr",
"Y ",
"Zr",
"Nb",
"Mo",
"Tc",
"Ru",
"Rh",
"Pd",
"Ag",
"Cd",
"In",
"Sn",
"Sb",
"Te",
"I ",
"Xe",
"Cs",
"Ba",
"La",
"Ce",
"Pr",
"Nd",
"Pm",
"Sm",
"Eu",
"Gd",
"Tb",
"Dy",
"Ho",
"Er",
"Tm",
"Yb",
"Lu",
"Hf",
"Ta",
"W ",
"Re",
"Os",
"Ir",
"Pt",
"Au",
"Hg",
"Ti",
"Pb",
"Bi",
"Po",
"At",
"Rn",
"Fr",
"Ra",
"Ac",
"Th",
"Pa",
"U ",
"Np",
"Pu",
"Am",
"Cm",
"Bk",
"Cf",
"Es",
"Fm",
"Md",
"No",
"Lr",
dice_ghost_label,
)
"H ", "He",
"Li","Be", "B ","C ","N ","O ","F ","Ne",
"Na","Mg", "Al","Si","P ","S ","Cl","Ar",
"K ","Ca","Sc","Ti","V ","Cr","Mn","Fe","Co","Ni","Cu","Zn","Ga","Ge","As","Se","Br","Kr",
"Rb","Sr","Y ","Zr","Nb","Mo","Tc","Ru","Rh","Pd","Ag","Cd","In","Sn","Sb","Te","I ","Xe",
"Cs","Ba",
"La","Ce","Pr","Nd","Pm","Sm","Eu","Gd","Tb","Dy","Ho","Er","Tm","Yb","Lu",
"Hf","Ta","W ","Re","Os","Ir","Pt","Au","Hg","Ti","Pb","Bi","Po","At","Rn",
"Fr","Ra",
"Ac","Th","Pa","U ","Np","Pu","Am","Cm","Bk","Cf","Es","Fm","Md","No","Lr",
dice_ghost_label )
#### Tuple of atom masses
atommass = ( 0.0,
1.0079, 4.0026,
6.9410,9.0122, 10.811,12.011,14.007,15.999,18.998,20.180,
22.990,24.305, 26.982,28.086,30.974,32.065,35.453,39.948,
39.098,40.078,44.956,47.867,50.942,51.996,54.938,55.845,58.933,58.693,63.546,65.409,69.723,72.640,74.922,78.960,79.904,83.798,
85.468,87.620,88.906,91.224,92.906,95.940,98.000,101.07,102.91,106.42,107.87,112.41,114.82,118.71,121.76,127.60,126.90,131.29,
132.91,137.33,
138.91,140.12,140.91,144.24,145.00,150.36,151.96,157.25,158.93,162.50,164.93,167.26,168.93,173.04,174.97,
178.49,180.95,183.84,186.21,190.23,192.22,195.08,196.97,200.59,204.38,207.20,208.98,209.00,210.00,222.00,
223.00,226.00,
227.00,232.04,231.04,238.03,237.00,244.00,243.00,247.00,247.00,251.00,252.00,257.00,258.00,259.00,262.00,
0.000 )
#### Tuple of atom masses
atommass = (
0.0,
1.0079,
4.0026,
6.9410,
9.0122,
10.811,
12.011,
14.007,
15.999,
18.998,
20.180,
22.990,
24.305,
26.982,
28.086,
30.974,
32.065,
35.453,
39.948,
39.098,
40.078,
44.956,
47.867,
50.942,
51.996,
54.938,
55.845,
58.933,
58.693,
63.546,
65.409,
69.723,
72.640,
74.922,
78.960,
79.904,
83.798,
85.468,
87.620,
88.906,
91.224,
92.906,
95.940,
98.000,
101.07,
102.91,
106.42,
107.87,
112.41,
114.82,
118.71,
121.76,
127.60,
126.90,
131.29,
132.91,
137.33,
138.91,
140.12,
140.91,
144.24,
145.00,
150.36,
151.96,
157.25,
158.93,
162.50,
164.93,
167.26,
168.93,
173.04,
174.97,
178.49,
180.95,
183.84,
186.21,
190.23,
192.22,
195.08,
196.97,
200.59,
204.38,
207.20,
208.98,
209.00,
210.00,
222.00,
223.00,
226.00,
227.00,
232.04,
231.04,
238.03,
237.00,
244.00,
243.00,
247.00,
247.00,
251.00,
252.00,
257.00,
258.00,
259.00,
262.00,
0.000,
)
#### Number of the ghost atom
ghost_number = len(atomsymb) - 1
ghost_number = len(atomsymb) - 1