Bug? MATLAB MEX changes the kind of the default logical

Question

When interfacing a piece of Fortran 2003 (or above) code with MATLAB by MEX, I am surprised to find that MEX changes the kind of the default logical. This is fatal, because a piece of perfectly compilable Fortran code may fail to be mexified due to a mismatch of types, which did happen in my project.

Here is a minimal working example.

Name the following code as "test_kind.F", compile it by mex test_kind.F in MATLAB, and then run test_kind in MATLAB. This will produce a plain text file named fort.99, which contains two numbers "4" and then "8" as the result of the WRITE instructions.

! test_kind.F
! Tested by MATLAB 9.8.0.1323502 (R2020a) with GNU Fortran (Ubuntu 9.3.0-17ubuntu1~20.04) 9.3.0

#include "fintrf.h"

      subroutine mexFunction(nlhs, plhs, nrhs, prhs)

      use ieee_arithmetic, only : ieee_is_nan
      implicit none
      mwPointer plhs(*), prhs(*)
      integer nlhs, nrhs

      write(99, *) kind(ieee_is_nan(1.0))  ! This prints a number in fort.99
      write(99, *) kind(.false.)  ! A benchmark, which should print the same number in fort.99
      close(99)

      end subroutine mexFunction

I thought the two printed numbers should always be equal to each other, although the specific value depends on the compiler and needs not to be 4 or 8. (As stressed by Doctor Fortran @SteveLionel, the Fortran standard imposes no relation between these kind numbers and the number of bytes used to represent the data. See Steve's blog Doctor Fortran in “It Takes All KINDs” for a nice elaboration on this topic.)

[Update: The above speculation about kind(ieee_is_nan(1.0)=kind(.false.) turns out wrong even though it is what the Fortran 2018 standard suggests. With certain compiler options, kind(ieee_is_nan(1.0) and kind(.false.) can actually differ from each other and hence violate the specifications in the Fortran standard. See the answer by @francescalus and my summary at the end of this question.]

[Update: in the original question, I flipped "4" and "8", so that I thought kind(ieee_is_nan(1.0)) was changed; indeed, it is kind(.false.) that is altered from 4 to 8, while kind(ieee_is_nan(1.0)) remains 4 always. So everything is fully explained by the nice answer from @francescalus, thank you!]

For comparison, here is the same code without interfaced with MATLAB, which prints "4" and "4" on the screen when compiled with gfortran. With the nagfor compiler, the numbers remain equal although become 3.

! test_kind.f90
! Tested by GNU Fortran (Ubuntu 9.3.0-17ubuntu1~20.04) 9.3.0

program test_kind

use ieee_arithmetic, only: ieee_is_nan
implicit none

write(*, *) kind(ieee_is_nan(1.0))  ! This prints a number on STDOUT (screen)
write(*, *) kind(.false.)   ! A benchmark, which should print the same number on STDOUT (screen)

end program test_kind

For your reference, below is the section on ieee_is_nan in the Fortran 2018 standard. It specifies that the ieee_is_nan returns a "default logical", which I guess should be identical to the type of the intrinsic .true. or .false. constant --- or did I get it wrongly?

17.11.13 IEEE_IS_NAN (X)
1 Description. Whether a value is an IEEE NaN.
2 Class. Elemental function.
3 Argument. X shall be of type real.
4 Restriction. IEEE_IS_NAN (X) shall not be invoked if IEEE_SUPPORT_NAN (X) has the value false.
5 Result Characteristics. Default logical.
6 Result Value. The result has the value true if the value of X is an IEEE NaN; otherwise, it has the value false.

It seems extraordinary to me that MEX can change the type of the default logical without taking care of ieee_is_nan. Maybe there exists an option of MEX that can rectify this behavior, but why should it be the default in the first place?

I tried the same code on more machines with other versions of MATLAB and Fortran compiler. The results are the same.

1. MATLAB 9.7.0.1319299 (R2019b) Update 5 with GNU Fortran (GCC) 8.3.1 20191121 (Red Hat 8.3.1-5):

   kind(ieee_is_nan(1.0))= 4, kind(.false.) = 8

   The same compiler without interfacing with MATLAB:

   kind(ieee_is_nan(1.0)) = 4 = kind(.false.)

2. MATLAB 9.5.0.1049112 (R2018b) Update 3 with GNU Fortran (Ubuntu 9.3.0-17ubuntu1~20.04) 9.3.0:

   kind(ieee_is_nan(1.0)) = 4, kind(.false.) = 8

   The same compiler without interfacing with MATLAB:

   kind(ieee_is_nan(1.0)) = 4 = kind(.false.)

3. (On Windows 10) MATLAB 9.5.0.944444 (R2018b) with Intel(R) Visual Fortran Intel (R) 64 Compiler, Version 19.1.1.216 Build 20200306

   kind(ieee_is_nan(1.0)) = 4, kind(.false.) = 8

   The same compiler without interfacing with MATLAB:

   kind(ieee_is_nan(1.0)) = 4 = kind(.false.)

---

A summary after accepting @francescalus' answer

1. It turns out that the mismatch between kind(.false.) and kind(ieee_is_nan(1.0)) comes from the gfortran option -fdefault-integer-8, which is adopted by MEX as the default. This option enforces gfortran to use 64-bit integer and 64-bit logical as the default kind, but does not change the returned kind of ieee_is_nan, even though the Fortran standard specifies that ieee_is_nan should return the default logical kind. This may be because ieee_is_nan is not really an intrinsic procedure but only a procedure from the intrinsic module ieee_arithmetic.

2. Note that ifort (version ifort (IFORT) 2021.2.0 20210228) and nagfor (NAG Fortran Compiler Release 7.0 (Yurakucho) Build 7036) behave also in the way described above, the corresponding option being -i8 for both of them. So the compiler vendors agree that it is fine to break the consistency of some intrinsic modules when certain options are enforced. This is surprising to me. Fortunately, flang (under clang 7.1.0) follows the Fortran standard even with -fdefault-integer-8 imposed, keeping kind(is_ieee_nan(1.0)) == kind(.false.) --- so it is not a mission impossible.

3. The NAG compiler nagfor warns about ieee_is_nan when -i8 is adopted, pointing out that they are incompatible; gfortran and ifort, however, keep absolutely silent even if you invoke them with -Wall -Wexta and -warn all respectively. This is even more surprising to me.

4. Given these facts, I decide not to use ieee_is_nan but to implement my own is_nan and stay vigilant about (or away from) all the procedures provided by intrinsic modules. Otherwise, my package will fail to compile due to a mismatch of types if the users choose to enforce the 64-bit integer as the default (without thinking about the logical kind; why should they?). More seriously, MATLAB has already made such a choice for all its users without telling them.

5. I am happy to see that my question leads to interesting discussions. Since some compilers seem to be in need of improvements regarding the problems discovered here (you may have different opinions), I made a post on Fortran Discourse about this topic. I hope that it will draw more attention from the community and someone will patch at least the open-source compilers like gfortran.

Thank you very much for any comments or criticism.

Cammiecammy

They should be the same - as you say MATLAB is probably using some horrible compiler flags that changes the types in your code, and so breaks the Fortran. – Pall

As possibly pointless pedantry (just in case it's relevant): ieee_is_nan is not a (standard) intrinsic procedure. – Morelos

Thank you for pointing out this. It is a procedure that is specified in Fortran standard 2018. Does this make it qualified as an intrinsic procedure? Maybe I am wrong with the terminology. – Cammiecammy

It's a procedure in an intrinsic module, but such procedures are themselves not intrinsic. Turns out that that's actually relevant (if my answer is correct). – Morelos

@IanBush It turns out that they may not be the same under some compiler options. – Cammiecammy

Answer 1

By default MEX compiles with the gfortran option -fdefault-integer-8. The way gfortran handles this results in what you see.

Consider the non-MEX program

  use, intrinsic :: ieee_arithmetic, only : ieee_is_nan, ieee_support_nan
  implicit none

  if (.not.ieee_support_nan(1.)) error stop "Lack of support"
  print*, KIND(ieee_is_nan(1.)), KIND(.TRUE.)

end

compiled with/without -fdefault-integer-8.

This option makes default integer and logical variables 8 bytes wide (and have kind parameter 8).

This all makes sense. However, it seems that gfortran doesn't change the function result kind parameter for functions in the intrinsic module ieee_arithmetic with this option (this is reasonable in some way: it's just taking what the "pre-compiled" module files say, which is "return is logical(4)".

(If things seem confusing when you see that it does seem to use the correct kinds in another intrinsic module like iso_c_binding, note that this second intrinsic module isn't supported by a supplied module file. You will see the same behaviour for IEEE modules with the OpenMP intrinsic module which is also provided as a file.)

As a workaround, you can use LOGICAL(ieee_is_nan(1.)) to make things match again.

---

Compiler options like -fdefault-integer-8 are not toys. They should be used only in the rarest of rare circumstances. With gfortran, if you use this option by itself you are telling the compiler "Go ahead and do whatever you like: I don't care whether or not you treat my program as the Fortran standard specifies." In this case, this is a choice that MATLAB has made for you.

In particular, if you use a compiler option to change kind parameters of things (such as default kinds, or the kind numbering scheme as seen in the NAG compiler for example), you must compile all constituent parts of your program with that option to be sure of consistency. As seen in this case, that's not always easy, making these horribly dangerous options to have as defaults. Even intrinsic modules may be included in the list of things which need to be considered for consistency.

Morelos

This analysis is very informative. Thank you very much. 1. So the inconsistency comes from gfortran, due to the option that MEX uses. 2. However, why in my very first example, kind(ieee_is_nan(1.0))=8 while kind(.false.)=4? It should be reversed as in the other examples. – Cammiecammy

I was going to suggest that this probably isn't the desired behaviour of gfortran, so things may well change over versions, etc. But if you don't see such variation then that's less important to note. (I blame MATLAB for using inexcusably using this option.) – Morelos

" Even intrinsic modules may be included in the list of things which need to be considered for consistency." I totally agree. It is indeed astonishing to me. The NAG compiler nagfor warns about ieee_is_nan when -i8 is adopted, pointing out that they are incompatible; gfortran and ifort, however, keep absolutely silent even if you invoke them with -Wall -Wexta and -warn all respectively. This is surprising. – Cammiecammy

I'd seriously consider filing bug reports (QoI issue at the least) in these cases. I have sympathy for compiler vendors in this regard but a warning is perhaps the least that can be done. I suspect Matlab has gone down this road far too long for this reliance to be fixed, but again documentation can be changed. – Morelos

So you’re saying there’s a compiler option that breaks things, but you don’t blame the compiler maker for creating this option? :D – Navarre

@CrisLuengo, sounds silly when put that (tongue-in-cheek) way, but a shockingly high proportion of my professional life involves support for decisions made thirty years ago which make new developments tricky. An entirely different way of implementing (some) intrinsic modules just to make it easier for users to shoot themselves in the foot? If you're using the IEEE modules you should be coding everything in a "modern way", perhaps. – Morelos

@CrisLuengo, as the person that created the names for these options and made the options work (for some definition "work"), I cordially invite you to join the one remaining gfortran volunteer to "fix" these options. My proposed fix was rejected years ago during code review. These options should be used only when porting code from one precision to another; and, definitely, never used in production code. – Palacio

Answer 2

(Edit - I misread the statement. See comment below.)

Certainly, these two numbers depend on the compiler --- they are 3 for the nagfor compiler, but I thought they should always be equal to each other.

NO! There is no implicit correspondence between the kind numbers of one type and the kind numbers of another type! Do not be fooled by the common use of byte sizes for kind numbers - this is just a convention some compilers adopted to look more familiar to programmers used to the integer*4 extension.

It would be perfectly valid for a compiler to have kinds 14, 27 and 830 for LOGICAL and 2,9 and 13 for REAL. The only things you can count on are that:

• Default integer, default real and default logical each occupy one "numeric storage unit"
• Double precision and default complex occupy two numeric storage units

(See Fortran 2018 19.5.3.2 Storage Sequence)

Of course, as already pointed out, if you use compiler options to change the default kind for one type but not another, you break this rule.

For further reading, see my blog post, Doctor Fortran in “It Takes All KINDs”

Cowley

Hi Dr. Fortran @Steve Lionel, it is a great honour to see you answering my question. However, I think I did not make myself clear enough in that sentence. I know very well (thanks to your writings) that the value of the kind has absolutely no relation with the number of bytes representing the data. That's why I said these values are compiler-dependent. In that sentence, I meant that kind(.false.) and kind(ieee_is_nan(1.0)) should be equal, although the specific value depends on the compiler (can be 830 as you said). Nevertheless, this claim has been disproved. Thank you very much. – Cammiecammy

I agree with your points, but I don't see too much relevance I'm afraid. In particular: "There is no implicit correspondence between the kind numbers of one type and the kind numbers of another type!" is true, but there's only the logical type being queried here? – Morelos

OK, I see I didn't understand the point being made well enough. Apologies. I would expect both results to be default logical. – Cowley

Thank you @SteveLionel. No, it was me who did not make the point clear. My English needs improvement. By the way, as a guru of Fortran and ifort, what is your opinion on the fact that kind(ieee_is_nan(1.0)) /= kind(.false.) when ifort is invoked with -i8? In adition, the compiler is silent about it even with -warn all. Than you very much! – Cammiecammy

@Nuno, I am not seeing that behavior with a current version of ifort. Feel free to ask about it at community.intel.com/t5/Intel-Fortran-Compiler/bd-p/… and be sure to provide a small but complete test case. – Cowley