www.digitalmars.com         C & C++   DMDScript  

digitalmars.D - maybe a floating point issue?

reply czylabsonasa <nobody dev.null> writes:
Hi.

I'm trying to solve a programming contest 
[problem](https://open.kattis.com/problems/anthony) in D. The 
almost identical solution in C++ got accepted, but I'm unable to 
find the issue in the D code... Maybe some D expert catch some 
trivial error in my code?
Here is the D source:

```d
// Anthony_and_Cora_anthony

import std;

alias i32=int;
alias i64=long;
alias f64=double;
alias f128=real;
alias fT=f64;

void main(){
    i32 A,C; readf("%d %d\n",A,C);
    auto p=new fT[](A+C);
    for(i32 i=1;i<=A+C-1;i++){
       readf("%f\n",p[i]);
    }
    auto win=new fT[][](A+1,C+1);
    for(i32 i=0;i<=A;i++){
       for(i32 j=0;j<=C;j++){
          win[i][j]=-1.0;
       }
    }
    for(i32 i=1;i<=A;i++){
       win[i][0]=1.0;
    }
    for(i32 j=0;j<=C;j++){ // j=0 ?
       win[0][j]=0.0;
    }

    void Trav(i32 a,i32 c,i32 k){
       if(win[a][c]<-0.5){
          Trav(a-1,c,k+1);
          Trav(a,c-1,k+1);
          win[a][c]=p[k]*win[a][c-1]+(1-p[k])*win[a-1][c];
       }
    }

    Trav(A,C,1);
    writef("%.7f\n",win[A][C]);
}
```

The approach is a top-down one with memoization. Tried both 
double and real, with no success. What is "surprising" the 
following C++ code is accepted by the judge system:

```c++
// Anthony_and_Cora_anthony
#include <bits/stdc++.h>
using namespace std;

using i32=int;
using f64=double;
using f128=long double;
using fT=f64;

vector<vector<fT>> win;
vector<fT> p;
void Trav(i32 a,i32 c,i32 k){
    if(win[a][c]<-0.5){
       Trav(a-1,c,k+1);
       Trav(a,c-1,k+1);
       win[a][c]=p[k]*win[a][c-1]+(1-p[k])*win[a-1][c];
    }
}

int main(){
    i32 A,C; cin>>A>>C;
    p.resize(A+C+1);
    for(i32 i=1;i<=A+C-1;i++){
       cin>>p[i];
    }
    win.resize(A+1);
    for(i32 i=0;i<=A;i++){
       win[i].resize(C+1);
    }
    for(i32 i=0;i<=A;i++){
       for(i32 j=0;j<=C;j++){
          win[i][j]=-1.0;
       }
    }
    for(i32 i=1;i<=A;i++){
       win[i][0]=1.0;
    }
    for(i32 j=0;j<=C;j++){
       win[0][j]=0.0;
    }

    Trav(A,C,1);
    cout<<setprecision(7)<<fixed<<win[A][C]<<'\n';
}

```

The D code was tested against the C++ code w/ randomly generated 
cases, using dmd,ldc2 and gdc-12 compilers, but even `diff` did 
not give any difference...

cheers.
Sep 18
next sibling parent reply Serg Gini <kornburn yandex.ru> writes:
On Friday, 19 September 2025 at 06:23:02 UTC, czylabsonasa wrote:

 Hi.

 I'm trying to solve a programming contest 
 [problem](https://open.kattis.com/problems/anthony) in D. The


Let's see:

    auto p=new fT[](A+C);


In D code you allocating p vector of size A+C


    p.resize(A+C+1);


In C++ code you are allocating p vector of size A+C+1

Maybe on some corner cases of the logic Trav function it can 
affect the result
Sep 19
parent reply czylabsonasa <nobody dev.null> writes:
On Friday, 19 September 2025 at 07:04:16 UTC, Serg Gini wrote:

 On Friday, 19 September 2025 at 06:23:02 UTC, czylabsonasa 
 wrote:

 Hi.

 I'm trying to solve a programming contest 
 [problem](https://open.kattis.com/problems/anthony) in D. The


 Let's see:

    auto p=new fT[](A+C);


 In D code you allocating p vector of size A+C


    p.resize(A+C+1);


 In C++ code you are allocating p vector of size A+C+1

 Maybe on some corner cases of the logic Trav function it can 
 affect the result


Thanks, Serg, good catch, usually i set the limits a little bit 
larger than the requirements, but this time i missed for the D 
code.
To be sure (the devil isn't sleeping...) tried w/ `A+C+1` sized 
allocation for the D, but no success, and `A+C` sized for the C++ 
- still accepted. (only the 1,...,A+C-1 indices are accessed).

cheers
Sep 19
parent reply czylabsonasa <nobody dev.null> writes:
... and even the following `julia` code is OK:

```jl


const _DBG_=false

fT=Float64

function solveIt()
    rL()=readline()
    rT(S,T=Int)=parse(T,S)

    A,C=rT.(rL()|>split)
    p=fill(fT(0.0),A+C-1)
    for i in 1:A+C-1
       p[i]=rT(rL(),fT)
    end

    for i in 1:A
       win[1+i,1+0]=1.0
    end
    for j in 0:C
       win[1+0,1+j]=0.0
    end
    function Trav(a,c,k)
       if win[1+a,1+c]<-0.5
          Trav(a-1,c,k+1)
          Trav(a,c-1,k+1)
          win[1+a,1+c]=p[k]*win[1+a,1+c-1]+(1-p[k])*win[1+a-1,1+c]
       end
    end

    Trav(A,C,1)

    println(round(win[1+A,1+C];digits=7))
end

solveIt()


```
(sorry, but no syntax highlighting here for julia)

As i see, the kattis site is using the `-O2` switch for compiling 
the source, as usual.
In the old times - 10 years ago - i experienced striking 
differences w/ gcc/g++ w/
floating point computations depending on optimization is enabled 
or not. Because eventually the 3 codes are the same, i am 
supposing that this the root cause of the issue.
Sep 19
parent reply Serg Gini <kornburn yandex.ru> writes:
On Friday, 19 September 2025 at 08:31:07 UTC, czylabsonasa wrote:

 ... and even the following `julia` code is OK:
 floating point computations depending on optimization is 
 enabled or not. Because eventually the 3 codes are the same, i 
 am supposing that this the root cause of the issue.


And you said you made some random tests - and results are the 
same?
even for 7 digits precision?

It's weird
Sep 19
parent czylabsonasa <nobody dev.null> writes:
 And you said you made some random tests - and results are the 
 same?
 even for 7 digits precision?

 It's weird


Yes, i tested the C++ and D code with a thousand (literally) of 
random instances and got no difference, even comparing the 
outputs as strings.
(For compiling g++ -O2 and gdc-12 -O2 were used in accordance 
with the kattis site). On my system (lmde/faye) only gdc-12 is 
available, and i don't want to "mess" it with a different 
compiler.
Sep 19
prev sibling next sibling parent reply Denis Feklushkin <ddd ddd.dddd> writes:
On Friday, 19 September 2025 at 06:23:02

 ```d
 // Anthony_and_Cora_anthony

 import std;

 alias i32=int;
 alias i64=long;
 alias f64=double;
 alias f128=real; // may be is not same as long double in C++


```
Sep 19
parent czylabsonasa <nobody dev.null> writes:
On Friday, 19 September 2025 at 10:06:05 UTC, Denis Feklushkin 
wrote:

 On Friday, 19 September 2025 at 06:23:02

 ```d
 // Anthony_and_Cora_anthony

 import std;

 alias i32=int;
 alias i64=long;
 alias f64=double;
 alias f128=real; // may be is not same as long double in C++


 ```


Denis, the kattis site probably using linux+gdc-14, and i tested 
my linux+gdc-12 compiled code gives ca. 2^(-16445) for the 
smallest pos. real, which shows that it is not double, at least 
on my platform. But i got accepted using double with C++ and 
Float64 with julia, so
the issue is hiding is somewhere else.
Sep 19
prev sibling next sibling parent reply Dennis <dkorpel gmail.com> writes:
On Friday, 19 September 2025 at 06:23:02 UTC, czylabsonasa wrote:

 I'm trying to solve a programming contest 
 [problem](https://open.kattis.com/problems/anthony) in D. The 
 almost identical solution in C++ got accepted, but I'm unable 
 to find the issue in the D code...


One difference between the C++ and D code is that your D code 
expects each line to end with a newline character:

```D
    for(i32 i=1;i<=A+C-1;i++){
       readf("%f\n",p[i]);
    }
```

So if your input file is this:

```
1 1
0.500000

```

Your C++ and D programs both output 0.500000. But if your input 
file lacks the final empty line:

```
1 1
0.500000
```

The C++ version still works, but D throws an exception:

```
std/format/spec.d(569): parseToFormatSpec: Cannot find character '
' in the input string.
```

Perhaps change the D code to this:

```D
    for(i32 i=1;i<=A+C-1;i++){
       p[i] = readln().to!float;
    }
```
Sep 19
parent czylabsonasa <nobody dev.null> writes:
On Friday, 19 September 2025 at 10:25:37 UTC, Dennis wrote:

 On Friday, 19 September 2025 at 06:23:02 UTC, czylabsonasa 
 wrote:

 I'm trying to solve a programming contest 
 [problem](https://open.kattis.com/problems/anthony) in D. The 
 almost identical solution in C++ got accepted, but I'm unable 
 to find the issue in the D code...


 One difference between the C++ and D code is that your D code 
 expects each line to end with a newline character:

 ```D
    for(i32 i=1;i<=A+C-1;i++){
       readf("%f\n",p[i]);
    }
 ```

 So if your input file is this:

 ```
 1 1
 0.500000

 ```

 Your C++ and D programs both output 0.500000. But if your input 
 file lacks the final empty line:

 ```
 1 1
 0.500000
 ```

 The C++ version still works, but D throws an exception:

 ```
 std/format/spec.d(569): parseToFormatSpec: Cannot find 
 character '
 ' in the input string.
 ```

 Perhaps change the D code to this:

 ```D
    for(i32 i=1;i<=A+C-1;i++){
       p[i] = readln().to!float;
    }
 ```


Denis, thanks for your time, but if there would such an input 
file, then it would resulted in runtime error. Anyway, i tested 
with readf(" %d %d",...)+read(" %f",...) and no change it is 
wrong answer.
Sep 19
prev sibling next sibling parent reply matheus <matheus gmail.com> writes:
On Friday, 19 September 2025 at 06:23:02 UTC, czylabsonasa wrote:

 Hi.

 I'm trying to solve a programming contest 
 [problem](https://open.kattis.com/problems/anthony) in D...


Hi, have you tried checking the Assembly code between both 
version with godbolt.org ?

Matheus.
Sep 19
parent czylabsonasa <nobody dev.null> writes:
 Hi, have you tried checking the Assembly code between both 
 version with godbolt.org ?

 Matheus.


Matheus, i'm using programming langs as a black box :-), i know 
almost nothing about assembly, but thanks for the suggestion!
Sep 19
prev sibling parent reply Timon Gehr <timon.gehr gmx.ch> writes:
On 9/19/25 08:23, czylabsonasa wrote:

 Hi.
 
 I'm trying to solve a programming contest [problem](https:// 
 open.kattis.com/problems/anthony) in D. The almost identical solution in 
 C++ got accepted, but I'm unable to find the issue in the D code... 
 Maybe some D expert catch some trivial error in my code?
 Here is the D source:
 
 ```d
 // Anthony_and_Cora_anthony
 
 import std;
 
 alias i32=int;
 alias i64=long;
 alias f64=double;
 alias f128=real;
 alias fT=f64;
 
 void main(){
     i32 A,C; readf("%d %d\n",A,C);
     auto p=new fT[](A+C);
     for(i32 i=1;i<=A+C-1;i++){
        readf("%f\n",p[i]);
     }
     auto win=new fT[][](A+1,C+1);
     for(i32 i=0;i<=A;i++){
        for(i32 j=0;j<=C;j++){
           win[i][j]=-1.0;
        }
     }
     for(i32 i=1;i<=A;i++){
        win[i][0]=1.0;
     }
     for(i32 j=0;j<=C;j++){ // j=0 ?
        win[0][j]=0.0;
     }
 
     void Trav(i32 a,i32 c,i32 k){
        if(win[a][c]<-0.5){
           Trav(a-1,c,k+1);
           Trav(a,c-1,k+1);
           win[a][c]=p[k]*win[a][c-1]+(1-p[k])*win[a-1][c];
        }
     }
 
     Trav(A,C,1);
     writef("%.7f\n",win[A][C]);
 }
 ```
 
 The approach is a top-down one with memoization. Tried both double and 
 real, with no success. What is "surprising" the following C++ code is 
 accepted by the judge system:
 
 ```c++
 // Anthony_and_Cora_anthony
 #include <bits/stdc++.h>
 using namespace std;
 
 using i32=int;
 using f64=double;
 using f128=long double;
 using fT=f64;
 
 vector<vector<fT>> win;
 vector<fT> p;
 void Trav(i32 a,i32 c,i32 k){
     if(win[a][c]<-0.5){
        Trav(a-1,c,k+1);
        Trav(a,c-1,k+1);
        win[a][c]=p[k]*win[a][c-1]+(1-p[k])*win[a-1][c];
     }
 }
 
 int main(){
     i32 A,C; cin>>A>>C;
     p.resize(A+C+1);
     for(i32 i=1;i<=A+C-1;i++){
        cin>>p[i];
     }
     win.resize(A+1);
     for(i32 i=0;i<=A;i++){
        win[i].resize(C+1);
     }
     for(i32 i=0;i<=A;i++){
        for(i32 j=0;j<=C;j++){
           win[i][j]=-1.0;
        }
     }
     for(i32 i=1;i<=A;i++){
        win[i][0]=1.0;
     }
     for(i32 j=0;j<=C;j++){
        win[0][j]=0.0;
     }
 
     Trav(A,C,1);
     cout<<setprecision(7)<<fixed<<win[A][C]<<'\n';
 }
 
 ```
 
 The D code was tested against the C++ code w/ randomly generated cases, 
 using dmd,ldc2 and gdc-12 compilers, but even `diff` did not give any 
 difference...
 
 cheers.
 
 



FWIW: if we improve the ability of the optimizer to do alias analysis, 
it passes:


```d
import std;

alias i32=int;
alias i64=long;
alias f64=double;
alias f128=real;
alias fT=f64;

void main(){
    i32 A,C; readf("%d %d\n",A,C);

    auto p=new fT[](A+C);
    for(i32 i=1;i<=A+C-1;i++){
       readf("%f\n",p[i]);
    }

    const i32 stride=C+1;
    auto win=new fT[]((A+1)*(C+1));

    for(i32 i=0;i<=A;i++){
       for(i32 j=0;j<=C;j++){
          win[i*stride+j]=-1.0;
       }
    }
    for(i32 i=1;i<=A;i++){
       win[i*stride+0]=1.0;
    }
    for(i32 j=0;j<=C;j++){
       win[0*stride+j]=0.0;
    }

    void Trav(i32 a,i32 c,i32 k){<
       const idx=a*stride+c;
       if(win[idx]<-0.5){
          Trav(a-1,c,k+1);
          Trav(a,c-1,k+1);
          win[idx]=p[k]*win[a*stride+(c-1)]+(1-p[k])*win[(a-1)*stride+c];
       }
    }

    Trav(A,C,1);
    writef("%.7f\n",win[A*stride+C]);
}
```

Not sure why this happens with GDC on -O2, I think LDC would refrain 
from changing results by default.
Sep 19
parent reply czylabsonasa <nobody dev.null> writes:
 Not sure why this happens with GDC on -O2, I think LDC would 
 refrain from changing results by default.


Timon, good finding! Then, with D, on the kattis site we should 
begin with the second step :-) I'm not happy about it, bcos there 
are plenty of interesting geometry/probability tasks on that 
site, where floating point numbers and (naturally) higher dim. 
arrays involved,
instead of focusing the problem, we must fight with 
language/compiler.
Sep 19
parent Sergey <kornburn yandex.ru> writes:
On Friday, 19 September 2025 at 20:33:31 UTC, czylabsonasa wrote:

 Not sure why this happens with GDC on -O2, I think LDC would 
 refrain from changing results by default.


 Timon, good finding! Then, with D, on the kattis site we should 
 begin with the second step :-) I'm not happy about it, bcos 
 there are plenty of interesting geometry/probability tasks on 
 that site, where floating point numbers and (naturally) higher 
 dim. arrays involved,
 instead of focusing the problem, we must fight with 
 language/compiler.


Maybe we should submit a ticket with request to change the 
compiler to ldc?
I saw they are using many GNU languages (C,C++,Fortran,Modula-2).

But they also have LLVM langs (Rust, Crystal).
So it should not be a huge problem to install LDC instead of GDC 
in their environment.
Sep 19