View Issue Details
ID  Project  Category  View Status  Date Submitted  Last Update 

0033167  FPC  RTL  public  20180213 19:30  20191209 13:09 
Reporter  wp  Assigned To  
Priority  normal  Severity  minor  Reproducibility  have not tried 
Status  new  Resolution  reopened  
Product Version  3.1.1  Product Build  
Target Version  Fixed in Version  3.1.1  
Summary  0033167: Calculation errors of floating point "mod" operation  
Description  The unit math of fpctrunk implements the "mod" operator for floating point numbers. This is coded as operator mod(const a,b:float) c:float;inline; begin c:= ab * Int(a/b); end; However, since the integer part of ratio a/b may be off by 1 the result can be off by one divisor if the result is expected to be close to 0. Example: a = 7.7 b = 1.1 Simple math tells that a mod b should be 0 here because 7.7/1.1 = 7.0. But due to the finite precision of floating point values the ratio actually turns out to be 6.9999999. Together with the Int() operation the result of the mod operation is off by 1 divisor. However, since the remainder of a division is always less than the divisor this error case can be detected by checking c and b for equality with the function SameValue which takes machine precision into account. Furthermore, the function abs should be used to cover the case of negative values. operator mod(const a,b:float) c:float;inline; begin c:= ab * Int(a/b); if SameValue(abs(c), abs(b)) then c := 0.0; end;  
Steps To Reproduce  Run this project: program Project1; uses Math; operator mod(const a,b:Double) c:Double; inline; begin c:= ab * Int(a/b); //if SameValue(abs(c), abs(b)) then c := 0; // Uncomment to see effect of fix end; procedure Test(Dividend, Divisor: double); var res: Double; begin WriteLn('Dividend: ', Dividend:10:5); WriteLn('Divisor : ', Divisor:10:5); WriteLn('Dividend/Divisor = ', Dividend/Divisor); res := Dividend mod Divisor; WriteLn('Result: ', res:10:5); WriteLn; end; begin Test(7.7, 1.1); Test(7.7, 1.1); Test(7.7, 1.1); Test(7.7, 1.1); { Test(7.2, 1.1); Test(7.2, 1.1); Test(7.2, 1.1); Test(7.2, 1.1); } ReadLn; end.  
Additional Information  I am not sure if negative values are handled by the fpc code correctly, anyway. I compared with Wolfram Alpha output: 7.2 mod 1.1 = 0.6 (Wolfram Alpha), same here 7.2 mod 1.1 = 0.5 (Wolfram Alpha), but 0.6 here 7.2 mod 1.1 = 0.5 (Wolfram Alpha), but 0.6 here 7.2 mod 1.1 = 0.6 (Wolfram Alpha), same here To me it looks as if there is no clear definition how to handle negative values: https://en.wikipedia.org/wiki/Modulo_operation  
Tags  No tags attached.  
Fixed in Revision  38332  
FPCOldBugId  
FPCTarget  
Attached Files 


Calling Floor() instead of Int() yields the same signs as Wolfram Alpha, and the abs() can be skipped as well, i.e: operator mod(const a,b:float) c:float; inline; begin c:= ab * Floor(a/b); if SameValue(c, b) then c := 0.0; end; 

Integer case: 72 mod 11=6 72 mod 11=6 So I think we shouldn't follow Wofram Alpha. 

I checked with Mathematica, Wolfram rounds towards negative infinity. That is why there is a difference when the result is negative: In[13]:= 72/11 // N Out[13]= 6.54545 In[14]:= Quotient[72, 11] Out[14]= 7 IMHO, we should follow the currently active FPU rounding mode. 

How does the rounding mode influence the sign? 

It would, if the definition was c:= ab * Round(a/b); But I have since found out from 0030744 that the current Trunc/Roundtozero is libccompatible, so that's better. 

I have written that code and tested against GNU C and Java indeed. 

Please, in this discussion about rounding, don't forget the main concern of this report: the error that "mod" may be off by one divisor if the expected result is 0. 

Fixed as proposed. 

Works only when the operands are "float" or "extended". "double" still yields the old error.  Run the test project given under "Steps to reproduce" (remove the local operator mod). 

