How to limit calculation precision?

Dear Programmers,

I am working with a computational commercial software package which does the calculation up to 6 digits. I would like to limit matlab calculation as well to 6 digits for all of the calculations executed by an m-file script. Is there any command or option setting for this purpose?

Thanks in advance

Iman

Answer 1

[James Tursa]

You can get about that precision (a little less) by using the half data type:

https://www.mathworks.com/help/fixedpoint/ref/half.html?s_tid=srchtitle

Answer 2

[Bruno Luong]

Perhaps all you need is cast all your input data/parameters/constants related your program to SINGLE.

The relative precision is

>> eps('single')

ans =

  single

  1.1921e-07

so somwhere between 6 and 7 digits.

I really doubt any SW works on exactly 6 digits. 99.999% SW works on IEEE floating point precision, single or double.

[John D'Errico]

Good point that single is at least close to the desired precision. The virtue being that single will at least be highly efficient compared to any other choice like HPF or VPA.

HPF is truly a decimal tool, so 6 digits is 6 digits.

Answer 3

[John D'Errico]

There are only two standard choices in MATLAB for floating point computations, and two ways to get round the problem.

That is, you can use single or double for floaing point numbers. But even there, you need to be careful, as MATLAB defaults to a double. You cannot set some flag in the top line of your code that says "From here on, ALL new numbers created, and all precision in your computations will be done in K number of digits". You cannot even set single precision as the default.

Having said that, there are two ways out. You can use the symbolic toolbox.

digits 6

X = vpa(pi)
X =
3.14159

X^2 - 3
ans =
6.8696

Or, you can use my HPF toolbox, available for free download from the file exchange.

https://www.mathworks.com/matlabcentral/fileexchange/36534-hpf-a-big-decimal-class

% set the default precision to 6 decimal digits,
% with no hidden digits to minimize floating point
% errors. As well, set the storage mode to be blocks
% of 3 decimal digits.
DefaultNumberOfDigits 6 0
DefaultDecimalBase 3

X = hpf('pi')
X =
3.14159

X^2 - 3
ans =
6.86960

Answer 4

[Rik]

You could round the end result. This will of course not replicate any rounding errors.

In general calculations will run at the precision of the data type (IEEE single or double, or an integer type). So if you insist, you will have to create a custom data type that will replicate all rounding errors arising from working with 6 digit precision. I suspect the easiest would be to use int64, multiply the input by 1e6, and divide the output by 1e6 (but casting to double/single before that division).

[John D'Errico]

This might appear to be not a bad idea, except that some computations are not defined for int64. Just picking the first function that comes to mind, I see:

normcdf(1)
ans =
         0.841344746068543
         
normcdf(int64(1))
Error using NaN
Input following 'like' is not a single or double array.

And while this is a strange error, it does fail.

There are other important reasons why this idea will fail. Consider any function that is nonlinear.

f = @(x) x.^2;

f(pi)
ans =
          9.86960440108936
          
double(f(int64(pi*1e6))/1e6)
ans =
     9869607

As you can see, there is an extra factor of 1e6 in there.

This suggests your idea works only on fully linear computations. A matrix multiply, perhaps. Or, does it?

X = 1./magic(3)
X =
                     0.125                         1         0.166666666666667
         0.333333333333333                       0.2         0.142857142857143
                      0.25         0.111111111111111                       0.5
>> Xi = int64(1e6*X)
Xi =
  3×3 int64 matrix
    125000   1000000    166667
    333333    200000    142857
    250000    111111    500000
>> Xi*Xi
Error using  * 
MTIMES (*) is not fully supported for integer classes. At least one argument must be scalar. 

Drat. Even mtimes fails to run for int64. But if it did, then you would still have seen an error in the conversion back into doubles, because again, we would see an extra factor of 1e6 in the result.

And finally, int64 overflows way sooner than a double overflows.

intmax('int64')
ans =
  int64
   9223372036854775807
   
realmax('double')
ans =
     1.79769313486232e+308

Since we start out by multiplying by 1e6, that gets us 1e6 closer to the overflow point for int64.

[Rik]

Hmm. So my suggested solution is even worse than I initially thought. Thanks for the thorough explanation. I guess the only side-effect-free solution (short of re-implementing most of Matlab with a custom data type) is to round the end result and live with the fact that this might contain rounding errors in the wrong direction.

I have never worked with vpa, but since you can specify a precision, it may be a better choice than normal rounding of the end result.