Linear Algebra and the C Language/a0lt

__NOTOC__ Install and compile this file in your working directory. /* ------------------------------------ */ /* Save as : c05d.c */ /* ------------------------------------ */ #include "v_a.h" /* ------------------------------------ */ /* ------------------------------------ */ #define RB R5 #define CB C3 /* ------------------------------------ */ /* ------------------------------------ */ int main(void) { double b[RB*(CB)]={ +0.000000000000, -2.000000000000, -0.500000000000, +0.000000000000, -1.666666666667, +0.000000000000, +1.000000000000, +0.000000000000, +0.000000000000, +0.000000000000, +1.000000000000, +0.000000000000, +0.000000000000, +0.000000000000, +1.000000000000 }; double x[RB*(C1)]={ -1, 2, -3, -3, -1, }; double **B = ca_A_mR(b, i_mR(RB,CB)); double **BT = transpose_mR(B, i_mR(CB,RB)); double **BTB = mul_mR(BT,B, i_mR(CB,CB)); double **invBTB = inv_mR(BTB, i_mR(CB,CB)); double **invBTB_BT = mul_mR(invBTB,BT, i_mR(CB,RB)); double **B_invBTB_BT = mul_mR(B,invBTB_BT, i_mR(RB,RB)); double **Id = eye_mR( i_mR(RB,RB)); double **V = sub_mR(Id,B_invBTB_BT, i_mR(RB,RB)); /* Id-(B inv(BT B) BT) */ double **X = ca_A_mR(x, i_mR(RB,C1)); double **VX = mul_mR(V,X, i_mR(RB,C1)); clrscrn(); printf(" B is a basis for the orthogonal complement of AT: \n\n" " Find a transformation matrix for \n" " a projection onto R%d: \n\n" " Proj(x) = [Id-(B inv(BTB)BT)] \n\n" " B: See the previous c file",RB); p_mR(B,S5,P4,C7); stop(); clrscrn(); printf(" BT:"); p_mR(BT,S5,P4,C7); printf(" BT B:"); p_mR(BTB,S5,P4,C7); printf(" inv(BTB):"); p_mR(invBTB,S5,P4,C7); stop(); clrscrn(); printf(" inv(BTB) BT:"); p_mR(invBTB_BT,S5,P4,C7); printf(" B inv(BTB)BT:"); p_mR(B_invBTB_BT,S5,P4,C7); printf(" V = Id - (B inv(BT B) BT)"); p_mR(V,S5,P4,C7); stop(); clrscrn(); printf(" V is transformation matrix for \n" " a projection onto a R%d subspace \n\n" " V:",RB); p_mR(V,S5,P4,C7); printf(" Proj(x) = [Id-(B inv(BTB)BT)] x \n\n" " Proj(x) = V x Verify the result with the first C file"); p_mR(VX,S5,P4,C7); stop(); f_mR(B); f_mR(BT); f_mR(BTB); f_mR(invBTB); f_mR(invBTB_BT); f_mR(V); /* Id-(B inv(BT B) BT) */ f_mR(X); f_mR(VX); return 0; } /* ------------------------------------ */ /* ------------------------------------ */ Screen output example: B is a basis for the orthogonal complement of AT: Find a transformation matrix for a projection onto R5: Proj(x) = [Id-(B inv(BTB)BT)] B: See the previous c file +0.0000 -2.0000 -0.5000 +0.0000 -1.6667 +0.0000 +1.0000 +0.0000 +0.0000 +0.0000 +1.0000 +0.0000 +0.0000 +0.0000 +1.0000 Press return to continue. BT: +0.0000 +0.0000 +1.0000 +0.0000 +0.0000 -2.0000 -1.6667 +0.0000 +1.0000 +0.0000 -0.5000 +0.0000 +0.0000 +0.0000 +1.0000 BT B: +1.0000 +0.0000 +0.0000 +0.0000 +7.7778 +1.0000 +0.0000 +1.0000 +1.2500 inv(BTB): +1.0000 -0.0000 -0.0000 -0.0000 +0.1433 -0.1146 -0.0000 -0.1146 +0.8917 Press return to continue. inv(BTB) BT: +0.0000 +0.0000 +1.0000 +0.0000 +0.0000 -0.2293 -0.2389 +0.0000 +0.1433 -0.1146 -0.2166 +0.1911 +0.0000 -0.1146 +0.8917 B inv(BTB)BT: +0.5669 +0.3822 +0.0000 -0.2293 -0.2166 +0.3822 +0.3981 +0.0000 -0.2389 +0.1911 +0.0000 +0.0000 +1.0000 +0.0000 +0.0000 -0.2293 -0.2389 +0.0000 +0.1433 -0.1146 -0.2166 +0.1911 +0.0000 -0.1146 +0.8917 V = Id - (B inv(BT B) BT) +0.4331 -0.3822 +0.0000 +0.2293 +0.2166 -0.3822 +0.6019 +0.0000 +0.2389 -0.1911 +0.0000 +0.0000 +0.0000 +0.0000 +0.0000 +0.2293 +0.2389 +0.0000 +0.8567 +0.1146 +0.2166 -0.1911 +0.0000 +0.1146 +0.1083 Press return to continue. V is transformation matrix for a projection onto a R5 subspace V: +0.4331 -0.3822 +0.0000 +0.2293 +0.2166 -0.3822 +0.6019 +0.0000 +0.2389 -0.1911 +0.0000 +0.0000 +0.0000 +0.0000 +0.0000 +0.2293 +0.2389 +0.0000 +0.8567 +0.1146 +0.2166 -0.1911 +0.0000 +0.1146 +0.1083 Proj(x) = [Id-(B inv(BTB)BT)] x Proj(x) = V x Verify the result with the first C file -2.1019 +1.0605 +0.0000 -2.4363 -1.0510 Press return to continue.