# cosine function returning long double (high resolution)

Discussion in 'C Programming' started by ratcharit@gmail.com, Mar 5, 2008.

1. ### Guest

Currently using cosine function in math.h

Currently I get:

1 = cos(1e^-7)

Is there another way for cos to return value of high accuracy say:

0.999999 = cos(1e^-7)

, Mar 5, 2008

2. ### jacob naviaGuest

wrote:
> Currently using cosine function in math.h
>
> Currently I get:
>
> 1 = cos(1e^-7)
>
> Is there another way for cos to return value of high accuracy say:
>
> 0.999999 = cos(1e^-7)
>

If you use the lcc-win compiler you get

cos(1e-7)=
0.999999999999995000000000000004166666666666665277777777777778025793650793650766093474426807762228769173214
using the extra precision floating point (qfloat).

In standard C the biggest floating point type is long double.

No way to do that in standard C.

--
jacob navia
jacob at jacob point remcomp point fr
logiciels/informatique
http://www.cs.virginia.edu/~lcc-win32

jacob navia, Mar 5, 2008

3. ### jacob naviaGuest

jacob navia wrote:
> wrote:
>> Currently using cosine function in math.h
>>
>> Currently I get:
>>
>> 1 = cos(1e^-7)
>>
>> Is there another way for cos to return value of high accuracy say:
>>
>> 0.999999 = cos(1e^-7)
>>

>
> If you use the lcc-win compiler you get
>
> cos(1e-7)=
> 0.999999999999995000000000000004166666666666665277777777777778025793650793650766093474426807762228769173214
>
> using the extra precision floating point (qfloat).
>
>
> In standard C the biggest floating point type is long double.
>
> No way to do that in standard C.
>

P.S.

WAIT

This is .00000000000000499999999999999, a number that is bigger than
DBL_EPSILON.

MAYBE your problem is a printing problem?

Which format did you use to print it?

--
jacob navia
jacob at jacob point remcomp point fr
logiciels/informatique
http://www.cs.virginia.edu/~lcc-win32

jacob navia, Mar 5, 2008
4. ### jacob naviaGuest

wrote:
> Currently using cosine function in math.h
>
> Currently I get:
>
> 1 = cos(1e^-7)
>
> Is there another way for cos to return value of high accuracy say:
>
> 0.999999 = cos(1e^-7)
>

This program:
#include <math.h>
#include <stdio.h>
int main(void)
{
printf("%40.17Lf\n",cosl(1e-7));
}

prints
0.99999999999999500

--
jacob navia
jacob at jacob point remcomp point fr
logiciels/informatique
http://www.cs.virginia.edu/~lcc-win32

jacob navia, Mar 5, 2008
5. ### Richard TobinGuest

In article <>,
<> wrote:

>Currently I get:
>
>1 = cos(1e^-7)
>
>Is there another way for cos to return value of high accuracy say:
>
>0.999999 = cos(1e^-7)

cos(1e-7) is roughly 1 - 1e-14 / 2, which is closer to 1 than to
0.999999. It should be distinguishable from 1 as an IEEE double,
but not as a float.

-- Richard

--
:wq

Richard Tobin, Mar 5, 2008
6. ### jacob naviaGuest

Richard Tobin wrote:
>
> cos(1e-7) is roughly 1 - 1e-14 / 2, which is closer to 1 than to
> 0.999999. It should be distinguishable from 1 as an IEEE double,
> but not as a float.
>

And surely as a long double. I think the OP has a
printing precision problem, i.e. has written
printf("%g\n",cos(1e-7));

--
jacob navia
jacob at jacob point remcomp point fr
logiciels/informatique
http://www.cs.virginia.edu/~lcc-win32

jacob navia, Mar 5, 2008
7. ### Micah CowanGuest

jacob navia <> writes:

> jacob navia wrote:
>> wrote:
>>> Currently using cosine function in math.h
>>>
>>> Currently I get:
>>>
>>> 1 = cos(1e^-7)
>>>
>>> Is there another way for cos to return value of high accuracy say:
>>>
>>> 0.999999 = cos(1e^-7)
>>>

