aboutsummaryrefslogtreecommitdiff
path: root/src/base/math/sqrt.c
diff options
context:
space:
mode:
authorNicholas Noll <nbnoll@eml.cc>2021-12-05 16:16:21 -0800
committerNicholas Noll <nbnoll@eml.cc>2021-12-05 16:16:21 -0800
commit07e77936d535e58b0aeb4f2a11400c1050556739 (patch)
treefb50fad6436ecbc159505bee73a92cd289b99a07 /src/base/math/sqrt.c
parentb48327d357e0818d1a6ae2a064cfa7d1567e1242 (diff)
Feat: added math library
Used Plan9's libc as starting point. This cleans up dangling references due to loss of libc.
Diffstat (limited to 'src/base/math/sqrt.c')
-rw-r--r--src/base/math/sqrt.c54
1 files changed, 54 insertions, 0 deletions
diff --git a/src/base/math/sqrt.c b/src/base/math/sqrt.c
new file mode 100644
index 0000000..a94d69f
--- /dev/null
+++ b/src/base/math/sqrt.c
@@ -0,0 +1,54 @@
+/*
+ sqrt returns the square root of its floating
+ point argument. Newton's method.
+
+ calls frexp
+*/
+
+#include <u.h>
+#include <base.h>
+
+double
+math·sqrt(double arg)
+{
+ double x, temp;
+ int exp, i;
+
+ if(arg <= 0) {
+ if(arg < 0)
+ return math·NaN();
+ return 0;
+ }
+ if(math·isInf(arg, 1))
+ return arg;
+ x = math·frexp(arg, &exp);
+ while(x < 0.5) {
+ x *= 2;
+ exp--;
+ }
+ /*
+ * NOTE
+ * this wont work on 1's comp
+ */
+ if(exp & 1) {
+ x *= 2;
+ exp--;
+ }
+ temp = 0.5 * (1.0+x);
+
+ while(exp > 60) {
+ temp *= (1L<<30);
+ exp -= 60;
+ }
+ while(exp < -60) {
+ temp /= (1L<<30);
+ exp += 60;
+ }
+ if(exp >= 0)
+ temp *= 1L << (exp/2);
+ else
+ temp /= 1L << (-exp/2);
+ for(i=0; i<=4; i++)
+ temp = 0.5*(temp + arg/temp);
+ return temp;
+}