(A, B, maxiter=200)
| 243 | |
| 244 | # author: A. Schoefl |
| 245 | def TOAR (A, B, maxiter=200): |
| 246 | |
| 247 | r = A.CreateVector() |
| 248 | |
| 249 | tmp_np = np.zeros((len(r), 2))#, dtype=complex) |
| 250 | r.SetRandom() |
| 251 | r.FV().imag = 0 |
| 252 | r /= Norm(r) |
| 253 | tmp_np[:,0] = r.FV().real |
| 254 | r.SetRandom() |
| 255 | r.FV().imag = 0 |
| 256 | r /= Norm(r) |
| 257 | tmp_np[:,1] = r.FV().real |
| 258 | |
| 259 | # tmp_np[:,1] = tmp_np[:,0] # just for testing |
| 260 | |
| 261 | |
| 262 | Q_np, X_np, perm = la.qr(tmp_np, pivoting=True, mode="economic") |
| 263 | # print(Q_np) |
| 264 | # print(Q_np.shape, X_np.shape, X_np) |
| 265 | |
| 266 | Q = MultiVector(r,0) |
| 267 | # r.FV().NumPy()[:] = Q_np[:,0] |
| 268 | |
| 269 | r.FV().real[:] = Q_np[:,0] |
| 270 | r.FV().imag = 0 |
| 271 | Q.Append(r) |
| 272 | # assign rank |
| 273 | if np.isclose(X_np[1,1], 0): |
| 274 | eta = 1 |
| 275 | else: |
| 276 | eta = 2 |
| 277 | |
| 278 | # r.FV().NumPy()[:] = Q_np[:,1] |
| 279 | r.FV().real[:] = Q_np[:,1] |
| 280 | r.FV().imag = 0 |
| 281 | Q.Append(r) |
| 282 | |
| 283 | # print (Q_np) |
| 284 | # print (Q[0], Q[1]) |
| 285 | print (InnerProduct(Q,Q)) |
| 286 | |
| 287 | gamma = np.linalg.norm(tmp_np, ord='fro') |
| 288 | |
| 289 | U1 = Matrix(eta,1, True) |
| 290 | U1.NumPy()[:,0] = X_np[:eta,1]/gamma |
| 291 | |
| 292 | U2 = Matrix(eta,1, True) |
| 293 | U2.NumPy()[:,0] = X_np[:eta,0]/gamma |
| 294 | print ("X_nb = ", X_np[:eta,:]) |
| 295 | print ("U1 = ", U1) |
| 296 | print ("U2 = ", U2) |
| 297 | |
| 298 | H = Matrix(maxiter, maxiter-1, True) |
| 299 | |
| 300 | # TODO: would be more efficient in C++ |
| 301 | for j in range(maxiter-1): |
| 302 |
nothing calls this directly
no test coverage detected