Hide keyboard shortcuts

Hot-keys on this page

r m x p   toggle line displays

j k   next/prev highlighted chunk

0   (zero) top of page

1   (one) first highlighted chunk

1import numpy as np 

2 

3# NB! This module was ported from a 4 year old CamposASE2 module. 

4 

5"""Bayesian Error Estimation 

6 

7For details, see: "Bayesian Error Estimation in Density Functional 

8Theory", J. J. Mortensen, K. Kaasbjerg, S. L. Frederiksen, 

9J. K. Norskov, J. P. Sethna, K. W. Jacobsen, Phys. Rev. Lett. 95, 

10216401 (2005).""" 

11 

12# T 

13# cost(c) = cost0 + 0.5 * (c - c0) H (c - c0) 

14# 

15 

16# Cost function minimum value: 

17cost0 = 3.4660625596 

18 

19# Best fit parameters: 

20c0 = np.array([1.000787451, 0.1926284063, 1.896191546]) 

21 

22# Hessian: 

23# H = np.array([[ 1.770035168e+03, -3.732470432e+02, -2.105836167e+02], 

24# [-3.732470432e+02, 1.188857209e+02, 6.054102443e+01], 

25# [-2.105836167e+02, 6.054102443e+01, 3.211200293e+01]]) 

26# 

27# 0.5 * np * T = cost0 (np=3: number of parameters) 

28T = cost0 * 2 / 3 

29 

30 

31def make_ensemble(N=1000, rng=np.random): 

32 M = np.array([(0.066, -0.812, 1.996), 

33 (0.055, 0.206, 0.082), 

34 (-0.034, 0.007, 0.004)]) 

35 alpha = rng.normal(0.0, 1.0, (N, 3)) 

36 return c0 + np.dot(alpha, M) 

37 

38 

39def get_ensemble_energies(atoms, ensemble): 

40 if hasattr(atoms, 'get_calculator'): 

41 coefs = atoms.calc.get_ensemble_coefficients() 

42 else: 

43 coefs = atoms 

44 return coefs[0] + np.dot(ensemble, coefs[1:])