regtests.py 39.7 KB
Newer Older
1 2 3 4 5 6 7 8
import os
import unittest
import logging
import numpy as np
from cp2kparser import CP2KParser
from nomadcore.unit_conversion.unit_conversion import convert_unit


9 10
def get_result(folder, metaname=None):
    """Get the results from the calculation in the given folder. By default goes through different
11 12 13

    Args:
        folder: The folder relative to the directory of this script where the
14
            parsed calculation resides.
15 16
        metaname(str): Optional quantity to return. If not specified, returns
            the full dictionary of results.
17 18
    """
    dirname = os.path.dirname(__file__)
19
    filename = os.path.join("cp2k_{}".format(VERSION), dirname, folder, "unittest.out")
20
    parser = CP2KParser(filename, None, debug=True, log_level=logging.CRITICAL)
21
    results = parser.parse()
22

23 24
    if metaname is None:
        return results
25
    else:
26
        return results[metaname]
27 28


29
class TestErrors(unittest.TestCase):
30
    """Test error situations which may occur during the parsing.
31 32
    """
    def test_no_file(self):
33 34 35 36
        """File is not no present.
        """
        with self.assertRaises(IOError):
            get_result("errors/no_file", "XC_functional")
37 38

    def test_invalid_file(self):
39 40 41 42
        """Main file is invalid.
        """
        with self.assertRaises(RuntimeError):
            get_result("errors/invalid_file", "XC_functional")
43 44

    def test_invalid_run_type(self):
45 46 47 48 49
        """Unrecognized run type.
        """
        with self.assertRaises(KeyError):
            get_result("errors/invalid_run_type", "XC_functional")

50

51 52 53
class TestUnknownInput(unittest.TestCase):
    """Tests for cases where unknown information is encountered in the parsing.
    """
54
    def test_unknown_version(self):
55 56
        """Test how a new version is handled.
        """
57 58
        get_result("errors/unknown_version", "XC_functional")

59
    def test_unknown_input_keyword(self):
60 61
        """Test how an unknown input keyword is handled.
        """
62 63 64
        get_result("errors/unknown_input_keyword", "XC_functional")

    def test_unknown_input_section(self):
65 66
        """Test unknown input file section.
        """
67 68 69
        get_result("errors/unknown_input_section", "XC_functional")

    def test_unknown_input_section_parameter(self):
70 71
        """
        """
72 73
        get_result("errors/unknown_input_section_parameter", "XC_functional")

74

75
class TestXCFunctional(unittest.TestCase):
76 77
    """Tests that the XC functionals can be properly parsed.
    """
78 79
    def test_pade(self):
        xc = get_result("XC_functional/pade", "XC_functional")
80
        self.assertEqual(xc, "1*LDA_XC_TETER93")
81 82 83

    def test_lda(self):
        xc = get_result("XC_functional/lda", "XC_functional")
84
        self.assertEqual(xc, "1*LDA_XC_TETER93")
85 86 87

    def test_blyp(self):
        xc = get_result("XC_functional/blyp", "XC_functional")
88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108
        self.assertEqual(xc, "1*GGA_C_LYP+1*GGA_X_B88")

    def test_b3lyp(self):
        xc = get_result("XC_functional/b3lyp", "XC_functional")
        self.assertEqual(xc, "1*HYB_GGA_XC_B3LYP")

    def test_olyp(self):
        xc = get_result("XC_functional/olyp", "XC_functional")
        self.assertEqual(xc, "1*GGA_C_LYP+1*GGA_X_OPTX")

    def test_hcth120(self):
        xc = get_result("XC_functional/hcth120", "XC_functional")
        self.assertEqual(xc, "1*GGA_XC_HCTH_120")

    def test_pbe0(self):
        xc = get_result("XC_functional/pbe0", "XC_functional")
        self.assertEqual(xc, "1*HYB_GGA_XC_PBEH")

    def test_pbe(self):
        xc = get_result("XC_functional/pbe", "XC_functional")
        self.assertEqual(xc, "1*GGA_C_PBE+1*GGA_X_PBE")
109 110


111 112 113 114 115 116 117 118 119 120 121 122 123
class TestSCFConvergence(unittest.TestCase):
    """Tests whether the convergence status and number of SCF step can be
    parsed correctly.
    """
    def test_converged(self):
        result = get_result("convergence/converged", "single_configuration_calculation_converged")
        self.assertTrue(result)

    def test_non_converged(self):
        result = get_result("convergence/non_converged", "single_configuration_calculation_converged")
        self.assertFalse(result)


124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146
class TestForceFiles(unittest.TestCase):
    """Tests that different force files that can be output, can actually be
    found and parsed.
    """
    def test_single_point(self):

        result = get_result("force_file/single_point", "atom_forces")
        expected_result = convert_unit(
            np.array([
                [0.00000000, 0.00000000, 0.00000000],
                [0.00000000, 0.00000001, 0.00000001],
                [0.00000001, 0.00000001, 0.00000000],
                [0.00000001, 0.00000000, 0.00000001],
                [-0.00000001, -0.00000001, -0.00000001],
                [-0.00000001, -0.00000001, -0.00000001],
                [-0.00000001, -0.00000001, -0.00000001],
                [-0.00000001, -0.00000001, -0.00000001],
            ]),
            "forceAu"
        )
        self.assertTrue(np.array_equal(result, expected_result))


