Re: Harrison 1967 Continued

JohnNichols · ‎01-29-2026

I enjoyed your comments. I started a new thread as the old one was so long I could not find anything.

I can just remember the end of those old days with the early programmers, but we dealt with those issues probably up until the early 1990s when memory was removed as a "real" issue.

Harrison is from Sydney University, they developed a lot of steel structures software, which was according to an old academic friend incorporated into a 1980's version of the Australian steel code. He said this with some interesting comments.

Harrison and others of the time pointed out issues that we are now grappling with in Structural Reliability, it is just a long slow process. The graph shows the amplification from the elastic to the second order in results, these are significant for Harrison's Simple Frame 1.

Ularc is from a similar vintage at UCB, it clearly states no second order, but I use it for Structural Reliability. Time to move up a notch, as Harrison noted it should be easy to add plastic to ELSO, or fix ULARC. So this is just the best place in the world to ask questions about Fortran, and I unlike a lot of people have only every used CDC or Intel.

mecej4 · ‎01-30-2026

The zip file that you attached contains versions of the Fortran source file(s) for three obsolete platforms: (i) CDC mainframes (ii) DEC VAX supermini-s and (iii) 16 bit IBM PC. I tried compiling and running those using (i) IFX 2024.2 and (ii) GFortran 13.4. I found that for one reason or another it was either impossible to generate an EXE that gave an error-free run, or to produce an EXE at all. For the few cases when I could run and display program output, the numbers in the results were quite inaccurate, or even absurd.

Please consider providing just one source file that can be used to generate reliable output using current compilers. It would also be helpful to have reference output files for one or two test problems to verify correct functioning of the code.

I recognize that the program is 59 years old, and the world has changed considerably. Nevertheless, with a little bit of work it may remain useful for another decade or two.

JohnNichols · ‎02-03-2026

This is as close to the original UCB 72 code as currently exists.

The input and output file are in the zip as forum does not allow .in and .out files.

The problem comes from ULARCA.pdf, which is a copy from UCB, I had an original copy of the book which I cannot find, which was readable.

The PDF has some pages in the problems that cannot be read, but the author gives us the source and I got the Hodge source over the weekend through Interlibrary loan, so I can recreate some of those pages.

A Sudhakar died in the 1990's, I corresponded with his wife. He did not write the program alone, some of the code is common to AXISHELL, which is a 1968 program, written by a engineer who was at Halyburton in Houston, 20 years ago. Axishell is a great program and I also talked the author, AXISHELL designed 4000 tonne coal bins in the early 1980's and reduced the steel mass from hand design a lot.

To get the code running on the latest ONEAPI in 64 bit I had to fix two small errors, the lingen subroutine variable prop did not match the calling routine type and an input file read line uses N twice, I swapped the second n for (k=1,NJST).

I have used this program many times, it has been updated in several stages, I do not touch the stuff that A Sudhakar wrote unless I absolutely have to change it, so the common is left as is. He wrote it that way and it works, who am I to change it.

This ran this morning on VS 2022 the latest and the latest ONEAPI. I have not touched this one for several decades, but it is the original.

Now the revised latest will calculate Beta values. This is the latest MC version, it still uses the common's.

Finally you never trust a FEM code, so you always do it with two programs, such as STRAND7 or the Harrison code to ensure the answers are approximate and the answers are an approximation.

mecej4 · ‎02-04-2026

I am puzzled by this post. You attached two files, which appear to be nearly identical. ULARC.FOR includes a main program, whereas ularcsub.for contains no main program nor Subroutine INPT. Compiling and running ULARC.FOR and typing in "ularc.in" as the name of the input file causes an I/O error: "Invalid decimal character L was detected".

You mention "Monte Carlo", and your red screen-copy does, too. I do not know what relevance the phrase has, since I have been under the impression that the ULARC program does an elastic analysis of a plane frame under static loads.

JohnNichols · ‎02-05-2026

The setup in a solution file for the two Fortran text files is :

This running 64 bit in debug mode on an intel compiler, I do not use anything but the Intel Compiler, which flows from the Microsoft 3.13 or 3.03 compiler last seen in the late 1980's and travelling through several iterations is now oneapi. It is just habit, but I have no desire to move off a program that generally works and is well supported here.

If I build it then this is the output, it is latest compiler and I use IFX, because IFORT gives me an annoying message.

It starts and asks for the input file, I have used \ as the stop mechanism for the input file since Microsoft Fortran. if this is unsupported in other languages I am sorry, the original program used other input methods for UCB Main computer,

I supplied ularc.dat and ularc.out, this version has the "new" on the data file so it stops you overwriting old output files, change it to unknown if wanted.

It runs in 64 bit mode, the output is as expected for three load combinations the plastic iterations are shown by the call force repeats.

This is the structure in metres, so it is a large warehouse, the nodes are

COMPLETE NODE COORDINATES

NODE	X-ORD	Y-ORD

