Complete the implement `opBinary` between two Fixed values

[mw66]

right now, we can only divide by isIntegral or isFloatingPoint. https://github.com/m3m0ry/fixedpoint/blob/master/fixedpoint/fixed.d#L256 https://github.com/m3m0ry/fixedpoint/blob/master/fixedpoint/fixed.d#L268

we should be able to divide money amounts of the same type, e.g. Bill Gates is 123456.789 x richer than me.

The return type is double.

(same issue here: https://github.com/qznc/d-money/issues/9 )

[m3m0ry]

(Your examples are really hard to understand. Instead of "Bill Gates" could you provide a code example: Fixed("123456.789") * Fixed("0.001") -> ...)

Your request is: Complete the implement opBinary between two Fixed values. +- of same scaling should work. */ of same scaling and */+- of different scaling are not yet implemented.

I have been thinking about it for some time now.

Why are only +- of the same scaling implemented? Because they are trivial. There is no change of the scaling needed, i.e. we are not losing precision adding or subtracting those two numbers. You might argue that there is also no checking of an overflow and you are probably partially right. However, an overflow could be expected, whereas slowly loosing precision is hard to find-out, track or debug. (also you can now use BigInt as a template for the value, so you can never have overflow. This remains to be tested.)

Implementation of remaining opBinary functions: +- of two Fixed values of different scaling: This should return a Fixed value of the higher scaling of the two: Fixed("0.1") + Fixed("0.01") == Fixed("0.11") rather then == Fixed("0.1").

*/ of two Fixed values, independent of scaling might or might not need more precision. In the worst case this is this.scaling * other.scaling, which can be quite large pretty quickly and is thus not feasible. One solution would be to take the more precise value and loose precision. Your proposal is to simply return double, which might loose precision as well (although significantly small) and introduce and maybe unexpected type in the middle of an expression.

Honestly i don't know what the better approach is. I think your proposal could be better then mine.

[m3m0ry]

I'll look into it tonight.

[mw66]

(Your examples are really hard to understand. Instead of "Bill Gates" could you provide a code example: Fixed("123456.789") * Fixed("0.001") -> ...)

Oh, I may not make myself clear, it's not about rounding, but is about provide such functionality: There are use cases where we need this, e.g. after calculating the fixedpoint value using this library, the value needs to feed into another library for some math operation (or for drawing a chart, i.e value => pixel height), which takes D native type double. In this charting example, I don't care the precision. But I need this functionality to be able to cast to double.

the implementation I added to the other package: https://github.com/qznc/d-money/blob/master/source/money.d#L189

        return (amount.to!double) / (rhs.amount.to!double);

in your package, it will be:

        return (this.value.to!double) / (rhs.value.to!double);

[m3m0ry]

I think it is weird that Fixed / Fixed will produce a double and the computation is in double as well. Fixed op Fixed has to produce an outcome of type Fixed. I understand that you might like it, but others would not.

If you need to compute with doubles, you can do it as well: Fixed!3(4.5).to!double / Fixed!4(5.6).to!double.

I have also pushed to master with now working */ for Fixed objects of the same type. They just keep the same precision, i.e. loose information.

Btw. does your method work? With return (amount.to!double) / (rhs.amount.to!double); you are loosing all information about scaling, or "dec_places". What does currency!"eur"(100) / currency!"eur"(0.01) produce?

[mw66]

/+dub.sdl:
dependency "money" version="~>2.5.0"
+/

import std.stdio;
import money;

alias EUR = currency!("EUR");

void main() {
  EUR x = EUR( 1.23456 );
  EUR y = EUR("1.23456");
  //writeln(x, x.sizeof);
  //writeln(y, x.sizeof);
  //assert(x == y);

  auto n = currency!"eur"(100) / currency!"eur"(0.01);
  writeln(n);
}

$ dub m.d
10000

[mw66]

Fixed op Fixed has to produce an outcome of type Fixed.

I see what you mean.

Right, the other library has 'unit' attached to it, after division of the same currency, the result is unit-less.

Your library has no 'unit', so it's fine to return Fixed.

[m3m0ry]

Ah, I have not seen that part. Sorry for confusion.

I hope my answer and implementation at least of the Fixed operator satisfied you.

[mw66]

Thanks.

[mw66]

FYI, I do encounter a problem due to Fixed / Fixed return Fixed, instead of double:

As I'm draw chart, need to calculate pixels (int), so I need to convert a (Fixed / Fixed) => int, but for my domain, I need to take the ceiling, the two number are:

alias decimal = Fixed!4;

a = decimal("43.5001")
b = decimal( "4.3500")

if we return double:

a / b = 10.00002298850574712643, and
ceil(a / b) = 11

Since we return Fixed,

a / b = 10
ceil (a/b) = 10

(11 - 10) = 1 of my drawing unit (x-number-of-pixels), so the graph doesn't look right.

I'm wondering if we really need to make the / return Fixed. In mathematics, the decimal defined in Fixed is actually an integer, integer operations are closed for + - *, but not /. Only rationals are closed under /.

https://en.wikipedia.org/wiki/Integer#Algebraic_properties

Algebraic properties

Like the natural numbers, Z is closed under the operations of addition and multiplication, that is, the sum and product of any two integers is an integer. However, with the inclusion of the negative natural numbers, and, importantly, 0, Z (unlike the natural numbers) is also closed under subtraction. The integers form a unital ring which is the most basic one, in the following sense: for any unital ring, there is a unique ring homomorphism from the integers into this ring. This universal property, namely to be an initial object in the category of rings, characterizes the ring Z.

Z is not closed under division, since the quotient of two integers (e.g., 1 divided by 2), need not be an integer. Although the natural numbers are closed under exponentiation, the integers are not (since the result can be a fraction when the exponent is negative).

[m3m0ry]

(Now this is an awesome issue report :) )

