# qr с выбором

 ``` 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108``` ```__author__ = 'Lenovo' import pprint from pprint import pprint from math import sqrt import numpy as np np.set_printoptions(suppress=True, linewidth=160, precision=10) import scipy from scipy import linalg import pprint A = scipy.array([[1,2,3], [1,5,6], [1,8,9], [1,11,12]]) Q, R, P = scipy.linalg.qr(A, pivoting=True) print "A:" pprint.pprint(A) print "Q:" pprint.pprint(Q) print "R:" pprint.pprint(R) print "P:" pprint.pprint(P) #################################################################################### def sign(x): return cmp(x, 0) def norm(x): return sum(i ** 2 for i in x) ** 0.5 def house(x): n = len(x) mu = norm(x) v = x if mu != 0: beta = x[0] + sign(x[0]) * mu v = v / beta v[0] = 1 return v def rowHouse(a, v): m = len(v) beta = -2.0 / sum(i ** 2 for i in v) vm = np.array([v]).transpose() w = beta * np.dot(a.transpose(), vm) p = np.dot(vm, w.transpose()) pa = a + p return pa def swap(a, r, k): temp = a[:, k].copy() a[:, k] = a[:, r] a[:, r] = temp def swapElems(c, r, k): temp = c[k] c[k] = c[r] c[r] = temp def qrColumnPivoting(permutation, a): m = len(a) n = len(a[0]) c = np.zeros(n) for i in range(0, n): c[i] = sum(j**2 for j in a[:,i]) r = 0 t = c.max() k = c.argmax() while t > 0: permutation[r] = k permutation[k] = r swap(a, r, k) swapElems(c, r, k) v = house(a[r:m, r]) a[r:m, r:n] = rowHouse(a[r:m, r:n], v) a[r+1:m, r] = v[1:m] for i in range(r+1, n): c[i] -= a[r, i]**2 if r < n-1: t = c[r+1:].max() k = c[r+1:].argmax() else: t = 0 r += 1 return permutation, a, r def fod(a): return qt, u a = np.array([[1,2,3], [1,5,6], [1,8,9], [1,11,12]], dtype=float) p, a, rank = qrColumnPivoting(np.arange(len(a[0])), a) qt, u = fod(a) print "A:" pprint.pprint(a) print "P:" pprint.pprint(p) ```