int main() { array double complex a[2][4] = {1, 5, -7, 4, 3, 3, 2, -13}; /* a[m][n] */ array double complex a1[2][4] = {complex(1,2), 5, -7, 4, 3, 3, 2, -13}; /* a[m][n] */ int m = 2, n = 4; int lqr = min(m,n); array double complex q[m][lqr]; // in the case of m<n array double complex r[lqr][n] ; qrdecomp(a, q, r); printf("q^H*q = \n%f\n", conj(transpose(q))*q); printf("q = \n%f\n", q); printf("r = \n%f\n", r); printf("q*r = \n%f\n", q*r); qrdecomp(a1, q, r); printf("q^H*q = \n%f\n", conj(transpose(q))*q); printf("q = \n%f\n", q); printf("r = \n%f\n", r); printf("q*r = \n%f\n", q*r); m = 4; n = 2; lqr = min(m,n); array double complex a2[4][2] = {complex(1,2), complex(3,5), complex(-7,2), complex(3,4), complex(3,0), complex(2,3), complex(2,0), complex(2,-13)}; // a[m][n] array double complex q2[m][lqr]; // in the case of m>=n array double complex r2[lqr][n] ; qrdecomp(a2, q2, r2); printf("q2^H*q2 = \n%f\n", conj(transpose(q2))*q2); printf("q2 = \n%f\n", q2); printf("r2 = \n%f\n", r2); printf("q2*r2 = \n%f\n", q2*r2); }
q^H*q = complex(1.000000,0.000000) complex(0.000000,0.000000) complex(0.000000,0.000000) complex(1.000000,0.000000) q = complex(-0.316228,0.000000) complex(-0.948683,0.000000) complex(-0.948683,-0.000000) complex(0.316228,0.000000) r = complex(-3.162278,0.000000) complex(-4.427189,0.000000) complex(0.316228,0.000000) complex(11.067972,0.000000) complex(0.000000,0.000000) complex(-3.794733,0.000000) complex(7.273239,0.000000) complex(-7.905694,0.000000) q*r = complex(1.000000,0.000000) complex(5.000000,0.000000) complex(-7.000000,0.000000) complex(4.000000,0.000000) complex(3.000000,0.000000) complex(3.000000,0.000000) complex(2.000000,0.000000) complex(-13.000000,0.000000) q^H*q = complex(1.000000,0.000000) complex(-0.000000,0.000000) complex(-0.000000,-0.000000) complex(1.000000,0.000000) q = complex(-0.267261,-0.534522) complex(0.717137,-0.358569) complex(-0.801784,0.000000) complex(-0.000000,0.597614) r = complex(-3.741657,0.000000) complex(-3.741657,2.672612) complex(0.267261,-3.741657) complex(9.354143,2.138090) complex(0.000000,0.000000) complex(3.585686,0.000000) complex(-5.019960,-3.705209) complex(2.868549,9.203260) q*r = complex(1.000000,2.000000) complex(5.000000,0.000000) complex(-7.000000,0.000000) complex(4.000000,-0.000000) complex(3.000000,-0.000000) complex(3.000000,0.000000) complex(2.000000,0.000000) complex(-13.000000,-0.000000) q2^H*q2 = complex(1.000000,0.000000) complex(0.000000,0.000000) complex(0.000000,-0.000000) complex(1.000000,0.000000) q2 = complex(-0.118678,-0.237356) complex(-0.097267,-0.380224) complex(0.830747,-0.237356) complex(-0.175866,0.098249) complex(-0.356034,0.000000) complex(-0.110039,-0.362539) complex(-0.237356,0.000000) complex(-0.119864,0.804660) r2 = complex(-8.426150,0.000000) complex(-1.186782,6.171265) complex(0.000000,0.000000) complex(-14.335517,0.000000) q2*r2 = complex(1.000000,2.000000) complex(3.000000,5.000000) complex(-7.000000,2.000000) complex(3.000000,4.000000) complex(3.000000,-0.000000) complex(2.000000,3.000000) complex(2.000000,0.000000) complex(2.000000,-13.000000)