HideyoshiNakazone/DicePlayer
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

SetGlobas/MolHandling(fixes)/DicePlayer Translation

Browse Source 
Fixed Numerus functions in the SetGlobals file and MolHandling, translated DicePlayer and tested till, such as inertia_tensor, add_type, center_of_mass_to_origin. standard_orientation; and removed chglevel keyword

Signed-off-by: Vitor Hideyoshi <vitor.h.n.batista@gmail.com>
This commit is contained in:
Vitor Hideyoshi
2021-10-24 00:55:11 -03:00
committed by Vitor Hideyoshi
parent 316acf4761
commit 6017f272b3
5 changed files with 64 additions and 77 deletions
Download Patch File Download Diff File Expand all files Collapse all files

68
DPpack/MolHandling.py
View File

@@ -28,7 +28,7 @@ class System:
def add_type(self,nmols, m): def add_type(self,nmols, m):
self.molecule.append(m) self.molecule.append(m)
self.nmols = nmols self.nmols.append(nmols)
# Função que calcula a distância entre dois centros de massa # Função que calcula a distância entre dois centros de massa
# e por se tratar de uma função de dois atomos não deve ser # e por se tratar de uma função de dois atomos não deve ser
@@ -208,21 +208,17 @@ class Molecule:
def center_of_mass(self): def center_of_mass(self):
com = np.zeros(3) self.com = np.zeros(3)
total_mass = 0.0
for atom in self.atom: for atom in self.atom:
total_mass += atom.mass self.com += atom.mass * np.array([atom.rx, atom.ry, atom.rz])
com += atom.mass * np.array([atom.rx, atom.ry, atom.rz])
com = com / total_mass self.com = self.com / self.total_mass
self.com = com
def center_of_mass_to_origin(self): def center_of_mass_to_origin(self):
com = self.center_of_mass() self.center_of_mass()
for atom in self.atom: for atom in self.atom:
@@ -260,6 +256,30 @@ class Molecule:
return np.array(distances).reshape(dim, dim) return np.array(distances).reshape(dim, dim)
def inertia_tensor(self):
self.center_of_mass()
Ixx = Ixy = Ixz = Iyy = Iyz = Izz = 0.0
for atom in self.atom:
#### Obtain the displacement from the center of mass
dx = atom.rx - self.com[0]
dy = atom.ry - self.com[1]
dz = atom.rz - self.com[2]
#### Update the diagonal components of the tensor
Ixx += atom.mass * (dy**2 + dz**2)
Iyy += atom.mass * (dz**2 + dx**2)
Izz += atom.mass * (dx**2 + dy**2)
#### Update the off-diagonal components of the tensor
Ixy += atom.mass * dx * dy * -1
Ixz += atom.mass * dx * dz * -1
Iyz += atom.mass * dy * dz * -1
return np.array([[Ixx, Ixy, Ixz],
[Ixy, Iyy, Iyz],
[Ixz, Iyz, Izz]])
def eixos(self): def eixos(self):
eixos = np.zeros(3) eixos = np.zeros(3)
@@ -288,36 +308,13 @@ class Molecule:
return eixos return eixos
def inertia_tensor(self):
Ixx = Ixy = Ixz = Iyy = Iyz = Izz = 0.0
for atom in self.atom:
#### Obtain the displacement from the center of mass
dx = atom.rx - self.com[0]
dy = atom.ry - self.com[1]
dz = atom.rz - self.com[2]
#### Update the diagonal components of the tensor
Ixx += atom.mass * (dy**2 + dz**2)
Iyy += atom.mass * (dz**2 + dx**2)
Izz += atom.mass * (dx**2 + dy**2)
#### Update the off-diagonal components of the tensor
Ixy += atom.mass * dx * dy * -1
Ixz += atom.mass * dx * dz * -1
Iyz += atom.mass * dy * dz * -1
return np.array([ [Ixx, Ixy, Ixz],
[Ixy, Iyy, Iyz],
[Ixz, Iyz, Izz] ])
def principal_axes(self): def principal_axes(self):
try: try:
evals, evecs = linalg.eigh(self.inertia_tensor()) evals, evecs = linalg.eigh(self.inertia_tensor())
except: except:
sys.exit("Error: diagonalization of inertia tensor did not converge") sys.exit("Error: diagonalization of inertia tensor did not converge")
return evals, evecs return evals, evecs
def read_position(self): def read_position(self):
@@ -366,7 +363,6 @@ class Molecule:
def standard_orientation(self): def standard_orientation(self):
self.center_of_mass_to_origin() self.center_of_mass_to_origin()
tensor = self.inertia_tensor()
evals, evecs = self.principal_axes() evals, evecs = self.principal_axes()
if round(linalg.det(evecs)) == -1: if round(linalg.det(evecs)) == -1:
@@ -435,6 +431,4 @@ class Atom:
self.chg = chg # Double self.chg = chg # Double
self.eps = eps # Double self.eps = eps # Double
self.sig = sig # Double self.sig = sig # Double
self.mass = atommass[self.na] # Double self.mass = atommass[self.na] # Double