Every calculator I have here tells me that Wolfram Alpha is wrong. I used the following and got all the same answers with my fx7400g, Windows calculor, and other ti's etc. I suppose I could dig out my HP, Too. Function myMod(x,y:Currency):Currency; Begin Result := XY*Int(X/Y); end; I don't know what to say, are we so sure Wolfram Alpha is really that correct? 

Note to wp: Float maps to the biggest available float size, i.e the highest available precision:  double on x86_64  extended on i38632  double on ARM and AARCH64 This is the way the math unit handles float. I will run tests against GNU C and Java again to see if there is a incompatibility introduced... I am not sure about the fmod function overloads that were already in math before my patch. Compatibility with the behavior of GNU C(++) Delphi and Java is more important to me than Wolfram/Mathematica. On first sight the solution if SameValue(abs(c), abs(b)) then c := 0.0; should work for all floats except when you pass mixed precision values. In the latter case there need to be overloaded operators for samevalue to work. That's because of the rounding of the least significant bits for different sizes. Java is always double, C++ either single or double. In the case of passing a single and a double the same value test will fail because the single is expanded to double with zero's, not e.g nine's 

Note that in my original patch the results are correct for the *mathematically* correct fact that float mods are always positive, so I used *floor*. However that is not *computationally* correct since in most languages negative values for modulo are allowed. So maybe store sign and use floor() * sign(). I will spend some time on it today an put it on the mailing list to prevent bloat here. Since leftmost determines sign for modulo (that is an existing rule) what I am referring to is: operator mod(const a,b:float) c:float;inline; begin c:= ab * floor(a/b) * Sign(a); end; 

FModfunction has same bug. 

I tested the similar approach for double and single and at least for single it fails (see my attachement). So I think it is better just to limit the mod operator to float or even remove it again altogether, this seems to be a very slippery slope. @delfion: FMod shall be probably delphi compatible, so no need to change it. 

tw33167.pp (1,684 bytes)
{$ifndef FPUNONE} uses math; operator mod(const a,b:single) c:single;inline; begin c:= ab * Int(a/b); if SameValue(abs(c),abs(b)) then c:=0.0; end; { procedure testextended; var a,b,c : extended; begin a:=7.7; b:=1.1; c:=a mod b; if not(SameValue(c,0.0)) then begin writeln(c); halt(1); end; a:=7.2; b:=1.1; c:=a mod b; if not(SameValue(c,0.6)) then begin writeln(c); halt(1); end; a:=7.2; b:=1.1; c:=a mod b; if not(SameValue(c,0.6)) then begin writeln(c); halt(1); end; writeln('ok'); end; procedure testdouble; var a,b,c : double; begin a:=7.7; b:=1.1; c:=a mod b; if not(SameValue(c,0.0)) then begin writeln(c); halt(1); end; a:=7.2; b:=1.1; c:=a mod b; if not(SameValue(c,0.6)) then begin writeln(c); halt(1); end; a:=7.2; b:=1.1; c:=a mod b; if not(SameValue(c,0.6)) then begin writeln(c); halt(1); end; writeln('ok'); end; } procedure testsingle; var a,b,c : single; begin a:=7.7; b:=1.1; c:=a mod b; if not(SameValue(c,0.0)) then begin writeln(c); halt(1); end; a:=7.2; b:=1.1; c:=a mod b; if not(SameValue(c,0.6)) then begin writeln(c); halt(1); end; a:=7.2; b:=1.1; c:=a mod b; if not(SameValue(c,0.6)) then begin writeln(c); halt(1); end; writeln('ok'); end; begin testsingle; // testdouble; // testextended; {$else FPUNONE} begin {$endif FPUNONE} end. 

Florian, If you think so, remove it. It is your child. I merely tried (and tested) a common feature for common languages. And it works as advertised: try the java or c++ code. It is compatible. I already explained that mod isn't mod... But then you must also remove the fmod's which behave exactly the same. 

@delphion: it is not a bug. It is a matter of interpretation consistent with other languages. 

Since it is compatible with major languages I suggest to leave it in. Both mod and fmod. 

