math.cpp 4.5 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170
  1. /*
  2. * Copyright (c) 2020, the SerenityOS developers.
  3. *
  4. * SPDX-License-Identifier: BSD-2-Clause
  5. */
  6. // TODO: Is there a compiler flag we can use to still get these math functions? (and compile with -nostdlib)
  7. // FIXME: What do we do with this regarding copyright stuff?
  8. // code for math functions taken from here:
  9. // https://gitlab.incom.co/CM-Shield/u-boot/commit/aa7839b39c2ee77f9ab8c393c56b8d812507dbb7
  10. // https://github.com/zayac/qemu-arm/blob/master/qemu/roms/ipxe/src/libgcc/__udivmoddi4.c
  11. // https://code.woboq.org/llvm/compiler-rt/lib/builtins/divdi3.c.html
  12. #include "math.h"
  13. extern "C" {
  14. union overlay64 {
  15. u64 longw;
  16. struct {
  17. u32 lower;
  18. u32 higher;
  19. } words;
  20. };
  21. u64 __ashldi3(u64 num, unsigned int shift)
  22. {
  23. union overlay64 output;
  24. output.longw = num;
  25. if (shift >= 32) {
  26. output.words.higher = output.words.lower << (shift - 32);
  27. output.words.lower = 0;
  28. } else {
  29. if (!shift)
  30. return num;
  31. output.words.higher = (output.words.higher << shift) | (output.words.lower >> (32 - shift));
  32. output.words.lower = output.words.lower << shift;
  33. }
  34. return output.longw;
  35. }
  36. u64 __lshrdi3(u64 num, unsigned int shift)
  37. {
  38. union overlay64 output;
  39. output.longw = num;
  40. if (shift >= 32) {
  41. output.words.lower = output.words.higher >> (shift - 32);
  42. output.words.higher = 0;
  43. } else {
  44. if (!shift)
  45. return num;
  46. output.words.lower = output.words.lower >> shift | (output.words.higher << (32 - shift));
  47. output.words.higher = output.words.higher >> shift;
  48. }
  49. return output.longw;
  50. }
  51. #define MAX_32BIT_UINT ((((u64)1) << 32) - 1)
  52. static u64 _64bit_divide(u64 dividend, u64 divider, u64* rem_p)
  53. {
  54. u64 result = 0;
  55. /*
  56. * If divider is zero - let the rest of the system care about the
  57. * exception.
  58. */
  59. if (!divider)
  60. return 1 / (u32)divider;
  61. /* As an optimization, let's not use 64 bit division unless we must. */
  62. if (dividend <= MAX_32BIT_UINT) {
  63. if (divider > MAX_32BIT_UINT) {
  64. result = 0;
  65. if (rem_p)
  66. *rem_p = divider;
  67. } else {
  68. result = (u32)dividend / (u32)divider;
  69. if (rem_p)
  70. *rem_p = (u32)dividend % (u32)divider;
  71. }
  72. return result;
  73. }
  74. while (divider <= dividend) {
  75. u64 locald = divider;
  76. u64 limit = __lshrdi3(dividend, 1);
  77. int shifts = 0;
  78. while (locald <= limit) {
  79. shifts++;
  80. locald = locald + locald;
  81. }
  82. result |= __ashldi3(1, shifts);
  83. dividend -= locald;
  84. }
  85. if (rem_p)
  86. *rem_p = dividend;
  87. return result;
  88. }
  89. u64 __udivdi3(u64 num, u64 den)
  90. {
  91. return _64bit_divide(num, den, nullptr);
  92. }
  93. u64 __umoddi3(u64 num, u64 den)
  94. {
  95. u64 v = 0;
  96. _64bit_divide(num, den, &v);
  97. return v;
  98. }
  99. uint64_t __udivmoddi4(uint64_t num, uint64_t den, uint64_t* rem_p)
  100. {
  101. uint64_t quot = 0, qbit = 1;
  102. if (den == 0) {
  103. return 1 / ((unsigned)den); /* Intentional divide by zero, without
  104. triggering a compiler warning which
  105. would abort the build */
  106. }
  107. /* Left-justify denominator and count shift */
  108. while ((int64_t)den >= 0) {
  109. den <<= 1;
  110. qbit <<= 1;
  111. }
  112. while (qbit) {
  113. if (den <= num) {
  114. num -= den;
  115. quot += qbit;
  116. }
  117. den >>= 1;
  118. qbit >>= 1;
  119. }
  120. if (rem_p)
  121. *rem_p = num;
  122. return quot;
  123. }
  124. int64_t __divdi3(int64_t a, int64_t b)
  125. {
  126. const int bits_in_dword_m1 = (int)(sizeof(int64_t) * 8) - 1;
  127. int64_t s_a = a >> bits_in_dword_m1; // s_a = a < 0 ? -1 : 0
  128. int64_t s_b = b >> bits_in_dword_m1; // s_b = b < 0 ? -1 : 0
  129. a = (a ^ s_a) - s_a; // negate if s_a == -1
  130. b = (b ^ s_b) - s_b; // negate if s_b == -1
  131. s_a ^= s_b; // sign of quotient
  132. return (__udivmoddi4(a, b, (uint64_t*)0) ^ s_a) - s_a; // negate if s_a == -1
  133. }
  134. int64_t __moddi3(int64_t a, int64_t b)
  135. {
  136. const int bits_in_dword_m1 = (int)(sizeof(int64_t) * 8) - 1;
  137. int64_t s = b >> bits_in_dword_m1; // s = b < 0 ? -1 : 0
  138. b = (b ^ s) - s; // negate if s == -1
  139. s = a >> bits_in_dword_m1; // s = a < 0 ? -1 : 0
  140. a = (a ^ s) - s; // negate if s == -1
  141. uint64_t r = 0;
  142. __udivmoddi4(a, b, &r);
  143. return ((int64_t)r ^ s) - s; // negate if s == -1
  144. }
  145. }