[unix-history] / .ref-BSD-4_3_Net_2 / usr / src / usr.bin / gcc / cc1 / math-68881.h

/******************************************************************\
*								   *
*  <math-68881.h>		last modified: 18 May 1989.	   *
*								   *
*  Copyright (C) 1989 by Matthew Self.				   *
*  You may freely distribute verbatim copies of this software	   *
*  provided that this copyright notice is retained in all copies.  *
*  You may distribute modifications to this software under the     *
*  conditions above if you also clearly note such modifications    *
*  with their author and date.			   	     	   *
*								   *
*  Note:  errno is not set to EDOM when domain errors occur for    *
*  most of these functions.  Rather, it is assumed that the	   *
*  68881's OPERR exception will be enabled and handled		   *
*  appropriately by the	operating system.  Similarly, overflow	   *
*  and underflow do not set errno to ERANGE.			   *
*								   *
*  Send bugs to Matthew Self (self@bayes.arc.nasa.gov).		   *
*								   *
\******************************************************************/

/* Modified by Richard Stallman, November 1990, to initialize HUGE_VAL
   specially on a Sun.
   December 1989, add parens around `&' in pow.  */

#include <errno.h>

#ifndef HUGE_VAL
#ifdef __sun__
/* The Sun assembler fails to handle the hex constant in the usual defn.  */
#define HUGE_VAL							\
({									\
  static union { int i[2]; double d; } u = { {0x7ff00000, 0} };		\
  u.d;									\
})
#else
#define HUGE_VAL							\
({									\
  double huge_val;							\
									\
  __asm ("fmove%.d #0x7ff0000000000000,%0"	/* Infinity */		\
	 : "=f" (huge_val)						\
	 : /* no inputs */);						\
  huge_val;								\
})
#endif
#endif

__inline static const double sin (double x)
{
  double value;

  __asm ("fsin%.x %1,%0"
	 : "=f" (value)
	 : "f" (x));
  return value;
}

__inline static const double cos (double x)
{
  double value;

  __asm ("fcos%.x %1,%0"
	 : "=f" (value)
	 : "f" (x));
  return value;
}

__inline static const double tan (double x)
{
  double value;

  __asm ("ftan%.x %1,%0"
	 : "=f" (value)
	 : "f" (x));
  return value;
}

__inline static const double asin (double x)
{
  double value;

  __asm ("fasin%.x %1,%0"
	 : "=f" (value)
	 : "f" (x));
  return value;
}

__inline static const double acos (double x)
{
  double value;

  __asm ("facos%.x %1,%0"
	 : "=f" (value)
	 : "f" (x));
  return value;
}

__inline static const double atan (double x)
{
  double value;

  __asm ("fatan%.x %1,%0"
	 : "=f" (value)
	 : "f" (x));
  return value;
}

__inline static const double atan2 (double y, double x)
{
  double pi, pi_over_2;

  __asm ("fmovecr%.x %#0,%0"		/* extended precision pi */
	 : "=f" (pi)
	 : /* no inputs */ );
  __asm ("fscale%.b %#-1,%0"		/* no loss of accuracy */
	 : "=f" (pi_over_2)
	 : "0" (pi));
  if (x > 0)
    {
      if (y > 0)
	{
	  if (x > y)
	    return atan (y / x);
	  else
	    return pi_over_2 - atan (x / y);
	}
      else
	{
	  if (x > -y)
	    return atan (y / x);
	  else
	    return - pi_over_2 - atan (x / y);
	}
    }
  else
    {
      if (y > 0)
	{
	  if (-x > y)
	    return pi + atan (y / x);
	  else
	    return pi_over_2 - atan (x / y);
	}
      else
	{
	  if (-x > -y)
	    return - pi + atan (y / x);
	  else if (y < 0)
	    return - pi_over_2 - atan (x / y);
	  else
	    {
	      double value;

	      errno = EDOM;
	      __asm ("fmove%.d %#0rnan,%0" 	/* quiet NaN */
		     : "=f" (value)
		     : /* no inputs */);
	      return value;
	    }
	}
    }
}

__inline static const double sinh (double x)
{
  double value;

  __asm ("fsinh%.x %1,%0"
	 : "=f" (value)
	 : "f" (x));
  return value;
}

__inline static const double cosh (double x)
{
  double value;

  __asm ("fcosh%.x %1,%0"
	 : "=f" (value)
	 : "f" (x));
  return value;
}

__inline static const double tanh (double x)
{
  double value;

  __asm ("ftanh%.x %1,%0"
	 : "=f" (value)
	 : "f" (x));
  return value;
}