>>
>> If you use the lcc-win compiler you get
>>
>> cos(1e-7)=
>> 0.999999999999995000000000000004166666666666665277777777777778025793650793650766093474426807762228769173214
>>
>> using the extra precision floating point (qfloat).
>>
>>
>> In standard C the biggest floating point type is long double.
>>
>> No way to do that in standard C.
>>

>
> P.S.
>
> WAIT
>
> This is .00000000000000499999999999999, a number that is bigger than
> DBL_EPSILON.
>
> MAYBE your problem is a printing problem?
>
> Which format did you use to print it?

That's rather what I was thinking.

\$ cat foo.c
#include <stdio.h>

int main(void)
{
double foo = 0.999999999;
printf("foo = %f\n", foo);
printf("foo = %.9f\n", foo);
return 0;
}
\$ make foo
cc foo.c -o foo
\$ ./foo
foo = 1.000000
foo = 0.999999999

--
Micah J. Cowan
Programmer, musician, typesetting enthusiast, gamer...
http://micah.cowan.name/

Micah Cowan, Mar 5, 2008
8. ### Micah CowanGuest

jacob navia <> writes:

> Richard Tobin wrote:
>>
>> cos(1e-7) is roughly 1 - 1e-14 / 2, which is closer to 1 than to
>> 0.999999. It should be distinguishable from 1 as an IEEE double,
>> but not as a float.
>>

>
> And surely as a long double. I think the OP has a
> printing precision problem, i.e. has written
> printf("%g\n",cos(1e-7));

Note that the C standard only requires accuracy up to 10 decimal
digits (including for long double). That's why Richard took care to
qualify his statement as "IEEE double".

--
Micah J. Cowan
Programmer, musician, typesetting enthusiast, gamer...
http://micah.cowan.name/

Micah Cowan, Mar 5, 2008
9. ### Martin AmbuhlGuest

wrote:
> Currently using cosine function in math.h
>
> Currently I get:
>
> 1 = cos(1e^-7)
>
> Is there another way for cos to return value of high accuracy say:
>
> 0.999999 = cos(1e^-7)

That is a grossly inaccurate value. Higher precision is not high accuracy:

#include <stdio.h>
#include <math.h>
#include <float.h>

int main(void)
{
int i;
double x = 1.0e-7, y = cos(x);
printf("With default precisions, cos(%g) is shown as:\n"
" (%%f) %f, (%%g) %g, (%%e) %e\n\n", x, y, y, y);
printf("With various precision up to DBL_DIG (%d), cos(%g) is\n",
DBL_DIG, x);
for (i = 0; i <= DBL_DIG; i++)
printf(" precision = %d, %.*f, %.*g, %.*e\n",
i, i, y, i, y, i, y);

x = acos(0.999999);
y = cos(x);
printf("\nWith default precisions, cos(%g) is shown as:\n"
" (%%f) %f, (%%g) %g, (%%e) %e\n\n", x, y, y, y);
printf("With various precision up to DBL_DIG (%d), cos(%g) is\n",
DBL_DIG, x);
for (i = 0; i <= DBL_DIG; i++)
printf(" precision = %d, %.*f, %.*g, %.*e\n",
i, i, y, i, y, i, y);

return 0;
}

With default precisions, cos(1e-07) is shown as:
(%f) 1.000000, (%g) 1, (%e) 1.000000e+00

With various precision up to DBL_DIG (15), cos(1e-07) is
precision = 0, 1, 1, 1e+00
precision = 1, 1.0, 1, 1.0e+00
precision = 2, 1.00, 1, 1.00e+00
precision = 3, 1.000, 1, 1.000e+00
precision = 4, 1.0000, 1, 1.0000e+00
precision = 5, 1.00000, 1, 1.00000e+00
precision = 6, 1.000000, 1, 1.000000e+00
precision = 7, 1.0000000, 1, 1.0000000e+00
precision = 8, 1.00000000, 1, 1.00000000e+00
precision = 9, 1.000000000, 1, 1.000000000e+00
precision = 10, 1.0000000000, 1, 1.0000000000e+00
precision = 11, 1.00000000000, 1, 1.00000000000e+00
precision = 12, 1.000000000000, 1, 1.000000000000e+00
precision = 13, 1.0000000000000, 1, 1.0000000000000e+00
precision = 14, 1.00000000000000, 1, 9.99999999999995e-01
precision = 15, 0.999999999999995, 0.999999999999995,
9.999999999999950e-01