147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170
class TestSelfInteractionCorrectionMethod(unittest.TestCase):
    """Tests that the self-interaction correction can be properly parsed.
    """
    def test_no(self):
        sic = get_result("sic/no", "self_interaction_correction_method")
        self.assertEqual(sic, "")

    def test_ad(self):
        sic = get_result("sic/ad", "self_interaction_correction_method")
        self.assertEqual(sic, "SIC_AD")

    def test_explicit_orbitals(self):
        sic = get_result("sic/explicit_orbitals", "self_interaction_correction_method")
        self.assertEqual(sic, "SIC_EXPLICIT_ORBITALS")

    def test_mauri_spz(self):
        sic = get_result("sic/mauri_spz", "self_interaction_correction_method")
        self.assertEqual(sic, "SIC_MAURI_SPZ")

    def test_mauri_us(self):
        sic = get_result("sic/mauri_us", "self_interaction_correction_method")
        self.assertEqual(sic, "SIC_MAURI_US")


171 172 173 174 175
class TestStressTensorMethods(unittest.TestCase):
    """Tests that the stress tensor can be properly parsed for different
    calculation methods.
    """
    def test_none(self):
176
        get_result("stress_tensor/none", "section_stress_tensor")
177 178

    def test_analytical(self):
179
        results = get_result("stress_tensor/analytical")
180 181 182 183 184
        method = results["stress_tensor_method"]
        results["stress_tensor"]
        self.assertEqual(method, "Analytical")

    def test_numerical(self):
185
        results = get_result("stress_tensor/numerical")
186 187 188 189 190
        method = results["stress_tensor_method"]
        results["stress_tensor"]
        self.assertEqual(method, "Numerical")

    def test_diagonal_analytical(self):
191
        results = get_result("stress_tensor/diagonal_analytical")
192 193 194 195 196
        method = results["stress_tensor_method"]
        results["stress_tensor"]
        self.assertEqual(method, "Diagonal analytical")

    def test_diagonal_numerical(self):
197
        results = get_result("stress_tensor/diagonal_numerical")
198 199 200 201 202
        method = results["stress_tensor_method"]
        results["stress_tensor"]
        self.assertEqual(method, "Diagonal numerical")


203 204 205
class TestConfigurationPeriodicDimensions(unittest.TestCase):
    """Tests that the self-interaction correction can be properly parsed.
    """
206
    folder = "configuration_periodic_dimensions/"
207 208

    def test_default(self):
209
        result = get_result(self.folder+"default", "configuration_periodic_dimensions")
210 211 212
        self.assertTrue(np.array_equal(result, np.array((True, True, True))))

    def test_none(self):
213
        result = get_result(self.folder+"none", "configuration_periodic_dimensions")
214 215 216
        self.assertTrue(np.array_equal(result, np.array((False, False, False))))

    def test_x(self):
217
        result = get_result(self.folder+"x", "configuration_periodic_dimensions")
218 219 220
        self.assertTrue(np.array_equal(result, np.array((True, False, False))))

    def test_y(self):
221
        result = get_result(self.folder+"y", "configuration_periodic_dimensions")
222 223 224
        self.assertTrue(np.array_equal(result, np.array((False, True, False))))

    def test_z(self):
225
        result = get_result(self.folder+"z", "configuration_periodic_dimensions")
226 227 228
        self.assertTrue(np.array_equal(result, np.array((False, False, True))))

    def test_xy(self):
229
        result = get_result(self.folder+"xy", "configuration_periodic_dimensions")
230 231 232
        self.assertTrue(np.array_equal(result, np.array((True, True, False))))

    def test_xyz(self):
233
        result = get_result(self.folder+"xyz", "configuration_periodic_dimensions")
234 235 236
        self.assertTrue(np.array_equal(result, np.array((True, True, True))))

    def test_xz(self):
237
        result = get_result(self.folder+"xz", "configuration_periodic_dimensions")
238 239 240
        self.assertTrue(np.array_equal(result, np.array((True, False, True))))

    def test_yz(self):
241
        result = get_result(self.folder+"yz", "configuration_periodic_dimensions")
242 243 244
        self.assertTrue(np.array_equal(result, np.array((False, True, True))))


245 246 247
class TestEnergyForce(unittest.TestCase):
    """Tests for a CP2K calculation with RUN_TYPE ENERGY_FORCE.
    """
248 249
    @classmethod
    def setUpClass(cls):
250
        cls.results = get_result("energy_force")
251
        # cls.results.print_summary()
252

253
    def test_energy_total_scf_iteration(self):
254
        result = self.results["energy_total_scf_iteration"]
255
        expected_result = convert_unit(np.array(-32.2320848878), "hartree")
256
        self.assertTrue(np.array_equal(result[0], expected_result))
257

258 259 260 261 262 263 264 265
    def test_program_name(self):
        result = self.results["program_name"]
        self.assertEqual(result, "CP2K")

    def test_program_compilation_host(self):
        result = self.results["program_compilation_host"]
        self.assertEqual(result, "lenovo700")

266 267 268 269
    def test_scf_max_iteration(self):
        result = self.results["scf_max_iteration"]
        self.assertEqual(result, 300)

270 271 272 273
    def test_basis_set(self):
        section_basis_set = self.results["section_basis_set"][0]

        # Basis name
274
        name = section_basis_set["basis_set_name"]
275 276 277
        self.assertEqual(name, "DZVP-GTH-PADE_PW_150.0")

        # Basis kind
278
        kind = section_basis_set["basis_set_kind"]
279 280 281
        self.assertEqual(kind, "wavefunction")

        # Cell dependent basis mapping
282
        cell_basis_mapping = section_basis_set["mapping_section_basis_set_cell_dependent"]
