Search code examples
fortran77

Minimal and maximal magnitude in Fortran

I'm trying to rewrite minpack Fortran77 library to Java (for my own needs), so I met this in minpack.f source code:

integer mcheps(4)
integer minmag(4)
integer maxmag(4)
double precision dmach(3)
equivalence (dmach(1),mcheps(1))
equivalence (dmach(2),minmag(1))
equivalence (dmach(3),maxmag(1))
...
data dmach(1) /2.22044604926d-16/
data dmach(2) /2.22507385852d-308/
data dmach(3) /1.79769313485d+308/
dpmpar = dmach(i)
return

What are minmag and maxmag functions, and why dmach(2) and dmach(3) have these values? There is an explanation in comments: 

c         dpmpar(1) = b**(1 - t), the machine precision,
c         dpmpar(2) = b**(emin - 1), the smallest magnitude,
c         dpmpar(3) = b**emax*(1 - b**(-t)), the largest magnitude.

What is smallest and largest magnitude? There must be a way to count these values on runtime; machine constants in source code is a bad style.

EDIT: I suppose that static fields Double.MIN_VALUE and Double.MAX_VALUE are those values I looked for.

Solution

  • minmag and maxmag (and mcheps too) are not functions, they are declared to be rank 1 integer arrays with 4 elements each. Likewise dmach is a rank 1 3 element array of double precision values. It is very likely, but not certain, that each integer value occupies 4 bytes and each d-p value 8 bytes. Bear this in mind as the answer progresses.

    So an expression such as mcheps(1) is not a function call but a reference to the 1st element of an array.

    equivalence is an old FORTRAN feature, now deprecated both by language standards and by software engineering practices. A statement such as 

    equivalence (dmach(1),mcheps(1))

    states that the first element of dmach is located, in memory, at the same address as the first element of mcheps. By implication, this also means that the 24 bytes of dmach occupy the same addresses as the 16 bytes of mcheps, and another 8 bytes too. I'll leave you to draw a picture of what is going on. Note that it is conceivable that the code originally (and perhaps still) uses 8 byte integers so that the elements of the equivalenced arrays match 1:1.

    Note that equivalence gives, essentially, more than one name, and more than one interpretation, to the same memory locations. mcheps(1) is the name of an integer stored in 4 bytes of memory which form part of the storage for dmach(1). Equivalencing used to be used to implement all sorts of 'clever' tricks back in the days when every byte was precious.

    Then the data statements assign values to the elements of dmach. To me those values look to be just what the comment tells us they are.

    EDIT: The comment indicates that those magnitudes are the smallest and largest representable double precision numbers on the platform for which the code was last compiled. I think that in Java they are probably called doubles. I don't know Java so don't know what facilities it has for returning the value of the largest and smallest doubles, if you don't know this either hit the 'net or ask another SO question -- to which you'll probably get responses along the lines of search the net.

    Most of this you should be able to ignore entirely. As you write, a better approach would be to find out those values at run-time by enquiry using intrinsic functions. Fortran 90 (and later) have such functions, I imagine Java has too but that's your domain not mine.