Linear Algebra and the C Language/a0hs


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