ICE in GFN-FF implementation (#218)

- Build both GCC7 and GCC8 on Travis CI
- Workaround https://gcc.gnu.org/bugzilla/show_bug.cgi?id=84412

Co-authored-by: Sebastian Ehlert <28669218+awvwgk@users.noreply.github.com>

.travis.yml

@@ -8,18 +8,24 @@ addons:
       packages:
          - libblas-dev
          - liblapack-dev
+         - gfortran-7
          - gfortran-8
          - python-pip
 
-env: FC=gfortran-8 CC=gcc-8
-
 jobs:
    include:
       - name: meson-build
+        env: FC=gfortran-8 CC=gcc-8
+        script:
+        - meson build --buildtype release --prefix $PWD/install --warnlevel 0 -Dla_backend=netlib
+        - OMP_NUM_THREADS=2 ninja -C build test install
+      - name: meson-build
+        env: FC=gfortran-7 CC=gcc-7
         script:
         - meson build --buildtype release --prefix $PWD/install --warnlevel 0 -Dla_backend=netlib
         - OMP_NUM_THREADS=2 ninja -C build test install
       - name: CMake-build
+        env: FC=gfortran-8 CC=gcc-8
         script:
         - mkdir _build
         - cd _build

src/CMakeLists.txt

@@ -66,7 +66,7 @@ list(APPEND srcs
   "${dir}/eqrot.f90"
   "${dir}/esp.f"
   "${dir}/expire.f90"
-  "${dir}/filetools.f90"
+  "${dir}/filetools.F90"
   "${dir}/fixparam.f90"
   "${dir}/foden.f90"
   "${dir}/fragment.f90"

src/filetools.f90 → src/filetools.F90

@@ -15,6 +15,12 @@
 ! You should have received a copy of the GNU Lesser General Public License
 ! along with xtb.  If not, see <https://www.gnu.org/licenses/>.
 
+#ifdef __GNUC__
+#if GCC_VERSION < 70500
+#define HAS_GCC_BUG_84412
+#endif
+#endif
+
 subroutine print_filelist(iunit)
    use xtb_mctc_filetools
    use xtb_readin
@@ -211,8 +217,20 @@ subroutine close_file(unit)
    implicit none
    integer,intent(in) :: unit
    logical :: opened
+   integer :: i
    !$omp critical (io)
+   ! GCC7 cannot inquire on files with negative unit, fixed in GCC7.5.0
+   ! see https://gcc.gnu.org/bugzilla/show_bug.cgi?id=84412
+#ifdef HAS_GCC_BUG_84412
+   opened = .false.
+   do i = 1, nfiles
+      if (unit == filelist(i)%unit) then
+         opened = opened .or. filelist(i)%open
+      end if
+   end do
+#else
    inquire(unit=unit,opened=opened)
+#endif
    if (opened) then
       close(unit)
       call pop_file(unit)
@@ -225,8 +243,20 @@ subroutine remove_file(unit)
    implicit none
    integer,intent(in) :: unit
    logical :: opened
+   integer :: i
    !$omp critical (io)
+   ! GCC7 cannot inquire on files with negative unit, fixed in GCC7.5.0
+   ! see https://gcc.gnu.org/bugzilla/show_bug.cgi?id=84412
+#ifdef HAS_GCC_BUG_84412
+   opened = .false.
+   do i = 1, nfiles
+      if (unit == filelist(i)%unit) then
+         opened = opened .or. filelist(i)%open
+      end if
+   end do
+#else
    inquire(unit=unit,opened=opened)
+#endif
    if (opened) then
       close(unit,status='delete')
       call pop_file(unit,'d')

src/gfnff/gfnff_ini.f90

@@ -836,7 +836,7 @@ integer function itabrow6(i)
       pisea= 0
 
       if(pr) then
-         write(env%unit,'(10x,"iterarive Hueckel run to get P ...")')
+         write(env%unit,'(10x,"iterative Hueckel run to get P ...")')
       endif
       do pis=1,picount ! loop over pi systems
       npi   =0
@@ -992,7 +992,7 @@ integer function itabrow6(i)
                dum=dum+topo%qa(k)
             endif
          enddo
-         write(env%unit,'(5x,i2,10x,i4,10x,f8.3)')i,m,dum
+         write(env%unit,'(5x,i3,10x,i4,10x,f8.3)')i,m,dum
       enddo
       write(env%unit,*)
       endif

src/gfnff/gfnff_setup.f90

@@ -122,6 +122,7 @@ subroutine gfnff_input(env, mol, topo)
   logical           :: ex
   character(len=80) :: atmp
   character(len=80) :: s(10)
+  integer, allocatable :: rn(:)
 
   if (.not.allocated(topo%nb))       allocate( topo%nb(20,mol%n), source = 0 )
   if (.not.allocated(topo%qfrag))    allocate( topo%qfrag(mol%n), source = 0.0d0 )
@@ -134,19 +135,22 @@ subroutine gfnff_input(env, mol, topo)
   case(fileType%pdb)
     ini = .true.
     ifrag=0
-    associate(rn => mol%pdb%residue_number, qatom => mol%pdb%charge)
-      do iresidue = minval(rn),maxval(rn)
-        if (any(iresidue .eq. rn)) then
-          ifrag=ifrag+1
-          where(iresidue .eq. rn) topo%fraglist = ifrag
-        end if
-      end do
-      topo%nfrag = maxval(topo%fraglist)
-      do iatom=1,mol%n
-        topo%qfrag(topo%fraglist(iatom)) = topo%qfrag(topo%fraglist(iatom)) + dble(qatom(iatom))
-      end do
-      topo%qpdb = qatom
-    end associate
+    allocate(rn(mol%n))
+    rn(:) = mol%pdb%residue_number
+    do iresidue = minval(rn),maxval(rn)
+      if (any(iresidue .eq. rn)) then
+        ifrag=ifrag+1
+        where(iresidue .eq. rn) topo%fraglist = ifrag
+      end if
+    end do
+    deallocate(rn)
+    topo%nfrag = maxval(topo%fraglist)
+    if (.not.allocated(topo%qpdb)) allocate(topo%qpdb(mol%n))
+    do iatom=1,mol%n
+      topo%qfrag(topo%fraglist(iatom)) = topo%qfrag(topo%fraglist(iatom)) &
+        & + dble(mol%pdb(iatom)%charge)
+      topo%qpdb(iatom) = mol%pdb(iatom)%charge
+    end do
     ichrg=idint(sum(topo%qfrag(1:topo%nfrag)))
     write(env%unit,'(10x,"charge from pdb residues: ",i0)') ichrg
   !--------------------------------------------------------------------

src/meson.build

@@ -64,7 +64,7 @@ srcs += files(
   'eqrot.f90',
   'esp.f',
   'expire.f90',
-  'filetools.f90',
+  'filetools.F90',
   'fixparam.f90',
   'foden.f90',
   'fragment.f90',

src/type/CMakeLists.txt

@@ -36,7 +36,7 @@ list(APPEND srcs
   "${dir}/options.f90"
   "${dir}/param.f90"
   "${dir}/pcem.f90"
-  "${dir}/reader.f90"
+  "${dir}/reader.F90"
   "${dir}/setvar.f90"
   "${dir}/solvent.f90"
   "${dir}/timer.f90"

src/type/meson.build

@@ -34,7 +34,7 @@ srcs += files(
   'options.f90',
   'param.f90',
   'pcem.f90',
-  'reader.f90',
+  'reader.F90',
   'setvar.f90',
   'solvent.f90',
   'timer.f90',

src/type/reader.f90 → src/type/reader.F90

@@ -15,6 +15,12 @@
 ! You should have received a copy of the GNU Lesser General Public License
 ! along with xtb.  If not, see <https://www.gnu.org/licenses/>.
 
+#ifdef __GNUC__
+#if GCC_VERSION < 70500
+#define HAS_GCC_BUG_84412
+#endif
+#endif
+
 !> Wrapper to read from Fortran units
 module xtb_type_reader
    implicit none
@@ -70,12 +76,16 @@ subroutine openUnit(self, unit)
    !> IO unit connect to
    integer, intent(in) :: unit
 
+#ifdef HAS_GCC_BUG_84412
+   self%unit = unit
+#else
    logical :: opened
 
    inquire(opened=opened, unit=unit)
    if (opened) then
       self%unit = unit
    end if
+#endif
 
 end subroutine openUnit
 
@@ -155,12 +165,16 @@ subroutine closeUnit(self)
    !> Instance of reader
    class(TReader), intent(inout) :: self
 
+#ifdef HAS_GCC_BUG_84412
+   close(self%unit)
+#else
    logical :: opened
 
    inquire(opened=opened, unit=self%unit)
    if (opened) then
       close(self%unit)
    end if
+#endif
 
 end subroutine closeUnit