A few years ago mecej4 on the MKL Forum helped me to fix - actually he fixed it and I just watched in awe, a bit like Shock and Awe in Fortran a program for the analysis of water supply systems. In mecej4's fix he used interfaces.
So in thinking about the Fortran Manuals - I dragged out Fortan 95/2003 explained. The Metcalf book is more confusing than the Sonic Menu for hamburgers, but I digress. I could not make head nor tail of interfaces in the book.
In the oldest iteration of the fixed program, we have a module called
Network that contains a lot of variable declarations, it ends then there are a series of subroutines, which are in the same file, but not in the module.
In the main program unit there are :
program Magni use Network implicit none interface subroutine PumpHead( H1, H4, a0, a1, a2, a3, HRatio, retQ, flow0, pumpQ) use Network real(kind=dp) H1 real(kind=dp) H4 real(kind=dp) a0 real(kind=dp) a1 real(kind=dp) a2 real(kind=dp) a3 real(kind=dp) HRatio real(kind=dp) retQ real(kind=dp) flow0 real(kind=dp) pumpQ end subroutine PumpHead end interface
a series of interface calls -- this program still compiles and run fine on the latest
A later version of the program has the subroutines within the module
MODULE network INTEGER, PARAMETER :: dp = selected_real_kind(15, 307) REAL (KIND=dp) :: g = 9.806, pi = 3.14159265D0 TYPE snode INTEGER id INTEGER :: tank = -1 INTEGER :: zone = 0 REAL (KIND=dp) x REAL (KIND=dp) y REAL (KIND=dp) e REAL (KIND=dp) c0 REAL (KIND=dp) qbase REAL (KIND=dp) cs REAL (KIND=dp) h REAL (KIND=dp) c REAL (KIND=dp) inddemand REAL (KIND=dp) parkdemand REAL (KIND=dp) leakdemand REAL (KIND=dp) resdemand REAL (KIND=dp) dumc REAL (KIND=dp) q END TYPE TYPE spipe INTEGER id INTEGER n1 INTEGER n2 INTEGER :: zone = 0 INTEGER pumpstatus INTEGER nodei, nodej REAL (KIND=dp) :: a0 = 0.0 REAL (KIND=dp) a1 REAL (KIND=dp) a2 REAL (KIND=dp) a3 REAL (KIND=dp) length REAL (KIND=dp) diameter REAL (KIND=dp) c REAL (KIND=dp) wlrf !pipe wall roughness, m REAL (KIND=dp) e REAL (KIND=dp) pumpq REAL (KIND=dp) deltaq REAL (KIND=dp) q REAL (KIND=dp) h0 REAL (KIND=dp) h1 REAL (KIND=dp) h2 REAL (KIND=dp) h3 REAL (KIND=dp) h4 REAL (KIND=dp) qt REAL (KIND=dp) v REAL (KIND=dp) :: k = 0.0 REAL (KIND=dp) :: k1 REAL (KIND=dp) :: RE REAL (KIND=dp) :: fr integer :: cv = 0 integer :: prv = 0 REAL (KIND=dp) :: prvcoeffic = 1.0 REAL (KIND=dp) :: prvzeropressure REAL (KIND=dp) kvis !kinematic viscosity, m^2/s REAL (KIND=dp) tc !water temperature, celsius END TYPE TYPE spump INTEGER :: pump = -1 INTEGER :: pipe = 0 INTEGER :: ontime = 0 INTEGER :: offtime = 0 INTEGER :: tank = -1 REAL (KIND=dp) :: onlevel = 0.0 REAL (KIND=dp) :: offlevel = 0.0 END TYPE TYPE stank INTEGER id REAL (KIND=dp) area REAL (KIND=dp) c REAL (KIND=dp) h0 REAL (KIND=dp) hmax REAL (KIND=dp) hmin REAL (KIND=dp) hnext REAL (KIND=dp) k1 INTEGER node REAL (KIND=dp) v REAL (KIND=dp) centreX REAL (KIND=dp) centerY REAL (KIND=dp) radius END TYPE TYPE constants REAL (KIND=dp) :: resdemfact = 0.15 REAL (KIND=dp) :: inddemfact = 0.15 REAL (KIND=dp) :: parkfactor = 540.0 REAL (KIND=dp) :: demand = 0.0 REAL (KIND=dp) :: flow0 = 0.01 REAL (KIND=dp) :: gpmpercfs = 448.4 REAL (KIND=dp) :: htol = 0.01 REAL (KIND=dp) :: huge_int = 999999999 REAL (KIND=dp) :: lperet = 0.15 REAL (KIND=dp) :: lperft3 = 28.3 REAL (KIND=dp) :: metricratio = 1000.0 REAL (KIND=dp) :: limit = 1000000 INTEGER :: maxhrs = 24 INTEGER :: maxid INTEGER :: maxiter = 2000 INTEGER :: maxnodes INTEGER :: maxpipes INTEGER :: maxpumps INTEGER :: maxsegs INTEGER :: maxtanks INTEGER :: maxtexts = 50 INTEGER :: nodedataccount = 12 INTEGER :: runtime = 1 REAL (KIND=dp) :: qtol = 0.01 REAL (KIND=dp) :: ttol = 5.0 REAL (KIND=dp) :: unknown = -999999 CHARACTER *2 :: units = 'SI' END TYPE !----------------------------------------------------- ! ! ! !----------------------------------------------------- CONTAINS SUBROUTINE timingline(sw) IMPLICIT NONE INTEGER :: sw INTEGER :: tmpday, tmpmonth, tmpyear INTEGER :: tmphour, tmpminute, tmpsecond, tmphund, values(8) CHARACTER *2 :: mer CALL date_and_time(values=values) tmpyear = values(1) tmpmonth = values(2) tmpday = values(3) tmphour = values(5) tmpminute = values(6) tmpsecond = values(7) tmphund = values(8) IF (tmphour.GT.12) THEN mer = 'pm' tmphour = tmphour - 12 ELSE mer = 'am' END IF WRITE (sw, 100) tmpday, tmpmonth, tmpyear WRITE (sw, 110) tmphour, tmpminute, tmpsecond, tmphund, mer WRITE (*, 100) tmpday, tmpmonth, tmpyear WRITE (*, 110) tmphour, tmpminute, tmpsecond, tmphund, mer WRITE (sw, 120) RETURN 100 FORMAT (' Analysis Date : ', I2, '/', I2.2, '/', I4.4) 110 FORMAT (' Time : ', I2, ':', I2.2, ':', I2.2, ':', I3.3, & ' ', A) 120 FORMAT (' ') END SUBROUTINE !------------------------------------------------------------------- !
And there are no inferfaces:
I was wondering which is the preferred standard and really if the module contains the subroutines, is there any need for interfaces?
I talked about this in Doctor Fortran Gets Explicit - Again! In general, putting Fortran procedures in modules (or making them internal procedures) is preferred over writing interface blocks. INTERFACE is still needed for generic declarations, and it's a way of declaring user-defined derived type I/O, but the actual routines themselves should be in modules. If you're writing interface blocks for Fortran routines, you're doing it wrong. (A possible exception is if you're "updating" an old program and want a quick-and-dirty way of creating explicit interfaces for external procedures.)
Thanks for the note. It prompted me to go back and look at the older code that I have not touched since 2014.
mecej4 at the time showed me how to use PARDISO and also helped a lot to fix the water supply analysis program - which is based on the Gradient Algorithm - I had an early 1987 program based on the Hardy Cross method that had convergence problems with 2000 node water supply models. mecej4 suggested Fortran as being faster than C# - as always he was correct as the graph shows. I am sure he showed me how to use modules properly as well.
For anyone who is not using PARDISO to invert, you are really wasting a lot computer time. - Note that it is a log scale for time not linear.
Teh water supply program final version uses your suggested subroutines in the modules. it is quite fun to look at the original PASCAL and C code from the EPANET - and then the version rewritten in C# (ok some of us have no choice sometimes) and then the Fortran version. As you can see the Fortran can handle bigger problems faster. Those are dense arrays as well not sparse - fully populated with random numbers.
A real pity I got sidetracked onto other problems, this is still the fastest water supply analysis program that I am aware of, and I am happy to give it away.
Dear Mr John Nichols
I am also in the field of Water Supply network analysis and I am happy to see that you have done lot of work in this field.
I am having a network of pipes in AutoCAD from which I am able to make a matrix in AutoLisp which I want to solve using Gradient Method. Can you please explain how can I get matrix inverse in Autolisp ?
Thanks a lot
Interesting problem, I have not written an inverter in LISP, should not be that hard. Take any FORTRAN inversion routine -- there are a lot and recode it in LISP -- it will be very slow as it is not complied.
My data is stored in AutoCAD and I use LISP routines to export the data, the issue is having the node numbers in a set location for each node and being able to recover them accurately. Pipes also need to be connected, which in AutoCAD is not hard, but you have to make sure it happens.
If I can write an AI program in AutoLISP then you can do inversion.
If you want the LISP code and FORTRAN code for what I have send me a note.