I am still stubborn and any operation between two Fixed needs to produce a Fixed type as well. This is not about algebra, this is about D. Any operation between two int results into an int, yes loosing information when dividing.

And this is also exactly what I meant by "loosing precision". If you multiply both numbers with 10**4 to get the integer, do the operation and then divide with 10**4 (because that what Fixed is doing for you), you will get the same result. That is what I wanted to achieve.

However I understand your point of you and problem. And it got me thinking, that I return Fixed, when I have operation with a floating point number, e.g.: Fixed and double. I am not sure if this is the smartest way to do. If I would return the floating point number, you could write something like this: ceil(1.0 * a / b) == 11 or ceil(a.to!double /b) == 11. Would that be satisfactory? It is certainly shorter then ceil(a.to!double / b.to!double) == 11.

[m3m0ry]

Another solution would be to provide a function mul and div, which would keep the information, e.g.: Fixed!2.div(Fixed!2) -> Fixed!4.

[m3m0ry]

The current commit and version has these changes. Operation with floating point numbers will be a floating point number.

[mw66]

I think if we step back: the raw operation to make use of the stored internal representation (without lose unnecessary precision) is:

  return ceil(cast(double)(a.value) / cast(double)(b.value));

To achieve this, there are several ways:

1) let the user do it explicitly outside the library, just as the above code, fortunately, Fixed is a struct, and value can be accessed by user.

2) the whole operation is provided internally by the package Fixed, which means Fixed need to provide all these ceil, floor on those intermediate result (cast(double)(a.value) / cast(double)(b.value) internally, which is not a good design, and also unncessary.

3) provide a template function div(T), let the user explicitly pass in the return type the user want to use, e.g.

T div(T)(typeof(this) b) {
  return (cast(T)(this.value) / cast(T)(b.value));  // directly use the value, to preserve as much precision as possible
}

and the user can call like this:

  Fixed!4 a, b;
  int n = a.div!int(b);
  double d = a.div!double(b);

This is most flexible solution I think.

[m3m0ry]

I have added an arithmetic module. What do you think?

[mw66]

You directly cast lhs & rhs will lost precision, by the factor

https://github.com/m3m0ry/fixedpoint/blob/master/fixedpoint/fixed.d#L231

what I'm proposing is: "rhs.value.to!T" ~ op ~ "lhs.value.to!T"; then (/ factor) for +/-; and (/ factor / factor) for *; but nothing need to do for "/" to preserver as much precision as possible.

Added comment in your PR.

[mw66]

https://github.com/m3m0ry/fixedpoint/commit/f5c4177c9a9629fb2ac8c8d33e5c51b7cad4f618#comments

[m3m0ry]

Yes but that only works for fixed numbers of same scaling.

I think I have now a solution even for operations between two different Fixed types.

[mw66]

Yes but that only works for fixed numbers of same scaling.

I think I have now a solution even for operations between two different Fixed types.

We should try our best to maximize the precision, for different scaling: e.g. Fixed!4 / Fixed!2, we can do:

  f4.value / f2.value / 100  
  // instead of cast(double)f4 / cast(double)f2; 

i.e. always work with the value, which keep the precision.

[m3m0ry]

I don't understand. Until now, you wanted to cast first and then do the operation. I don't see how it increases precision.

