Linear Algebra and the C Language/a0b4
Install and compile this file in your working directory.

/* ------------------------------------ */
/* Save as : c00b.c */
/* ------------------------------------ */
#include "v_a.h"
/* ------------------------------------ */
/* ------------------------------------ */
#define RA R4
#define CA C3
#define Cb C1
/* ------------------------------------ */
/* ------------------------------------ */
int main(void)
{
double ta[RA*(CA+Cb)]={
// x1 x3 x4
+1, +0, +0, // A
+0, +1, +0, // B
+0, -1, +1, // C
-1, +0, -1 // D
};
double tb[RA*(CA+Cb)]={
+20+10, // A
+60 -10 -30, // B
+20, // C
-100+30, // D
};
double **A = ca_A_mR(ta,i_mR(RA,CA));
double **b = ca_A_mR(tb,i_mR(RA,Cb));
double **Q = i_mR(RA,CA);
double **Q_T = i_mR(CA,RA);
double **R = i_mR(CA,CA);
double **invR = i_mR(CA,CA);
double **invR_Q_T = i_mR(CA,RA);
double **x = i_mR(CA,C1);
clrscrn();
printf(" Copy/Paste into the octave windows \n\n");
p_Octave_mR(A,"a",P0);
printf(" [Q, R] = qr (a,0) \n\n");
QR_mR(A,Q,R);
printf(" Q :");
p_mR(Q, S10,P4, C10);
printf(" R :");
p_mR(R, S10,P4, C10);
stop();
clrscrn();
transpose_mR(Q,Q_T);
printf(" Q_T :");
pE_mR(Q_T,S9,P5, C3);
inv_mR(R,invR);
printf(" invR :");
pE_mR(invR,S9,P5, C6);
stop();
clrscrn();
printf(" Solving this system yields a unique\n"
" least squares solution, namely \n\n");
mul_mR(invR,Q_T,invR_Q_T);
mul_mR(invR_Q_T,b,x);
printf(" x = invR * Q_T * b :");
p_mR(x,S9,P5 ,C6);
stop();
f_mR(A);
f_mR(b);
f_mR(Q);
f_mR(Q_T);
f_mR(R);
f_mR(invR);
f_mR(x);
return 0;
}
/* ------------------------------------ */
/* ------------------------------------ */

Screen output example:
Copy/Paste into the octave windows
a=[
+1,+0,+0;
+0,+1,+0;
+0,-1,+1;
-1,+0,-1]
[Q, R] = qr (a,0)
Q :
+0.7071 +0.0000 -0.5000
+0.0000 +0.7071 +0.5000
+0.0000 -0.7071 +0.5000
-0.7071 +0.0000 -0.5000
R :
+1.4142 +0.0000 +0.7071
+0.0000 +1.4142 -0.7071
+0.0000 +0.0000 +1.0000
Press return to continue.
Q_T :
+7.07107e-01 +0.00000e+00 +0.00000e+00
+0.00000e+00 +7.07107e-01 -7.07107e-01
-5.00000e-01 +5.00000e-01 +5.00000e-01
-7.07107e-01
+0.00000e+00
-5.00000e-01
invR :
+7.07107e-01 -0.00000e+00 -5.00000e-01
-0.00000e+00 +7.07107e-01 +5.00000e-01
-0.00000e+00 -0.00000e+00 +1.00000e+00
Press return to continue.
Solving this system yields a unique
least squares solution, namely
x = invR * Q_T * b :
+30.00000
+20.00000
+40.00000
Press return to continue.