283 284 285
        self.assertEqual(cell_basis_mapping, 0)

        # # Atom centered basis mapping
286
        atom_basis_mapping = section_basis_set["mapping_section_basis_set_atom_centered"]
287 288
        self.assertTrue(np.array_equal(atom_basis_mapping, np.array(8*[0])))

289 290
    def test_scf_threshold_energy_change(self):
        result = self.results["scf_threshold_energy_change"]
291
        self.assertEqual(result, convert_unit(1.00E-07, "hartree"))
292 293 294 295 296

    def test_number_of_spin_channels(self):
        result = self.results["number_of_spin_channels"]
        self.assertEqual(result, 1)

297 298 299 300
    def test_electronic_structure_method(self):
        result = self.results["electronic_structure_method"]
        self.assertEqual(result, "DFT")

301 302 303 304 305 306 307 308 309 310
    def test_energy_change_scf_iteration(self):
        energy_change = self.results["energy_change_scf_iteration"]
        expected_result = convert_unit(np.array(-3.22E+01), "hartree")
        self.assertTrue(np.array_equal(energy_change[0], expected_result))

    def test_energy_XC_scf_iteration(self):
        result = self.results["energy_XC_scf_iteration"]
        expected_result = convert_unit(np.array(-9.4555961214), "hartree")
        self.assertTrue(np.array_equal(result[0], expected_result))

311
    def test_energy_total(self):
312
        result = self.results["energy_total"]
313
        expected_result = convert_unit(np.array(-31.297885372811074), "hartree")
314
        self.assertTrue(np.array_equal(result, expected_result))
315

316 317
    def test_electronic_kinetic_energy(self):
        result = self.results["electronic_kinetic_energy"]
318
        expected_result = convert_unit(np.array(13.31525592466419), "hartree")
319 320
        self.assertTrue(np.array_equal(result, expected_result))

321
    def test_atom_forces(self):
322
        result = self.results["atom_forces"]
323 324 325 326 327 328 329 330 331 332 333
        expected_result = convert_unit(
            np.array([
                [0.00000000, 0.00000000, 0.00000000],
                [0.00000000, 0.00000001, 0.00000001],
                [0.00000001, 0.00000001, 0.00000000],
                [0.00000001, 0.00000000, 0.00000001],
                [-0.00000001, -0.00000001, -0.00000001],
                [-0.00000001, -0.00000001, -0.00000001],
                [-0.00000001, -0.00000001, -0.00000001],
                [-0.00000001, -0.00000001, -0.00000001],
            ]),
334
            "forceAu"
335
        )
336
        self.assertTrue(np.array_equal(result, expected_result))
337

338
    def test_atom_label(self):
339
        atom_labels = self.results["atom_labels"]
340
        expected_labels = np.array(8*["Si"])
Lauri Himanen's avatar
Lauri Himanen committed
341 342 343
        self.assertTrue(np.array_equal(atom_labels, expected_labels))

    def test_simulation_cell(self):
Lauri Himanen's avatar
Lauri Himanen committed
344
        cell = self.results["simulation_cell"]
Lauri Himanen's avatar
Lauri Himanen committed
345 346 347 348 349 350 351 352
        n_vectors = cell.shape[0]
        n_dim = cell.shape[1]
        self.assertEqual(n_vectors, 3)
        self.assertEqual(n_dim, 3)
        expected_cell = convert_unit(np.array([[5.431, 0, 0], [0, 5.431, 0], [0, 0, 5.431]]), "angstrom")
        self.assertTrue(np.array_equal(cell, expected_cell))

    def test_number_of_atoms(self):
Lauri Himanen's avatar
Lauri Himanen committed
353
        n_atoms = self.results["number_of_atoms"]
Lauri Himanen's avatar
Lauri Himanen committed
354 355 356
        self.assertEqual(n_atoms, 8)

    def test_atom_position(self):
357
        atom_position = self.results["atom_positions"]
Lauri Himanen's avatar
Lauri Himanen committed
358 359 360
        expected_position = convert_unit(np.array([4.073023, 4.073023, 1.357674]), "angstrom")
        self.assertTrue(np.array_equal(atom_position[-1, :], expected_position))

361
    def test_x_cp2k_filenames(self):
362
        input_filename = self.results["x_cp2k_input_filename"]
363 364 365
        expected_input = "si_bulk8.inp"
        self.assertTrue(input_filename, expected_input)

366
        bs_filename = self.results["x_cp2k_basis_set_filename"]
367 368 369
        expected_bs = "../BASIS_SET"
        self.assertEqual(bs_filename, expected_bs)

370
        geminal_filename = self.results["x_cp2k_geminal_filename"]
371 372 373
        expected_geminal = "BASIS_GEMINAL"
        self.assertEqual(geminal_filename, expected_geminal)

374
        potential_filename = self.results["x_cp2k_potential_filename"]
375 376 377
        expected_potential = "../GTH_POTENTIALS"
        self.assertEqual(potential_filename, expected_potential)

378
        mm_potential_filename = self.results["x_cp2k_mm_potential_filename"]
379 380 381
        expected_mm_potential = "MM_POTENTIAL"
        self.assertEqual(mm_potential_filename, expected_mm_potential)

382
        coordinate_filename = self.results["x_cp2k_coordinate_filename"]
383 384
        expected_coordinate = "__STD_INPUT__"
        self.assertEqual(coordinate_filename, expected_coordinate)
385

386
    def test_target_multiplicity(self):
387
        multiplicity = self.results["spin_target_multiplicity"]
