1 /* 2 * Copyright (c) 1987 Regents of the University of California. 3 * All rights reserved. 4 * 5 * Redistribution and use in source and binary forms are permitted 6 * provided that this notice is preserved and that due credit is given 7 * to the University of California at Berkeley. The name of the University 8 * may not be used to endorse or promote products derived from this 9 * software without specific prior written permission. This software 10 * is provided ``as is'' without express or implied warranty. 11 * 12 * All recipients should regard themselves as participants in an ongoing 13 * research project and hence should feel obligated to report their 14 * experiences (good or bad) with these elementary function codes, using 15 * the sendbug(8) program, to the authors. 16 */ 17 18 #ifndef lint 19 static char sccsid[] = "@(#)sincos.c 5.2 (Berkeley) 04/29/88"; 20 #endif /* not lint */ 21 22 #include "trig.h" 23 double 24 sin(x) 25 double x; 26 { 27 double a,c,z; 28 29 if(!finite(x)) /* sin(NaN) and sin(INF) must be NaN */ 30 return x-x; 31 x=drem(x,PI2); /* reduce x into [-PI,PI] */ 32 a=copysign(x,one); 33 if (a >= PIo4) { 34 if(a >= PI3o4) /* ... in [3PI/4,PI] */ 35 x = copysign((a = PI-a),x); 36 else { /* ... in [PI/4,3PI/4] */ 37 a = PIo2-a; /* rtn. sign(x)*C(PI/2-|x|) */ 38 z = a*a; 39 c = cos__C(z); 40 z *= half; 41 a = (z >= thresh ? half-((z-half)-c) : one-(z-c)); 42 return copysign(a,x); 43 } 44 } 45 46 if (a < small) { /* rtn. S(x) */ 47 big+a; 48 return x; 49 } 50 return x+x*sin__S(x*x); 51 } 52 53 double 54 cos(x) 55 double x; 56 { 57 double a,c,z,s = 1.0; 58 59 if(!finite(x)) /* cos(NaN) and cos(INF) must be NaN */ 60 return x-x; 61 x=drem(x,PI2); /* reduce x into [-PI,PI] */ 62 a=copysign(x,one); 63 if (a >= PIo4) { 64 if (a >= PI3o4) { /* ... in [3PI/4,PI] */ 65 a = PI-a; 66 s = negone; 67 } 68 else { /* ... in [PI/4,3PI/4] */ 69 a = PIo2-a; 70 return a+a*sin__S(a*a); /* rtn. S(PI/2-|x|) */ 71 } 72 } 73 if (a < small) { 74 big+a; 75 return s; /* rtn. s*C(a) */ 76 } 77 z = a*a; 78 c = cos__C(z); 79 z *= half; 80 a = (z >= thresh ? half-((z-half)-c) : one-(z-c)); 81 return copysign(a,s); 82 } 83