[mw66]

every operation may lost precision; so 1) we want use as less op as possible, by remove the unnecessary ops. 2) always use the original value first.

fixed.value holds the original value, let's consider "/".

cast(double)(a.value) / cast(double)(b.value)    // 2 cast + 1 div = 3 op in total

for each cast(double)(fixed), internally it need to do:

 cast(double)(fixed.value) / cast(double)factor    // 2 cast + 1 div = 3 op in total

so:

  cast(double)(a) / cast(double)(b)   // (2 x 3 internal op) + 1 div = 7 op in total

Let's try the following example: cat precision.d

/+dub.sdl:
dependency "fixedpoint" version="~>0.1.0"
+/

import std.format;
import std.stdio;
import fixedpoint.fixed : Fixed;

void main() {
  alias decimal = Fixed!8;

  auto a = decimal("43.50011234");
  auto b = decimal( "4.3500");

  double d1 = (cast(double)a) / (cast(double)b);

  double d2 = (cast(double)(a.value)) / (cast(double)(b.value));

  writeln(format("%.20f", d1));
  writeln(format("%.20f", d2));
}

$ dub precision.d
10.00002582528735750600
10.00002582528735572964

The result differs, I would trust the 2nd result, since it has the least ops.

[m3m0ry]

Ok, now I understand.

But I am not sure what part of the code you mean. The divTo uses a.value.to!double. (Notice, that I didn't do a version bump yet)

[mw66]

e.g.

https://github.com/m3m0ry/fixedpoint/blob/master/fixedpoint/fixed.d#L296

        return Fixed.make(mixin("((value * factor" ~ op ~ "rhs.value * factor) / factor).round.to!V"));

5 ops before passing to Fixed.make(...)

BTW: value * factor itself is a bug, I think you mean:

        return Fixed.make(mixin("((value / factor" ~ op ~ "rhs.value / factor) * factor).round.to!V"));

while this can be simplified to:

        return Fixed.make(mixin("((value" ~ op ~ "rhs.value) * factor).round.to!V"));

3 ops before passing to Fixed.make(...)

Numerical stability is a big topic, if you are interested:

https://en.wikipedia.org/wiki/Round-off_error https://en.wikipedia.org/wiki/Numerical_analysis#Numerical_stability_and_well-posed_problems https://en.wikipedia.org/wiki/Numerical_stability

If you are interested.

[mw66]

Also I think

https://github.com/m3m0ry/fixedpoint/blob/master/fixedpoint/fixed.d#L86

    /// Direct construction of a Fixed struct
    static pure nothrow Fixed make(const V v)
    {
        Fixed fixed;
        fixed.value = v;
        return fixed;
    }

better decl this func as private, since it directly assign to value, the outside client code can easily make mistakes.

[mw66]

very interesting example:

https://en.wikipedia.org/wiki/Numerical_analysis#Numerical_stability_and_well-posed_problems

https://wikimedia.org/api/rest_v1/media/math/render/svg/bd01a68a0d14cee9d71920b6a63169d72c1301f1

Numerical stability is affected by the number of the significant digits the machine keeps on, if a machine is used that keeps only the four most significant decimal digits, a good example on loss of significance can be given by these two equivalent functions {\displaystyle f(x)=x\left({\sqrt {x+1}}-{\sqrt {x}}\right){\text{ and }}g(x)={\frac {x}{{\sqrt {x+1}}+{\sqrt {x}}}}.} f(x)=x\left(\sqrt{x+1}-\sqrt{x}\right) \text{ and } g(x)=\frac{x}{\sqrt{x+1}+\sqrt{x}}. Comparing the results of {\displaystyle f(500)=500\left({\sqrt {501}}-{\sqrt {500}}\right)=500\left(22.38-22.36\right)=500(0.02)=10} f(500)=500 \left(\sqrt{501}-\sqrt{500} \right)=500 \left(22.38-22.36 \right)=500(0.02)=10 and {\displaystyle {\begin{alignedat}{3}g(500)&={\frac {500}{{\sqrt {501}}+{\sqrt {500}}}}\&={\frac {500}{22.38+22.36}}\&={\frac {500}{44.74}}=11.17\end{alignedat}}} \begin{alignat}{3}g(500)&=\frac{500}{\sqrt{501}+\sqrt{500}}\ &=\frac{500}{22.38+22.36}\ &=\frac{500}{44.74}=11.17 \end{alignat}

[m3m0ry]

Yeah, I was not happy with that operator either, but at that day, I was too tired to bother with it.

The solution is actually:

 return Fixed.make(mixin("((value * factor") ~ op ~ "rhs.value).round.to!V"));

