Comparison element by element.

The problem

The reason that your equality test fails is almost certainly because one or both of the variables is not exactly 60. There are many resources on the internet which discuss these issues--like this link.

Floating point numbers typically have approximately 15.5 digits of precision. You can change the formatting of the print statement to show more with the format statement like this:

// Give each number room for 16 digits
// and print out 14 digits after the decimal point
format /rd 16,14;

print w_grid2[N_wage];
print (N_wageOloop1)*(Wage_Truncation/N_wageOloop1);

This should show you that they are different.

The solution

To resolve this problem, you can either truncate the numbers to be compared (if they represent integers), or you can use a fuzzy comparison operator which only checks to see if a certain number of digits are equal. For example:

a = 60.000000001;
b = 60.000000002;

print (a == b);

0.000000

// Compare after zeroing all digits
// after the decimial point
print (trunc(a) == trunc(b));

// Use a fuzzy comparison to check
// only the first 6 digits after the decimal point
print feqmt(a, b, 1e-6);

1.000000
1.000000

Thank you for your answer, indeed the was a 4:

print w_grid2[N_wage];
print (N_wageOloop1)*(Wage_Truncation/N_wageOloop1);

60.00000000000004
60.00000000000000

What is super weird is that the vector with that 4 is generated in the following routine w_grid2[N_wage] is the last element of that vector and the last number should be equal to 60. so I don't really know how the number 60.00000000000004 is generated:

Wage_Truncation = 60;// Truncation for the wage distribution
w_0 = 2 ;// Initial value of wage grid

N_wage=50;
w_grid = seqa(w_0,(Wage_Truncation-w_0)/(N_Wage-1),N_Wage)' ;

w_grid2=w_grid;

The range for double precision floating point numbers, which is used by GAUSS and most scientific computing software, has a range from approximately -1e300 to +1e300. With only 64-bits available (a bit can only be 1 or 0), it is not possible to provide 15.5 digits of precision for every possible number in this space. Not all numbers can be represented exactly. This means that some small rounding error is inherent in floating point math.

Your code starts with 2.0 and adds 1.1836734693877551 to each element, 50 times. There is some rounding that occurs in this calculation, that is why you end up with 60.00000000000004.

I repeated this series of additions in C, GAUSS, Matlab, Python, and R. Each time my computer returned: 60.0000000000000355.

Here is a link to the most famous discussion of floating point arithmetic. You do not need to have a full and deep understanding, but having some awareness of these issues will allow you to avoid some issues and make your work easier.

Thank for your answer!