388 389 390 391 392 393
        self.assertEqual(multiplicity, 1)

    def test_total_charge(self):
        charge = self.results["total_charge"]
        self.assertEqual(charge, 0)

394 395
    def test_section_basis_set_atom_centered(self):
        basis = self.results["section_basis_set_atom_centered"][0]
396 397
        name = basis["basis_set_atom_centered_short_name"]
        number = basis["basis_set_atom_number"]
398 399 400
        self.assertEquals(name, "DZVP-GTH-PADE")
        self.assertEquals(number, 14)

Lauri Himanen's avatar
Lauri Himanen committed
401 402
    def test_section_basis_set_cell_dependent(self):
        basis = self.results["section_basis_set_cell_dependent"][0]
403
        cutoff = basis["basis_set_planewave_cutoff"]
Lauri Himanen's avatar
Lauri Himanen committed
404 405
        self.assertEquals(cutoff, convert_unit(300.0, "hartree"))

406 407
    def test_section_method_atom_kind(self):
        kind = self.results["section_method_atom_kind"][0]
408 409
        self.assertEqual(kind["method_atom_kind_atom_number"], 14)
        self.assertEqual(kind["method_atom_kind_label"], "Si")
410 411 412

    def test_section_method_basis_set(self):
        kind = self.results["section_method_basis_set"][0]
413 414
        self.assertEqual(kind["method_basis_set_kind"], "wavefunction")
        self.assertEqual(kind["number_of_basis_sets_atom_centered"], 1)
415
        self.assertTrue(np.array_equal(kind["mapping_section_method_basis_set_atom_centered"], np.array([[0, 0]])))
416

417 418 419 420 421
    def test_single_configuration_calculation_converged(self):
        result = self.results["single_configuration_calculation_converged"]
        self.assertTrue(result)

    def test_scf_dft_number_of_iterations(self):
422
        result = self.results["number_of_scf_iterations"]
423 424 425 426 427 428 429
        self.assertEqual(result, 10)

    def test_single_configuration_to_calculation_method_ref(self):
        result = self.results["single_configuration_to_calculation_method_ref"]
        self.assertEqual(result, 0)

    def test_single_configuration_calculation_to_system_description_ref(self):
430
        result = self.results["single_configuration_calculation_to_system_ref"]
431 432
        self.assertEqual(result, 0)

433
    def test_stress_tensor(self):
434
        result = self.results["stress_tensor"]
435 436
        expected_result = convert_unit(
            np.array([
437 438 439
                [7.77640934, -0.00000098, -0.00000099],
                [-0.00000098, 7.77640935, -0.00000101],
                [-0.00000099, -0.00000101, 7.77640935],
440 441 442 443 444
            ]),
            "GPa"
        )
        self.assertTrue(np.array_equal(result, expected_result))

445 446
    def test_stress_tensor_eigenvalues(self):
        result = self.results["x_cp2k_stress_tensor_eigenvalues"]
447
        expected_result = convert_unit(np.array([7.77640735, 7.77641033, 7.77641036]), "GPa")
448 449 450 451 452
        self.assertTrue(np.array_equal(result, expected_result))

    def test_stress_tensor_eigenvectors(self):
        result = self.results["x_cp2k_stress_tensor_eigenvectors"]
        expected_result = np.array([
453 454 455
            [0.57490332, -0.79965737, -0.17330395],
            [0.57753686, 0.54662171, -0.60634634],
            [0.57960102, 0.24850110, 0.77608624],
456 457 458 459 460
        ])
        self.assertTrue(np.array_equal(result, expected_result))

    def test_stress_tensor_determinant(self):
        result = self.results["x_cp2k_stress_tensor_determinant"]
461
        expected_result = convert_unit(4.70259243E+02, "GPa^3")
462 463 464 465
        self.assertTrue(np.array_equal(result, expected_result))

    def test_stress_tensor_one_third_of_trace(self):
        result = self.results["x_cp2k_stress_tensor_one_third_of_trace"]
466
        expected_result = convert_unit(7.77640934E+00, "GPa")
467 468
        self.assertTrue(np.array_equal(result, expected_result))

469 470
    def test_program_basis_set_type(self):
        result = self.results["program_basis_set_type"]
471
        self.assertEqual(result, "gaussians")
472 473 474


class TestPreprocessor(unittest.TestCase):
475 476 477
    """Tests that the parser can read input files with preprocessor
    declarations, such as variables.
    """
478
    def test_include(self):
479
        result = get_result("input_preprocessing/include", "x_cp2k_input_GLOBAL.PRINT_LEVEL")
480 481 482
        self.assertEqual(result, "LOW")

    def test_variable(self):
483
        result = get_result("input_preprocessing/variable", "x_cp2k_input_GLOBAL.PROJECT_NAME")
484 485
        self.assertEqual(result, "variable_test")

486
    def test_variable_multiple(self):
487
        result = get_result("input_preprocessing/variable_multiple", "x_cp2k_input_FORCE_EVAL.DFT.MGRID.CUTOFF")
488
        self.assertEqual(result, "50")
489

Lauri Himanen's avatar
Lauri Himanen committed
490
    def test_comments(self):
491
        result = get_result("input_preprocessing/comments", "x_cp2k_input_FORCE_EVAL.DFT.MGRID.CUTOFF")
492
        self.assertEqual(result, "120")
Lauri Himanen's avatar
Lauri Himanen committed
493

494
    def test_tabseparator(self):