__inline static const double atanh (double x)
{
  double value;

  __asm ("fatanh%.x %1,%0"
	 : "=f" (value)
	 : "f" (x));
  return value;
}

__inline static const double exp (double x)
{
  double value;

  __asm ("fetox%.x %1,%0"
	 : "=f" (value)
	 : "f" (x));
  return value;
}

__inline static const double expm1 (double x)
{
  double value;

  __asm ("fetoxm1%.x %1,%0"
	 : "=f" (value)
	 : "f" (x));
  return value;
}

__inline static const double log (double x)
{
  double value;

  __asm ("flogn%.x %1,%0"
	 : "=f" (value)
	 : "f" (x));
  return value;
}

__inline static const double log1p (double x)
{
  double value;

  __asm ("flognp1%.x %1,%0"
	 : "=f" (value)
	 : "f" (x));
  return value;
}

__inline static const double log10 (double x)
{
  double value;

  __asm ("flog10%.x %1,%0"
	 : "=f" (value)
	 : "f" (x));
  return value;
}

__inline static const double sqrt (double x)
{
  double value;

  __asm ("fsqrt%.x %1,%0"
	 : "=f" (value)
	 : "f" (x));
  return value;
}

__inline static const double pow (const double x, const double y)
{
  if (x > 0)
    return exp (y * log (x));
  else if (x == 0)
    {
      if (y > 0)
	return 0.0;
      else
	{
	  double value;

	  errno = EDOM;
	  __asm ("fmove%.d %#0rnan,%0"		/* quiet NaN */
		 : "=f" (value)
		 : /* no inputs */);
	  return value;
	}
    }
  else
    {
      double temp;

      __asm ("fintrz%.x %1,%0"
	     : "=f" (temp)			/* integer-valued float */
	     : "f" (y));
      if (y == temp)
        {
	  int i = (int) y;
	  
	  if ((i & 1) == 0)			/* even */
	    return exp (y * log (x));
	  else
	    return - exp (y * log (x));
        }
      else
        {
	  double value;

	  errno = EDOM;
	  __asm ("fmove%.d %#0rnan,%0"		/* quiet NaN */
		 : "=f" (value)
		 : /* no inputs */);
	  return value;
        }
    }
}

__inline static const double fabs (double x)
{
  double value;

  __asm ("fabs%.x %1,%0"
	 : "=f" (value)
	 : "f" (x));
  return value;
}

__inline static const double ceil (double x)
{
  int rounding_mode, round_up;
  double value;

  __asm volatile ("fmove%.l fpcr,%0"
		  : "=dm" (rounding_mode)
		  : /* no inputs */ );
  round_up = rounding_mode | 0x30;
  __asm volatile ("fmove%.l %0,fpcr"
		  : /* no outputs */
		  : "dmi" (round_up));
  __asm volatile ("fint%.x %1,%0"
		  : "=f" (value)
		  : "f" (x));
  __asm volatile ("fmove%.l %0,fpcr"
		  : /* no outputs */
		  : "dmi" (rounding_mode));
  return value;
}

__inline static const double floor (double x)
{
  int rounding_mode, round_down;
  double value;

  __asm volatile ("fmove%.l fpcr,%0"
		  : "=dm" (rounding_mode)
		  : /* no inputs */ );
  round_down = (rounding_mode & ~0x10)
		| 0x20;
  __asm volatile ("fmove%.l %0,fpcr"
		  : /* no outputs */
		  : "dmi" (round_down));
  __asm volatile ("fint%.x %1,%0"
		  : "=f" (value)
		  : "f" (x));
  __asm volatile ("fmove%.l %0,fpcr"
		  : /* no outputs */
		  : "dmi" (rounding_mode));
  return value;
}

__inline static const double rint (double x)
{
  int rounding_mode, round_nearest;
  double value;

  __asm volatile ("fmove%.l fpcr,%0"
		  : "=dm" (rounding_mode)
		  : /* no inputs */ );
  round_nearest = rounding_mode & ~0x30;
  __asm volatile ("fmove%.l %0,fpcr"
		  : /* no outputs */
		  : "dmi" (round_nearest));
  __asm volatile ("fint%.x %1,%0"
		  : "=f" (value)
		  : "f" (x));
  __asm volatile ("fmove%.l %0,fpcr"
		  : /* no outputs */
		  : "dmi" (rounding_mode));
  return value;
}

__inline static const double fmod (double x, double y)
{
  double value;

  __asm ("fmod%.x %2,%0"
	 : "=f" (value)
	 : "0" (x),
	   "f" (y));
  return value;
}

