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 

2from scipy.linalg import null_space 

3import itertools 

4 

5 

6class Constraints: 

7 """ Class for handling linear constraints with right hand side equal to zero. 

8 

9 Parameters 

10 ---------- 

11 n_params 

12 number of parameters in model 

13 

14 Example 

15 ------- 

16 The following example demonstrates fitting of a cluster expansion under the 

17 constraint that parameter 2 and parameter 4 should be equal:: 

18 

19 >>> from icet.tools import Constraints 

20 >>> from icet.fitting import Optimizer 

21 >>> import numpy as np 

22 

23 >>> # Set up random sensing matrix and target "energies" 

24 >>> n_params = 10 

25 >>> A = np.random.random((10, n_params)) 

26 >>> y = np.random.random(10) 

27 

28 >>> # Define constraints 

29 >>> c = Constraints(n_params=n_params) 

30 >>> M = np.zeros((1, n_params)) 

31 >>> M[0, [2, 4]] = 1 

32 >>> c.add_constraint(M) 

33 

34 >>> # Do the actual fit and finally extract parameters 

35 >>> A_constrained = c.transform(A) 

36 >>> opt = Optimizer((A_constrained, y), fit_method='ridge') 

37 >>> opt.train() 

38 >>> parameters = c.inverse_transform(opt.parameters) 

39 

40 """ 

41 

42 def __init__(self, n_params: int): 

43 self.M = np.empty((0, n_params)) 

44 self.constraint_vectors = np.eye(n_params) 

45 

46 def transform(self, A: np.ndarray) -> np.ndarray: 

47 """ Transform array to constrained parameter space 

48 

49 Parameters 

50 ---------- 

51 A 

52 array to be transformed 

53 """ 

54 return A.dot(self.constraint_vectors) 

55 

56 def inverse_transform(self, A: np.ndarray) -> np.ndarray: 

57 """ Inverse transform array from constrained parameter space 

58 to unconstrained space 

59 

60 Parameters 

61 ---------- 

62 A 

63 array to be inversed transformed 

64 """ 

65 return self.constraint_vectors.dot(A) 

66 

67 def add_constraint(self, M: np.ndarray) -> None: 

68 """ Add a constraint matrix and resolve for the constraint space 

69 

70 Parameters 

71 ---------- 

72 M 

73 Constraint matrix with each constraint as a row. Can (but need not be) 

74 cluster vectors. 

75 """ 

76 M = np.array(M) 

77 self.M = np.vstack((self.M, M)) 

78 self.constraint_vectors = null_space(self.M) 

79 

80 

81def get_mixing_energy_constraints(cluster_space) -> Constraints: 

82 """ 

83 A cluster expansion of *mixing energy* should ideally predict zero energy 

84 for concentration 0 and 1. This function constructs a :class:`Constraints` 

85 object that enforces that condition during fitting. 

86 

87 Parameters 

88 ---------- 

89 cluster_space : ClusterSpace 

90 Cluster space corresponding to cluster expansion for which constraints 

91 should be imposed 

92 

93 Example 

94 ------- 

95 This example demonstrates how to constrain the mixing energy to zero 

96 at the pure phases in a toy example with random cluster vectors and 

97 random target energies:: 

98 

99 >>> from icet.tools import get_mixing_energy_constraints 

100 >>> from icet.fitting import Optimizer 

101 >>> from icet import ClusterSpace 

102 >>> from ase.build import bulk 

103 >>> import numpy as np 

104 

105 >>> # Set up cluster space along with random sensing matrix and target "energies" 

106 >>> prim = bulk('Au') 

107 >>> cs = ClusterSpace(prim, cutoffs=[6.0, 5.0], chemical_symbols=['Au', 'Ag']) 

108 >>> n_params = len(cs) 

109 >>> A = np.random.random((10, len(cs))) 

110 >>> y = np.random.random(10) 

111 

112 >>> # Define constraints 

113 >>> c = get_mixing_energy_constraints(cs) 

114 

115 >>> # Do the actual fit and finally extract parameters 

116 >>> A_constrained = c.transform(A) 

117 >>> opt = Optimizer((A_constrained, y), fit_method='ridge') 

118 >>> opt.train() 

119 >>> parameters = c.inverse_transform(opt.parameters) 

120 

121 Warning 

122 ------- 

123 Constraining the energy of one structure is always done at the expense of the 

124 fit quality of the others. Always expect that your :term:`CV` scores will increase 

125 somewhat when using this function. 

126 """ 

127 M = [] 

128 

129 prim = cluster_space.primitive_structure.copy() 

130 sublattices = cluster_space.get_sublattices(prim) 

131 chemical_symbols = [subl.chemical_symbols for subl in sublattices] 

132 

133 # Loop over all combinations of pure phases 

134 for symbols in itertools.product(*chemical_symbols): 

135 structure = prim.copy() 

136 for subl, symbol in zip(sublattices, symbols): 

137 for atom_index in subl.indices: 

138 structure[atom_index].symbol = symbol 

139 

140 # Add constraint for this pure phase 

141 M.append(cluster_space.get_cluster_vector(structure)) 

142 

143 c = Constraints(n_params=len(cluster_space)) 

144 c.add_constraint(M) 

145 return c