Counter Example: Take 2 Fixed!2 numbers: 5.23 / 7.10 which is 0.73661971831 according to my calculator. (523 * 10^^2) / 710 == 73 (because everything is int), which is 0.73 in "human readable" format. Error is: 0.0066197183099

(523 / 710) * 10^^2 == 0 (because everything is int) Error is: 0.7366197183099

Integers except for divisions and overflows have no problems with stability.

Regarding the private/public think. I do lot of things in Python, so maybe I am affected by that, but I don't see a reason, why to make it private. It could be useful. Also useful for me, if I ever expand on that arithmetic module.

Thank you for forcing me into thinking about it! I really appreciate your reviews.

[mw66]

@schveiguy FYI:

I still have problems, which I haven't fully debugged, but I suspect the issue could caused by: Fixed / Fixed => Fixed.

https://github.com/m3m0ry/fixedpoint/blob/master/fixedpoint/fixed.d#L294

As I noted above :

the decimal defined in Fixed is actually an integer (long value), integer operations are closed for + - *, but not /. Only rationals are closed under /.

https://en.wikipedia.org/wiki/Integer#Algebraic_properties

Algebraic properties

Like the natural numbers, Z is closed under the operations of addition and multiplication, that is, the sum and product of any two integers is an integer. However, with the inclusion of the negative natural numbers, and, importantly, 0, Z (unlike the natural numbers) is also closed under subtraction. The integers form a unital ring which is the most basic one, in the following sense: for any unital ring, there is a unique ring homomorphism from the integers into this ring. This universal property, namely to be an initial object in the category of rings, characterizes the ring Z.

Z is not closed under division, since the quotient of two integers (e.g., 1 divided by 2), need not be an integer. Although the natural numbers are closed under exponentiation, the integers are not (since the result can be a fraction when the exponent is negative).

[mw66]

@schveiguy this is a long discussion, you can walk thru it when you have time.

[mw66]

BTW, I know divTo!double is provided, but it's just inconvenient to use.

https://github.com/m3m0ry/fixedpoint/blob/master/fixedpoint/arithmetic.d#L35

[mw66]

One more issue:

https://github.com/m3m0ry/fixedpoint/blob/master/fixedpoint/fixed.d#L258

    T opBinary(string op, T)(const T rhs) const
            if ((op == "+" || op == "-" || op == "*" || op == "/") && isFloatingPoint!T)
    {
        return mixin("this.to!T" ~ op ~ "rhs");
    }

This computation order is not optimal, suppose T is float or double: when this.to!T is calculated first, it may lose some precision already, then that intermediate result is op-ed with rhs (which is the user specified value), will lose more precision than this order:

( this.value `op` rhs ) / this.factor 

i.e. the general principle is:

let the original values of two operands op (meet) first, then divide by the const (un-significant value), as I noted previously:

https://github.com/m3m0ry/fixedpoint/issues/5#issuecomment-639995618

There are maybe other places in the code which have this order problem, please give it a check when you have time.

Numerical stability is very subtle to get it right.

[m3m0ry]

On thursday I have a day off. I'll look into it then. I know i have right now a "working" Fixed struct, and there is room for improvement.

[schveiguy]

I feel like multiplication and division between Fixed should result in another Fixed. I specifically punted on multiplication and division between Fixed because it's unclear what people may want as a result.

Division is difficult, because the actual result may not be representable

The way I would handle it is to make the division between the two integers, and then figure out the factor. The rounding mode could be specified for the last digit either as a parameter to the call or a parameter to the template. I would default to truncation since these are in essence integer types.

Note that overflows are going to come easily in these operations.

[m3m0ry]

I still try to keep Fixed similar to integral types. So any operation between two Fixed results into Fixed. If two Fixed with different template come along, i will take the higher/bigger version. Similar to int op long -> long. Yes truncation might happen, like with integral types. I'll think on Thursday about my operators, if by restructuring the operation would result into less likely overflow.

Note that I do have tests, but for a solid library i would need WAY more.

[m3m0ry]

One more issue:

https://github.com/m3m0ry/fixedpoint/blob/master/fixedpoint/fixed.d#L258

    T opBinary(string op, T)(const T rhs) const
            if ((op == "+" || op == "-" || op == "*" || op == "/") && isFloatingPoint!T)
    {
        return mixin("this.to!T" ~ op ~ "rhs");
    }