__inline static const double drem (double x, double y)
{
  double value;

  __asm ("frem%.x %2,%0"
	 : "=f" (value)
	 : "0" (x),
	   "f" (y));
  return value;
}

__inline static const double scalb (double x, int n)
{
  double value;

  __asm ("fscale%.l %2,%0"
	 : "=f" (value)
	 : "0" (x),
	   "dmi" (n));
  return value;
}

__inline static double logb (double x)
{
  double exponent;

  __asm ("fgetexp%.x %1,%0"
	 : "=f" (exponent)
	 : "f" (x));
  return exponent;
}

__inline static const double ldexp (double x, int n)
{
  double value;

  __asm ("fscale%.l %2,%0"
	 : "=f" (value)
	 : "0" (x),
	   "dmi" (n));
  return value;
}

__inline static double frexp (double x, int *exp)
{
  double float_exponent;
  int int_exponent;
  double mantissa;

  __asm ("fgetexp%.x %1,%0"
	 : "=f" (float_exponent) 	/* integer-valued float */
	 : "f" (x));
  int_exponent = (int) float_exponent;
  __asm ("fgetman%.x %1,%0"
	 : "=f" (mantissa)		/* 1.0 <= mantissa < 2.0 */
	 : "f" (x));
  if (mantissa != 0)
    {
      __asm ("fscale%.b %#-1,%0"
	     : "=f" (mantissa)		/* mantissa /= 2.0 */
	     : "0" (mantissa));
      int_exponent += 1;
    }
  *exp = int_exponent;
  return mantissa;
}

