Skip to content

Commit

Permalink
finished tabulated potential norm and testing
Browse files Browse the repository at this point in the history
  • Loading branch information
jrudz committed Sep 24, 2024
1 parent 7647721 commit 0af786c
Show file tree
Hide file tree
Showing 3 changed files with 157 additions and 100 deletions.
66 changes: 39 additions & 27 deletions src/nomad_simulations/schema_packages/force_field.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -251,49 +251,61 @@ def normalize(self, archive, logger) -> None:
- np.min(self.energies.to('kJ').magnitude * MOL)
)
if np.all([x < tol for x in energies_diff]):
print('consistent energies')
logger.warning(
f'Forces were generated from energies in {self.name},'
f' but appear to be roughly consistent, with rtol={tol}. '
f'Tabulated forces were generated from energies in {self.name},'
f'with consistency errors less than tol={tol}. '
)
else:
# logger.warning('inconsistent energies')
print('inconsistent energies')
logger.warning(
f'Unable to derive tabulated forces from energies in {self.name},'
f'consistency errors were greater than tol={tol}.'
)
self.forces = None
# print(energies.to('kJ').magnitude * MOL)
# print(
# self.energies.to('kJ').magnitude * MOL
# - np.min(self.energies.to('kJ').magnitude * MOL)
# )
# print(energies_diff)
except Exception:
pass
except ValueError as e:
logger.warning(
f'Unable to derive tabulated forces from energies in {self.name},'
f'Unkown error occurred in derivation: {e}'
)