32
DPpack/SetGlobals.py
View File

@@ -1,6 +1,7 @@
import os, sys import os, sys
import shutil import shutil
import textwrap import textwrap
import types
from DPpack.MolHandling import * from DPpack.MolHandling import *
from DPpack.PTable import * from DPpack.PTable import *
@@ -179,8 +180,8 @@ class Internal:
except ValueError: except ValueError:
sys.exit(err) sys.exit(err)
elif key in ('level', 'chglevel'): elif key == 'level':
setattr(self.gaussian, key, value[0]) setattr(self.gaussian, key, value)
elif key in ('gmiddle', 'gbottom'): elif key in ('gmiddle', 'gbottom'):
setattr(self.gaussian, key, value[0]) setattr(self.gaussian, key, value[0])
@@ -243,9 +244,6 @@ class Internal:
if self.gaussian.pop != "chelpg" and (self.player.ghosts == "yes" or self.player.lps == "yes"): if self.gaussian.pop != "chelpg" and (self.player.ghosts == "yes" or self.player.lps == "yes"):
sys.exit("Error: ghost atoms or lone pairs only available with 'pop = chelpg')") sys.exit("Error: ghost atoms or lone pairs only available with 'pop = chelpg')")
if self.gaussian.chglevel == None:
self.gaussian.chglevel = self.gaussian.level
## Check only if QM program is Molcas: ## Check only if QM program is Molcas:
# if self.player.qmprog == "molcas": # if self.player.qmprog == "molcas":
@@ -636,7 +634,7 @@ class Internal:
pass pass
### I still have to talk with Herbet about
def populate_asec_vdw(self, cycle): def populate_asec_vdw(self, cycle):
asec_charges = [] # (rx, ry, rz, chg) asec_charges = [] # (rx, ry, rz, chg)
@@ -651,17 +649,17 @@ class Internal:
nsitesref = len(self.system.molecule[0]) + len(ghost_atoms) + len(lp_atoms) nsitesref = len(self.system.molecule[0]) + len(ghost_atoms) + len(lp_atoms)
nsites_total = dice['nmol'][0] * nsitesref nsites_total = self.dice.nmol[0] * nsitesref
for i in range(1, len(dice['nmol'])): for i in range(1, len(self.dice.nmol)):
nsites_total += dice['nmol'][i] * len(molecules[i]) nsites_total += self.dice.nmol[i] * len(self.system.molecule[i])
thickness = [] thickness = []
picked_mols = [] picked_mols = []
for proc in range(1, player['nprocs'] + 1): ## Run over folders for proc in range(1, self.player.nprocs + 1): ## Run over folders
path = "step{:02d}".format(cycle) + os.sep + "p{:02d}".format(proc) path = "step{:02d}".format(cycle) + os.sep + "p{:02d}".format(proc)
file = path + os.sep + dice['outname'] + ".xyz" file = path + os.sep + self.dice.outname + ".xyz"
if not os.path.isfile(file): if not os.path.isfile(file):
sys.exit("Error: cannot find file {}".format(file)) sys.exit("Error: cannot find file {}".format(file))
try: try:
@@ -677,23 +675,23 @@ class Internal:
box = xyzfile.pop(0).split()[-3:] box = xyzfile.pop(0).split()[-3:]
box = [ float(box[0]), float(box[1]), float(box[2]) ] box = [ float(box[0]), float(box[1]), float(box[2]) ]
sizes = sizes_of_molecule(molecules[0]) sizes = self.system.molecule[0].sizes_of_molecule()
thickness.append( min([ (box[0] - sizes[0])/2, (box[1] - sizes[1])/2, thickness.append( min([ (box[0] - sizes[0])/2, (box[1] - sizes[1])/2,
(box[2] - sizes[2])/2 ]) ) (box[2] - sizes[2])/2 ]) )
xyzfile = xyzfile[nsitesref:] ## Skip the first (reference) molecule xyzfile = xyzfile[nsitesref:] ## Skip the first (reference) molecule
mol_count = 0 mol_count = 0
for type in range(len(dice['nmol'])): ## Run over types of molecules for type in range(len(self.dice.nmol)): ## Run over types of molecules
if type == 0: if type == 0:
nmols = dice['nmol'][0] - 1 nmols = self.dice.nmol[0] - 1
else: else:
nmols = dice['nmol'][type] nmols = self.dice.nmol[type]
for mol in range(nmols): ## Run over molecules of each type for mol in range(nmols): ## Run over molecules of each type
new_molecule = [] new_molecule = []
for site in range(len(molecules[type])): ## Run over sites of each molecule for site in range(len(self.system.molecule[types].atoms)): ## Run over sites of each molecule
new_molecule.append({}) new_molecule.append({})
line = xyzfile.pop(0).split() line = xyzfile.pop(0).split()
@@ -831,8 +829,6 @@ class Internal:
self.gmiddle = None # In each case, if a filename is given, its content will be self.gmiddle = None # In each case, if a filename is given, its content will be
self.gbottom = None # inserted in the gaussian input self.gbottom = None # inserted in the gaussian input
self.pop = "chelpg" self.pop = "chelpg"
self.chglevel = None
self.level = None self.level = None
# class Molcas: # class Molcas:

