mirror of
https://github.com/shopspring/decimal.git
synced 2024-11-23 21:00:49 +01:00
2436 lines
66 KiB
Go
2436 lines
66 KiB
Go
// Package decimal implements an arbitrary precision fixed-point decimal.
|
|
//
|
|
// The zero-value of a Decimal is 0, as you would expect.
|
|
//
|
|
// The best way to create a new Decimal is to use decimal.NewFromString, ex:
|
|
//
|
|
// n, err := decimal.NewFromString("-123.4567")
|
|
// n.String() // output: "-123.4567"
|
|
//
|
|
// To use Decimal as part of a struct:
|
|
//
|
|
// type StructName struct {
|
|
// Number Decimal
|
|
// }
|
|
//
|
|
// Note: This can "only" represent numbers with a maximum of 2^31 digits after the decimal point.
|
|
package decimal
|
|
|
|
import (
|
|
"database/sql/driver"
|
|
"encoding/binary"
|
|
"fmt"
|
|
"math"
|
|
"math/big"
|
|
"regexp"
|
|
"strconv"
|
|
"strings"
|
|
)
|
|
|
|
// DivisionPrecision is the number of decimal places in the result when it
|
|
// doesn't divide exactly.
|
|
//
|
|
// Example:
|
|
//
|
|
// d1 := decimal.NewFromFloat(2).Div(decimal.NewFromFloat(3))
|
|
// d1.String() // output: "0.6666666666666667"
|
|
// d2 := decimal.NewFromFloat(2).Div(decimal.NewFromFloat(30000))
|
|
// d2.String() // output: "0.0000666666666667"
|
|
// d3 := decimal.NewFromFloat(20000).Div(decimal.NewFromFloat(3))
|
|
// d3.String() // output: "6666.6666666666666667"
|
|
// decimal.DivisionPrecision = 3
|
|
// d4 := decimal.NewFromFloat(2).Div(decimal.NewFromFloat(3))
|
|
// d4.String() // output: "0.667"
|
|
var DivisionPrecision = 16
|
|
|
|
// PowPrecisionNegativeExponent specifies the maximum precision of the result (digits after decimal point)
|
|
// when calculating decimal power. Only used for cases where the exponent is a negative number.
|
|
// This constant applies to Pow, PowInt32 and PowBigInt methods, PowWithPrecision method is not constrained by it.
|
|
//
|
|
// Example:
|
|
//
|
|
// d1, err := decimal.NewFromFloat(15.2).PowInt32(-2)
|
|
// d1.String() // output: "0.0043282548476454"
|
|
//
|
|
// decimal.PowPrecisionNegativeExponent = 24
|
|
// d2, err := decimal.NewFromFloat(15.2).PowInt32(-2)
|
|
// d2.String() // output: "0.004328254847645429362881"
|
|
var PowPrecisionNegativeExponent = 16
|
|
|
|
// MarshalJSONWithoutQuotes should be set to true if you want the decimal to
|
|
// be JSON marshaled as a number, instead of as a string.
|
|
// WARNING: this is dangerous for decimals with many digits, since many JSON
|
|
// unmarshallers (ex: Javascript's) will unmarshal JSON numbers to IEEE 754
|
|
// double-precision floating point numbers, which means you can potentially
|
|
// silently lose precision.
|
|
var MarshalJSONWithoutQuotes = false
|
|
|
|
// ExpMaxIterations specifies the maximum number of iterations needed to calculate
|
|
// precise natural exponent value using ExpHullAbrham method.
|
|
var ExpMaxIterations = 1000
|
|
|
|
// Zero constant, to make computations faster.
|
|
// Zero should never be compared with == or != directly, please use decimal.Equal or decimal.Cmp instead.
|
|
var Zero = New(0, 1)
|
|
|
|
var zeroInt = big.NewInt(0)
|
|
var oneInt = big.NewInt(1)
|
|
var twoInt = big.NewInt(2)
|
|
var fourInt = big.NewInt(4)
|
|
var fiveInt = big.NewInt(5)
|
|
var tenInt = big.NewInt(10)
|
|
var twentyInt = big.NewInt(20)
|
|
|
|
var factorials = []Decimal{New(1, 0)}
|
|
|
|
// Decimal represents a fixed-point decimal. It is immutable.
|
|
// number = value * 10 ^ exp
|
|
type Decimal struct {
|
|
value *big.Int
|
|
|
|
// NOTE(vadim): this must be an int32, because we cast it to float64 during
|
|
// calculations. If exp is 64 bit, we might lose precision.
|
|
// If we cared about being able to represent every possible decimal, we
|
|
// could make exp a *big.Int but it would hurt performance and numbers
|
|
// like that are unrealistic.
|
|
exp int32
|
|
}
|
|
|
|
// New returns a new fixed-point decimal, value * 10 ^ exp.
|
|
func New(value int64, exp int32) Decimal {
|
|
return Decimal{
|
|
value: big.NewInt(value),
|
|
exp: exp,
|
|
}
|
|
}
|
|
|
|
// NewFromInt converts an int64 to Decimal.
|
|
//
|
|
// Example:
|
|
//
|
|
// NewFromInt(123).String() // output: "123"
|
|
// NewFromInt(-10).String() // output: "-10"
|
|
func NewFromInt(value int64) Decimal {
|
|
return Decimal{
|
|
value: big.NewInt(value),
|
|
exp: 0,
|
|
}
|
|
}
|
|
|
|
// NewFromInt32 converts an int32 to Decimal.
|
|
//
|
|
// Example:
|
|
//
|
|
// NewFromInt(123).String() // output: "123"
|
|
// NewFromInt(-10).String() // output: "-10"
|
|
func NewFromInt32(value int32) Decimal {
|
|
return Decimal{
|
|
value: big.NewInt(int64(value)),
|
|
exp: 0,
|
|
}
|
|
}
|
|
|
|
// NewFromUint64 converts an uint64 to Decimal.
|
|
//
|
|
// Example:
|
|
//
|
|
// NewFromUint64(123).String() // output: "123"
|
|
func NewFromUint64(value uint64) Decimal {
|
|
return Decimal{
|
|
value: new(big.Int).SetUint64(value),
|
|
exp: 0,
|
|
}
|
|
}
|
|
|
|
// NewFromBigInt returns a new Decimal from a big.Int, value * 10 ^ exp
|
|
func NewFromBigInt(value *big.Int, exp int32) Decimal {
|
|
return Decimal{
|
|
value: new(big.Int).Set(value),
|
|
exp: exp,
|
|
}
|
|
}
|
|
|
|
// NewFromBigRat returns a new Decimal from a big.Rat. The numerator and
|
|
// denominator are divided and rounded to the given precision.
|
|
//
|
|
// Example:
|
|
//
|
|
// d1 := NewFromBigRat(big.NewRat(0, 1), 0) // output: "0"
|
|
// d2 := NewFromBigRat(big.NewRat(4, 5), 1) // output: "0.8"
|
|
// d3 := NewFromBigRat(big.NewRat(1000, 3), 3) // output: "333.333"
|
|
// d4 := NewFromBigRat(big.NewRat(2, 7), 4) // output: "0.2857"
|
|
func NewFromBigRat(value *big.Rat, precision int32) Decimal {
|
|
return Decimal{
|
|
value: new(big.Int).Set(value.Num()),
|
|
exp: 0,
|
|
}.DivRound(Decimal{
|
|
value: new(big.Int).Set(value.Denom()),
|
|
exp: 0,
|
|
}, precision)
|
|
}
|
|
|
|
// NewFromString returns a new Decimal from a string representation.
|
|
// Trailing zeroes are not trimmed.
|
|
//
|
|
// Example:
|
|
//
|
|
// d, err := NewFromString("-123.45")
|
|
// d2, err := NewFromString(".0001")
|
|
// d3, err := NewFromString("1.47000")
|
|
func NewFromString(value string) (Decimal, error) {
|
|
originalInput := value
|
|
var intString string
|
|
var exp int64
|
|
|
|
// Check if number is using scientific notation
|
|
eIndex := strings.IndexAny(value, "Ee")
|
|
if eIndex != -1 {
|
|
expInt, err := strconv.ParseInt(value[eIndex+1:], 10, 32)
|
|
if err != nil {
|
|
if e, ok := err.(*strconv.NumError); ok && e.Err == strconv.ErrRange {
|
|
return Decimal{}, fmt.Errorf("can't convert %s to decimal: fractional part too long", value)
|
|
}
|
|
return Decimal{}, fmt.Errorf("can't convert %s to decimal: exponent is not numeric", value)
|
|
}
|
|
value = value[:eIndex]
|
|
exp = expInt
|
|
}
|
|
|
|
pIndex := -1
|
|
vLen := len(value)
|
|
for i := 0; i < vLen; i++ {
|
|
if value[i] == '.' {
|
|
if pIndex > -1 {
|
|
return Decimal{}, fmt.Errorf("can't convert %s to decimal: too many .s", value)
|
|
}
|
|
pIndex = i
|
|
}
|
|
}
|
|
|
|
if pIndex == -1 {
|
|
// There is no decimal point, we can just parse the original string as
|
|
// an int
|
|
intString = value
|
|
} else {
|
|
if pIndex+1 < vLen {
|
|
intString = value[:pIndex] + value[pIndex+1:]
|
|
} else {
|
|
intString = value[:pIndex]
|
|
}
|
|
expInt := -len(value[pIndex+1:])
|
|
exp += int64(expInt)
|
|
}
|
|
|
|
var dValue *big.Int
|
|
// strconv.ParseInt is faster than new(big.Int).SetString so this is just a shortcut for strings we know won't overflow
|
|
if len(intString) <= 18 {
|
|
parsed64, err := strconv.ParseInt(intString, 10, 64)
|
|
if err != nil {
|
|
return Decimal{}, fmt.Errorf("can't convert %s to decimal", value)
|
|
}
|
|
dValue = big.NewInt(parsed64)
|
|
} else {
|
|
dValue = new(big.Int)
|
|
_, ok := dValue.SetString(intString, 10)
|
|
if !ok {
|
|
return Decimal{}, fmt.Errorf("can't convert %s to decimal", value)
|
|
}
|
|
}
|
|
|
|
if exp < math.MinInt32 || exp > math.MaxInt32 {
|
|
// NOTE(vadim): I doubt a string could realistically be this long
|
|
return Decimal{}, fmt.Errorf("can't convert %s to decimal: fractional part too long", originalInput)
|
|
}
|
|
|
|
return Decimal{
|
|
value: dValue,
|
|
exp: int32(exp),
|
|
}, nil
|
|
}
|
|
|
|
// NewFromFormattedString returns a new Decimal from a formatted string representation.
|
|
// The second argument - replRegexp, is a regular expression that is used to find characters that should be
|
|
// removed from given decimal string representation. All matched characters will be replaced with an empty string.
|
|
//
|
|
// Example:
|
|
//
|
|
// r := regexp.MustCompile("[$,]")
|
|
// d1, err := NewFromFormattedString("$5,125.99", r)
|
|
//
|
|
// r2 := regexp.MustCompile("[_]")
|
|
// d2, err := NewFromFormattedString("1_000_000", r2)
|
|
//
|
|
// r3 := regexp.MustCompile("[USD\\s]")
|
|
// d3, err := NewFromFormattedString("5000 USD", r3)
|
|
func NewFromFormattedString(value string, replRegexp *regexp.Regexp) (Decimal, error) {
|
|
parsedValue := replRegexp.ReplaceAllString(value, "")
|
|
d, err := NewFromString(parsedValue)
|
|
if err != nil {
|
|
return Decimal{}, err
|
|
}
|
|
return d, nil
|
|
}
|
|
|
|
// RequireFromString returns a new Decimal from a string representation
|
|
// or panics if NewFromString had returned an error.
|
|
//
|
|
// Example:
|
|
//
|
|
// d := RequireFromString("-123.45")
|
|
// d2 := RequireFromString(".0001")
|
|
func RequireFromString(value string) Decimal {
|
|
dec, err := NewFromString(value)
|
|
if err != nil {
|
|
panic(err)
|
|
}
|
|
return dec
|
|
}
|
|
|
|
// NewFromFloat converts a float64 to Decimal.
|
|
//
|
|
// The converted number will contain the number of significant digits that can be
|
|
// represented in a float with reliable roundtrip.
|
|
// This is typically 15 digits, but may be more in some cases.
|
|
// See https://www.exploringbinary.com/decimal-precision-of-binary-floating-point-numbers/ for more information.
|
|
//
|
|
// For slightly faster conversion, use NewFromFloatWithExponent where you can specify the precision in absolute terms.
|
|
//
|
|
// NOTE: this will panic on NaN, +/-inf
|
|
func NewFromFloat(value float64) Decimal {
|
|
if value == 0 {
|
|
return New(0, 0)
|
|
}
|
|
return newFromFloat(value, math.Float64bits(value), &float64info)
|
|
}
|
|
|
|
// NewFromFloat32 converts a float32 to Decimal.
|
|
//
|
|
// The converted number will contain the number of significant digits that can be
|
|
// represented in a float with reliable roundtrip.
|
|
// This is typically 6-8 digits depending on the input.
|
|
// See https://www.exploringbinary.com/decimal-precision-of-binary-floating-point-numbers/ for more information.
|
|
//
|
|
// For slightly faster conversion, use NewFromFloatWithExponent where you can specify the precision in absolute terms.
|
|
//
|
|
// NOTE: this will panic on NaN, +/-inf
|
|
func NewFromFloat32(value float32) Decimal {
|
|
if value == 0 {
|
|
return New(0, 0)
|
|
}
|
|
// XOR is workaround for https://github.com/golang/go/issues/26285
|
|
a := math.Float32bits(value) ^ 0x80808080
|
|
return newFromFloat(float64(value), uint64(a)^0x80808080, &float32info)
|
|
}
|
|
|
|
func newFromFloat(val float64, bits uint64, flt *floatInfo) Decimal {
|
|
if math.IsNaN(val) || math.IsInf(val, 0) {
|
|
panic(fmt.Sprintf("Cannot create a Decimal from %v", val))
|
|
}
|
|
exp := int(bits>>flt.mantbits) & (1<<flt.expbits - 1)
|
|
mant := bits & (uint64(1)<<flt.mantbits - 1)
|
|
|
|
switch exp {
|
|
case 0:
|
|
// denormalized
|
|
exp++
|
|
|
|
default:
|
|
// add implicit top bit
|
|
mant |= uint64(1) << flt.mantbits
|
|
}
|
|
exp += flt.bias
|
|
|
|
var d decimal
|
|
d.Assign(mant)
|
|
d.Shift(exp - int(flt.mantbits))
|
|
d.neg = bits>>(flt.expbits+flt.mantbits) != 0
|
|
|
|
roundShortest(&d, mant, exp, flt)
|
|
// If less than 19 digits, we can do calculation in an int64.
|
|
if d.nd < 19 {
|
|
tmp := int64(0)
|
|
m := int64(1)
|
|
for i := d.nd - 1; i >= 0; i-- {
|
|
tmp += m * int64(d.d[i]-'0')
|
|
m *= 10
|
|
}
|
|
if d.neg {
|
|
tmp *= -1
|
|
}
|
|
return Decimal{value: big.NewInt(tmp), exp: int32(d.dp) - int32(d.nd)}
|
|
}
|
|
dValue := new(big.Int)
|
|
dValue, ok := dValue.SetString(string(d.d[:d.nd]), 10)
|
|
if ok {
|
|
return Decimal{value: dValue, exp: int32(d.dp) - int32(d.nd)}
|
|
}
|
|
|
|
return NewFromFloatWithExponent(val, int32(d.dp)-int32(d.nd))
|
|
}
|
|
|
|
// NewFromFloatWithExponent converts a float64 to Decimal, with an arbitrary
|
|
// number of fractional digits.
|
|
//
|
|
// Example:
|
|
//
|
|
// NewFromFloatWithExponent(123.456, -2).String() // output: "123.46"
|
|
func NewFromFloatWithExponent(value float64, exp int32) Decimal {
|
|
if math.IsNaN(value) || math.IsInf(value, 0) {
|
|
panic(fmt.Sprintf("Cannot create a Decimal from %v", value))
|
|
}
|
|
|
|
bits := math.Float64bits(value)
|
|
mant := bits & (1<<52 - 1)
|
|
exp2 := int32((bits >> 52) & (1<<11 - 1))
|
|
sign := bits >> 63
|
|
|
|
if exp2 == 0 {
|
|
// specials
|
|
if mant == 0 {
|
|
return Decimal{}
|
|
}
|
|
// subnormal
|
|
exp2++
|
|
} else {
|
|
// normal
|
|
mant |= 1 << 52
|
|
}
|
|
|
|
exp2 -= 1023 + 52
|
|
|
|
// normalizing base-2 values
|
|
for mant&1 == 0 {
|
|
mant = mant >> 1
|
|
exp2++
|
|
}
|
|
|
|
// maximum number of fractional base-10 digits to represent 2^N exactly cannot be more than -N if N<0
|
|
if exp < 0 && exp < exp2 {
|
|
if exp2 < 0 {
|
|
exp = exp2
|
|
} else {
|
|
exp = 0
|
|
}
|
|
}
|
|
|
|
// representing 10^M * 2^N as 5^M * 2^(M+N)
|
|
exp2 -= exp
|
|
|
|
temp := big.NewInt(1)
|
|
dMant := big.NewInt(int64(mant))
|
|
|
|
// applying 5^M
|
|
if exp > 0 {
|
|
temp = temp.SetInt64(int64(exp))
|
|
temp = temp.Exp(fiveInt, temp, nil)
|
|
} else if exp < 0 {
|
|
temp = temp.SetInt64(-int64(exp))
|
|
temp = temp.Exp(fiveInt, temp, nil)
|
|
dMant = dMant.Mul(dMant, temp)
|
|
temp = temp.SetUint64(1)
|
|
}
|
|
|
|
// applying 2^(M+N)
|
|
if exp2 > 0 {
|
|
dMant = dMant.Lsh(dMant, uint(exp2))
|
|
} else if exp2 < 0 {
|
|
temp = temp.Lsh(temp, uint(-exp2))
|
|
}
|
|
|
|
// rounding and downscaling
|
|
if exp > 0 || exp2 < 0 {
|
|
halfDown := new(big.Int).Rsh(temp, 1)
|
|
dMant = dMant.Add(dMant, halfDown)
|
|
dMant = dMant.Quo(dMant, temp)
|
|
}
|
|
|
|
if sign == 1 {
|
|
dMant = dMant.Neg(dMant)
|
|
}
|
|
|
|
return Decimal{
|
|
value: dMant,
|
|
exp: exp,
|
|
}
|
|
}
|
|
|
|
// Copy returns a copy of decimal with the same value and exponent, but a different pointer to value.
|
|
func (d Decimal) Copy() Decimal {
|
|
d.ensureInitialized()
|
|
return Decimal{
|
|
value: new(big.Int).Set(d.value),
|
|
exp: d.exp,
|
|
}
|
|
}
|
|
|
|
// rescale returns a rescaled version of the decimal. Returned
|
|
// decimal may be less precise if the given exponent is bigger
|
|
// than the initial exponent of the Decimal.
|
|
// NOTE: this will truncate, NOT round
|
|
//
|
|
// Example:
|
|
//
|
|
// d := New(12345, -4)
|
|
// d2 := d.rescale(-1)
|
|
// d3 := d2.rescale(-4)
|
|
// println(d1)
|
|
// println(d2)
|
|
// println(d3)
|
|
//
|
|
// Output:
|
|
//
|
|
// 1.2345
|
|
// 1.2
|
|
// 1.2000
|
|
func (d Decimal) rescale(exp int32) Decimal {
|
|
d.ensureInitialized()
|
|
|
|
if d.exp == exp {
|
|
return Decimal{
|
|
new(big.Int).Set(d.value),
|
|
d.exp,
|
|
}
|
|
}
|
|
|
|
// NOTE(vadim): must convert exps to float64 before - to prevent overflow
|
|
diff := math.Abs(float64(exp) - float64(d.exp))
|
|
value := new(big.Int).Set(d.value)
|
|
|
|
expScale := new(big.Int).Exp(tenInt, big.NewInt(int64(diff)), nil)
|
|
if exp > d.exp {
|
|
value = value.Quo(value, expScale)
|
|
} else if exp < d.exp {
|
|
value = value.Mul(value, expScale)
|
|
}
|
|
|
|
return Decimal{
|
|
value: value,
|
|
exp: exp,
|
|
}
|
|
}
|
|
|
|
// Abs returns the absolute value of the decimal.
|
|
func (d Decimal) Abs() Decimal {
|
|
if !d.IsNegative() {
|
|
return d
|
|
}
|
|
d.ensureInitialized()
|
|
d2Value := new(big.Int).Abs(d.value)
|
|
return Decimal{
|
|
value: d2Value,
|
|
exp: d.exp,
|
|
}
|
|
}
|
|
|
|
// Add returns d + d2.
|
|
func (d Decimal) Add(d2 Decimal) Decimal {
|
|
rd, rd2 := RescalePair(d, d2)
|
|
|
|
d3Value := new(big.Int).Add(rd.value, rd2.value)
|
|
return Decimal{
|
|
value: d3Value,
|
|
exp: rd.exp,
|
|
}
|
|
}
|
|
|
|
// Sub returns d - d2.
|
|
func (d Decimal) Sub(d2 Decimal) Decimal {
|
|
rd, rd2 := RescalePair(d, d2)
|
|
|
|
d3Value := new(big.Int).Sub(rd.value, rd2.value)
|
|
return Decimal{
|
|
value: d3Value,
|
|
exp: rd.exp,
|
|
}
|
|
}
|
|
|
|
// Neg returns -d.
|
|
func (d Decimal) Neg() Decimal {
|
|
d.ensureInitialized()
|
|
val := new(big.Int).Neg(d.value)
|
|
return Decimal{
|
|
value: val,
|
|
exp: d.exp,
|
|
}
|
|
}
|
|
|
|
// Mul returns d * d2.
|
|
func (d Decimal) Mul(d2 Decimal) Decimal {
|
|
d.ensureInitialized()
|
|
d2.ensureInitialized()
|
|
|
|
expInt64 := int64(d.exp) + int64(d2.exp)
|
|
if expInt64 > math.MaxInt32 || expInt64 < math.MinInt32 {
|
|
// NOTE(vadim): better to panic than give incorrect results, as
|
|
// Decimals are usually used for money
|
|
panic(fmt.Sprintf("exponent %v overflows an int32!", expInt64))
|
|
}
|
|
|
|
d3Value := new(big.Int).Mul(d.value, d2.value)
|
|
return Decimal{
|
|
value: d3Value,
|
|
exp: int32(expInt64),
|
|
}
|
|
}
|
|
|
|
// Shift shifts the decimal in base 10.
|
|
// It shifts left when shift is positive and right if shift is negative.
|
|
// In simpler terms, the given value for shift is added to the exponent
|
|
// of the decimal.
|
|
func (d Decimal) Shift(shift int32) Decimal {
|
|
d.ensureInitialized()
|
|
return Decimal{
|
|
value: new(big.Int).Set(d.value),
|
|
exp: d.exp + shift,
|
|
}
|
|
}
|
|
|
|
// Div returns d / d2. If it doesn't divide exactly, the result will have
|
|
// DivisionPrecision digits after the decimal point.
|
|
func (d Decimal) Div(d2 Decimal) Decimal {
|
|
return d.DivRound(d2, int32(DivisionPrecision))
|
|
}
|
|
|
|
// QuoRem does division with remainder
|
|
// d.QuoRem(d2,precision) returns quotient q and remainder r such that
|
|
//
|
|
// d = d2 * q + r, q an integer multiple of 10^(-precision)
|
|
// 0 <= r < abs(d2) * 10 ^(-precision) if d>=0
|
|
// 0 >= r > -abs(d2) * 10 ^(-precision) if d<0
|
|
//
|
|
// Note that precision<0 is allowed as input.
|
|
func (d Decimal) QuoRem(d2 Decimal, precision int32) (Decimal, Decimal) {
|
|
d.ensureInitialized()
|
|
d2.ensureInitialized()
|
|
if d2.value.Sign() == 0 {
|
|
panic("decimal division by 0")
|
|
}
|
|
scale := -precision
|
|
e := int64(d.exp) - int64(d2.exp) - int64(scale)
|
|
if e > math.MaxInt32 || e < math.MinInt32 {
|
|
panic("overflow in decimal QuoRem")
|
|
}
|
|
var aa, bb, expo big.Int
|
|
var scalerest int32
|
|
// d = a 10^ea
|
|
// d2 = b 10^eb
|
|
if e < 0 {
|
|
aa = *d.value
|
|
expo.SetInt64(-e)
|
|
bb.Exp(tenInt, &expo, nil)
|
|
bb.Mul(d2.value, &bb)
|
|
scalerest = d.exp
|
|
// now aa = a
|
|
// bb = b 10^(scale + eb - ea)
|
|
} else {
|
|
expo.SetInt64(e)
|
|
aa.Exp(tenInt, &expo, nil)
|
|
aa.Mul(d.value, &aa)
|
|
bb = *d2.value
|
|
scalerest = scale + d2.exp
|
|
// now aa = a ^ (ea - eb - scale)
|
|
// bb = b
|
|
}
|
|
var q, r big.Int
|
|
q.QuoRem(&aa, &bb, &r)
|
|
dq := Decimal{value: &q, exp: scale}
|
|
dr := Decimal{value: &r, exp: scalerest}
|
|
return dq, dr
|
|
}
|
|
|
|
// DivRound divides and rounds to a given precision
|
|
// i.e. to an integer multiple of 10^(-precision)
|
|
//
|
|
// for a positive quotient digit 5 is rounded up, away from 0
|
|
// if the quotient is negative then digit 5 is rounded down, away from 0
|
|
//
|
|
// Note that precision<0 is allowed as input.
|
|
func (d Decimal) DivRound(d2 Decimal, precision int32) Decimal {
|
|
// QuoRem already checks initialization
|
|
q, r := d.QuoRem(d2, precision)
|
|
// the actual rounding decision is based on comparing r*10^precision and d2/2
|
|
// instead compare 2 r 10 ^precision and d2
|
|
var rv2 big.Int
|
|
rv2.Abs(r.value)
|
|
rv2.Lsh(&rv2, 1)
|
|
// now rv2 = abs(r.value) * 2
|
|
r2 := Decimal{value: &rv2, exp: r.exp + precision}
|
|
// r2 is now 2 * r * 10 ^ precision
|
|
var c = r2.Cmp(d2.Abs())
|
|
|
|
if c < 0 {
|
|
return q
|
|
}
|
|
|
|
if d.value.Sign()*d2.value.Sign() < 0 {
|
|
return q.Sub(New(1, -precision))
|
|
}
|
|
|
|
return q.Add(New(1, -precision))
|
|
}
|
|
|
|
// Mod returns d % d2.
|
|
func (d Decimal) Mod(d2 Decimal) Decimal {
|
|
_, r := d.QuoRem(d2, 0)
|
|
return r
|
|
}
|
|
|
|
// Pow returns d to the power of d2.
|
|
// When exponent is negative the returned decimal will have maximum precision of PowPrecisionNegativeExponent places after decimal point.
|
|
//
|
|
// Pow returns 0 (zero-value of Decimal) instead of error for power operation edge cases, to handle those edge cases use PowWithPrecision
|
|
// Edge cases not handled by Pow:
|
|
// - 0 ** 0 => undefined value
|
|
// - 0 ** y, where y < 0 => infinity
|
|
// - x ** y, where x < 0 and y is non-integer decimal => imaginary value
|
|
//
|
|
// Example:
|
|
//
|
|
// d1 := decimal.NewFromFloat(4.0)
|
|
// d2 := decimal.NewFromFloat(4.0)
|
|
// res1 := d1.Pow(d2)
|
|
// res1.String() // output: "256"
|
|
//
|
|
// d3 := decimal.NewFromFloat(5.0)
|
|
// d4 := decimal.NewFromFloat(5.73)
|
|
// res2 := d3.Pow(d4)
|
|
// res2.String() // output: "10118.08037125"
|
|
func (d Decimal) Pow(d2 Decimal) Decimal {
|
|
baseSign := d.Sign()
|
|
expSign := d2.Sign()
|
|
|
|
if baseSign == 0 {
|
|
if expSign == 0 {
|
|
return Decimal{}
|
|
}
|
|
if expSign == 1 {
|
|
return Decimal{zeroInt, 0}
|
|
}
|
|
if expSign == -1 {
|
|
return Decimal{}
|
|
}
|
|
}
|
|
|
|
if expSign == 0 {
|
|
return Decimal{oneInt, 0}
|
|
}
|
|
|
|
// TODO: optimize extraction of fractional part
|
|
one := Decimal{oneInt, 0}
|
|
expIntPart, expFracPart := d2.QuoRem(one, 0)
|
|
|
|
if baseSign == -1 && !expFracPart.IsZero() {
|
|
return Decimal{}
|
|
}
|
|
|
|
intPartPow, _ := d.PowBigInt(expIntPart.value)
|
|
|
|
// if exponent is an integer we don't need to calculate d1**frac(d2)
|
|
if expFracPart.value.Sign() == 0 {
|
|
return intPartPow
|
|
}
|
|
|
|
// TODO: optimize NumDigits for more performant precision adjustment
|
|
digitsBase := d.NumDigits()
|
|
digitsExponent := d2.NumDigits()
|
|
|
|
precision := digitsBase
|
|
|
|
if digitsExponent > precision {
|
|
precision += digitsExponent
|
|
}
|
|
|
|
precision += 6
|
|
|
|
// Calculate x ** frac(y), where
|
|
// x ** frac(y) = exp(ln(x ** frac(y)) = exp(ln(x) * frac(y))
|
|
fracPartPow, err := d.Abs().Ln(-d.exp + int32(precision))
|
|
if err != nil {
|
|
return Decimal{}
|
|
}
|
|
|
|
fracPartPow = fracPartPow.Mul(expFracPart)
|
|
|
|
fracPartPow, err = fracPartPow.ExpTaylor(-d.exp + int32(precision))
|
|
if err != nil {
|
|
return Decimal{}
|
|
}
|
|
|
|
// Join integer and fractional part,
|
|
// base ** (expBase + expFrac) = base ** expBase * base ** expFrac
|
|
res := intPartPow.Mul(fracPartPow)
|
|
|
|
return res
|
|
}
|
|
|
|
// PowWithPrecision returns d to the power of d2.
|
|
// Precision parameter specifies minimum precision of the result (digits after decimal point).
|
|
// Returned decimal is not rounded to 'precision' places after decimal point.
|
|
//
|
|
// PowWithPrecision returns error when:
|
|
// - 0 ** 0 => undefined value
|
|
// - 0 ** y, where y < 0 => infinity
|
|
// - x ** y, where x < 0 and y is non-integer decimal => imaginary value
|
|
//
|
|
// Example:
|
|
//
|
|
// d1 := decimal.NewFromFloat(4.0)
|
|
// d2 := decimal.NewFromFloat(4.0)
|
|
// res1, err := d1.PowWithPrecision(d2, 2)
|
|
// res1.String() // output: "256"
|
|
//
|
|
// d3 := decimal.NewFromFloat(5.0)
|
|
// d4 := decimal.NewFromFloat(5.73)
|
|
// res2, err := d3.PowWithPrecision(d4, 5)
|
|
// res2.String() // output: "10118.080371595015625"
|
|
//
|
|
// d5 := decimal.NewFromFloat(-3.0)
|
|
// d6 := decimal.NewFromFloat(-6.0)
|
|
// res3, err := d5.PowWithPrecision(d6, 10)
|
|
// res3.String() // output: "0.0013717421"
|
|
func (d Decimal) PowWithPrecision(d2 Decimal, precision int32) (Decimal, error) {
|
|
baseSign := d.Sign()
|
|
expSign := d2.Sign()
|
|
|
|
if baseSign == 0 {
|
|
if expSign == 0 {
|
|
return Decimal{}, fmt.Errorf("cannot represent undefined value of 0**0")
|
|
}
|
|
if expSign == 1 {
|
|
return Decimal{zeroInt, 0}, nil
|
|
}
|
|
if expSign == -1 {
|
|
return Decimal{}, fmt.Errorf("cannot represent infinity value of 0 ** y, where y < 0")
|
|
}
|
|
}
|
|
|
|
if expSign == 0 {
|
|
return Decimal{oneInt, 0}, nil
|
|
}
|
|
|
|
// TODO: optimize extraction of fractional part
|
|
one := Decimal{oneInt, 0}
|
|
expIntPart, expFracPart := d2.QuoRem(one, 0)
|
|
|
|
if baseSign == -1 && !expFracPart.IsZero() {
|
|
return Decimal{}, fmt.Errorf("cannot represent imaginary value of x ** y, where x < 0 and y is non-integer decimal")
|
|
}
|
|
|
|
intPartPow, _ := d.powBigIntWithPrecision(expIntPart.value, precision)
|
|
|
|
// if exponent is an integer we don't need to calculate d1**frac(d2)
|
|
if expFracPart.value.Sign() == 0 {
|
|
return intPartPow, nil
|
|
}
|
|
|
|
// TODO: optimize NumDigits for more performant precision adjustment
|
|
digitsBase := d.NumDigits()
|
|
digitsExponent := d2.NumDigits()
|
|
|
|
if int32(digitsBase) > precision {
|
|
precision = int32(digitsBase)
|
|
}
|
|
if int32(digitsExponent) > precision {
|
|
precision += int32(digitsExponent)
|
|
}
|
|
// increase precision by 10 to compensate for errors in further calculations
|
|
precision += 10
|
|
|
|
// Calculate x ** frac(y), where
|
|
// x ** frac(y) = exp(ln(x ** frac(y)) = exp(ln(x) * frac(y))
|
|
fracPartPow, err := d.Abs().Ln(precision)
|
|
if err != nil {
|
|
return Decimal{}, err
|
|
}
|
|
|
|
fracPartPow = fracPartPow.Mul(expFracPart)
|
|
|
|
fracPartPow, err = fracPartPow.ExpTaylor(precision)
|
|
if err != nil {
|
|
return Decimal{}, err
|
|
}
|
|
|
|
// Join integer and fractional part,
|
|
// base ** (expBase + expFrac) = base ** expBase * base ** expFrac
|
|
res := intPartPow.Mul(fracPartPow)
|
|
|
|
return res, nil
|
|
}
|
|
|
|
// PowInt32 returns d to the power of exp, where exp is int32.
|
|
// Only returns error when d and exp is 0, thus result is undefined.
|
|
//
|
|
// When exponent is negative the returned decimal will have maximum precision of PowPrecisionNegativeExponent places after decimal point.
|
|
//
|
|
// Example:
|
|
//
|
|
// d1, err := decimal.NewFromFloat(4.0).PowInt32(4)
|
|
// d1.String() // output: "256"
|
|
//
|
|
// d2, err := decimal.NewFromFloat(3.13).PowInt32(5)
|
|
// d2.String() // output: "300.4150512793"
|
|
func (d Decimal) PowInt32(exp int32) (Decimal, error) {
|
|
if d.IsZero() && exp == 0 {
|
|
return Decimal{}, fmt.Errorf("cannot represent undefined value of 0**0")
|
|
}
|
|
|
|
isExpNeg := exp < 0
|
|
exp = abs(exp)
|
|
|
|
n, result := d, New(1, 0)
|
|
|
|
for exp > 0 {
|
|
if exp%2 == 1 {
|
|
result = result.Mul(n)
|
|
}
|
|
exp /= 2
|
|
|
|
if exp > 0 {
|
|
n = n.Mul(n)
|
|
}
|
|
}
|
|
|
|
if isExpNeg {
|
|
return New(1, 0).DivRound(result, int32(PowPrecisionNegativeExponent)), nil
|
|
}
|
|
|
|
return result, nil
|
|
}
|
|
|
|
// PowBigInt returns d to the power of exp, where exp is big.Int.
|
|
// Only returns error when d and exp is 0, thus result is undefined.
|
|
//
|
|
// When exponent is negative the returned decimal will have maximum precision of PowPrecisionNegativeExponent places after decimal point.
|
|
//
|
|
// Example:
|
|
//
|
|
// d1, err := decimal.NewFromFloat(3.0).PowBigInt(big.NewInt(3))
|
|
// d1.String() // output: "27"
|
|
//
|
|
// d2, err := decimal.NewFromFloat(629.25).PowBigInt(big.NewInt(5))
|
|
// d2.String() // output: "98654323103449.5673828125"
|
|
func (d Decimal) PowBigInt(exp *big.Int) (Decimal, error) {
|
|
return d.powBigIntWithPrecision(exp, int32(PowPrecisionNegativeExponent))
|
|
}
|
|
|
|
func (d Decimal) powBigIntWithPrecision(exp *big.Int, precision int32) (Decimal, error) {
|
|
if d.IsZero() && exp.Sign() == 0 {
|
|
return Decimal{}, fmt.Errorf("cannot represent undefined value of 0**0")
|
|
}
|
|
|
|
tmpExp := new(big.Int).Set(exp)
|
|
isExpNeg := exp.Sign() < 0
|
|
|
|
if isExpNeg {
|
|
tmpExp.Abs(tmpExp)
|
|
}
|
|
|
|
n, result := d, New(1, 0)
|
|
|
|
for tmpExp.Sign() > 0 {
|
|
if tmpExp.Bit(0) == 1 {
|
|
result = result.Mul(n)
|
|
}
|
|
tmpExp.Rsh(tmpExp, 1)
|
|
|
|
if tmpExp.Sign() > 0 {
|
|
n = n.Mul(n)
|
|
}
|
|
}
|
|
|
|
if isExpNeg {
|
|
return New(1, 0).DivRound(result, precision), nil
|
|
}
|
|
|
|
return result, nil
|
|
}
|
|
|
|
// ExpHullAbrham calculates the natural exponent of decimal (e to the power of d) using Hull-Abraham algorithm.
|
|
// OverallPrecision argument specifies the overall precision of the result (integer part + decimal part).
|
|
//
|
|
// ExpHullAbrham is faster than ExpTaylor for small precision values, but it is much slower for large precision values.
|
|
//
|
|
// Example:
|
|
//
|
|
// NewFromFloat(26.1).ExpHullAbrham(2).String() // output: "220000000000"
|
|
// NewFromFloat(26.1).ExpHullAbrham(20).String() // output: "216314672147.05767284"
|
|
func (d Decimal) ExpHullAbrham(overallPrecision uint32) (Decimal, error) {
|
|
// Algorithm based on Variable precision exponential function.
|
|
// ACM Transactions on Mathematical Software by T. E. Hull & A. Abrham.
|
|
if d.IsZero() {
|
|
return Decimal{oneInt, 0}, nil
|
|
}
|
|
|
|
currentPrecision := overallPrecision
|
|
|
|
// Algorithm does not work if currentPrecision * 23 < |x|.
|
|
// Precision is automatically increased in such cases, so the value can be calculated precisely.
|
|
// If newly calculated precision is higher than ExpMaxIterations the currentPrecision will not be changed.
|
|
f := d.Abs().InexactFloat64()
|
|
if ncp := f / 23; ncp > float64(currentPrecision) && ncp < float64(ExpMaxIterations) {
|
|
currentPrecision = uint32(math.Ceil(ncp))
|
|
}
|
|
|
|
// fail if abs(d) beyond an over/underflow threshold
|
|
overflowThreshold := New(23*int64(currentPrecision), 0)
|
|
if d.Abs().Cmp(overflowThreshold) > 0 {
|
|
return Decimal{}, fmt.Errorf("over/underflow threshold, exp(x) cannot be calculated precisely")
|
|
}
|
|
|
|
// Return 1 if abs(d) small enough; this also avoids later over/underflow
|
|
overflowThreshold2 := New(9, -int32(currentPrecision)-1)
|
|
if d.Abs().Cmp(overflowThreshold2) <= 0 {
|
|
return Decimal{oneInt, d.exp}, nil
|
|
}
|
|
|
|
// t is the smallest integer >= 0 such that the corresponding abs(d/k) < 1
|
|
t := d.exp + int32(d.NumDigits()) // Add d.NumDigits because the paper assumes that d.value [0.1, 1)
|
|
|
|
if t < 0 {
|
|
t = 0
|
|
}
|
|
|
|
k := New(1, t) // reduction factor
|
|
r := Decimal{new(big.Int).Set(d.value), d.exp - t} // reduced argument
|
|
p := int32(currentPrecision) + t + 2 // precision for calculating the sum
|
|
|
|
// Determine n, the number of therms for calculating sum
|
|
// use first Newton step (1.435p - 1.182) / log10(p/abs(r))
|
|
// for solving appropriate equation, along with directed
|
|
// roundings and simple rational bound for log10(p/abs(r))
|
|
rf := r.Abs().InexactFloat64()
|
|
pf := float64(p)
|
|
nf := math.Ceil((1.453*pf - 1.182) / math.Log10(pf/rf))
|
|
if nf > float64(ExpMaxIterations) || math.IsNaN(nf) {
|
|
return Decimal{}, fmt.Errorf("exact value cannot be calculated in <=ExpMaxIterations iterations")
|
|
}
|
|
n := int64(nf)
|
|
|
|
tmp := New(0, 0)
|
|
sum := New(1, 0)
|
|
one := New(1, 0)
|
|
for i := n - 1; i > 0; i-- {
|
|
tmp.value.SetInt64(i)
|
|
sum = sum.Mul(r.DivRound(tmp, p))
|
|
sum = sum.Add(one)
|
|
}
|
|
|
|
ki := k.IntPart()
|
|
res := New(1, 0)
|
|
for i := ki; i > 0; i-- {
|
|
res = res.Mul(sum)
|
|
}
|
|
|
|
resNumDigits := int32(res.NumDigits())
|
|
|
|
var roundDigits int32
|
|
if resNumDigits > abs(res.exp) {
|
|
roundDigits = int32(currentPrecision) - resNumDigits - res.exp
|
|
} else {
|
|
roundDigits = int32(currentPrecision)
|
|
}
|
|
|
|
res = res.Round(roundDigits)
|
|
|
|
return res, nil
|
|
}
|
|
|
|
// ExpTaylor calculates the natural exponent of decimal (e to the power of d) using Taylor series expansion.
|
|
// Precision argument specifies how precise the result must be (number of digits after decimal point).
|
|
// Negative precision is allowed.
|
|
//
|
|
// ExpTaylor is much faster for large precision values than ExpHullAbrham.
|
|
//
|
|
// Example:
|
|
//
|
|
// d, err := NewFromFloat(26.1).ExpTaylor(2).String()
|
|
// d.String() // output: "216314672147.06"
|
|
//
|
|
// NewFromFloat(26.1).ExpTaylor(20).String()
|
|
// d.String() // output: "216314672147.05767284062928674083"
|
|
//
|
|
// NewFromFloat(26.1).ExpTaylor(-10).String()
|
|
// d.String() // output: "220000000000"
|
|
func (d Decimal) ExpTaylor(precision int32) (Decimal, error) {
|
|
// Note(mwoss): Implementation can be optimized by exclusively using big.Int API only
|
|
if d.IsZero() {
|
|
return Decimal{oneInt, 0}.Round(precision), nil
|
|
}
|
|
|
|
var epsilon Decimal
|
|
var divPrecision int32
|
|
if precision < 0 {
|
|
epsilon = New(1, -1)
|
|
divPrecision = 8
|
|
} else {
|
|
epsilon = New(1, -precision-1)
|
|
divPrecision = precision + 1
|
|
}
|
|
|
|
decAbs := d.Abs()
|
|
pow := d.Abs()
|
|
factorial := New(1, 0)
|
|
|
|
result := New(1, 0)
|
|
|
|
for i := int64(1); ; {
|
|
step := pow.DivRound(factorial, divPrecision)
|
|
result = result.Add(step)
|
|
|
|
// Stop Taylor series when current step is smaller than epsilon
|
|
if step.Cmp(epsilon) < 0 {
|
|
break
|
|
}
|
|
|
|
pow = pow.Mul(decAbs)
|
|
|
|
i++
|
|
|
|
// Calculate next factorial number or retrieve cached value
|
|
if len(factorials) >= int(i) && !factorials[i-1].IsZero() {
|
|
factorial = factorials[i-1]
|
|
} else {
|
|
// To avoid any race conditions, firstly the zero value is appended to a slice to create
|
|
// a spot for newly calculated factorial. After that, the zero value is replaced by calculated
|
|
// factorial using the index notation.
|
|
factorial = factorials[i-2].Mul(New(i, 0))
|
|
factorials = append(factorials, Zero)
|
|
factorials[i-1] = factorial
|
|
}
|
|
}
|
|
|
|
if d.Sign() < 0 {
|
|
result = New(1, 0).DivRound(result, precision+1)
|
|
}
|
|
|
|
result = result.Round(precision)
|
|
return result, nil
|
|
}
|
|
|
|
// Ln calculates natural logarithm of d.
|
|
// Precision argument specifies how precise the result must be (number of digits after decimal point).
|
|
// Negative precision is allowed.
|
|
//
|
|
// Example:
|
|
//
|
|
// d1, err := NewFromFloat(13.3).Ln(2)
|
|
// d1.String() // output: "2.59"
|
|
//
|
|
// d2, err := NewFromFloat(579.161).Ln(10)
|
|
// d2.String() // output: "6.3615805046"
|
|
func (d Decimal) Ln(precision int32) (Decimal, error) {
|
|
// Algorithm based on The Use of Iteration Methods for Approximating the Natural Logarithm,
|
|
// James F. Epperson, The American Mathematical Monthly, Vol. 96, No. 9, November 1989, pp. 831-835.
|
|
if d.IsNegative() {
|
|
return Decimal{}, fmt.Errorf("cannot calculate natural logarithm for negative decimals")
|
|
}
|
|
|
|
if d.IsZero() {
|
|
return Decimal{}, fmt.Errorf("cannot represent natural logarithm of 0, result: -infinity")
|
|
}
|
|
|
|
calcPrecision := precision + 2
|
|
z := d.Copy()
|
|
|
|
var comp1, comp3, comp2, comp4, reduceAdjust Decimal
|
|
comp1 = z.Sub(Decimal{oneInt, 0})
|
|
comp3 = Decimal{oneInt, -1}
|
|
|
|
// for decimal in range [0.9, 1.1] where ln(d) is close to 0
|
|
usePowerSeries := false
|
|
|
|
if comp1.Abs().Cmp(comp3) <= 0 {
|
|
usePowerSeries = true
|
|
} else {
|
|
// reduce input decimal to range [0.1, 1)
|
|
expDelta := int32(z.NumDigits()) + z.exp
|
|
z.exp -= expDelta
|
|
|
|
// Input decimal was reduced by factor of 10^expDelta, thus we will need to add
|
|
// ln(10^expDelta) = expDelta * ln(10)
|
|
// to the result to compensate that
|
|
ln10 := ln10.withPrecision(calcPrecision)
|
|
reduceAdjust = NewFromInt32(expDelta)
|
|
reduceAdjust = reduceAdjust.Mul(ln10)
|
|
|
|
comp1 = z.Sub(Decimal{oneInt, 0})
|
|
|
|
if comp1.Abs().Cmp(comp3) <= 0 {
|
|
usePowerSeries = true
|
|
} else {
|
|
// initial estimate using floats
|
|
zFloat := z.InexactFloat64()
|
|
comp1 = NewFromFloat(math.Log(zFloat))
|
|
}
|
|
}
|
|
|
|
epsilon := Decimal{oneInt, -calcPrecision}
|
|
|
|
if usePowerSeries {
|
|
// Power Series - https://en.wikipedia.org/wiki/Logarithm#Power_series
|
|
// Calculating n-th term of formula: ln(z+1) = 2 sum [ 1 / (2n+1) * (z / (z+2))^(2n+1) ]
|
|
// until the difference between current and next term is smaller than epsilon.
|
|
// Coverage quite fast for decimals close to 1.0
|
|
|
|
// z + 2
|
|
comp2 = comp1.Add(Decimal{twoInt, 0})
|
|
// z / (z + 2)
|
|
comp3 = comp1.DivRound(comp2, calcPrecision)
|
|
// 2 * (z / (z + 2))
|
|
comp1 = comp3.Add(comp3)
|
|
comp2 = comp1.Copy()
|
|
|
|
for n := 1; ; n++ {
|
|
// 2 * (z / (z+2))^(2n+1)
|
|
comp2 = comp2.Mul(comp3).Mul(comp3)
|
|
|
|
// 1 / (2n+1) * 2 * (z / (z+2))^(2n+1)
|
|
comp4 = NewFromInt(int64(2*n + 1))
|
|
comp4 = comp2.DivRound(comp4, calcPrecision)
|
|
|
|
// comp1 = 2 sum [ 1 / (2n+1) * (z / (z+2))^(2n+1) ]
|
|
comp1 = comp1.Add(comp4)
|
|
|
|
if comp4.Abs().Cmp(epsilon) <= 0 {
|
|
break
|
|
}
|
|
}
|
|
} else {
|
|
// Halley's Iteration.
|
|
// Calculating n-th term of formula: a_(n+1) = a_n - 2 * (exp(a_n) - z) / (exp(a_n) + z),
|
|
// until the difference between current and next term is smaller than epsilon
|
|
var prevStep Decimal
|
|
maxIters := calcPrecision*2 + 10
|
|
|
|
for i := int32(0); i < maxIters; i++ {
|
|
// exp(a_n)
|
|
comp3, _ = comp1.ExpTaylor(calcPrecision)
|
|
// exp(a_n) - z
|
|
comp2 = comp3.Sub(z)
|
|
// 2 * (exp(a_n) - z)
|
|
comp2 = comp2.Add(comp2)
|
|
// exp(a_n) + z
|
|
comp4 = comp3.Add(z)
|
|
// 2 * (exp(a_n) - z) / (exp(a_n) + z)
|
|
comp3 = comp2.DivRound(comp4, calcPrecision)
|
|
// comp1 = a_(n+1) = a_n - 2 * (exp(a_n) - z) / (exp(a_n) + z)
|
|
comp1 = comp1.Sub(comp3)
|
|
|
|
if prevStep.Add(comp3).IsZero() {
|
|
// If iteration steps oscillate we should return early and prevent an infinity loop
|
|
// NOTE(mwoss): This should be quite a rare case, returning error is not necessary
|
|
break
|
|
}
|
|
|
|
if comp3.Abs().Cmp(epsilon) <= 0 {
|
|
break
|
|
}
|
|
|
|
prevStep = comp3
|
|
}
|
|
}
|
|
|
|
comp1 = comp1.Add(reduceAdjust)
|
|
|
|
return comp1.Round(precision), nil
|
|
}
|
|
|
|
// NumDigits returns the number of digits of the decimal coefficient (d.Value)
|
|
func (d Decimal) NumDigits() int {
|
|
if d.value == nil {
|
|
return 1
|
|
}
|
|
|
|
if d.value.IsInt64() {
|
|
i64 := d.value.Int64()
|
|
// restrict fast path to integers with exact conversion to float64
|
|
if i64 <= (1<<53) && i64 >= -(1<<53) {
|
|
if i64 == 0 {
|
|
return 1
|
|
}
|
|
return int(math.Log10(math.Abs(float64(i64)))) + 1
|
|
}
|
|
}
|
|
|
|
estimatedNumDigits := int(float64(d.value.BitLen()) / math.Log2(10))
|
|
|
|
// estimatedNumDigits (lg10) may be off by 1, need to verify
|
|
digitsBigInt := big.NewInt(int64(estimatedNumDigits))
|
|
errorCorrectionUnit := digitsBigInt.Exp(tenInt, digitsBigInt, nil)
|
|
|
|
if d.value.CmpAbs(errorCorrectionUnit) >= 0 {
|
|
return estimatedNumDigits + 1
|
|
}
|
|
|
|
return estimatedNumDigits
|
|
}
|
|
|
|
// IsInteger returns true when decimal can be represented as an integer value, otherwise, it returns false.
|
|
func (d Decimal) IsInteger() bool {
|
|
// The most typical case, all decimal with exponent higher or equal 0 can be represented as integer
|
|
if d.exp >= 0 {
|
|
return true
|
|
}
|
|
// When the exponent is negative we have to check every number after the decimal place
|
|
// If all of them are zeroes, we are sure that given decimal can be represented as an integer
|
|
var r big.Int
|
|
q := new(big.Int).Set(d.value)
|
|
for z := abs(d.exp); z > 0; z-- {
|
|
q.QuoRem(q, tenInt, &r)
|
|
if r.Cmp(zeroInt) != 0 {
|
|
return false
|
|
}
|
|
}
|
|
return true
|
|
}
|
|
|
|
// Abs calculates absolute value of any int32. Used for calculating absolute value of decimal's exponent.
|
|
func abs(n int32) int32 {
|
|
if n < 0 {
|
|
return -n
|
|
}
|
|
return n
|
|
}
|
|
|
|
// Cmp compares the numbers represented by d and d2 and returns:
|
|
//
|
|
// -1 if d < d2
|
|
// 0 if d == d2
|
|
// +1 if d > d2
|
|
func (d Decimal) Cmp(d2 Decimal) int {
|
|
d.ensureInitialized()
|
|
d2.ensureInitialized()
|
|
|
|
if d.exp == d2.exp {
|
|
return d.value.Cmp(d2.value)
|
|
}
|
|
|
|
rd, rd2 := RescalePair(d, d2)
|
|
|
|
return rd.value.Cmp(rd2.value)
|
|
}
|
|
|
|
// Compare compares the numbers represented by d and d2 and returns:
|
|
//
|
|
// -1 if d < d2
|
|
// 0 if d == d2
|
|
// +1 if d > d2
|
|
func (d Decimal) Compare(d2 Decimal) int {
|
|
return d.Cmp(d2)
|
|
}
|
|
|
|
// Equal returns whether the numbers represented by d and d2 are equal.
|
|
func (d Decimal) Equal(d2 Decimal) bool {
|
|
return d.Cmp(d2) == 0
|
|
}
|
|
|
|
// Deprecated: Equals is deprecated, please use Equal method instead.
|
|
func (d Decimal) Equals(d2 Decimal) bool {
|
|
return d.Equal(d2)
|
|
}
|
|
|
|
// GreaterThan (GT) returns true when d is greater than d2.
|
|
func (d Decimal) GreaterThan(d2 Decimal) bool {
|
|
return d.Cmp(d2) == 1
|
|
}
|
|
|
|
// GreaterThanOrEqual (GTE) returns true when d is greater than or equal to d2.
|
|
func (d Decimal) GreaterThanOrEqual(d2 Decimal) bool {
|
|
cmp := d.Cmp(d2)
|
|
return cmp == 1 || cmp == 0
|
|
}
|
|
|
|
// LessThan (LT) returns true when d is less than d2.
|
|
func (d Decimal) LessThan(d2 Decimal) bool {
|
|
return d.Cmp(d2) == -1
|
|
}
|
|
|
|
// LessThanOrEqual (LTE) returns true when d is less than or equal to d2.
|
|
func (d Decimal) LessThanOrEqual(d2 Decimal) bool {
|
|
cmp := d.Cmp(d2)
|
|
return cmp == -1 || cmp == 0
|
|
}
|
|
|
|
// Sign returns:
|
|
//
|
|
// -1 if d < 0
|
|
// 0 if d == 0
|
|
// +1 if d > 0
|
|
func (d Decimal) Sign() int {
|
|
if d.value == nil {
|
|
return 0
|
|
}
|
|
return d.value.Sign()
|
|
}
|
|
|
|
// IsPositive return
|
|
//
|
|
// true if d > 0
|
|
// false if d == 0
|
|
// false if d < 0
|
|
func (d Decimal) IsPositive() bool {
|
|
return d.Sign() == 1
|
|
}
|
|
|
|
// IsNegative return
|
|
//
|
|
// true if d < 0
|
|
// false if d == 0
|
|
// false if d > 0
|
|
func (d Decimal) IsNegative() bool {
|
|
return d.Sign() == -1
|
|
}
|
|
|
|
// IsZero return
|
|
//
|
|
// true if d == 0
|
|
// false if d > 0
|
|
// false if d < 0
|
|
func (d Decimal) IsZero() bool {
|
|
return d.Sign() == 0
|
|
}
|
|
|
|
// Exponent returns the exponent, or scale component of the decimal.
|
|
func (d Decimal) Exponent() int32 {
|
|
return d.exp
|
|
}
|
|
|
|
// Coefficient returns the coefficient of the decimal. It is scaled by 10^Exponent()
|
|
func (d Decimal) Coefficient() *big.Int {
|
|
d.ensureInitialized()
|
|
// we copy the coefficient so that mutating the result does not mutate the Decimal.
|
|
return new(big.Int).Set(d.value)
|
|
}
|
|
|
|
// CoefficientInt64 returns the coefficient of the decimal as int64. It is scaled by 10^Exponent()
|
|
// If coefficient cannot be represented in an int64, the result will be undefined.
|
|
func (d Decimal) CoefficientInt64() int64 {
|
|
d.ensureInitialized()
|
|
return d.value.Int64()
|
|
}
|
|
|
|
// IntPart returns the integer component of the decimal.
|
|
func (d Decimal) IntPart() int64 {
|
|
scaledD := d.rescale(0)
|
|
return scaledD.value.Int64()
|
|
}
|
|
|
|
// BigInt returns integer component of the decimal as a BigInt.
|
|
func (d Decimal) BigInt() *big.Int {
|
|
scaledD := d.rescale(0)
|
|
return scaledD.value
|
|
}
|
|
|
|
// BigFloat returns decimal as BigFloat.
|
|
// Be aware that casting decimal to BigFloat might cause a loss of precision.
|
|
func (d Decimal) BigFloat() *big.Float {
|
|
f := &big.Float{}
|
|
f.SetString(d.String())
|
|
return f
|
|
}
|
|
|
|
// Rat returns a rational number representation of the decimal.
|
|
func (d Decimal) Rat() *big.Rat {
|
|
d.ensureInitialized()
|
|
if d.exp <= 0 {
|
|
// NOTE(vadim): must negate after casting to prevent int32 overflow
|
|
denom := new(big.Int).Exp(tenInt, big.NewInt(-int64(d.exp)), nil)
|
|
return new(big.Rat).SetFrac(d.value, denom)
|
|
}
|
|
|
|
mul := new(big.Int).Exp(tenInt, big.NewInt(int64(d.exp)), nil)
|
|
num := new(big.Int).Mul(d.value, mul)
|
|
return new(big.Rat).SetFrac(num, oneInt)
|
|
}
|
|
|
|
// Float64 returns the nearest float64 value for d and a bool indicating
|
|
// whether f represents d exactly.
|
|
// For more details, see the documentation for big.Rat.Float64
|
|
func (d Decimal) Float64() (f float64, exact bool) {
|
|
return d.Rat().Float64()
|
|
}
|
|
|
|
// InexactFloat64 returns the nearest float64 value for d.
|
|
// It doesn't indicate if the returned value represents d exactly.
|
|
func (d Decimal) InexactFloat64() float64 {
|
|
f, _ := d.Float64()
|
|
return f
|
|
}
|
|
|
|
// String returns the string representation of the decimal
|
|
// with the fixed point.
|
|
//
|
|
// Example:
|
|
//
|
|
// d := New(-12345, -3)
|
|
// println(d.String())
|
|
//
|
|
// Output:
|
|
//
|
|
// -12.345
|
|
func (d Decimal) String() string {
|
|
return d.string(true)
|
|
}
|
|
|
|
// StringFixed returns a rounded fixed-point string with places digits after
|
|
// the decimal point.
|
|
//
|
|
// Example:
|
|
//
|
|
// NewFromFloat(0).StringFixed(2) // output: "0.00"
|
|
// NewFromFloat(0).StringFixed(0) // output: "0"
|
|
// NewFromFloat(5.45).StringFixed(0) // output: "5"
|
|
// NewFromFloat(5.45).StringFixed(1) // output: "5.5"
|
|
// NewFromFloat(5.45).StringFixed(2) // output: "5.45"
|
|
// NewFromFloat(5.45).StringFixed(3) // output: "5.450"
|
|
// NewFromFloat(545).StringFixed(-1) // output: "550"
|
|
func (d Decimal) StringFixed(places int32) string {
|
|
rounded := d.Round(places)
|
|
return rounded.string(false)
|
|
}
|
|
|
|
// StringFixedBank returns a banker rounded fixed-point string with places digits
|
|
// after the decimal point.
|
|
//
|
|
// Example:
|
|
//
|
|
// NewFromFloat(0).StringFixedBank(2) // output: "0.00"
|
|
// NewFromFloat(0).StringFixedBank(0) // output: "0"
|
|
// NewFromFloat(5.45).StringFixedBank(0) // output: "5"
|
|
// NewFromFloat(5.45).StringFixedBank(1) // output: "5.4"
|
|
// NewFromFloat(5.45).StringFixedBank(2) // output: "5.45"
|
|
// NewFromFloat(5.45).StringFixedBank(3) // output: "5.450"
|
|
// NewFromFloat(545).StringFixedBank(-1) // output: "540"
|
|
func (d Decimal) StringFixedBank(places int32) string {
|
|
rounded := d.RoundBank(places)
|
|
return rounded.string(false)
|
|
}
|
|
|
|
// StringFixedCash returns a Swedish/Cash rounded fixed-point string. For
|
|
// more details see the documentation at function RoundCash.
|
|
func (d Decimal) StringFixedCash(interval uint8) string {
|
|
rounded := d.RoundCash(interval)
|
|
return rounded.string(false)
|
|
}
|
|
|
|
// Round rounds the decimal to places decimal places (half away from zero).
|
|
// If places < 0, it will round the integer part to the nearest 10^(-places).
|
|
//
|
|
// Example:
|
|
//
|
|
// NewFromFloat(5.45).Round(1).String() // output: "5.5"
|
|
// NewFromFloat(545).Round(-1).String() // output: "550"
|
|
func (d Decimal) Round(places int32) Decimal {
|
|
if d.exp == -places {
|
|
return d
|
|
}
|
|
// truncate to places + 1
|
|
ret := d.rescale(-places - 1)
|
|
|
|
// add sign(d) * 0.5
|
|
if ret.value.Sign() < 0 {
|
|
ret.value.Sub(ret.value, fiveInt)
|
|
} else {
|
|
ret.value.Add(ret.value, fiveInt)
|
|
}
|
|
|
|
// floor for positive numbers, ceil for negative numbers
|
|
_, m := ret.value.DivMod(ret.value, tenInt, new(big.Int))
|
|
ret.exp++
|
|
if ret.value.Sign() < 0 && m.Cmp(zeroInt) != 0 {
|
|
ret.value.Add(ret.value, oneInt)
|
|
}
|
|
|
|
return ret
|
|
}
|
|
|
|
// RoundHalfTowardZero rounds the decimal to places decimal places (half toward zero).
|
|
// If places < 0, it will round the integer part to the nearest 10^(-places).
|
|
//
|
|
// Example:
|
|
//
|
|
// NewFromFloat(5.45).RoundHalfTowardZero(1).String() // output: "5.4"
|
|
// NewFromFloat(545).RoundHalfTowardZero(-1).String() // output: "540"
|
|
func (d Decimal) RoundHalfTowardZero(places int32) Decimal {
|
|
if d.exp == -places {
|
|
return d
|
|
}
|
|
// truncate to places + 1
|
|
ret := d.rescale(-places - 1)
|
|
|
|
// add sign(d) * 0.4
|
|
if ret.value.Sign() < 0 {
|
|
ret.value.Sub(ret.value, fourInt)
|
|
} else {
|
|
ret.value.Add(ret.value, fourInt)
|
|
}
|
|
|
|
// floor for positive numbers, ceil for negative numbers
|
|
_, m := ret.value.DivMod(ret.value, tenInt, new(big.Int))
|
|
ret.exp++
|
|
if ret.value.Sign() < 0 && m.Cmp(zeroInt) != 0 {
|
|
ret.value.Add(ret.value, oneInt)
|
|
}
|
|
|
|
return ret
|
|
}
|
|
|
|
// RoundHalfUp rounds the decimal half towards +infinity.
|
|
//
|
|
// Example:
|
|
//
|
|
// NewFromFloat(545).RoundHalfUp(-2).String() // output: "500"
|
|
// NewFromFloat(500).RoundHalfUp(-2).String() // output: "500"
|
|
// NewFromFloat(1.1001).RoundHalfUp(2).String() // output: "1.10"
|
|
// NewFromFloat(-1.454).RoundHalfUp(1).String() // output: "-1.4"
|
|
// NewFromFloat(-1.464).RoundHalfUp(1).String() // output: "-1.5"
|
|
func (d Decimal) RoundHalfUp(places int32) Decimal {
|
|
if d.exp == -places {
|
|
return d
|
|
}
|
|
// truncate to places + 1
|
|
ret := d.rescale(-places - 1)
|
|
|
|
// add sign(d) * 0.5 if sign(d) >= 0 else sign(d) * 0.4
|
|
if ret.value.Sign() < 0 {
|
|
ret.value.Sub(ret.value, fourInt)
|
|
} else {
|
|
ret.value.Add(ret.value, fiveInt)
|
|
}
|
|
|
|
// floor for positive numbers, ceil for negative numbers
|
|
_, m := ret.value.DivMod(ret.value, tenInt, new(big.Int))
|
|
ret.exp++
|
|
if ret.value.Sign() < 0 && m.Cmp(zeroInt) != 0 {
|
|
ret.value.Add(ret.value, oneInt)
|
|
}
|
|
|
|
return ret
|
|
}
|
|
|
|
// RoundHalfDown rounds the decimal half towards -infinity.
|
|
//
|
|
// Example:
|
|
//
|
|
// NewFromFloat(550).RoundHalfDown(-2).String() // output: "500"
|
|
// NewFromFloat(560).RoundHalfDown(-2).String() // output: "600"
|
|
// NewFromFloat(1.1001).RoundHalfDown(2).String() // output: "1.11"
|
|
// NewFromFloat(-1.454).RoundHalfDown(1).String() // output: "-1.5"
|
|
// NewFromFloat(-1.444).RoundHalfDown(1).String() // output: "-1.4"
|
|
func (d Decimal) RoundHalfDown(places int32) Decimal {
|
|
if d.exp == -places {
|
|
return d
|
|
}
|
|
// truncate to places + 1
|
|
ret := d.rescale(-places - 1)
|
|
|
|
// add sign(d) * 0.5 if sign(d) < 0 else sign(d) * 0.4
|
|
if ret.value.Sign() < 0 {
|
|
ret.value.Sub(ret.value, fiveInt)
|
|
} else {
|
|
ret.value.Add(ret.value, fourInt)
|
|
}
|
|
|
|
// floor for positive numbers, ceil for negative numbers
|
|
_, m := ret.value.DivMod(ret.value, tenInt, new(big.Int))
|
|
ret.exp++
|
|
if ret.value.Sign() < 0 && m.Cmp(zeroInt) != 0 {
|
|
ret.value.Add(ret.value, oneInt)
|
|
}
|
|
|
|
return ret
|
|
}
|
|
|
|
// RoundCeil rounds the decimal towards +infinity.
|
|
//
|
|
// Example:
|
|
//
|
|
// NewFromFloat(545).RoundCeil(-2).String() // output: "600"
|
|
// NewFromFloat(500).RoundCeil(-2).String() // output: "500"
|
|
// NewFromFloat(1.1001).RoundCeil(2).String() // output: "1.11"
|
|
// NewFromFloat(-1.454).RoundCeil(1).String() // output: "-1.4"
|
|
func (d Decimal) RoundCeil(places int32) Decimal {
|
|
if d.exp >= -places {
|
|
return d
|
|
}
|
|
|
|
rescaled := d.rescale(-places)
|
|
if d.Equal(rescaled) {
|
|
return d
|
|
}
|
|
|
|
if d.value.Sign() > 0 {
|
|
rescaled.value.Add(rescaled.value, oneInt)
|
|
}
|
|
|
|
return rescaled
|
|
}
|
|
|
|
// RoundFloor rounds the decimal towards -infinity.
|
|
//
|
|
// Example:
|
|
//
|
|
// NewFromFloat(545).RoundFloor(-2).String() // output: "500"
|
|
// NewFromFloat(-500).RoundFloor(-2).String() // output: "-500"
|
|
// NewFromFloat(1.1001).RoundFloor(2).String() // output: "1.1"
|
|
// NewFromFloat(-1.454).RoundFloor(1).String() // output: "-1.5"
|
|
func (d Decimal) RoundFloor(places int32) Decimal {
|
|
if d.exp >= -places {
|
|
return d
|
|
}
|
|
|
|
rescaled := d.rescale(-places)
|
|
if d.Equal(rescaled) {
|
|
return d
|
|
}
|
|
|
|
if d.value.Sign() < 0 {
|
|
rescaled.value.Sub(rescaled.value, oneInt)
|
|
}
|
|
|
|
return rescaled
|
|
}
|
|
|
|
// RoundUp rounds the decimal away from zero.
|
|
//
|
|
// Example:
|
|
//
|
|
// NewFromFloat(545).RoundUp(-2).String() // output: "600"
|
|
// NewFromFloat(500).RoundUp(-2).String() // output: "500"
|
|
// NewFromFloat(1.1001).RoundUp(2).String() // output: "1.11"
|
|
// NewFromFloat(-1.454).RoundUp(1).String() // output: "-1.5"
|
|
func (d Decimal) RoundUp(places int32) Decimal {
|
|
if d.exp >= -places {
|
|
return d
|
|
}
|
|
|
|
rescaled := d.rescale(-places)
|
|
if d.Equal(rescaled) {
|
|
return d
|
|
}
|
|
|
|
if d.value.Sign() > 0 {
|
|
rescaled.value.Add(rescaled.value, oneInt)
|
|
} else if d.value.Sign() < 0 {
|
|
rescaled.value.Sub(rescaled.value, oneInt)
|
|
}
|
|
|
|
return rescaled
|
|
}
|
|
|
|
// RoundDown rounds the decimal towards zero.
|
|
//
|
|
// Example:
|
|
//
|
|
// NewFromFloat(545).RoundDown(-2).String() // output: "500"
|
|
// NewFromFloat(-500).RoundDown(-2).String() // output: "-500"
|
|
// NewFromFloat(1.1001).RoundDown(2).String() // output: "1.1"
|
|
// NewFromFloat(-1.454).RoundDown(1).String() // output: "-1.4"
|
|
func (d Decimal) RoundDown(places int32) Decimal {
|
|
if d.exp >= -places {
|
|
return d
|
|
}
|
|
|
|
rescaled := d.rescale(-places)
|
|
if d.Equal(rescaled) {
|
|
return d
|
|
}
|
|
return rescaled
|
|
}
|
|
|
|
// RoundBank rounds the decimal to places decimal places.
|
|
// If the final digit to round is equidistant from the nearest two integers the
|
|
// rounded value is taken as the even number
|
|
//
|
|
// If places < 0, it will round the integer part to the nearest 10^(-places).
|
|
//
|
|
// Examples:
|
|
//
|
|
// NewFromFloat(5.45).RoundBank(1).String() // output: "5.4"
|
|
// NewFromFloat(545).RoundBank(-1).String() // output: "540"
|
|
// NewFromFloat(5.46).RoundBank(1).String() // output: "5.5"
|
|
// NewFromFloat(546).RoundBank(-1).String() // output: "550"
|
|
// NewFromFloat(5.55).RoundBank(1).String() // output: "5.6"
|
|
// NewFromFloat(555).RoundBank(-1).String() // output: "560"
|
|
func (d Decimal) RoundBank(places int32) Decimal {
|
|
|
|
round := d.Round(places)
|
|
remainder := d.Sub(round).Abs()
|
|
|
|
half := New(5, -places-1)
|
|
if remainder.Cmp(half) == 0 && round.value.Bit(0) != 0 {
|
|
if round.value.Sign() < 0 {
|
|
round.value.Add(round.value, oneInt)
|
|
} else {
|
|
round.value.Sub(round.value, oneInt)
|
|
}
|
|
}
|
|
|
|
return round
|
|
}
|
|
|
|
// RoundCash aka Cash/Penny/öre rounding rounds decimal to a specific
|
|
// interval. The amount payable for a cash transaction is rounded to the nearest
|
|
// multiple of the minimum currency unit available. The following intervals are
|
|
// available: 5, 10, 25, 50 and 100; any other number throws a panic.
|
|
//
|
|
// 5: 5 cent rounding 3.43 => 3.45
|
|
// 10: 10 cent rounding 3.45 => 3.50 (5 gets rounded up)
|
|
// 25: 25 cent rounding 3.41 => 3.50
|
|
// 50: 50 cent rounding 3.75 => 4.00
|
|
// 100: 100 cent rounding 3.50 => 4.00
|
|
//
|
|
// For more details: https://en.wikipedia.org/wiki/Cash_rounding
|
|
func (d Decimal) RoundCash(interval uint8) Decimal {
|
|
var iVal *big.Int
|
|
switch interval {
|
|
case 5:
|
|
iVal = twentyInt
|
|
case 10:
|
|
iVal = tenInt
|
|
case 25:
|
|
iVal = fourInt
|
|
case 50:
|
|
iVal = twoInt
|
|
case 100:
|
|
iVal = oneInt
|
|
default:
|
|
panic(fmt.Sprintf("Decimal does not support this Cash rounding interval `%d`. Supported: 5, 10, 25, 50, 100", interval))
|
|
}
|
|
dVal := Decimal{
|
|
value: iVal,
|
|
}
|
|
|
|
// TODO: optimize those calculations to reduce the high allocations (~29 allocs).
|
|
return d.Mul(dVal).Round(0).Div(dVal).Truncate(2)
|
|
}
|
|
|
|
// Floor returns the nearest integer value less than or equal to d.
|
|
func (d Decimal) Floor() Decimal {
|
|
d.ensureInitialized()
|
|
|
|
if d.exp >= 0 {
|
|
return d
|
|
}
|
|
|
|
exp := big.NewInt(10)
|
|
|
|
// NOTE(vadim): must negate after casting to prevent int32 overflow
|
|
exp.Exp(exp, big.NewInt(-int64(d.exp)), nil)
|
|
|
|
z := new(big.Int).Div(d.value, exp)
|
|
return Decimal{value: z, exp: 0}
|
|
}
|
|
|
|
// Ceil returns the nearest integer value greater than or equal to d.
|
|
func (d Decimal) Ceil() Decimal {
|
|
d.ensureInitialized()
|
|
|
|
if d.exp >= 0 {
|
|
return d
|
|
}
|
|
|
|
exp := big.NewInt(10)
|
|
|
|
// NOTE(vadim): must negate after casting to prevent int32 overflow
|
|
exp.Exp(exp, big.NewInt(-int64(d.exp)), nil)
|
|
|
|
z, m := new(big.Int).DivMod(d.value, exp, new(big.Int))
|
|
if m.Cmp(zeroInt) != 0 {
|
|
z.Add(z, oneInt)
|
|
}
|
|
return Decimal{value: z, exp: 0}
|
|
}
|
|
|
|
// Truncate truncates off digits from the number, without rounding.
|
|
//
|
|
// NOTE: precision is the last digit that will not be truncated (must be >= 0).
|
|
//
|
|
// Example:
|
|
//
|
|
// decimal.NewFromString("123.456").Truncate(2).String() // "123.45"
|
|
func (d Decimal) Truncate(precision int32) Decimal {
|
|
d.ensureInitialized()
|
|
if precision >= 0 && -precision > d.exp {
|
|
return d.rescale(-precision)
|
|
}
|
|
return d
|
|
}
|
|
|
|
// UnmarshalJSON implements the json.Unmarshaler interface.
|
|
func (d *Decimal) UnmarshalJSON(decimalBytes []byte) error {
|
|
if string(decimalBytes) == "null" {
|
|
return nil
|
|
}
|
|
|
|
str, err := unquoteIfQuoted(decimalBytes)
|
|
if err != nil {
|
|
return fmt.Errorf("error decoding string '%s': %s", decimalBytes, err)
|
|
}
|
|
|
|
decimal, err := NewFromString(str)
|
|
*d = decimal
|
|
if err != nil {
|
|
return fmt.Errorf("error decoding string '%s': %s", str, err)
|
|
}
|
|
return nil
|
|
}
|
|
|
|
// MarshalJSON implements the json.Marshaler interface.
|
|
func (d Decimal) MarshalJSON() ([]byte, error) {
|
|
var str string
|
|
if MarshalJSONWithoutQuotes {
|
|
str = d.String()
|
|
} else {
|
|
str = "\"" + d.String() + "\""
|
|
}
|
|
return []byte(str), nil
|
|
}
|
|
|
|
// UnmarshalBinary implements the encoding.BinaryUnmarshaler interface. As a string representation
|
|
// is already used when encoding to text, this method stores that string as []byte
|
|
func (d *Decimal) UnmarshalBinary(data []byte) error {
|
|
// Verify we have at least 4 bytes for the exponent. The GOB encoded value
|
|
// may be empty.
|
|
if len(data) < 4 {
|
|
return fmt.Errorf("error decoding binary %v: expected at least 4 bytes, got %d", data, len(data))
|
|
}
|
|
|
|
// Extract the exponent
|
|
d.exp = int32(binary.BigEndian.Uint32(data[:4]))
|
|
|
|
// Extract the value
|
|
d.value = new(big.Int)
|
|
if err := d.value.GobDecode(data[4:]); err != nil {
|
|
return fmt.Errorf("error decoding binary %v: %s", data, err)
|
|
}
|
|
|
|
return nil
|
|
}
|
|
|
|
// MarshalBinary implements the encoding.BinaryMarshaler interface.
|
|
func (d Decimal) MarshalBinary() (data []byte, err error) {
|
|
// exp is written first, but encode value first to know output size
|
|
var valueData []byte
|
|
if valueData, err = d.value.GobEncode(); err != nil {
|
|
return nil, err
|
|
}
|
|
|
|
// Write the exponent in front, since it's a fixed size
|
|
expData := make([]byte, 4, len(valueData)+4)
|
|
binary.BigEndian.PutUint32(expData, uint32(d.exp))
|
|
|
|
// Return the byte array
|
|
return append(expData, valueData...), nil
|
|
}
|
|
|
|
// Scan implements the sql.Scanner interface for database deserialization.
|
|
func (d *Decimal) Scan(value interface{}) error {
|
|
// first try to see if the data is stored in database as a Numeric datatype
|
|
switch v := value.(type) {
|
|
|
|
case float32:
|
|
*d = NewFromFloat(float64(v))
|
|
return nil
|
|
|
|
case float64:
|
|
// numeric in sqlite3 sends us float64
|
|
*d = NewFromFloat(v)
|
|
return nil
|
|
|
|
case int64:
|
|
// at least in sqlite3 when the value is 0 in db, the data is sent
|
|
// to us as an int64 instead of a float64 ...
|
|
*d = New(v, 0)
|
|
return nil
|
|
|
|
case uint64:
|
|
// while clickhouse may send 0 in db as uint64
|
|
*d = NewFromUint64(v)
|
|
return nil
|
|
|
|
default:
|
|
// default is trying to interpret value stored as string
|
|
str, err := unquoteIfQuoted(v)
|
|
if err != nil {
|
|
return err
|
|
}
|
|
*d, err = NewFromString(str)
|
|
return err
|
|
}
|
|
}
|
|
|
|
// Value implements the driver.Valuer interface for database serialization.
|
|
func (d Decimal) Value() (driver.Value, error) {
|
|
return d.String(), nil
|
|
}
|
|
|
|
// UnmarshalText implements the encoding.TextUnmarshaler interface for XML
|
|
// deserialization.
|
|
func (d *Decimal) UnmarshalText(text []byte) error {
|
|
str := string(text)
|
|
|
|
dec, err := NewFromString(str)
|
|
*d = dec
|
|
if err != nil {
|
|
return fmt.Errorf("error decoding string '%s': %s", str, err)
|
|
}
|
|
|
|
return nil
|
|
}
|
|
|
|
// MarshalText implements the encoding.TextMarshaler interface for XML
|
|
// serialization.
|
|
func (d Decimal) MarshalText() (text []byte, err error) {
|
|
return []byte(d.String()), nil
|
|
}
|
|
|
|
// GobEncode implements the gob.GobEncoder interface for gob serialization.
|
|
func (d Decimal) GobEncode() ([]byte, error) {
|
|
return d.MarshalBinary()
|
|
}
|
|
|
|
// GobDecode implements the gob.GobDecoder interface for gob serialization.
|
|
func (d *Decimal) GobDecode(data []byte) error {
|
|
return d.UnmarshalBinary(data)
|
|
}
|
|
|
|
// StringScaled first scales the decimal then calls .String() on it.
|
|
//
|
|
// Deprecated: buggy and unintuitive. Use StringFixed instead.
|
|
func (d Decimal) StringScaled(exp int32) string {
|
|
return d.rescale(exp).String()
|
|
}
|
|
|
|
func (d Decimal) string(trimTrailingZeros bool) string {
|
|
if d.exp >= 0 {
|
|
return d.rescale(0).value.String()
|
|
}
|
|
|
|
abs := new(big.Int).Abs(d.value)
|
|
str := abs.String()
|
|
|
|
var intPart, fractionalPart string
|
|
|
|
// NOTE(vadim): this cast to int will cause bugs if d.exp == INT_MIN
|
|
// and you are on a 32-bit machine. Won't fix this super-edge case.
|
|
dExpInt := int(d.exp)
|
|
if len(str) > -dExpInt {
|
|
intPart = str[:len(str)+dExpInt]
|
|
fractionalPart = str[len(str)+dExpInt:]
|
|
} else {
|
|
intPart = "0"
|
|
|
|
num0s := -dExpInt - len(str)
|
|
fractionalPart = strings.Repeat("0", num0s) + str
|
|
}
|
|
|
|
if trimTrailingZeros {
|
|
i := len(fractionalPart) - 1
|
|
for ; i >= 0; i-- {
|
|
if fractionalPart[i] != '0' {
|
|
break
|
|
}
|
|
}
|
|
fractionalPart = fractionalPart[:i+1]
|
|
}
|
|
|
|
number := intPart
|
|
if len(fractionalPart) > 0 {
|
|
number += "." + fractionalPart
|
|
}
|
|
|
|
if d.value.Sign() < 0 {
|
|
return "-" + number
|
|
}
|
|
|
|
return number
|
|
}
|
|
|
|
func (d *Decimal) ensureInitialized() {
|
|
if d.value == nil {
|
|
d.value = new(big.Int)
|
|
}
|
|
}
|
|
|
|
// Min returns the smallest Decimal that was passed in the arguments.
|
|
//
|
|
// To call this function with an array, you must do:
|
|
//
|
|
// Min(arr[0], arr[1:]...)
|
|
//
|
|
// This makes it harder to accidentally call Min with 0 arguments.
|
|
func Min(first Decimal, rest ...Decimal) Decimal {
|
|
ans := first
|
|
for _, item := range rest {
|
|
if item.Cmp(ans) < 0 {
|
|
ans = item
|
|
}
|
|
}
|
|
return ans
|
|
}
|
|
|
|
// Max returns the largest Decimal that was passed in the arguments.
|
|
//
|
|
// To call this function with an array, you must do:
|
|
//
|
|
// Max(arr[0], arr[1:]...)
|
|
//
|
|
// This makes it harder to accidentally call Max with 0 arguments.
|
|
func Max(first Decimal, rest ...Decimal) Decimal {
|
|
ans := first
|
|
for _, item := range rest {
|
|
if item.Cmp(ans) > 0 {
|
|
ans = item
|
|
}
|
|
}
|
|
return ans
|
|
}
|
|
|
|
// Sum returns the combined total of the provided first and rest Decimals
|
|
func Sum(first Decimal, rest ...Decimal) Decimal {
|
|
total := first
|
|
for _, item := range rest {
|
|
total = total.Add(item)
|
|
}
|
|
|
|
return total
|
|
}
|
|
|
|
// Avg returns the average value of the provided first and rest Decimals
|
|
func Avg(first Decimal, rest ...Decimal) Decimal {
|
|
count := New(int64(len(rest)+1), 0)
|
|
sum := Sum(first, rest...)
|
|
return sum.Div(count)
|
|
}
|
|
|
|
// RescalePair rescales two decimals to common exponential value (minimal exp of both decimals)
|
|
func RescalePair(d1 Decimal, d2 Decimal) (Decimal, Decimal) {
|
|
d1.ensureInitialized()
|
|
d2.ensureInitialized()
|
|
|
|
if d1.exp < d2.exp {
|
|
return d1, d2.rescale(d1.exp)
|
|
} else if d1.exp > d2.exp {
|
|
return d1.rescale(d2.exp), d2
|
|
}
|
|
|
|
return d1, d2
|
|
}
|
|
|
|
func unquoteIfQuoted(value interface{}) (string, error) {
|
|
var bytes []byte
|
|
|
|
switch v := value.(type) {
|
|
case string:
|
|
bytes = []byte(v)
|
|
case []byte:
|
|
bytes = v
|
|
default:
|
|
return "", fmt.Errorf("could not convert value '%+v' to byte array of type '%T'", value, value)
|
|
}
|
|
|
|
// If the amount is quoted, strip the quotes
|
|
if len(bytes) > 2 && bytes[0] == '"' && bytes[len(bytes)-1] == '"' {
|
|
bytes = bytes[1 : len(bytes)-1]
|
|
}
|
|
return string(bytes), nil
|
|
}
|
|
|
|
// NullDecimal represents a nullable decimal with compatibility for
|
|
// scanning null values from the database.
|
|
type NullDecimal struct {
|
|
Decimal Decimal
|
|
Valid bool
|
|
}
|
|
|
|
func NewNullDecimal(d Decimal) NullDecimal {
|
|
return NullDecimal{
|
|
Decimal: d,
|
|
Valid: true,
|
|
}
|
|
}
|
|
|
|
// Scan implements the sql.Scanner interface for database deserialization.
|
|
func (d *NullDecimal) Scan(value interface{}) error {
|
|
if value == nil {
|
|
d.Valid = false
|
|
return nil
|
|
}
|
|
d.Valid = true
|
|
return d.Decimal.Scan(value)
|
|
}
|
|
|
|
// Value implements the driver.Valuer interface for database serialization.
|
|
func (d NullDecimal) Value() (driver.Value, error) {
|
|
if !d.Valid {
|
|
return nil, nil
|
|
}
|
|
return d.Decimal.Value()
|
|
}
|
|
|
|
// UnmarshalJSON implements the json.Unmarshaler interface.
|
|
func (d *NullDecimal) UnmarshalJSON(decimalBytes []byte) error {
|
|
if string(decimalBytes) == "null" {
|
|
d.Valid = false
|
|
return nil
|
|
}
|
|
d.Valid = true
|
|
return d.Decimal.UnmarshalJSON(decimalBytes)
|
|
}
|
|
|
|
// MarshalJSON implements the json.Marshaler interface.
|
|
func (d NullDecimal) MarshalJSON() ([]byte, error) {
|
|
if !d.Valid {
|
|
return []byte("null"), nil
|
|
}
|
|
return d.Decimal.MarshalJSON()
|
|
}
|
|
|
|
// UnmarshalText implements the encoding.TextUnmarshaler interface for XML
|
|
// deserialization
|
|
func (d *NullDecimal) UnmarshalText(text []byte) error {
|
|
str := string(text)
|
|
|
|
// check for empty XML or XML without body e.g., <tag></tag>
|
|
if str == "" {
|
|
d.Valid = false
|
|
return nil
|
|
}
|
|
if err := d.Decimal.UnmarshalText(text); err != nil {
|
|
d.Valid = false
|
|
return err
|
|
}
|
|
d.Valid = true
|
|
return nil
|
|
}
|
|
|
|
// MarshalText implements the encoding.TextMarshaler interface for XML
|
|
// serialization.
|
|
func (d NullDecimal) MarshalText() (text []byte, err error) {
|
|
if !d.Valid {
|
|
return []byte{}, nil
|
|
}
|
|
return d.Decimal.MarshalText()
|
|
}
|
|
|
|
// Trig functions
|
|
|
|
// Atan returns the arctangent, in radians, of x.
|
|
func (d Decimal) Atan() Decimal {
|
|
if d.Equal(NewFromFloat(0.0)) {
|
|
return d
|
|
}
|
|
if d.GreaterThan(NewFromFloat(0.0)) {
|
|
return d.satan()
|
|
}
|
|
return d.Neg().satan().Neg()
|
|
}
|
|
|
|
func (d Decimal) xatan() Decimal {
|
|
P0 := NewFromFloat(-8.750608600031904122785e-01)
|
|
P1 := NewFromFloat(-1.615753718733365076637e+01)
|
|
P2 := NewFromFloat(-7.500855792314704667340e+01)
|
|
P3 := NewFromFloat(-1.228866684490136173410e+02)
|
|
P4 := NewFromFloat(-6.485021904942025371773e+01)
|
|
Q0 := NewFromFloat(2.485846490142306297962e+01)
|
|
Q1 := NewFromFloat(1.650270098316988542046e+02)
|
|
Q2 := NewFromFloat(4.328810604912902668951e+02)
|
|
Q3 := NewFromFloat(4.853903996359136964868e+02)
|
|
Q4 := NewFromFloat(1.945506571482613964425e+02)
|
|
z := d.Mul(d)
|
|
b1 := P0.Mul(z).Add(P1).Mul(z).Add(P2).Mul(z).Add(P3).Mul(z).Add(P4).Mul(z)
|
|
b2 := z.Add(Q0).Mul(z).Add(Q1).Mul(z).Add(Q2).Mul(z).Add(Q3).Mul(z).Add(Q4)
|
|
z = b1.Div(b2)
|
|
z = d.Mul(z).Add(d)
|
|
return z
|
|
}
|
|
|
|
// satan reduces its argument (known to be positive)
|
|
// to the range [0, 0.66] and calls xatan.
|
|
func (d Decimal) satan() Decimal {
|
|
Morebits := NewFromFloat(6.123233995736765886130e-17) // pi/2 = PIO2 + Morebits
|
|
Tan3pio8 := NewFromFloat(2.41421356237309504880) // tan(3*pi/8)
|
|
pi := NewFromFloat(3.14159265358979323846264338327950288419716939937510582097494459)
|
|
|
|
if d.LessThanOrEqual(NewFromFloat(0.66)) {
|
|
return d.xatan()
|
|
}
|
|
if d.GreaterThan(Tan3pio8) {
|
|
return pi.Div(NewFromFloat(2.0)).Sub(NewFromFloat(1.0).Div(d).xatan()).Add(Morebits)
|
|
}
|
|
return pi.Div(NewFromFloat(4.0)).Add((d.Sub(NewFromFloat(1.0)).Div(d.Add(NewFromFloat(1.0)))).xatan()).Add(NewFromFloat(0.5).Mul(Morebits))
|
|
}
|
|
|
|
// sin coefficients
|
|
var _sin = [...]Decimal{
|
|
NewFromFloat(1.58962301576546568060e-10), // 0x3de5d8fd1fd19ccd
|
|
NewFromFloat(-2.50507477628578072866e-8), // 0xbe5ae5e5a9291f5d
|
|
NewFromFloat(2.75573136213857245213e-6), // 0x3ec71de3567d48a1
|
|
NewFromFloat(-1.98412698295895385996e-4), // 0xbf2a01a019bfdf03
|
|
NewFromFloat(8.33333333332211858878e-3), // 0x3f8111111110f7d0
|
|
NewFromFloat(-1.66666666666666307295e-1), // 0xbfc5555555555548
|
|
}
|
|
|
|
// Sin returns the sine of the radian argument x.
|
|
func (d Decimal) Sin() Decimal {
|
|
PI4A := NewFromFloat(7.85398125648498535156e-1) // 0x3fe921fb40000000, Pi/4 split into three parts
|
|
PI4B := NewFromFloat(3.77489470793079817668e-8) // 0x3e64442d00000000,
|
|
PI4C := NewFromFloat(2.69515142907905952645e-15) // 0x3ce8469898cc5170,
|
|
M4PI := NewFromFloat(1.273239544735162542821171882678754627704620361328125) // 4/pi
|
|
|
|
if d.Equal(NewFromFloat(0.0)) {
|
|
return d
|
|
}
|
|
// make argument positive but save the sign
|
|
sign := false
|
|
if d.LessThan(NewFromFloat(0.0)) {
|
|
d = d.Neg()
|
|
sign = true
|
|
}
|
|
|
|
j := d.Mul(M4PI).IntPart() // integer part of x/(Pi/4), as integer for tests on the phase angle
|
|
y := NewFromFloat(float64(j)) // integer part of x/(Pi/4), as float
|
|
|
|
// map zeros to origin
|
|
if j&1 == 1 {
|
|
j++
|
|
y = y.Add(NewFromFloat(1.0))
|
|
}
|
|
j &= 7 // octant modulo 2Pi radians (360 degrees)
|
|
// reflect in x axis
|
|
if j > 3 {
|
|
sign = !sign
|
|
j -= 4
|
|
}
|
|
z := d.Sub(y.Mul(PI4A)).Sub(y.Mul(PI4B)).Sub(y.Mul(PI4C)) // Extended precision modular arithmetic
|
|
zz := z.Mul(z)
|
|
|
|
if j == 1 || j == 2 {
|
|
w := zz.Mul(zz).Mul(_cos[0].Mul(zz).Add(_cos[1]).Mul(zz).Add(_cos[2]).Mul(zz).Add(_cos[3]).Mul(zz).Add(_cos[4]).Mul(zz).Add(_cos[5]))
|
|
y = NewFromFloat(1.0).Sub(NewFromFloat(0.5).Mul(zz)).Add(w)
|
|
} else {
|
|
y = z.Add(z.Mul(zz).Mul(_sin[0].Mul(zz).Add(_sin[1]).Mul(zz).Add(_sin[2]).Mul(zz).Add(_sin[3]).Mul(zz).Add(_sin[4]).Mul(zz).Add(_sin[5])))
|
|
}
|
|
if sign {
|
|
y = y.Neg()
|
|
}
|
|
return y
|
|
}
|
|
|
|
// cos coefficients
|
|
var _cos = [...]Decimal{
|
|
NewFromFloat(-1.13585365213876817300e-11), // 0xbda8fa49a0861a9b
|
|
NewFromFloat(2.08757008419747316778e-9), // 0x3e21ee9d7b4e3f05
|
|
NewFromFloat(-2.75573141792967388112e-7), // 0xbe927e4f7eac4bc6
|
|
NewFromFloat(2.48015872888517045348e-5), // 0x3efa01a019c844f5
|
|
NewFromFloat(-1.38888888888730564116e-3), // 0xbf56c16c16c14f91
|
|
NewFromFloat(4.16666666666665929218e-2), // 0x3fa555555555554b
|
|
}
|
|
|
|
// Cos returns the cosine of the radian argument x.
|
|
func (d Decimal) Cos() Decimal {
|
|
|
|
PI4A := NewFromFloat(7.85398125648498535156e-1) // 0x3fe921fb40000000, Pi/4 split into three parts
|
|
PI4B := NewFromFloat(3.77489470793079817668e-8) // 0x3e64442d00000000,
|
|
PI4C := NewFromFloat(2.69515142907905952645e-15) // 0x3ce8469898cc5170,
|
|
M4PI := NewFromFloat(1.273239544735162542821171882678754627704620361328125) // 4/pi
|
|
|
|
// make argument positive
|
|
sign := false
|
|
if d.LessThan(NewFromFloat(0.0)) {
|
|
d = d.Neg()
|
|
}
|
|
|
|
j := d.Mul(M4PI).IntPart() // integer part of x/(Pi/4), as integer for tests on the phase angle
|
|
y := NewFromFloat(float64(j)) // integer part of x/(Pi/4), as float
|
|
|
|
// map zeros to origin
|
|
if j&1 == 1 {
|
|
j++
|
|
y = y.Add(NewFromFloat(1.0))
|
|
}
|
|
j &= 7 // octant modulo 2Pi radians (360 degrees)
|
|
// reflect in x axis
|
|
if j > 3 {
|
|
sign = !sign
|
|
j -= 4
|
|
}
|
|
if j > 1 {
|
|
sign = !sign
|
|
}
|
|
|
|
z := d.Sub(y.Mul(PI4A)).Sub(y.Mul(PI4B)).Sub(y.Mul(PI4C)) // Extended precision modular arithmetic
|
|
zz := z.Mul(z)
|
|
|
|
if j == 1 || j == 2 {
|
|
y = z.Add(z.Mul(zz).Mul(_sin[0].Mul(zz).Add(_sin[1]).Mul(zz).Add(_sin[2]).Mul(zz).Add(_sin[3]).Mul(zz).Add(_sin[4]).Mul(zz).Add(_sin[5])))
|
|
} else {
|
|
w := zz.Mul(zz).Mul(_cos[0].Mul(zz).Add(_cos[1]).Mul(zz).Add(_cos[2]).Mul(zz).Add(_cos[3]).Mul(zz).Add(_cos[4]).Mul(zz).Add(_cos[5]))
|
|
y = NewFromFloat(1.0).Sub(NewFromFloat(0.5).Mul(zz)).Add(w)
|
|
}
|
|
if sign {
|
|
y = y.Neg()
|
|
}
|
|
return y
|
|
}
|
|
|
|
var _tanP = [...]Decimal{
|
|
NewFromFloat(-1.30936939181383777646e+4), // 0xc0c992d8d24f3f38
|
|
NewFromFloat(1.15351664838587416140e+6), // 0x413199eca5fc9ddd
|
|
NewFromFloat(-1.79565251976484877988e+7), // 0xc1711fead3299176
|
|
}
|
|
var _tanQ = [...]Decimal{
|
|
NewFromFloat(1.00000000000000000000e+0),
|
|
NewFromFloat(1.36812963470692954678e+4), //0x40cab8a5eeb36572
|
|
NewFromFloat(-1.32089234440210967447e+6), //0xc13427bc582abc96
|
|
NewFromFloat(2.50083801823357915839e+7), //0x4177d98fc2ead8ef
|
|
NewFromFloat(-5.38695755929454629881e+7), //0xc189afe03cbe5a31
|
|
}
|
|
|
|
// Tan returns the tangent of the radian argument x.
|
|
func (d Decimal) Tan() Decimal {
|
|
|
|
PI4A := NewFromFloat(7.85398125648498535156e-1) // 0x3fe921fb40000000, Pi/4 split into three parts
|
|
PI4B := NewFromFloat(3.77489470793079817668e-8) // 0x3e64442d00000000,
|
|
PI4C := NewFromFloat(2.69515142907905952645e-15) // 0x3ce8469898cc5170,
|
|
M4PI := NewFromFloat(1.273239544735162542821171882678754627704620361328125) // 4/pi
|
|
|
|
if d.Equal(NewFromFloat(0.0)) {
|
|
return d
|
|
}
|
|
|
|
// make argument positive but save the sign
|
|
sign := false
|
|
if d.LessThan(NewFromFloat(0.0)) {
|
|
d = d.Neg()
|
|
sign = true
|
|
}
|
|
|
|
j := d.Mul(M4PI).IntPart() // integer part of x/(Pi/4), as integer for tests on the phase angle
|
|
y := NewFromFloat(float64(j)) // integer part of x/(Pi/4), as float
|
|
|
|
// map zeros to origin
|
|
if j&1 == 1 {
|
|
j++
|
|
y = y.Add(NewFromFloat(1.0))
|
|
}
|
|
|
|
z := d.Sub(y.Mul(PI4A)).Sub(y.Mul(PI4B)).Sub(y.Mul(PI4C)) // Extended precision modular arithmetic
|
|
zz := z.Mul(z)
|
|
|
|
if zz.GreaterThan(NewFromFloat(1e-14)) {
|
|
w := zz.Mul(_tanP[0].Mul(zz).Add(_tanP[1]).Mul(zz).Add(_tanP[2]))
|
|
x := zz.Add(_tanQ[1]).Mul(zz).Add(_tanQ[2]).Mul(zz).Add(_tanQ[3]).Mul(zz).Add(_tanQ[4])
|
|
y = z.Add(z.Mul(w.Div(x)))
|
|
} else {
|
|
y = z
|
|
}
|
|
if j&2 == 2 {
|
|
y = NewFromFloat(-1.0).Div(y)
|
|
}
|
|
if sign {
|
|
y = y.Neg()
|
|
}
|
|
return y
|
|
}
|