495
        result = get_result("input_preprocessing/tabseparator", "x_cp2k_input_FORCE_EVAL.DFT.MGRID.CUTOFF")
496
        self.assertEqual(result, "120")
497

498

499
class TestGeoOpt(unittest.TestCase):
500 501
    """Tests that geometry optimizations are correctly parsed.
    """
502 503
    @classmethod
    def setUpClass(cls):
504
        cls.results = get_result("geo_opt/cg")
505 506 507 508 509

    def test_geometry_optimization_converged(self):
        result = self.results["geometry_optimization_converged"]
        self.assertTrue(result)

510
    def test_number_of_frames_in_sequence(self):
511
        result = self.results["number_of_frames_in_sequence"]
512 513
        self.assertEqual(result, 7)

514 515 516 517 518 519
    def test_frame_sequence_to_sampling_ref(self):
        result = self.results["frame_sequence_to_sampling_ref"]
        self.assertEqual(result, 0)

    def test_frame_sequence_local_frames_ref(self):
        result = self.results["frame_sequence_local_frames_ref"]
520
        expected_result = np.array([0, 1, 2, 3, 4, 5, 6])
521 522
        self.assertTrue(np.array_equal(result, expected_result))

523 524 525 526 527 528 529 530 531 532 533 534 535 536 537 538 539 540 541 542 543 544 545 546 547 548 549 550 551 552 553 554 555 556
    def test_sampling_method(self):
        result = self.results["sampling_method"]
        self.assertEqual(result, "geometry_optimization")

    def test_geometry_optimization_method(self):
        result = self.results["geometry_optimization_method"]
        self.assertEqual(result, "conjugate_gradient")

    def test_geometry_optimization_geometry_change(self):
        result = self.results["geometry_optimization_geometry_change"]
        expected_result = convert_unit(
            0.0010000000,
            "bohr"
        )
        self.assertEqual(result, expected_result)

    def test_geometry_optimization_threshold_force(self):
        result = self.results["geometry_optimization_threshold_force"]
        expected_result = convert_unit(
            0.0010000000,
            "bohr^-1*hartree"
        )
        self.assertEqual(result, expected_result)

    def test_frame_sequence_potential_energy(self):
        result = self.results["frame_sequence_potential_energy"]
        expected_result = convert_unit(
            np.array([
                -17.1534159246,
                -17.1941015290,
                -17.2092321965,
                -17.2097667733,
                -17.2097743028,
                -17.2097743229,
557
                -17.20977820662248,
558 559 560 561 562
            ]),
            "hartree"
        )
        self.assertTrue(np.array_equal(result, expected_result))

563 564 565 566 567 568 569 570 571 572 573 574 575 576 577 578 579 580 581
    def test_atom_positions(self):
        result = self.results["atom_positions"]
        expected_start = convert_unit(
            np.array([
                [12.2353220000, 1.3766420000, 10.8698800000],
                [12.4175775999, 2.2362362573, 11.2616216864],
                [11.9271436933, 1.5723516602, 10.0115134757],
            ]),
            "angstrom"
        )

        expected_end = convert_unit(
            np.array([
                [12.2353220000, 1.3766420000, 10.8698800000],
                [12.4958164689, 2.2307248873, 11.3354322515],
                [11.9975558616, 1.5748085240, 10.0062792262],
            ]),
            "angstrom"
        )
582 583
        result_start = result[0, :, :]
        result_end = result[-1, :, :]
584 585
        self.assertTrue(np.array_equal(result_start, expected_start))
        self.assertTrue(np.array_equal(result_end, expected_end))
586 587


588
class TestGeoOptTrajFormats(unittest.TestCase):
589 590
    """Different trajectory formats in geometry optimization.
    """
591
    def test_xyz(self):
592

593
        result = get_result("geo_opt/geometry_formats/xyz", "atom_positions")
594 595 596 597 598 599 600 601 602 603 604 605 606 607 608 609
        expected_start = convert_unit(
            np.array([
                [12.2353220000, 1.3766420000, 10.8698800000],
                [12.4175624065, 2.2362390825, 11.2616392180],
                [11.9271777126, 1.5723402996, 10.0115089094],
            ]),
            "angstrom"
        )
        expected_end = convert_unit(
            np.array([
                [12.2353220000, 1.3766420000, 10.8698800000],
                [12.4957995882, 2.2307218433, 11.3354453867],
                [11.9975764125, 1.5747996320, 10.0062529540],
            ]),
            "angstrom"
        )
610 611
        result_start = result[0, :, :]
        result_end = result[-1, :, :]
612 613 614 615
        self.assertTrue(np.array_equal(result_start, expected_start))
        self.assertTrue(np.array_equal(result_end, expected_end))

    def test_pdb(self):
616
        result = get_result("geo_opt/geometry_formats/pdb", "atom_positions")
617 618 619 620 621 622 623 624 625 626 627 628 629 630 631 632
        expected_start = convert_unit(
            np.array([
                [12.235, 1.377, 10.870],
                [12.418, 2.236, 11.262],
                [11.927, 1.572, 10.012],
            ]),
            "angstrom"
        )
        expected_end = convert_unit(
            np.array([
                [12.235, 1.377, 10.870],
                [12.496, 2.231, 11.335],
                [11.998, 1.575, 10.006],
            ]),
            "angstrom"
        )
633 634
        result_start = result[0, :, :]
        result_end = result[-1, :, :]
635 636 637 638
        self.assertTrue(np.array_equal(result_start, expected_start))
        self.assertTrue(np.array_equal(result_end, expected_end))

    def test_dcd(self):
