# Path: molmod/test/common.py
# class BaseTestCase(unittest.TestCase):
#     def assertArraysEqual(self, a, b):
#         self.assertEqual(a.shape, b.shape, "The array shapes do not match.")
#         assert (a==b).all(), "The array values do not match."
# 
#     def assertArrayConstant(self, arr, const):
#         assert (arr==const).all(), "Some/All array values do not match the constant."
# 
#     def assertArraysAlmostEqual(self, a, b, threshold=1e-6, mean=False, doabs=False):
#         self.assertEqual(a.shape, b.shape, "The array shapes do not match.")
#         if a.shape == (0,):
#             return
#         if mean:
#             abserr = abs(a-b).mean()
#             oom = 0.5*(abs(a).mean()+abs(b).mean())
#         else:
#             abserr = abs(a-b).max()
#             oom = 0.5*(abs(a).max()+abs(b).max())
#         if doabs:
#             assert abserr <= threshold, "The absolute error is larger than the given threshold: %5.3e > %5.3e" % (abserr, threshold)
#         else:
#             relerr = abserr/oom
#             assert relerr <= threshold, "The relative error is larger than the given threshold: %5.3e > %5.3e" % (relerr, threshold)
# 
#     def assertArrayAlmostConstant(self, arr, const, relerr_threshold):
#         error = abs(arr-const).max()
#         oom = const
#         relerr = error/oom
#         assert relerr <= relerr_threshold, "The relative error is larger than the given threshold: %5.3e > %5.3e" % (relerr, relerr_threshold)
# 
#     def assertArrayAlmostZero(self, arr, abserr_threshold):
#         abserr = abs(arr).max()
#         assert abserr <= abserr_threshold, "The absolute error is larger than the given threshold: %5.3e > %5.3e" % (abserr, abserr_threshold)


from builtins import range
from molmod.test.common import BaseTestCase
from molmod import *
import unittest
import numpy as np
        #print hessian
        #print ((numerical_hessian - hessian)*diagonal_mask)
        error = sum(((numerical_hessian - hessian)*diagonal_mask).ravel()**2)
        reference = sum((numerical_hessian*diagonal_mask).ravel()**2)
        self.assertAlmostEqual(error, 0.0, 3, "2b) The diagonal blocks of the analytical hessian are incorrect: % 12.8f / %12.8f" % (error, reference))

        error = sum(((numerical_hessian - hessian)*off_diagonal_mask).ravel()**2)
        reference = sum((numerical_hessian*off_diagonal_mask).ravel()**2)
        self.assertAlmostEqual(error, 0.0, 3, "2b) The off-diagonal blocks of the analytical hessian are incorrect: % 12.8f / %12.8f" % (error, reference))

        # 2c) test the analytical hessian in another way
        numerical_hessian = np.zeros(hessian.shape, float)
        for atom in range(len(coordinates)):
            for index in range(3):
                delta_coordinates = coordinates.copy()
                delta_coordinates[atom,index] += delta
                ff.update_coordinates(delta_coordinates)
                numerical_hessian[atom,index,:,:] = (ff.gradient() - gradient) / delta

        #print  numerical_hessian

        error = sum(((numerical_hessian - hessian)*diagonal_mask).ravel()**2)
        reference = sum((numerical_hessian*diagonal_mask).ravel()**2)
        self.assertAlmostEqual(error, 0.0, 3, "2b) The diagonal blocks of the analytical hessian are incorrect: % 12.8f / %12.8f" % (error, reference))

        error = sum(((numerical_hessian - hessian)*off_diagonal_mask).ravel()**2)
        reference = sum((numerical_hessian*off_diagonal_mask).ravel()**2)
        self.assertAlmostEqual(error, 0.0, 3, "2b) The off-diagonal blocks of the analytical hessian are incorrect: % 12.8f / %12.8f" % (error, reference))