1	0	0
2	0	2.167
3	0	4.333
4	0	6.5
5	1.97	6.82
6	3.411	7.055
7	4.852	7.29
8	6.294	7.525
9	7.735	7.76
10	9.176	7.995
11	10.617	8.23
12	12.059	8.465
13	13.5	8.7
14	14.941	8.465
15	16.382	8.23
16	17.824	7.995
17	19.265	7.76
18	20.706	7.525
19	22.147	7.29
20	23.589	7.055
21	25.03	6.82
22	27	6.5
23	27	4.333
24	27	2.167
25	27	0

One of the failure modes is at node 22 as is expected, it has a lambda as is tradition at 1.594, it needs to be above 1.696 for safety, so the design would need to be modified.

It is a simple system with a simple failure mode. Ularc uses rotational stiffness at the restraints which is a great idea.

This is the UCB version, since then it has been upgraded to allow for stress analysis, it does not allow for stability functions as for ELSO, so the deflections are linear.

If you are designing for coal mines etc, one is not limited by floor deflections, a person driving one of five D10's bull dozers on a suspended slab over a reclaim tunnel is not really worrying about floor movement for a school.

The last iteration has a Monte Carlo analysis routine in Intel Fortran with MKL for the RNG which is a pain as you need to take care with the seeds, the output is the red shot. But this method allows one to calculate a Beta for the range of loads and not just one single standard load. Of course you need a FEM that calculates frequency and as a check on the ULARC. Which is the 3D Harrison program from several years ago, or Strand7 or Abacus in extremis.

It depends on the money and the available time, if you have a weekend for a large bridge one uses ULARC.

I hope this explains this program, I was just interested to show the old UCB code as it was, and now runs with limited changes.

witwald · ‎02-11-2026

@JohnNichols, in the the ULARC72.FOR code, there is a line that outputs the displacements in millimetres. To do that, it multiplies all three elements of array data by 1000 to convert from metres to millimetres. That includes the ROTATION parameter.

Is that correct? I'm unsure as to why the ROTATION parameter would need to be multiplied by 1000.

witwald · ‎02-14-2026

I've taken the ULARC72.for and ULARCSUB.for code and amalgamated it into one file, and created the ULARC72_refactor.f program. This new version has been refactored as follows:

Now IMPLICIT NONE is used throughout.
All Hollerith edit descriptors have been removed and replaced with character strings.
All FORMAT statements have been removed, allowing the removal of their associated labels.
All GOTO statements have been removed, resulting in a more linear program flow through the code. As necessary, the GOTO statements have been replaced by:
- IF blocks
- DO WHILE blocks
- EXIT statements
- CYCLE statements
- ENDDO statements
The program now defaults to expecting 'ularc.in' and 'ularc.out' as the names of the input and output files.
The value Pi/180 is now computed instead of being hardcoded.
The variable IBUST has been renamed ICOLL to better reflect its attachment to identifying collapse states.
The formatting of the output has been reorganised to make it tighter and less voluminous, as well as making different sections easier to identify when looking at the results.

The new program is attached. It has been tested against the sample output file provided earlier, and its results match up with the old one. Hopefully the test data exercised a fair bit of the code.

I trust that ULARC_refactor.f may prove useful to somebody in the future.

JohnNichols · ‎02-15-2026

Really neat.

    Module Base
    
    implicit none

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

    INTEGER, PARAMETER :: sw = 16, sr = 1, st = 11, sa = 3, sm = 4, sc = 12, sRT = 25
    Integer, parameter :: slog = 9, sOUT = 10, sCAD = 21, sINI = 13
    Integer, parameter :: array_size = 300
    integer, parameter :: nt = 120
    integer, parameter :: mt = 200          !   Number of members    
    integer, parameter :: mn = 2    
    integer, parameter :: ml = 3
    integer, parameter :: mn_12 = 12
    integer, parameter :: mn_4 = 4
    integer, parameter :: count1 = 16385
    integer, parameter :: nFA = 8

    REAL (KIND=dp) :: g = 9.806, pi = 3.14159265D0
    
    REAL (KIND=dp) :: delta = 0.000000000001d0
    REAL (KIND=dp) :: h = 0.0005d0
    REAL (KIND=dp) :: alpha = 0.16666d0
    REAL (KIND=dp) :: beta = 0.25d0
    REAL (KIND=dp) :: gamma = 0.5d0
    REAL (KIND=dp) :: ONE = 1.0d0
    REAL (KIND=dp) :: ZERO = 0.0d0
    REAL (KIND=dp) :: MM = 1000.0d0
    
    
    contains
    
    
    
    end module Base

My first step in upgrading any program is to include BASE.f90, as above. @mecej4 showed me this trick.

Ularc is now an advanced Monte Carlo Analysis machine looking at statistical variations in the properties and loads.

If you use Fortran one does not run down the MATLAB rabbit hole.

