# Fun with 40 year old code

**URL:** <https://fortran-lang.discourse.group/t/fun-with-40-year-old-code/3726>\
**Category:** Uncategorized\
**Created:** [June 11, 2022, 5:59pm UTC](https://fortran-lang.discourse.group/t/fun-with-40-year-old-code/3726 "2022-06-11T17:59:37Z")\
**Posts on this page:** 1\
**Showing post:** 3

<div class="post-metadata">

**Author:** ![mecej4](https://yyz2.discourse-cdn.com/free1/user_avatar/fortran-lang.discourse.group/mecej4/32/855_2.png) [@mecej4](https://fortran-lang.discourse.group/u/mecej4)\
**Post date:** [June 12, 2022, 12:00am UTC](https://fortran-lang.discourse.group/t/fun-with-40-year-old-code/3726/3 "2022-06-12T00:00:34Z")

</div>

> [@jacobwilliams](#):
>
> So why didn’t they just do that? If this really was the intent, then in modern Fortran I believe we can just replace it with:
> 
> ```auto
> nbits = digits(1.0_wp)
> 
> ```

That would not work on an [old IBM mainframe, which uses base 16, not base 2](https://en.wikipedia.org/wiki/IBM_hexadecimal_floating-point). You can check in the file i1mach.f, and you will find a small number of mostly obsolete machines for which I1mach(14) = 14, rather than the value 53 for IEEE double precision. For the menagerie of machines of the 70s and 80s, they had to convert between bits, nybbles, bytes, digits, etc.

Next, about the rounding of log10(2d0). We may speculate that some fudging was needed to get the answer to come out equal to a known result for the common case where I1mach(14) = 53. My speculation: the programmer tried at first

```auto
ANIB = D1MACH(5)*K/0.30103000D0

```

but that gave ANIB = 52.99999923659105…, and that made NIB = 52 (truncation to next lower integer), which was not “right”. Using 0.30102000 gives ANIB = 53.00175991, which leads to NIB = 53, which is what was wanted. I don’t understand why the ANINT intrinsic was not used.

Those were the days!

Additional comments:  
In this context, we may take note of a hidden thorn that awaits us in the Netlib file d1mach.f. Once upon a time, in order to avoid the messy equivalence statements in that file, I replaced the executable statements by

```auto
      SELECT CASE (I)
      CASE (1)
         D1MACH = TINY(1d0)
      CASE (2)
         D1MACH = HUGE(1d0)
      CASE (3)
         D1MACH = EPSILON(1d0)/2
      CASE (4)
         D1MACH = EPSILON(1d0)
      CASE (5)
         D1MACH = LOG10(2d0)
      END SELECT

```

and that gave mysterious errors in my programs.

The problem is that LOG10 is declared as an integer array, and many compilers allowed 2d0 as an index into that array (allowing a double precision index as an extension to the language), so LOG10(2) is equivalenced to the upper 32 bits of DMACH(5). One fix is to replace the variable LOG10 by ILOG10, so that it does not clash with the intrinsic function LOG10. Alternatively, one could remove all the now-superfluous local variables from the routine.

---

_[View the full topic](https://fortran-lang.discourse.group/t/fun-with-40-year-old-code/3726)._