639
        result = get_result("geo_opt/geometry_formats/dcd", "atom_positions")
640 641 642 643
        frames = result.shape[0]
        self.assertEqual(frames, 7)


644
class TestGeoOptOptimizers(unittest.TestCase):
645 646
    """Different optimization methods in gemoetry optimization.
    """
647 648 649 650 651 652 653 654
    def test_bfgs(self):
        result = get_result("geo_opt/bfgs", "geometry_optimization_method")
        self.assertEqual(result, "bfgs")

    def test_lbfgs(self):
        result = get_result("geo_opt/lbfgs", "geometry_optimization_method")
        self.assertEqual(result, "bfgs")

655 656

class TestGeoOptTrajectory(unittest.TestCase):
657 658 659
    """Tests that the print settings for geometry optimization are handled
    correctly.
    """
660 661 662 663
    def test_each_and_add_last(self):
        """Test that the EACH and ADD_LAST settings affect the parsing
        correctly.
        """
664
        results = get_result("geo_opt/each")
665 666 667 668 669

        single_conf = results["section_single_configuration_calculation"]
        systems = results["section_system"]

        i_conf = 0
670
        for calc in single_conf:
671
            system_index = calc["single_configuration_calculation_to_system_ref"]
672 673 674
            system = systems[system_index]

            if i_conf == 0 or i_conf == 2 or i_conf == 4:
675
                with self.assertRaises(KeyError):
676
                    pos = system["atom_positions"]
677
            else:
678
                pos = system["atom_positions"]
679 680 681 682 683 684 685 686 687 688 689 690 691 692 693 694 695 696 697 698 699 700 701 702 703 704 705 706 707 708 709 710 711 712 713 714 715 716 717 718 719 720 721 722
                if i_conf == 1:
                    expected_pos = convert_unit(
                        np.array([
                            [12.2353220000, 1.3766420000, 10.8698800000],
                            [12.4618486015, 2.2314871691, 11.3335607388],
                            [11.9990227122, 1.5776813026, 10.0384213366],
                        ]),
                        "angstrom"
                    )
                    self.assertTrue(np.array_equal(pos, expected_pos))
                if i_conf == 3:
                    expected_pos = convert_unit(
                        np.array([
                            [12.2353220000, 1.3766420000, 10.8698800000],
                            [12.4962705528, 2.2308411983, 11.3355758433],
                            [11.9975151486, 1.5746309898, 10.0054430868],
                        ]),
                        "angstrom"
                    )
                    self.assertTrue(np.array_equal(pos, expected_pos))
                if i_conf == 5:
                    expected_pos = convert_unit(
                        np.array([
                            [12.2353220000, 1.3766420000, 10.8698800000],
                            [12.4958168364, 2.2307249171, 11.3354322532],
                            [11.9975556812, 1.5748088251, 10.0062793864],
                        ]),
                        "angstrom"
                    )
                    self.assertTrue(np.array_equal(pos, expected_pos))

                if i_conf == 6:
                    expected_pos = convert_unit(
                        np.array([
                            [12.2353220000, 1.3766420000, 10.8698800000],
                            [12.4958164689, 2.2307248873, 11.3354322515],
                            [11.9975558616, 1.5748085240, 10.0062792262],
                        ]),
                        "angstrom"
                    )
                    self.assertTrue(np.array_equal(pos, expected_pos))

            i_conf += 1

Lauri Himanen's avatar
Lauri Himanen committed
723 724

class TestMD(unittest.TestCase):
725 726
    """Molecular dynamics tests.
    """
Lauri Himanen's avatar
Lauri Himanen committed
727 728
    @classmethod
    def setUpClass(cls):
729
        cls.results = get_result("md/nve")
Lauri Himanen's avatar
Lauri Himanen committed
730 731 732 733 734 735 736 737 738
        cls.temp = convert_unit(
            np.array([
                300.000000000,
                275.075405378,
                235.091633019,
                202.752506973,
                192.266488819,
                201.629598676,
                218.299664775,
739
                230.324748557,
Lauri Himanen's avatar
Lauri Himanen committed
740 741 742 743 744 745 746 747 748 749 750 751 752 753 754 755 756 757 758 759 760 761 762 763 764 765 766 767 768 769 770 771 772 773 774 775 776 777 778 779 780 781 782 783 784 785 786 787 788 789 790 791 792 793
                232.691881533,
                226.146979313,
                213.165337396,
            ]),
            "K"
        )
        cls.cons = convert_unit(
            np.array([
                -34.323271136,
                -34.323245645,
                -34.323206964,
                -34.323183380,
                -34.323187747,
                -34.323208962,
                -34.323227533,
                -34.323233583,
                -34.323230715,
                -34.323227013,
                -34.323224123,
            ]),
            "hartree"
        )
        cls.pot = convert_unit(
            np.array([
                -34.330396471,
                -34.329778993,
                -34.328790653,
                -34.327998978,
                -34.327754290,
                -34.327997890,
                -34.328412394,
                -34.328704052,
                -34.328757407,
                -34.328598255,
                -34.328287038,
            ]),
            "hartree"
        )
        cls.kin = convert_unit(
            np.array([
                0.007125335,
                0.006533348,
                0.005583688,
                0.004815598,
                0.004566544,
                0.004788928,
                0.005184860,
                0.005470470,
                0.005526692,
                0.005371243,
                0.005062914,
            ]),
            "hartree"
        )
Lauri Himanen's avatar
Lauri Himanen committed
794