# Path: molmod/test/common.py
# class BaseTestCase(unittest.TestCase):
#     def assertArraysEqual(self, a, b):
#         self.assertEqual(a.shape, b.shape, "The array shapes do not match.")
#         assert (a==b).all(), "The array values do not match."
# 
#     def assertArrayConstant(self, arr, const):
#         assert (arr==const).all(), "Some/All array values do not match the constant."
# 
#     def assertArraysAlmostEqual(self, a, b, threshold=1e-6, mean=False, doabs=False):
#         self.assertEqual(a.shape, b.shape, "The array shapes do not match.")
#         if a.shape == (0,):
#             return
#         if mean:
#             abserr = abs(a-b).mean()
#             oom = 0.5*(abs(a).mean()+abs(b).mean())
#         else:
#             abserr = abs(a-b).max()
#             oom = 0.5*(abs(a).max()+abs(b).max())
#         if doabs:
#             assert abserr <= threshold, "The absolute error is larger than the given threshold: %5.3e > %5.3e" % (abserr, threshold)
#         else:
#             relerr = abserr/oom
#             assert relerr <= threshold, "The relative error is larger than the given threshold: %5.3e > %5.3e" % (relerr, threshold)
# 
#     def assertArrayAlmostConstant(self, arr, const, relerr_threshold):
#         error = abs(arr-const).max()
#         oom = const
#         relerr = error/oom
#         assert relerr <= relerr_threshold, "The relative error is larger than the given threshold: %5.3e > %5.3e" % (relerr, relerr_threshold)
# 
#     def assertArrayAlmostZero(self, arr, abserr_threshold):
#         abserr = abs(arr).max()
#         assert abserr <= abserr_threshold, "The absolute error is larger than the given threshold: %5.3e > %5.3e" % (abserr, abserr_threshold)


# Path: molmod/test/test_pairff.py
from builtins import range
from molmod.test.common import BaseTestCase
from molmod import *
import unittest
import numpy as np

        #print hessian
        #print ((numerical_hessian - hessian)*diagonal_mask)
        error = sum(((numerical_hessian - hessian)*diagonal_mask).ravel()**2)
        reference = sum((numerical_hessian*diagonal_mask).ravel()**2)
        self.assertAlmostEqual(error, 0.0, 3, "2b) The diagonal blocks of the analytical hessian are incorrect: % 12.8f / %12.8f" % (error, reference))

        error = sum(((numerical_hessian - hessian)*off_diagonal_mask).ravel()**2)
        reference = sum((numerical_hessian*off_diagonal_mask).ravel()**2)
        self.assertAlmostEqual(error, 0.0, 3, "2b) The off-diagonal blocks of the analytical hessian are incorrect: % 12.8f / %12.8f" % (error, reference))

        # 2c) test the analytical hessian in another way
        numerical_hessian = np.zeros(hessian.shape, float)
        for atom in range(len(coordinates)):
            for index in range(3):
                delta_coordinates = coordinates.copy()
                delta_coordinates[atom,index] += delta
                ff.update_coordinates(delta_coordinates)
                numerical_hessian[atom,index,:,:] = (ff.gradient() - gradient) / delta

        #print  numerical_hessian

        error = sum(((numerical_hessian - hessian)*diagonal_mask).ravel()**2)
        reference = sum((numerical_hessian*diagonal_mask).ravel()**2)
        self.assertAlmostEqual(error, 0.0, 3, "2b) The diagonal blocks of the analytical hessian are incorrect: % 12.8f / %12.8f" % (error, reference))

        error = sum(((numerical_hessian - hessian)*off_diagonal_mask).ravel()**2)
        reference = sum((numerical_hessian*off_diagonal_mask).ravel()**2)
        self.assertAlmostEqual(error, 0.0, 3, "2b) The off-diagonal blocks of the analytical hessian are incorrect: % 12.8f / %12.8f" % (error, reference))



class CoulombFFTestCase(BaseTestCase):