__inline static double modf (double x, double *ip)
{
  double temp;

  __asm ("fintrz%.x %1,%0"
	 : "=f" (temp)			/* integer-valued float */
	 : "f" (x));
  *ip = temp;
  return x - temp;
}
Commit	Line	Data
d0880e41 C	1	/******************************************************************\
	2	* *
	3	* <math-68881.h> last modified: 18 May 1989. *
	4	* *
	5	* Copyright (C) 1989 by Matthew Self. *
	6	* You may freely distribute verbatim copies of this software *
	7	* provided that this copyright notice is retained in all copies. *
	8	* You may distribute modifications to this software under the *
	9	* conditions above if you also clearly note such modifications *
	10	* with their author and date. *
	11	* *
	12	* Note: errno is not set to EDOM when domain errors occur for *
	13	* most of these functions. Rather, it is assumed that the *
	14	* 68881's OPERR exception will be enabled and handled *
	15	* appropriately by the operating system. Similarly, overflow *
	16	* and underflow do not set errno to ERANGE. *
	17	* *
	18	* Send bugs to Matthew Self (self@bayes.arc.nasa.gov). *
	19	* *
	20	\******************************************************************/
	21
	22	/* Modified by Richard Stallman, November 1990, to initialize HUGE_VAL
	23	specially on a Sun.
	24	December 1989, add parens around `&' in pow. */
	25
	26	#include <errno.h>
	27
	28	#ifndef HUGE_VAL
	29	#ifdef __sun__
	30	/* The Sun assembler fails to handle the hex constant in the usual defn. */
	31	#define HUGE_VAL \
	32	({ \
	33	static union { int i[2]; double d; } u = { {0x7ff00000, 0} }; \
	34	u.d; \
	35	})
	36	#else
	37	#define HUGE_VAL \
	38	({ \
	39	double huge_val; \
	40	\
	41	__asm ("fmove%.d #0x7ff0000000000000,%0" /* Infinity */ \
	42	: "=f" (huge_val) \
	43	: /* no inputs */); \
	44	huge_val; \
	45	})
	46	#endif
	47	#endif
	48
	49	__inline static const double sin (double x)
	50	{
	51	double value;
	52
	53	__asm ("fsin%.x %1,%0"
	54	: "=f" (value)
	55	: "f" (x));
	56	return value;
	57	}
	58
	59	__inline static const double cos (double x)
	60	{
	61	double value;
	62
	63	__asm ("fcos%.x %1,%0"
	64	: "=f" (value)
65	: "f" (x));
66	return value;
67	}
68
69	__inline static const double tan (double x)
70	{
71	double value;
72
73	__asm ("ftan%.x %1,%0"
74	: "=f" (value)
75	: "f" (x));
76	return value;
77	}
78
79	__inline static const double asin (double x)
80	{
81	double value;
82
83	__asm ("fasin%.x %1,%0"
84	: "=f" (value)
85	: "f" (x));
86	return value;
87	}
88
89	__inline static const double acos (double x)
90	{
91	double value;
92
93	__asm ("facos%.x %1,%0"
94	: "=f" (value)
95	: "f" (x));
96	return value;
97	}
98
99	__inline static const double atan (double x)
100	{
101	double value;
102
103	__asm ("fatan%.x %1,%0"
104	: "=f" (value)
105	: "f" (x));
106	return value;
107	}
108
109	__inline static const double atan2 (double y, double x)
110	{
111	double pi, pi_over_2;
112
113	__asm ("fmovecr%.x %#0,%0" /* extended precision pi */
114	: "=f" (pi)
115	: /* no inputs */ );
116	__asm ("fscale%.b %#-1,%0" /* no loss of accuracy */
117	: "=f" (pi_over_2)
118	: "0" (pi));
119	if (x > 0)
120	{
121	if (y > 0)
122	{
123	if (x > y)
124	return atan (y / x);
125	else
126	return pi_over_2 - atan (x / y);
127	}
128	else
129	{
130	if (x > -y)
131	return atan (y / x);
132	else
133	return - pi_over_2 - atan (x / y);
134	}
135	}
136	else
137	{
138	if (y > 0)
139	{
140	if (-x > y)
141	return pi + atan (y / x);
142	else
143	return pi_over_2 - atan (x / y);
144	}
145	else
146	{
147	if (-x > -y)
148	return - pi + atan (y / x);
149	else if (y < 0)
150	return - pi_over_2 - atan (x / y);
151	else
152	{
153	double value;
154
155	errno = EDOM;
156	__asm ("fmove%.d %#0rnan,%0" /* quiet NaN */
157	: "=f" (value)
158	: /* no inputs */);
159	return value;
160	}
161	}
162	}
163	}
164
165	__inline static const double sinh (double x)
166	{
167	double value;
168
169	__asm ("fsinh%.x %1,%0"
170	: "=f" (value)
171	: "f" (x));
172	return value;
173	}
174
175	__inline static const double cosh (double x)
176	{
177	double value;
178
179	__asm ("fcosh%.x %1,%0"
180	: "=f" (value)
181	: "f" (x));
182	return value;
183	}
184
185	__inline static const double tanh (double x)
186	{
187	double value;
188
189	__asm ("ftanh%.x %1,%0"
190	: "=f" (value)
191	: "f" (x));
192	return value;
193	}
194
195	__inline static const double atanh (double x)
196	{
197	double value;
198
199	__asm ("fatanh%.x %1,%0"
200	: "=f" (value)
201	: "f" (x));
202	return value;
203	}
204
205	__inline static const double exp (double x)
206	{
207	double value;
208
209	__asm ("fetox%.x %1,%0"
210	: "=f" (value)
211	: "f" (x));
212	return value;
213	}
214
215	__inline static const double expm1 (double x)
216	{
217	double value;
218
219	__asm ("fetoxm1%.x %1,%0"
220	: "=f" (value)
221	: "f" (x));
222	return value;
223	}
224
225	__inline static const double log (double x)
226	{
227	double value;
228
229	__asm ("flogn%.x %1,%0"
230	: "=f" (value)
231	: "f" (x));
232	return value;
233	}
234
235	__inline static const double log1p (double x)
236	{
237	double value;
238
239	__asm ("flognp1%.x %1,%0"
240	: "=f" (value)
241	: "f" (x));
242	return value;
243	}
244
245	__inline static const double log10 (double x)
246	{
247	double value;
248
249	__asm ("flog10%.x %1,%0"
250	: "=f" (value)
251	: "f" (x));
252	return value;
253	}
254
255	__inline static const double sqrt (double x)
256	{
257	double value;
258
259	__asm ("fsqrt%.x %1,%0"
260	: "=f" (value)
261	: "f" (x));
262	return value;
263	}
264
265	__inline static const double pow (const double x, const double y)
266	{
267	if (x > 0)
268	return exp (y * log (x));
269	else if (x == 0)
270	{
271	if (y > 0)
272	return 0.0;
273	else
274	{
275	double value;
276
277	errno = EDOM;
278	__asm ("fmove%.d %#0rnan,%0" /* quiet NaN */
279	: "=f" (value)
280	: /* no inputs */);
281	return value;
282	}
283	}
284	else
285	{
286	double temp;
287
288	__asm ("fintrz%.x %1,%0"
289	: "=f" (temp) /* integer-valued float */
290	: "f" (y));
291	if (y == temp)
292	{
293	int i = (int) y;
294
295	if ((i & 1) == 0) /* even */
296	return exp (y * log (x));
297	else
298	return - exp (y * log (x));
299	}
300	else
301	{
302	double value;
303
304	errno = EDOM;
305	__asm ("fmove%.d %#0rnan,%0" /* quiet NaN */
306	: "=f" (value)
307	: /* no inputs */);
308	return value;
309	}
310	}
311	}
312
313	__inline static const double fabs (double x)
314	{
315	double value;
316
317	__asm ("fabs%.x %1,%0"
318	: "=f" (value)
319	: "f" (x));
320	return value;
321	}
322
323	__inline static const double ceil (double x)
324	{
325	int rounding_mode, round_up;
326	double value;
327
328	__asm volatile ("fmove%.l fpcr,%0"
329	: "=dm" (rounding_mode)
330	: /* no inputs */ );
331	round_up = rounding_mode \| 0x30;
332	__asm volatile ("fmove%.l %0,fpcr"
333	: /* no outputs */
334	: "dmi" (round_up));
335	__asm volatile ("fint%.x %1,%0"
336	: "=f" (value)
337	: "f" (x));
338	__asm volatile ("fmove%.l %0,fpcr"
339	: /* no outputs */
340	: "dmi" (rounding_mode));
341	return value;
342	}
343
344	__inline static const double floor (double x)
345	{
346	int rounding_mode, round_down;
347	double value;
348
349	__asm volatile ("fmove%.l fpcr,%0"
350	: "=dm" (rounding_mode)
351	: /* no inputs */ );
352	round_down = (rounding_mode & ~0x10)
353	\| 0x20;
354	__asm volatile ("fmove%.l %0,fpcr"
355	: /* no outputs */
356	: "dmi" (round_down));
357	__asm volatile ("fint%.x %1,%0"
358	: "=f" (value)
359	: "f" (x));
360	__asm volatile ("fmove%.l %0,fpcr"
361	: /* no outputs */
362	: "dmi" (rounding_mode));
363	return value;
364	}
365
366	__inline static const double rint (double x)
367	{
368	int rounding_mode, round_nearest;
369	double value;
370
371	__asm volatile ("fmove%.l fpcr,%0"
372	: "=dm" (rounding_mode)
373	: /* no inputs */ );
374	round_nearest = rounding_mode & ~0x30;
375	__asm volatile ("fmove%.l %0,fpcr"
376	: /* no outputs */
377	: "dmi" (round_nearest));
378	__asm volatile ("fint%.x %1,%0"
379	: "=f" (value)
380	: "f" (x));
381	__asm volatile ("fmove%.l %0,fpcr"
382	: /* no outputs */
383	: "dmi" (rounding_mode));
384	return value;
385	}
386
387	__inline static const double fmod (double x, double y)
388	{
389	double value;
390
391	__asm ("fmod%.x %2,%0"
392	: "=f" (value)
393	: "0" (x),
394	"f" (y));
395	return value;
396	}
397
398	__inline static const double drem (double x, double y)
399	{
400	double value;
401
402	__asm ("frem%.x %2,%0"
403	: "=f" (value)
404	: "0" (x),
405	"f" (y));
406	return value;
407	}
408
409	__inline static const double scalb (double x, int n)
410	{
411	double value;
412
413	__asm ("fscale%.l %2,%0"
414	: "=f" (value)
415	: "0" (x),
416	"dmi" (n));
417	return value;
418	}
419
420	__inline static double logb (double x)
421	{
422	double exponent;
423
424	__asm ("fgetexp%.x %1,%0"
425	: "=f" (exponent)
426	: "f" (x));
427	return exponent;
428	}
429
430	__inline static const double ldexp (double x, int n)
431	{
432	double value;
433
434	__asm ("fscale%.l %2,%0"
435	: "=f" (value)
436	: "0" (x),
437	"dmi" (n));
438	return value;
439	}
440
441	__inline static double frexp (double x, int *exp)
442	{
443	double float_exponent;
444	int int_exponent;
445	double mantissa;
446
447	__asm ("fgetexp%.x %1,%0"
448	: "=f" (float_exponent) /* integer-valued float */
449	: "f" (x));
450	int_exponent = (int) float_exponent;
451	__asm ("fgetman%.x %1,%0"
452	: "=f" (mantissa) /* 1.0 <= mantissa < 2.0 */
453	: "f" (x));
454	if (mantissa != 0)
455	{
456	__asm ("fscale%.b %#-1,%0"
457	: "=f" (mantissa) /* mantissa /= 2.0 */
458	: "0" (mantissa));
459	int_exponent += 1;
460	}
461	*exp = int_exponent;
462	return mantissa;
463	}
464
465	__inline static double modf (double x, double *ip)
466	{
467	double temp;
468
469	__asm ("fintrz%.x %1,%0"
470	: "=f" (temp) /* integer-valued float */
471	: "f" (x));
472	*ip = temp;
473	return x - temp;
474	}
475