7
diceplayer.py
View File

@@ -95,6 +95,7 @@ if __name__ == '__main__':
#### Bring the molecules to standard orientation and prints info about them #### Bring the molecules to standard orientation and prints info about them
for i in range(len(internal.system.molecule)): for i in range(len(internal.system.molecule)):
internal.outfile.write("\nMolecule type {}:\n\n".format(i + 1)) internal.outfile.write("\nMolecule type {}:\n\n".format(i + 1))
internal.system.molecule[i].print_mol_info(internal.outfile) internal.system.molecule[i].print_mol_info(internal.outfile)
internal.outfile.write(" Translating and rotating molecule to standard orientation...") internal.outfile.write(" Translating and rotating molecule to standard orientation...")
@@ -134,9 +135,9 @@ if __name__ == '__main__':
# #Molcas.read_forces("grad_hessian.dat") # #Molcas.read_forces("grad_hessian.dat")
# #Molcas.read_hessian("grad_hessian.dat") # #Molcas.read_hessian("grad_hessian.dat")
# #### ####
# #### Start the iterative process #### Start the iterative process
# #### ####
# for cycle in range(internal.player.cyc, internal.player.cyc + internal.player.maxcyc): # for cycle in range(internal.player.cyc, internal.player.cyc + internal.player.maxcyc):

25
run.log
View File

