123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483 |
- /*
- * Copyright 2012-15 Advanced Micro Devices, Inc.
- *
- * Permission is hereby granted, free of charge, to any person obtaining a
- * copy of this software and associated documentation files (the "Software"),
- * to deal in the Software without restriction, including without limitation
- * the rights to use, copy, modify, merge, publish, distribute, sublicense,
- * and/or sell copies of the Software, and to permit persons to whom the
- * Software is furnished to do so, subject to the following conditions:
- *
- * The above copyright notice and this permission notice shall be included in
- * all copies or substantial portions of the Software.
- *
- * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
- * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
- * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
- * THE COPYRIGHT HOLDER(S) OR AUTHOR(S) BE LIABLE FOR ANY CLAIM, DAMAGES OR
- * OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE,
- * ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
- * OTHER DEALINGS IN THE SOFTWARE.
- *
- * Authors: AMD
- *
- */
- #include "dm_services.h"
- #include "include/fixed31_32.h"
- static inline unsigned long long abs_i64(
- long long arg)
- {
- if (arg > 0)
- return (unsigned long long)arg;
- else
- return (unsigned long long)(-arg);
- }
- /*
- * @brief
- * result = dividend / divisor
- * *remainder = dividend % divisor
- */
- static inline unsigned long long complete_integer_division_u64(
- unsigned long long dividend,
- unsigned long long divisor,
- unsigned long long *remainder)
- {
- unsigned long long result;
- ASSERT(divisor);
- result = div64_u64_rem(dividend, divisor, remainder);
- return result;
- }
- #define FRACTIONAL_PART_MASK \
- ((1ULL << FIXED31_32_BITS_PER_FRACTIONAL_PART) - 1)
- #define GET_INTEGER_PART(x) \
- ((x) >> FIXED31_32_BITS_PER_FRACTIONAL_PART)
- #define GET_FRACTIONAL_PART(x) \
- (FRACTIONAL_PART_MASK & (x))
- struct fixed31_32 dc_fixpt_from_fraction(long long numerator, long long denominator)
- {
- struct fixed31_32 res;
- bool arg1_negative = numerator < 0;
- bool arg2_negative = denominator < 0;
- unsigned long long arg1_value = arg1_negative ? -numerator : numerator;
- unsigned long long arg2_value = arg2_negative ? -denominator : denominator;
- unsigned long long remainder;
- /* determine integer part */
- unsigned long long res_value = complete_integer_division_u64(
- arg1_value, arg2_value, &remainder);
- ASSERT(res_value <= LONG_MAX);
- /* determine fractional part */
- {
- unsigned int i = FIXED31_32_BITS_PER_FRACTIONAL_PART;
- do {
- remainder <<= 1;
- res_value <<= 1;
- if (remainder >= arg2_value) {
- res_value |= 1;
- remainder -= arg2_value;
- }
- } while (--i != 0);
- }
- /* round up LSB */
- {
- unsigned long long summand = (remainder << 1) >= arg2_value;
- ASSERT(res_value <= LLONG_MAX - summand);
- res_value += summand;
- }
- res.value = (long long)res_value;
- if (arg1_negative ^ arg2_negative)
- res.value = -res.value;
- return res;
- }
- struct fixed31_32 dc_fixpt_mul(struct fixed31_32 arg1, struct fixed31_32 arg2)
- {
- struct fixed31_32 res;
- bool arg1_negative = arg1.value < 0;
- bool arg2_negative = arg2.value < 0;
- unsigned long long arg1_value = arg1_negative ? -arg1.value : arg1.value;
- unsigned long long arg2_value = arg2_negative ? -arg2.value : arg2.value;
- unsigned long long arg1_int = GET_INTEGER_PART(arg1_value);
- unsigned long long arg2_int = GET_INTEGER_PART(arg2_value);
- unsigned long long arg1_fra = GET_FRACTIONAL_PART(arg1_value);
- unsigned long long arg2_fra = GET_FRACTIONAL_PART(arg2_value);
- unsigned long long tmp;
- res.value = arg1_int * arg2_int;
- ASSERT(res.value <= LONG_MAX);
- res.value <<= FIXED31_32_BITS_PER_FRACTIONAL_PART;
- tmp = arg1_int * arg2_fra;
- ASSERT(tmp <= (unsigned long long)(LLONG_MAX - res.value));
- res.value += tmp;
- tmp = arg2_int * arg1_fra;
- ASSERT(tmp <= (unsigned long long)(LLONG_MAX - res.value));
- res.value += tmp;
- tmp = arg1_fra * arg2_fra;
- tmp = (tmp >> FIXED31_32_BITS_PER_FRACTIONAL_PART) +
- (tmp >= (unsigned long long)dc_fixpt_half.value);
- ASSERT(tmp <= (unsigned long long)(LLONG_MAX - res.value));
- res.value += tmp;
- if (arg1_negative ^ arg2_negative)
- res.value = -res.value;
- return res;
- }
- struct fixed31_32 dc_fixpt_sqr(struct fixed31_32 arg)
- {
- struct fixed31_32 res;
- unsigned long long arg_value = abs_i64(arg.value);
- unsigned long long arg_int = GET_INTEGER_PART(arg_value);
- unsigned long long arg_fra = GET_FRACTIONAL_PART(arg_value);
- unsigned long long tmp;
- res.value = arg_int * arg_int;
- ASSERT(res.value <= LONG_MAX);
- res.value <<= FIXED31_32_BITS_PER_FRACTIONAL_PART;
- tmp = arg_int * arg_fra;
- ASSERT(tmp <= (unsigned long long)(LLONG_MAX - res.value));
- res.value += tmp;
- ASSERT(tmp <= (unsigned long long)(LLONG_MAX - res.value));
- res.value += tmp;
- tmp = arg_fra * arg_fra;
- tmp = (tmp >> FIXED31_32_BITS_PER_FRACTIONAL_PART) +
- (tmp >= (unsigned long long)dc_fixpt_half.value);
- ASSERT(tmp <= (unsigned long long)(LLONG_MAX - res.value));
- res.value += tmp;
- return res;
- }
- struct fixed31_32 dc_fixpt_recip(struct fixed31_32 arg)
- {
- /*
- * @note
- * Good idea to use Newton's method
- */
- ASSERT(arg.value);
- return dc_fixpt_from_fraction(
- dc_fixpt_one.value,
- arg.value);
- }
- struct fixed31_32 dc_fixpt_sinc(struct fixed31_32 arg)
- {
- struct fixed31_32 square;
- struct fixed31_32 res = dc_fixpt_one;
- int n = 27;
- struct fixed31_32 arg_norm = arg;
- if (dc_fixpt_le(
- dc_fixpt_two_pi,
- dc_fixpt_abs(arg))) {
- arg_norm = dc_fixpt_sub(
- arg_norm,
- dc_fixpt_mul_int(
- dc_fixpt_two_pi,
- (int)div64_s64(
- arg_norm.value,
- dc_fixpt_two_pi.value)));
- }
- square = dc_fixpt_sqr(arg_norm);
- do {
- res = dc_fixpt_sub(
- dc_fixpt_one,
- dc_fixpt_div_int(
- dc_fixpt_mul(
- square,
- res),
- n * (n - 1)));
- n -= 2;
- } while (n > 2);
- if (arg.value != arg_norm.value)
- res = dc_fixpt_div(
- dc_fixpt_mul(res, arg_norm),
- arg);
- return res;
- }
- struct fixed31_32 dc_fixpt_sin(struct fixed31_32 arg)
- {
- return dc_fixpt_mul(
- arg,
- dc_fixpt_sinc(arg));
- }
- struct fixed31_32 dc_fixpt_cos(struct fixed31_32 arg)
- {
- /* TODO implement argument normalization */
- const struct fixed31_32 square = dc_fixpt_sqr(arg);
- struct fixed31_32 res = dc_fixpt_one;
- int n = 26;
- do {
- res = dc_fixpt_sub(
- dc_fixpt_one,
- dc_fixpt_div_int(
- dc_fixpt_mul(
- square,
- res),
- n * (n - 1)));
- n -= 2;
- } while (n != 0);
- return res;
- }
- /*
- * @brief
- * result = exp(arg),
- * where abs(arg) < 1
- *
- * Calculated as Taylor series.
- */
- static struct fixed31_32 fixed31_32_exp_from_taylor_series(struct fixed31_32 arg)
- {
- unsigned int n = 9;
- struct fixed31_32 res = dc_fixpt_from_fraction(
- n + 2,
- n + 1);
- /* TODO find correct res */
- ASSERT(dc_fixpt_lt(arg, dc_fixpt_one));
- do
- res = dc_fixpt_add(
- dc_fixpt_one,
- dc_fixpt_div_int(
- dc_fixpt_mul(
- arg,
- res),
- n));
- while (--n != 1);
- return dc_fixpt_add(
- dc_fixpt_one,
- dc_fixpt_mul(
- arg,
- res));
- }
- struct fixed31_32 dc_fixpt_exp(struct fixed31_32 arg)
- {
- /*
- * @brief
- * Main equation is:
- * exp(x) = exp(r + m * ln(2)) = (1 << m) * exp(r),
- * where m = round(x / ln(2)), r = x - m * ln(2)
- */
- if (dc_fixpt_le(
- dc_fixpt_ln2_div_2,
- dc_fixpt_abs(arg))) {
- int m = dc_fixpt_round(
- dc_fixpt_div(
- arg,
- dc_fixpt_ln2));
- struct fixed31_32 r = dc_fixpt_sub(
- arg,
- dc_fixpt_mul_int(
- dc_fixpt_ln2,
- m));
- ASSERT(m != 0);
- ASSERT(dc_fixpt_lt(
- dc_fixpt_abs(r),
- dc_fixpt_one));
- if (m > 0)
- return dc_fixpt_shl(
- fixed31_32_exp_from_taylor_series(r),
- (unsigned char)m);
- else
- return dc_fixpt_div_int(
- fixed31_32_exp_from_taylor_series(r),
- 1LL << -m);
- } else if (arg.value != 0)
- return fixed31_32_exp_from_taylor_series(arg);
- else
- return dc_fixpt_one;
- }
- struct fixed31_32 dc_fixpt_log(struct fixed31_32 arg)
- {
- struct fixed31_32 res = dc_fixpt_neg(dc_fixpt_one);
- /* TODO improve 1st estimation */
- struct fixed31_32 error;
- ASSERT(arg.value > 0);
- /* TODO if arg is negative, return NaN */
- /* TODO if arg is zero, return -INF */
- do {
- struct fixed31_32 res1 = dc_fixpt_add(
- dc_fixpt_sub(
- res,
- dc_fixpt_one),
- dc_fixpt_div(
- arg,
- dc_fixpt_exp(res)));
- error = dc_fixpt_sub(
- res,
- res1);
- res = res1;
- /* TODO determine max_allowed_error based on quality of exp() */
- } while (abs_i64(error.value) > 100ULL);
- return res;
- }
- /* this function is a generic helper to translate fixed point value to
- * specified integer format that will consist of integer_bits integer part and
- * fractional_bits fractional part. For example it is used in
- * dc_fixpt_u2d19 to receive 2 bits integer part and 19 bits fractional
- * part in 32 bits. It is used in hw programming (scaler)
- */
- static inline unsigned int ux_dy(
- long long value,
- unsigned int integer_bits,
- unsigned int fractional_bits)
- {
- /* 1. create mask of integer part */
- unsigned int result = (1 << integer_bits) - 1;
- /* 2. mask out fractional part */
- unsigned int fractional_part = FRACTIONAL_PART_MASK & value;
- /* 3. shrink fixed point integer part to be of integer_bits width*/
- result &= GET_INTEGER_PART(value);
- /* 4. make space for fractional part to be filled in after integer */
- result <<= fractional_bits;
- /* 5. shrink fixed point fractional part to of fractional_bits width*/
- fractional_part >>= FIXED31_32_BITS_PER_FRACTIONAL_PART - fractional_bits;
- /* 6. merge the result */
- return result | fractional_part;
- }
- static inline unsigned int clamp_ux_dy(
- long long value,
- unsigned int integer_bits,
- unsigned int fractional_bits,
- unsigned int min_clamp)
- {
- unsigned int truncated_val = ux_dy(value, integer_bits, fractional_bits);
- if (value >= (1LL << (integer_bits + FIXED31_32_BITS_PER_FRACTIONAL_PART)))
- return (1 << (integer_bits + fractional_bits)) - 1;
- else if (truncated_val > min_clamp)
- return truncated_val;
- else
- return min_clamp;
- }
- unsigned int dc_fixpt_u3d19(struct fixed31_32 arg)
- {
- return ux_dy(arg.value, 3, 19);
- }
- unsigned int dc_fixpt_u2d19(struct fixed31_32 arg)
- {
- return ux_dy(arg.value, 2, 19);
- }
- unsigned int dc_fixpt_u0d19(struct fixed31_32 arg)
- {
- return ux_dy(arg.value, 0, 19);
- }
- unsigned int dc_fixpt_clamp_u0d14(struct fixed31_32 arg)
- {
- return clamp_ux_dy(arg.value, 0, 14, 1);
- }
- unsigned int dc_fixpt_clamp_u0d10(struct fixed31_32 arg)
- {
- return clamp_ux_dy(arg.value, 0, 10, 1);
- }
- int dc_fixpt_s4d19(struct fixed31_32 arg)
- {
- if (arg.value < 0)
- return -(int)ux_dy(dc_fixpt_abs(arg).value, 4, 19);
- else
- return ux_dy(arg.value, 4, 19);
- }
|