## Re: [Mingw-w64-public] [RFC] sinf and cosf for ARM

 Re: [Mingw-w64-public] [RFC] sinf and cosf for ARM From: André Hentschel - 2014-06-20 19:26:07 ```-----BEGIN PGP SIGNED MESSAGE----- Hash: SHA1 Am 19.06.2014 22:36, schrieb Pavel: > Hi André, Hi, wow, thanks for having such a close look. > I have some comments to the sinf implementation. First of all, > calculating the factorials in separate function is quite inefficient, > since you repeat multiplication of what has already been multiplied. yeah, i guess that's even true when making it inline... > Slightly better would be to declare local variable and calculate the > factorial in the for loop: > > float fact = 1; > float y = 1; > float res = y; > for(int i = 1; i < 6; i++) > { > y *= x; > fact *= i; > res += y/fact; > } i think this could make the code more unreadable > this is rather an example for exp function, sine would require a > modification. And since you only calculate 6 values, it would be even > better to precalculate them and hardcode them in a table. yes, that's my best bet i guess, too. > However, calculating goniometric functions from Taylor expansion is > quite inefficient and slow in general. In your algorithm, you use 6 > values and power of 11, which itself introduces some error. The accuracy > of your implementation is 4 decimal digits at pi/2, where I believe the > deviation is the highest. > > I suggest the following formula for calculating sine at (0, pi/2) > interval: > > sin(x) = (0.99996736*x - 0.09779183*x^2 - 0.12839086*x^3 + > 0.01826627*x^4) / (1 - 0.09809942*x + 0.03936879*x^2) > > This is accurate at 6 decimal digits in the whole interval, moreover it > matches the values at 0, pi/12, pi/6, pi/4, pi/3, 5*pi/12 and pi/2 > exactly. The highest deviation magnitude is about 7.8e-7. I have derived > the coefficients using R. could you please attach some R code for this? My plan for now is to make those functions work, not spending too much time optimizing precision and performance. After that me and others could start improving both. So you can try with your formula later, hopefully also for other functions. > cosine can be then calculated as sin(pi/2 - x) at the (0, pi/2) > interval. > > So I hope it might help you. It does, thanks! If you wanna continue looking into that, i put my code here: https://github.com/AndreRH/mingw-w64/commits/master Note that you sometimes get in trouble when fetching from it, as i might do forced pushes when needed. -----BEGIN PGP SIGNATURE----- Version: GnuPG v1 Comment: Using GnuPG with Thunderbird - http://www.enigmail.net/ iQGcBAEBAgAGBQJTpIrCAAoJEGm5GZTakYsseOML/3Lry/5hzegZBCC6LJAk397F qx8E59Y2lr6cSRgsp0HudkErh84+YALkb5tFxX24yHQwPs6sbxL9JClpReXJJNYH val9Vkj8lV0BiG2kSlPE1nPqqTUbHUbkCRc2mfylTzT1v2cLWNgkZhM+eUqtp1Xp j0tsyBCR0L6LU5SOFaqbpgnMn4A2DO/fN3GyPwZDubQ04iZXw/6toCWiN8RzRr8/ 8BBHCdMtJghS9l8LoEjOTD/ffSy6Pz7CEre7iuFRZA+LUcfncrstaMebAjjkui9p BBgtBQUas9gl6GtovmdD1jwFUJZyr7pezqOQTq2dBI9jmartYPogqc/+t1DBnJdf qpYigBC2NFzqSnOMN2E6U7FVCuQIS+IIzMNVKcqiw22V1JJsUNKQL9jsWcHT6ahr fmeMkdJyFjIdPDKo6gT6/7Jw0UzxDeFhx9iUFM79ieF/cjdNFqv72hyVNoNSSNZQ +xGDU05ahsSiuwLenh99yHS9iR76VCbbyl5n/TBX9w== =h95E -----END PGP SIGNATURE----- ```

 [Mingw-w64-public] [RFC] sinf and cosf for ARM From: André Hentschel - 2014-06-16 19:17:07 Attachments: tmp.diff tmp.diff.sig ```-----BEGIN PGP SIGNED MESSAGE----- Hash: SHA1 Please review this initial implementation for the remaining math functions I'm unsure about the folder name and the license for taylor_private.h, it mostly contains defines/typedefs/consts from FreeBSD... -----BEGIN PGP SIGNATURE----- Version: GnuPG v1 Comment: Using GnuPG with Thunderbird - http://www.enigmail.net/ iQGcBAEBAgAGBQJTn0KlAAoJEGm5GZTakYssggYL/iNc6+WuWR8Nx+k06Zjy4yxl CCB3/c9AL/uqSZdUaqHr0HGTZsy6ny8iLXUy082XnY5wE/87QInWPZ7hZjJ/5WKc 3pay4H1FufDKfrHJRBkcWD1z6PDiYNWboTKkHVOrNLcA6rQV8iOjFocIdepQjiF7 nIJDx+/i+YQ2VDNt3KgHhYnuR3smZFC2vhKSOR61XVx72v+CZ78IBE072ONvJyfM UxKa/FNEzKx5yF+ebC60nX15KldS+xBGcnsSsq8rdQpRN8dog0Z/vfQepV3yrzLV mx3s7nR3wU3b4j6rnxwnAkIAqfjkZk6JbPVB6I+rx+zKekkg966ZkNZX/TQeh2BM WyX+jlnRpwolLRPXN77Ywt/FQ6UJ8S8AoCaIWh7rpPtBADUh6z/mZxhH1F+h0Zs/ rTSE6w+dGWj3Szh1K/boqSdM6RdBg8gqrDr7SzNFy2XKgvuRjp2O/nUT3SB/nJHh nNoDvEt4LerCYnplOa/rfFDBd9kLAQiSacVw1iCPow== =4G6y -----END PGP SIGNATURE----- ```
 Re: [Mingw-w64-public] [RFC] sinf and cosf for ARM From: Kai Tietz - 2014-06-16 19:51:32 ```2014-06-16 21:16 GMT+02:00 André Hentschel : > -----BEGIN PGP SIGNED MESSAGE----- > Hash: SHA1 > > Please review this initial implementation for the remaining math functions > I'm unsure about the folder name and the license for taylor_private.h, it mostly contains > defines/typedefs/consts from FreeBSD... Well, license is BSD. That's actuallly ok. What is about the accuracy of the sinf/cosf implementation? It might be of interest to replace for x64 the current x87 implementation by soft-float, too. Kai ```
 Re: [Mingw-w64-public] [RFC] sinf and cosf for ARM From: André Hentschel - 2014-06-16 21:04:22 ```-----BEGIN PGP SIGNED MESSAGE----- Hash: SHA1 Am 16.06.2014 21:51, schrieb Kai Tietz: > 2014-06-16 21:16 GMT+02:00 André Hentschel : >> -----BEGIN PGP SIGNED MESSAGE----- >> Hash: SHA1 >> >> Please review this initial implementation for the remaining math functions >> I'm unsure about the folder name and the license for taylor_private.h, it mostly contains >> defines/typedefs/consts from FreeBSD... > > Well, license is BSD. That's actuallly ok. ok, true > What is about the accuracy of the sinf/cosf implementation? > It might be of interest to replace for x64 the current x87 > implementation by soft-float, too. As discussed on IRC i'll continue the plan of bringing the algorithms and build system changes in and later we see what we can do for amd64 -----BEGIN PGP SIGNATURE----- Version: GnuPG v1 Comment: Using GnuPG with Thunderbird - http://www.enigmail.net/ iQGcBAEBAgAGBQJTn1vLAAoJEGm5GZTakYssUFML/1WN9ff0SZLmMy5j8pgPIO3h lHIwEPL8Bau+aHr/TNlah7ihnKToXgqGHlOkP6y757ZAdr1t5AZoOKmGILVFJ8xz 0xBSWp3LNCHtivBvz4x9AFLiTQCvXuHkpFgl6Kcw/sec3n6++ZN4OxOg8InoKgXs JTz5ylLcwLiwTftWa+2Twx90q+TVw61AbOCb7KK7EiA2cDZCBNf7b+6tBIH5BcdR xKjJFsuW+Wqdp4+IvHsI2zKm3Dx2IfJPn1q/Ht1RCJE/Ly5SMYzBKAtzrgJt63tj RkOltTUADuM29uKoR2IvAFRNcRy7fbkMLtcqd0u8yJXXMfWeyrXVQqKm0PjNQVYc 4JFqhfrW3ra7tCfPZmC9d5dN8EGj45nLSJGw7kTpuDKcAj95Nr2NmlI5/2IYe03E uvOGynAGh4dNVJH6JELYFJ1B7KJSKHaBylv7Dvajg18vrDlQcoi8lsWSMXFloPRB xYTl//LvfKHGpxBmv8AHxR/G/ZCS8FNVHb6yxJXpZA== =vomx -----END PGP SIGNATURE----- ```
 Re: [Mingw-w64-public] [RFC] sinf and cosf for ARM From: Pavel - 2014-06-19 20:55:22 ```Hi André, I have some comments to the sinf implementation. First of all, calculating the factorials in separate function is quite inefficient, since you repeat multiplication of what has already been multiplied. Slightly better would be to declare local variable and calculate the factorial in the for loop: float fact = 1; float y = 1; float res = y; for(int i = 1; i < 6; i++) { y *= x; fact *= i; res += y/fact; } this is rather an example for exp function, sine would require a modification. And since you only calculate 6 values, it would be even better to precalculate them and hardcode them in a table. However, calculating goniometric functions from Taylor expansion is quite inefficient and slow in general. In your algorithm, you use 6 values and power of 11, which itself introduces some error. The accuracy of your implementation is 4 decimal digits at pi/2, where I believe the deviation is the highest. I suggest the following formula for calculating sine at (0, pi/2) interval: sin(x) = (0.99996736*x - 0.09779183*x^2 - 0.12839086*x^3 + 0.01826627*x^4) / (1 - 0.09809942*x + 0.03936879*x^2) This is accurate at 6 decimal digits in the whole interval, moreover it matches the values at 0, pi/12, pi/6, pi/4, pi/3, 5*pi/12 and pi/2 exactly. The highest deviation magnitude is about 7.8e-7. I have derived the coefficients using R. cosine can be then calculated as sin(pi/2 - x) at the (0, pi/2) interval. So I hope it might help you. Pavel On Mon, 2014-06-16 at 23:04 +0200, André Hentschel wrote: > -----BEGIN PGP SIGNED MESSAGE----- > Hash: SHA1 > > Am 16.06.2014 21:51, schrieb Kai Tietz: > > 2014-06-16 21:16 GMT+02:00 André Hentschel : > >> -----BEGIN PGP SIGNED MESSAGE----- > >> Hash: SHA1 > >> > >> Please review this initial implementation for the remaining math functions > >> I'm unsure about the folder name and the license for taylor_private.h, it mostly contains > >> defines/typedefs/consts from FreeBSD... > > > > Well, license is BSD. That's actuallly ok. > > ok, true > > > What is about the accuracy of the sinf/cosf implementation? > > It might be of interest to replace for x64 the current x87 > > implementation by soft-float, too. > > As discussed on IRC i'll continue the plan of bringing the algorithms and build system changes in and later we see what we can do for amd64 > > -----BEGIN PGP SIGNATURE----- > Version: GnuPG v1 > Comment: Using GnuPG with Thunderbird - http://www.enigmail.net/ > > iQGcBAEBAgAGBQJTn1vLAAoJEGm5GZTakYssUFML/1WN9ff0SZLmMy5j8pgPIO3h > lHIwEPL8Bau+aHr/TNlah7ihnKToXgqGHlOkP6y757ZAdr1t5AZoOKmGILVFJ8xz > 0xBSWp3LNCHtivBvz4x9AFLiTQCvXuHkpFgl6Kcw/sec3n6++ZN4OxOg8InoKgXs > JTz5ylLcwLiwTftWa+2Twx90q+TVw61AbOCb7KK7EiA2cDZCBNf7b+6tBIH5BcdR > xKjJFsuW+Wqdp4+IvHsI2zKm3Dx2IfJPn1q/Ht1RCJE/Ly5SMYzBKAtzrgJt63tj > RkOltTUADuM29uKoR2IvAFRNcRy7fbkMLtcqd0u8yJXXMfWeyrXVQqKm0PjNQVYc > 4JFqhfrW3ra7tCfPZmC9d5dN8EGj45nLSJGw7kTpuDKcAj95Nr2NmlI5/2IYe03E > uvOGynAGh4dNVJH6JELYFJ1B7KJSKHaBylv7Dvajg18vrDlQcoi8lsWSMXFloPRB > xYTl//LvfKHGpxBmv8AHxR/G/ZCS8FNVHb6yxJXpZA== > =vomx > -----END PGP SIGNATURE----- > > ------------------------------------------------------------------------------ > HPCC Systems Open Source Big Data Platform from LexisNexis Risk Solutions > Find What Matters Most in Your Big Data with HPCC Systems > Open Source. Fast. Scalable. Simple. Ideal for Dirty Data. > Leverages Graph Analysis for Fast Processing & Easy Data Exploration > http://p.sf.net/sfu/hpccsystems > _______________________________________________ > Mingw-w64-public mailing list > Mingw-w64-public@... > https://lists.sourceforge.net/lists/listinfo/mingw-w64-public ```
 Re: [Mingw-w64-public] [RFC] sinf and cosf for ARM From: André Hentschel - 2014-06-20 19:26:07 ```-----BEGIN PGP SIGNED MESSAGE----- Hash: SHA1 Am 19.06.2014 22:36, schrieb Pavel: > Hi André, Hi, wow, thanks for having such a close look. > I have some comments to the sinf implementation. First of all, > calculating the factorials in separate function is quite inefficient, > since you repeat multiplication of what has already been multiplied. yeah, i guess that's even true when making it inline... > Slightly better would be to declare local variable and calculate the > factorial in the for loop: > > float fact = 1; > float y = 1; > float res = y; > for(int i = 1; i < 6; i++) > { > y *= x; > fact *= i; > res += y/fact; > } i think this could make the code more unreadable > this is rather an example for exp function, sine would require a > modification. And since you only calculate 6 values, it would be even > better to precalculate them and hardcode them in a table. yes, that's my best bet i guess, too. > However, calculating goniometric functions from Taylor expansion is > quite inefficient and slow in general. In your algorithm, you use 6 > values and power of 11, which itself introduces some error. The accuracy > of your implementation is 4 decimal digits at pi/2, where I believe the > deviation is the highest. > > I suggest the following formula for calculating sine at (0, pi/2) > interval: > > sin(x) = (0.99996736*x - 0.09779183*x^2 - 0.12839086*x^3 + > 0.01826627*x^4) / (1 - 0.09809942*x + 0.03936879*x^2) > > This is accurate at 6 decimal digits in the whole interval, moreover it > matches the values at 0, pi/12, pi/6, pi/4, pi/3, 5*pi/12 and pi/2 > exactly. The highest deviation magnitude is about 7.8e-7. I have derived > the coefficients using R. could you please attach some R code for this? My plan for now is to make those functions work, not spending too much time optimizing precision and performance. After that me and others could start improving both. So you can try with your formula later, hopefully also for other functions. > cosine can be then calculated as sin(pi/2 - x) at the (0, pi/2) > interval. > > So I hope it might help you. It does, thanks! If you wanna continue looking into that, i put my code here: https://github.com/AndreRH/mingw-w64/commits/master Note that you sometimes get in trouble when fetching from it, as i might do forced pushes when needed. -----BEGIN PGP SIGNATURE----- Version: GnuPG v1 Comment: Using GnuPG with Thunderbird - http://www.enigmail.net/ iQGcBAEBAgAGBQJTpIrCAAoJEGm5GZTakYsseOML/3Lry/5hzegZBCC6LJAk397F qx8E59Y2lr6cSRgsp0HudkErh84+YALkb5tFxX24yHQwPs6sbxL9JClpReXJJNYH val9Vkj8lV0BiG2kSlPE1nPqqTUbHUbkCRc2mfylTzT1v2cLWNgkZhM+eUqtp1Xp j0tsyBCR0L6LU5SOFaqbpgnMn4A2DO/fN3GyPwZDubQ04iZXw/6toCWiN8RzRr8/ 8BBHCdMtJghS9l8LoEjOTD/ffSy6Pz7CEre7iuFRZA+LUcfncrstaMebAjjkui9p BBgtBQUas9gl6GtovmdD1jwFUJZyr7pezqOQTq2dBI9jmartYPogqc/+t1DBnJdf qpYigBC2NFzqSnOMN2E6U7FVCuQIS+IIzMNVKcqiw22V1JJsUNKQL9jsWcHT6ahr fmeMkdJyFjIdPDKo6gT6/7Jw0UzxDeFhx9iUFM79ieF/cjdNFqv72hyVNoNSSNZQ +xGDU05ahsSiuwLenh99yHS9iR76VCbbyl5n/TBX9w== =h95E -----END PGP SIGNATURE----- ```
 Re: [Mingw-w64-public] [RFC] sinf and cosf for ARM From: Pavel - 2014-06-20 20:49:34 ```> > I suggest the following formula for calculating sine at (0, pi/2) > > interval: > > > > sin(x) = (0.99996736*x - 0.09779183*x^2 - 0.12839086*x^3 + > > 0.01826627*x^4) / (1 - 0.09809942*x + 0.03936879*x^2) > > > could you please attach some R code for this? > My plan for now is to make those functions work, not spending too much time optimizing precision and performance. > After that me and others could start improving both. So you can try with your formula later, hopefully also for other functions. Well, the idea behind the formula is quite simple. The aim is to find rational function, which match sine at predefined number of points. Thus: sin(x)*(a0+a1*x+a2*x^2+...) = b0+b1*x+b2*x^2+... from x=0 it follows that b0=0 and if we create linear system from the predefined points, we find that the matrix has less rank than the number of rows/columns. So we can choose one more parameter and the natural choice is a0=1. Now let's choose the points we know the values of sine, the R command is: x0 <- c(pi/12, pi/6, pi/4, pi/3, 5*pi/12, pi/2) the corresponding sine values are: b <- c(sqrt((1 - sqrt(3)/2)/2), 1/2, sqrt(2)/2, sqrt(3)/2, sqrt((1 + sqrt(3)/2)/2), 1) if we are lazy, we can exploit R, and also we can do it if we want more points: b <- sin(x0) Then, if we decide that the upper polynomial (b0,b1...) will be of degree 4 and the lower one (a0,a1...) of degree 2 - let me make a short step aside - we have 6 known values, so we can estimate 6 coefficients. We know a0 and b0, so the sum of the polynomials degrees should be 6. I tried all combinations, 6-0, 5-1, etc., 4-2 seems to be the best - then the matrix can be constructed as: A <- cbind(x0, x0^2, x0^3, x0^4, -b*x0, -b*(x0^2)) and the coefficients can be calculated as: x1 <- solve(A, b) first four elements of x1 are the bi coefficients, the last two are the ai coefficients. You can create your own function mysin <- function(x){tst <- (x1[1]*x + x1[2]*x^2 + x1[3]*x^3 + x1[4]*x^4)/(1 + x1[5]*x + x1[6]*x^2) tst} and make a test like: mysin(pi*c(1:10)/20) - sin(pi*c(1:10)/20) this will show match to 6 decimal places for all values. ---- So, what I want to say, that if looking for sine approximation at interval (0, pi/2), which is of course correct idea, then for given number of coefficients, there is always better approximation that Taylor expansion with the same number of coefficients. Even the simple polynomial as of the situation 6-0 (so only bx coefficients) is better than Taylor expansion. It is even almost as good as the rational function 4-2. I have also made a quick research of how the sine is being numerically calculated. It is quite difficult to find something particular, but my impression is, that the most used method nowadays is a table of predefined values, and the rest is somehow interpolated, probably by a linear function. So that's all from my side :-) Pavel ```

JavaScript is required for this form.

No, thanks