> it is not a bug. It is a matter of interpretation consistent with other languages. What are you talking about, the rounding error, or the sign inconsistency? I don't care so much about the sign, but rounding must be correct for all data types, i.e. 7.7 mod 1.1 must be 0 independently of whether the input values are provided a single, double or extended. Otherwise this operator is useless and must be removed. 

Any news on this? I see that current trunk still exhibits the issue reported. 

Function MyFMod(A,B:Double):Double; Var cA,cB:Currency; Begin IF A < 0 Then cA := A else cA := A; If B < 0 Then cB := B else cB := B; Result := cACb*Int(cA/cB); If B < 0 Then Result := Result; End; This yields Zero for your 7.7, 1.1 test/ procedure TForm1.Button1Click(Sender: TObject); begin Caption := MyFMod(7.7,1.1).ToString; end; Testing with other values posted here, they are same except for the two 0.5, with my results its actually 0.5999999999 and its displayed as 0.6 The trick is to move this all into integers and perform the math that way.. I used Currency because of that, I suppose two side by side Int64 per digit can be processed to yield a long fraction it needed. 

To me the compatibility with the GCC compiler suite and Java is essential. I wrote that  original  code based on math.fmod(), nothing more (and suggested sign() as a possible solution). Any other solution would be to adapt fmod too and file reports against the GCC compiler suite. I would keep mod and fmod asis. Modulo is already slow enough and does not need more comparisons. 
Date Modified  Username  Field  Change 

20180213 19:30  wp  New Issue  
20180213 19:51  wp  Note Added: 0106373  
20180213 20:00  wp  Note Edited: 0106373  View Revisions 
20180213 22:00  Florian  Note Added: 0106374  
20180213 22:38  Martok  Note Added: 0106376  
20180213 23:19  Florian  Note Added: 0106377  
20180214 08:32  Martok  Note Added: 0106383  
20180214 08:36  Thaddy de Koning  Note Added: 0106384  
20180214 08:38  Thaddy de Koning  Note Edited: 0106384  View Revisions 
20180214 09:07  wp  Note Added: 0106385  
20180224 19:23  Florian  Fixed in Revision  => 38332 
20180224 19:23  Florian  Note Added: 0106591  
20180224 19:23  Florian  Status  new => resolved 
20180224 19:23  Florian  Fixed in Version  => 3.1.1 
20180224 19:23  Florian  Resolution  open => fixed 
20180224 19:23  Florian  Assigned To  => Florian 
20180226 01:14  wp  Note Added: 0106635  
20180226 01:14  wp  Status  resolved => feedback 
20180226 01:14  wp  Resolution  fixed => reopened 
20180226 03:40  jamie philbrook  Note Added: 0106636  
20180226 03:43  jamie philbrook  Note Edited: 0106636  View Revisions 
20180226 07:46  Thaddy de Koning  Note Added: 0106637  
20180226 07:50  Thaddy de Koning  Note Edited: 0106637  View Revisions 
20180226 07:55  Thaddy de Koning  Note Edited: 0106637  View Revisions 
20180226 08:00  Thaddy de Koning  Note Edited: 0106637  View Revisions 
20180226 08:06  Thaddy de Koning  Note Added: 0106638  
20180226 08:16  Thaddy de Koning  Note Edited: 0106638  View Revisions 
20180226 10:28  delfion  Note Added: 0106642  
20180314 21:43  Florian  Note Added: 0107114  
20180314 21:44  Florian  File Added: tw33167.pp  
20180314 21:44  Florian  Assigned To  Florian => 
20180314 21:50  Thaddy de Koning  Note Added: 0107115  
20180314 21:52  Thaddy de Koning  Note Added: 0107116  
20180314 21:53  Thaddy de Koning  Note Edited: 0107115  View Revisions 
20180314 21:56  Thaddy de Koning  Note Added: 0107117  
20180315 09:07  wp  Note Added: 0107123  
20180315 09:07  wp  Status  feedback => new 
20191206 10:46  wp  Note Added: 0119656  
20191207 03:04  jamie philbrook  Note Added: 0119671  
20191209 13:05  Thaddy de Koning  Note Added: 0119707  
20191209 13:08  Thaddy de Koning  Note Edited: 0119707  View Revisions 
20191209 13:09  Thaddy de Koning  Note Edited: 0119707  View Revisions 