With default precisions, cos(0.00141421) is shown as:
(%f) 0.999999, (%g) 0.999999, (%e) 9.999990e-01

With various precision up to DBL_DIG (15), cos(0.00141421) is
precision = 0, 1, 1, 1e+00
precision = 1, 1.0, 1, 1.0e+00
precision = 2, 1.00, 1, 1.00e+00
precision = 3, 1.000, 1, 1.000e+00
precision = 4, 1.0000, 1, 1.0000e+00
precision = 5, 1.00000, 1, 9.99999e-01
precision = 6, 0.999999, 0.999999, 9.999990e-01
precision = 7, 0.9999990, 0.999999, 9.9999900e-01
precision = 8, 0.99999900, 0.999999, 9.99999000e-01
precision = 9, 0.999999000, 0.999999, 9.999990000e-01
precision = 10, 0.9999990000, 0.999999, 9.9999900000e-01
precision = 11, 0.99999900000, 0.999999, 9.99999000000e-01
precision = 12, 0.999999000000, 0.999999, 9.999990000000e-01
precision = 13, 0.9999990000000, 0.999999, 9.9999900000000e-01
precision = 14, 0.99999900000000, 0.999999, 9.99999000000000e-01
precision = 15, 0.999999000000000, 0.999999, 9.999990000000000e-01

Martin Ambuhl, Mar 6, 2008
10. ### CBFalconerGuest

wrote:
>
> Currently using cosine function in math.h. Currently I get:
>
> 1 = cos(1e^-7)
>
> Is there another way for cos to return value of high accuracy say:
>
> 0.999999 = cos(1e^-7)

Most math function packages evaluate cos (or sin) via fairly long
approximation formula, involving several multiplications and
additions, maybe divisions. If the designer is smart he designs
sin to:

double sin(double x) {
if (x < MAGICVALUE) return x;
else return sinfunction(x);
}

and cos(y) just returns sin(PI/2 - y). Suitable mods for quadrants
needed. Note that MAGICVALUE depends on the actual approximation
etc. used.

--
[mail]: Chuck F (cbfalconer at maineline dot net)
[page]: <http://cbfalconer.home.att.net>

--
Posted via a free Usenet account from http://www.teranews.com

CBFalconer, Mar 6, 2008
11. ### peteGuest

wrote:
>
> Currently using cosine function in math.h
>
> Currently I get:
>
> 1 = cos(1e^-7)
>
> Is there another way for cos to return value of high accuracy say:
>
> 0.999999 = cos(1e^-7)

/* BEGIN new.c output */

LDBL_EPSILON is 2.220446e-016

fs_cosl(1.0 / 10000000) - 1.0 is -4.996004e-015
fs_cosl(1.0 / 1000000) - 1.0 is -5.000445e-013
fs_cosl(1.0 / 100000) - 1.0 is -5.000000e-011
fs_cosl(1.0 / 10000) - 1.0 is -5.000000e-009
fs_cosl(1.0 / 100) - 1.0 is -4.999958e-005

/* END new.c output */

/* BEGIN new.c */

#include <stdio.h>
#include <math.h>
#include <float.h>

#include "fs_cosl.h"

#define DEGREES_TO_RADIANS(A) ((A) * atan(1) / 45)

