Resistance to modernization

1

What do you think are the main reasons to avoid modernization/refactoring of legacy Fortran codes if we exclude obvious reasons such as cost/man-hours needed?

By modernization I mean rewriting as free-form, explicit typing, removal of obsolescent control-flow constructs and non-standard features, introduction of new intrinsics (e.g. replacing MACH with epsilon) and other changes which do not require a change in the routine interfaces. In other words, non-breaking changes.

I was thinking about this after going through some of the open LAPACK issues:

The main reasons I could identify (not just limited to LAPACK) are:

  • no (perceived) need; the code works as is - “if it ain’t broken, don’t fix it”
  • massive reformatting effort for potentially little gain beyond perhaps improved readability
  • modernization (as in updating to free-form) doesn’t actually bring new features
  • losing the option to use f2c as a portable means of distributing libraries given that C compilers are ubiquitous
  • using a restricted language subset can make it easier to port code to other frameworks (e.g. Java or .NET)
  • retain support for legacy compilers and platforms
  • high risk of introducing new bugs
  • maintaining consistency with existing code-base and project guidelines
5 Likes
2

All of the above. If you are already maintaining the code and want to introduce modern features, fine, but I think it is counterproductive to rewrite the code just because it’s old-style. This is one of the strengths of Fortran.

12 Likes
3

I generally agree with the “if it ain’t broke, don’t fix it” for legacy code. However, in the case of LAPACK/BLAS, there is a concrete advantage to changing the dummy array arguments to assumed shape. As it is now, with explicit shape or assumed size dummy arguments, then many actual arguments must be associated with the dummy arguments with copy-in/copy-out rather than the simpler/cheaper address association. That is expensive and arguably unnecessary. So my argument basically is that the current LAPACK is a little broken as is, so maybe it does need to be fixed.

However, this change would require explicit interfaces (either through modules or with interface blocks), something that is not required now, and it would make it a little more difficult to interface the library to other languages, such as C and python and so on. But this might be a blessing in disguise. If LAPACK is deemed to be important enough, then this might be the mechanism to finally get standardized definitions of the fortran array (and pointer) metadata. If those could be standardized somehow across all fortran compilers, then interfacing generic modern fortran (not just a specific compiler) to other languages would be much easier.

6 Likes
4

IMHO code rewrites should be the last resort, and it should happen only if your legacy code is broken beyond repair and built upon a really bad architecture that holds you back when implementing new features, otherwise writing higher-level interfaces is the best option.
There is always a risk-reward involved, so rewriting for the sake of it may simply not be worth it, since you are investing your time with no feasible return other than a programmer in the future (even yourself) will now easily understand such code, which is only good in cases where there is actually a good amount of space for improvement in the first place.

2 Likes
5

I just wanted to point out that the idea of writing higher-level interface (e.g. LAPACK95) does not solve the problem associated with copy-in/copy-out. While the user code can call the LAPACK95 routine correctly, the copy-in/copy-out then occurs between the LAPACK95 and the LAPACK routines. So the interface just pushes the problem down one level. To actually fix this problem requires that all of LAPACK be rewritten to use assumed shape arrays, every routine at every level. This would be a big change.

3 Likes
6

As a sidenote, the motivation for this thread came from refactoring the LAPACK function lsame. The Reference LAPACK version is cluttered by additional else-if branches for processors that use other encodings than ASCII such as IBM mainframes using EBCDIC or the Prime OS @MarDie posted about.

lsame is really just a “trivial” support routine used to compare two characters. It is used in the routines which can perform operations on non-transpose, transpose or conjugate forms of matrices, i.e. 'N', 'T', to check which character was provided at the call site and perform the corresponding operation.

The original and F95 code are given in the drop-down arrows below. Observe that gfortran generates the same assembly code when using a reasonable optimization flag. Based on the character encoding it infers from the statement ZCODE = ichar('Z') it can simply eliminate the dead branches.

Original 
      LOGICAL FUNCTION LSAME( CA, CB )
*
*  -- LAPACK auxiliary routine --
*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
*
*     .. Scalar Arguments ..
      CHARACTER          CA, CB
*     ..
*
* =====================================================================
*
*     .. Intrinsic Functions ..
      INTRINSIC          ICHAR
*     ..
*     .. Local Scalars ..
      INTEGER            INTA, INTB, ZCODE
*     ..
*     .. Executable Statements ..
*
*     Test if the characters are equal
*
      LSAME = CA.EQ.CB
      IF( LSAME )
     $   RETURN
*
*     Now test for equivalence if both characters are alphabetic.
*
      ZCODE = ICHAR( 'Z' )
