1# -*- coding: utf-8 -*- 2# 3# Copyright (c) 2018, the cclib development team 4# 5# This file is part of cclib (http://cclib.github.io) and is distributed under 6# the terms of the BSD 3-Clause License. 7 8"""Test logfiles with vibration output in cclib""" 9 10import os 11import unittest 12 13from skip import skipForParser, skipForLogfile 14 15 16__filedir__ = os.path.realpath(os.path.dirname(__file__)) 17 18 19class GenericIRTest(unittest.TestCase): 20 """Generic vibrational frequency unittest""" 21 22 # Unit tests should normally give this value for the largest IR intensity. 23 max_IR_intensity = 100 24 25 # Unit tests may give these values for the largest force constant and reduced mass, respectively. 26 max_force_constant = 10.0 27 max_reduced_mass = 6.9 28 29 # reference zero-point correction from Gaussian 16 dvb_ir.out 30 zpve = 0.1771 31 32 def setUp(self): 33 """Initialize the number of vibrational frequencies on a per molecule basis""" 34 self.numvib = 3*len(self.data.atomnos) - 6 35 36 def testbasics(self): 37 """Are basic attributes correct?""" 38 self.assertEqual(self.data.natom, 20) 39 40 def testvibdisps(self): 41 """Are the dimensions of vibdisps consistent with numvib x N x 3""" 42 self.assertEqual(len(self.data.vibfreqs), self.numvib) 43 self.assertEqual(self.data.vibdisps.shape, 44 (self.numvib, len(self.data.atomnos), 3)) 45 46 def testlengths(self): 47 """Are the lengths of vibfreqs and vibirs (and if present, vibsyms, vibfconnsts and vibrmasses) correct?""" 48 self.assertEqual(len(self.data.vibfreqs), self.numvib) 49 if hasattr(self.data, 'vibirs'): 50 self.assertEqual(len(self.data.vibirs), self.numvib) 51 if hasattr(self.data, 'vibsyms'): 52 self.assertEqual(len(self.data.vibsyms), self.numvib) 53 if hasattr(self.data, 'vibfconsts'): 54 self.assertEqual(len(self.data.vibfconsts), self.numvib) 55 if hasattr(self.data, 'vibrmasses'): 56 self.assertEqual(len(self.data.vibrmasses), self.numvib) 57 58 def testfreqval(self): 59 """Is the highest freq value 3630 +/- 200 wavenumber?""" 60 self.assertAlmostEqual(max(self.data.vibfreqs), 3630, delta=200) 61 62 @skipForParser('Psi4', 'Psi cannot print IR intensities') 63 def testirintens(self): 64 """Is the maximum IR intensity 100 +/- 10 km/mol?""" 65 self.assertAlmostEqual(max(self.data.vibirs), self.max_IR_intensity, delta=10) 66 67 @skipForParser('ADF', 'ADF cannot print force constants') 68 @skipForParser('DALTON', 'DALTON cannot print force constants') 69 @skipForParser('GAMESS', 'GAMESS-US cannot print force constants') 70 @skipForParser('GAMESSUK', 'GAMESS-UK cannot print force constants') 71 @skipForParser('Molcas', 'Molcas cannot print force constants') 72 @skipForParser('Molpro', 'Molpro cannot print force constants') 73 @skipForParser('ORCA', 'ORCA cannot print force constants') 74 @skipForParser('Turbomole', 'Turbomole cannot print force constants') 75 @skipForLogfile('Jaguar/Jaguar4.2', 'Data file does not contain force constants') 76 @skipForLogfile('Psi4/Psi4-1.0', 'Data file contains vibrational info with cartesian coordinates') 77 def testvibfconsts(self): 78 """Is the maximum force constant 10. +/- 0.1 mDyn/angstrom?""" 79 self.assertAlmostEqual(max(self.data.vibfconsts), self.max_force_constant, delta=0.1) 80 81 @skipForParser('ADF', 'ADF cannot print reduced masses') 82 @skipForParser('DALTON', 'DALTON cannot print reduced masses') 83 @skipForParser('GAMESSUK', 'GAMESSUK cannot print reduced masses') 84 @skipForParser('Molpro', 'Molpro cannot print reduced masses') 85 @skipForParser('ORCA', 'ORCA cannot print reduced masses') 86 @skipForLogfile('GAMESS/PCGAMESS', 'Data file does not contain reduced masses') 87 @skipForLogfile('Psi4/Psi4-1.0', 'Data file does not contain reduced masses') 88 def testvibrmasses(self): 89 """Is the maximum reduced mass 6.9 +/- 0.1 daltons?""" 90 self.assertAlmostEqual(max(self.data.vibrmasses), self.max_reduced_mass, delta=0.1) 91 92 @skipForParser('Psi3', 'not implemented yet') 93 def testzeropointcorrection(self): 94 """Is the zero-point correction correct?""" 95 self.assertAlmostEqual(self.data.zpve, self.zpve, delta=1.0e-3) 96 97 98class ADFIRTest(GenericIRTest): 99 """Customized vibrational frequency unittest""" 100 101 # ??? 102 def testzeropointcorrection(self): 103 """Is the zero-point correction correct?""" 104 self.assertAlmostEqual(self.data.zpve, self.zpve, delta=1.0e-2) 105 106 107class FireflyIRTest(GenericIRTest): 108 """Customized vibrational frequency unittest""" 109 110 max_IR_intensity = 135 111 # ??? 112 zpve = 0.1935 113 114 115class GaussianIRTest(GenericIRTest): 116 """Customized vibrational frequency unittest""" 117 118 def testvibsyms(self): 119 """Is the length of vibsyms correct?""" 120 self.assertEqual(len(self.data.vibsyms), self.numvib) 121 122 def testzeropointcorrection(self): 123 # reference zero-point correction from dvb_ir.out 124 zpve = 0.1771 125 """Is the zero-point correction correct?""" 126 self.assertAlmostEqual(self.data.zpve, zpve, delta=0.001) 127 128 entropy_places = 6 129 enthalpy_places = 3 130 freeenergy_places = 3 131 132 def testtemperature(self): 133 """Is the temperature 298.15 K?""" 134 self.assertAlmostEqual(298.15, self.data.temperature) 135 136 def testpressure(self): 137 """Is the pressure 1 atm?""" 138 self.assertAlmostEqual(1, self.data.pressure) 139 140 def testentropy(self): 141 """Is the entropy reasonable""" 142 self.assertAlmostEqual(0.0001462623335480945, self.data.entropy, self.entropy_places) 143 144 def testenthalpy(self): 145 """Is the enthalpy reasonable""" 146 self.assertAlmostEqual(-382.12130688525264, self.data.enthalpy, self.enthalpy_places) 147 148 def testfreeenergy(self): 149 """Is the freeenergy reasonable""" 150 self.assertAlmostEqual(-382.164915, self.data.freeenergy, self.freeenergy_places) 151 152 def testfreeenergyconsistency(self): 153 """Does G = H - TS hold""" 154 self.assertAlmostEqual(self.data.enthalpy - self.data.temperature * self.data.entropy, self.data.freeenergy, self.freeenergy_places) 155 156class JaguarIRTest(GenericIRTest): 157 """Customized vibrational frequency unittest""" 158 159 # Jagur outputs vibrational info with cartesian coordinates 160 max_force_constant = 3.7 161 max_reduced_mass = 2.3 162 163 def testvibsyms(self): 164 """Is the length of vibsyms correct?""" 165 self.assertEqual(len(self.data.vibsyms), self.numvib) 166 167 168class MolcasIRTest(GenericIRTest): 169 """Customized vibrational frequency unittest""" 170 171 max_IR_intensity = 65 172 zpve = 0.1783 173 174 entropy_places = 6 175 enthalpy_places = 3 176 freeenergy_places = 3 177 178 def testtemperature(self): 179 """Is the temperature 298.15 K?""" 180 self.assertAlmostEqual(298.15, self.data.temperature) 181 182 def testpressure(self): 183 """Is the pressure 1 atm?""" 184 self.assertAlmostEqual(1, self.data.pressure) 185 186 def testentropy(self): 187 """Is the entropy reasonable""" 188 self.assertAlmostEqual(0.00013403320476271246, self.data.entropy, self.entropy_places) 189 190 def testenthalpy(self): 191 """Is the enthalpy reasonable""" 192 self.assertAlmostEqual(-382.11385, self.data.enthalpy, self.enthalpy_places) 193 194 def testfreeenergy(self): 195 """Is the freeenergy reasonable""" 196 self.assertAlmostEqual(-382.153812, self.data.freeenergy, self.freeenergy_places) 197 198 def testfreeenergyconsistency(self): 199 """Does G = H - TS hold""" 200 self.assertAlmostEqual(self.data.enthalpy - self.data.temperature * self.data.entropy, self.data.freeenergy, self.freeenergy_places) 201 202 203class OrcaIRTest(GenericIRTest): 204 """Customized vibrational frequency unittest""" 205 206 # ORCA has a bug in the intensities for version < 4.0 207 max_IR_intensity = 215 208 zpve = 0.1921 209 210 enthalpy_places = 3 211 entropy_places = 6 212 freeenergy_places = 3 213 214 def testtemperature(self): 215 """Is the temperature 298.15 K?""" 216 self.assertAlmostEqual(298.15, self.data.temperature) 217 218 def testpressure(self): 219 """Is the pressure 1 atm?""" 220 self.assertAlmostEqual(1, self.data.pressure) 221 222 def testenthalpy(self): 223 """Is the enthalpy reasonable""" 224 self.assertAlmostEqual(-381.85224835, self.data.enthalpy, self.enthalpy_places) 225 226 def testentropy(self): 227 """Is the entropy reasonable""" 228 self.assertAlmostEqual(0.00012080325339594164, self.data.entropy, self.entropy_places) 229 230 def testfreeenergy(self): 231 """Is the freeenergy reasonable""" 232 self.assertAlmostEqual(-381.88826585, self.data.freeenergy, self.freeenergy_places) 233 234 def testfreeenergyconsistency(self): 235 """Does G = H - TS hold""" 236 self.assertAlmostEqual(self.data.enthalpy - self.data.temperature * self.data.entropy, self.data.freeenergy, self.freeenergy_places) 237 238 239class QChemIRTest(GenericIRTest): 240 """Customized vibrational frequency unittest""" 241 242 enthalpy_places = 3 243 entropy_places = 6 244 freeenergy_places = 3 245 246 def testtemperature(self): 247 """Is the temperature 298.15 K?""" 248 self.assertEqual(298.15, self.data.temperature) 249 250 def testpressure(self): 251 """Is the pressure 1 atm?""" 252 self.assertAlmostEqual(1, self.data.pressure) 253 254 def testenthalpy(self): 255 """Is the enthalpy reasonable""" 256 self.assertAlmostEqual(0.1871270552135131, self.data.enthalpy, self.enthalpy_places) 257 258 def testentropy(self): 259 """Is the entropy reasonable""" 260 self.assertAlmostEqual(0.00014667348271900577, self.data.entropy, self.entropy_places) 261 262 def testfreeenergy(self): 263 """Is the freeenergy reasonable""" 264 self.assertAlmostEqual(0.14339635634084155, self.data.freeenergy, self.freeenergy_places) 265 266 # Molecular mass of DVB in mD. 267 molecularmass = 130078.25 268 269 def testatommasses(self): 270 """Do the atom masses sum up to the molecular mass (130078.25+-0.1mD)?""" 271 mm = 1000*sum(self.data.atommasses) 272 self.assertAlmostEqual(mm, 130078.25, delta=0.1, msg = "Molecule mass: %f not 130078 +- 0.1mD" % mm) 273 274 def testhessian(self): 275 """Do the frequencies from the Hessian match the printed frequencies?""" 276 277 def testfreeenergyconsistency(self): 278 """Does G = H - TS hold""" 279 self.assertAlmostEqual(self.data.enthalpy - self.data.temperature * self.data.entropy, self.data.freeenergy, self.freeenergy_places) 280 281 282class GamessIRTest(GenericIRTest): 283 """Customized vibrational frequency unittest""" 284 # Molecular mass of DVB in mD. 285 molecularmass = 130078.25 286 enthalpy_places = 3 287 entropy_places = 6 288 freeenergy_places = 3 289 290 def testatommasses(self): 291 """Do the atom masses sum up to the molecular mass (130078.25+-0.1mD)?""" 292 mm = 1000*sum(self.data.atommasses) 293 self.assertAlmostEqual(mm, 130078.25, delta=0.1, msg = "Molecule mass: %f not 130078 +- 0.1mD" % mm) 294 295 def testtemperature(self): 296 """Is the temperature 298.15 K?""" 297 self.assertAlmostEqual(298.15, self.data.temperature) 298 299 def testpressure(self): 300 """Is the pressure 1 atm?""" 301 self.assertAlmostEqual(1, self.data.pressure) 302 303 def testenthalpy(self): 304 """Is the enthalpy reasonable""" 305 self.assertAlmostEqual(-381.86372805188300, self.data.enthalpy, self.enthalpy_places) 306 307 def testentropy(self): 308 """Is the entropy reasonable""" 309 self.assertAlmostEqual(0.00014875961938, self.data.entropy, self.entropy_places) 310 311 def testfreeenergy(self): 312 """Is the freeenergy reasonable""" 313 self.assertAlmostEqual(-381.90808120060200, self.data.freeenergy, self.freeenergy_places) 314 315 def testfreeenergyconsistency(self): 316 """Does G = H - TS hold""" 317 self.assertAlmostEqual(self.data.enthalpy - self.data.temperature * self.data.entropy, self.data.freeenergy, self.freeenergy_places) 318 319 320class Psi4IRTest(GenericIRTest): 321 """Customized vibrational frequency unittest""" 322 323 # RHF is used for Psi4 IR test data instead of B3LYP 324 max_force_constant = 9.37 325 zpve = 0.1917 326 327 328class TurbomoleIRTest(GenericIRTest): 329 """Customized vibrational frequency unittest""" 330 331 # ??? 332 zpve = 0.1725 333 334 335class GenericIRimgTest(unittest.TestCase): 336 """Generic imaginary vibrational frequency unittest""" 337 338 def setUp(self): 339 """Initialize the number of vibrational frequencies on a per molecule basis""" 340 self.numvib = 3*len(self.data.atomnos) - 6 341 342 def testvibdisps(self): 343 """Are the dimensions of vibdisps consistent with numvib x N x 3""" 344 self.assertEqual(self.data.vibdisps.shape, 345 (self.numvib, len(self.data.atomnos), 3)) 346 347 def testlengths(self): 348 """Are the lengths of vibfreqs and vibirs correct?""" 349 self.assertEqual(len(self.data.vibfreqs), self.numvib) 350 self.assertEqual(len(self.data.vibirs), self.numvib) 351 352 def testfreqval(self): 353 """Is the lowest freq value negative?""" 354 self.assertTrue(self.data.vibfreqs[0] < 0) 355 356## def testmaxvibdisps(self): 357## """What is the maximum value of displacement for a H vs a C?""" 358## Cvibdisps = compress(self.data.atomnos==6, self.data.vibdisps, 1) 359## Hvibdisps = compress(self.data.atomnos==1, self.data.vibdisps, 1) 360## self.assertEqual(max(abs(Cvibdisps).flat), 1.0) 361 362 363class GenericRamanTest(unittest.TestCase): 364 """Generic Raman unittest""" 365 366 # This value is in amu. 367 max_raman_intensity = 575 368 369 def setUp(self): 370 """Initialize the number of vibrational frequencies on a per molecule basis""" 371 self.numvib = 3*len(self.data.atomnos) - 6 372 373 def testlengths(self): 374 """Is the length of vibramans correct?""" 375 self.assertEqual(len(self.data.vibramans), self.numvib) 376 377 # The tolerance for this number has been increased, since ORCA 378 # failed to make it inside +/-5, but it would be nice in the future 379 # to determine is it's not too much work whether this is due to 380 # algorithmic differences, or to differences in the input basis set 381 # or coordinates. The first would be OK, but in the second case the 382 # unit test jobs should be made more comparable. With cclib, we first 383 # of all want to succeed in parsing, but would also like to remain 384 # as comparable between programs as possible (for these tests). 385 # Note also that this value is adjusted for Gaussian and DALTON - why? 386 def testramanintens(self): 387 """Is the maximum Raman intensity correct?""" 388 self.assertAlmostEqual(max(self.data.vibramans), self.max_raman_intensity, delta=8) 389 390 # We used to test this, but it seems to vary wildly between 391 # programs... perhaps we could use it if we knew why... 392 #self.assertInside(self.data.vibramans[1], 2.6872, 0.0001) 393 394 def testvibdisps(self): 395 """Is the length and value of vibdisps correct?""" 396 assert hasattr(self.data, "vibdisps") 397 assert len(self.data.vibdisps) == 54 398 399 400class DALTONRamanTest(GenericRamanTest): 401 """Customized Raman unittest""" 402 403 max_raman_intensity = 745 404 405 406class GaussianRamanTest(GenericRamanTest): 407 """Customized Raman unittest""" 408 409 max_raman_intensity = 1066 410 411 412class OrcaRamanTest(GenericRamanTest): 413 """Customized Raman unittest""" 414 415 max_raman_intensity = 1048 416 417 418class QChemRamanTest(GenericRamanTest): 419 """Customized Raman unittest""" 420 421 max_raman_intensity = 588 422 423if __name__=="__main__": 424 425 import sys 426 sys.path.insert(1, os.path.join(__filedir__, "..")) 427 428 from test_data import DataSuite 429 suite = DataSuite(['vib']) 430 suite.testall() 431