Well done splice the main brace.

DavidWhite · ‎02-16-2026

Far better to set PI=ACOS(-1D0)

witwald · ‎02-16-2026

@DavidWhite. I'm curious as to why PI=ACOS(-1D0) is far better than PI=4.0*ATAN(1.0D0)? Can you elaborate? Cheers.

DavidWhite · ‎02-16-2026

I was commenting on the use of an approximate, hardcoded value.

Whether ACOS or ATAN is better/faster, I don't know. I prefer the ACOS version because it is a single call, and no extra calcs.

DavidWhite · ‎02-16-2026

I think you win on this one ... from Stackoverflow

JohnNichols · ‎02-17-2026

The earliest version of base.f90 that I can find is from 2016, at that time @mecej4 helped me mightily to get Harrison's 3D code running with a eigenvalue solver and inverter.

He showed me how to use base.f90. The then version has now been fixed for the new pi. At that time I had not known about the use of selected_real_kind and there is no way I could have solved the matrix input form on my own. So it was for line 3 that base was created as it saves a lot of copying of dp.

On Sunday I was running Strand7 which uses a Fortran motor and the 267 natural frequency set for a cable stayed bridge took 69 seconds to solve, in 1996 it took a week.

    Module Base

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

    INTEGER, PARAMETER :: sw = 2, srA = 1, st = 11, sa = 3, sm = 4, sc = 12, sRT = 25
    Integer, parameter :: slog = 9, sOUT = 10, sCAD = 21, sINI = 13, sWULF = 18
    Integer, parameter :: array_size = 300
    integer, parameter :: nt = 120
    integer, parameter :: mt = 200          !   Number of members    
    integer, parameter :: mn = 2    
    integer, parameter :: ml = 3
    integer, parameter :: mn_12 = 12
    integer, parameter :: mn_4 = 4
    integer, parameter :: count1 = 20000
    integer, parameter :: nFA = 8
    integer :: readFlag = 0

    REAL (KIND=dp) :: g = 9.806, pi = 4.0d0(arctan(1.0d0)
    
    REAL (KIND=dp) :: delta = 0.000000000001d0
    REAL (KIND=dp) :: h = 0.0005d0
    REAL (KIND=dp) :: alpha = 0.16666d0
    REAL (KIND=dp) :: beta = 0.25d0
    REAL (KIND=dp) :: gamma = 0.5d0
    REAL (KIND=dp) :: ONE = 1.0d0
    REAL (KIND=dp) :: ZERO = 0.0d0
    REAL (KIND=dp) :: DAMPRATIO = 100.0d0
    
    
    contains
    
    
    
    end module Base

. At that

andrew_4619 · ‎02-17-2026

Well that is a moot point as surely these should all be parameters so are evaluated at compile time! One would not want the value of pi to change!0 Also I think FPATAN was an x87 instruction so that is obsolete hardware. On X64 chips is seems acos and atan are not native instructions.

JohnNichols · ‎02-17-2026

Of course, thank you.

DavidWhite · ‎02-18-2026

I found somewhere that ACOS and ATAN are both determined from polynomials in X64. However, the ACOS needs more terms than ATAN, so that is preferrable.

However, if PI is stored as a PARAMETER, then this will be executed at compile time, so you'll only see the difference then, not at run time.

andrew_4619 · ‎02-15-2026

That looks much better all those numbered line hurt your eyes, good riddance! I am Curious why you chose to keep the fixed format, something I have not used for many years. It seems strange to refactor and not to use free form which in this case would take minutes! 80 character lines are far too limiting to make readable code particularly so when the first 6 are wasted!

JohnNichols · ‎02-15-2026

If you look at how 3 people played with the Fortran code in the Harrison program in the last month, one would see a lot of variation. I stole some of both other programmers code as they had useful features. Getting rid of arithmetic if's is like entering a nightmare.

The problem with lines longer than 80 you cannot see them in the editor, so even with long lines you continue.

andrew_4619 · ‎02-15-2026

"The problem with lines longer than 80 you cannot see them in the editor,"

I don't get that, what editor? Every editor I use is just fine. Lines up to 132 are just fine maybe not so much on 1980s screens.

JohnNichols · ‎02-16-2026

On the three editors I use, I can see 155 characters of text.

So whilst moving from .for to .f90 is nice, it depends on the available time.

So yes f90 is much better.

I grant your point.

witwald · ‎02-16-2026

@Andrew_44619, thanks for your interest in my efforts. You're right about all those numbered lines—they are very distracting noise when one is used to more modern versions of Fortran.

I concentrated on getting rid of those GOTOs and FORMATs and the (horrible) Holleriths, which was a large part of the work. By keeping most of the original style, it seemed to make comparisons with the original code a little bit easier for me when I was checking the refactoring steps. The next stage, which would be to first convert to Fortran 90 style free format lines, is a relatively easy one now for anyone who wants to go down that particular path.