*
*     Use 'Z' rather than 'A' so that ASCII can be detected on Prime
*     machines, on which ICHAR returns a value with bit 8 set.
*     ICHAR('A') on Prime machines returns 193 which is the same as
*     ICHAR('A') on an EBCDIC machine.
*
      INTA = ICHAR( CA )
      INTB = ICHAR( CB )
*
      IF( ZCODE.EQ.90 .OR. ZCODE.EQ.122 ) THEN
*
*        ASCII is assumed - ZCODE is the ASCII code of either lower or
*        upper case 'Z'.
*
         IF( INTA.GE.97 .AND. INTA.LE.122 ) INTA = INTA - 32
         IF( INTB.GE.97 .AND. INTB.LE.122 ) INTB = INTB - 32
*
      ELSE IF( ZCODE.EQ.233 .OR. ZCODE.EQ.169 ) THEN
*
*        EBCDIC is assumed - ZCODE is the EBCDIC code of either lower or
*        upper case 'Z'.
*
         IF( INTA.GE.129 .AND. INTA.LE.137 .OR.
     $       INTA.GE.145 .AND. INTA.LE.153 .OR.
     $       INTA.GE.162 .AND. INTA.LE.169 ) INTA = INTA + 64
         IF( INTB.GE.129 .AND. INTB.LE.137 .OR.
     $       INTB.GE.145 .AND. INTB.LE.153 .OR.
     $       INTB.GE.162 .AND. INTB.LE.169 ) INTB = INTB + 64
*
      ELSE IF( ZCODE.EQ.218 .OR. ZCODE.EQ.250 ) THEN
*
*        ASCII is assumed, on Prime machines - ZCODE is the ASCII code
*        plus 128 of either lower or upper case 'Z'.
*
         IF( INTA.GE.225 .AND. INTA.LE.250 ) INTA = INTA - 32
         IF( INTB.GE.225 .AND. INTB.LE.250 ) INTB = INTB - 32
      END IF
      LSAME = INTA.EQ.INTB
*
*     RETURN
*
*     End of LSAME
*
      END

Output from compiler explorer, with gfortran-12 -ffixed-form -O2

lsame_:
        movzx   eax, BYTE PTR [rdi]
        movzx   edx, BYTE PTR [rsi]
        mov     ecx, 1
        cmp     al, dl
        je      .L1
        lea     esi, [rax-97]
        lea     ecx, [rax-32]
        cmp     esi, 26
        lea     esi, [rdx-97]
        cmovb   eax, ecx
        cmp     esi, 26
        lea     ecx, [rdx-32]
        cmovb   edx, ecx
        xor     ecx, ecx
        cmp     eax, edx
        sete    cl
.L1:
        mov     eax, ecx
        ret
Refactored 
      pure logical function lsame( ca, cb )
         character(len=1), intent(in) :: ca, cb

         integer, parameter :: diff_case = iachar('a')-iachar('A')
         integer :: ia, ib

         lsame = ca == cb
         if (lsame) return 

         ! case-insensitive test
         ia = iachar(ca)
         ib = iachar(cb)

         ! Capital alphabetic characters are between 65 and 90
         ! in the ASCII character encoding
         if ( 64 < ia .and. ia < 91) ia = ia + diff_case
         if ( 64 < ib .and. ib < 91) ib = ib + diff_case

         lsame = ia == ib

      end

Output from compiler explorer, with gfortran-12 -ffixed-form -f95 -O2

lsame_:
        movzx   eax, BYTE PTR [rdi]
        movzx   edx, BYTE PTR [rsi]
        mov     ecx, 1
        cmp     al, dl
        je      .L1
        lea     esi, [rax-65]
        lea     ecx, [rax+32]
        cmp     esi, 26
        lea     esi, [rdx-65]
        cmovb   eax, ecx
        cmp     esi, 26
        lea     ecx, [rdx+32]
        cmovb   edx, ecx
        xor     ecx, ecx
        cmp     eax, edx
        sete    cl
.L1:
        mov     eax, ecx
        ret

From a functional point of view, you could say the refactoring brings little to no change. However the impact on readability and intent is non-negligible IMO.

