Switching between single and double precision in Fortran
I want to write fortran 90 or fortran 2003 code and be able to select between either single/double precision at compile time.
In C++ it is extremely easy to do this using template meta-programming. Code:
template<typename T> However I have very little idea on how to implement such a generic function in Fortran. I can select between single and double precision for basic arithmetic using selected_real_kind method and using real(xp) to define the required precision. Code:
! Define single & double precision Code:
! Constant fractions So, I tried to define generic interfaces to mathematical functions like : Code:
!==-------------------------------------------------------------------------==! Code:
pure function sqrt_dp(x) I can use sqrt_gn(9.0_xp) to get the correct square root value for my preferred precision and the code works just fine. But I have some questions.... QUESTIONS : + Is there a better way to do this? <-- MOST IMPORTANT REQUEST + Is there any performance loss due to using this indirection through the interface and then to the pure functions? |
While Fortran doesn't have, for some good reasons, the same form of polymorphism used by C++, you are clearly misinterpreting how selected_xxx_kind types are meant to be used. Once you have the xp kind you can have on input to any function/subroutine that knows it something like:
REAL(xp), INTENT(IN) :: input Because xp is just a placeholder for some already premapped type. Of course, this only works for basic types and not, say, for derived ones (still, they can be parameterized) Also, you should not, in any case, use something like dsqrt. sqrt would have worked just fine independently from the declaration issue To give you an example: MODULE mysqrt_m USE myprecision_m !where xp is defined FUNCTION mysqrt(input) RESULT(res) REAL(xp), INTENT(IN) :: input REAL(xp) :: res res = sqrt(input) ENDFUNCTION mysqrt ENDMODULE mysqrt_m |
Why not just use the compiler flag to specify precision? Most FORTRAN compilers have a flag (For intel you can use -r8). Of course, if you are building a module that you want to be portable and usable for either single or double precision then the correct answer is to create single and double precision routines and provide the explicit interface. With regard to performance issues, I have never seen any significant problems, even dealing with complex-valued functions and routines. The interface allows the compiler to know ahead of time which routine to use based on the data type being sent. So there is no type conversion hit to performance (or so I've been told).
|
agd gives good suggestion. I have always done it like this. Makes the code more readable also. Declare floats as "real :: x", set values like
x = 1.0 Use sqrt (not dsqrt), etc. Then while compiling, say with gfortran, use -fdefault-real-8 if you want double precision version. |
Unless there is a very specific business case, honestly, I would kind avoid approaches based on compiler flags. If switching precision is what you want, selected_xxx_kind are there also for this, and you should always add _xp to the reals so defined.
In my opinion code should always speak for itself. Of course, if you still need to write subroutines for multiple types, you still need an interface etc., but that is a different issue |
Quote:
This is explained in this stackoverflow answer : https://stackoverflow.com/a/24070364 |
Quote:
In the first case it is possible to use both the compiler flags and the selected kinds approaches. Still, as I wrote, I would prefer the latter. In the second case, of course, there is no chance to use the two previous approaches and you need a different route. The tedious, formal one is the interface (which is a non solution, because you actually write a lot of identical code + boiler plate stuff for the interface). The meta-language approach is another one. Personally, I am an integralist of formalisms, and the fact that you are using some meta language to solve an issue with your selected language, well... I know it is used a lot, but I feel pain every time I see such solutions. Honestly, if a library is what you want to write, and you don't want to use Fortran interfaces and stuff, there is clear indication that you have selected the wrong language in the first place, at least for now (may change in the future, who knows). I don't write code for others to be used in a library, so I might be skewed, but whenever I find myself in the need of the same subroutine for different types (and honestly, this only happened for linked lists and sort/search stuff) I first take the chance to analyze the thing for every different data type, in order to understand if I can rearrange stuff differently, reuse memory in some cases, etc. Then I just put the different version in a module and keep the names different and clear. Later in time, when everything is stabilized, I decide if the boilerplate code is actually worth the effort just for the sake of using the same name for calling the function/subroutine on different types... the answer is typically not, but there is certainly bias due to my use case. |
Quote:
Use modern Fortran and just use sqrt(), which will take either a single or double precision real and return a value of the same kind. Second point: Don't use selected_real_kind(). Use the REAL32/REAL64-Types from the ISO_FORTRAN_ENV module Code:
use, intrinsic :: ISO_FORTRAN_ENV, only : INT8, INT16, INT32, INT64, & ! Standardized Fortran 2008 kind definitions |
Dear Sayan:
If you just want to switch between single and double precision float points during compile time, you can declare them as REAL(KIND=4)'s (single precision) in the code, and then use the gfortran compiler switch -freal-4-real-8 to force them into REAL(KIND=8)'s (double precision). You can even convert them into "long double" with the -freal-4-real-16 switch if you want. There used to be an -r8 switch which will force all REAL's into DOUBLE PRECISION's in g77, but it seems that it's not in gfortran anymore. You might want to check for portability first, but I suspect there should be compiler options in other Fortran products that are equivalent. Gerry. P.S. - here is a list of gfortran compiler options |
Quote:
|
Quote:
Um, yes, they are? Page 400: https://j3-fortran.org/doc/year/10/10-007r1.pdf |
Quote:
Code:
If the processor supports no kind of a particular size, that constant shall be equal to −2 if the processor supports a kind with larger size and −1 otherwise. |
|
All times are GMT -4. The time now is 07:04. |