The NASA/NRL Horizontal Wind Model code HWM14 consists of two Fortran-90 files, two text data files and one unformatted Fortran data file. Unfortunately, the code contains several errors. Some subroutine arguments and some local variables are never used after declaration. Several two-dimensional arrays are used without initialization. Perhaps, zero-initialization was assumed. After spending a few hours on the code, I am puzzled as to how to reproduce the result files that are provided at the NRL site (I can match the output with some compilers, but not reliably so.) Any suggestions welcome.
2 Likes
subroutine alfbasis(nmax,mmax,theta,P,V,W) subroutine alfbasis(nmax,mmax,theta,P,V,W)
implicit none implicit none
integer(4), intent(in) :: nmax, mmax integer(4), intent(in) :: nmax, mmax
real(8), intent(in) :: theta real(8), intent(in) :: theta
real(8), intent(out) :: P(0:nmax,0:mmax) real(8), intent(out) :: P(0:nmax,0:mmax)
real(8), intent(out) :: V(0:nmax,0:mmax) real(8), intent(out) :: V(0:nmax,0:mmax)
real(8), intent(out) :: W(0:nmax,0:mmax) real(8), intent(out) :: W(0:nmax,0:mmax)
integer(8) :: n, m integer(8) :: n, m
real(8) :: x, y real(8) :: x, y
real(8), parameter :: p00 = 0.70710678118654746d real(8), parameter :: p00 = 0.70710678118654746d
> real(8), parameter :: zero = 0.d0
P(0,0) = p00 P(0,0) = p00
> do n=0,min(nmax,mmax-1)
> W(n,n+1) = zero
> end do
> W(0:nmax,0) = zero
> V(0,0) = zero
x = dcos(theta) x = dcos(theta)
y = dsin(theta) y = dsin(theta)
do m = 1, mmax do m = 1, mmax
W(m,m) = cm(m) * P(m-1,m-1) W(m,m) = cm(m) * P(m-1,m-1)
P(m,m) = y * en(m) * W(m,m) P(m,m) = y * en(m) * W(m,m)
do n = m+1, nmax do n = m+1, nmax
W(n,m) = anm(n,m) * x * W(n-1,m) - bnm(n,m) * W(n,m) = anm(n,m) * x * W(n-1,m) - bnm(n,m) *
P(n,m) = y * en(n) * W(n,m) P(n,m) = y * en(n) * W(n,m)
V(n,m) = narr(n) * x * W(n,m) - dnm(n,m) * W( V(n,m) = narr(n) * x * W(n,m) - dnm(n,m) * W(
W(n-2,m) = marr(m) * W(n-2,m) W(n-2,m) = marr(m) * W(n-2,m)
It seems to be that INTENT(OUT) arrays are not defined in all cases. The above extra assignments seem to do the trick.
1 Like
Thanks, @themos. Sorry, I have trouble copying and pasting your two-column code sections into a compilable Fortran source file. For example, your last line of code has
W(n-2,m) = marr(m) * W(n-2,m)
repeated on one line, with blanks in between. Am I to insert a semicolon between the two copies? Or fold the two pieces into two successive lines? I should be very grateful if you can just attach the source file to your reply.
Thank you!
I am also a bit worried that âThe above extra assignmentsâ may do more than just âdo the trickâ.
1 Like
subroutine alfbasis(nmax,mmax,theta,P,V,W)
implicit none
integer(4), intent(in) :: nmax, mmax
real(8), intent(in) :: theta
real(8), intent(out) :: P(0:nmax,0:mmax)
real(8), intent(out) :: V(0:nmax,0:mmax)
real(8), intent(out) :: W(0:nmax,0:mmax)
integer(8) :: n, m
real(8) :: x, y
real(8), parameter :: p00 = 0.70710678118654746d0
real(8), parameter :: zero = 0.d0
P(0,0) = p00
do n=0,min(nmax,mmax-1)
W(n,n+1) = zero
end do
W(0:nmax,0) = zero
V(0,0) = zero
x = dcos(theta)
y = dsin(theta)
do m = 1, mmax
W(m,m) = cm(m) * P(m-1,m-1)
P(m,m) = y * en(m) * W(m,m)
do n = m+1, nmax
W(n,m) = anm(n,m) * x * W(n-1,m) - bnm(n,m) * W(n-2,m)
P(n,m) = y * en(n) * W(n,m)
V(n,m) = narr(n) * x * W(n,m) - dnm(n,m) * W(n-1,m)
W(n-2,m) = marr(m) * W(n-2,m)
enddo
W(nmax-1,m) = marr(m) * W(nmax-1,m)
W(nmax,m) = marr(m) * W(nmax,m)
V(m,m) = x * W(m,m)
enddo
P(1,0) = anm(1,0) * x * P(0,0)
V(1,0) = -P(1,1)
do n = 2, nmax
P(n,0) = anm(n,0) * x * P(n-1,0) - bnm(n,0) * P(n-2,0)
V(n,0) = -P(n,1)
enddo
return
end subroutine alfbasis
1 Like
The only extra work is assigning parts of W and V to zero.
1 Like
Thanks, @Themos.
I substituted your code for Subroutine Alfbasis in file HWM14.f90. The program ran without any errors, but the output did not agree with the check outputs posted at the HWM14 site. At this point I think that I may have to ask the authors of the code to help me understand what the code may be doing. Earlier, I had found that parts of arrays gvbar( : , : ) and dwbar( : , : ) were being used before definition.
By trial and error, I found that setting dwbar(1:nmax,0)=0 in dwm07b, setting gvbar(0,0)=0 in gd2qd, and setting dwbar(1:nmax,0)=0 after allocating dwbar in subroutine initdwm gave runs without error from uninitialized variables and also gave output that agreed with the reference outputs.
Another puzzling error: Subroutine DWM07 contains this strange source line:
real(4), parameter :: talt=125.0 !, twidth=5.0
Note the â!â, which makes half the line a comment, leaves the local variable âtwidthâ undeclared and uninitialized. A few lines later in the same subroutine, this variable is used:
dw = dw / (1 + exp(-(alt - talt)/twidth))
Try running the code under valgrind to find more cases of uninitialized variables.
Thanks for the suggestions. After initializing dwbar(1:nmax,0)=0, and gvbar(0,0)=0, I eliminated some
subroutine arguments such as STL,F107A and F107 from a number of subroutine calls because the subroutine ignored the values passed. The program now runs and delivers results that agree with the provided outputs.
2 Likes
i do not like these definitions, as I would prefer
real(8), parameter :: zero = 0
real(8), parameter :: half = 0.5
real(8), parameter :: p00 = sqrt (half)
There is no precision difference between 0 and 0.d0 for this usage and sqrt (half) is a more stable presentation of precision.
It is surprising (as the files are Fortran 90 vintage) that
a) uninitialised variable have peristed, as zero-initialization of variables is now not common, and
b) an unformatted Fortran data file is not portable. A text file would be more informative.
Most legacy codes have passed through a number of different compilers with different limitations, which should have eliminated these features.
4 Likes
A feature I dislike is declarations beginning real(8) that mean different things to different compilers. Last time I used nagfor one would need to set an option to allow it. That is presumably why we may write real(dp) with dp as an integer constant = kind(1d0) or selected_real_kind(15) or âŚ
4 Likes
The NAG option that you mention is â-kind=byteâ, and I have found myself surprised to relearn that âinteger(4)â and âreal(8)â are not portable across compilers.
1 Like
I have always used this module since you showed me how + unitialized variables are a pain in old code, merry xmas
Module Base
INTEGER, PARAMETER :: dp = selected_real_kind(15, 307)
INTEGER, PARAMETER :: sw = 2 ! Output file
INTEGER, PARAMETER :: srA = 15 ! output.txt file
INTEGER, PARAMETER :: srB = 16 ! output.txt file
INTEGER, PARAMETER :: st = 14
INTEGER, PARAMETER :: sCAD = 12
INTEGER, PARAMETER :: sa = 3 ! Output file
INTEGER, PARAMETER :: smWrite = 4
INTEGER, PARAMETER :: si = 1
Integer, parameter :: slog = 9 ! Log file
Integer, parameter :: nta = 200 ! Log file
Integer, parameter :: outNode = 63 ! Log file
Integer, parameter :: inNode = 0 ! Log file
integer, parameter :: nt1 = 2000
integer, parameter :: mt1 = 2000 ! Number of members
integer, parameter :: mn1 = 2
integer, parameter :: ml1 = 3
integer, parameter :: ml2 = 4
integer, parameter :: ml30 = 9000
integer, parameter :: limit = 9000
REAL (KIND=dp), PARAMETER :: gr = 9.806_DP, pi = 3.14159265_DP ! Standard parameters
REAL (KIND=dp), PARAMETER :: delta = 0.001_DP ! Error number of checking for zero
REAL (KIND=dp), PARAMETER :: deltafreq = 0.00001_DP ! 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
REAL (KIND=dp), PARAMETER :: THREE = 3.0_DP
REAL (KIND=dp), PARAMETER :: THIRTY = 30.0_DP
REAL (KIND=dp), PARAMETER :: FOUR = 4.0_DP
@mecej4 you should post the working code somewhere, if you can, so that itâs easier for others to run it.
1 Like
OndĹej, Of course, you are correct. Here is a link to a file sharing site where anyone with the link can access the file .
The forum does not allow a file with the â.zipâ extension to be uploaded. What is the approved way of providing a .zip or .tgz file? I could, obviously, rename the zip file with a .jpg extension, but that might upset someone because of the deception.
The zip contains the two .f90 source files, two text data files and a .bin unformatted data file (compatible with Gfortran, NAG Fortran, Intel and Absoft)
I suspect that there are more instances of uninitialized variables in the code. Without some users of the code participating in the process, I am afraid that I am at the point of fixing things that âainât brokenâ.
Nevertheless, I am surprised that despite the several instances of uninitialized variables in it, the HWM14 code outputs results that are not only reasonable but quite accurate. A comparison of the âGfortran resultsâ placed over a decade ago on the NRL site and results from the NRL source files compiled by me using Gfortran 13.4 today show just 12 differences among 900 numbers. Two examples of these differences: 26.175 vs 26.176; 155.577 vs 155.576. That there are so few differences indicates that the consequence of having uninitialized variables is not catastrophe, as one would imagine, and that assumed zero-initialization, as Fortran programmers expected and compilers delivered for decades, was not often wrong. Even then, getting good results from a faulty program is disconcerting.
Perfect, thanks. The standard way to share code is to post it on GitHub or similar sites (GitLab, Codeberg, etc.), there are many options.
1 Like
Here is an update. After some stumbling, I was able to download the Fortran source files and data files for HWM07, the older version of HWM14. Somewhat unbelievably, the output results from HWM14 and HWM07 are nearly identical. The older Fortran source code, too, uses uninitialized variables, but their presence does not affect the output unless one insists on checking for them! Here is a link to the files for HWM07
2 Likes
I think that I now know where the undefined variables originate. In file HWM14.f90, Subroutine DWM07 is declared with seven arguments, the last of which is
`real(4),intent(out) :: DW(2)â
I think that the intent should be (in out). According to the rules of Fortran, an argument with intent(out) technically becomes âundefinedâ at subroutine entry, although the compiler need not do anything to make that happen. Not many compilers can detect this kind of error. As a result, the error can persist for decades, as we have learned from HWM07 and HWM14.
FYI: hereâs an archive of this model on GitHub: GitHub - jacobwilliams/HWM14: Unofficial mirror of Horizontal Wind Model 2014
I donât remember if I made any changes or not. But we could work on fixing the bugs there if we want.
2 Likes
Thanks for responding, and I am keen to do what I can to fix the error. It is nice that you have an unofficial mirror set up, since I ran into several problems accessing the files at the official site . Attempting to download one of the data files using the wget utility results in he message âforbiddenâ â strange for a public access govt. web site.
Just now, I checked the file src/hwm14.f90 in your mirror site, and it has the incorrect intent(out) on line 752, in subroutine DWM07.
Likewise, the arguments P,V,W in Subroutine Alfbasis in your alf.f90 should be declared intent(in out), rather than (out), as they are now.
Assuming you can agree with my analysis of the error and the proposed correction, do you know how to get the correction propagated to the NRL site?
It was a post made in 2012 by another user of the older HWM93 on the Silverfrost Fortran user forum that first brought this issue to my attention.
Added on Jan 11 2026: As I noted earlier, some of the intent declarations in the source files on your Github for HWM14 have incorrect intent declarations. Today, I took out all the intent declarations, recompiled, linked and ran the test problem again. The output results to STDOUT were byte-for-byte identical with the results that I recei ved using the code with the intents left âas-wasâ. One implication is that the code generated by the compiler is not affected by what the intent declarations in the source codes say.
@JohnCampbell, Here is evidence from one of the NRL source files to support your comments regarding named constants. The file hwm14.f90 contains the following lines
266: real(8), parameter :: pi=3.1415926535897932
528: real(8),parameter :: twoPi = 2.0d0*3.1415926535897932384626433832795d0
1158: real(8), parameter :: pi = 3.1415926535897932d0
1406: real(4), parameter :: pi=3.141592653590```
Not only are there three values of Pi given, but they differ in the number of significant digits as well as kind.