795 796 797 798 799
    def test_number_of_atoms(self):
        result = self.results["number_of_atoms"]
        expected_result = np.array(11*[6])
        self.assertTrue(np.array_equal(result, expected_result))

800 801 802 803 804 805 806 807 808 809 810 811 812
    def test_simulation_cell(self):
        result = self.results["simulation_cell"]
        self.assertEqual(len(result), 11)
        expected_start = convert_unit(
            np.array([
                [9.853, 0, 0],
                [0, 9.853, 0],
                [0, 0, 9.853],
            ]),
            "angstrom"
        )
        self.assertTrue(np.array_equal(result[0], expected_start))

Lauri Himanen's avatar
Lauri Himanen committed
813
    def test_ensemble_type(self):
814
        result = self.results["ensemble_type"]
Lauri Himanen's avatar
Lauri Himanen committed
815 816
        self.assertEqual(result, "NVE")

817 818 819 820
    def test_sampling_method(self):
        result = self.results["sampling_method"]
        self.assertEqual(result, "molecular_dynamics")

Lauri Himanen's avatar
Lauri Himanen committed
821 822 823 824
    def test_number_of_frames_in_sequence(self):
        result = self.results["number_of_frames_in_sequence"]
        self.assertEqual(result, 11)

825 826 827 828 829 830 831 832 833 834 835 836 837 838 839 840 841 842 843 844 845 846 847 848
    def test_atom_positions(self):
        result = self.results["atom_positions"]
        expected_start = convert_unit(
            np.array([
                [2.2803980000, 9.1465390000, 5.0886960000],
                [1.2517030000, 2.4062610000, 7.7699080000],
                [1.7620190000, 9.8204290000, 5.5284540000],
                [3.0959870000, 9.1070880000, 5.5881860000],
                [0.5541290000, 2.9826340000, 8.0820240000],
                [1.7712570000, 2.9547790000, 7.1821810000],
            ]),
            "angstrom"
        )
        expected_end = convert_unit(
            np.array([
                [2.2916014875, 9.1431763260, 5.0868100688],
                [1.2366834078, 2.4077552776, 7.7630044423],
                [1.6909790671, 9.8235337924, 5.5042564094],
                [3.1130341664, 9.0372111810, 5.6100739746],
                [0.5652070478, 3.0441761067, 8.1734257299],
                [1.8669280879, 2.9877213524, 7.2364955946],
            ]),
            "angstrom"
        )
849 850
        self.assertTrue(np.array_equal(result[0, :], expected_start))
        self.assertTrue(np.array_equal(result[-1, :], expected_end))
851

852 853 854 855 856 857 858 859 860 861 862
    def test_atom_velocities(self):
        result = self.results["atom_velocities"]
        expected_start = convert_unit(
            np.array([
                [0.0000299284, 0.0000082360, -0.0000216368],
                [-0.0001665963, 0.0001143863, -0.0000622640],
                [-0.0005732926, -0.0003112611, -0.0007149779],
                [0.0013083605, -0.0009262219, 0.0006258560],
                [0.0012002313, -0.0003701042, 0.0002810523],
                [0.0002340810, -0.0003388418, 0.0011398583],
            ]),
863
            "bohr*(hbar/hartree)^-1"
864 865 866 867 868 869 870 871 872 873
        )
        expected_end = convert_unit(
            np.array([
                [0.0001600263, -0.0000383308, 0.0000153662],
                [-0.0001269381, -0.0000005151, -0.0000726214],
                [0.0000177093, -0.0003257814, -0.0000257852],
                [-0.0015067045, -0.0001700489, -0.0003651605],
                [0.0000307926, 0.0006886719, 0.0008431321],
                [0.0007424681, 0.0003614127, 0.0005749089],
            ]),
874
            "bohr*(hbar/hartree)^-1"
875 876
        )

877 878
        self.assertTrue(np.array_equal(result[0, :], expected_start))
        self.assertTrue(np.array_equal(result[-1, :], expected_end))
879

880 881
    def test_frame_sequence_potential_energy(self):
        result = self.results["frame_sequence_potential_energy"]
Lauri Himanen's avatar
Lauri Himanen committed
882 883 884 885 886 887 888 889 890 891 892 893 894 895 896 897
        self.assertTrue(np.array_equal(result, self.pot))

    def test_frame_sequence_kinetic_energy(self):
        result = self.results["frame_sequence_kinetic_energy"]
        self.assertTrue(np.array_equal(result, self.kin))

    def test_frame_sequence_conserved_quantity(self):
        result = self.results["frame_sequence_conserved_quantity"]
        self.assertTrue(np.array_equal(result, self.cons))

    def test_frame_sequence_temperature(self):
        result = self.results["frame_sequence_temperature"]
        self.assertTrue(np.array_equal(result, self.temp))

    def test_frame_sequence_time(self):
        result = self.results["frame_sequence_time"]
898 899
        expected_result = convert_unit(
            np.array([
Lauri Himanen's avatar
Lauri Himanen committed
900 901 902 903 904 905 906 907 908 909 910
                0.000000,
                0.500000,
                1.000000,
                1.500000,
                2.000000,
                2.500000,
                3.000000,
                3.500000,
                4.000000,
                4.500000,
                5.000000,
911
            ]),
Lauri Himanen's avatar
Lauri Himanen committed
912
            "fs"
913 914 915
        )
        self.assertTrue(np.array_equal(result, expected_result))

