From d5216879d463f1f65d83bb675285e3694f52a108 Mon Sep 17 00:00:00 2001
From: whitequark <whitequark@whitequark.org>
Date: Sun, 20 Dec 2015 04:06:03 +0800
Subject: [PATCH] runtime: implement round() (fixes #198).

From https://github.com/bpowers/musl/blob/master/src/math/round.c.
---
 artiq/runtime/ksupport.c | 33 +++++++++++++++++++++++++++++++++
 1 file changed, 33 insertions(+)

diff --git a/artiq/runtime/ksupport.c b/artiq/runtime/ksupport.c
index b1bcc984f..2f775cf4b 100644
--- a/artiq/runtime/ksupport.c
+++ b/artiq/runtime/ksupport.c
@@ -17,6 +17,8 @@
 #include "dds.h"
 #include "rtio.h"
 
+double round(double x);
+
 void ksupport_abort(void);
 
 int64_t now;
@@ -82,6 +84,9 @@ static const struct symbol runtime_exports[] = {
     {"__moddi3", &__moddi3},
     {"__powidf2", &__powidf2},
 
+    /* libm */
+    {"round", &round},
+
     /* exceptions */
     {"_Unwind_Resume", &_Unwind_Resume},
     {"__artiq_personality", &__artiq_personality},
@@ -118,6 +123,34 @@ static const struct symbol runtime_exports[] = {
     {NULL, NULL}
 };
 
+double round(double x)
+{
+    union {double f; uint64_t i;} u = {x};
+    int e = u.i >> 52 & 0x7ff;
+    double y;
+
+    if (e >= 0x3ff+52)
+        return x;
+    if (u.i >> 63)
+        x = -x;
+    if (e < 0x3ff-1) {
+        /* we don't do it in ARTIQ */
+        /* raise inexact if x!=0 */
+        // FORCE_EVAL(x + 0x1p52);
+        return 0*u.f;
+    }
+    y = (double)(x + 0x1p52) - 0x1p52 - x;
+    if (y > 0.5)
+        y = y + x - 1;
+    else if (y <= -0.5)
+        y = y + x + 1;
+    else
+        y = y + x;
+    if (u.i >> 63)
+        y = -y;
+    return y;
+}
+
 /* called by libunwind */
 int fprintf(FILE *stream, const char *fmt, ...)
 {