int main(void)
{
puts("/* BEGIN new.c output */\n");
printf("LDBL_EPSILON is %Le\n\n",
LDBL_EPSILON);
printf("fs_cosl(1.0 / 10000000) - 1.0 is %Le\n",
fs_cosl(1.0 / 10000000) - 1.0);
printf("fs_cosl(1.0 / 1000000) - 1.0 is %Le\n",
fs_cosl(1.0 / 1000000) - 1.0);
printf("fs_cosl(1.0 / 100000) - 1.0 is %Le\n",
fs_cosl(1.0 / 100000) - 1.0);
printf("fs_cosl(1.0 / 10000) - 1.0 is %Le\n",
fs_cosl(1.0 / 10000) - 1.0);
printf("fs_cosl(1.0 / 100) - 1.0 is %Le\n\n",
fs_cosl(1.0 / 100) - 1.0);
puts("/* END new.c output */");
return 0;
}

/* END new.c */

/* BEGIN fs_cosl.h */
/*
** Portable freestanding code.
*/
#ifndef H_FS_COSL_H
#define H_FS_COSL_H

long double fs_cosl(long double x);
long double fs_sqrtl(long double x);
long double fs_fmodl(long double x, long double y);

#endif

/* END fs_cosl.h */

/* BEGIN fs_cosl.c */

#include <float.h>

#include "fs_cosl.h"
/*
** pi == (atan(1.0 / 3) + atan(1.0 / 2)) * 4
*/
static long double fs_pil(void);

long double fs_cosl(long double x)
{
long unsigned n;
int negative, sine;
long double a, b, c;
static long double pi, two_pi, half_pi, third_pi;
static int initialized;

if (0 > x) {
x = -x;
}
if (LDBL_MAX >= x) {
if (!initialized) {
initialized = 1;
pi = fs_pil();
two_pi = 2 * pi;
half_pi = pi / 2;
third_pi = pi / 3;
}
if (x > two_pi) {
x = fs_fmodl(x, two_pi);
}
if (x > pi) {
x = two_pi - x;
}
if (x > half_pi) {
x = pi - x;
negative = 1;
} else {
negative = 0;
}
if (x > third_pi) {
x = half_pi - x;
sine = 1;
} else {
sine = 0;
}
c = x * x;
x = n = 0;
a = 1;
do {
b = a;
a *= c;
a /= ++n;
a /= ++n;
b -= a;
a *= c;
a /= ++n;
a /= ++n;
x += b;
} while (b > LDBL_EPSILON / 2);
if (sine) {
x = fs_sqrtl((1 - x) * (1 + x));
}
if (negative) {
x = -x;
}
} else {
x = -LDBL_MAX;
}
return x;
}

long double fs_sqrtl(long double x)
{
long int n;
long double a, b;

if (x > 0 && LDBL_MAX >= x) {
for (n = 0; x > 2; x /= 4) {
++n;
}
while (0.5 > x) {
--n;
x *= 4;
}
a = x;
b = (1 + x) / 2;
do {
x = b;
b = (a / x + x) / 2;
} while (x > b);
while (n > 0) {
x *= 2;
--n;
}
while (0 > n) {
x /= 2;
++n;
}
} else {
if (x != 0) {
x = LDBL_MAX;
}
}
return x;
}

long double fs_fmodl(long double x, long double y)
{
long double a, b;
const long double c = x;

if (0 > c) {
x = -x;
}
if (0 > y) {
y = -y;
}
if (y != 0 && LDBL_MAX >= y && LDBL_MAX >= x) {
while (x >= y) {
a = x / 2;
b = y;
while (a >= b) {
b *= 2;
}
x -= b;
}
} else {
x = 0;
}
return 0 > c ? -x : x;
}

static long double fs_pil(void)
{
long unsigned n;
long double a, b;
static long double p;
static int initialized;

if (!initialized) {
initialized = 1;
n = 1;
a = 3;
do {
a /= 9;
b = a / n;
n += 2;
a /= 9;
b -= a / n;
n += 2;
p += b;
} while (b > LDBL_EPSILON / 4);
n = 1;
a = 2;
do {
a /= 4;
b = a / n;
n += 2;
a /= 4;
b -= a / n;
n += 2;
p += b;
} while (b > LDBL_EPSILON / 2);
p *= 4;
}
return p;
}

/* END fs_cosl.c */

--
pete

pete, Mar 6, 2008