/*************************************************************** __ST_COS.C This file contains source code of functions for COS constants operations. PART OF : MPU - library . USAGE : Internal only . NOTE : NONE . Copyright (C) 2000 - 2024 by Andrew V.Kosteltsev. All Rights Reserved. ***************************************************************/ #ifdef HAVE_CONFIG_H #include #endif #include /* errno(3) */ #include /* strcpy(3) */ #include /* bzero(3) */ #include #include #include #include #include #include #include #include #include #include #include /*************************************************************** Кодировка имен файлов: Трехзначное десятичное число, представляющее количество 128-и битных слов, из которых состоят вещественные числа размещенные в массивах: размер чисел в битах кодировка -------------------- --------- 128 001 256 002 512 004 1024 008 2048 016 4096 032 8192 064 16384 128 32768 256 65536 512 (это предел); ПРИМЕРЫ: ------- ei_cos_001_emu32lsb.dfn - 128-бит, ei_cos_512_emu32lsb.dfn - 65536-бит. ***************************************************************/ #if MPU_MATH_FN_LIMIT >= 128 #if MPU_WORD_ORDER_BIG_ENDIAN == 0 #include #else #include #endif #endif /* MPU_MATH_FN_LIMIT >= 128 */ #if MPU_MATH_FN_LIMIT >= 256 #if MPU_WORD_ORDER_BIG_ENDIAN == 0 #include #else #include #endif #endif /* MPU_MATH_FN_LIMIT >= 256 */ #if MPU_MATH_FN_LIMIT >= 512 #if MPU_WORD_ORDER_BIG_ENDIAN == 0 #include #else #include #endif #endif /* MPU_MATH_FN_LIMIT >= 512 */ #if MPU_MATH_FN_LIMIT >= 1024 #if MPU_WORD_ORDER_BIG_ENDIAN == 0 #include #else #include #endif #endif /* MPU_MATH_FN_LIMIT >= 1024 */ #if MPU_MATH_FN_LIMIT >= 2048 #if MPU_WORD_ORDER_BIG_ENDIAN == 0 #include #else #include #endif #endif /* MPU_MATH_FN_LIMIT >= 2048 */ #if MPU_MATH_FN_LIMIT >= 4096 #if MPU_WORD_ORDER_BIG_ENDIAN == 0 #include #else #include #endif #endif /* MPU_MATH_FN_LIMIT >= 4096 */ #if MPU_MATH_FN_LIMIT >= 8192 #if MPU_WORD_ORDER_BIG_ENDIAN == 0 #include #else #include #endif #endif /* MPU_MATH_FN_LIMIT >= 8192 */ #if MPU_MATH_FN_LIMIT >= 16384 #if MPU_WORD_ORDER_BIG_ENDIAN == 0 #include #else #include #endif #endif /* MPU_MATH_FN_LIMIT >= 16384 */ #if MPU_MATH_FN_LIMIT >= 32768 #if MPU_WORD_ORDER_BIG_ENDIAN == 0 #include #else #include #endif #endif /* MPU_MATH_FN_LIMIT >= 32768 */ #if MPU_MATH_FN_LIMIT >= 65536 #if MPU_WORD_ORDER_BIG_ENDIAN == 0 #include #else #include #endif #endif /* MPU_MATH_FN_LIMIT >= 65536 */ static int _get_n_cos( int nb ) { int rc = 0; if( nb < NBR_32 || nb > MPU_MATH_FN_LIMIT ) { /* error: Invalid size of operand(s) */ __real_error_no = __R_ESIZE__; __STIND; /* Set REAL ind-produsing operation Flag */ return( rc ); } switch( nb ) { #if MPU_MATH_FN_LIMIT >= 128 case NBR_32 : case NBR_64 : case NBR_128 : rc = N_COS_C128; break; #endif /* MPU_MATH_FN_LIMIT >= 128 */ #if MPU_MATH_FN_LIMIT >= 256 case NBR_256 : rc = N_COS_C256; break; #endif /* MPU_MATH_FN_LIMIT >= 256 */ #if MPU_MATH_FN_LIMIT >= 512 case NBR_512 : rc = N_COS_C512; break; #endif /* MPU_MATH_FN_LIMIT >= 512 */ #if MPU_MATH_FN_LIMIT >= 1024 case NBR_1024 : rc = N_COS_C1024; break; #endif /* MPU_MATH_FN_LIMIT >= 1024 */ #if MPU_MATH_FN_LIMIT >= 2048 case NBR_2048 : rc = N_COS_C2048; break; #endif /* MPU_MATH_FN_LIMIT >= 2048 */ #if MPU_MATH_FN_LIMIT >= 4096 case NBR_4096 : rc = N_COS_C4096; break; #endif /* MPU_MATH_FN_LIMIT >= 4096 */ #if MPU_MATH_FN_LIMIT >= 8192 case NBR_8192 : rc = N_COS_C8192; break; #endif /* MPU_MATH_FN_LIMIT >= 8192 */ #if MPU_MATH_FN_LIMIT >= 16384 case NBR_16384: rc = N_COS_C16384; break; #endif /* MPU_MATH_FN_LIMIT >= 16384 */ default: { /* error: Invalid size of operand(s) */ __real_error_no = __R_ESIZE__; __STIND; /* Set REAL ind-produsing operation Flag */ break; } } /* End of switch( nb ) */ return( rc ); } /* End of _get_n_cos() */ static EMUSHORT *_get_cos__C_ptr( int nb ) { EMUSHORT *rc = (EMUSHORT *)NULL; if( nb < NBR_32 || nb > MPU_MATH_FN_LIMIT ) { /* error: Invalid size of operand(s) */ __real_error_no = __R_ESIZE__; __STIND; /* Set REAL ind-produsing operation Flag */ return( rc ); } switch( nb ) { #if MPU_MATH_FN_LIMIT >= 128 case NBR_32 : case NBR_64 : case NBR_128 : rc = (EMUSHORT *)&_ei_cos__C_128_[0][0]; break; #endif /* MPU_MATH_FN_LIMIT >= 128 */ #if MPU_MATH_FN_LIMIT >= 256 case NBR_256 : rc = (EMUSHORT *)&_ei_cos__C_256_[0][0]; break; #endif /* MPU_MATH_FN_LIMIT >= 256 */ #if MPU_MATH_FN_LIMIT >= 512 case NBR_512 : rc = (EMUSHORT *)&_ei_cos__C_512_[0][0]; break; #endif /* MPU_MATH_FN_LIMIT >= 512 */ #if MPU_MATH_FN_LIMIT >= 1024 case NBR_1024 : rc = (EMUSHORT *)&_ei_cos__C_1024_[0][0]; break; #endif /* MPU_MATH_FN_LIMIT >= 1024 */ #if MPU_MATH_FN_LIMIT >= 2048 case NBR_2048 : rc = (EMUSHORT *)&_ei_cos__C_2048_[0][0]; break; #endif /* MPU_MATH_FN_LIMIT >= 2048 */ #if MPU_MATH_FN_LIMIT >= 4096 case NBR_4096 : rc = (EMUSHORT *)&_ei_cos__C_4096_[0][0]; break; #endif /* MPU_MATH_FN_LIMIT >= 4096 */ #if MPU_MATH_FN_LIMIT >= 8192 case NBR_8192 : rc = (EMUSHORT *)&_ei_cos__C_8192_[0][0]; break; #endif /* MPU_MATH_FN_LIMIT >= 8192 */ #if MPU_MATH_FN_LIMIT >= 16384 case NBR_16384: rc = (EMUSHORT *)&_ei_cos__C_16384_[0][0]; break; #endif /* MPU_MATH_FN_LIMIT >= 16384 */ default: { /* error: Invalid size of operand(s) */ __real_error_no = __R_ESIZE__; __STIND; /* Set REAL ind-produsing operation Flag */ break; } } /* End of switch( nb ) */ return( rc ); } /* End of _get_cos__C_ptr() */ void ei_cos__C( EMUSHORT *eiy, EMUSHORT *eix, int nb ) /*************************************************************** Description : ei_cos__C() Работает с internal e-type data struct. EIX*EIX Concepts : RETURN in EIY = cos( EIX*k ) - 1 + --------- 2 on [-PI/4, PI/4], where k = pi/PI, PI is the rounted value of pi in machine precision. METHOD : Let z = x*x. Create a polinomial approximation to cos(k*x)-1+z/2 = z*z*(C0 + C1*z^1 + C2*z^2 + ... ). Then ei_cos__C(x*x) = z*z*(C0 + C1*z^1 + C2*z^2 + ... ). ACCURACY : In the absence of rounding error, the appriximation has absolute error less than EPSILON [see: FLOATP.H]. NOTE : The coefficient C's are obtained by a special Remez algorithm. =============================================== В INTERNET я нашел следующие алгоритмы: 409 cacm 355 356 14 5 May 1971 e2 A60 ------------------------------------------ H. Schmitt; Discrete {Chebychev} Curve Fit approximation;Chebyshev approximation; Chebyshev curve fitting; +Chebyshev polynomial; curve approximation;curve fitting; exchange algorithm; +polynomial approximation;Remez algorithm; 501 toms 95 97 2 1 March 1976 e2 F K2 --------------------------------------------- J. C. Simpson; {FORTRAN} Translation of Algorithm 409 Discrete {Chebyshev} Curve Fit approximation;polynomial approximation; exchange algorithm; +Chebyshev approximation; polynomial approximation; R,toms,95,4,1,March,1978,F. Futrell; последний из которых я перевел на "С", затем на язык операций повышенной разрядности. =============================================== Use Global Variable: Use Functions : internal_np( nb ); | real.c ei_copy( eiy, eix, nb ); | real.c _gen_zero( eic, nb ); | real.c ei_add( eic,eia,eib,nb ); | real.c ei_mul( eic,eia,eib,nb ); | real.c Parameters : EMUSHORT *eiy; - указатель на internal e-type data struct. TARGET; EMUSHORT *eix; - указатель на internal e-type data struct. SOURCE; int nb; - количество бит в external e-type data struct. Return : [void] ***************************************************************/ { EMUSHORT *x = NULL, *y = NULL; EMUSHORT *p; int np; int i, n = _get_n_cos( nb ); np = internal_np( nb ); /*** Allocate memory for x, y . *****************************/ x = (EMUSHORT *)__mpu_sbrk( (int)(np*SIZE_OF_EMUSHORT) ); if( !x ) { /* fatal error */ return; } y = (EMUSHORT *)__mpu_sbrk( (int)(np*SIZE_OF_EMUSHORT) ); if( !y ) { /* fatal error */ /* FREE x *****************/ __mpu_sbrk( -(int)(np*SIZE_OF_EMUSHORT) ); /**************************/ return; } /************************************************************/ _gen_zero( y, nb ); ei_copy( x, eix, nb ); p = _get_cos__C_ptr( nb ); for( i = 0; i < n; i++ ) { ei_add( y, y, p, nb ); ei_mul( y, y, x, nb ); p += np; } ei_mul( y, y, x, nb ); ei_copy( eiy, y, nb ); /* FREE x *****************/ /* FREE y *****************/ __mpu_sbrk( -(int)(2*np*SIZE_OF_EMUSHORT) ); /**************************/ } /* End of ei_cos__C() */ /*************************************************************** Hide internal symbols: ***************************************************************/ __mpu_hidden_decl(ei_cos__C); /* End of hide internal symbols. ***************************************************************/