Orange Pi5 kernel

^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300   1) /*
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300   2) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300   3)   fp_trig.c: floating-point math routines for the Linux-m68k
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300   4)   floating point emulator.
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300   5) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300   6)   Copyright (c) 1998-1999 David Huggins-Daines / Roman Zippel.
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300   7) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300   8)   I hereby give permission, free of charge, to copy, modify, and
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300   9)   redistribute this software, in source or binary form, provided that
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  10)   the above copyright notice and the following disclaimer are included
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  11)   in all such copies.
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  12) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  13)   THIS SOFTWARE IS PROVIDED "AS IS", WITH ABSOLUTELY NO WARRANTY, REAL
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  14)   OR IMPLIED.
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  15) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  16) */
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  17) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  18) #include "fp_emu.h"
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  19) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  20) static const struct fp_ext fp_one =
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  21) {
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  22) 	.exp = 0x3fff,
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  23) };
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  24) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  25) extern struct fp_ext *fp_fadd(struct fp_ext *dest, const struct fp_ext *src);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  26) extern struct fp_ext *fp_fdiv(struct fp_ext *dest, const struct fp_ext *src);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  27) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  28) struct fp_ext *
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  29) fp_fsqrt(struct fp_ext *dest, struct fp_ext *src)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  30) {
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  31) 	struct fp_ext tmp, src2;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  32) 	int i, exp;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  33) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  34) 	dprint(PINSTR, "fsqrt\n");
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  35) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  36) 	fp_monadic_check(dest, src);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  37) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  38) 	if (IS_ZERO(dest))
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  39) 		return dest;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  40) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  41) 	if (dest->sign) {
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  42) 		fp_set_nan(dest);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  43) 		return dest;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  44) 	}
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  45) 	if (IS_INF(dest))
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  46) 		return dest;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  47) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  48) 	/*
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  49) 	 *		 sqrt(m) * 2^(p)	, if e = 2*p
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  50) 	 * sqrt(m*2^e) =
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  51) 	 *		 sqrt(2*m) * 2^(p)	, if e = 2*p + 1
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  52) 	 *
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  53) 	 * So we use the last bit of the exponent to decide whether to
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  54) 	 * use the m or 2*m.
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  55) 	 *
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  56) 	 * Since only the fractional part of the mantissa is stored and
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  57) 	 * the integer part is assumed to be one, we place a 1 or 2 into
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  58) 	 * the fixed point representation.
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  59) 	 */
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  60) 	exp = dest->exp;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  61) 	dest->exp = 0x3FFF;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  62) 	if (!(exp & 1))		/* lowest bit of exponent is set */
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  63) 		dest->exp++;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  64) 	fp_copy_ext(&src2, dest);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  65) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  66) 	/*
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  67) 	 * The taylor row around a for sqrt(x) is:
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  68) 	 *	sqrt(x) = sqrt(a) + 1/(2*sqrt(a))*(x-a) + R
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  69) 	 * With a=1 this gives:
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  70) 	 *	sqrt(x) = 1 + 1/2*(x-1)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  71) 	 *		= 1/2*(1+x)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  72) 	 */
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  73) 	fp_fadd(dest, &fp_one);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  74) 	dest->exp--;		/* * 1/2 */
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  75) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  76) 	/*
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  77) 	 * We now apply the newton rule to the function
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  78) 	 *	f(x) := x^2 - r
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  79) 	 * which has a null point on x = sqrt(r).
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  80) 	 *
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  81) 	 * It gives:
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  82) 	 *	x' := x - f(x)/f'(x)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  83) 	 *	    = x - (x^2 -r)/(2*x)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  84) 	 *	    = x - (x - r/x)/2
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  85) 	 *          = (2*x - x + r/x)/2
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  86) 	 *	    = (x + r/x)/2
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  87) 	 */
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  88) 	for (i = 0; i < 9; i++) {
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  89) 		fp_copy_ext(&tmp, &src2);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  90) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  91) 		fp_fdiv(&tmp, dest);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  92) 		fp_fadd(dest, &tmp);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  93) 		dest->exp--;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  94) 	}
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  95) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  96) 	dest->exp += (exp - 0x3FFF) / 2;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  97) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  98) 	return dest;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300  99) }
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 100) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 101) struct fp_ext *
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 102) fp_fetoxm1(struct fp_ext *dest, struct fp_ext *src)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 103) {
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 104) 	uprint("fetoxm1\n");
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 105) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 106) 	fp_monadic_check(dest, src);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 107) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 108) 	return dest;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 109) }
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 110) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 111) struct fp_ext *
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 112) fp_fetox(struct fp_ext *dest, struct fp_ext *src)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 113) {
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 114) 	uprint("fetox\n");
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 115) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 116) 	fp_monadic_check(dest, src);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 117) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 118) 	return dest;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 119) }
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 120) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 121) struct fp_ext *
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 122) fp_ftwotox(struct fp_ext *dest, struct fp_ext *src)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 123) {
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 124) 	uprint("ftwotox\n");
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 125) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 126) 	fp_monadic_check(dest, src);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 127) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 128) 	return dest;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 129) }
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 130) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 131) struct fp_ext *
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 132) fp_ftentox(struct fp_ext *dest, struct fp_ext *src)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 133) {
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 134) 	uprint("ftentox\n");
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 135) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 136) 	fp_monadic_check(dest, src);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 137) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 138) 	return dest;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 139) }
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 140) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 141) struct fp_ext *
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 142) fp_flogn(struct fp_ext *dest, struct fp_ext *src)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 143) {
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 144) 	uprint("flogn\n");
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 145) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 146) 	fp_monadic_check(dest, src);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 147) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 148) 	return dest;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 149) }
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 150) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 151) struct fp_ext *
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 152) fp_flognp1(struct fp_ext *dest, struct fp_ext *src)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 153) {
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 154) 	uprint("flognp1\n");
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 155) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 156) 	fp_monadic_check(dest, src);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 157) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 158) 	return dest;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 159) }
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 160) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 161) struct fp_ext *
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 162) fp_flog10(struct fp_ext *dest, struct fp_ext *src)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 163) {
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 164) 	uprint("flog10\n");
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 165) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 166) 	fp_monadic_check(dest, src);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 167) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 168) 	return dest;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 169) }
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 170) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 171) struct fp_ext *
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 172) fp_flog2(struct fp_ext *dest, struct fp_ext *src)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 173) {
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 174) 	uprint("flog2\n");
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 175) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 176) 	fp_monadic_check(dest, src);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 177) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 178) 	return dest;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 179) }
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 180) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 181) struct fp_ext *
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 182) fp_fgetexp(struct fp_ext *dest, struct fp_ext *src)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 183) {
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 184) 	dprint(PINSTR, "fgetexp\n");
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 185) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 186) 	fp_monadic_check(dest, src);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 187) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 188) 	if (IS_INF(dest)) {
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 189) 		fp_set_nan(dest);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 190) 		return dest;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 191) 	}
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 192) 	if (IS_ZERO(dest))
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 193) 		return dest;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 194) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 195) 	fp_conv_long2ext(dest, (int)dest->exp - 0x3FFF);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 196) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 197) 	fp_normalize_ext(dest);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 198) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 199) 	return dest;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 200) }
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 201) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 202) struct fp_ext *
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 203) fp_fgetman(struct fp_ext *dest, struct fp_ext *src)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 204) {
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 205) 	dprint(PINSTR, "fgetman\n");
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 206) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 207) 	fp_monadic_check(dest, src);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 208) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 209) 	if (IS_ZERO(dest))
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 210) 		return dest;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 211) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 212) 	if (IS_INF(dest))
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 213) 		return dest;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 214) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 215) 	dest->exp = 0x3FFF;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 216) 
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 217) 	return dest;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 218) }
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 219)