Amusingly, the C interfaces to this function (which from a modern Fortran perspective shouldn’t even be part of the public interface) are different:

  • Intel MKL: int lsame( const char* ca, const char* cb, int lca, int lcb );
  • BLIS: int PASTEF770(lsame)(const char *ca, const char *cb, int ca_len, int cb_len)
  • OpenBLAS: int NAME(char *A, char *B){
  • Reference LAPACK: lapack_logical LAPACK_lsame_base( char* ca, char* cb, lapack_int lca, lapack_int lcb

OpenBLAS also maintains several assembly specializations for DEC Alpha, ia64, power, sparc and x86. Btw, OpenBLAS also seems to have kicked out EBCDIC and Prime OS. Now, to add further to my confusion as a Fortran programmer, in C you also need a version which can accept literals by value:

lapack_logical LAPACKE_lsame( char ca,  char cb )
{
    return (lapack_logical) LAPACK_lsame( &ca, &cb, 1, 1 );
}

Since the function is used a lot in LAPACK, it tripped up the GCC static analyzer, leading the developers to add __attribute__((const)), which is kind of similar to what simple procedures will provide in Fortran 202x but without requiring a preprocessor.

TL;DR Doing a case-insensitive comparison led to a big language mess.

It makes me afraid of what idiosyncrasies/baggage can be found in the rest of the interfaces.

2 Likes
7

Typically we only see the arguments from Big Corps with massive legacy codes that the public never sees. Maybe they don’t want to modernize them. But, for the public-facing codes, it does great harm to the perception of Fortran when people see these codes written in a style or dialect of Fortran that was obsolete decades ago. LAPACK is one of the worst offenders. So my comments below are from the point of view of somebody in the Fortran advocacy camp, not the legacy camp:

Old code frequently has hidden bugs that are obscured by all the line numbers, common blocks, gotos, etc. that are necessary to do anything complicated. With modern-style code it can be much easier to detect bugs.

Improved readability is not nothing. It’s actually very important. Nobody wants to work on FORTRAN 77 spaghetti code. Even the Big Corps must know that they are and will continue to have trouble hiring people to do that.

Totally false.

The sooner the world forgets about f2c the better. This abomination has held back modernization of more than one Fortran library. Fortran has standard C interfacing now. People should use it. Also can remove subtle errors.

Nobody is writing a restricted subset to Python so that it can be converted to another language. Why should we do that? Show what Fortran and do when used correctly and maybe more people will want to use it.

People need to move on. We don’t need to be using RATFOR anymore.

Modernization can remove bugs and make maintenance easier. Fixed to free format conversion can be done nearly automatically in a lot of codes.

People need to update their guidelines if they say to write code like it’s 1977 and to use compilers that haven’t had bug fixes in decades.

9 Likes
8

I think the f2c developers held the view that Fortran is the abomination:

/* f2c.h  --  Standard Fortran to C header file */

/**  barf  [ba:rf]  2.  "He suggested using FORTRAN, and everybody barfed."

	- From The Shogakukan DICTIONARY OF NEW ENGLISH (Second edition) */
3 Likes
9

Speaking of Python there are many efforts disjoint efforts to compile a subset of the language, Numba being the most famous one. There is more than a dozen of Python compilers at this point (e.g. see the list compiled by @certik).

I also think libraries shouldn’t limit themselves to a subset just for the sake of porting to other languages. It sounds like an anti-pattern. The proper solution would be to build better Fortran to .NET CIL compilers. Unfortunately it’s easier said than done.

5 Likes
10

OK, OK, I should have said “nobody is writing Python 1.0 code in order to make it easier to convert to other languages”. :slight_smile:

4 Likes
11

In this specific case, if I’m not reading it wrong:

Unless LSAME is not doing what is supposed to be doing, there is no point in changing its code. As the docs say:

LSAME returns .TRUE. if CA is the same letter as CB
regardless of case.

If the function was named as something like same_letter or is_same_letter I guess that confusion wouldn’t be made and you wouldn’t resort to trying to read the code (or the docs, that may lie). The only scenario that I can see someone revisiting this piece of code is if a new encoding standard takes over and the code stop working (the risk exists, it may occur in a far future, but the odds are too small to consider).

I can see the reasoning behind not rewriting code that is working already (especially ones with such small scope). That’s why risk analysis and planning are so important, so we can move forward and not revisit things that are already solved.

That being said, Modern Fortran is a tool to reduce a lot of risks: missing old school Fortran programmers in the future, keep using tested and reliable software while allowing it to be extended among other stuff.

2 Likes
12

The greatest loss of code I have seen is caused by not being ready for change. If a code is not well documented, lacks automated test cases and uses many obsolete features it will die. Just running it till it dies can work up to a point, just like with an old car you decide to no longer invest in and just run till it quits working. So if you just leave a code as-is because it works you are sealing it’s doom. If you invest continuously in the code at a low level then when it does need changed (and sooner or later everything does) then it will be affordable to add the new feature or port it to a new programming environment. If you leave it as-is long enough it will be easier to write a new one that to decipher and upgrade the old one and all the value of the original code will have been wasted.

12 Likes
13

At least for BLAS, it may be possible to build an assumed-shape variant based upon the future C++ stdBLAS and mdspan container in combination the the C-Fortran array descriptor available in ISO_Fortran_binding.h.

Using GEMV interface from the Fortran 95 Thin BLAS as an example,

interface gemv
   subroutine cpp_dgemv(alpha,a,op_a,b,beta,c) bind(c)
      import c_double, blas_key
      real(c_double), intent(in), optional :: alpha, beta
      real(c_double), intent(in) :: a(:,:), b(:)
      integer(blas_key), optional :: op_a          ! blas_trans, or blas_notrans
      real(c_double), intent(inout) :: c(:)
   end subroutine
   ! ... other real kinds
end interface

I think it will be possible to define a family of routines for different kinds in C++ as shown in the draft below. (Disclaimer: the code below contains future language features and is meant only for demonstration purposes; I have not attempted any verification or validation whatsoever).

#include <mdspan> // P0009r16 - MDSPAN
#include <linalg> // P1673R9 - A free function linear algebra interfaces based on the BLAS

#include <ISO_Fortran_binding.h>

enum blas_key
{
    blas_notrans,
    blas_trans
};

template<typename T>
auto as_mdspan_2d(CFI_cdesc_t *obj) {

    using ext_t = std::extents<std::size_t,std::dynamic_extent,std::dynamic_extent>;
    ext_t extent{obj->dim[0].extent,obj->dim[1].extent};

    using mat_t = std::mdspan<std::size_t,ext_t,std::layout_left>;
    mat_t a{static_cast<T *>(obj->base_addr),extent};

    return a;
};

template<typename T>
auto as_mdspan(CFI_cdesc_t *obj) {

    using ext_t = std::extents<size_t,std::dynamic_extent,1>;
    ext_t extent{obj->dim[0].extent};

    using vec_t = std::mdspan<size_t,ext_t>;
    vec_t a{static_cast<T *>(obj->base_addr),extent};

    return a;
};

template<typename T>
void cpp_Tgemv(const T *alpha, const CFI_cdesc_t *A, const blas_key *op_a,
    const CFI_cdesc_t *b, const T *beta, CFI_cdesc_t *c) {

    using std::linalg::scaled;
    using std::linalg::transposed;
    using std::linalg::matrix_vector_product;

    T _alpha = (T) 1.0
    if (alpha) _alpha = *alpha;

    blas_key _op_a = blas_notrans;
    if (op_a) _op_a = *op_a;

    T _beta = (T) 0.0;
    if (beta) _beta = *beta;

    // assume a, b, c contiguous for now; for non-contiguous or strided
    // array submdspan is needed

    auto _A = as_mdspan_2d<T>(A);
    auto _b =    as_mdspan<T>(b);
    auto _c =    as_mdspan<T>(c);

    matrix_vector_product(
        scaled(_alpha, _op_a == blas_notrans ? _A : transposed(_A)),
        _b,scaled(_beta,_c),_c);
}

/* ... instantiate function templates for different 
       real kinds: float, c_double std::complex, ...
*/

I’m not sure if there is any communication between the C++ and Fortran committees on this topic.

(For brevity, handling of optional arguments can be deferred to a small function template.)

2 Likes
14

one thing missing from this conversation so far is a labor of examining the assumption that the algorithms of these old code bases are good algorithms that are worth modernizing rather than just rewriting the program from scratch. of the author was a good programmer, they wrote a program for a computer that is very different from a modern one. it was written for a computer without vectorization, where memory accesses and if statements are completely free, fma doesn’t exist and multithreading doesn’t exist. they also don’t know about the last couple decades of theoretical advances in methods to solve the problem. porting this code will make it easier to read, but you need more than a port to fix the fact that it’s probably a bad algorithm implemented poorly for modern architectures.

5 Likes
15

The reasons listed are entirely bogus.

Fortranners specifically and anyone in the scientific computing generally will be far better off following the axiom - refactor each and every codebase that has coding practices which you assiduously will not adopt in any new code.

Computational approach now has advanced enough to be seen as the veritable third leg on which rests the scientific advancement. Fortran plays a major role in this and with suitable extension of the language standard, Fortran can make a bigger and more profound contribution to science and technology.

In the scientific pursuit then which is also aided by broader technical and engineering computing wherein lies LAPACK, if there is anything more inane than having to teach implicit none to all of posterity, it is to misuse the argument of “if it ain’t broke, don’t fix it”:

  • if a codebase relies on FORTRAN I style implicit mapping, it is absolutely “broke”!,
  • if a codebase uses external procedures without explicit interfaces and ignores rank considerations with actual/effective and received arguments, it is “broke”,
  • if a codebase employs assigned or arithmetic GO TO, it is indeed “broke”,
  • if a codebase branches to an END IF or uses a non-block DO construct, it is “broke”,
  • if a codebase uses assignment on declaration without explicit SAVE, it is “broke”
  • … (you get the idea)
3 Likes
16

LFortran will very soon be able to fix most of these basic modernizations almost for free:

  • fixed-form → free-form
  • implicit typing
  • implicit save
  • transformation to C and C++ of any modern Fortran (so no need for f2c)

Just by transforming the internal representation (ASR) back to Fortran.

And the following with an optional ASR->ASR transformation pass first:

  • goto (at least the most common cases)
  • implicit interfaces
  • Moving global subroutines/functions to modules
  • it will also be able to transform your modern code to a stricter subset to not lose support for legacy compilers and platforms

The workflow will be that you compile your code and make sure it works. Then you tell LFortran to rewrite it, and it will work, modulo compiler bugs, that we will fix as we encounter them.

And I can imagine if things are successful, we’ll get more people contributing and we can write more such ASR->ASR passes to do more complicated refactorings. But even the above, which technically is relatively straightforward, will be a huge help.

Once it is as easy as just adding a compiler option to a compiler that you already use, I think projects like Lapack will be much more open to do this.

15 Likes
17

I used two versions of optimization TOMS algorithm. Modernized one was done by experienced scientist and Fortran 95 programmer Allan Miller, who followed aforementioned approach and used his tool to make modernization. Finally, i quickly discovered the old one in f77 worked perfectly while the new one failed giving hard to trace down error. I gave up after few days and used the original implementation.
The second example concerns potentially efficient whole array constructs like where. They were never been as efficient as explicit loops since the latest were easily recognized patterns for advanced automatic compiler optimizations. My take is: be careful with modernization

2 Likes
18

No offense, but this is an “other than that, how was the play Mrs. Lincoln” situation. HPC codes have the unique challenge on top of general man-hour concerns that you’re dealing with highly specialized code that, in most cases, require a PhD in a highly specialized field to understand, and the twin terrors of the pay structure and incentive structure of academia further restrict the labor force.

It’s hard to refactor a code when you keep losing all your workforce that’s skilled enough to effectively refactor code to FAANG.

Avoiding this elephant in the room, the other issue is that when any sufficiently battletested programmer hears the word “modernization”, they reach for their gun. Fortran is a special case here, since we still have FORTRAN to deal with, but today’s “modern code practices” have a funny habit of being tomorrow’s “code anti-patterns”. To use an example relevant to F20XX, I’m old enough to remember back when exceptions were the Future of Error Handling, now everyone avoids them like the plague even in cases where they’re perfectly valid.

8 Likes
19

The lesson for Fortran people here is not to just run some fixed2free script, dump the result on a website and never touch it again. What was the bug? Probably something that can be easily fixed. If it was on GitHub or GitLab we could have bug reports, versions, releases, etc. New algorithms and methods might be incorporated over time. Maybe new language features that allow for better or faster code, etc…

2 Likes
20

Legacy code was used by… generations of programmers. FOOPACK packages passed the test of time again and again. Yes, it is fixed form. Yes, it is spaghetti. Yes, I would not start a new project in Fortran 77 or earlier. Yes, trying to read and understand legacy code becomes terrifying very quickly. But it works.
Have you ever tried to compile huge packages like Scilab? It requires Fortran installed and it does compile legacy Fortran code, which becomes a Scilab command in the end. Matlab’s syntax for pretty much everything related is suspiciously similar. I can think of other packages doing the same. Legacy code is everywhere.

Modernizing such code is enormous work, extremely prone to introducing bugs, and feels like reinventing the wheel - a wheel that is not broken, I might add. This is Fortran, after all. Notorious for backwards compatibility.
What I do when I have to use legacy code is to create a Fortran 90+ module which defines “human friendly” versions of the legacy procedures. Not only that but I “strip” the legacy code to the part needed for the task at hand. For example, ODEPACK needs some of SLATEC; I strip both ODEPACK and SLATEC procedures to the absolute minimum for my task. Legacy code is located in a subdirectory and it is “seen” privately by the module only. Everything else in the project just uses the module’s public modern procedures and it is completely unaware of the “offending” legacy code. It takes some time to implement such a module, but it pays off. You do it once, you use it forever. And it is way quicker and safer than modernizing the legacy code itself.

I would find useful the exact opposite instead. But I can think of reasons it can’t be done.
I never used f2c (except out of curiosity), and never will.

4 Likes