organise

1 files changed, 402 insertions, 0 deletions

diff --git a/2/homework/2.92.c b/2/homework/2.92.c
new file mode 100644
index 0000000..bf1096e
--- /dev/null
+++ b/2/homework/2.92.c
@@ -0,0 +1,402 @@
+#include <assert.h>
+#include <stdio.h>
+#include <stdint.h>
+#include <limits.h>
+#include <math.h>
+#include "../code/data/show-bytes.c"
+
+typedef unsigned float_bits;
+
+#define M ~0
+
+unsigned i2u(int i) {
+    union {
+        unsigned u;
+        int i;
+    } ui;
+    ui.i = i;
+    return ui.u;
+}
+
+float u2f(unsigned int my_u) {
+    assert(sizeof(unsigned int) == sizeof(float));
+    union {
+        unsigned int u;
+        float f;
+    } my_union;
+    /* Apparently you cannot print the bytes of a float directly? like this printf("%x\n", f); */
+    my_union.u = my_u;
+    return my_union.f;
+}
+
+float fpwr2(int x)
+{
+    /* Result exponent and fraction */
+    unsigned exp, frac;
+    unsigned u;
+
+    int k = 8;
+    int n = 23;
+
+    int bias = (1<<(k-1))-1;
+    int max_e = ((1<<k) - 2) - bias;
+
+    if (x < 1-bias-n) {
+        /*  Too small. Return 0.0 */
+        exp = 0;
+        frac = 0;
+    } else if (x < 1-bias) {
+        /* Denormalized, [(1-bias-n), (1-bias)]  [-139, -126] */
+        exp = 0;
+        frac = 1<<(x+1-bias-n);
+    } else if (x<=max_e) {
+        /* Normalized */
+        exp = x+bias;
+        frac = 0;
+    } else {
+        /* Too big. Return +oo */
+        exp = (1<<8)-1;
+        frac = 0;
+    }
+
+    /* Pack exp and frac into 32 bits */
+    u = (exp << 23) | frac;
+    /* Return as float */
+    return u2f(u);
+}
+
+/*
+  int and unsigned and their operations are allowed.
+ */
+float_bits float_denorm_zero(float_bits f) {
+    unsigned sign = f>>31;
+    unsigned exp =  f>>23 & 0xFF;
+    unsigned frac = f     & 0x7FFFFF;
+    if (exp == 0) {
+        frac = 0;
+    }
+    return (sign << 31) | (exp << 23) | frac;
+}
+
+float_bits float_negate(float_bits f) {
+    unsigned is_valid_exp =  (f>>23) ^ 0xFF;
+    unsigned frac = f     & 0x7FFFFF;
+    if (!is_valid_exp && frac) return f;
+    return (1<<31) ^ f;
+}
+
+float_bits float_absval(float_bits f) {
+    unsigned is_valid_exp =  (f>>23) ^ 0xFF;
+    unsigned frac = f     & 0x7FFFFF;
+    if (!is_valid_exp && frac) return f;
+    return ~(1<<31) & f;
+}
+
+float_bits float_twice(float_bits f) {
+    /*
+     */
+    unsigned sign = f>>31;
+    unsigned exp = (f>>23) & 0xFF;
+    unsigned frac = f     & 0x7FFFFF;
+    unsigned is_valid_exp =  exp ^ 0xFF;
+    if (!is_valid_exp) return f;
+
+    if (!exp) {
+        /*
+          Denormalized because we have a zero exponent
+          f<<1
+         */
+        return (sign<<31) | (frac<<1);
+    } else {
+        /*
+          Normalized value.
+          = (2^k(1+f))*2
+          = 2^(k+1)(1+f)
+         */
+        unsigned twice_exp = exp + 1;
+        if (twice_exp == 0xFF) {
+            frac = 0;
+        }
+        return (sign<<31) | (twice_exp<<23) | frac;
+    }
+}
+
+float_bits float_half(float_bits f) {
+    unsigned sign = f>>31;
+    unsigned exp = (f>>23) & 0xFF;
+    unsigned frac = f     & 0x7FFFFF;
+    unsigned is_valid_exp =  exp ^ 0xFF;
+    if (!is_valid_exp) return f;
+
+    if (!exp) {
+        /*
+          Denormalized because we have a zero exponent
+          f>>1
+         */
+        if (frac&1 && frac^1) {
+            /*
+              Round to nearest even?
+            */
+            frac >>= 1;
+            frac += frac&1;
+        } else {
+            frac >>= 1;
+        }
+        return (sign<<31) | frac;
+    } else {
+        /*
+          Normalized value.
+         */
+        unsigned half_exp = exp;
+        /* Can we halve with exp? */
+        if ((half_exp - 1) == 0) {
+            /* Is it the value in the middle between norm and denorm? */
+            if ((frac ^ 0x7FFFFF) == 0) {
+                frac = 0;
+            } else {
+                /* Otherwise, normally divide the (1+f) by two */
+                half_exp -= 1;
+                if (frac&1 && frac^1) {
+                    frac >>= 1;
+                    frac += frac&1;
+                } else {
+                    frac >>= 1;
+                }
+                frac |= 1<<22;
+            }
+        } else {
+            half_exp -= 1;
+        }
+        return (sign<<31) | (half_exp<<23) | frac;
+    }
+}
+
+/*
+  Compute (int) f.
+  If conversion causes overflow or f is NaN, return 0x80 00 00 00
+ */
+int float_f2i(float_bits f) {
+    int invalid = 0x80000000;
+    unsigned sign = f>>31;
+    unsigned exp = (f>>23) & 0xFF;
+    unsigned frac = f     & 0x7FFFFF;
+    unsigned is_valid_exp =  exp ^ 0xFF;
+    if (!is_valid_exp) return invalid;
+
+    unsigned bias = (1<<7) - 1;
+    if (exp < bias) return 0;
+    exp -= bias;
+    /*
+      Smallest odd that cannot be represented exactly
+      2^(n+1) + 1
+    */
+    if (exp < 24) {
+        int result = (1<<exp) + (frac>>(23 - exp));
+        if (sign) {
+            result = ~result + 1;
+        }
+        return result;
+    }
+    if (exp < 31) {
+        int result = (1<<exp) + (frac<<(exp - 23));
+        if (sign) {
+            result = ~result + 1;
+        }
+        return result;
+    }
+    return invalid;
+}
+
+float_bits float_i2f(int i) {
+    if (i==0) return 0;
+    /*
+      Need to calculate exp and frac based on integer, so we will take apart the integer into exp and frac using the fundamenal property of division.
+      p = q + r, where q and r are integers and q is the result of integer division of p
+     */
+    unsigned sign = i>>31;
+    int exp, frac, k, n, bias;
+
+    k = 8;
+    n = 23;
+    bias = (1<<(k-1))-1;
+    exp = bias;
+    frac = 0;
+
+    if (i==INT_MIN) {
+        /* Special case since ~i+1=i, if i=INT_MIN */
+        return (sign<<31) | ((exp+31)<<23) | 0;
+    }
+
+    /* Get the greatest power of 2, p2 and remainder */
+    int p2, rem;
+    int u = (i & ~sign) | ((~i + 1) & sign);
+    int u_orig = u;
+    int b16, b8, b4, b2, b1, b0;
+
+    b16 = (!!(u>>16))<<4;
+    u >>= b16;
+    b8 = (!!(u>>8))<<3;
+    u >>= b8;
+    b4 = (!!(u>>4))<<2;
+    u >>= b4;
+    b2 = (!!(u>>2))<<1;
+    u >>= b2;
+    b1 = (!!(u>>1));
+
+    p2 = b16 + b8 + b4 + b2 + b1;
+    exp += p2;
+    rem = u_orig - (1<<p2);
+
+    /* How much greater than precision are we */
+    int prec = (p2-23);
+    int prec2 = (1<<prec);
+
+    if (rem == 0) {
+        frac = 0;
+    } else {
+        if (p2 > 23) {
+            /* Outside of float32 precision */
+
+            /* Special case, Is remainder greater than half precision removed from highest power? */
+            if ((1<<p2) - (prec2>>1) <= rem) {
+                /* Round up into next power level */
+                frac = 0;
+                exp += 1;
+            } else {
+                /* Encode remainder as fraction with prec denom */
+                frac = rem>>prec;
+                if (prec>1) {
+                    /* Rounding is non binary, need to check if round to even applies or round up */
+                    int r = rem%prec2;
+                    if (r > (prec2>>1)) {
+                        frac += 1;
+                    } else if (r==(prec2>>1)) {
+                        frac += frac&1;
+                    }
+                } else {
+                    /* Binary rounding is easy, even rem don't need it and odd rem need it if the result is odd */
+                    frac += rem&1 && frac&1;
+                }
+            }
+        } else {
+            /* Within precision, just encode remainder as fraction with prec denom */
+            frac = (rem<<(23-p2));
+        }
+    }
+    return (sign<<31) | (exp<<23) | frac;
+    /*
+      1 1001 1101 [0]+
+      1 1001 1100 [1]+
+     */
+    /*
+      Denormalized
+       V = M * 2^E
+       M = f < 1-e
+       E = 1 - bias
+       <=> V < 1
+
+       Normalized
+       M = 1 + f
+       V = 1 = 2^E *(1+f)    <==>    e>=bias, e=bias <=> V=(2^0)(1+f)
+       E = e-bias
+       bias=127
+
+       smallest int
+       -(2^31)
+       biggest int
+       (2^31)-1
+       1 <> 0 0111 1111 [0]+
+       2 <> 0 1000 0000 [0]+
+       3 <> 0 1000 0000 1[0]+
+       4 <> 0 1000 0001 [0]+
+       5 <> 0 1000 0001 01[0]+
+       6 <> 0 1000 0001 1[0]+
+       7 <> 0 1000 0001 11[0]+
+       8 <> 0 1000 0010 000[0]+
+       9 <> 0 1000 0010 001[0]+
+       10 <> 0 1000 0010 010[0]+
+       11 <> 0 1000 0010 011[0]+
+     */
+    return 1;
+}
+
+int main(void) {
+    /* 2.97 */
+    int j;
+    for (j=INT_MIN; j<INT_MAX; j++) {
+        float v = u2f(float_i2f(j));
+        int result = (((float) j) == v);
+            /* show_float((float) j); */
+            /* show_float(v); */
+        if (!result) {
+            printf("want: %f, got: %f\n", ((float) j), v);
+            show_float((float) j);
+            show_float(v);
+        }
+        assert(result);
+    }
+    return 0;
+
+    /* 2.96 */
+    unsigned i;
+    for (i=0; i<M; i++) {
+        int v = float_f2i(i);
+        int result = isnan(u2f(i)) || (((int) u2f(i)) == v);
+        if (!result) {
+            printf("want: %d, got: %d\n", ((int) u2f(i)), v);
+            show_float(u2f(i));
+        }
+        assert(result);
+    }
+    return 0;
+
+    /* 2.95 */
+    for (i=0; i<M; i++) {
+        float v = u2f(float_half(i));
+        int result = isnan(u2f(i)) || (0.5f*u2f(i) == v);
+        if (!result) {
+            printf("want: %f, got: %f\n", 0.5f*u2f(i), v);
+            show_float(u2f(1));
+            show_float(u2f(i));
+            show_float(0.5f*u2f(i));
+            show_float(v);
+        }
+        assert(result);
+    }
+    return 0;
+
+    /* 2.94 */
+    for (i=0; i<M; i++) {
+        float v = u2f(float_twice(i));
+        int result = isnan(u2f(i)) || (2*u2f(i) == v);
+        if (!result) {
+            printf("want: %f, got: %f\n", 2*u2f(i), v);
+            show_float(2*u2f(i));
+            show_float(v);
+        }
+        assert(result);
+    }
+    return 0;
+
+    /* 2.93 */
+    for (i=0; i<M; i++) {
+        float v = u2f(float_absval(i));
+        int result = isnan(v) || (fabsf(u2f(i)) == v);
+        if (!result) {
+            printf("%f\n", fabsf(u2f(i)));
+        }
+        assert(result);
+    }
+    return 0;
+
+    /* 2.92 */
+    for (i=0; i<M; i++) {
+        float v = u2f(float_negate(i));
+        int result = isnan(v) || (-u2f(i) == v);
+        if (!result) {
+            printf("%f\n", u2f(i));
+        }
+        assert(result);
+    }
+    return 0;
+}