Intel® Fortran Compiler
Build applications that can scale for the future with optimized code designed for Intel® Xeon® and compatible processors.
The Intel sign-in experience has changed to support enhanced security controls. If you sign in, click here for more information.


Valued Contributor II

Many years ago, with a lot of help from the good people on this site, one in particular,  we got a structural analysis package running based on code from two good sources. 

I was looking at a unique engineering problem in Europe that for the first time allows me to use the full range of the package, beams and plate elements in a close approximation to reality.  

I have a test file that I am slowly working on to set up a sound model, I run the model as I add a few nodes so I can check the engineering data.  Very easy to make a mistake and blasted hard to find if you try all at once.  

I have been using IFORT, along with Pardiso and the Eigensolver.  I get the following answers, Screenshot 2022-08-08 100243.png

If I make a change to IFX I get the following answers, 

Screenshot 2022-08-08 100154.png

The difference is not exactly 2, but about 1.9 +- 0.2 in the ratio.  

When I compile with IFX I get 



I realize IFX is a work in progress, I have no idea how to find out, what these errors mean in terms of line numbers,  I still have a day of work to set up the final model,.

I am merely pointing out the problem.  

I have a TB or so of measured data, so I will shortly know the real measurements. 


0 Kudos
41 Replies
Valued Contributor II

The 32-bit MKL is a separate download.  It is the third file that is now added to the HPC download, - I think. 

Each square local matrix is say (6 by 6) it is placed in the large stiffness matrix at offsets, so you end up with a mainly banded matrix.  

The beams program from Harrison did not have plates, the problem is plates are useful for modelling some bridge decks.  

    ! PURPOSE Form stiffness of 3-node, 18-dof flat triangle shell
    ! AUTHOR C. A. Felippa, October 1996
    !    Converted to Matlab from Mathemadea, by S. Gohnon, N. Ledford,
    !      L. Liu, Apri] 2009
    ! KEYWORDS Kirchhoffthin shell isotropic homogeneous finite element triangle stiffness matrix

I converted some elements of it back to Fortran and C Felippa sent me some earlier Fortran code.  

I added beams and plates, it potentially can solve some problems well, others not, but for this bridge, it is theoretically ideal.  

From FEAST you can get the equation numbers with the eigenvalues, so you have some idea of what is affected by the frequency.  

We have frequencies that are temperature dependent and some that are not, so that is a long potential research area. 

@mecej4 developed the good type stuff used in the program.  

I have a more advanced beam version, but this is the version that runs the plates and beams, the answers have not been checked except against the small Felippa problem.  I now have a solvable problem and some good data, so it is worth looking to see if it works.  

Does that make sense?   Everything is designed to run small, not huge FEM. 



Black Belt

JohnNichols, I looked at the source files that you attached, but there are some problems.

  • You use many lines longer than 132 characters, even in the fixed-form files. This is allowed by Intel Fortran, but becomes an impediment when I want to use some utilities for checking the source for potential errors.
  • Some of the source files contain modules with CONTAINed subprograms, but many of those do not have the required form END SUBROUTINE.
  • In the routine that assembles file names, you treat arrays such as character(1) :: fname(14)  as if they are equivalent to character(14) :: fname.
  • In subroutine SM3MHFF, you have a "DO 2500" statement with no index variable or limits. I think that you need to add "j = 1, 3".
  • In subroutine READBEAMS, the variable DONE is used before it is initialized.
  • The data file is opened first on unit 9 in subroutine FILES; subsequently, with this connection still present, an attempt is made to open the same file on unit 14, in subroutine ReadBeams.
Valued Contributor II

@mecej4, are you using FLINT? 



Black Belt

No, I do not have Flint (I tried it several years ago).

Instead, I use a number of tools, which together serve better to catch bugs than any one of them. That is why I prefer that the program under study be in standard format, so that I do not have to spend time bringing the code into that format in order to be able to apply the tools.


Valued Contributor II

I have just amended borr so it has no lines longer than 132 characters, and made sure all of the subroutines have a proper ending, and now working on the rest.  



