Linear Algebra and the C Language/a0hs

__NOTOC__ The Maclaurin series of the exponential function is : Exp(A) = Ide + A + (1/2) A2 + (1/6) A3 + ... + 1/(n-1)! A(n-1) Install and compile this file in your working directory. /* ------------------------------------ */ /* Save as : c00a.c */ /* ------------------------------------ */ #include "v_a.h" /* ------------------------------------ */ #define RCA RC4 /* ------------------------------------ */ /* ------------------------------------ */ double **X_Exp_mR( double **A, double **ExpA ) { double **A_n = i_RC_mR(A[R_SIZE][C0],A[C_SIZE][C0]); double **sA = i_RC_mR(A[R_SIZE][C0],A[C_SIZE][C0]); double **NewExpA = i_RC_mR(A[R_SIZE][C0],A[C_SIZE][C0]); int n = 25; int i = 0; for(i = 0; i <= n; i++) { pow_mR(i, A, A_n); smul_mR(1./factorial(i), A_n, sA); add_mR(ExpA,sA,NewExpA); c_mR(NewExpA,ExpA); } f_mR(A_n); f_mR(sA); f_mR(NewExpA); return (ExpA); } /* ------------------------------------ */ /* ------------------------------------ */ void fun(void) { double **A = rE_mR(i_mR(RCA,RCA),999,+1.E-3); double **ExpA = i_mR(RCA,RCA); clrscrn(); printf(" Copy/Paste into the octave window.\n"); p_Octave_mR(A, "A", P4); printf(" expm (A)\n\n"); Exp_mR(A,ExpA); printf(" ExpA:"); p_mR(ExpA, S8,P4,C6); f_mR(A); f_mR(ExpA); } /* ------------------------------------ */ int main(void) { time_t t; srand(time(&t)); do { fun(); } while(stop_w()); return 0; } /* ------------------------------------ */ /* ------------------------------------ */ Screen output example: Copy/Paste into the octave window. A=[ +0.5580,+0.7720,-0.5640,-0.3210; -0.5240,+0.7820,-0.5920,+0.9560; -0.6750,-0.1450,-0.8890,+0.2210; +0.0060,+0.4320,+0.0260,+0.3950] expm (A) ExpA: +1.6507 +1.5186 -0.7929 +0.0607 -0.8313 +2.2963 -0.4757 +1.8280 -0.5796 -0.4233 +0.6154 +0.1441 -0.1706 +0.7971 -0.1038 +1.8480 Press return to continue Press X return to stop