This computation order is not optimal, suppose T is float or double: when this.to!T is calculated first, it may lose some precision already, then that intermediate result is op-ed with rhs (which is the user specified value), will lose more precision than this order:

( this.value `op` rhs ) / this.factor 

i.e. the general principle is:

let the original values of two operands op (meet) first, then divide by the const (un-significant value), as I noted previously:

#5 (comment)

There are maybe other places in the code which have this order problem, please give it a check when you have time.

Numerical stability is very subtle to get it right.

First of all, this would only work for */. *- need:

( this.value `op` rhs*factor ) / factor 

Let's look at */ operators. Case double. Current version:

this.to!double * rhs
-> ((this.value.to!double / factor).to!double) * rhs
-> (int.to!double / int).to!double * double 
                             \\ the to!double outside does nothing, since inside is already a double
                             \\ int needs to be cast implicitely
-> (int.to!double / int.to!double) * double
----> (value.to!double / factor.to!double) * rhs

Your version:

( this.value * rhs ) / this.factor
-> ( int * double ) / int                \\ implicit cast to double
-> (int.to!double * double) / int.to!double
-----> (value.to!double * rhs) / factor.to!double

All operations are between doubles and are */. Without knowing the numbers, I cannot tell which version is more stable.

[schveiguy]

Multiplication and division are simple as long as you are fine with integer truncation:

Fixed!(scaling + F) opBinary(string op : "*", int F)(Fixed!F other) {
    return Fixed!(scaling + F).make(value * rhs.value);
}

Fixed!(scaling - F) opBinary(string op : "/", int F)(Fixed!F other) {
    return Fixed!(scaling - F).make(value / rhs.value);
}

But the reason I originally left it out is because this may not be what you want (to adjust the scaling factor).

I think fixedpoint should avoid all vagaries of floating point, that should be something you opt into by asking for it specifically.

[m3m0ry]

The current implementation of binary operations between two different Fixeds looks, which of the two operands have the bigger scaling, and chooses this as the common type. It then casts the smaller one to the bigger one and does the operation. This is similar to int + long. Int also gets casted to the bigger one, i.e. long. This will be the default, since it cannot surprise anyone. Yes truncation might happen.

(Automatically) Adjusting the scaling factor should be a function, which would be in the arithmetic module. I'll implement it (some day :smile: )

[schveiguy]

int * long is not the same -- no data is lost.

1 1L => 1L Fixed!2("0.01") Fixed!2("0.01") => Fixed!2(0) Correct:

Fixed!2("0.01") * Fixed!2("0.01") => Fixed!4("0.0001")

[m3m0ry]

That is not an analogical example. The equivalent would be something like: short.max + int.max. As far as I know, the result will be an int and it will overflow, producing something negative. It won't be magically a long so no information is lost.

Also Fixed!2("10.00") * Fixed!2("10.00") => Fixed!4("10.0000") Doing this a few times will result into an overflow pretty quickly.

The questions is, what should Fixed provide as opBinary.

1. Choose the biggest scaling from the operands (current version)
2. Ensure that no information is lost, e.g. add scaling during a multiplication and risk an overflow
3. Something more complicated (choose scaling during runtime based on values)

I thought 1. is the most simple, straight forward one. Something which we are used to.

[schveiguy]

short.max + int.max is not equivalent.

A better example:

Fixed!2("10.01") * Fixed!2("10.01")

The true answer is 100.2001. But Fixed!2 will only represent 100.2.

short.max + int.max is going to result in something COMPLETELY unrelated (and I think actually will be negative). It's not just losing data from the result of the operation, like the 0.0001.

If someone is using fixed, they possibly prefer not to lose data. Otherwise they would use doubles.

It's also possible they only ever want to use Fixed!2, which means they don't care about the decimals dropping off. I'd consider this last case as equivalent to short + short, which results in int, and you have to cast to declare "I don't care about the possible problems".

Doing this a few times will result into an overflow pretty quickly.

Yes, which is why it's hard to determine which is more important. The solution could be adjusting the scale manually after the math is done.

I thought 1. is the most simple, straight forward one. Something which we are used to.

IMO, the most straightforward one is what I wrote. It's one multiplication or division. Not only that, but you don't have to worry about an overflow poisoning the result (if the result is representable by the type, it will be valid).

We may want to consider providing ways to specify the behavior (like std.experimental.checkedint), so we can support all use cases.

[m3m0ry]

I like the idea with hooks. I'll look into it.

[m3m0ry]

I have implemented the hooks for binary operations. However one of the hooks is not yet tested (the one, with your proposal @schveiguy )