Valued Contributor II


  • In the routine that assembles file names, you treat arrays such as character(1) :: fname(14)  as if they are equivalent to character(14) :: fname.

This code was written in 1988 when I developed a sewer design and a pipe network analysis package for waste sewers and storm drains. At the time there was only AutoCAD to draw, it was slow and tedious, almost as slow as by hand.  

As I developed the code, there were only people at the Uni who coded, so it was quite alone, I found that it helped if I used a single name of the files for each job and just changed the extension.  We used job numbers.  So I developed a set of routines to open and close the files just from the job number.  

VFILE is the job number with a maximum length of 9, remember MS-DOS  3.0 MS Fortran 3 and large floppy drives.  XFILE is the extension.  

FILECRE is passed  Character*14 and assumes it is character*1 fname(14). I use Fname to assemble the vfile and xfile into a full file and return them to FOPEN to open the file using FNAME.  

It has always worked.  Other than to add implicit none, it has not changed since 1988.  



C      ****************************************************************
C      ***************************************************************
      Implicit None

      Integer I,J
      Character*1 Xfile(4)


      character dh
      character dh1



      IF (J .EQ. 1 ) THEN
          IF (.NOT. EXISTS) THEN
              !CALL CLS
90                FORMAT( ////////,' Unable to find file => ',(A))
                  STOP ' Please create the data files.'
          ELSEIF(J .EQ. 2) THEN

      END subroutine FOPEN

C      ****************************************************************
C      ***************************************************************
      Implicit None

      File_CODE = 1

      IF(I .EQ. 1) THEN

      END subroutine CLOSE_FILES

C       ***************************************************************
C       ***************************************************************
      Implicit None



      character dh
      character dh1


      DO 10 I=1,9

          IF (vFILE(I) .NE. DH1) THEN
          ELSEIF (vFILE(I) .EQ. DH1) THEN
              IF (I .EQ. 1) THEN
              ELSEIF (I .EQ. 2) THEN
              ELSEIF (I .EQ. 3) THEN
              ELSEIF (I .EQ. 4) THEN
              ELSEIF (I .EQ. 5) THEN
              ELSEIF (I .EQ. 6) THEN
              ELSEIF (I .EQ. 7) THEN
              ELSEIF (I .EQ.  THEN
              ELSEIF (I .EQ. 9) THEN
30            FORMAT(A4)
              GOTO 40


20    FORMAT(A1)


      end subroutine FILCRE


Valued Contributor II

J=1,3 fixed, thanks.

This routine is not used, I am only interested in a full plate element with 18 of freedom.   But IMHO, line 53 should throw an error, there is no control structure.  The program does not use it, but it compiles ok.  



 subroutine    SM3MHFF (x, y, dm, f, ls, sm, m, status)

    implicit none
    !                              ARGUMENTS
    character*(*)  status
    integer       ls(9), m
    double precision x(3), y(3), dm(3,3), f, sm(m,m)
    !                  LOCAL   VARIABLES
    double precision xc(3), yc(3), dxc(3), dyc(3), hh(3,9)
    double precision sqh(3,3), qx(3,3), qy(3,3), r(3,3)
    double precision  area, area2, lambda, cj, sj, csj
    double precision e11, e12, e13, e22, e23, e33, jxx, jxy, jyy
    double precision det, gamma, ggg, mu, mux, muy, mumu, tau
    double precision sum, s1, s2, s3, s4, s5, s6, x0, y0
    integer i, j, k, l
    status =   ' '
    area2 =    (y(2)-y(1))*(x(1)-x(3)) - (x(2)-x(1))*(y(1)-y(3))
    if (area2 .le. 0.0)      then
        status = 'SM3MHFF: Negative area'
        if (area2 .eq. 0.0)   status = 'SM3MHFF: Zero area'
    end if
    if (f .eq. 0.0)         return
    x0 = (x(1)+x(2)+x(3))/3.0
    y0 = (y(1)+y(2)+y(3))/3.0
    area = 0.5*area2
    lambda  = 1.0/sqrt(area)
    xc(1) =  lambda * (x(1)-x0)
    xc(2) = lambda * (x(2)-x0)
    xc(3) = lambda * (x(3)-x0)
    yc(1) = lambda * (y(1)-y0)
    yc(2) = lambda * (y(2)-y0)
    yc(3) = lambda * (y(3)-y0)
    dxc(1) = xc(3) - xc(2)
    dxc(2) = xc(1) - xc(3)
    dxc(3) = xc(2) - xc(1)
    dyc(1) = yc(3) - yc(2)
    dyc(2) = yc(1) - yc(3)
    dyc(3) = yc(2) - yc(l)
    e11 = dm(1,1) * f
    e22 = dm(2,2) * f
    e33 = dm(3,3) * f
    e12 = dm(1,2) * f
    e13 = dm(1,3) * f
    e23 = dm(2,3) * f
    jxx = -2.*(xc(1)*xc(2)+xc(2)*xc(3)+xc(3)*xc(1))/3.0
    jxy = (xc(1)*yc(1)+xc(2)*yc(2)+xc(3)*yc(3))/3.0
    jyy = -2.*(yc(1)*yc(2)+yc(2)*yc(3)+yc(3)*yc(1))/3.0
    do 2500 
        mux =-3.0*xc(j)/2.0
        muy = -3.0*yc(j)/2.0
        mumu = mux**2 + muy**2
        mu = sqrt(mumu)
        gamma =    2.0/mu
        tau = 0.25D0*(dxc(j)**2+dyc(j)**2-gamma**2)
        ggg =(mumu-3.0*tau)*gamma*lambda/24.
        CJ = mux/mu
        sj = muy/mu
        r(1,j) =    -lambda * (cj*xc(l)+sj*yc(1)) + ggg
        r(2,j) =    -lambda * (cj*xc(2)+sj*yc(2)) + ggg
        r(3,j) =    -lambda * (cj*xc(3)+sj*yc(3)) + ggg
        csj =        cj*sj
        qx(j,1) =   -0.5*csj*cj
        qx(j,2) =   -0.5*sj**3
        qx(j,3) =   -csj*sj
        qy(j,1) =    0.5*cj**3
        qy(j,2) =    0.5*csj*sj
        qy(j,3) =    csj*cj
        s1 =      e11*qx(j,l) + e12*qx(j,2) + e13*qx(j,3)
        s2 =      e12*qx(j,l) + e22*qx(j,2) + e23*qx(j,3)
        s3 =      e13*qx(j,l) + e23*qx(j,2) + e33*qx(j,3)
        s4 =      e11*qy(j,l) + e12*qy(j,2) + e13*qy(j,3)
        s5 =      e12*qy(j,l) + e22*qy(j,2) + e23*qy(j,3)
        s6 =      e13*qy(j,l) + e23*qy(j,2) + e33*qy(j,3)
        do 2400 i = 1,j
            sqh(i,j) = jxx * (qx(1,1)*s1+qx(i,2)*s2+qx(i,3)*s3) + jxy&
            * (qx(i,1)*s4+qx(i,2)*s5+qx(i,3)*s6+ qy(i,1)*s1+qy(i, 2) *s2+qy(i, 3) *s3)&
                + jyy * (qy(i,1)*s4+qy(i,2)*s5+qy(i,3)*s6)
2400    continue
2500 continue




Black Belt

The whole subroutine FILCRE is rather superfluous. Instead of 


you can simply write

fname = trim(vfile) // trim(xfile)

and discard the entire subroutine FILCRE .

There is a major error in the program. The arguments in the call

call LinSolve(count, nnz, error, atemp, jatemp, iatemp, b, X)

do not match the argument list in

subroutine linsolve(n, error, nnz, a, ja, ia, b, x)

If this mismatch actually exists in the code that you ran on the various problems, it follows that any conclusions regarding the calculated results are invalid.

Blunders of this nature can be caught and fixed if you use a compiler with good error-checking capabilities.

Valued Contributor II

It has only been compiled on Intel Fortran Compilers since 2006 or thereabouts.  It is a nasty error, thanks, but why did IFORT not flag the error?  

What do you consider to be a compiler with good error-checking? 

The trim - thanks 

I wrote FILECRE in 1988, it is a personal friend.  

Valued Contributor II

Answer - The compiler will not detect two integers swapped in the call statement.  I added check routines. 



@JohnNichols we use NAG Fortran on our team for checking.  Most would agree that it's the gold standard in parsing and diagnostics of the Fortran Standard(s).   

Valued Contributor II

1.  @Ron_Green   =====   My daughter thought I was dying from laughing so hard.  She said no mouth to mouth, it was time I was gone.    She is her father's daughter, but she speaks Chinese.  

2.  I tried to change my ID NAME to NAG, but that part of your website does not work, I tried thrice (三次).  Chinese write the low numbers sideways.  

3. I tried to install FLINT and it throws an error, I asked the FLINT people and they said - we have not seen that error before, you can have a 30-day license if you help with the error.  

4. NAG does not tell you a price, so I have no idea, but generally, that means I cannot afford it.  

5. Why do all these Fortran websites, always have a Japan option? Is Fortran popular in Japan?

6.  .done.   -- @mecej4, I ran the code and the logical variable done is created as a false. Not sure if this is standard, I added a line to make it false before I use it. 

7.   I use two files for special things, base.f90 is the place I store all the things you showed me including this stuff.  It is used in several programs, so there are routines unrelated to this program in there, hence the duplicate open.  Only one is used and I added a close so it is complete.    

8. The same problem can occur with Scotia.  When I started to write sewer and drain software, which made us a lot of money in the consulting office as we were ahead of the wave, I had a library for common files, called Scotia.  This is the name of the ship my family sailed from England to Australia in 1849.  These routines can be very old school. 

9. I know I have to be careful using base and scotia.  

Module Base

    INTEGER, PARAMETER :: dp = selected_real_kind(15, 307)
    INTEGER, PARAMETER :: dpA = selected_real_kind(15, 307)

    INTEGER, PARAMETER :: sw = 2                    !   Output file
    INTEGER, PARAMETER :: srA = 15                  !   output.txt file
    INTEGER, PARAMETER :: st = 14
    INTEGER, PARAMETER :: sa = 3                    !   Output file
    INTEGER, PARAMETER :: smWrite = 4
    INTEGER, PARAMETER :: si = 1
    Integer, parameter :: slog = 9                  !   Log file
    Integer, parameter :: nta = 100                  !   Log file
    Integer, parameter :: outNode = 63                  !   Log file
    Integer, parameter :: inNode = 0                  !   Log file

    REAL (KIND=dp), PARAMETER :: gr = 9.806, pi = 3.14159265D0  !   Standard parameters
    REAL (KIND=dp), PARAMETER :: delta = 0.0001                !   Error number of checking for zero
    REAL (KIND=dp), PARAMETER :: ZERO = 0.0_DP
    REAL (KIND=dp), PARAMETER :: ONE = 1.0_DP
    REAL (KIND=dp), PARAMETER :: TWO = 2.0_DP



Black Belt

I notice that you read in from the data file a variable "MaxNodeNumber".

Are (is) the nodes referenced as (0:MaxNodeNumber) or (1:MaxNodeNumber)?

If former, then your (1-based) allocations need to be bumped by 1


Jim Dempsey

New Contributor II



Some areas where converting old programs like Ularc to 64-bit could be:

  • When changing from COMMON array sliceing to ALLOCATE, the past code approach could go out of bounds on array subscripts and address other parts of the shared COMMON array with little damage. Now with ALLOCATE this can be more problematic, especially with a subscript=0, which can be a commonly missed error in old codes. You may need to check for this !
  • When converting to 64-bit, typical 32-bit subscripts can now overflow. This was not checked in old code, as this was not anticipated.
  • You were previously experincing errors of 1.e-8.  This may be where 32-bit reals are now being assumed for real constants, where they were previously assumed 64-bit, especially for old reference results ?
  • Mixing of character data in real arrays can change with compilers, producing unexpected IF results. Ularc et al had a lot of this.

However, it is impressive how these old codes still function so well. I recently did an eigen-vector extraction for 300 vectors from 450,000 equations and no sign of your 12^2 problem ! Those old algorithms still work well in 64-bit.

Valued Contributor II

Dear All:

Thank you for your time and patience.  

John: Ularc has one advantage, it can very quickly tell you the safety of a steel structure, from an engineering perspective it is really good.  It is very old, I have updated it and I still use it, one uses it with caution on a limited number of problems.  So may I enquire as to the eigensolver you used?  It told me that the designers of the bridge stopped the cross support at the point of maximum likely moment and plastic failure.  Utility trumps design.  Powell's group developed wonderful engineering concepts, which we still use, it is also available for 25 per year.  They are all gone, but the shadow lives on. 

Jim:  Yes, I agree, that is one of the errors, I spent two days reaching that conclusion, you spent a few hours. I looked for the original paper to remember the plate model. Interesting it is only found on one place on the internet, the FORTRAN INTEL FORUM, I placed a copy there in 2009, when I added plates to the Harrison program, after @mecej4 very kindly helped me, (he did the lion's share I merely tagged along), change the old Harrison code to a typed system. 

If you stick with beams this program will solve a lot of interesting problems, quickly.  It is not perfect, but it is not $3000 a copy.   The problem is the old sources for the code and the papers that were written are disappearing as they die and their owner's websites are removed, the internet is not a library, it is a window and the window moves. 

@mecej4 , I often have an internal debate late at night as to whether, you, Jim or Steve are the best of programmers and people one can stupid questions.  I do not have your skills in looking at code, I accepted my failing a long time ago, but I have the second best thing, a person I can ask questions, who always answers, sometimes using interesting words, but who is patient and kind, as are Steve and Jim.   So thank you for the comments, you have helped me with this program over many years, I added plates in 2009, and then all of the subsequent bridge problems were such that they were of no use. So I have not used that code since, nor was it vetted,  as with all code in engineering, you start to use it and then you solve the problems. 

If you look at the first iron bridge it is a slender masonry bridge with everything in compression, he had no idea how to design it, he just built it copying older masonry stuff.  The English accepted that bridges failed, as they increased the mass of the trains, finally sense prevailed and they started to use standards.  

I was once at a conference, where a very old distinguished engineer commented on the use of 3 as a factor of safety,  he asked why in a roundabout fashion. I did not have the heart to tell him the real story of the British lawyer, who asked 7 of the world's famous engineers in about 1850, what factor of safety to adopt, added the numbers, divided by 7 and came up with 3. You can still read the paper, it is on the UK government website somewhere.  

Thank you, I will correct the errors.  

Finally, people ask challenging questions,  but if you are going to ask challenging questions this is the place to do it, by far for Fortran. 

I have no idea what the gfortran compiler is, or where to find it, I have been an MS ---> to Intel Fortran human since 1986, and I am not going to change now.  Grammarly wanted to change GFORTRAN to gfortran, and then to Fortran.





New Contributor I

hi john,


i just wanted to comment on gfortran. i wouldn't bother messing with it. intel fortran produces exe files that run 2-3x faster in my testing. i had been comparing them for about 12 years before i finally stopped. there are other issues with gfortran that are really annoying as well. one is that the documentation does not tell you that many of the features only work on linux. you find out the hard way as a windows user. the other problem is it use to be very easy to download and install gfortran on windows. not any more. it's now very complicated and a rolling release. so stuff that worked will suddenly stop working etc. one instance for me was the static option. i have to have that working to distribute my code. i finally got fed up and found intel fortran was now free. since then, i just use intel fortran. for a time, i also used the portland group fortran compiler. it didn't create exe files that ran as fast as intel fortran exe files. but they ran faster than gfortran created exe files. the nice thing with the portland group compiler was it didn't have a preference for amd or intel cpus. so it was easier to understand how to compile code. however, nvidia bought them and killed the program. now they use something like ifx. however, they don't have a windows version available. so again, intel fortran is definitely the way to go. the only real issue i have with intel fortran is it's pretty hard to figure out how to set the compiler options so that the exe files you create will run equally well on amd and intel cpus. i think i finally figured that out though.



Valued Contributor II


I have spent an enjoyable day buried deep inside the program.   Like anything you leave for a long time, it takes a while to get back into it.  

I do not use anything but Intel chips, life is too short for playing with some things.   

I just got four new Intel NUCS, Steve was correct, the core i3's are out, relatively cheap and incredibly easy to set up, one hour from opening the box, assembling and running VS.  The craftsmanship is superb.  Pity they are going to be stuck in a 250 yr old boat in a museum.  

Thanks, mecej4, as always you are great. 

No idea how that happened.  


Loading the submatrices into the main matrix has been fun, your comments helped a lot, again thanks.  

I still have to check the beam additions, something is still funny.  






Black Belt

I just noticed another source of trouble in subroutine Read_Beams, in file Elements.f90:


    CHARACTER *4        af                  !   Read input line type - TRIA
    CHARACTER *120      iline               !   Input line
    EQUIVALENCE         (af, iline)         !   Simple dual read method
    READ (si, '(A)', ERR=1000, END=400) iline    ! Read line
    READ(iline,*, err = 400, end = 1000)af, IV, MCON1, &
       MCON2,GAMMAA,beamID,Frame, X1,Y1,Z1,X2,Y2,Z2,h,j

The Fortran standard prohibits I/O to/from an internal file when the file overlaps with one or more items in the I/O list. In your code, af overlaps part of iline because of the equivalence declaration.

Valued Contributor II

It works quite well, I am sure I learned how to do that step on this site.  

It can be changed.

The code at the moment is useful for engineering, it is not perfect, but Jim and you have helped make it a lot better, I leave flaws, but the main element is that the answers are checked against standard answers.  That is the process at the moment. 

Dr. Golmon published a paper on the use of this method and provides some samples.  In reading the paper, there are some interesting questions, and STRAND7 does not get her answers, close but enough to wonder.  An 18 dof plate is a useful tool, and the problem is people build weird bridges and they want to model them and it has to be fast and cheap to work, otherwise, it is not done.  

QA is critical. 

Chris Niggeler, CSS, has kindly offered me a month-long license on FLINT.  I warned him I would talk about it here.  We shall see what it picks up, once I get it properly installed.  

Why do arrays in IFORT initialize to zero and scalars do not? I zeroed out all the variables before I use them, which popped up a few errors.  

I found some elements of SMSHELL on a git hub site, it is a huge geotechnical FEM.  I have yet to look at it.  

This would not be possible without all I have learnt on this site.  

PURPOSE Form stiffness of 3-node, 18-dof flat triangle shell
! AUTHOR C. A. Felippa, October 1996

This work should not be lost, at least the code is safely available here, with extensive comments.  

It is translated to Mathematica to make teaching easy for the professor, if you remember Professors are lazy -- they would say highly efficient.  The academic system is not a good method to produce useful long-term results, one day it may be fixed, but this is what happens when you use a 15th-century method in the 21st century.  PS the new hellfire missile, almost certainly parts of it designed with Intel Fortran, is just a stone from a modern trébuchet. 

My students in 2002, built a trébuchet and used bowling balls to demolish a dunny. (Australian vernacular).  A trébuchet is a dangerous tool. 


Black Belt Retired Employee

NAG Fortran, though NAG is a UK company, is developed in Japan (the lead developer is the Fortran standard's editor). Yes, Fortran is popular in Japan. And yes, it is paid software. NAG Fortran's diagnostics are second to none, and I often wish Intel Fortran's messages were a tenth as good, though I understand why they aren't.

Valued Contributor II

>>>though I understand why they aren't.<<<


may I gently enquire as to why? 


I am happy with Intel Fortran, nothing is perfect, but it is good software and I enjoy working with it.  

I have been through the program using VS debugger, looking at the output at all stages, I now have IFORT and IFX giving the same answers, the Gods be praised.  

The code is extremely sensitive to the ordering of the plates and nodes.  If you do not follow the pattern the program just explodes. 

Dr Golmon, who wrote the standard paper on the code in Mathematica, now works at Draper.  I sent her a polite request to talk about the code.  I can approximate their standard sample, but it is not exact.  

The frequency results are very sensitive to Poisson's ratio as I discovered when comparing to STRAND7 results.