diff options
author | jack <[email protected]> | 2024-04-22 10:06:09 -0600 |
---|---|---|
committer | jack <[email protected]> | 2024-04-22 10:06:09 -0600 |
commit | 0ceb7b89f95b52e07f4fcbed0622f1e9758560b3 (patch) | |
tree | 9baf56dc27154046387a03fd0ff60207f500490b | |
parent | 90971c7e02768082b092020cfced5f05566f594a (diff) | |
download | qmk_me-0ceb7b89f95b52e07f4fcbed0622f1e9758560b3.tar.gz qmk_me-0ceb7b89f95b52e07f4fcbed0622f1e9758560b3.zip |
add keymaps to link script, continue mathlib
-rw-r--r-- | keymaps/onekey/stm32f411/config.h | 9 | ||||
-rw-r--r-- | keymaps/onekey/stm32f411/qfplib/functions.c | 1026 | ||||
-rw-r--r-- | keymaps/onekey/stm32f411/qfplib/qfplib.s | 825 | ||||
-rw-r--r-- | keymaps/onekey/stm32f411/rules.mk | 3 | ||||
-rwxr-xr-x | util/link.sh | 16 |
5 files changed, 426 insertions, 1453 deletions
diff --git a/keymaps/onekey/stm32f411/config.h b/keymaps/onekey/stm32f411/config.h new file mode 100644 index 0000000..84ac83f --- /dev/null +++ b/keymaps/onekey/stm32f411/config.h @@ -0,0 +1,9 @@ +// Copyright 2024 jack (@waffle87) +// SPDX-License-Identifier: GPL-2.0-or-later +#pragma once + +#undef ADC_PIN +#undef BACKLIGHT_PIN +#undef MATRIX_ROW_PINS +#undef MATRIX_COL_PINS +#define DIRECT_PINS {{A0}} diff --git a/keymaps/onekey/stm32f411/qfplib/functions.c b/keymaps/onekey/stm32f411/qfplib/functions.c deleted file mode 100644 index 7a2b051..0000000 --- a/keymaps/onekey/stm32f411/qfplib/functions.c +++ /dev/null @@ -1,1026 +0,0 @@ -// Copyright 2018 @lupyuen -// SPDX-License-Identifier: GPL-2.0-or-later -#include "qfplib.h" -#include <math.h> -#include <stdint.h> -#define M_LN10 2.30258509299404568402 // Natural log of 10 -#define M_PI_2 1.57079632679489661923 // Pi divided by 2 - -// We count the number of times each function was called. So we can check -// whether our unit tests cover all functions. -enum float_usage_index { - FIRST_USAGE_INDEX, - USAGE_AEABI_DDIV, - USAGE_AEABI_DMUL, - USAGE_AEABI_DADD, - USAGE_AEABI_DSUB, - USAGE_AEABI_DCMPEQ, - USAGE_AEABI_DCMPLT, - USAGE_AEABI_DCMPLE, - USAGE_AEABI_DCMPGE, - USAGE_AEABI_DCMPGT, - USAGE_AEABI_DCMPUN, - USAGE_AEABI_FDIV, - USAGE_AEABI_FADD, - USAGE_AEABI_FSUB, - USAGE_AEABI_FMUL, - USAGE_AEABI_D2IZ, - USAGE_AEABI_D2UIZ, - USAGE_ADDSF3, - USAGE_SUBSF3, - USAGE_MULSF3, - USAGE_DIVSF3, - USAGE_SQRT, - USAGE_LOG, - USAGE_EXP, - USAGE_LOG2, - USAGE_LOG10, - USAGE_POW, - USAGE_LDEXP, - USAGE_SIN, - USAGE_COS, - USAGE_TAN, - USAGE_ATAN2, - USAGE_ATAN, - USAGE_ASIN, - USAGE_ACOS, - USAGE_SINH, - USAGE_COSH, - USAGE_TANH, - USAGE_ASINH, - USAGE_ACOSH, - USAGE_ATANH, - USAGE_TRUNC, - USAGE_FLOOR, - USAGE_CEIL, - USAGE_FMOD, - USAGE_FABS, - LAST_FLOAT_USAGE_INDEX -}; -static uint8_t float_usage[LAST_FLOAT_USAGE_INDEX]; -uint8_t *get_float_usage(uint16_t *size) { - *size = (uint16_t)LAST_FLOAT_USAGE_INDEX - 1; - return float_usage; -} - -// Run-time ABI for the ARM Architecture. The function names are wrapped via -// "-Wl,-wrap,__aeabi..." in newlib/CMakeLists.txt. See -// http://infocenter.arm.com/help/topic/com.arm.doc.ihi0043d/IHI0043D_rtabi.pdf - -/////////////////////////////////////////////////////////////////////////////// -// Table 2, Standard double precision floating-point arithmetic helper -// functions - -// double-precision division, n / d -double __wrap___aeabi_ddiv(double n, double d) { - float_usage[USAGE_AEABI_DDIV]++; - return qfp_fdiv_fast(n, d); -} -// Unit Tests: -// aeabi_ddiv(2205.1969, 270.8886) = 8.140604292687105 -// aeabi_ddiv(-2205.1969, 270.8886) = -8.140604292687105 -// aeabi_ddiv(2205.1969, -270.8886) = -8.140604292687105 -// aeabi_ddiv(-2205.1969, -270.8886) = 8.140604292687105 - -// double-precision multiplication -double __wrap___aeabi_dmul(double x, double y) { - float_usage[USAGE_AEABI_DMUL]++; - return qfp_fmul(x, y); -} -// Unit Tests: -// aeabi_dmul(2205.1969, 270.8886) = 597362.70096534 -// aeabi_dmul(-2205.1969, 270.8886) = -597362.70096534 -// aeabi_dmul(2205.1969, -270.8886) = -597362.70096534 -// aeabi_dmul(-2205.1969, -270.8886) = 597362.70096534 - -double __wrap___aeabi_dadd(double a, double b) { - float_usage[USAGE_AEABI_DADD]++; - return qfp_fadd(a, b); -} -// Unit Tests: -// aeabi_dadd(2205.1969, 270.8886) = 2476.0855 -// aeabi_dadd(-2205.1969, 270.8886) = -1934.3083 -// aeabi_dadd(2205.1969, -270.8886) = 1934.3083 -// aeabi_dadd(-2205.1969, -270.8886) = -2476.0855 - -double __wrap___aeabi_dsub(double a, double b) { - float_usage[USAGE_AEABI_DSUB]++; - return qfp_fsub(a, b); -} -// Unit Tests: -// aeabi_dsub(2205.1969, 270.8886) = 1934.3083 -// aeabi_dsub(-2205.1969, 270.8886) = -2476.0855 -// aeabi_dsub(2205.1969, -270.8886) = 2476.0855 -// aeabi_dsub(-2205.1969, -270.8886) = -1934.3083 - -/////////////////////////////////////////////////////////////////////////////// -// Table 3, double precision floating-point comparison helper functions - -// qfp_fcmp(r0, r1): -// equal? return 0 -// r0 > r1? return +1 -// r0 < r1: return -1 - -// result (1, 0) denotes (=, ?<>) [2], use for C == and != -int __wrap___aeabi_dcmpeq(double x, double y) { - float_usage[USAGE_AEABI_DCMPEQ]++; - return (qfp_fcmp(x, y) == 0) // x == y - ? 1 - : 0; -} -// Unit Tests: -// aeabi_dcmpeq(2205.196, 2205.196) = 1 -// aeabi_dcmpeq(2205.196, 2205.195) = 0 -// aeabi_dcmpeq(2205.196, 2205.197) = 0 -// aeabi_dcmpeq(2205.196, 0) = 0 -// aeabi_dcmpeq(-2205.196, -2205.196) = 1 -// aeabi_dcmpeq(-2205.196, -2205.195) = 0 -// aeabi_dcmpeq(-2205.196, -2205.197) = 0 -// aeabi_dcmpeq(-2205.196, 0) = 0 - -// result (1, 0) denotes (<, ?>=) [2], use for C < -int __wrap___aeabi_dcmplt(double x, double y) { - float_usage[USAGE_AEABI_DCMPLT]++; - return (qfp_fcmp(x, y) < 0) // x < y - ? 1 - : 0; -} -// Unit Tests: -// aeabi_dcmplt(2205.196, 2205.196) = 0 -// aeabi_dcmplt(2205.196, 2205.195) = 0 -// aeabi_dcmplt(2205.196, 2205.197) = 1 -// aeabi_dcmplt(2205.196, 0) = 0 -// aeabi_dcmplt(-2205.196, -2205.196) = 0 -// aeabi_dcmplt(-2205.196, -2205.195) = 1 -// aeabi_dcmplt(-2205.196, -2205.197) = 0 -// aeabi_dcmplt(-2205.196, 0) = 1 - -// result (1, 0) denotes (<=, ?>) [2], use for C <= -int __wrap___aeabi_dcmple(double x, double y) { - float_usage[USAGE_AEABI_DCMPLE]++; - return (qfp_fcmp(x, y) > 0) // x > y - ? 0 - : 1; -} -// Unit Tests: -// aeabi_dcmple(2205.196, 2205.196) = 1 -// aeabi_dcmple(2205.196, 2205.195) = 0 -// aeabi_dcmple(2205.196, 2205.197) = 1 -// aeabi_dcmple(2205.196, 0) = 0 -// aeabi_dcmple(-2205.196, -2205.196) = 1 -// aeabi_dcmple(-2205.196, -2205.195) = 1 -// aeabi_dcmple(-2205.196, -2205.197) = 0 -// aeabi_dcmple(-2205.196, 0) = 1 - -// result (1, 0) denotes (>=, ?<) [2], use for C >= -int __wrap___aeabi_dcmpge(double x, double y) { - float_usage[USAGE_AEABI_DCMPGE]++; - return (qfp_fcmp(x, y) < 0) // x < y - ? 0 - : 1; -} -// Unit Tests: -// aeabi_dcmpge(2205.196, 2205.196) = 1 -// aeabi_dcmpge(2205.196, 2205.195) = 1 -// aeabi_dcmpge(2205.196, 2205.197) = 0 -// aeabi_dcmpge(2205.196, 0) = 1 -// aeabi_dcmpge(-2205.196, -2205.196) = 1 -// aeabi_dcmpge(-2205.196, -2205.195) = 0 -// aeabi_dcmpge(-2205.196, -2205.197) = 1 -// aeabi_dcmpge(-2205.196, 0) = 0 - -// result (1, 0) denotes (>, ?<=) [2], use for C > -int __wrap___aeabi_dcmpgt(double x, double y) { - float_usage[USAGE_AEABI_DCMPGT]++; - return (qfp_fcmp(x, y) > 0) // x > y - ? 1 - : 0; -} -// Unit Tests: -// aeabi_dcmpgt(2205.196, 2205.196) = 0 -// aeabi_dcmpgt(2205.196, 2205.195) = 1 -// aeabi_dcmpgt(2205.196, 2205.197) = 0 -// aeabi_dcmpgt(2205.196, 0) = 1 -// aeabi_dcmpgt(-2205.196, -2205.196) = 0 -// aeabi_dcmpgt(-2205.196, -2205.195) = 0 -// aeabi_dcmpgt(-2205.196, -2205.197) = 1 -// aeabi_dcmpgt(-2205.196, 0) = 0 - -// result (1, 0) denotes (?, <=>) [2], use for C99 isunordered() -int __wrap___aeabi_dcmpun(double x, double y) { - float_usage[USAGE_AEABI_DCMPUN]++; - return (qfp_fcmp(x, y) == 0) // x == y - ? 0 - : 1; -} -// Unit Tests: -// aeabi_dcmpun(2205.196, 2205.196) = 0 -// aeabi_dcmpun(2205.196, 2205.195) = 1 -// aeabi_dcmpun(2205.196, 2205.197) = 1 -// aeabi_dcmpun(2205.196, 0) = 1 -// aeabi_dcmpun(-2205.196, -2205.196) = 0 -// aeabi_dcmpun(-2205.196, -2205.195) = 1 -// aeabi_dcmpun(-2205.196, -2205.197) = 1 -// aeabi_dcmpun(-2205.196, 0) = 1 - -/////////////////////////////////////////////////////////////////////////////// -// Table 4, Standard single precision floating-point arithmetic helper -// functions - -// single-precision division, n / d -float __wrap___aeabi_fdiv(float n, float d) { - float_usage[USAGE_AEABI_FDIV]++; - return qfp_fdiv_fast(n, d); -} -// Unit Tests: -// aeabi_fdiv(2205.1969, 270.8886) = 8.140604292687105 -// aeabi_fdiv(-2205.1969, 270.8886) = -8.140604292687105 -// aeabi_fdiv(2205.1969, -270.8886) = -8.140604292687105 -// aeabi_fdiv(-2205.1969, -270.8886) = 8.140604292687105 - -float __wrap___aeabi_fadd(float a, float b) { - float_usage[USAGE_AEABI_FADD]++; - return qfp_fadd(a, b); -} -// Unit Tests: -// aeabi_fadd(2205.1969, 270.8886) = 2476.0855 -// aeabi_fadd(-2205.1969, 270.8886) = -1934.3083 -// aeabi_fadd(2205.1969, -270.8886) = 1934.3083 -// aeabi_fadd(-2205.1969, -270.8886) = -2476.0855 - -float __wrap___aeabi_fsub(float a, float b) { - float_usage[USAGE_AEABI_FSUB]++; - return qfp_fsub(a, b); -} -// Unit Tests: -// aeabi_fsub(2205.1969, 270.8886) = 1934.3083 -// aeabi_fsub(-2205.1969, 270.8886) = -2476.0855 -// aeabi_fsub(2205.1969, -270.8886) = 2476.0855 -// aeabi_fsub(-2205.1969, -270.8886) = -1934.3083 - -float __wrap___aeabi_fmul(float a, float b) { - float_usage[USAGE_AEABI_FMUL]++; - return qfp_fmul(a, b); -} -// Unit Tests: -// aeabi_fmul(2205.1969, 270.8886) = 597362.70096534 -// aeabi_fmul(-2205.1969, 270.8886) = -597362.70096534 -// aeabi_fmul(2205.1969, -270.8886) = -597362.70096534 -// aeabi_fmul(-2205.1969, -270.8886) = 597362.70096534 - -/////////////////////////////////////////////////////////////////////////////// -// Table 6, Standard floating-point to integer conversions - -// double to integer C-style conversion. "z" means round towards 0. -int __wrap___aeabi_d2iz(double x) { - float_usage[USAGE_AEABI_D2IZ]++; - if (qfp_fcmp(x, 0) == 0) { - return 0; - } - // qfp_float2int() works like floor(). If x is negative, we add 1 to the - // result. - int xfloored = qfp_float2int(x); - if (xfloored < 0) { - return xfloored + 1; - } - return xfloored; -} -// Unit Tests: -// aeabi_d2iz(0) = 0 -// aeabi_d2iz(2205.1969) = 2205 -// aeabi_d2iz(-2205.1969) = -2205 - -// double to unsigned C-style conversion. "z" means round towards 0. -unsigned __wrap___aeabi_d2uiz(double x) { - float_usage[USAGE_AEABI_D2UIZ]++; - if (qfp_fcmp(x, 0) == 0) { - return 0; - } - if (qfp_fcmp(x, 0) < 0) { - return 0; - } - return qfp_float2uint(x); -} -// Unit Tests: -// aeabi_d2iz(0) = 0 -// aeabi_d2uiz(2205.1969) = 2205 -// aeabi_d2uiz(-2205.1969) = 0 - -/////////////////////////////////////////////////////////////////////////////// -// GNU C Library Routines for floating point emulation -// From https://gcc.gnu.org/onlinedocs/gccint/Soft-float-library-routines.html - -float __wrap___addsf3(float a, float b) { - float_usage[USAGE_ADDSF3]++; - return qfp_fadd(a, b); -} -// Unit Tests: -// addsf3(2205.1969, 270.8886) = 2476.0855 -// addsf3(-2205.1969, 270.8886) = -1934.3083 -// addsf3(2205.1969, -270.8886) = 1934.3083 -// addsf3(-2205.1969, -270.8886) = -2476.0855 - -float __wrap___subsf3(float a, float b) { - float_usage[USAGE_SUBSF3]++; - return qfp_fsub(a, b); -} -// Unit Tests: -// subsf3(2205.1969, 270.8886) = 1934.3083 -// subsf3(-2205.1969, 270.8886) = -2476.0855 -// subsf3(2205.1969, -270.8886) = 2476.0855 -// subsf3(-2205.1969, -270.8886) = -1934.3083 - -float __wrap___mulsf3(float a, float b) { - float_usage[USAGE_MULSF3]++; - return qfp_fmul(a, b); -} -// Unit Tests: -// mulsf3(2205.1969, 270.8886) = 597362.70096534 -// mulsf3(-2205.1969, 270.8886) = -597362.70096534 -// mulsf3(2205.1969, -270.8886) = -597362.70096534 -// mulsf3(-2205.1969, -270.8886) = 597362.70096534 - -float __wrap___divsf3(float a, float b) { - float_usage[USAGE_DIVSF3]++; - return qfp_fdiv_fast(a, b); -} -// Unit Tests: -// divsf3(2205.1969, 270.8886) = 8.140604292687105 -// divsf3(-2205.1969, 270.8886) = -8.140604292687105 -// divsf3(2205.1969, -270.8886) = -8.140604292687105 -// divsf3(-2205.1969, -270.8886) = 8.140604292687105 - -/////////////////////////////////////////////////////////////////////////////// -// <math.h> Functions - -// CMakeFiles/STM32_BLUE_PILL.dir/pxtapp/base/core.cpp.o: In function -// `Math_::sqrt(pxt::TValueStruct*)': /src/pxtapp/base/core.cpp:925: undefined -// reference to `sqrt' - -double sqrt(double x) { - float_usage[USAGE_SQRT]++; - return qfp_fsqrt_fast(x); -} -// Unit Tests: -// sqrt(100) = 10.000000 -// sqrt(2) = 1.414214 -// sqrt(-0) = -0.000000 -//// sqrt(-1.0) = -nan - -// CMakeFiles/STM32_BLUE_PILL.dir/pxtapp/base/core.cpp.o: In function -// `Math_::log(pxt::TValueStruct*)': /src/pxtapp/base/core.cpp:901: undefined -// reference to `log' - -double log(double x) { - float_usage[USAGE_LOG]++; - return qfp_fln(x); -} -// Unit Tests: -// log(1) = 0.000000 -// log(2) = _M_LN2 -// log(10) = M_LN10 -//// log(+Inf) = inf -//// log(0) = -inf - -double exp(double x) { - float_usage[USAGE_EXP]++; - return qfp_fexp(x); -} -// Unit Tests: -// exp(1) = 2.718282 -// exp(_M_LN2) = 2 -// exp(M_LN10) = 10 -// exp(-0) = 1.000000 -//// exp(-Inf) = 0.000000 - -// log2(x) = ln(x) / ln(2) -#undef log2 -double log2(double x) { - float_usage[USAGE_LOG2]++; - return qfp_fmul(qfp_fln(x), - 1.0 / _M_LN2 // Constant - ); -} -// Unit Tests: -// log2(65536) = 16.000000 -// log2(0.125) = -3.000000 -// log2(527) = 9.041659 -// log2(1) = 0.000000 -//// log2(+Inf) = inf -//// log2(0) = -inf - -// CMakeFiles/STM32_BLUE_PILL.dir/pxtapp/base/core.cpp.o: In function -// `Math_::log10(pxt::TValueStruct*)': /src/pxtapp/base/core.cpp:904: undefined -// reference to `log10' - -// log10(x) = ln(x) / ln(10) -// e.g. log10(1000) = ln(1000) / ln(10) = 3 -double log10(double x) { - float_usage[USAGE_LOG10]++; - return qfp_fmul(qfp_fln(x), - 1.0 / M_LN10 // Constant - ); -} -// Unit Tests: -// log10(1000) = 3.000000 -// log10(0.001) = -3.000000 -// base-5 logarithm of 125 = 3.000000 -// log10(1) = 0.000000 -//// log10(+Inf) = inf -//// log10(0) = -inf - -// pow(b, x) = pow(e, log(b) * x) = exp(log(b) * x) -// e.g. pow(10, 3) = exp(log(10) * 3) = 1000 -double pow(double b, double x) { - float_usage[USAGE_POW]++; - // If b and x are both negative and x is odd, then return -pow(-b, x). - // e.g. pow(-2, -3) = -pow(2, -3) - if (qfp_fcmp(b, 0) < 0 && qfp_fcmp(x, 0) < 0) { - int xfloored = qfp_float2int(-x); - if (xfloored % 2 == 1) { - return -qfp_fexp(qfp_fmul(qfp_fln(-b), x)); - } - } - return qfp_fexp(qfp_fmul(qfp_fln(b), x)); -} -// Unit Tests: -// pow(2, 10) = 1024.000000 -// pow(2, 0.5) = 1.414214 -// pow(-2, -3) = -0.125000 -// pow(-1, NAN) = nan -// pow(+1, NAN) = 1.000000 -//// pow(INFINITY, 2) = inf -//// pow(INFINITY, -1) = 0.000000 -//// pow(-1, 1/3) = -nan - -// ldexp(x, ex) = x * pow(2, ex) -// = x * exp(log(2) * ex) -double ldexp(double x, int ex) { - float_usage[USAGE_LDEXP]++; - return qfp_fmul(x, qfp_fexp(qfp_fmul(_M_LN2, ex))); -} -// Unit Tests: -// ldexp(7, -4) = 0.437500 -// ldexp(-0, 10) = -0.000000 -//// ldexp(-Inf, -1) = -inf -//// ldexp(1, 1024) = inf - -// CMakeFiles/STM32_BLUE_PILL.dir/pxtapp/base/core.cpp.o: In function -// `Math_::sin(pxt::TValueStruct*)': /src/pxtapp/base/core.cpp:910: undefined -// reference to `sin' - -double sin(double x) { - float_usage[USAGE_SIN]++; - return qfp_fsin(x); -} -// Unit Tests: -// sin(pi/6) = 0.500000 -// sin(pi/2) = 1.000000 -// sin(-3*pi/4) = -0.707107 -// sin(+0) = 0.000000 -// sin(-0) = -0.000000 -//// sin(INFINITY) = -nan - -// CMakeFiles/STM32_BLUE_PILL.dir/pxtapp/base/core.cpp.o: In function -// `Math_::cos(pxt::TValueStruct*)': /src/pxtapp/base/core.cpp:913: undefined -// reference to `cos' - -double cos(double x) { - float_usage[USAGE_COS]++; - // If x is Pi/2, return 0. - if (qfp_fcmp(x, M_PI_2) == 0) { - return 0; - } - return qfp_fcos(x); -} -// Unit Tests: -// cos(pi/3) = 0.500000 -// cos(pi/2) = 0.000000 -// cos(-3*pi/4) = -0.707107 -// cos(+0) = 1.000000 -// cos(-0) = 1.000000 -//// cos(INFINITY) = -nan - -// CMakeFiles/STM32_BLUE_PILL.dir/pxtapp/base/core.cpp.o: In function -// `Math_::tan(pxt::TValueStruct*)': /src/pxtapp/base/core.cpp:907: undefined -// reference to `tan' - -double tan(double x) { - float_usage[USAGE_TAN]++; - return qfp_ftan(x); -} -// Unit Tests: -// tan (pi/4) = +1.000000 -// tan(3*pi/4) = -1.000000 -// tan(5*pi/4) = +1.000000 -// tan(7*pi/4) = -1.000000 -// tan(+0) = 0.000000 -// tan(-0) = -0.000000 -//// tan(INFINITY) = -nan - -// CMakeFiles/STM32_BLUE_PILL.dir/pxtapp/base/core.cpp.o: In function -// `Math_::atan(pxt::TValueStruct*)': /src/pxtapp/base/core.cpp:916: undefined -// reference to `atan' - -double atan2(double y, double x) { - float_usage[USAGE_ATAN2]++; - // If x and y are both 0, return 0. - if (qfp_fcmp(y, 0) == 0 && qfp_fcmp(x, 0) == 0) { - return 0; - } - return qfp_fatan2(y, x); -} -// Unit Tests: -// atan2(+1,+1) = 0.785398 -// atan2(+1,-1) = 2.356194 -// atan2(-1,-1) = -2.356194 -// atan2(-1,+1) = -0.785398 -// atan2(0, 0) = 0.000000 -//// atan2(0, -0)=3.141593 -// atan2(7, 0) = 1.570796 -// atan2(7, -0)=1.570796 - -//// TODO: Confirm -double atan(double y_over_x) { - float_usage[USAGE_ATAN]++; - // If the argument is NaN, NaN is returned - if (isnan(y_over_x)) { - return NAN; - } - - // If the argument is ±0, it is returned unmodified - if (qfp_fcmp(y_over_x, 0) == 0) { - return y_over_x; - } - - // If the argument is +∞, +π/2 is returned - if (isinf(y_over_x) && qfp_fcmp(y_over_x, 0) > 0) { - return M_PI_2; - } - - // If the argument is -∞, -π/2 is returned - if (isinf(y_over_x) && qfp_fcmp(y_over_x, 0) < 0) { - return M_PI_2; - } - - // Must call atan2 instead of qfp_fatan2 in case the values are 0. - return atan2(y_over_x, 1); -} -// Unit Tests: -// atan(1) = 0.785398 -//// atan(Inf) = 1.570796 -// atan(-0.0) = -0.000000 -// atan(+0.0) = +0.000000 - -// CMakeFiles/STM32_BLUE_PILL.dir/pxtapp/base/core.cpp.o: In function -// `Math_::asin(pxt::TValueStruct*)': /src/pxtapp/base/core.cpp:919: undefined -// reference to `asin' - -// arcsin(x) = arctan( x / sqrt( 1 - x^2 ) ) -// = arctan2( x , sqrt( 1 - (x*x) ) ) -double asin(double x) { - float_usage[USAGE_ASIN]++; - // If the argument is NaN, NaN is returned - if (isnan(x)) { - return NAN; - } - - // If the argument is ±0, it is returned unmodified - if (qfp_fcmp(x, 0) == 0) { - return x; - } - - // If |arg| > 1, a domain error occurs and NaN is returned. - if (qfp_fcmp(x, 1) > 0) { - return NAN; - } - if (qfp_fcmp(x, -1) < 0) { - return NAN; - } - - // Must call atan2 instead of qfp_fatan2 in case the values are 0. - return atan2(x, qfp_fsqrt_fast(qfp_fsub(1, qfp_fmul(x, x)))); -} -// Unit Tests: -// asin( 1.0) = +1.570796 -// asin(-0.5) = -0.523599 -// asin(0.0) = 0.000000 -// asin(-0.0)=-0.000000 -//// asin(1.1) = nan - -// CMakeFiles/STM32_BLUE_PILL.dir/pxtapp/base/core.cpp.o: In function -// `Math_::acos(pxt::TValueStruct*)': /src/pxtapp/base/core.cpp:922: undefined -// reference to `acos' - -// arccos(x) = 2 * arctan( -// sqrt( 1 - x^2 ) / -// ( 1 + x ) -// ) where -1 < x <= 1 -double acos(double x) { - float_usage[USAGE_ACOS]++; - // if the argument is NaN, NaN is returned - if (isnan(x)) { - return NAN; - } - - // If the argument is +1, the value +0 is returned. - if (qfp_fcmp(x, 1) == 0) { - return 0; - } - - // If |arg| > 1, a domain error occurs and NaN is returned. - if (qfp_fcmp(x, 1) > 0) { - return NAN; - } - if (qfp_fcmp(x, -1) < 0) { - return NAN; - } - - // If x is -1, return Pi. - if (qfp_fcmp(x, -1) == 0) { - return M_PI_2 * 2.0; - } - - // Must call atan2 instead of qfp_fatan2 in case the values are 0. - return qfp_fmul(2.0f, atan2(qfp_fsqrt_fast(qfp_fsub(1.0f, qfp_fmul(x, x))), - qfp_fadd(1.0f, x))); -} -// Unit Tests: -// acos(-1) = 3.141593 -// acos(-0.9) = 2.690565 -// acos(0.0) = 1.570796 -// acos(0.5) = 1.047198 -// acos(1) = 0.000000 -//// acos(1.1) = nan - -// From https://en.wikipedia.org/wiki/Hyperbolic_function - -// Computes hyperbolic sine of arg. -// 0.5 * { e^x - e^{-x} } -double sinh(double x) { - float_usage[USAGE_SINH]++; - // if the argument is NaN, NaN is returned - if (isnan(x)) { - return NAN; - } - - // If the argument is ±0, it is returned unmodified - if (qfp_fcmp(x, 0) == 0) { - return x; - } - - // If the argument is ±∞, it is returned unmodified - if (isinf(x)) { - return x; - } - - return qfp_fmul(0.5f, qfp_fsub(qfp_fexp(x), qfp_fexp(-x))); -} -// Unit Tests: -// sinh(1) = 1.175201 -// sinh(-1)=-1.175201 -// sinh(+0) = 0.000000 -// sinh(-0)=-0.000000 -//// sinh(710.5) = inf - -// Computes the hyperbolic cosine of arg. -// 0.5 * { e^x + e^{-x} } -double cosh(double x) { - float_usage[USAGE_COSH]++; - // if the argument is NaN, NaN is returned - if (isnan(x)) { - return NAN; - } - - // If the argument is ±0, 1 is returned - if (qfp_fcmp(x, 0) == 0) { - return 1; - } - - // If the argument is ±∞, return +∞ - if (isinf(x)) { - return x; - } // TODO: Handle -∞ - - return qfp_fmul(0.5f, qfp_fadd(qfp_fexp(x), qfp_fexp(-x))); -} -// Unit Tests: -// cosh(1) = 1.543081 -// cosh(-1)= 1.543081 -// cosh(+0) = 1.000000 -// cosh(-0) = 1.000000 -//// cosh(710.5) = inf - -// Computes the hyperbolic tangent of arg. -// { e^{2x} - 1 } / { e^{2x} + 1 } -double tanh(double x) { - float_usage[USAGE_TANH]++; - // if the argument is NaN, NaN is returned - if (isnan(x)) { - return NAN; - } - - // If the argument is ±0, ±0 is returned - if (qfp_fcmp(x, 0) == 0) { - return x; - } - - // If the argument is ±∞, return ±1 - if (isinf(x)) { - return 1; - } // TODO: Handle -1 - - float e2x = qfp_fexp( // e^{2x} - qfp_fmul(2.0f, x)); - return qfp_fdiv_fast(qfp_fsub(e2x, 1.0f), qfp_fadd(e2x, 1.0f)); -} -// Unit Tests: -// tanh(1) = 0.761594 -// tanh(-1) = -0.761594 -// tanh(+0) = 0.000000 -// tanh(-0) = -0.000000 - -// From https://en.wikipedia.org/wiki/Inverse_hyperbolic_functions - -// Computes the inverse hyperbolic sine of arg. -// ln ( x + sqrt{ x^2 + 1 } ) -double asinh(double x) { - float_usage[USAGE_ASINH]++; - // if the argument is NaN, NaN is returned - if (isnan(x)) { - return NAN; - } - - // If the argument is ±0, it is returned unmodified - if (qfp_fcmp(x, 0) == 0) { - return x; - } - - // If the argument is ±∞, it is returned unmodified - if (isinf(x)) { - return x; - } - - return qfp_fln(qfp_fadd(x, - qfp_fsqrt_fast( // sqrt{ x^2 + 1 } - qfp_fadd(qfp_fmul(x, x), 1.0f)))); -} -// Unit Tests: -// asinh(1) = 0.881374 -// asinh(-1) = -0.881374 -// asinh(+0) = 0.000000 -// asinh(-0) = -0.000000 - -// Computes the inverse hyperbolic cosine of arg. -// ln ( x + sqrt{ x^2 - 1 } ) -double acosh(double x) { - float_usage[USAGE_ACOSH]++; - // if the argument is NaN, NaN is returned - if (isnan(x)) { - return NAN; - } - - // If the argument is less than 1, NaN is returned - if (qfp_fcmp(x, 1) < 0) { - return NAN; - } - - // If the argument is 1, 0 is returned - if (qfp_fcmp(x, 1) == 0) { - return 0; - } - - // If the argument is +∞, it is returned unmodified - if (isinf(x)) { - return x; - } - - return qfp_fln(qfp_fadd(x, - qfp_fsqrt_fast( // sqrt{ x^2 - 1 } - qfp_fsub(qfp_fmul(x, x), 1.0f)))); -} -// Unit Tests: -// acosh(1) = 0.000000 -// acosh(10) = 2.993223 -//// acosh(DBL_MAX) = 710.475860 -//// acosh(Inf) = inf -//// acosh(0.5) = -nan - -// Computes the inverse hyperbolic tangent of arg. -// 0.5 * ln ( {1+x} / {1-x} ) -double atanh(double x) { - float_usage[USAGE_ATANH]++; - // if the argument is NaN, NaN is returned - if (isnan(x)) { - return NAN; - } - - // If the argument is 0, it is returned unmodified - if (qfp_fcmp(x, 0) == 0) { - return x; - } - - // If the argument is ±1, ±∞ is returned - if (qfp_fcmp(x, 1) == 0) { - return INFINITY; - } - if (qfp_fcmp(x, -1) == 0) { - return -INFINITY; - } - - // if |arg| > 1, NaN is returned - if (qfp_fcmp(x, 1) > 0) { - return NAN; - } - if (qfp_fcmp(x, -1) < 0) { - return NAN; - } - - return qfp_fmul(0.5f, - qfp_fln( // ln ( {1+x} / {1-x} ) - qfp_fdiv_fast(qfp_fadd(1.0f, x), qfp_fsub(1.0f, x)))); -} -// Unit Tests: -// atanh(0) = 0.000000 -// atanh(-0) = -0.000000 -// atanh(0.9) = 1.472219 -//// atanh(-1) = -inf - -// CMakeFiles/STM32_BLUE_PILL.dir/pxtapp/base/core.cpp.o: In function -// `Math_::trunc(pxt::TValueStruct*)': /src/pxtapp/base/core.cpp:934: undefined -// reference to `trunc' - -// Computes the nearest integer not greater in magnitude than x. -// TODO: Warn if number is out of 32-bit int range. -double trunc(double x) { - float_usage[USAGE_TRUNC]++; - // If arg is NaN, NaN is returned - if (isnan(x)) { - return NAN; - } - - // If arg is ±∞, it is returned, unmodified - if (isinf(x)) { - return x; - } - - // If arg is ±0, it is returned, unmodified - if (qfp_fcmp(x, 0) == 0) { - return x; - } - - // qfp_float2int() behaves like floor(). If negative, add one. - int floored = qfp_float2int(x); - if (floored < 0) { - return floored + 1; - } - return floored; -} -// Unit Tests: -// trunc(+2.7) = +2.0 -// trunc(-2.7) = -2.0 -// trunc(-0.0) = -0.0 -// trunc(2205.1969) = 2205.000000 -// trunc(-270.8886) = -270.000000 - -// CMakeFiles/STM32_BLUE_PILL.dir/pxtapp/base/core.cpp.o: In function -// `Math_::floor(pxt::TValueStruct*)': /src/pxtapp/base/core.cpp:928: undefined -// reference to `floor' - -// Computes the largest integer value not greater than arg. -// TODO: Warn if number is out of 32-bit int range. -double floor(double x) { - float_usage[USAGE_FLOOR]++; - // If arg is NaN, NaN is returned - if (isnan(x)) { - return NAN; - } - - // If arg is ±∞, it is returned, unmodified - if (isinf(x)) { - return x; - } - - // If arg is ±0, it is returned, unmodified - if (qfp_fcmp(x, 0) == 0) { - return x; - } - - // qfp_float2int() behaves like floor(). - int floored = qfp_float2int(x); - return floored; -} -// Unit Tests: -// floor(+2.7) = +2.0 -// floor(-2.7) = -3.0 -// floor(-0.0) = -0.0 -// floor(2205.1969) = 2205.000000 -// floor(-270.8886) = -271.000000 - -// CMakeFiles/STM32_BLUE_PILL.dir/pxtapp/base/core.cpp.o: In function -// `Math_::ceil(pxt::TValueStruct*)': /src/pxtapp/base/core.cpp:931: undefined -// reference to `ceil' - -// Computes the smallest integer value not less than arg. -// TODO: Warn if number is out of 32-bit int range. -double ceil(double x) { - float_usage[USAGE_CEIL]++; - // If arg is NaN, NaN is returned - if (isnan(x)) { - return NAN; - } - - // If arg is ±∞, it is returned, unmodified - if (isinf(x)) { - return x; - } - - // If arg is ±0, it is returned, unmodified - if (qfp_fcmp(x, 0) == 0) { - return x; - } - - // qfp_float2int() behaves like floor(). Always add one. - int floored = qfp_float2int(x); - return floored + 1; -} -// Unit Tests: -// ceil(+2.4) = +3.0 -// ceil(-2.4) = -2.0 -// ceil(-0.0) = -0.0 -// ceil(2205.1969) = 2206.000000 -// ceil(-270.8886) = -270.000000 - -// Computes the floating-point remainder of the division operation x/y -// i.e. x - n*y, where n is x/y with its fractional part truncated. -double fmod(double x, double y) { - float_usage[USAGE_FMOD]++; - // If either argument is NaN, NaN is returned - if (isnan(x) || isnan(y)) { - return NAN; - } - - // If x is ±0 and y is not zero, ±0 is returned - if (qfp_fcmp(x, 0) == 0 && qfp_fcmp(y, 0) != 0) { - return 0; - } - - // If x is ±∞ and y is not NaN, NaN is returned and FE_INVALID is raised - if (isinf(x) && !isnan(y)) { - return NAN; - } - - // If y is ±0 and x is not NaN, NaN is returned and FE_INVALID is raised - if (qfp_fcmp(y, 0) == 0 && !isnan(x)) { - return NAN; - } - - // If y is ±∞ and x is finite, x is returned. - if (isinf(y) && !isinf(x)) { - return x; - } - - // From https://en.cppreference.com/w/c/numeric/math/fmod - double xabs = fabs(x); - double yabs = fabs(y); - // Was: double result = remainder(fabs(x), (y = fabs(y))); - double n = trunc(qfp_fdiv_fast(xabs, yabs)); - float result = qfp_fsub(xabs, qfp_fmul(n, yabs)); // x - n*y, always positive - - // Was: if (signbit(result)) result += y; - if (qfp_fcmp(result, 0) < 0) { - result = qfp_fadd(result, yabs); - } - - // Composes a floating point value with the magnitude of result and the sign - // of x. Was: return copysign(result, x); - return (qfp_fcmp(x, 0) < 0) ? -result : result; -} -// Unit Tests: -// fmod(+5.1, +3.0) = 2.1 -// fmod(-5.1, +3.0) = -2.1 -// fmod(+5.1, -3.0) = 2.1 -// fmod(-5.1, -3.0) = -2.1 -// fmod(+0.0, 1.0) = 0.0 -// fmod(-0.0, 1.0) = -0.0 -// fmod(+5.1, Inf) = 5.1 -//// fmod(+5.1, 0) = nan - -// Computes the absolute value of a floating point value arg. -double fabs(double x) { - float_usage[USAGE_FABS]++; - if (isnan(x) || isinf(x) || qfp_fcmp(x, 0) == 0) { - return x; - } - if (qfp_fcmp(x, 0) < 0) { - return -x; - } - return x; -} -// Unit Tests: -// fabs(+3) = 3.000000 -// fabs(-3) = 3.000000 -// fabs(-0) = 0.000000 -// fabs(2205.1969) = 2205.1969 -// fabs(-270.8886) = 270.8886 -//// fabs(-Inf) = inf diff --git a/keymaps/onekey/stm32f411/qfplib/qfplib.s b/keymaps/onekey/stm32f411/qfplib/qfplib.s index f21812a..e035472 100644 --- a/keymaps/onekey/stm32f411/qfplib/qfplib.s +++ b/keymaps/onekey/stm32f411/qfplib/qfplib.s @@ -1,26 +1,5 @@ -@ Copyright 2015-2020 Mark Owen -@ http://www.quinapalus.com -@ E-mail: [email protected] -@ -@ This file is free software: you can redistribute it and/or modify -@ it under the terms of version 2 of the GNU General Public License -@ as published by the Free Software Foundation. -@ -@ This file is distributed in the hope that it will be useful, -@ but WITHOUT ANY WARRANTY; without even the implied warranty of -@ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the -@ GNU General Public License for more details. -@ -@ You should have received a copy of the GNU General Public License -@ along with this file. If not, see <http://www.gnu.org/licenses/> or -@ write to the Free Software Foundation, Inc., 51 Franklin Street, -@ Fifth Floor, Boston, MA 02110-1301, USA. - -@ Thanks to Li Ling for optimisation suggestions and discussion. - [email protected] include_faster,0 @ include fast divide and square root? [email protected] include_conversions,1 @ include float <-> fixed point conversion functions? [email protected] include_scientific,1 @ include trignometic, exponential etc. functions? +// Copyright 2015 Mark Owen ([email protected]) +// SPDX-License-Identifier: GPL-2.0-or-later .ifndef include_faster .equ include_faster,1 @@ -42,8 +21,6 @@ .cpu cortex-m0 .thumb -@ exported symbols - .global qfp_fadd .global qfp_fsub .global qfp_fmul @@ -74,15 +51,15 @@ .global qfp_fsqrt_fast .endif -@ exchange r0<->r1, r2<->r3 +// exchange r0<->r1, r2<->r3 xchxy: push {r0,r2,r14} mov r0,r1 mov r2,r3 pop {r1,r3,r15} -@ IEEE single precision floats in r0,r1-> mantissae in r1,r0 exponents in r3,r2 *respectively* -@ trashes r4 +// IEEE single precision floats in r0,r1-> mantissae in r1,r0 exponents in r3,r2 *respectively* +// trashes r4 unpackxy: push {r14} bl unpackx @@ -90,192 +67,192 @@ unpackxy: pop {r4} mov r14,r4 -@ IEEE single in r0-> signed (two's complemennt) mantissa in r0 9Q23 (24 significant bits), signed exponent (bias removed) in r2 -@ trashes r4; zero, denormal -> mantissa=+/-1, exponent=-380; Inf, NaN -> mantissa=+/-1, exponent=+640 +// IEEE single in r0-> signed (two's complemennt) mantissa in r0 9Q23 (24 significant bits), signed exponent (bias removed) in r2 +// trashes r4; zero, denormal -> mantissa=+/-1, exponent=-380; Inf, NaN -> mantissa=+/-1, exponent=+640 unpackx: - lsrs r2,r0,#23 @ save exponent and sign - lsls r0,#9 @ extract mantissa + lsrs r2,r0,#23 // save exponent and sign + lsls r0,#9 // extract mantissa lsrs r0,#9 movs r4,#1 lsls r4,#23 - orrs r0,r4 @ reinstate implied leading 1 - cmp r2,#255 @ test sign bit - uxtb r2,r2 @ clear it - bls 1f @ branch on positive - rsbs r0,#0 @ negate mantissa + orrs r0,r4 // reinstate implied leading 1 + cmp r2,#255 // test sign bit + uxtb r2,r2 // clear it + bls 1f // branch on positive + rsbs r0,#0 // negate mantissa 1: subs r2,#1 - cmp r2,#254 @ zero/denormal/Inf/NaN? + cmp r2,#254 // zero/denormal/Inf/NaN? bhs 2f - subs r2,#126 @ remove exponent bias: can now be -126..+127 + subs r2,#126 // remove exponent bias: can now be -126..+127 bx r14 -2: @ here with special-case values +2: // here with special-case values cmp r0,#0 - mov r0,r4 @ set mantissa to +1 + mov r0,r4 // set mantissa to +1 bpl 3f - rsbs r0,#0 @ zero/denormal/Inf/NaN: mantissa=+/-1 + rsbs r0,#0 // zero/denormal/Inf/NaN: mantissa=+/-1 3: - subs r2,#126 @ zero/denormal: exponent -> -127; Inf, NaN: exponent -> 128 - lsls r2,#2 @ zero/denormal: exponent -> -508; Inf, NaN: exponent -> 512 - adds r2,#128 @ zero/denormal: exponent -> -380; Inf, NaN: exponent -> 640 + subs r2,#126 // zero/denormal: exponent -> -127; Inf, NaN: exponent -> 128 + lsls r2,#2 // zero/denormal: exponent -> -508; Inf, NaN: exponent -> 512 + adds r2,#128 // zero/denormal: exponent -> -380; Inf, NaN: exponent -> 640 bx r14 -@ normalise and pack signed mantissa in r0 nominally 3Q29, signed exponent in r2-> IEEE single in r0 -@ trashes r4, preserves r1,r3 -@ r5: "sticky bits", must be zero iff all result bits below r0 are zero for correct rounding +// normalise and pack signed mantissa in r0 nominally 3Q29, signed exponent in r2-> IEEE single in r0 +// trashes r4, preserves r1,r3 +// r5: "sticky bits", must be zero iff all result bits below r0 are zero for correct rounding packx: - lsrs r4,r0,#31 @ save sign bit - lsls r4,r4,#31 @ sign now in b31 - bpl 2f @ skip if positive + lsrs r4,r0,#31 // save sign bit + lsls r4,r4,#31 // sign now in b31 + bpl 2f // skip if positive cmp r5,#0 beq 11f - adds r0,#1 @ fiddle carry in to following rsb if sticky bits are non-zero + adds r0,#1 // fiddle carry in to following rsb if sticky bits are non-zero 11: - rsbs r0,#0 @ can now treat r0 as unsigned + rsbs r0,#0 // can now treat r0 as unsigned packx0: - bmi 3f @ catch r0=0x80000000 case + bmi 3f // catch r0=0x80000000 case 2: - subs r2,#1 @ normalisation loop + subs r2,#1 // normalisation loop adds r0,r0 - beq 1f @ zero? special case - bpl 2b @ normalise so leading "1" in bit 31 + beq 1f // zero? special case + bpl 2b // normalise so leading "1" in bit 31 3: - adds r2,#129 @ (mis-)offset exponent - bne 12f @ special case: highest denormal can round to lowest normal - adds r0,#0x80 @ in special case, need to add 256 to r0 for rounding - bcs 4f @ tripped carry? then have leading 1 in C as required + adds r2,#129 // (mis-)offset exponent + bne 12f // special case: highest denormal can round to lowest normal + adds r0,#0x80 // in special case, need to add 256 to r0 for rounding + bcs 4f // tripped carry? then have leading 1 in C as required 12: - adds r0,#0x80 @ rounding - bcs 4f @ tripped carry? then have leading 1 in C as required (and result is even so can ignore sticky bits) + adds r0,#0x80 // rounding + bcs 4f // tripped carry? then have leading 1 in C as required (and result is even so can ignore sticky bits) cmp r5,#0 - beq 7f @ sticky bits zero? + beq 7f // sticky bits zero? 8: - lsls r0,#1 @ remove leading 1 + lsls r0,#1 // remove leading 1 9: - subs r2,#1 @ compensate exponent on this path + subs r2,#1 // compensate exponent on this path 4: cmp r2,#254 - bge 5f @ overflow? - adds r2,#1 @ correct exponent offset - ble 10f @ denormal/underflow? - lsrs r0,#9 @ align mantissa - lsls r2,#23 @ align exponent - orrs r0,r2 @ assemble exponent and mantissa + bge 5f // overflow? + adds r2,#1 // correct exponent offset + ble 10f // denormal/underflow? + lsrs r0,#9 // align mantissa + lsls r2,#23 // align exponent + orrs r0,r2 // assemble exponent and mantissa 6: - orrs r0,r4 @ apply sign + orrs r0,r4 // apply sign 1: bx r14 5: - movs r0,#0xff @ create infinity + movs r0,#0xff // create infinity lsls r0,#23 b 6b 10: - movs r0,#0 @ create zero + movs r0,#0 // create zero bx r14 -7: @ sticky bit rounding case - lsls r5,r0,#24 @ check bottom 8 bits of r0 - bne 8b @ in rounding-tie case? - lsrs r0,#9 @ ensure even result +7: // sticky bit rounding case + lsls r5,r0,#24 // check bottom 8 bits of r0 + bne 8b // in rounding-tie case? + lsrs r0,#9 // ensure even result lsls r0,#10 b 9b -@ unpack two arguments (r0,r1) and shift one down to have common exponent, returned in r2; note that arguments are exchanged -@ sticky bits shifted off bottom of smaller argument in r5 -@ following code is unnecessarily general for fadd, but is shared with atan2 +// unpack two arguments (r0,r1) and shift one down to have common exponent, returned in r2; note that arguments are exchanged +// sticky bits shifted off bottom of smaller argument in r5 +// following code is unnecessarily general for fadd, but is shared with atan2 unpackxyalign: push {r14} bl unpackxy - lsls r0,r0,#6 @ Q29 - lsls r1,r1,#6 @ Q29 - subs r4,r2,r3 @ calculate shift - bge 1f @ x>=y? - mov r2,r3 @ no: take common exponent from y - mov r5,r0 @ potential sticky bits from x - rsbs r4,#0 @ make shift positive + lsls r0,r0,#6 // Q29 + lsls r1,r1,#6 // Q29 + subs r4,r2,r3 // calculate shift + bge 1f // x>=y? + mov r2,r3 // no: take common exponent from y + mov r5,r0 // potential sticky bits from x + rsbs r4,#0 // make shift positive asrs r0,r4 cmp r4,#32 blo 2f - movs r0,#0 @ large shift, so all bits are sticky and result is zero + movs r0,#0 // large shift, so all bits are sticky and result is zero pop {r15} 1: - mov r5,r1 @ common exponent from x; potential sticky bits from y + mov r5,r1 // common exponent from x; potential sticky bits from y asrs r1,r4 cmp r4,#32 blo 2f - movs r1,#0 @ large shift, so all bits are sticky and result is zero + movs r1,#0 // large shift, so all bits are sticky and result is zero pop {r15} 2: rsbs r4,#0 adds r4,#32 - lsls r5,r4 @ extract sticky bits + lsls r5,r4 // extract sticky bits pop {r15} .thumb_func qfp_fsub: - movs r2,#1 @ subtract: flip sign bit of second argument and fall through to fadd + movs r2,#1 // subtract: flip sign bit of second argument and fall through to fadd lsls r2,#31 eors r1,r2 .thumb_func qfp_fadd: push {r4,r5,r14} bl unpackxyalign - adds r0,r1 @ do addition - bne 2f @ not in Inf-Inf case? + adds r0,r1 // do addition + bne 2f // not in Inf-Inf case? cmp r2,#200 blt 2f movs r0,#1 - lsls r0,#29 @ for Inf-Inf, set mantissa to +1 to prevent zero result + lsls r0,#29 // for Inf-Inf, set mantissa to +1 to prevent zero result 2: -packret: @ common return point: "pack and return" +packret: // common return point: "pack and return" bl packx pop {r4,r5,r15} -@ signed multiply r0 1Q23 by r1 4Q23, result in r0 7Q25, sticky bits in r5 -@ trashes r3,r4 +// signed multiply r0 1Q23 by r1 4Q23, result in r0 7Q25, sticky bits in r5 +// trashes r3,r4 mul0: - uxth r3,r0 @ Q23 - asrs r4,r1,#16 @ Q7 - muls r3,r4 @ L*H, Q30 signed - asrs r4,r0,#16 @ Q7 - uxth r5,r1 @ Q23 - muls r4,r5 @ H*L, Q30 signed - adds r3,r4 @ sum of middle partial products + uxth r3,r0 // Q23 + asrs r4,r1,#16 // Q7 + muls r3,r4 // L*H, Q30 signed + asrs r4,r0,#16 // Q7 + uxth r5,r1 // Q23 + muls r4,r5 // H*L, Q30 signed + adds r3,r4 // sum of middle partial products uxth r4,r0 - muls r4,r5 @ L*L, Q46 unsigned - lsls r5,r4,#16 @ initialise sticky bits from low half of low partial product - lsrs r4,#16 @ Q25 - adds r3,r4 @ add high half of low partial product to sum of middle partial products - @ (cannot generate carry by limits on input arguments) - asrs r0,#16 @ Q7 - asrs r1,#16 @ Q7 - muls r0,r1 @ H*H, Q14 signed - lsls r0,#11 @ high partial product Q25 - lsls r1,r3,#27 @ sticky - orrs r5,r1 @ collect further sticky bits - asrs r1,r3,#5 @ middle partial products Q25 - adds r0,r1 @ final result + muls r4,r5 // L*L, Q46 unsigned + lsls r5,r4,#16 // initialise sticky bits from low half of low partial product + lsrs r4,#16 // Q25 + adds r3,r4 // add high half of low partial product to sum of middle partial products + // (cannot generate carry by limits on input arguments) + asrs r0,#16 // Q7 + asrs r1,#16 // Q7 + muls r0,r1 // H*H, Q14 signed + lsls r0,#11 // high partial product Q25 + lsls r1,r3,#27 // sticky + orrs r5,r1 // collect further sticky bits + asrs r1,r3,#5 // middle partial products Q25 + adds r0,r1 // final result bx r14 .thumb_func qfp_fcmp: - movs r2,#1 @ initialise result - lsls r3,r2,#31 @ r3=0x80000000 - tst r0,r3 @ check sign of first argument + movs r2,#1 // initialise result + lsls r3,r2,#31 // r3=0x80000000 + tst r0,r3 // check sign of first argument beq 1f - subs r0,r3,r0 @ convert to 2's complement form for direct comparison + subs r0,r3,r0 // convert to 2's complement form for direct comparison 1: - tst r1,r3 @ repeat for second argument + tst r1,r3 // repeat for second argument beq 2f subs r1,r3,r1 2: - subs r0,r1 @ perform comparison - beq 4f @ equal? return 0 - bgt 3f @ r0>r1? return +1 - rsbs r2,#0 @ r0<r1: return -1 + subs r0,r1 // perform comparison + beq 4f // equal? return 0 + bgt 3f // r0>r1? return +1 + rsbs r2,#0 // r0<r1: return -1 3: mov r0,r2 4: @@ -285,8 +262,8 @@ qfp_fcmp: qfp_fmul: push {r4,r5,r14} bl unpackxy - add r2,r3 @ sum exponents - adds r2,#4 @ adjust exponent for pack which expects Q29 + add r2,r3 // sum exponents + adds r2,#4 // adjust exponent for pack which expects Q29 fmul0: bl mul0 b packret @@ -296,56 +273,56 @@ qfp_fdiv: push {r4,r5,r14} fdiv_n: bl unpackxy - movs r5,#1 @ result cannot fall exactly half-way between two representable numbers (exercise for reader - note that - @ we do not handle denormals) so there will always be sticky bits - cmp r0,#0 @ check divisor sign + movs r5,#1 // result cannot fall exactly half-way between two representable numbers (exercise for reader - note that + // we do not handle denormals) so there will always be sticky bits + cmp r0,#0 // check divisor sign bpl 1f - rsbs r0,#0 @ ensure divisor positive - rsbs r1,#0 @ preserve sign of result + rsbs r0,#0 // ensure divisor positive + rsbs r1,#0 // preserve sign of result 1: movs r4,#0 - cmp r1,#0 @ check sign of dividend + cmp r1,#0 // check sign of dividend bpl 2f - rsbs r1,#0 @ result will be negative - mvns r4,r4 @ save sign as 0x00000000 or 0xffffffff -2: @ now do unsigned division on unpacked values {r1,r3}/{r0,r2} - cmp r3,#200 @ dividend is an infinity? return it + rsbs r1,#0 // result will be negative + mvns r4,r4 // save sign as 0x00000000 or 0xffffffff +2: // now do unsigned division on unpacked values {r1,r3}/{r0,r2} + cmp r3,#200 // dividend is an infinity? return it bge 3f rsbs r2,#0 - cmp r2,#200 @ divisor is zero? return infinity + cmp r2,#200 // divisor is zero? return infinity bge 3f - adds r2,r3 @ difference of exponents - movs r3,#0x40 @ termination marker (calculate enough bits to do rounding correctly) -2: @ division loop - subs r1,r0 @ trial subtraction + adds r2,r3 // difference of exponents + movs r3,#0x40 // termination marker (calculate enough bits to do rounding correctly) +2: // division loop + subs r1,r0 // trial subtraction bcs 1f - add r1,r0 @ restore if failed + add r1,r0 // restore if failed 1: - adcs r3,r3 @ shift in result bit - add r1,r1 @ shift up dividend - bcc 2b @ loop until marker appears in carry - lsls r0,r3,#4 @ align for packing + adcs r3,r3 // shift in result bit + add r1,r1 // shift up dividend + bcc 2b // loop until marker appears in carry + lsls r0,r3,#4 // align for packing 4: - eors r0,r4 @ apply sign + eors r0,r4 // apply sign b packret infret: movs r4,#0 3: mov r0,r1 - movs r2,#255 @ return infinity + movs r2,#255 // return infinity b 4b .if include_faster -@ The fast divide routine uses an initial approximation to the reciprocal of the divisor based on the top four bits of the mantissa -@ followed by three Newton-Raphson iterations, resulting in about 27 bits of accuracy. This reciprocal is then multiplied by -@ the dividend. -@ The (fixed-point) reciprocal calculation is carefully implemented to preserve the necessary accuracy throughout. In places -@ the implied binary point is not within the stored value. For example where "Q47" is shown below it means that the least -@ significant bit of the value has significance 2^-47 and hence the most significant bit has significance 2^-16. In these -@ cases the value is known to be very close to an integer (usually 1) and so the bits of greater significance do not need -@ to be stored. -@ The reciprocal calculation has been tested for all possible input mantissa values. +// The fast divide routine uses an initial approximation to the reciprocal of the divisor based on the top four bits of the mantissa +// followed by three Newton-Raphson iterations, resulting in about 27 bits of accuracy. This reciprocal is then multiplied by +// the dividend. +// The (fixed-point) reciprocal calculation is carefully implemented to preserve the necessary accuracy throughout. In places +// the implied binary point is not within the stored value. For example where "Q47" is shown below it means that the least +// significant bit of the value has significance 2^-47 and hence the most significant bit has significance 2^-16. In these +// cases the value is known to be very close to an integer (usually 1) and so the bits of greater significance do not need +// to be stored. +// The reciprocal calculation has been tested for all possible input mantissa values. .thumb_func qfp_fdiv_fast: push {r4,r5,r14} @@ -354,115 +331,115 @@ fdiv_fast_n: cmp r0,#0 bpl 1f rsbs r0,#0 - rsbs r1,#0 @ ensure divisor positive + rsbs r1,#0 // ensure divisor positive 1: cmp r3,#200 - bge infret @ dividend is an infinity? return it + bge infret // dividend is an infinity? return it rsbs r2,#0 - cmp r2,#200 @ divisor is zero? - bge infret @ return infinity - adds r2,r3 @ difference of exponents - adr r4,rcpapp-8 @ the first 8 elements of the table are never accessed because of the mantissa's leading 1 - lsrs r3,r0,#20 @ y Q23; y>>20 Q7 - ldrb r4,[r4,r3] @ m=rcpapp[(y>>20)&7]; // Q8, .5<m<1 - - lsls r3,r4,#2 @ m<<2 // Q10 first Newton-Raphson iteration - muls r3,r0 @ s=y*(m<<2); // Q33 - lsls r4,#8 @ m<<=8; // Q16 - asrs r3,#21 @ s>>=21; // Q12 - muls r3,r4 @ s*=m; // Q28 - asrs r3,#12 @ s>>=12; // Q16 - subs r4,r3 @ m=m-s; // Q16 - - mov r3,r4 @ s=y*m // Q39 second Newton-Raphson iteration - muls r4,r0 @ ... - asrs r4,#16 @ s>>=16; // Q23 - muls r4,r3 @ s*=m; // Q39 - lsls r3,#8 @ m<<=8; // Q24 - asrs r4,#15 @ s>>=15; // Q24 - subs r3,r4 @ m=m-s; // Q24 - - lsls r4,r3,#7 @ \/ s=y*m; // Q47 third Newton-Raphson iteration - muls r3,r0 @ /\ m<<=7; // Q31 - asrs r3,#15 @ s>>=15; // Q32 - lsrs r0,r4,#16 @ s*=(m>>16); // Q47 - muls r3,r0 @ ... - asrs r3,#16 @ s>>=16; // Q31 - subs r0,r4,r3 @ m=m-s; // Q31 + cmp r2,#200 // divisor is zero? + bge infret // return infinity + adds r2,r3 // difference of exponents + adr r4,rcpapp-8 // the first 8 elements of the table are never accessed because of the mantissa's leading 1 + lsrs r3,r0,#20 // y Q23; y>>20 Q7 + ldrb r4,[r4,r3] // m=rcpapp[(y>>20)&7]; // Q8, .5<m<1 + + lsls r3,r4,#2 // m<<2 // Q10 first Newton-Raphson iteration + muls r3,r0 // s=y*(m<<2); // Q33 + lsls r4,#8 // m<<=8; // Q16 + asrs r3,#21 // s>>=21; // Q12 + muls r3,r4 // s*=m; // Q28 + asrs r3,#12 // s>>=12; // Q16 + subs r4,r3 // m=m-s; // Q16 + + mov r3,r4 // s=y*m // Q39 second Newton-Raphson iteration + muls r4,r0 // ... + asrs r4,#16 // s>>=16; // Q23 + muls r4,r3 // s*=m; // Q39 + lsls r3,#8 // m<<=8; // Q24 + asrs r4,#15 // s>>=15; // Q24 + subs r3,r4 // m=m-s; // Q24 + + lsls r4,r3,#7 // \/ s=y*m; // Q47 third Newton-Raphson iteration + muls r3,r0 // /\ m<<=7; // Q31 + asrs r3,#15 // s>>=15; // Q32 + lsrs r0,r4,#16 // s*=(m>>16); // Q47 + muls r3,r0 // ... + asrs r3,#16 // s>>=16; // Q31 + subs r0,r4,r3 // m=m-s; // Q31 div0: - adds r0,#7 @ rounding; reduce systematic error - lsrs r0,#4 @ Q27 - b fmul0 @ drop into multiplication code to calculate result - -@ The fast square root routine uses an initial approximation to the reciprocal of the square root of the argument based -@ on the top four bits of the mantissa (possibly shifted one place to make the exponent even). It then performs three -@ Newton-Raphson iterations, resulting in about 28-29 bits of accuracy. This reciprocal is then multiplied by -@ the original argument to produce the result. -@ Again, the fixed-point calculation is carefully implemented to preserve accuracy, and similar comments to those -@ made above on the fast division routine apply. -@ The reciprocal square root calculation has been tested for all possible (possibly shifted) input mantissa values. + adds r0,#7 // rounding; reduce systematic error + lsrs r0,#4 // Q27 + b fmul0 // drop into multiplication code to calculate result + +// The fast square root routine uses an initial approximation to the reciprocal of the square root of the argument based +// on the top four bits of the mantissa (possibly shifted one place to make the exponent even). It then performs three +// Newton-Raphson iterations, resulting in about 28-29 bits of accuracy. This reciprocal is then multiplied by +// the original argument to produce the result. +// Again, the fixed-point calculation is carefully implemented to preserve accuracy, and similar comments to those +// made above on the fast division routine apply. +// The reciprocal square root calculation has been tested for all possible (possibly shifted) input mantissa values. .thumb_func qfp_fsqrt_fast: push {r4,r5,r14} bl unpackx movs r1,r0 - bmi infret @ negative? return -Inf - asrs r0,r2,#1 @ check LSB of exponent + bmi infret // negative? return -Inf + asrs r0,r2,#1 // check LSB of exponent bcc 1f - lsls r1,#1 @ was odd: double mantissa; mantissa y now 1..4 Q23 + lsls r1,#1 // was odd: double mantissa; mantissa y now 1..4 Q23 1: - adds r2,#4 @ correction for packing - adr r4,rsqrtapp-4@ first four table entries are never accessed because of the mantissa's leading 1 - lsrs r3,r1,#21 @ y>>21 Q2 - ldrb r4,[r4,r3] @ initial approximation to reciprocal square root m Q8 - - lsrs r0,r1,#7 @ y>>7 // Q16 first Newton-Raphson iteration - muls r0,r4 @ m*y - muls r0,r4 @ s=m*y*y // Q32 - asrs r0,#12 @ s>>12 - muls r0,r4 @ m*s // Q28 - asrs r0,#13 @ m*s // Q15 - lsls r4,#8 @ m // Q16 - subs r4,r0 @ m=(m<<8)-(s>>13) // Q16-Q15/2 -> Q16 - - mov r0,r4 @ // second Newton-Raphson iteration - muls r0,r0 @ u=m*m // Q32 - lsrs r0,#16 @ u>>16 // Q16 - lsrs r3,r1,#7 @ y>>7 // Q16 - muls r0,r3 @ s=u*(y>>7) // Q32 - asrs r0,#12 @ s>>12 // Q20 - muls r0,r4 @ s*m // Q36 - asrs r0,#21 @ s*m // Q15 - subs r4,r0 @ m=m-s // Q16-Q15/2 - - mov r0,r4 @ // third Newton-Raphson iteration - muls r0,r0 @ u=m*m // Q32 - lsrs r3,r0,#12 @ now multiply u and y in two parts: u>>12 - muls r3,r1 @ first partial product (u>>12)*y Q43 + adds r2,#4 // correction for packing + adr r4,rsqrtapp-4// first four table entries are never accessed because of the mantissa's leading 1 + lsrs r3,r1,#21 // y>>21 Q2 + ldrb r4,[r4,r3] // initial approximation to reciprocal square root m Q8 + + lsrs r0,r1,#7 // y>>7 // Q16 first Newton-Raphson iteration + muls r0,r4 // m*y + muls r0,r4 // s=m*y*y // Q32 + asrs r0,#12 // s>>12 + muls r0,r4 // m*s // Q28 + asrs r0,#13 // m*s // Q15 + lsls r4,#8 // m // Q16 + subs r4,r0 // m=(m<<8)-(s>>13) // Q16-Q15/2 -> Q16 + + mov r0,r4 // // second Newton-Raphson iteration + muls r0,r0 // u=m*m // Q32 + lsrs r0,#16 // u>>16 // Q16 + lsrs r3,r1,#7 // y>>7 // Q16 + muls r0,r3 // s=u*(y>>7) // Q32 + asrs r0,#12 // s>>12 // Q20 + muls r0,r4 // s*m // Q36 + asrs r0,#21 // s*m // Q15 + subs r4,r0 // m=m-s // Q16-Q15/2 + + mov r0,r4 // // third Newton-Raphson iteration + muls r0,r0 // u=m*m // Q32 + lsrs r3,r0,#12 // now multiply u and y in two parts: u>>12 + muls r3,r1 // first partial product (u>>12)*y Q43 lsls r0,#20 - lsrs r0,#20 @ u&0xfff - lsrs r5,r1,#12 @ y>>12 - muls r0,r5 @ second partial product (u&0xfff)*(y>>12) Q43 - add r0,r3 @ s=u*y // Q43 - asrs r0,#15 @ s>>15 // Q28 - muls r0,r4 @ (s>>15)*m // Q44 - lsls r4,#13 @ m<<13 // Q29 - asrs r0,#16 @ s>>16 // Q28 - subs r0,r4,r0 @ // Q29-Q28/2 - - asrs r2,#1 @ halve exponent - bcc div0 @ was y shifted? + lsrs r0,#20 // u&0xfff + lsrs r5,r1,#12 // y>>12 + muls r0,r5 // second partial product (u&0xfff)*(y>>12) Q43 + add r0,r3 // s=u*y // Q43 + asrs r0,#15 // s>>15 // Q28 + muls r0,r4 // (s>>15)*m // Q44 + lsls r4,#13 // m<<13 // Q29 + asrs r0,#16 // s>>16 // Q28 + subs r0,r4,r0 // // Q29-Q28/2 + + asrs r2,#1 // halve exponent + bcc div0 // was y shifted? lsrs r0,#1 - lsls r1,#1 @ shift y back - b div0 @ round and complete with multiplication + lsls r1,#1 // shift y back + b div0 // round and complete with multiplication .align 2 -@ round(2^15./[136:16:248]) +// round(2^15./[136:16:248]) rcpapp: .byte 0xf1,0xd8,0xc3,0xb2, 0xa4,0x98,0x8d,0x84 -@ round(sqrt(2^22./[72:16:248])) +// round(sqrt(2^22./[72:16:248])) rsqrtapp: .byte 0xf1,0xda,0xc9,0xbb, 0xb0,0xa6,0x9e,0x97, 0x91,0x8b,0x86,0x82 @@ -470,86 +447,86 @@ rsqrtapp: .if include_conversions -@ convert float to signed int, rounding towards -Inf, clamping +// convert float to signed int, rounding towards -Inf, clamping .thumb_func qfp_float2int: - movs r1,#0 @ fall through + movs r1,#0 // fall through -@ convert float in r0 to signed fixed point in r0, clamping +// convert float in r0 to signed fixed point in r0, clamping .thumb_func qfp_float2fix: push {r4,r14} bl unpackx - add r2,r1 @ incorporate binary point position into exponent - subs r2,#23 @ r2 is now amount of left shift required - blt 1f @ requires right shift? - cmp r2,#7 @ overflow? + add r2,r1 // incorporate binary point position into exponent + subs r2,#23 // r2 is now amount of left shift required + blt 1f // requires right shift? + cmp r2,#7 // overflow? ble 4f -3: @ overflow - asrs r1,r0,#31 @ +ve:0 -ve:0xffffffff - mvns r1,r1 @ +ve:0xffffffff -ve:0 +3: // overflow + asrs r1,r0,#31 // +ve:0 -ve:0xffffffff + mvns r1,r1 // +ve:0xffffffff -ve:0 movs r0,#1 lsls r0,#31 5: - eors r0,r1 @ +ve:0x7fffffff -ve:0x80000000 (unsigned path: 0xffffffff) + eors r0,r1 // +ve:0x7fffffff -ve:0x80000000 (unsigned path: 0xffffffff) pop {r4,r15} 1: - rsbs r2,#0 @ right shift for r0, >0 + rsbs r2,#0 // right shift for r0, >0 cmp r2,#32 - blt 2f @ more than 32 bits of right shift? + blt 2f // more than 32 bits of right shift? movs r2,#32 2: asrs r0,r0,r2 pop {r4,r15} -@ unsigned version +// unsigned version .thumb_func qfp_float2uint: - movs r1,#0 @ fall through + movs r1,#0 // fall through .thumb_func qfp_float2ufix: push {r4,r14} bl unpackx - add r2,r1 @ incorporate binary point position into exponent + add r2,r1 // incorporate binary point position into exponent movs r1,r0 - bmi 5b @ negative? return zero - subs r2,#23 @ r2 is now amount of left shift required - blt 1b @ requires right shift? - mvns r1,r0 @ ready to return 0xffffffff - cmp r2,#8 @ overflow? + bmi 5b // negative? return zero + subs r2,#23 // r2 is now amount of left shift required + blt 1b // requires right shift? + mvns r1,r0 // ready to return 0xffffffff + cmp r2,#8 // overflow? bgt 5b 4: - lsls r0,r0,r2 @ result fits, left shifted + lsls r0,r0,r2 // result fits, left shifted pop {r4,r15} -@ convert signed int to float, rounding +// convert signed int to float, rounding .thumb_func qfp_int2float: - movs r1,#0 @ fall through + movs r1,#0 // fall through -@ convert signed fix to float, rounding; number of r0 bits after point in r1 +// convert signed fix to float, rounding; number of r0 bits after point in r1 .thumb_func qfp_fix2float: push {r4,r5,r14} 1: movs r2,#29 - subs r2,r1 @ fix exponent -packretns: @ pack and return, sticky bits=0 + subs r2,r1 // fix exponent +packretns: // pack and return, sticky bits=0 movs r5,#0 b packret -@ unsigned version +// unsigned version .thumb_func qfp_uint2float: - movs r1,#0 @ fall through + movs r1,#0 // fall through .thumb_func qfp_ufix2float: push {r4,r5,r14} cmp r0,#0 - bge 1b @ treat <2^31 as signed + bge 1b // treat <2^31 as signed movs r2,#30 - subs r2,r1 @ fix exponent - lsls r5,r0,#31 @ one sticky bit + subs r2,r1 // fix exponent + lsls r5,r0,#31 // one sticky bit lsrs r0,#1 b packret @@ -557,207 +534,207 @@ qfp_ufix2float: .if include_scientific -@ All the scientific functions are implemented using the CORDIC algorithm. For notation, -@ details not explained in the comments below, and a good overall survey see -@ "50 Years of CORDIC: Algorithms, Architectures, and Applications" by Meher et al., -@ IEEE Transactions on Circuits and Systems Part I, Volume 56 Issue 9. - -@ Register use: -@ r0: x -@ r1: y -@ r2: z/omega -@ r3: coefficient pointer -@ r4,r12: m -@ r5: i (shift) - -cordic_start: @ initialisation - movs r5,#0 @ initial shift=0 +// All the scientific functions are implemented using the CORDIC algorithm. For notation, +// details not explained in the comments below, and a good overall survey see +// "50 Years of CORDIC: Algorithms, Architectures, and Applications" by Meher et al., +// IEEE Transactions on Circuits and Systems Part I, Volume 56 Issue 9. + +// Register use: +// r0: x +// r1: y +// r2: z/omega +// r3: coefficient pointer +// r4,r12: m +// r5: i (shift) + +cordic_start: // initialisation + movs r5,#0 // initial shift=0 mov r12,r4 b 5f -cordic_vstep: @ one step of algorithm in vector mode - cmp r1,#0 @ check sign of y +cordic_vstep: // one step of algorithm in vector mode + cmp r1,#0 // check sign of y bgt 4f b 1f -cordic_rstep: @ one step of algorithm in rotation mode - cmp r2,#0 @ check sign of angle +cordic_rstep: // one step of algorithm in rotation mode + cmp r2,#0 // check sign of angle bge 1f 4: - subs r1,r6 @ negative rotation: y=y-(x>>i) + subs r1,r6 // negative rotation: y=y-(x>>i) rsbs r7,#0 - adds r2,r4 @ accumulate angle + adds r2,r4 // accumulate angle b 2f 1: - adds r1,r6 @ positive rotation: y=y+(x>>i) - subs r2,r4 @ accumulate angle + adds r1,r6 // positive rotation: y=y+(x>>i) + subs r2,r4 // accumulate angle 2: mov r4,r12 - muls r7,r4 @ apply sign from m - subs r0,r7 @ finish rotation: x=x{+/-}(y>>i) + muls r7,r4 // apply sign from m + subs r0,r7 // finish rotation: x=x{+/-}(y>>i) 5: - ldmia r3!,{r4} @ fetch next angle from table and bump pointer - lsrs r4,#1 @ repeated angle? + ldmia r3!,{r4} // fetch next angle from table and bump pointer + lsrs r4,#1 // repeated angle? bcs 3f - adds r5,#1 @ adjust shift if not + adds r5,#1 // adjust shift if not 3: mov r6,r0 - asrs r6,r5 @ x>>i + asrs r6,r5 // x>>i mov r7,r1 - asrs r7,r5 @ y>>i - lsrs r4,#1 @ shift end flag into carry + asrs r7,r5 // y>>i + lsrs r4,#1 // shift end flag into carry bx r14 -@ CORDIC rotation mode +// CORDIC rotation mode cordic_rot: push {r6,r7,r14} - bl cordic_start @ initialise + bl cordic_start // initialise 1: bl cordic_rstep - bcc 1b @ step until table finished - asrs r6,r0,#14 @ remaining small rotations can be linearised: see IV.B of paper referenced above + bcc 1b // step until table finished + asrs r6,r0,#14 // remaining small rotations can be linearised: see IV.B of paper referenced above asrs r7,r1,#14 asrs r2,#3 - muls r6,r2 @ all remaining CORDIC steps in a multiplication + muls r6,r2 // all remaining CORDIC steps in a multiplication muls r7,r2 mov r4,r12 muls r7,r4 asrs r6,#12 asrs r7,#12 - subs r0,r7 @ x=x{+/-}(yz>>k) - adds r1,r6 @ y=y+(xz>>k) + subs r0,r7 // x=x{+/-}(yz>>k) + adds r1,r6 // y=y+(xz>>k) cordic_exit: pop {r6,r7,r15} -@ CORDIC vector mode +// CORDIC vector mode cordic_vec: push {r6,r7,r14} - bl cordic_start @ initialise + bl cordic_start // initialise 1: bl cordic_vstep - bcc 1b @ step until table finished + bcc 1b // step until table finished 4: - cmp r1,#0 @ continue as in cordic_vstep but without using table; x is not affected as y is small - bgt 2f @ check sign of y - adds r1,r6 @ positive rotation: y=y+(x>>i) - subs r2,r4 @ accumulate angle + cmp r1,#0 // continue as in cordic_vstep but without using table; x is not affected as y is small + bgt 2f // check sign of y + adds r1,r6 // positive rotation: y=y+(x>>i) + subs r2,r4 // accumulate angle b 3f 2: - subs r1,r6 @ negative rotation: y=y-(x>>i) - adds r2,r4 @ accumulate angle + subs r1,r6 // negative rotation: y=y-(x>>i) + adds r2,r4 // accumulate angle 3: asrs r6,#1 - asrs r4,#1 @ next "table entry" + asrs r4,#1 // next "table entry" bne 4b b cordic_exit .thumb_func -qfp_fsin: @ calculate sin and cos using CORDIC rotation method +qfp_fsin: // calculate sin and cos using CORDIC rotation method push {r4,r5,r14} movs r1,#24 - bl qfp_float2fix @ range reduction by repeated subtraction/addition in fixed point + bl qfp_float2fix // range reduction by repeated subtraction/addition in fixed point ldr r4,pi_q29 - lsrs r4,#4 @ 2pi Q24 + lsrs r4,#4 // 2pi Q24 1: subs r0,r4 bge 1b 1: adds r0,r4 - bmi 1b @ now in range 0..2pi - lsls r2,r0,#2 @ z Q26 - lsls r5,r4,#1 @ pi Q26 (r4=pi/2 Q26) - ldr r0,=#0x136e9db4 @ initialise CORDIC x,y with scaling + bmi 1b // now in range 0..2pi + lsls r2,r0,#2 // z Q26 + lsls r5,r4,#1 // pi Q26 (r4=pi/2 Q26) + ldr r0,=#0x136e9db4 // initialise CORDIC x,y with scaling movs r1,#0 1: - cmp r2,r4 @ >pi/2? + cmp r2,r4 // >pi/2? blt 2f - subs r2,r5 @ reduce range to -pi/2..pi/2 - rsbs r0,#0 @ rotate vector by pi + subs r2,r5 // reduce range to -pi/2..pi/2 + rsbs r0,#0 // rotate vector by pi b 1b 2: - lsls r2,#3 @ Q29 - adr r3,tab_cc @ circular coefficients - movs r4,#1 @ m=1 + lsls r2,#3 // Q29 + adr r3,tab_cc // circular coefficients + movs r4,#1 // m=1 bl cordic_rot - adds r1,#9 @ fiddle factor to make sin(0)==0 - movs r2,#0 @ exponents to zero + adds r1,#9 // fiddle factor to make sin(0)==0 + movs r2,#0 // exponents to zero movs r3,#0 - movs r5,#0 @ no sticky bits - bl packx @ pack cosine + movs r5,#0 // no sticky bits + bl packx // pack cosine bl xchxy - b packretns @ pack sine + b packretns // pack sine .thumb_func qfp_fcos: push {r14} bl qfp_fsin - mov r0,r1 @ extract cosine result + mov r0,r1 // extract cosine result pop {r15} .thumb_func qfp_ftan: push {r4,r5,r14} - bl qfp_fsin @ sine in r0/r2, cosine in r1/r3 + bl qfp_fsin // sine in r0/r2, cosine in r1/r3 .if include_faster - b fdiv_fast_n @ sin/cos + b fdiv_fast_n // sin/cos .else b fdiv_n .endif .thumb_func -qfp_fexp: @ calculate cosh and sinh using rotation method; add to obtain exp +qfp_fexp: // calculate cosh and sinh using rotation method; add to obtain exp push {r4,r5,r14} movs r1,#24 - bl qfp_float2fix @ Q24: covers entire valid input range - asrs r1,r0,#16 @ Q8 - ldr r2,=#5909 @ log_2(e) Q12 - muls r1,r2 @ estimate exponent of result Q20 - asrs r1,#19 @ Q1 - adds r1,#1 @ rounding - asrs r1,#1 @ rounded estimate of exponent of result - push {r1} @ save for later - lsls r2,r0,#5 @ Q29 - ldr r0,=#0x162e42ff @ ln(2) Q29 - muls r1,r0 @ accurate contribution of estimated exponent - subs r2,r1 @ residual to be exponentiated, approximately -.5..+.5 Q29 - ldr r0,=#0x2c9e15ca @ initialise CORDIC x,y with scaling + bl qfp_float2fix // Q24: covers entire valid input range + asrs r1,r0,#16 // Q8 + ldr r2,=#5909 // log_2(e) Q12 + muls r1,r2 // estimate exponent of result Q20 + asrs r1,#19 // Q1 + adds r1,#1 // rounding + asrs r1,#1 // rounded estimate of exponent of result + push {r1} // save for later + lsls r2,r0,#5 // Q29 + ldr r0,=#0x162e42ff // ln(2) Q29 + muls r1,r0 // accurate contribution of estimated exponent + subs r2,r1 // residual to be exponentiated, approximately -.5..+.5 Q29 + ldr r0,=#0x2c9e15ca // initialise CORDIC x,y with scaling movs r1,#0 - adr r3,tab_ch @ hyperbolic coefficients - mvns r4,r1 @ m=-1 - bl cordic_rot @ calculate cosh and sinh - add r0,r1 @ exp=cosh+sinh - pop {r2} @ recover exponent - b packretns @ pack result + adr r3,tab_ch // hyperbolic coefficients + mvns r4,r1 // m=-1 + bl cordic_rot // calculate cosh and sinh + add r0,r1 // exp=cosh+sinh + pop {r2} // recover exponent + b packretns // pack result .thumb_func -qfp_fsqrt: @ calculate sqrt and ln using vector method +qfp_fsqrt: // calculate sqrt and ln using vector method push {r4,r5,r14} bl unpackx - movs r1,r0 @ -ve argument? - bmi 3f @ return -Inf, -Inf - ldr r1,=#0x0593C2B9 @ scale factor for CORDIC - bl mul0 @ Q29 - asrs r1,r2,#1 @ halve exponent + movs r1,r0 // -ve argument? + bmi 3f // return -Inf, -Inf + ldr r1,=#0x0593C2B9 // scale factor for CORDIC + bl mul0 // Q29 + asrs r1,r2,#1 // halve exponent bcc 1f - adds r1,#1 @ was odd: add 1 and shift mantissa + adds r1,#1 // was odd: add 1 and shift mantissa asrs r0,#1 1: - push {r1} @ save exponent/2 for later + push {r1} // save exponent/2 for later mov r1,r0 - ldr r3,=#0x0593C2B9 @ re-use constant + ldr r3,=#0x0593C2B9 // re-use constant lsls r3,#2 - adds r0,r3 @ "a+1" - subs r1,r3 @ "a-1" + adds r0,r3 // "a+1" + subs r1,r3 // "a-1" movs r2,#0 - adr r3,tab_ch @ hyperbolic coefficients - mvns r4,r2 @ m=-1 + adr r3,tab_ch // hyperbolic coefficients + mvns r4,r2 // m=-1 bl cordic_vec - mov r1,r2 @ keep ln result - pop {r2} @ retrieve exponent/2 + mov r1,r2 // keep ln result + pop {r2} // retrieve exponent/2 2: movs r3,r2 - b packretns @ pack sqrt result + b packretns // pack sqrt result 3: movs r2,#255 @@ -766,59 +743,59 @@ qfp_fsqrt: @ calculate sqrt and ln using vector method .thumb_func qfp_fln: push {r4,r5,r14} - bl qfp_fsqrt @ get unpacked ln in r1/r3; exponent has been halved - cmp r3,#70 @ ln(Inf)? - bgt 3b @ return Inf + bl qfp_fsqrt // get unpacked ln in r1/r3; exponent has been halved + cmp r3,#70 // ln(Inf)? + bgt 3b // return Inf rsbs r3,#0 cmp r3,#70 - bgt 1f @ ln(0)? return -Inf - ldr r0,=#0x0162e430 @ ln(4) Q24 - muls r0,r3 @ contribution from negated, halved exponent - adds r1,#8 @ round result of ln - asrs r1,#4 @ Q24 - subs r0,r1,r0 @ add in contribution from (negated) exponent - movs r2,#5 @ pack expects Q29 + bgt 1f // ln(0)? return -Inf + ldr r0,=#0x0162e430 // ln(4) Q24 + muls r0,r3 // contribution from negated, halved exponent + adds r1,#8 // round result of ln + asrs r1,#4 // Q24 + subs r0,r1,r0 // add in contribution from (negated) exponent + movs r2,#5 // pack expects Q29 b packretns 1: - mvns r0,r0 @ make result -Inf + mvns r0,r0 // make result -Inf b 3b .thumb_func qfp_fatan2: push {r4,r5,r14} - bl unpackxyalign @ convert to fixed point (ensure common exponent, which is discarded) + bl unpackxyalign // convert to fixed point (ensure common exponent, which is discarded) asrs r0,#1 asrs r1,#1 - movs r2,#0 @ initial angle - cmp r0,#0 @ x negative + movs r2,#0 // initial angle + cmp r0,#0 // x negative bge 5f - rsbs r0,#0 @ rotate to 1st/4th quadrants + rsbs r0,#0 // rotate to 1st/4th quadrants rsbs r1,#0 - ldr r2,pi_q29 @ pi Q29 + ldr r2,pi_q29 // pi Q29 5: - adr r3,tab_cc @ circular coefficients - movs r4,#1 @ m=1 - bl cordic_vec @ also produces magnitude (with scaling factor 1.646760119), which is discarded - mov r0,r2 @ result here is -pi/2..3pi/2 Q29 - ldr r2,pi_q29 @ pi Q29 - adds r4,r0,r2 @ attempt to fix -3pi/2..-pi case - bcs 6f @ -pi/2..0? leave result as is - subs r4,r0,r2 @ <pi? leave as is + adr r3,tab_cc // circular coefficients + movs r4,#1 // m=1 + bl cordic_vec // also produces magnitude (with scaling factor 1.646760119), which is discarded + mov r0,r2 // result here is -pi/2..3pi/2 Q29 + ldr r2,pi_q29 // pi Q29 + adds r4,r0,r2 // attempt to fix -3pi/2..-pi case + bcs 6f // -pi/2..0? leave result as is + subs r4,r0,r2 // <pi? leave as is bmi 6f - subs r0,r4,r2 @ >pi: take off 2pi + subs r0,r4,r2 // >pi: take off 2pi 6: - subs r0,#1 @ fiddle factor so atan2(0,1)==0 - movs r2,#0 @ exponent for pack + subs r0,#1 // fiddle factor so atan2(0,1)==0 + movs r2,#0 // exponent for pack b packretns .align 2 .ltorg -@ first entry in following table is pi Q29 +// first entry in following table is pi Q29 pi_q29: -@ circular CORDIC coefficients: atan(2^-i), b0=flag for preventing shift, b1=flag for end of table +// circular CORDIC coefficients: atan(2^-i), b0=flag for preventing shift, b1=flag for end of table tab_cc: -.word 0x1921fb54*4+1 @ no shift before first iteration +.word 0x1921fb54*4+1 // no shift before first iteration .word 0x0ed63383*4+0 .word 0x07d6dd7e*4+0 .word 0x03fab753*4+0 @@ -830,16 +807,16 @@ tab_cc: .word 0x000fffff*4+0 .word 0x0007ffff*4+0 .word 0x00040000*4+0 -.word 0x00020000*4+0+2 @ +2 marks end +.word 0x00020000*4+0+2 // +2 marks end -@ hyperbolic CORDIC coefficients: atanh(2^-i), flags as above +// hyperbolic CORDIC coefficients: atanh(2^-i), flags as above tab_ch: .word 0x1193ea7b*4+0 -.word 0x1193ea7b*4+1 @ repeat i=1 +.word 0x1193ea7b*4+1 // repeat i=1 .word 0x082c577d*4+0 .word 0x04056247*4+0 .word 0x0200ab11*4+0 -.word 0x0200ab11*4+1 @ repeat i=4 +.word 0x0200ab11*4+1 // repeat i=4 .word 0x01001559*4+0 .word 0x008002ab*4+0 .word 0x00400055*4+0 @@ -848,7 +825,7 @@ tab_ch: .word 0x00080001*4+0 .word 0x00040000*4+0 .word 0x00020000*4+0 -.word 0x00020000*4+1+2 @ repeat i=12 +.word 0x00020000*4+1+2 // repeat i=12 .endif diff --git a/keymaps/onekey/stm32f411/rules.mk b/keymaps/onekey/stm32f411/rules.mk index c070040..8b12a3d 100644 --- a/keymaps/onekey/stm32f411/rules.mk +++ b/keymaps/onekey/stm32f411/rules.mk @@ -1,2 +1,3 @@ CONSOLE_ENABLE = yes -SRC += $(KEYMAP_PATH)/qfplib/functions.c + +SRC += $(KEYMAP_PATH)/qfplib/*.s diff --git a/util/link.sh b/util/link.sh index 2851345..308feda 100755 --- a/util/link.sh +++ b/util/link.sh @@ -1,9 +1,8 @@ #!/bin/zsh -GITDIR="$HOME/git" QMKME="$HOME/git/qmk_me" QMKFIRMWARE="$HOME/git/qmk_firmware" -find $QMKFIRMWARE -type l -delete +find $QMKFIRMWARE -type l -not -path "$QMKFIRMWARE/.venv/*" -delete rm $QMKFIRMWARE/.git/info/exclude ln -sv $QMKME/waffle $QMKFIRMWARE/users/waffle @@ -15,3 +14,16 @@ for dir in "keyboards"/*; do ln -sv $src $dst echo "$dir" >> $QMKFIRMWARE/.git/info/exclude done + +link() { + src="$QMKME/keymaps/$1" + dst="$QMKFIRMWARE/keyboards/$2/keymaps/waffle" + ln -sv $QMKME/$1 $QMKFIRMWARE/keyboards/$2/keymaps/waffle + echo "keyboards/$2/keymaps/waffle" >> $QMKFIRMWARE/.git/info/exclude +} + +link crkbd crkbd +link f60 bioi/f60 +link ferris ferris +link kyria splitkb/kyria +link microdox boardsource/microdox |