<![CDATA[
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;
array double complex q[m][m]; // in the case of m<n
array double complex r[m][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;
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][m]; // in the case of m>=n
array double complex r2[m][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);
}
]]>
<![CDATA[
]]>