diff --git a/docs/pyerrors/correlators.html b/docs/pyerrors/correlators.html
index f4a2c4c3..02e101c0 100644
--- a/docs/pyerrors/correlators.html
+++ b/docs/pyerrors/correlators.html
@@ -1182,345 +1182,349 @@
969 content_string += "Description: " + self.tag + "\n"
970 if self.N != 1:
971 return content_string
- 972
- 973 if print_range[1]:
- 974 print_range[1] += 1
- 975 content_string += 'x0/a\tCorr(x0/a)\n------------------\n'
- 976 for i, sub_corr in enumerate(self.content[print_range[0]:print_range[1]]):
- 977 if sub_corr is None:
- 978 content_string += str(i + print_range[0]) + '\n'
- 979 else:
- 980 content_string += str(i + print_range[0])
- 981 for element in sub_corr:
- 982 content_string += '\t' + ' ' * int(element >= 0) + str(element)
- 983 content_string += '\n'
- 984 return content_string
- 985
- 986 def __str__(self):
- 987 return self.__repr__()
- 988
- 989 # We define the basic operations, that can be performed with correlators.
- 990 # While */+- get defined here, they only work for Corr*Obs and not Obs*Corr.
- 991 # This is because Obs*Corr checks Obs.__mul__ first and does not catch an exception.
- 992 # One could try and tell Obs to check if the y in __mul__ is a Corr and
- 993
- 994 def __add__(self, y):
- 995 if isinstance(y, Corr):
- 996 if ((self.N != y.N) or (self.T != y.T)):
- 997 raise Exception("Addition of Corrs with different shape")
- 998 newcontent = []
- 999 for t in range(self.T):
-1000 if _check_for_none(self, self.content[t]) or _check_for_none(y, y.content[t]):
-1001 newcontent.append(None)
-1002 else:
-1003 newcontent.append(self.content[t] + y.content[t])
-1004 return Corr(newcontent)
-1005
-1006 elif isinstance(y, (Obs, int, float, CObs)):
-1007 newcontent = []
-1008 for t in range(self.T):
-1009 if _check_for_none(self, self.content[t]):
-1010 newcontent.append(None)
-1011 else:
-1012 newcontent.append(self.content[t] + y)
-1013 return Corr(newcontent, prange=self.prange)
-1014 elif isinstance(y, np.ndarray):
-1015 if y.shape == (self.T,):
-1016 return Corr(list((np.array(self.content).T + y).T))
-1017 else:
-1018 raise ValueError("operands could not be broadcast together")
-1019 else:
-1020 raise TypeError("Corr + wrong type")
-1021
-1022 def __mul__(self, y):
-1023 if isinstance(y, Corr):
-1024 if not ((self.N == 1 or y.N == 1 or self.N == y.N) and self.T == y.T):
-1025 raise Exception("Multiplication of Corr object requires N=N or N=1 and T=T")
-1026 newcontent = []
-1027 for t in range(self.T):
-1028 if _check_for_none(self, self.content[t]) or _check_for_none(y, y.content[t]):
-1029 newcontent.append(None)
-1030 else:
-1031 newcontent.append(self.content[t] * y.content[t])
-1032 return Corr(newcontent)
-1033
-1034 elif isinstance(y, (Obs, int, float, CObs)):
-1035 newcontent = []
-1036 for t in range(self.T):
-1037 if _check_for_none(self, self.content[t]):
-1038 newcontent.append(None)
-1039 else:
-1040 newcontent.append(self.content[t] * y)
-1041 return Corr(newcontent, prange=self.prange)
-1042 elif isinstance(y, np.ndarray):
-1043 if y.shape == (self.T,):
-1044 return Corr(list((np.array(self.content).T * y).T))
-1045 else:
-1046 raise ValueError("operands could not be broadcast together")
-1047 else:
-1048 raise TypeError("Corr * wrong type")
-1049
-1050 def __truediv__(self, y):
-1051 if isinstance(y, Corr):
-1052 if not ((self.N == 1 or y.N == 1 or self.N == y.N) and self.T == y.T):
-1053 raise Exception("Multiplication of Corr object requires N=N or N=1 and T=T")
-1054 newcontent = []
-1055 for t in range(self.T):
-1056 if _check_for_none(self, self.content[t]) or _check_for_none(y, y.content[t]):
-1057 newcontent.append(None)
-1058 else:
-1059 newcontent.append(self.content[t] / y.content[t])
-1060 for t in range(self.T):
-1061 if _check_for_none(self, newcontent[t]):
-1062 continue
-1063 if np.isnan(np.sum(newcontent[t]).value):
-1064 newcontent[t] = None
-1065
-1066 if all([item is None for item in newcontent]):
-1067 raise Exception("Division returns completely undefined correlator")
-1068 return Corr(newcontent)
-1069
-1070 elif isinstance(y, (Obs, CObs)):
-1071 if isinstance(y, Obs):
-1072 if y.value == 0:
-1073 raise Exception('Division by zero will return undefined correlator')
-1074 if isinstance(y, CObs):
-1075 if y.is_zero():
-1076 raise Exception('Division by zero will return undefined correlator')
-1077
-1078 newcontent = []
-1079 for t in range(self.T):
-1080 if _check_for_none(self, self.content[t]):
-1081 newcontent.append(None)
-1082 else:
-1083 newcontent.append(self.content[t] / y)
-1084 return Corr(newcontent, prange=self.prange)
-1085
-1086 elif isinstance(y, (int, float)):
-1087 if y == 0:
-1088 raise Exception('Division by zero will return undefined correlator')
-1089 newcontent = []
-1090 for t in range(self.T):
-1091 if _check_for_none(self, self.content[t]):
-1092 newcontent.append(None)
-1093 else:
-1094 newcontent.append(self.content[t] / y)
-1095 return Corr(newcontent, prange=self.prange)
-1096 elif isinstance(y, np.ndarray):
-1097 if y.shape == (self.T,):
-1098 return Corr(list((np.array(self.content).T / y).T))
-1099 else:
-1100 raise ValueError("operands could not be broadcast together")
-1101 else:
-1102 raise TypeError('Corr / wrong type')
-1103
-1104 def __neg__(self):
-1105 newcontent = [None if _check_for_none(self, item) else -1. * item for item in self.content]
-1106 return Corr(newcontent, prange=self.prange)
-1107
-1108 def __sub__(self, y):
-1109 return self + (-y)
-1110
-1111 def __pow__(self, y):
-1112 if isinstance(y, (Obs, int, float, CObs)):
-1113 newcontent = [None if _check_for_none(self, item) else item**y for item in self.content]
-1114 return Corr(newcontent, prange=self.prange)
-1115 else:
-1116 raise TypeError('Type of exponent not supported')
-1117
-1118 def __abs__(self):
-1119 newcontent = [None if _check_for_none(self, item) else np.abs(item) for item in self.content]
-1120 return Corr(newcontent, prange=self.prange)
-1121
-1122 # The numpy functions:
-1123 def sqrt(self):
-1124 return self ** 0.5
-1125
-1126 def log(self):
-1127 newcontent = [None if _check_for_none(self, item) else np.log(item) for item in self.content]
-1128 return Corr(newcontent, prange=self.prange)
-1129
-1130 def exp(self):
-1131 newcontent = [None if _check_for_none(self, item) else np.exp(item) for item in self.content]
-1132 return Corr(newcontent, prange=self.prange)
-1133
-1134 def _apply_func_to_corr(self, func):
-1135 newcontent = [None if _check_for_none(self, item) else func(item) for item in self.content]
-1136 for t in range(self.T):
-1137 if _check_for_none(self, newcontent[t]):
-1138 continue
-1139 if np.isnan(np.sum(newcontent[t]).value):
-1140 newcontent[t] = None
-1141 if all([item is None for item in newcontent]):
-1142 raise Exception('Operation returns undefined correlator')
-1143 return Corr(newcontent)
-1144
-1145 def sin(self):
-1146 return self._apply_func_to_corr(np.sin)
-1147
-1148 def cos(self):
-1149 return self._apply_func_to_corr(np.cos)
-1150
-1151 def tan(self):
-1152 return self._apply_func_to_corr(np.tan)
-1153
-1154 def sinh(self):
-1155 return self._apply_func_to_corr(np.sinh)
-1156
-1157 def cosh(self):
-1158 return self._apply_func_to_corr(np.cosh)
-1159
-1160 def tanh(self):
-1161 return self._apply_func_to_corr(np.tanh)
-1162
-1163 def arcsin(self):
-1164 return self._apply_func_to_corr(np.arcsin)
-1165
-1166 def arccos(self):
-1167 return self._apply_func_to_corr(np.arccos)
-1168
-1169 def arctan(self):
-1170 return self._apply_func_to_corr(np.arctan)
-1171
-1172 def arcsinh(self):
-1173 return self._apply_func_to_corr(np.arcsinh)
-1174
-1175 def arccosh(self):
-1176 return self._apply_func_to_corr(np.arccosh)
-1177
-1178 def arctanh(self):
-1179 return self._apply_func_to_corr(np.arctanh)
-1180
-1181 # Right hand side operations (require tweak in main module to work)
-1182 def __radd__(self, y):
-1183 return self + y
+ 972 if isinstance(self[0], CObs):
+ 973 return content_string
+ 974
+ 975 if print_range[1]:
+ 976 print_range[1] += 1
+ 977 content_string += 'x0/a\tCorr(x0/a)\n------------------\n'
+ 978 for i, sub_corr in enumerate(self.content[print_range[0]:print_range[1]]):
+ 979 if sub_corr is None:
+ 980 content_string += str(i + print_range[0]) + '\n'
+ 981 else:
+ 982 content_string += str(i + print_range[0])
+ 983 for element in sub_corr:
+ 984 content_string += '\t' + ' ' * int(element >= 0) + str(element)
+ 985 content_string += '\n'
+ 986 return content_string
+ 987
+ 988 def __str__(self):
+ 989 return self.__repr__()
+ 990
+ 991 # We define the basic operations, that can be performed with correlators.
+ 992 # While */+- get defined here, they only work for Corr*Obs and not Obs*Corr.
+ 993 # This is because Obs*Corr checks Obs.__mul__ first and does not catch an exception.
+ 994 # One could try and tell Obs to check if the y in __mul__ is a Corr and
+ 995
+ 996 def __add__(self, y):
+ 997 if isinstance(y, Corr):
+ 998 if ((self.N != y.N) or (self.T != y.T)):
+ 999 raise Exception("Addition of Corrs with different shape")
+1000 newcontent = []
+1001 for t in range(self.T):
+1002 if _check_for_none(self, self.content[t]) or _check_for_none(y, y.content[t]):
+1003 newcontent.append(None)
+1004 else:
+1005 newcontent.append(self.content[t] + y.content[t])
+1006 return Corr(newcontent)
+1007
+1008 elif isinstance(y, (Obs, int, float, CObs)):
+1009 newcontent = []
+1010 for t in range(self.T):
+1011 if _check_for_none(self, self.content[t]):
+1012 newcontent.append(None)
+1013 else:
+1014 newcontent.append(self.content[t] + y)
+1015 return Corr(newcontent, prange=self.prange)
+1016 elif isinstance(y, np.ndarray):
+1017 if y.shape == (self.T,):
+1018 return Corr(list((np.array(self.content).T + y).T))
+1019 else:
+1020 raise ValueError("operands could not be broadcast together")
+1021 else:
+1022 raise TypeError("Corr + wrong type")
+1023
+1024 def __mul__(self, y):
+1025 if isinstance(y, Corr):
+1026 if not ((self.N == 1 or y.N == 1 or self.N == y.N) and self.T == y.T):
+1027 raise Exception("Multiplication of Corr object requires N=N or N=1 and T=T")
+1028 newcontent = []
+1029 for t in range(self.T):
+1030 if _check_for_none(self, self.content[t]) or _check_for_none(y, y.content[t]):
+1031 newcontent.append(None)
+1032 else:
+1033 newcontent.append(self.content[t] * y.content[t])
+1034 return Corr(newcontent)
+1035
+1036 elif isinstance(y, (Obs, int, float, CObs)):
+1037 newcontent = []
+1038 for t in range(self.T):
+1039 if _check_for_none(self, self.content[t]):
+1040 newcontent.append(None)
+1041 else:
+1042 newcontent.append(self.content[t] * y)
+1043 return Corr(newcontent, prange=self.prange)
+1044 elif isinstance(y, np.ndarray):
+1045 if y.shape == (self.T,):
+1046 return Corr(list((np.array(self.content).T * y).T))
+1047 else:
+1048 raise ValueError("operands could not be broadcast together")
+1049 else:
+1050 raise TypeError("Corr * wrong type")
+1051
+1052 def __truediv__(self, y):
+1053 if isinstance(y, Corr):
+1054 if not ((self.N == 1 or y.N == 1 or self.N == y.N) and self.T == y.T):
+1055 raise Exception("Multiplication of Corr object requires N=N or N=1 and T=T")
+1056 newcontent = []
+1057 for t in range(self.T):
+1058 if _check_for_none(self, self.content[t]) or _check_for_none(y, y.content[t]):
+1059 newcontent.append(None)
+1060 else:
+1061 newcontent.append(self.content[t] / y.content[t])
+1062 for t in range(self.T):
+1063 if _check_for_none(self, newcontent[t]):
+1064 continue
+1065 if np.isnan(np.sum(newcontent[t]).value):
+1066 newcontent[t] = None
+1067
+1068 if all([item is None for item in newcontent]):
+1069 raise Exception("Division returns completely undefined correlator")
+1070 return Corr(newcontent)
+1071
+1072 elif isinstance(y, (Obs, CObs)):
+1073 if isinstance(y, Obs):
+1074 if y.value == 0:
+1075 raise Exception('Division by zero will return undefined correlator')
+1076 if isinstance(y, CObs):
+1077 if y.is_zero():
+1078 raise Exception('Division by zero will return undefined correlator')
+1079
+1080 newcontent = []
+1081 for t in range(self.T):
+1082 if _check_for_none(self, self.content[t]):
+1083 newcontent.append(None)
+1084 else:
+1085 newcontent.append(self.content[t] / y)
+1086 return Corr(newcontent, prange=self.prange)
+1087
+1088 elif isinstance(y, (int, float)):
+1089 if y == 0:
+1090 raise Exception('Division by zero will return undefined correlator')
+1091 newcontent = []
+1092 for t in range(self.T):
+1093 if _check_for_none(self, self.content[t]):
+1094 newcontent.append(None)
+1095 else:
+1096 newcontent.append(self.content[t] / y)
+1097 return Corr(newcontent, prange=self.prange)
+1098 elif isinstance(y, np.ndarray):
+1099 if y.shape == (self.T,):
+1100 return Corr(list((np.array(self.content).T / y).T))
+1101 else:
+1102 raise ValueError("operands could not be broadcast together")
+1103 else:
+1104 raise TypeError('Corr / wrong type')
+1105
+1106 def __neg__(self):
+1107 newcontent = [None if _check_for_none(self, item) else -1. * item for item in self.content]
+1108 return Corr(newcontent, prange=self.prange)
+1109
+1110 def __sub__(self, y):
+1111 return self + (-y)
+1112
+1113 def __pow__(self, y):
+1114 if isinstance(y, (Obs, int, float, CObs)):
+1115 newcontent = [None if _check_for_none(self, item) else item**y for item in self.content]
+1116 return Corr(newcontent, prange=self.prange)
+1117 else:
+1118 raise TypeError('Type of exponent not supported')
+1119
+1120 def __abs__(self):
+1121 newcontent = [None if _check_for_none(self, item) else np.abs(item) for item in self.content]
+1122 return Corr(newcontent, prange=self.prange)
+1123
+1124 # The numpy functions:
+1125 def sqrt(self):
+1126 return self ** 0.5
+1127
+1128 def log(self):
+1129 newcontent = [None if _check_for_none(self, item) else np.log(item) for item in self.content]
+1130 return Corr(newcontent, prange=self.prange)
+1131
+1132 def exp(self):
+1133 newcontent = [None if _check_for_none(self, item) else np.exp(item) for item in self.content]
+1134 return Corr(newcontent, prange=self.prange)
+1135
+1136 def _apply_func_to_corr(self, func):
+1137 newcontent = [None if _check_for_none(self, item) else func(item) for item in self.content]
+1138 for t in range(self.T):
+1139 if _check_for_none(self, newcontent[t]):
+1140 continue
+1141 tmp_sum = np.sum(newcontent[t])
+1142 if hasattr(tmp_sum, "value"):
+1143 if np.isnan(tmp_sum.value):
+1144 newcontent[t] = None
+1145 if all([item is None for item in newcontent]):
+1146 raise Exception('Operation returns undefined correlator')
+1147 return Corr(newcontent)
+1148
+1149 def sin(self):
+1150 return self._apply_func_to_corr(np.sin)
+1151
+1152 def cos(self):
+1153 return self._apply_func_to_corr(np.cos)
+1154
+1155 def tan(self):
+1156 return self._apply_func_to_corr(np.tan)
+1157
+1158 def sinh(self):
+1159 return self._apply_func_to_corr(np.sinh)
+1160
+1161 def cosh(self):
+1162 return self._apply_func_to_corr(np.cosh)
+1163
+1164 def tanh(self):
+1165 return self._apply_func_to_corr(np.tanh)
+1166
+1167 def arcsin(self):
+1168 return self._apply_func_to_corr(np.arcsin)
+1169
+1170 def arccos(self):
+1171 return self._apply_func_to_corr(np.arccos)
+1172
+1173 def arctan(self):
+1174 return self._apply_func_to_corr(np.arctan)
+1175
+1176 def arcsinh(self):
+1177 return self._apply_func_to_corr(np.arcsinh)
+1178
+1179 def arccosh(self):
+1180 return self._apply_func_to_corr(np.arccosh)
+1181
+1182 def arctanh(self):
+1183 return self._apply_func_to_corr(np.arctanh)
1184
-1185 def __rsub__(self, y):
-1186 return -self + y
-1187
-1188 def __rmul__(self, y):
-1189 return self * y
-1190
-1191 def __rtruediv__(self, y):
-1192 return (self / y) ** (-1)
-1193
-1194 @property
-1195 def real(self):
-1196 def return_real(obs_OR_cobs):
-1197 if isinstance(obs_OR_cobs, CObs):
-1198 return obs_OR_cobs.real
-1199 else:
-1200 return obs_OR_cobs
-1201
-1202 return self._apply_func_to_corr(return_real)
-1203
-1204 @property
-1205 def imag(self):
-1206 def return_imag(obs_OR_cobs):
-1207 if isinstance(obs_OR_cobs, CObs):
-1208 return obs_OR_cobs.imag
-1209 else:
-1210 return obs_OR_cobs * 0 # So it stays the right type
-1211
-1212 return self._apply_func_to_corr(return_imag)
-1213
-1214 def prune(self, Ntrunc, tproj=3, t0proj=2, basematrix=None):
-1215 r''' Project large correlation matrix to lowest states
-1216
-1217 This method can be used to reduce the size of an (N x N) correlation matrix
-1218 to (Ntrunc x Ntrunc) by solving a GEVP at very early times where the noise
-1219 is still small.
+1185 # Right hand side operations (require tweak in main module to work)
+1186 def __radd__(self, y):
+1187 return self + y
+1188
+1189 def __rsub__(self, y):
+1190 return -self + y
+1191
+1192 def __rmul__(self, y):
+1193 return self * y
+1194
+1195 def __rtruediv__(self, y):
+1196 return (self / y) ** (-1)
+1197
+1198 @property
+1199 def real(self):
+1200 def return_real(obs_OR_cobs):
+1201 if isinstance(obs_OR_cobs.flatten()[0], CObs):
+1202 return np.vectorize(lambda x: x.real)(obs_OR_cobs)
+1203 else:
+1204 return obs_OR_cobs
+1205
+1206 return self._apply_func_to_corr(return_real)
+1207
+1208 @property
+1209 def imag(self):
+1210 def return_imag(obs_OR_cobs):
+1211 if isinstance(obs_OR_cobs.flatten()[0], CObs):
+1212 return np.vectorize(lambda x: x.imag)(obs_OR_cobs)
+1213 else:
+1214 return obs_OR_cobs * 0 # So it stays the right type
+1215
+1216 return self._apply_func_to_corr(return_imag)
+1217
+1218 def prune(self, Ntrunc, tproj=3, t0proj=2, basematrix=None):
+1219 r''' Project large correlation matrix to lowest states
1220
-1221 Parameters
-1222 ----------
-1223 Ntrunc: int
-1224 Rank of the target matrix.
-1225 tproj: int
-1226 Time where the eigenvectors are evaluated, corresponds to ts in the GEVP method.
-1227 The default value is 3.
-1228 t0proj: int
-1229 Time where the correlation matrix is inverted. Choosing t0proj=1 is strongly
-1230 discouraged for O(a) improved theories, since the correctness of the procedure
-1231 cannot be granted in this case. The default value is 2.
-1232 basematrix : Corr
-1233 Correlation matrix that is used to determine the eigenvectors of the
-1234 lowest states based on a GEVP. basematrix is taken to be the Corr itself if
-1235 is is not specified.
-1236
-1237 Notes
-1238 -----
-1239 We have the basematrix $C(t)$ and the target matrix $G(t)$. We start by solving
-1240 the GEVP $$C(t) v_n(t, t_0) = \lambda_n(t, t_0) C(t_0) v_n(t, t_0)$$ where $t \equiv t_\mathrm{proj}$
-1241 and $t_0 \equiv t_{0, \mathrm{proj}}$. The target matrix is projected onto the subspace of the
-1242 resulting eigenvectors $v_n, n=1,\dots,N_\mathrm{trunc}$ via
-1243 $$G^\prime_{i, j}(t) = (v_i, G(t) v_j)$$. This allows to reduce the size of a large
-1244 correlation matrix and to remove some noise that is added by irrelevant operators.
-1245 This may allow to use the GEVP on $G(t)$ at late times such that the theoretically motivated
-1246 bound $t_0 \leq t/2$ holds, since the condition number of $G(t)$ is decreased, compared to $C(t)$.
-1247 '''
-1248
-1249 if self.N == 1:
-1250 raise Exception('Method cannot be applied to one-dimensional correlators.')
-1251 if basematrix is None:
-1252 basematrix = self
-1253 if Ntrunc >= basematrix.N:
-1254 raise Exception('Cannot truncate using Ntrunc <= %d' % (basematrix.N))
-1255 if basematrix.N != self.N:
-1256 raise Exception('basematrix and targetmatrix have to be of the same size.')
-1257
-1258 evecs = basematrix.GEVP(t0proj, tproj, sort=None)[:Ntrunc]
-1259
-1260 tmpmat = np.empty((Ntrunc, Ntrunc), dtype=object)
-1261 rmat = []
-1262 for t in range(basematrix.T):
-1263 for i in range(Ntrunc):
-1264 for j in range(Ntrunc):
-1265 tmpmat[i][j] = evecs[i].T @ self[t] @ evecs[j]
-1266 rmat.append(np.copy(tmpmat))
-1267
-1268 newcontent = [None if (self.content[t] is None) else rmat[t] for t in range(self.T)]
-1269 return Corr(newcontent)
-1270
+1221 This method can be used to reduce the size of an (N x N) correlation matrix
+1222 to (Ntrunc x Ntrunc) by solving a GEVP at very early times where the noise
+1223 is still small.
+1224
+1225 Parameters
+1226 ----------
+1227 Ntrunc: int
+1228 Rank of the target matrix.
+1229 tproj: int
+1230 Time where the eigenvectors are evaluated, corresponds to ts in the GEVP method.
+1231 The default value is 3.
+1232 t0proj: int
+1233 Time where the correlation matrix is inverted. Choosing t0proj=1 is strongly
+1234 discouraged for O(a) improved theories, since the correctness of the procedure
+1235 cannot be granted in this case. The default value is 2.
+1236 basematrix : Corr
+1237 Correlation matrix that is used to determine the eigenvectors of the
+1238 lowest states based on a GEVP. basematrix is taken to be the Corr itself if
+1239 is is not specified.
+1240
+1241 Notes
+1242 -----
+1243 We have the basematrix $C(t)$ and the target matrix $G(t)$. We start by solving
+1244 the GEVP $$C(t) v_n(t, t_0) = \lambda_n(t, t_0) C(t_0) v_n(t, t_0)$$ where $t \equiv t_\mathrm{proj}$
+1245 and $t_0 \equiv t_{0, \mathrm{proj}}$. The target matrix is projected onto the subspace of the
+1246 resulting eigenvectors $v_n, n=1,\dots,N_\mathrm{trunc}$ via
+1247 $$G^\prime_{i, j}(t) = (v_i, G(t) v_j)$$. This allows to reduce the size of a large
+1248 correlation matrix and to remove some noise that is added by irrelevant operators.
+1249 This may allow to use the GEVP on $G(t)$ at late times such that the theoretically motivated
+1250 bound $t_0 \leq t/2$ holds, since the condition number of $G(t)$ is decreased, compared to $C(t)$.
+1251 '''
+1252
+1253 if self.N == 1:
+1254 raise Exception('Method cannot be applied to one-dimensional correlators.')
+1255 if basematrix is None:
+1256 basematrix = self
+1257 if Ntrunc >= basematrix.N:
+1258 raise Exception('Cannot truncate using Ntrunc <= %d' % (basematrix.N))
+1259 if basematrix.N != self.N:
+1260 raise Exception('basematrix and targetmatrix have to be of the same size.')
+1261
+1262 evecs = basematrix.GEVP(t0proj, tproj, sort=None)[:Ntrunc]
+1263
+1264 tmpmat = np.empty((Ntrunc, Ntrunc), dtype=object)
+1265 rmat = []
+1266 for t in range(basematrix.T):
+1267 for i in range(Ntrunc):
+1268 for j in range(Ntrunc):
+1269 tmpmat[i][j] = evecs[i].T @ self[t] @ evecs[j]
+1270 rmat.append(np.copy(tmpmat))
1271
-1272def _sort_vectors(vec_set, ts):
-1273 """Helper function used to find a set of Eigenvectors consistent over all timeslices"""
-1274 reference_sorting = np.array(vec_set[ts])
-1275 N = reference_sorting.shape[0]
-1276 sorted_vec_set = []
-1277 for t in range(len(vec_set)):
-1278 if vec_set[t] is None:
-1279 sorted_vec_set.append(None)
-1280 elif not t == ts:
-1281 perms = [list(o) for o in permutations([i for i in range(N)], N)]
-1282 best_score = 0
-1283 for perm in perms:
-1284 current_score = 1
-1285 for k in range(N):
-1286 new_sorting = reference_sorting.copy()
-1287 new_sorting[perm[k], :] = vec_set[t][k]
-1288 current_score *= abs(np.linalg.det(new_sorting))
-1289 if current_score > best_score:
-1290 best_score = current_score
-1291 best_perm = perm
-1292 sorted_vec_set.append([vec_set[t][k] for k in best_perm])
-1293 else:
-1294 sorted_vec_set.append(vec_set[t])
-1295
-1296 return sorted_vec_set
-1297
-1298
-1299def _check_for_none(corr, entry):
-1300 """Checks if entry for correlator corr is None"""
-1301 return len(list(filter(None, np.asarray(entry).flatten()))) < corr.N ** 2
+1272 newcontent = [None if (self.content[t] is None) else rmat[t] for t in range(self.T)]
+1273 return Corr(newcontent)
+1274
+1275
+1276def _sort_vectors(vec_set, ts):
+1277 """Helper function used to find a set of Eigenvectors consistent over all timeslices"""
+1278 reference_sorting = np.array(vec_set[ts])
+1279 N = reference_sorting.shape[0]
+1280 sorted_vec_set = []
+1281 for t in range(len(vec_set)):
+1282 if vec_set[t] is None:
+1283 sorted_vec_set.append(None)
+1284 elif not t == ts:
+1285 perms = [list(o) for o in permutations([i for i in range(N)], N)]
+1286 best_score = 0
+1287 for perm in perms:
+1288 current_score = 1
+1289 for k in range(N):
+1290 new_sorting = reference_sorting.copy()
+1291 new_sorting[perm[k], :] = vec_set[t][k]
+1292 current_score *= abs(np.linalg.det(new_sorting))
+1293 if current_score > best_score:
+1294 best_score = current_score
+1295 best_perm = perm
+1296 sorted_vec_set.append([vec_set[t][k] for k in best_perm])
+1297 else:
+1298 sorted_vec_set.append(vec_set[t])
+1299
+1300 return sorted_vec_set
+1301
1302
-1303
-1304def _GEVP_solver(Gt, G0):
-1305 """Helper function for solving the GEVP and sorting the eigenvectors.
+1303def _check_for_none(corr, entry):
+1304 """Checks if entry for correlator corr is None"""
+1305 return len(list(filter(None, np.asarray(entry).flatten()))) < corr.N ** 2
1306
-1307 The helper function assumes that both provided matrices are symmetric and
-1308 only processes the lower triangular part of both matrices. In case the matrices
-1309 are not symmetric the upper triangular parts are effectively discarded."""
-1310 return scipy.linalg.eigh(Gt, G0, lower=True)[1].T[::-1]
+1307
+1308def _GEVP_solver(Gt, G0):
+1309 """Helper function for solving the GEVP and sorting the eigenvectors.
+1310
+1311 The helper function assumes that both provided matrices are symmetric and
+1312 only processes the lower triangular part of both matrices. In case the matrices
+1313 are not symmetric the upper triangular parts are effectively discarded."""
+1314 return scipy.linalg.eigh(Gt, G0, lower=True)[1].T[::-1]
@@ -2495,304 +2499,308 @@
970 content_string += "Description: " + self.tag + "\n"
971 if self.N != 1:
972 return content_string
- 973
- 974 if print_range[1]:
- 975 print_range[1] += 1
- 976 content_string += 'x0/a\tCorr(x0/a)\n------------------\n'
- 977 for i, sub_corr in enumerate(self.content[print_range[0]:print_range[1]]):
- 978 if sub_corr is None:
- 979 content_string += str(i + print_range[0]) + '\n'
- 980 else:
- 981 content_string += str(i + print_range[0])
- 982 for element in sub_corr:
- 983 content_string += '\t' + ' ' * int(element >= 0) + str(element)
- 984 content_string += '\n'
- 985 return content_string
- 986
- 987 def __str__(self):
- 988 return self.__repr__()
- 989
- 990 # We define the basic operations, that can be performed with correlators.
- 991 # While */+- get defined here, they only work for Corr*Obs and not Obs*Corr.
- 992 # This is because Obs*Corr checks Obs.__mul__ first and does not catch an exception.
- 993 # One could try and tell Obs to check if the y in __mul__ is a Corr and
- 994
- 995 def __add__(self, y):
- 996 if isinstance(y, Corr):
- 997 if ((self.N != y.N) or (self.T != y.T)):
- 998 raise Exception("Addition of Corrs with different shape")
- 999 newcontent = []
-1000 for t in range(self.T):
-1001 if _check_for_none(self, self.content[t]) or _check_for_none(y, y.content[t]):
-1002 newcontent.append(None)
-1003 else:
-1004 newcontent.append(self.content[t] + y.content[t])
-1005 return Corr(newcontent)
-1006
-1007 elif isinstance(y, (Obs, int, float, CObs)):
-1008 newcontent = []
-1009 for t in range(self.T):
-1010 if _check_for_none(self, self.content[t]):
-1011 newcontent.append(None)
-1012 else:
-1013 newcontent.append(self.content[t] + y)
-1014 return Corr(newcontent, prange=self.prange)
-1015 elif isinstance(y, np.ndarray):
-1016 if y.shape == (self.T,):
-1017 return Corr(list((np.array(self.content).T + y).T))
-1018 else:
-1019 raise ValueError("operands could not be broadcast together")
-1020 else:
-1021 raise TypeError("Corr + wrong type")
-1022
-1023 def __mul__(self, y):
-1024 if isinstance(y, Corr):
-1025 if not ((self.N == 1 or y.N == 1 or self.N == y.N) and self.T == y.T):
-1026 raise Exception("Multiplication of Corr object requires N=N or N=1 and T=T")
-1027 newcontent = []
-1028 for t in range(self.T):
-1029 if _check_for_none(self, self.content[t]) or _check_for_none(y, y.content[t]):
-1030 newcontent.append(None)
-1031 else:
-1032 newcontent.append(self.content[t] * y.content[t])
-1033 return Corr(newcontent)
-1034
-1035 elif isinstance(y, (Obs, int, float, CObs)):
-1036 newcontent = []
-1037 for t in range(self.T):
-1038 if _check_for_none(self, self.content[t]):
-1039 newcontent.append(None)
-1040 else:
-1041 newcontent.append(self.content[t] * y)
-1042 return Corr(newcontent, prange=self.prange)
-1043 elif isinstance(y, np.ndarray):
-1044 if y.shape == (self.T,):
-1045 return Corr(list((np.array(self.content).T * y).T))
-1046 else:
-1047 raise ValueError("operands could not be broadcast together")
-1048 else:
-1049 raise TypeError("Corr * wrong type")
-1050
-1051 def __truediv__(self, y):
-1052 if isinstance(y, Corr):
-1053 if not ((self.N == 1 or y.N == 1 or self.N == y.N) and self.T == y.T):
-1054 raise Exception("Multiplication of Corr object requires N=N or N=1 and T=T")
-1055 newcontent = []
-1056 for t in range(self.T):
-1057 if _check_for_none(self, self.content[t]) or _check_for_none(y, y.content[t]):
-1058 newcontent.append(None)
-1059 else:
-1060 newcontent.append(self.content[t] / y.content[t])
-1061 for t in range(self.T):
-1062 if _check_for_none(self, newcontent[t]):
-1063 continue
-1064 if np.isnan(np.sum(newcontent[t]).value):
-1065 newcontent[t] = None
-1066
-1067 if all([item is None for item in newcontent]):
-1068 raise Exception("Division returns completely undefined correlator")
-1069 return Corr(newcontent)
-1070
-1071 elif isinstance(y, (Obs, CObs)):
-1072 if isinstance(y, Obs):
-1073 if y.value == 0:
-1074 raise Exception('Division by zero will return undefined correlator')
-1075 if isinstance(y, CObs):
-1076 if y.is_zero():
-1077 raise Exception('Division by zero will return undefined correlator')
-1078
-1079 newcontent = []
-1080 for t in range(self.T):
-1081 if _check_for_none(self, self.content[t]):
-1082 newcontent.append(None)
-1083 else:
-1084 newcontent.append(self.content[t] / y)
-1085 return Corr(newcontent, prange=self.prange)
-1086
-1087 elif isinstance(y, (int, float)):
-1088 if y == 0:
-1089 raise Exception('Division by zero will return undefined correlator')
-1090 newcontent = []
-1091 for t in range(self.T):
-1092 if _check_for_none(self, self.content[t]):
-1093 newcontent.append(None)
-1094 else:
-1095 newcontent.append(self.content[t] / y)
-1096 return Corr(newcontent, prange=self.prange)
-1097 elif isinstance(y, np.ndarray):
-1098 if y.shape == (self.T,):
-1099 return Corr(list((np.array(self.content).T / y).T))
-1100 else:
-1101 raise ValueError("operands could not be broadcast together")
-1102 else:
-1103 raise TypeError('Corr / wrong type')
-1104
-1105 def __neg__(self):
-1106 newcontent = [None if _check_for_none(self, item) else -1. * item for item in self.content]
-1107 return Corr(newcontent, prange=self.prange)
-1108
-1109 def __sub__(self, y):
-1110 return self + (-y)
-1111
-1112 def __pow__(self, y):
-1113 if isinstance(y, (Obs, int, float, CObs)):
-1114 newcontent = [None if _check_for_none(self, item) else item**y for item in self.content]
-1115 return Corr(newcontent, prange=self.prange)
-1116 else:
-1117 raise TypeError('Type of exponent not supported')
-1118
-1119 def __abs__(self):
-1120 newcontent = [None if _check_for_none(self, item) else np.abs(item) for item in self.content]
-1121 return Corr(newcontent, prange=self.prange)
-1122
-1123 # The numpy functions:
-1124 def sqrt(self):
-1125 return self ** 0.5
-1126
-1127 def log(self):
-1128 newcontent = [None if _check_for_none(self, item) else np.log(item) for item in self.content]
-1129 return Corr(newcontent, prange=self.prange)
-1130
-1131 def exp(self):
-1132 newcontent = [None if _check_for_none(self, item) else np.exp(item) for item in self.content]
-1133 return Corr(newcontent, prange=self.prange)
-1134
-1135 def _apply_func_to_corr(self, func):
-1136 newcontent = [None if _check_for_none(self, item) else func(item) for item in self.content]
-1137 for t in range(self.T):
-1138 if _check_for_none(self, newcontent[t]):
-1139 continue
-1140 if np.isnan(np.sum(newcontent[t]).value):
-1141 newcontent[t] = None
-1142 if all([item is None for item in newcontent]):
-1143 raise Exception('Operation returns undefined correlator')
-1144 return Corr(newcontent)
-1145
-1146 def sin(self):
-1147 return self._apply_func_to_corr(np.sin)
-1148
-1149 def cos(self):
-1150 return self._apply_func_to_corr(np.cos)
-1151
-1152 def tan(self):
-1153 return self._apply_func_to_corr(np.tan)
-1154
-1155 def sinh(self):
-1156 return self._apply_func_to_corr(np.sinh)
-1157
-1158 def cosh(self):
-1159 return self._apply_func_to_corr(np.cosh)
-1160
-1161 def tanh(self):
-1162 return self._apply_func_to_corr(np.tanh)
-1163
-1164 def arcsin(self):
-1165 return self._apply_func_to_corr(np.arcsin)
-1166
-1167 def arccos(self):
-1168 return self._apply_func_to_corr(np.arccos)
-1169
-1170 def arctan(self):
-1171 return self._apply_func_to_corr(np.arctan)
-1172
-1173 def arcsinh(self):
-1174 return self._apply_func_to_corr(np.arcsinh)
-1175
-1176 def arccosh(self):
-1177 return self._apply_func_to_corr(np.arccosh)
-1178
-1179 def arctanh(self):
-1180 return self._apply_func_to_corr(np.arctanh)
-1181
-1182 # Right hand side operations (require tweak in main module to work)
-1183 def __radd__(self, y):
-1184 return self + y
+ 973 if isinstance(self[0], CObs):
+ 974 return content_string
+ 975
+ 976 if print_range[1]:
+ 977 print_range[1] += 1
+ 978 content_string += 'x0/a\tCorr(x0/a)\n------------------\n'
+ 979 for i, sub_corr in enumerate(self.content[print_range[0]:print_range[1]]):
+ 980 if sub_corr is None:
+ 981 content_string += str(i + print_range[0]) + '\n'
+ 982 else:
+ 983 content_string += str(i + print_range[0])
+ 984 for element in sub_corr:
+ 985 content_string += '\t' + ' ' * int(element >= 0) + str(element)
+ 986 content_string += '\n'
+ 987 return content_string
+ 988
+ 989 def __str__(self):
+ 990 return self.__repr__()
+ 991
+ 992 # We define the basic operations, that can be performed with correlators.
+ 993 # While */+- get defined here, they only work for Corr*Obs and not Obs*Corr.
+ 994 # This is because Obs*Corr checks Obs.__mul__ first and does not catch an exception.
+ 995 # One could try and tell Obs to check if the y in __mul__ is a Corr and
+ 996
+ 997 def __add__(self, y):
+ 998 if isinstance(y, Corr):
+ 999 if ((self.N != y.N) or (self.T != y.T)):
+1000 raise Exception("Addition of Corrs with different shape")
+1001 newcontent = []
+1002 for t in range(self.T):
+1003 if _check_for_none(self, self.content[t]) or _check_for_none(y, y.content[t]):
+1004 newcontent.append(None)
+1005 else:
+1006 newcontent.append(self.content[t] + y.content[t])
+1007 return Corr(newcontent)
+1008
+1009 elif isinstance(y, (Obs, int, float, CObs)):
+1010 newcontent = []
+1011 for t in range(self.T):
+1012 if _check_for_none(self, self.content[t]):
+1013 newcontent.append(None)
+1014 else:
+1015 newcontent.append(self.content[t] + y)
+1016 return Corr(newcontent, prange=self.prange)
+1017 elif isinstance(y, np.ndarray):
+1018 if y.shape == (self.T,):
+1019 return Corr(list((np.array(self.content).T + y).T))
+1020 else:
+1021 raise ValueError("operands could not be broadcast together")
+1022 else:
+1023 raise TypeError("Corr + wrong type")
+1024
+1025 def __mul__(self, y):
+1026 if isinstance(y, Corr):
+1027 if not ((self.N == 1 or y.N == 1 or self.N == y.N) and self.T == y.T):
+1028 raise Exception("Multiplication of Corr object requires N=N or N=1 and T=T")
+1029 newcontent = []
+1030 for t in range(self.T):
+1031 if _check_for_none(self, self.content[t]) or _check_for_none(y, y.content[t]):
+1032 newcontent.append(None)
+1033 else:
+1034 newcontent.append(self.content[t] * y.content[t])
+1035 return Corr(newcontent)
+1036
+1037 elif isinstance(y, (Obs, int, float, CObs)):
+1038 newcontent = []
+1039 for t in range(self.T):
+1040 if _check_for_none(self, self.content[t]):
+1041 newcontent.append(None)
+1042 else:
+1043 newcontent.append(self.content[t] * y)
+1044 return Corr(newcontent, prange=self.prange)
+1045 elif isinstance(y, np.ndarray):
+1046 if y.shape == (self.T,):
+1047 return Corr(list((np.array(self.content).T * y).T))
+1048 else:
+1049 raise ValueError("operands could not be broadcast together")
+1050 else:
+1051 raise TypeError("Corr * wrong type")
+1052
+1053 def __truediv__(self, y):
+1054 if isinstance(y, Corr):
+1055 if not ((self.N == 1 or y.N == 1 or self.N == y.N) and self.T == y.T):
+1056 raise Exception("Multiplication of Corr object requires N=N or N=1 and T=T")
+1057 newcontent = []
+1058 for t in range(self.T):
+1059 if _check_for_none(self, self.content[t]) or _check_for_none(y, y.content[t]):
+1060 newcontent.append(None)
+1061 else:
+1062 newcontent.append(self.content[t] / y.content[t])
+1063 for t in range(self.T):
+1064 if _check_for_none(self, newcontent[t]):
+1065 continue
+1066 if np.isnan(np.sum(newcontent[t]).value):
+1067 newcontent[t] = None
+1068
+1069 if all([item is None for item in newcontent]):
+1070 raise Exception("Division returns completely undefined correlator")
+1071 return Corr(newcontent)
+1072
+1073 elif isinstance(y, (Obs, CObs)):
+1074 if isinstance(y, Obs):
+1075 if y.value == 0:
+1076 raise Exception('Division by zero will return undefined correlator')
+1077 if isinstance(y, CObs):
+1078 if y.is_zero():
+1079 raise Exception('Division by zero will return undefined correlator')
+1080
+1081 newcontent = []
+1082 for t in range(self.T):
+1083 if _check_for_none(self, self.content[t]):
+1084 newcontent.append(None)
+1085 else:
+1086 newcontent.append(self.content[t] / y)
+1087 return Corr(newcontent, prange=self.prange)
+1088
+1089 elif isinstance(y, (int, float)):
+1090 if y == 0:
+1091 raise Exception('Division by zero will return undefined correlator')
+1092 newcontent = []
+1093 for t in range(self.T):
+1094 if _check_for_none(self, self.content[t]):
+1095 newcontent.append(None)
+1096 else:
+1097 newcontent.append(self.content[t] / y)
+1098 return Corr(newcontent, prange=self.prange)
+1099 elif isinstance(y, np.ndarray):
+1100 if y.shape == (self.T,):
+1101 return Corr(list((np.array(self.content).T / y).T))
+1102 else:
+1103 raise ValueError("operands could not be broadcast together")
+1104 else:
+1105 raise TypeError('Corr / wrong type')
+1106
+1107 def __neg__(self):
+1108 newcontent = [None if _check_for_none(self, item) else -1. * item for item in self.content]
+1109 return Corr(newcontent, prange=self.prange)
+1110
+1111 def __sub__(self, y):
+1112 return self + (-y)
+1113
+1114 def __pow__(self, y):
+1115 if isinstance(y, (Obs, int, float, CObs)):
+1116 newcontent = [None if _check_for_none(self, item) else item**y for item in self.content]
+1117 return Corr(newcontent, prange=self.prange)
+1118 else:
+1119 raise TypeError('Type of exponent not supported')
+1120
+1121 def __abs__(self):
+1122 newcontent = [None if _check_for_none(self, item) else np.abs(item) for item in self.content]
+1123 return Corr(newcontent, prange=self.prange)
+1124
+1125 # The numpy functions:
+1126 def sqrt(self):
+1127 return self ** 0.5
+1128
+1129 def log(self):
+1130 newcontent = [None if _check_for_none(self, item) else np.log(item) for item in self.content]
+1131 return Corr(newcontent, prange=self.prange)
+1132
+1133 def exp(self):
+1134 newcontent = [None if _check_for_none(self, item) else np.exp(item) for item in self.content]
+1135 return Corr(newcontent, prange=self.prange)
+1136
+1137 def _apply_func_to_corr(self, func):
+1138 newcontent = [None if _check_for_none(self, item) else func(item) for item in self.content]
+1139 for t in range(self.T):
+1140 if _check_for_none(self, newcontent[t]):
+1141 continue
+1142 tmp_sum = np.sum(newcontent[t])
+1143 if hasattr(tmp_sum, "value"):
+1144 if np.isnan(tmp_sum.value):
+1145 newcontent[t] = None
+1146 if all([item is None for item in newcontent]):
+1147 raise Exception('Operation returns undefined correlator')
+1148 return Corr(newcontent)
+1149
+1150 def sin(self):
+1151 return self._apply_func_to_corr(np.sin)
+1152
+1153 def cos(self):
+1154 return self._apply_func_to_corr(np.cos)
+1155
+1156 def tan(self):
+1157 return self._apply_func_to_corr(np.tan)
+1158
+1159 def sinh(self):
+1160 return self._apply_func_to_corr(np.sinh)
+1161
+1162 def cosh(self):
+1163 return self._apply_func_to_corr(np.cosh)
+1164
+1165 def tanh(self):
+1166 return self._apply_func_to_corr(np.tanh)
+1167
+1168 def arcsin(self):
+1169 return self._apply_func_to_corr(np.arcsin)
+1170
+1171 def arccos(self):
+1172 return self._apply_func_to_corr(np.arccos)
+1173
+1174 def arctan(self):
+1175 return self._apply_func_to_corr(np.arctan)
+1176
+1177 def arcsinh(self):
+1178 return self._apply_func_to_corr(np.arcsinh)
+1179
+1180 def arccosh(self):
+1181 return self._apply_func_to_corr(np.arccosh)
+1182
+1183 def arctanh(self):
+1184 return self._apply_func_to_corr(np.arctanh)
1185
-1186 def __rsub__(self, y):
-1187 return -self + y
-1188
-1189 def __rmul__(self, y):
-1190 return self * y
-1191
-1192 def __rtruediv__(self, y):
-1193 return (self / y) ** (-1)
-1194
-1195 @property
-1196 def real(self):
-1197 def return_real(obs_OR_cobs):
-1198 if isinstance(obs_OR_cobs, CObs):
-1199 return obs_OR_cobs.real
-1200 else:
-1201 return obs_OR_cobs
-1202
-1203 return self._apply_func_to_corr(return_real)
-1204
-1205 @property
-1206 def imag(self):
-1207 def return_imag(obs_OR_cobs):
-1208 if isinstance(obs_OR_cobs, CObs):
-1209 return obs_OR_cobs.imag
-1210 else:
-1211 return obs_OR_cobs * 0 # So it stays the right type
-1212
-1213 return self._apply_func_to_corr(return_imag)
-1214
-1215 def prune(self, Ntrunc, tproj=3, t0proj=2, basematrix=None):
-1216 r''' Project large correlation matrix to lowest states
-1217
-1218 This method can be used to reduce the size of an (N x N) correlation matrix
-1219 to (Ntrunc x Ntrunc) by solving a GEVP at very early times where the noise
-1220 is still small.
+1186 # Right hand side operations (require tweak in main module to work)
+1187 def __radd__(self, y):
+1188 return self + y
+1189
+1190 def __rsub__(self, y):
+1191 return -self + y
+1192
+1193 def __rmul__(self, y):
+1194 return self * y
+1195
+1196 def __rtruediv__(self, y):
+1197 return (self / y) ** (-1)
+1198
+1199 @property
+1200 def real(self):
+1201 def return_real(obs_OR_cobs):
+1202 if isinstance(obs_OR_cobs.flatten()[0], CObs):
+1203 return np.vectorize(lambda x: x.real)(obs_OR_cobs)
+1204 else:
+1205 return obs_OR_cobs
+1206
+1207 return self._apply_func_to_corr(return_real)
+1208
+1209 @property
+1210 def imag(self):
+1211 def return_imag(obs_OR_cobs):
+1212 if isinstance(obs_OR_cobs.flatten()[0], CObs):
+1213 return np.vectorize(lambda x: x.imag)(obs_OR_cobs)
+1214 else:
+1215 return obs_OR_cobs * 0 # So it stays the right type
+1216
+1217 return self._apply_func_to_corr(return_imag)
+1218
+1219 def prune(self, Ntrunc, tproj=3, t0proj=2, basematrix=None):
+1220 r''' Project large correlation matrix to lowest states
1221
-1222 Parameters
-1223 ----------
-1224 Ntrunc: int
-1225 Rank of the target matrix.
-1226 tproj: int
-1227 Time where the eigenvectors are evaluated, corresponds to ts in the GEVP method.
-1228 The default value is 3.
-1229 t0proj: int
-1230 Time where the correlation matrix is inverted. Choosing t0proj=1 is strongly
-1231 discouraged for O(a) improved theories, since the correctness of the procedure
-1232 cannot be granted in this case. The default value is 2.
-1233 basematrix : Corr
-1234 Correlation matrix that is used to determine the eigenvectors of the
-1235 lowest states based on a GEVP. basematrix is taken to be the Corr itself if
-1236 is is not specified.
-1237
-1238 Notes
-1239 -----
-1240 We have the basematrix $C(t)$ and the target matrix $G(t)$. We start by solving
-1241 the GEVP $$C(t) v_n(t, t_0) = \lambda_n(t, t_0) C(t_0) v_n(t, t_0)$$ where $t \equiv t_\mathrm{proj}$
-1242 and $t_0 \equiv t_{0, \mathrm{proj}}$. The target matrix is projected onto the subspace of the
-1243 resulting eigenvectors $v_n, n=1,\dots,N_\mathrm{trunc}$ via
-1244 $$G^\prime_{i, j}(t) = (v_i, G(t) v_j)$$. This allows to reduce the size of a large
-1245 correlation matrix and to remove some noise that is added by irrelevant operators.
-1246 This may allow to use the GEVP on $G(t)$ at late times such that the theoretically motivated
-1247 bound $t_0 \leq t/2$ holds, since the condition number of $G(t)$ is decreased, compared to $C(t)$.
-1248 '''
-1249
-1250 if self.N == 1:
-1251 raise Exception('Method cannot be applied to one-dimensional correlators.')
-1252 if basematrix is None:
-1253 basematrix = self
-1254 if Ntrunc >= basematrix.N:
-1255 raise Exception('Cannot truncate using Ntrunc <= %d' % (basematrix.N))
-1256 if basematrix.N != self.N:
-1257 raise Exception('basematrix and targetmatrix have to be of the same size.')
-1258
-1259 evecs = basematrix.GEVP(t0proj, tproj, sort=None)[:Ntrunc]
-1260
-1261 tmpmat = np.empty((Ntrunc, Ntrunc), dtype=object)
-1262 rmat = []
-1263 for t in range(basematrix.T):
-1264 for i in range(Ntrunc):
-1265 for j in range(Ntrunc):
-1266 tmpmat[i][j] = evecs[i].T @ self[t] @ evecs[j]
-1267 rmat.append(np.copy(tmpmat))
-1268
-1269 newcontent = [None if (self.content[t] is None) else rmat[t] for t in range(self.T)]
-1270 return Corr(newcontent)
+1222 This method can be used to reduce the size of an (N x N) correlation matrix
+1223 to (Ntrunc x Ntrunc) by solving a GEVP at very early times where the noise
+1224 is still small.
+1225
+1226 Parameters
+1227 ----------
+1228 Ntrunc: int
+1229 Rank of the target matrix.
+1230 tproj: int
+1231 Time where the eigenvectors are evaluated, corresponds to ts in the GEVP method.
+1232 The default value is 3.
+1233 t0proj: int
+1234 Time where the correlation matrix is inverted. Choosing t0proj=1 is strongly
+1235 discouraged for O(a) improved theories, since the correctness of the procedure
+1236 cannot be granted in this case. The default value is 2.
+1237 basematrix : Corr
+1238 Correlation matrix that is used to determine the eigenvectors of the
+1239 lowest states based on a GEVP. basematrix is taken to be the Corr itself if
+1240 is is not specified.
+1241
+1242 Notes
+1243 -----
+1244 We have the basematrix $C(t)$ and the target matrix $G(t)$. We start by solving
+1245 the GEVP $$C(t) v_n(t, t_0) = \lambda_n(t, t_0) C(t_0) v_n(t, t_0)$$ where $t \equiv t_\mathrm{proj}$
+1246 and $t_0 \equiv t_{0, \mathrm{proj}}$. The target matrix is projected onto the subspace of the
+1247 resulting eigenvectors $v_n, n=1,\dots,N_\mathrm{trunc}$ via
+1248 $$G^\prime_{i, j}(t) = (v_i, G(t) v_j)$$. This allows to reduce the size of a large
+1249 correlation matrix and to remove some noise that is added by irrelevant operators.
+1250 This may allow to use the GEVP on $G(t)$ at late times such that the theoretically motivated
+1251 bound $t_0 \leq t/2$ holds, since the condition number of $G(t)$ is decreased, compared to $C(t)$.
+1252 '''
+1253
+1254 if self.N == 1:
+1255 raise Exception('Method cannot be applied to one-dimensional correlators.')
+1256 if basematrix is None:
+1257 basematrix = self
+1258 if Ntrunc >= basematrix.N:
+1259 raise Exception('Cannot truncate using Ntrunc <= %d' % (basematrix.N))
+1260 if basematrix.N != self.N:
+1261 raise Exception('basematrix and targetmatrix have to be of the same size.')
+1262
+1263 evecs = basematrix.GEVP(t0proj, tproj, sort=None)[:Ntrunc]
+1264
+1265 tmpmat = np.empty((Ntrunc, Ntrunc), dtype=object)
+1266 rmat = []
+1267 for t in range(basematrix.T):
+1268 for i in range(Ntrunc):
+1269 for j in range(Ntrunc):
+1270 tmpmat[i][j] = evecs[i].T @ self[t] @ evecs[j]
+1271 rmat.append(np.copy(tmpmat))
+1272
+1273 newcontent = [None if (self.content[t] is None) else rmat[t] for t in range(self.T)]
+1274 return Corr(newcontent)
@@ -4472,8 +4480,8 @@ specifies a custom path for the file (default '.')
- 1124 def sqrt(self):
-1125 return self ** 0.5
+ 1126 def sqrt(self):
+1127 return self ** 0.5
1127 def log(self):
-1128 newcontent = [None if _check_for_none(self, item) else np.log(item) for item in self.content]
-1129 return Corr(newcontent, prange=self.prange)
+ 1129 def log(self):
+1130 newcontent = [None if _check_for_none(self, item) else np.log(item) for item in self.content]
+1131 return Corr(newcontent, prange=self.prange)
1131 def exp(self):
-1132 newcontent = [None if _check_for_none(self, item) else np.exp(item) for item in self.content]
-1133 return Corr(newcontent, prange=self.prange)
+ 1133 def exp(self):
+1134 newcontent = [None if _check_for_none(self, item) else np.exp(item) for item in self.content]
+1135 return Corr(newcontent, prange=self.prange)
1146 def sin(self):
-1147 return self._apply_func_to_corr(np.sin)
+ 1150 def sin(self):
+1151 return self._apply_func_to_corr(np.sin)
1149 def cos(self):
-1150 return self._apply_func_to_corr(np.cos)
+ 1153 def cos(self):
+1154 return self._apply_func_to_corr(np.cos)
1152 def tan(self):
-1153 return self._apply_func_to_corr(np.tan)
+ 1156 def tan(self):
+1157 return self._apply_func_to_corr(np.tan)
1155 def sinh(self):
-1156 return self._apply_func_to_corr(np.sinh)
+ 1159 def sinh(self):
+1160 return self._apply_func_to_corr(np.sinh)
1158 def cosh(self):
-1159 return self._apply_func_to_corr(np.cosh)
+ 1162 def cosh(self):
+1163 return self._apply_func_to_corr(np.cosh)
1161 def tanh(self):
-1162 return self._apply_func_to_corr(np.tanh)
+ 1165 def tanh(self):
+1166 return self._apply_func_to_corr(np.tanh)
1164 def arcsin(self):
-1165 return self._apply_func_to_corr(np.arcsin)
+ 1168 def arcsin(self):
+1169 return self._apply_func_to_corr(np.arcsin)
1167 def arccos(self):
-1168 return self._apply_func_to_corr(np.arccos)
+ 1171 def arccos(self):
+1172 return self._apply_func_to_corr(np.arccos)
1170 def arctan(self):
-1171 return self._apply_func_to_corr(np.arctan)
+ 1174 def arctan(self):
+1175 return self._apply_func_to_corr(np.arctan)
1173 def arcsinh(self):
-1174 return self._apply_func_to_corr(np.arcsinh)
+ 1177 def arcsinh(self):
+1178 return self._apply_func_to_corr(np.arcsinh)
1176 def arccosh(self):
-1177 return self._apply_func_to_corr(np.arccosh)
+ 1180 def arccosh(self):
+1181 return self._apply_func_to_corr(np.arccosh)
1179 def arctanh(self):
-1180 return self._apply_func_to_corr(np.arctanh)
+ 1183 def arctanh(self):
+1184 return self._apply_func_to_corr(np.arctanh)
1215 def prune(self, Ntrunc, tproj=3, t0proj=2, basematrix=None):
-1216 r''' Project large correlation matrix to lowest states
-1217
-1218 This method can be used to reduce the size of an (N x N) correlation matrix
-1219 to (Ntrunc x Ntrunc) by solving a GEVP at very early times where the noise
-1220 is still small.
+ 1219 def prune(self, Ntrunc, tproj=3, t0proj=2, basematrix=None):
+1220 r''' Project large correlation matrix to lowest states
1221
-1222 Parameters
-1223 ----------
-1224 Ntrunc: int
-1225 Rank of the target matrix.
-1226 tproj: int
-1227 Time where the eigenvectors are evaluated, corresponds to ts in the GEVP method.
-1228 The default value is 3.
-1229 t0proj: int
-1230 Time where the correlation matrix is inverted. Choosing t0proj=1 is strongly
-1231 discouraged for O(a) improved theories, since the correctness of the procedure
-1232 cannot be granted in this case. The default value is 2.
-1233 basematrix : Corr
-1234 Correlation matrix that is used to determine the eigenvectors of the
-1235 lowest states based on a GEVP. basematrix is taken to be the Corr itself if
-1236 is is not specified.
-1237
-1238 Notes
-1239 -----
-1240 We have the basematrix $C(t)$ and the target matrix $G(t)$. We start by solving
-1241 the GEVP $$C(t) v_n(t, t_0) = \lambda_n(t, t_0) C(t_0) v_n(t, t_0)$$ where $t \equiv t_\mathrm{proj}$
-1242 and $t_0 \equiv t_{0, \mathrm{proj}}$. The target matrix is projected onto the subspace of the
-1243 resulting eigenvectors $v_n, n=1,\dots,N_\mathrm{trunc}$ via
-1244 $$G^\prime_{i, j}(t) = (v_i, G(t) v_j)$$. This allows to reduce the size of a large
-1245 correlation matrix and to remove some noise that is added by irrelevant operators.
-1246 This may allow to use the GEVP on $G(t)$ at late times such that the theoretically motivated
-1247 bound $t_0 \leq t/2$ holds, since the condition number of $G(t)$ is decreased, compared to $C(t)$.
-1248 '''
-1249
-1250 if self.N == 1:
-1251 raise Exception('Method cannot be applied to one-dimensional correlators.')
-1252 if basematrix is None:
-1253 basematrix = self
-1254 if Ntrunc >= basematrix.N:
-1255 raise Exception('Cannot truncate using Ntrunc <= %d' % (basematrix.N))
-1256 if basematrix.N != self.N:
-1257 raise Exception('basematrix and targetmatrix have to be of the same size.')
-1258
-1259 evecs = basematrix.GEVP(t0proj, tproj, sort=None)[:Ntrunc]
-1260
-1261 tmpmat = np.empty((Ntrunc, Ntrunc), dtype=object)
-1262 rmat = []
-1263 for t in range(basematrix.T):
-1264 for i in range(Ntrunc):
-1265 for j in range(Ntrunc):
-1266 tmpmat[i][j] = evecs[i].T @ self[t] @ evecs[j]
-1267 rmat.append(np.copy(tmpmat))
-1268
-1269 newcontent = [None if (self.content[t] is None) else rmat[t] for t in range(self.T)]
-1270 return Corr(newcontent)
+1222 This method can be used to reduce the size of an (N x N) correlation matrix
+1223 to (Ntrunc x Ntrunc) by solving a GEVP at very early times where the noise
+1224 is still small.
+1225
+1226 Parameters
+1227 ----------
+1228 Ntrunc: int
+1229 Rank of the target matrix.
+1230 tproj: int
+1231 Time where the eigenvectors are evaluated, corresponds to ts in the GEVP method.
+1232 The default value is 3.
+1233 t0proj: int
+1234 Time where the correlation matrix is inverted. Choosing t0proj=1 is strongly
+1235 discouraged for O(a) improved theories, since the correctness of the procedure
+1236 cannot be granted in this case. The default value is 2.
+1237 basematrix : Corr
+1238 Correlation matrix that is used to determine the eigenvectors of the
+1239 lowest states based on a GEVP. basematrix is taken to be the Corr itself if
+1240 is is not specified.
+1241
+1242 Notes
+1243 -----
+1244 We have the basematrix $C(t)$ and the target matrix $G(t)$. We start by solving
+1245 the GEVP $$C(t) v_n(t, t_0) = \lambda_n(t, t_0) C(t_0) v_n(t, t_0)$$ where $t \equiv t_\mathrm{proj}$
+1246 and $t_0 \equiv t_{0, \mathrm{proj}}$. The target matrix is projected onto the subspace of the
+1247 resulting eigenvectors $v_n, n=1,\dots,N_\mathrm{trunc}$ via
+1248 $$G^\prime_{i, j}(t) = (v_i, G(t) v_j)$$. This allows to reduce the size of a large
+1249 correlation matrix and to remove some noise that is added by irrelevant operators.
+1250 This may allow to use the GEVP on $G(t)$ at late times such that the theoretically motivated
+1251 bound $t_0 \leq t/2$ holds, since the condition number of $G(t)$ is decreased, compared to $C(t)$.
+1252 '''
+1253
+1254 if self.N == 1:
+1255 raise Exception('Method cannot be applied to one-dimensional correlators.')
+1256 if basematrix is None:
+1257 basematrix = self
+1258 if Ntrunc >= basematrix.N:
+1259 raise Exception('Cannot truncate using Ntrunc <= %d' % (basematrix.N))
+1260 if basematrix.N != self.N:
+1261 raise Exception('basematrix and targetmatrix have to be of the same size.')
+1262
+1263 evecs = basematrix.GEVP(t0proj, tproj, sort=None)[:Ntrunc]
+1264
+1265 tmpmat = np.empty((Ntrunc, Ntrunc), dtype=object)
+1266 rmat = []
+1267 for t in range(basematrix.T):
+1268 for i in range(Ntrunc):
+1269 for j in range(Ntrunc):
+1270 tmpmat[i][j] = evecs[i].T @ self[t] @ evecs[j]
+1271 rmat.append(np.copy(tmpmat))
+1272
+1273 newcontent = [None if (self.content[t] is None) else rmat[t] for t in range(self.T)]
+1274 return Corr(newcontent)