if self.forces is not None and self.energies is None:
print('in gen energies')
try:
# generated energies from forces numerically using spline
self.energies = self.compute_energies(
self.bins, self.forces, smoothing_factor=smoothing_factor
self.bins.magnitude,
self.forces.magnitude,
smoothing_factor=smoothing_factor,
)
# re-derive forces to check consistency of the energies
forces = self.compute_forces(
self.bins, self.energies, smoothing_factor=smoothing_factor
forces = (
self.compute_forces(
self.bins.magnitude,
self.energies.magnitude,
smoothing_factor=smoothing_factor,
)
* ureg.J
/ self.bins.units
)

forces_diff = forces.magnitude * 1000.0 * MOL - (
self.energies.to('kJ').magnitude * 1000.0 * MOL
- np.min(self.energies.to('kJ').magnitude * 1000.0 * MOL)
forces_diff = forces.to(f'kJ/{self.bins.units}').magnitude * MOL - (
self.forces.to(f'kJ/{self.bins.units}').magnitude * MOL
)
if np.all([x < tol for x in forces_diff]):
print('consistent forces')
# logger.warning(
# f'Energies were generated from forces in {self.name},'
# f'but appear to be consistent with rtol={rtol}. '
# )
logger.warning(
f'Tabulated energies were generated from forces in {self.name},'
f'with consistency errors less than tol={tol}. '
)
else:
logger.warning(
f'Unable to derive tabulated energies from forces in {self.name},'
f'consistency errors were greater than tol={tol}.'
)
self.energies = None
except Exception:
pass
except ValueError as e:
logger.warning(
f'Unable to derive tabulated energies from forces in {self.name},'
f'Unkown error occurred in derivation: {e}'
)


class BondPotential(Potential):
Expand Down Expand Up @@ -534,7 +546,7 @@ def normalize(self, archive, logger) -> None:
logger.warning('Incorrect functional form set for HarmonicAngle.')


class TabulatedAngle(BondPotential, TabulatedPotential):
class TabulatedAngle(AnglePotential, TabulatedPotential):
"""
Section containing information about a tabulated bond potential. The value of the potential and/or force
is stored for a set of corresponding bin distances.
Expand Down
2 changes: 0 additions & 2 deletions src/nomad_simulations/schema_packages/outputs.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -66,8 +66,6 @@ class Outputs(ArchiveSection):
information if the output `is_derived` from another output section or directly parsed from the simulation output files.
"""

# TODO add time quantities

normalizer_level = 2

model_system_ref = Quantity(
Expand Down
189 changes: 118 additions & 71 deletions tests/test_force_field.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -318,6 +318,28 @@ def populate_parameters(sec_potential, parameters):
7.19013405e-21,
8.38572215e-21,
],
'forces': [
1.3216958010639183e-10,
1.0784837919852664e-10,
8.347183948070384e-11,
5.903996095292366e-11,
3.521528791034704e-11,
1.0674227406164319e-11,
-1.3922171907975959e-11,
-3.857391003207359e-11,
-6.32809869661286e-11,
-8.804340271014125e-11,
-1.1218967953067358e-10,
-1.3706127725108878e-10,
-1.6198821378146133e-10,
-1.8697048912179138e-10,
-2.1200810327207885e-10,
-2.3710105623232373e-10,
-2.6156893683081195e-10,
-2.867710717674596e-10,
-3.1202854551406456e-10,
-3.373413580706269e-10,
],
'name': 'TabulatedBond',
'type': 'bond',
'functional_form': 'tabulated',
Expand Down Expand Up @@ -454,51 +476,73 @@ def populate_parameters(sec_potential, parameters):
'particle_indices': [[0, 1, 2], [3, 4, 5]],
'particle_labels': [['O', 'H', 'H'], ['O', 'H', 'H']],
'bins': [
7.600e-11,
7.970e-11,
8.340e-11,
8.710e-11,
9.070e-11,
9.440e-11,
9.810e-11,
1.018e-10,
1.055e-10,
1.092e-10,
1.128e-10,
1.165e-10,
1.202e-10,
1.239e-10,
1.276e-10,
1.313e-10,
1.349e-10,
1.386e-10,
1.423e-10,
1.460e-10,
92.99105014973262,
94.19727699638311,
95.40350384303362,
96.60973068968411,
97.81595810929241,
99.02218495594292,
100.22841180259343,
101.43463864924392,
102.64086549589443,
103.84709234254493,
105.05331976215322,
106.25954660880373,
107.46577345545423,
108.67200030210473,
109.87822714875524,
111.08445399540572,
112.29068141501403,
113.49690826166454,
114.70313510831502,
115.90936195496555,
],
'energies': [
1.32311751e-21,
8.81248069e-22,
5.28549577e-22,
2.65354139e-22,
9.18278089e-23,
8.30269520e-24,
1.47787975e-23,
1.11422170e-22,
2.98564919e-22,
5.76539155e-22,
9.45178822e-22,
1.40498208e-21,
1.95611499e-21,
2.59857754e-21,
3.33286791e-21,
4.15882003e-21,
5.07693206e-21,
6.08737007e-21,
7.19013405e-21,
8.38572215e-21,
2.483908793147618e-21,
1.9871270690593166e-21,
1.5455433097527105e-21,
1.1591575152954348e-21,
8.279695415419479e-22,
5.51979703171419e-22,
3.311878298023737e-22,
1.6559392146862525e-22,
5.519797819553163e-23,
0.0,
0.0,
5.519797819553239e-23,
1.6559392146862525e-22,
3.311878298023748e-22,
5.519797031714216e-22,
8.2796954154194335e-22,
1.1591575152954338e-21,
1.545543309752715e-21,
1.98712706905931e-21,
2.483908793147623e-21,
],
'name': 'TabulatedBond',
'type': 'bond',
'forces': [
4.347281042004185e-22,
3.8896725108092945e-22,
3.432063979614403e-22,
2.974455448419512e-22,
2.516846920122808e-22,
2.0592383889279167e-22,
1.6016298577330257e-22,
1.1440213265381345e-22,
6.864127953432432e-23,
2.288042641483519e-23,
-2.288042641483519e-23,
-6.864127953432432e-23,
-1.1440213265381345e-22,
-1.6016298577330257e-22,
-2.0592383889279167e-22,
-2.516846920122808e-22,
-2.974455448419512e-22,
-3.432063979614403e-22,
-3.8896725108092945e-22,
-4.347281042004185e-22,
],
'name': 'TabulatedAngle',
'type': 'angle',
'functional_form': 'tabulated',
}
data_tabulated_angle = (
Expand Down Expand Up @@ -531,30 +575,30 @@ def populate_parameters(sec_potential, parameters):
2.023,
]
* ureg.radian,
'energies': [
1.5,
1.20083102,
0.93490305,
0.70221607,
0.50277008,
0.3365651,
0.20360111,
0.10387812,
0.03739612,
0.00415512,
0.00415512,
0.03739612,
0.10387812,
0.20360111,
0.3365651,
0.50277008,
0.70221607,
0.93490305,
1.20083102,
1.5,
]
* ureg.kJ
/ MOL,
# 'energies': [ # ! These are the exact energies. We will test the auto-generation of these energies from the forces.
# 1.5,
# 1.20083102,
# 0.93490305,
# 0.70221607,
# 0.50277008,
# 0.3365651,
# 0.20360111,
# 0.10387812,
# 0.03739612,
# 0.00415512,
# 0.00415512,
# 0.03739612,
# 0.10387812,
# 0.20360111,
# 0.3365651,
# 0.50277008,
# 0.70221607,
# 0.93490305,
# 1.20083102,
# 1.5,
# ]
# * ureg.kJ
# / MOL,
'forces': [
15.0,
13.42105263,
Expand All @@ -576,7 +620,10 @@ def populate_parameters(sec_potential, parameters):
-11.84210526,
-13.42105263,
-15.0,
],
]
* ureg.kJ
/ MOL
/ ureg.radian,
},
results_tabulated_angle,
)
Expand All @@ -591,8 +638,8 @@ def populate_parameters(sec_potential, parameters):
# data_fene_bond,
data_tabulated_bond,
# data_custom_bond,
# data_harmonic_angle,
# data_tabulated_angle,
data_harmonic_angle,
data_tabulated_angle,
],
)
def test_potentials(
Expand Down Expand Up @@ -626,10 +673,10 @@ def test_potentials(
sec_FF.contributions[-1].normalize(EntryArchive, logger) # BoundLogger)

potential_dict = sec_FF.contributions[-1].m_to_dict()
potential_dict = {
potential_dict = { # ! The dev is required to add new results to the dictionary, this will not be caught by the test!
key: value for key, value in potential_dict.items() if key in results
}
print(potential_dict)
print(results)
assert 1 == 2
# assert 1 == 2
assert_dict_equal(potential_dict, results)

0 comments on commit 0af786c

Please sign in to comment.