Re: Elementary Functions

From: Richard Mathar (mathar_at_amer.strw.leidenuniv.nl)
Date: 10/08/04 

Date: Fri, 8 Oct 2004 22:16:14 +0000 (UTC)

In article <vabjp6w3vqd3@legacy>,
 derekroconnor@eircom.net ("Derek O'Connor") writes:
>Does anyone have a set of Fortran, C, Pascal, Matlab, or etc. programs
>that are simple implementations of the Elementary Functions (exp, sin
>cos, etc.)?
>
>I've searched around quite a lot and could find nothing. I have the
>LIBM implementations but these are too complicated to use in the class
>I'm teaching.
>..

After about 10 minutes of typing I'd propose the following (not tested!)
Something more advanced could use the Tchebyshev expansions tabulated
by Schonfelder in Math Comp 34 (149) (1980) 237-244

/*****************************************************************
 Naive attempt at a representation of the sine as it's Talor series
 sin(x)= x- x^3/3! +x^5/5! -x^7/7!..
***************************************************************/
double sinTayl(const double x)
{
        double resul=x ;
        double facto=resul ;
        for(int powe=3 ; powe< 100 ; powe += 2)
        {
                facto *= -x*x/((powe-1)*powe) ;
                resul += facto ;
                if( fabs(facto) < 1.e-13*fabs(resul) )
                        return resul ;
        }
        fprintf(stderr,"sinTayl does not converge for x=%e\n",x) ;
        return resul ;
}
/*****************************************************************
 Naive attempt at a representation of the sine as it's Talor series
 cos(x)= 1- x^2/2! +x^4/4! -x^6/6!..
***************************************************************/
double cosTayl(const double x)
{
        double resul=1. ;
        double facto=resul ;
        for(int powe=2 ; powe< 100 ; powe += 2)
        {
                facto *= -x*x/((powe-1)*powe) ;
                resul += facto ;
                if( fabs(facto) < 1.e-13*fabs(resul) )
                        return resul ;
        }
        fprintf(stderr,"cosTayl does not converge for x=%e\n",x) ;
        return resul ;
}
/*****************************************************************
 Naive attempt at a representation of the expon as it's Talor series
 exp(x)= 1+x+ x^2/2!+x^3/2!+x^4/4! +x^5/5!..
***************************************************************/
double expTayl(const double x)
{
        double resul=1. ;
        double facto=resul ;
        for(int powe=1 ; powe< 100 ; powe ++)
        {
                facto *= x/powe ;
                resul += facto ;
                if( fabs(facto) < 1.e-13*fabs(resul) )
                        return resul ;
        }
        fprintf(stderr,"expTayl does not converge for x=%e\n",x) ;
        return resul ;
}

Relevant Pages

  • Re: Elementary Functions
    ... >>LIBM implementations but these are too complicated to use in the ... > Naive attempt at a representation of the sine as it's Talor series ...
    (sci.math.num-analysis)
  • Re: Layout Hell-o
    ... The C standard requires that integral types use a "pure binary ... the representation match the "native" pure binary representation used ... C implementations essentially always do use the native ... integers if the CPU uses big-endian integers. ...
    (comp.lang.cobol)
  • Re: How to extract bytes from long?
    ... >> Implementations are allowed trap negative zero. ... zero" is meaningful in standard C (depending on the representation ... ...which is guaranteed to have an all-ones bit pattern, ...
    (comp.lang.c)
  • Re: what is trap representation?
    ... > of 'trap representation'. ... If the stored value of an object has such a representation ... Some implementations will issue some kind of exception, such as a SIGFPE, ... and therefore not part of standard C. ...
    (comp.lang.c.moderated)
  • Re: Google summer of code
    ... I wonder why implementations don't ... just wrap the C BSD sockets interface. ... The internal representation of an address wraps the C representation, ...
    (comp.lang.lisp)