HYF
HYF

Reputation: 173

How to change the rounding mode for floating point operations in MATLAB?

I want to change the rounding mode for floating point operations in MATLAB. According to IEEE 754-2008, there are 5 strategies for rounding:

Does MATLAB supports these 5 strategies? How to change the rounding mode for floating point operations in MATLAB?

Upvotes: 8

Views: 1775

Answers (2)

user2609605
user2609605

Reputation: 626

@HYF: I found that feature('setround', 0.5) leads to round to even. I checked the following: a=1+2^-52 This means that mantissa looks like so: 1.0...01 where the last 1 is still in the mantissa. The leading 1 is not stored in IIIE754 format. ( I checked that 1+2^-52 == 1 but not 1+2^-53 == 1)

Then I computed b = a + 2^-53. Without rounding it is 1.0...01|1 where the last digit is to be rounded away. I found the following true: b==1+2^-51 We have b == 1.0...010.

We have several rounding modes submodes of round to nearest: This can be round to inf, round away from 0 or round to even.

Next we check -b==-1-2^-51 to be true which excludes round to inf still allowing round away from 0 or round to even.

Then I checked 1.5==1.5+2^-53.
Of course 1.5 = 1.10...0 binary and and 1.5+2^-53 = 1.10...0|1 without rounding, the last digit to be rounded away. Rounding away from 0 would be 1.10...01 and rounding to even is 1.10...0. So the latter is the case.

Upvotes: 0

Andrew Janke
Andrew Janke

Reputation: 23848

Answer

Kind of. There is an undocumented feature('setround') function call that you can use to get or set the rounding mode used by Matlab.

So, it can be done, but you shouldn’t do it. :)

WARNING: This is an undocumented, unsupported feature! Use at your own peril!

This feature('setround') supports 4 of the 5 IEEE-754 rounding modes: there’s only one “nearest” mode, and I don't know if it’s “ties to even” or “ties away from zero”.

Supported modes:

  • feature('setround') – Get current rounding mode
  • feature('setround', 0.5) – Round toward nearest (don’t know if it’s ties to even or away from zero)
  • feature('setround', Inf) – Round up (towards +Inf)
  • feature('setround', 0) – Round toward zero
  • feature('setround', -Inf) – Round down (towards -Inf)

Note on testing: The IEEE-754 rounding mode does not affect round() and its relatives. Rather, it governs how arithmetic operations behave around the limits of floating-point precision.

Demonstration

%ROUNDINGEXAMPLE Demonstrates IEEE-754 Rounding Mode control
%
% This uses a completely undocumented and unsupported feature!
% Not for production use!

%% Setup
clear; clc

n = 2000;
X = ones(n)*1E-30; % matrix with n^2 elements
defaultRoundingMode = feature('setround'); % store default rounding mode

%%
feature('setround',0.5);
r1 = prettyPrint('Nearest', sum(X(:)));
%{
  sign   exponent                       mantissa
     0 01110110001 0011010101111100001010011001101001110101010000011110
     | \_________/ \__________________________________________________/
     |      |             ______________________|___________________________
     |      |            /                                                  \
(-1)^0 2^( 945 - 1023) 1.0011010101111100001010011001101001110101010000011110 = 4e-24
%}

%%
feature('setround',-Inf);
r2 = prettyPrint('To -Infinity', sum(X(:)));
%{
  sign   exponent                       mantissa
     0 01110110001 0011010101111100001010011001101001011100000111000110
     | \_________/ \__________________________________________________/
     |      |             ______________________|___________________________
     |      |            /                                                  \
(-1)^0 2^( 945 - 1023) 1.0011010101111100001010011001101001011100000111000110 = 4e-24
%}

%%
feature('setround',Inf);
r3 = prettyPrint('To Infinity', sum(X(:)));
%{
  sign   exponent                       mantissa
     0 01110110001 0011010101111100001010011001101010100011101100100001
     | \_________/ \__________________________________________________/
     |      |             ______________________|___________________________
     |      |            /                                                  \
(-1)^0 2^( 945 - 1023) 1.0011010101111100001010011001101010100011101100100001 = 4e-24
%}

%%
feature('setround',0);
r4 = prettyPrint('To zero', sum(X(:)));
%{
  sign   exponent                       mantissa
     0 01110110001 0011010101111100001010011001101001011100000111000110
     | \_________/ \__________________________________________________/
     |      |             ______________________|___________________________
     |      |            /                                                  \
(-1)^0 2^( 945 - 1023) 1.0011010101111100001010011001101001011100000111000110 = 4e-24
%}

%%
feature('setround',defaultRoundingMode);
r5 = prettyPrint('No accumulated roundoff error', 4e-24);
%{
  sign   exponent                       mantissa
     0 01110110001 0011010101111100001010011001101010001000111010100111
     | \_________/ \__________________________________________________/
     |      |             ______________________|___________________________
     |      |            /                                                  \
(-1)^0 2^( 945 - 1023) 1.0011010101111100001010011001101010001000111010100111 = 4e-24
%}

%% Helper function
function r = prettyPrint(s, r)
    fprintf('%s:\n%65.60f\n\n', s, r); 
end

I get:

Nearest:
   0.000000000000000000000003999999999966490758963870373537264729

To -Infinity:
   0.000000000000000000000003999999999789077070014108839608005726

To Infinity:
   0.000000000000000000000004000000000118618095059505975310731249

To zero:
   0.000000000000000000000003999999999789077070014108839608005726

No accumulated roundoff error:
   0.000000000000000000000003999999999999999694801998206811298525

Acknowledgments

Thanks to Ryan Klots at MathWorks Technical Support for setting me straight on this and providing the nice demo code!

Upvotes: 7

Related Questions