@@ -2,10 +2,10 @@
############# Welcome to DICEPLAYER version 1.0 ############# ############# Welcome to DICEPLAYER version 1.0 #############
########################################################################################## ##########################################################################################
Your python version is 3.8.8 (default, Apr 13 2021, 19:58:26) Your python version is 3.8.12 (default, Oct 7 2021, 05:40:48)
[GCC 7.3.0] [GCC 11.1.0]
Program started on Saturday, 25 Sep 2021 at 15:24:31 Program started on Sunday, 24 Oct 2021 at 00:46:13
Environment variables: Environment variables:
OMP_STACKSIZE = Not set OMP_STACKSIZE = Not set
@@ -50,7 +50,6 @@ title = Diceplayer run
GAUSSIAN variables being used in this run: GAUSSIAN variables being used in this run:
------------------------------------------------------------------------------------------ ------------------------------------------------------------------------------------------
chglevel = MP2/aug-cc-pVTZ
chgmult = 0 1 chgmult = 0 1
level = MP2/aug-cc-pVTZ level = MP2/aug-cc-pVTZ
pop = chelpg pop = chelpg
@@ -139,11 +138,11 @@ Molecule type 1:
Translating and rotating molecule to standard orientation... Done Translating and rotating molecule to standard orientation... Done
New values: New values:
Center of mass = ( -0.0000 , 0.0000 , 0.0000 ) Center of mass = ( -0.0000 , -0.0000 , -0.0000 )
Moments of inertia = 3.144054E+02 2.801666E+03 3.027366E+03 Moments of inertia = 3.144054E+02 2.801666E+03 3.027366E+03
Major principal axis = ( 1.000000 , 0.000000 , 0.000000 ) Major principal axis = ( 1.000000 , 0.000000 , 0.000000 )
Inter principal axis = ( -0.000000 , 1.000000 , -0.000000 ) Inter principal axis = ( -0.000000 , 1.000000 , -0.000000 )
Minor principal axis = ( -0.000000 , 0.000000 , 1.000000 ) Minor principal axis = ( 0.000000 , 0.000000 , 1.000000 )
Characteristic lengths = ( 11.97 , 5.27 , 2.99 ) Characteristic lengths = ( 11.97 , 5.27 , 2.99 )
Total mass = 226.28 au Total mass = 226.28 au
Total charge = -0.0000 e Total charge = -0.0000 e
@@ -165,14 +164,14 @@ Molecule type 2:
Translating and rotating molecule to standard orientation... Done Translating and rotating molecule to standard orientation... Done
New values: New values:
Center of mass = ( 1.6067 , -1.0929 , 0.0026 ) Center of mass = ( 0.0000 , -0.0000 , 0.0000 )
Moments of inertia = 1.561638E+02 6.739391E+02 8.270247E+02 Moments of inertia = 1.205279E+02 2.859254E+02 4.033761E+02
Major principal axis = ( 0.899747 , -0.436411 , 0.000541 ) Major principal axis = ( 1.000000 , -0.000000 , 0.000000 )
Inter principal axis = ( 0.436411 , 0.899747 , -0.001219 ) Inter principal axis = ( 0.000000 , 1.000000 , -0.000000 )
Minor principal axis = ( 0.000045 , 0.001333 , 0.999999 ) Minor principal axis = ( -0.000000 , 0.000000 , 1.000000 )
Characteristic lengths = ( 5.67 , 5.05 , 1.51 ) Characteristic lengths = ( 5.67 , 5.13 , 1.51 )
Total mass = 112.20 au Total mass = 112.20 au
Total charge = -0.0000 e Total charge = -0.0000 e
Dipole moment = ( 1.8148 , 1.4689 , -0.0016 ) Total = 2.3347 Debye Dipole moment = ( 1.7575 , 1.5369 , -0.0000 ) Total = 2.3347 Debye
========================================================================================== ==========================================================================================

9
run.log.backup
View File

@@ -2,10 +2,10 @@
############# Welcome to DICEPLAYER version 1.0 ############# ############# Welcome to DICEPLAYER version 1.0 #############
########################################################################################## ##########################################################################################
Your python version is 3.8.8 (default, Apr 13 2021, 19:58:26) Your python version is 3.8.12 (default, Oct 7 2021, 05:40:48)
[GCC 7.3.0] [GCC 11.1.0]
Program started on Saturday, 25 Sep 2021 at 15:23:44 Program started on Sunday, 24 Oct 2021 at 00:44:02
Environment variables: Environment variables:
OMP_STACKSIZE = Not set OMP_STACKSIZE = Not set
@@ -50,7 +50,6 @@ title = Diceplayer run
GAUSSIAN variables being used in this run: GAUSSIAN variables being used in this run:
------------------------------------------------------------------------------------------ ------------------------------------------------------------------------------------------
chglevel = MP2/aug-cc-pVTZ
chgmult = 0 1 chgmult = 0 1
level = MP2/aug-cc-pVTZ level = MP2/aug-cc-pVTZ
pop = chelpg pop = chelpg
@@ -176,5 +175,3 @@ Molecule type 2:
Dipole moment = ( 1.8148 , 1.4689 , -0.0016 ) Total = 2.3347 Debye Dipole moment = ( 1.8148 , 1.4689 , -0.0016 ) Total = 2.3347 Debye
========================================================================================== ==========================================================================================
Starting the iterative process.