[Wien] Compilation issues for lapw0 v.13 on Mac OS X (request for code review)

[Laurence Marks]

As mentioned in the other responses, so long as an array is large enough
there is no problem with having it dimensioned smaller than this in fortran
(with some caveats about passing values via subroutines).

The even case does for certain have a bug, and the relevant lines in
dergl_spline should be

        if(.not.even)then
                do J=1,n
                        rw(-J)=-r(j)
                        w1(-J)=-f(J)
                enddo
                rw(0)=0.D0
                w1(0)=0.D0
        else
                do J=1,n !-1   Bug here found by Kevin Jorissen
                        rw(-J+1)=-r(j)
                        w1(-J+1)= f(J)
                enddo
        endif

I think Peter may have already found this one and it is cured in the next
release, but I am not sure.



On Fri, Aug 15, 2014 at 5:20 PM, Kevin Jorissen <kevinjorissenpdx at gmail.com>
wrote:

>  Dear WIEN2k community,
>
>
>  I have identified some issues in v.13 of the SRC_lapw0 code that
> prevented me from compiling a working Mac OS X version.  I am also
> submitting a solution.  The details are below my signature.  After the code
> edits, I can compile and run lapw0 correctly using the recommended
> compilation settings (no additional flags like -O0 or -fp-model strict).
>
>  The issue is somewhat sloppy code, legal in Fortran but vulnerable to
> compiler trip-ups, in the "drho" section of lapw0.  It is quite possible
> that these issues do not occur on Linux, as memory is managed differently
> on the different OS's (e.g. which variables go on stack or heap depending
> on their allocation).  Nevertheless it seems useful to report.  They occur
> for me with both gfortran and ifort, indicating it's not just Intel
> goofiness.  I would be grateful for feedback from the developers or the
> community.  I'm happy to send source files - didn't include them now as I'm
> not sure the mailing list accepts attachments.
>
>  Finally, are total energies expected to be identical between wien2k_12
> and wien2k_13?  Or are small differences ~10^-4 or ~10^-5 normal?
>
>
>  Cheers,
>
>  Kevin Jorissen
>
>
>
>  Details:
> OS X 10.9.4
> ifort 14.0.1 or ifort 14.0.2 or gfortran 4.9.0
>
>  Test:
> * Standard "run_lapw" for HOPG graphite.
> * The test is passed if the calculation:  a/ does not crash ; b/ converges
> to approximately the same total energy (grep :ENE case.scf) as a wien2k-12
> calculation (compiled with the same compiler as wien2k-13 and with "safe"
> options (no optimization)).  Initialization is done with that same
> wien2k-12 version, i.e. I test only the "SCF" programs.
>
>  Outcome:
> * The test did not pass.  The wien2k-13 calculation fails almost always: A
> seclr4 error in lapw1 following NaN in the case.scf0 and case.output0, so
> the problem is in lapw0.  Case.scf0 points to the xcpot1 section.
> wien2k.at reports updates to the "drho" in v.13.  The crash occurs in
> varying SCF iteration for each run, indicating an "instability", i.e.
> numerical noise or memory issues.
> * Trying different configurations didn't help.  E.g. "safe" compilation
> options (-O0 -fp-model strict ...) ; e.g. swapping ifort for gfortran ;
> trying different fftw libraries ...
>
>  Solution (after tedious debugging ...):
> Several small code changes in drho.f and dergl.f to correct array size
> issues.
>
>  DRHO.F
> * The array "r" passed to dergl is size "nrad" (allocated elsewhere) while
> dergl expects an array of size "n=jri(jatom)".  E.g. 881 vs. 781 for my
> test.
> Solution (it could be done just by restricting call(...r(1:n)...) also):
>
>       real*8 ::  rlocal(n)
>
>        rlocal(:)=0.d0
>
>       do i=1,min(size(r),n)
>
>          rlocal(i) = r(i)
>
>       enddo
>
> !        call dergl(n,r,c,g,even,g0)
>
> !        call dergl(n,r,g,g2,.not.even,g0)
>
>         call dergl(n,rlocal,c,g,even,g0)
>
>         call dergl(n,rlocal,g,g2,.not.even,g0)
>
>  DERGL.F
> * In calls to Fornberg, the array "A(0:Ninp)" is made with a larger array
> "RW(I)".  Solution: replace with appropriate subarray of size Ninp+1:
>
> !                call  Fornberg(Ninp,Minp,RW(J),RW(J-N1),D)
>
>                 call  Fornberg(Ninp,Minp,RW(J),RW(J-N1:J-N1+Ninp),D)
>
>  !                call  Fornberg(Ninp,Minp,R(J),R(JL),D)
>
>                 call  Fornberg(Ninp,Minp,R(J),R(JL:JL+Ninp),D)
>
>
>   !      call  Fornberg(Ninp,Minp,0.0D0,RW(-N1),D)
>
>       call  Fornberg(Ninp,Minp,0.0D0,RW(-N1:-N1+Ninp),D)
>
>
>  * In subroutine dergl_spline, rw is not initialized correctly for the
> "even" case: the lower point rw(-n+1) is uninitialized and could be
> anything.  Solution:
>
>        else
>
>                 do J=1,n-1
>
>                         rw(-J+1)=-r(j)
>
>                         w1(-J+1)= f(J)
>
>                 enddo
>
>                 rw(-n+1)=-r(n) ! initialize this field!
>
>                 w1(-n+1)=f(n)  ! Technically this is redundant (see w1=0
> above).
>
>         endif
>
>
>  PS xcpot1.f has some unrelated issue with sigmamaxup and irsigmaup,
> which gets used without being initialized or something like that.  I don't
> remember the details now, and it doesn't seem to hinder the program, but it
> popped up during a runtime check at some point - I believe it was triggered
> by a WRITE statement containing r(irsigmaup) when irsigmaup can be
> uninitialized.
>
>
>
>  TEST SCRIPT:
>
>  ~/science/graphite13-opt-test% more ../MyWienTest
>
> #!/bin/tcsh -f
>
>
>  echo Running test in `pwd`
>
> echo WIENROOT is $WIENROOT
>
> rm -f *scf* *clm* *broy* *tmp* *vec* *def *error :log
>
> x dstart -d
>
> ~/science/wien2k-12-noopt/dstart dstart.def
>
> run -ec 0.000001
>
>
>  exit 0
>
>
>
>


-- 
Professor Laurence Marks
Department of Materials Science and Engineering
Northwestern University
www.numis.northwestern.edu
Corrosion in 4D: www.numis.northwestern.edu/MURI
Co-Editor, Acta Cryst A
"Research is to see what everybody else has seen, and to think what nobody
else has thought"
Albert Szent-Gyorgi
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://zeus.theochem.tuwien.ac.at/pipermail/wien/attachments/20140817/a56c6801/attachment.htm>