Lauri Himanen's avatar
Lauri Himanen committed
916 917 918 919 920 921 922 923 924 925 926 927 928 929 930 931 932 933 934 935 936 937
    def test_frame_sequence_potential_energy_stats(self):
        result = self.results["frame_sequence_potential_energy_stats"]
        expected_result = np.array([self.pot.mean(), self.pot.std()])
        self.assertTrue(np.array_equal(result, expected_result))

    def test_frame_sequence_kinetic_energy_stats(self):
        result = self.results["frame_sequence_kinetic_energy_stats"]
        expected_result = np.array([self.kin.mean(), self.kin.std()])
        self.assertTrue(np.array_equal(result, expected_result))

    def test_frame_sequence_conserved_quantity_stats(self):
        result = self.results["frame_sequence_conserved_quantity_stats"]
        expected_result = np.array([self.cons.mean(), self.cons.std()])
        self.assertTrue(np.array_equal(result, expected_result))

    def test_frame_sequence_temperature_stats(self):
        result = self.results["frame_sequence_temperature_stats"]
        expected_result = np.array([self.temp.mean(), self.temp.std()])
        self.assertTrue(np.array_equal(result, expected_result))


class TestMDEnsembles(unittest.TestCase):
938 939
    """Different ensembles in MD.
    """
Lauri Himanen's avatar
Lauri Himanen committed
940 941 942
    @classmethod
    def setUpClass(cls):
        cls.pressure = convert_unit(
943
            np.array([
944 945
                -0.192828092559E+04,
                -0.145371071470E+04,
Lauri Himanen's avatar
Lauri Himanen committed
946
                -0.210098903760E+03,
947 948
                0.167260570313E+04,
                0.395562042841E+04,
Lauri Himanen's avatar
Lauri Himanen committed
949
                0.630374855942E+04,
950 951
                0.836906136786E+04,
                0.983216022830E+04,
Lauri Himanen's avatar
Lauri Himanen committed
952
                0.104711540465E+05,
953
                0.102444821550E+05,
Lauri Himanen's avatar
Lauri Himanen committed
954
                0.931695792434E+04,
955
            ]),
Lauri Himanen's avatar
Lauri Himanen committed
956
            "bar"
957
        )
Lauri Himanen's avatar
Lauri Himanen committed
958 959

    def test_nvt(self):
960
        results = get_result("md/nvt")
Lauri Himanen's avatar
Lauri Himanen committed
961 962 963 964
        ensemble = results["ensemble_type"]
        self.assertEqual(ensemble, "NVT")

    def test_npt(self):
965
        results = get_result("md/npt")
Lauri Himanen's avatar
Lauri Himanen committed
966 967 968 969 970 971 972 973 974
        ensemble = results["ensemble_type"]
        self.assertEqual(ensemble, "NPT")

        pressure = results["frame_sequence_pressure"]
        self.assertTrue(np.array_equal(pressure, self.pressure))

        pressure_stats = results["frame_sequence_pressure_stats"]
        expected_pressure_stats = np.array([self.pressure.mean(), self.pressure.std()])
        self.assertTrue(np.array_equal(pressure_stats, expected_pressure_stats))
Lauri Himanen's avatar
Lauri Himanen committed
975

976 977 978 979 980 981 982 983 984 985 986 987 988 989 990 991 992 993 994 995 996 997 998 999 1000 1001 1002 1003 1004 1005 1006 1007 1008 1009 1010 1011
        simulation_cell = results["simulation_cell"]
        expected_cell_start = convert_unit(
            np.array(
                [[
                    6.0000000000,
                    0.0000000000,
                    0.0000000000,
                ], [
                    0.0000000000,
                    6.0000000000,
                    0.0000000000,
                ], [
                    0.0000000000,
                    0.0000000000,
                    6.0000000000,
                ]]),
            "angstrom"
        )
        expected_cell_end = convert_unit(
            np.array(
                [[
                    5.9960617905,
                    -0.0068118798,
                    -0.0102043036,
                ], [
                    -0.0068116027,
                    6.0225574669,
                    -0.0155044063,
                ], [
                    -0.0102048226,
                    -0.0155046726,
                    6.0083072343,
                ]]),
            "angstrom"
        )
        self.assertEqual(simulation_cell.shape[0], 11)
1012 1013
        self.assertTrue(np.array_equal(expected_cell_start, simulation_cell[0, :, :]))
        self.assertTrue(np.array_equal(expected_cell_end, simulation_cell[-1, :, :]))
Lauri Himanen's avatar
Lauri Himanen committed
1014

1015 1016

class TestElectronicStructureMethod(unittest.TestCase):
1017 1018
    """Tests that different methods are recognized correctly.
    """
1019
    def test_mp2(self):
1020
        results = get_result("electronic_structure_method/mp2")
1021 1022 1023 1024
        result = results["electronic_structure_method"]
        self.assertEqual(result, "MP2")

    def test_dft_plus_u(self):
1025
        results = get_result("electronic_structure_method/dft_plus_u")
1026 1027 1028 1029
        result = results["electronic_structure_method"]
        self.assertEqual(result, "DFT+U")

    def test_rpa(self):
1030
        results = get_result("electronic_structure_method/rpa")
1031 1032 1033 1034
        result = results["electronic_structure_method"]
        self.assertEqual(result, "RPA")


1035
if __name__ == '__main__':
1036

1037 1038 1039
    logger = logging.getLogger("cp2kparser")
    logger.setLevel(logging.ERROR)

1040 1041 1042 1043 1044 1045 1046 1047 1048 1049 1050 1051 1052 1053 1054 1055 1056 1057 1058 1059 1060 1061 1062
    VERSIONS = [</