Add datacomp to subroutine iabort

src/dft.F90

@@ -16,7 +16,7 @@
  subroutine formrfockexcor(fockdsum,fockd,energy,totalelec,cmo,atomvec,radpt,angpt, &
 & rad,ptweight,vao,vmo,xyzpt,rsqrd,transcmo,work,ndftatom, &
 & nproc,myrank,mpi_comm, &
-& datajob,datamol,databasis)
+& datajob,datamol,databasis,datacomp)
 !---------------------------------------------------------------------------------------
 !
 ! Driver of DFT Fock matrix formation from exchange-correlation functionals
@@ -33,11 +33,12 @@ subroutine formrfockexcor(fockdsum,fockd,energy,totalelec,cmo,atomvec,radpt,angp
 ! totalelec(Number of numerially integrated electrons)
 ! vao,vmo,xyzpt,work (work space)
 !
- use modtype, only : typejob, typemol, typebasis
+ use modtype, only : typejob, typemol, typebasis, typecomp
  implicit none
  type(typejob),intent(in) :: datajob
  type(typemol),intent(in) :: datamol
  type(typebasis),intent(in) :: databasis
+ type(typecomp),intent(in) :: datacomp
  integer,intent(in) :: ndftatom, nproc, myrank, mpi_comm
  integer :: ngridatom, iatom, irad, ileb, icount, ilebstart, jatom, imo
  real(8),parameter :: zero=0.0D+00, one=1.0D+00, two=2.0D+00
@@ -88,7 +89,7 @@ subroutine formrfockexcor(fockdsum,fockd,energy,totalelec,cmo,atomvec,radpt,angp
  enddo
 !
  call gridraomo(vao,vmo,transcmo,xyzpt,rsqrd,aocutoff,datajob%threshex,ndftatom, &
-& datamol,databasis)
+& datamol,databasis,datacomp)
 !
  rhoa= zero
  grhoa(1:3)= zero
@@ -214,7 +215,7 @@ subroutine formufockexcor(fockd1,fockd2,fockd3,energy,totalelec,cmoa,cmob,atomve
 & radpt,angpt,rad,ptweight,vao,vmoa,vmob,xyzpt,rsqrd, &
 & transcmoa,transcmob,work,ndftatom, &
 & nproc,myrank,mpi_comm, &
-& datajob,datamol,databasis)
+& datajob,datamol,databasis,datacomp)
 !---------------------------------------------------------------------------------------
 !
 ! Driver of unrestricted DFT Fock matrix formation from exchange-correlation functionals
@@ -225,11 +226,12 @@ subroutine formufockexcor(fockd1,fockd2,fockd3,energy,totalelec,cmoa,cmob,atomve
 ! fock2 (Beta Fock matrix)
 ! fock3 (Work space)
 !
- use modtype, only : typejob, typemol, typebasis
+ use modtype, only : typejob, typemol, typebasis, typecomp
  implicit none
  type(typejob),intent(in) :: datajob
  type(typemol),intent(in) :: datamol
  type(typebasis),intent(in) :: databasis
+ type(typecomp),intent(in) :: datacomp
  integer,intent(in) :: ndftatom, nproc, myrank, mpi_comm
  integer :: ngridatom, iatom, irad, ileb, icount, ilebstart, jatom, imo
  real(8),parameter :: zero=0.0D+00, one=1.0D+00, two=2.0D+00
@@ -284,7 +286,7 @@ subroutine formufockexcor(fockd1,fockd2,fockd3,energy,totalelec,cmoa,cmob,atomve
  enddo
 !
  call griduaomo(vao,vmoa,vmob,transcmoa,transcmob,xyzpt,rsqrd,aocutoff, &
-& datajob%threshex,ndftatom,datamol,databasis)
+& datajob%threshex,ndftatom,datamol,databasis,datacomp)
 !
  rhoa= zero
  grhoa(1:3)= zero
@@ -490,17 +492,19 @@ subroutine calcradpt(radpt,nrad)
 end
 
 
-!-----------------------------------
- subroutine calclebpt(angpt,nleb)
-!-----------------------------------
+!--------------------------------------------
+ subroutine calclebpt(angpt,nleb,datacomp)
+!--------------------------------------------
 !
 ! Calculate Levedev quadrature points and weights
 !
 ! In : nleb (number of Lebedev quadrature points)
 ! Out : angpt(1:3,*) (Levedev points)
 ! angpt(4 ,*) (Levedev weights)
 !
+ use modtype, only : typecomp
  implicit none
+ type(typecomp),intent(inout) :: datacomp
  integer,intent(in) :: nleb
  integer :: ileb
  real(8),parameter :: pi4=1.256637061435917D+01
@@ -555,7 +559,7 @@ subroutine calclebpt(angpt,nleb)
  call lebedev1454(angpt)
  case default
  write(*,'(" Error! Nleb=",i4," is not supported. ")')nleb
- call iabort
+ call iabort(datacomp)
  end select
 !
  do ileb= 1,nleb
@@ -648,15 +652,16 @@ subroutine calcgridweight(ptweight,rad,radpt,angpt,atomvec,surface,xyzpt,work,nr
 
 !----------------------------------------------------------------------------------
  subroutine gridraomo(vao,vmo,transcmo,xyzpt,rsqrd,aocutoff,threshex,ndftatom, &
-& datamol,databasis)
+& datamol,databasis,datacomp)
 !----------------------------------------------------------------------------------
 !
 ! Calculate closed-shell AO and MO values for a grid point
 !
- use modtype, only : typemol, typebasis
+ use modtype, only : typemol, typebasis, typecomp
  implicit none
  type(typemol),intent(in) :: datamol
  type(typebasis),intent(in) :: databasis
+ type(typecomp),intent(in) :: datacomp
  integer,intent(in) :: ndftatom
  integer :: ii, jj
  real(8),parameter :: zero=0.0D+00
@@ -667,7 +672,7 @@ subroutine gridraomo(vao,vmo,transcmo,xyzpt,rsqrd,aocutoff,threshex,ndftatom, &
 ! Calculate AO values for a grid point
 !
  vao(:,:)= zero
- call gridao(vao,xyzpt,rsqrd,threshex,ndftatom,databasis)
+ call gridao(vao,xyzpt,rsqrd,threshex,ndftatom,databasis,datacomp)
 !
 ! Calculate MO values for a grid point
 !
@@ -690,15 +695,16 @@ subroutine gridraomo(vao,vmo,transcmo,xyzpt,rsqrd,aocutoff,threshex,ndftatom, &
 
 !---------------------------------------------------------------------------------
  subroutine griduaomo(vao,vmoa,vmob,transcmoa,transcmob,xyzpt,rsqrd,aocutoff, &
-& threshex,ndftatom,datamol,databasis)
+& threshex,ndftatom,datamol,databasis,datacomp)
 !---------------------------------------------------------------------------------
 !
 ! Calculate open-shell AO and MO values for a grid point
 !
- use modtype, only : typemol, typebasis
+ use modtype, only : typemol, typebasis, typecomp
  implicit none
  type(typemol),intent(in) :: datamol
  type(typebasis),intent(in) :: databasis
+ type(typecomp),intent(in) :: datacomp
  integer,intent(in) :: ndftatom
  integer :: ii, jj
  real(8),parameter :: zero=0.0D+00
@@ -709,7 +715,7 @@ subroutine griduaomo(vao,vmoa,vmob,transcmoa,transcmob,xyzpt,rsqrd,aocutoff, &
 !
  vao(:,:)= zero
 !
- call gridao(vao,xyzpt,rsqrd,threshex,ndftatom,databasis)
+ call gridao(vao,xyzpt,rsqrd,threshex,ndftatom,databasis,datacomp)
 !
  vmoa(:,:)= zero
  vmob(:,:)= zero
@@ -735,15 +741,16 @@ subroutine griduaomo(vao,vmoa,vmob,transcmoa,transcmob,xyzpt,rsqrd,aocutoff, &
 end
 
 
-!--------------------------------------------------------------
- subroutine gridao(vao,xyzpt,rsqrd,threshex,ndftatom,databasis)
-!--------------------------------------------------------------
+!--------------------------------------------------------------------------
+ subroutine gridao(vao,xyzpt,rsqrd,threshex,ndftatom,databasis,datacomp)
+!--------------------------------------------------------------------------
 !
 ! Calculate AO values for a grid point
 !
- use modtype, only : typebasis
+ use modtype, only : typebasis, typecomp
  implicit none
  type(typebasis),intent(in) :: databasis
+ type(typecomp),intent(in) :: datacomp
  integer,intent(in) :: ndftatom
  integer :: icount, ish, numprim, iatom, iprim, nang, nbf, ilocbf, ii, jj
  real(8),parameter :: half=0.5D+00, one=1.0D+00, two=2.0D+00
@@ -1601,23 +1608,24 @@ subroutine gridao(vao,xyzpt,rsqrd,threshex,ndftatom,databasis)
  enddo
  case default
  write(*,'(" Error! Subroutine Gridao supports up to i functions.")')
- call iabort
+ call iabort(datacomp)
  end select
  enddo
 !
  return
 end
 
 
-!------------------------------------------------------------------
- subroutine gridd2ao(vg2ao,xyzpt,rsqrd,threshex,ndftatom,databasis)
-!------------------------------------------------------------------
+!------------------------------------------------------------------------------
+ subroutine gridd2ao(vg2ao,xyzpt,rsqrd,threshex,ndftatom,databasis,datacomp)
+!------------------------------------------------------------------------------
 !
 ! Calculate second derivatives of AO values for a grid point
 !
- use modtype, only : typebasis
+ use modtype, only : typebasis, typecomp
  implicit none
  type(typebasis),intent(in) :: databasis
+ type(typecomp),intent(in) :: datacomp
  integer,intent(in) :: ndftatom
  integer :: icount, ish, numprim, iatom, iprim, nang, nbf, iloc, ii, jj
  real(8),parameter :: zero=0.0D+00, half=0.5D+00, one=1.0D+00, two=2.0D+00, three=3.0D+00
@@ -2370,7 +2378,7 @@ subroutine gridd2ao(vg2ao,xyzpt,rsqrd,threshex,ndftatom,databasis)
  enddo
  case default
  write(*,'(" Error! Subroutine Gridgao supports up to h functions.")')
- call iabort
+ call iabort(datacomp)
  end select
  enddo
 !
@@ -2451,7 +2459,7 @@ subroutine gradrexcor(egrad,edftgrad,cmo,fulldmtrx,atomvec,surface,radpt,angpt,r
  enddo
 !
  call gridraomo(vao,vmo,transcmo,xyzpt,rsqrd,aocutoff,datajob%threshex,ndftatom, &
-& datamol,databasis)
+& datamol,databasis,datacomp)
 !
  rhoa= zero
  grhoa(1:3)= zero
@@ -2465,7 +2473,7 @@ subroutine gradrexcor(egrad,edftgrad,cmo,fulldmtrx,atomvec,surface,radpt,angpt,r
  grhoa(1)= grhoa(1)*two
  grhoa(2)= grhoa(2)*two
  grhoa(3)= grhoa(3)*two
- call gridd2ao(vao(1,5),xyzpt,rsqrd,datajob%threshex,ndftatom,databasis)
+ call gridd2ao(vao(1,5),xyzpt,rsqrd,datajob%threshex,ndftatom,databasis,datacomp)
 !
 ! Calculate weight derivative
 !
@@ -2581,7 +2589,7 @@ subroutine graduexcor(egrad,edftgrad,cmoa,cmob,fulldmtrx1,fulldmtrx2,atomvec,sur
  enddo
 !
  call griduaomo(vao,vmoa,vmob,transcmoa,transcmob,xyzpt,rsqrd,aocutoff, &
-& datajob%threshex,ndftatom,datamol,databasis)
+& datajob%threshex,ndftatom,datamol,databasis,datacomp)
 !
  rhoa= zero
  grhoa(1:3)= zero
@@ -2602,7 +2610,7 @@ subroutine graduexcor(egrad,edftgrad,cmoa,cmob,fulldmtrx1,fulldmtrx2,atomvec,sur
  if((rhoa+rhob) < rcutoff)cycle
  grhoa(1:3)= grhoa(1:3)*two
  grhob(1:3)= grhob(1:3)*two
- call gridd2ao(vao(1,5),xyzpt,rsqrd,datajob%threshex,ndftatom,databasis)
+ call gridd2ao(vao(1,5),xyzpt,rsqrd,datajob%threshex,ndftatom,databasis,datacomp)
 !
 ! Calculate weight derivative
 !
@@ -2710,7 +2718,7 @@ subroutine calcdgridweight(dweight,dpa,pa,uvec,atomvec,surface,rr,sphweight,kato
  else
  if(abs(g4) > threshg4) then
  if(datacomp%master) write(*,'(" Error! The value G4 in calcdergridweight is large.")')
- call iabort
+ call iabort(datacomp)
  endif
  endif
  if(abs(cutji) > threshcut) then
@@ -2721,7 +2729,7 @@ subroutine calcdgridweight(dweight,dpa,pa,uvec,atomvec,surface,rr,sphweight,kato
  else
  if(abs(g4) > threshg4) then
  if(datacomp%master) write(*,'(" Error! The value G4 in calcdergridweight is large.")')
- call iabort
+ call iabort(datacomp)
  endif
  endif
  enddo

src/ecpder.F90

@@ -12,19 +12,20 @@
 ! See the License for the specific language governing permissions and
 ! limitations under the License.
 !
-!---------------------------------------------------------------------------
- subroutine gradoneeiecp(egrad1,fulldmtrx,nproc,myrank,datamol,databasis)
-!---------------------------------------------------------------------------
+!------------------------------------------------------------------------------------
+ subroutine gradoneeiecp(egrad1,fulldmtrx,nproc,myrank,datamol,databasis,datacomp)
+!------------------------------------------------------------------------------------
 !
 ! Calculate ECP derivative terms and add them into energy gradient matrix
 !
 ! Inout : egrad1 (one electron gradient matrix)
 !
  use modecp, only : nterm1, nterm2
- use modtype, only : typemol, typebasis
+ use modtype, only : typemol, typebasis, typecomp
  implicit none
  type(typemol),intent(in) :: datamol
  type(typebasis),intent(in) :: databasis
+ type(typecomp),intent(in) :: datacomp
  integer,intent(in) :: nproc, myrank
  integer :: maxfunc(0:6)=(/1,3,6,10,15,21,28/)
  integer :: maxbasis, numtbasis, maxdim, llmax, maxecpdim, nsizecp1, nsizecp2, ish, jsh
@@ -38,11 +39,11 @@ subroutine gradoneeiecp(egrad1,fulldmtrx,nproc,myrank,datamol,databasis)
  llmax= maxval(databasis%maxangecp(1:datamol%natom))
  if(maxbasis >= 6) then
  write(*,'(" Error! This program supports up to g function in ecp derivative calculation.")')
- call iabort
+ call iabort(datacomp)
  endif
  if(llmax >= 5) then
  write(*,'(" Error! This program supports up to SPDFG core potentials in ecp calculation.")')
- call iabort
+ call iabort(datacomp)
  endif
  maxecpdim= max(maxbasis,llmax-1)
  numtbasis=(maxecpdim+1)*(maxecpdim+2)*(maxecpdim+3)/6
@@ -64,7 +65,7 @@ subroutine gradoneeiecp(egrad1,fulldmtrx,nproc,myrank,datamol,databasis)
  do jsh= 1,databasis%nshell
  call calcdintecp(egrad1,fulldmtrx,term1ecp,term2ecp,term0ecp,xyzintecp, &
 & label1ecp,label2ecp,num1ecp,num2ecp,numtbasis,ish,jsh,maxdim, &
-& datamol,databasis)
+& datamol,databasis,datacomp)
  enddo
 !$OMP enddo
  enddo
@@ -81,7 +82,7 @@ subroutine gradoneeiecp(egrad1,fulldmtrx,nproc,myrank,datamol,databasis)
 !---------------------------------------------------------------------------------------
  subroutine calcdintecp(egrad1,fulldmtrx,term1ecp,term2ecp,term0ecp,xyzintecp, &
 & label1ecp,label2ecp,num1ecp,num2ecp,numtbasis,ish,jsh,len1, &
-& datamol,databasis)
+& datamol,databasis,datacomp)
 !---------------------------------------------------------------------------------------
 !
 ! Driver of ECP derivative terms
@@ -92,10 +93,11 @@ subroutine calcdintecp(egrad1,fulldmtrx,term1ecp,term2ecp,term0ecp,xyzintecp, &
 !
  use modparam, only : mxprsh
  use modecp, only : nterm1, nterm2
- use modtype, only : typemol, typebasis
+ use modtype, only : typemol, typebasis, typecomp
  implicit none
  type(typemol),intent(in) :: datamol
  type(typebasis),intent(in) :: databasis
+ type(typecomp),intent(in) :: datacomp
  integer,intent(in) :: label1ecp(nterm1*9), label2ecp(nterm2*6), num1ecp(*), num2ecp(*)
  integer,intent(in) :: numtbasis, ish, jsh, len1
  integer :: nangij(2), nprimij(2), nbfij(2), iatom, jatom, katom, iloc, jloc, ilocbf, jlocbf
@@ -159,7 +161,7 @@ subroutine calcdintecp(egrad1,fulldmtrx,term1ecp,term2ecp,term0ecp,xyzintecp, &
 !
  if((nangij(1) > 6).or.(nangij(2) > 6))then
  write(*,'(" Error! This program supports up to h function in calcdintecp")')
- call iabort
+ call iabort(datacomp)
  endif
 !
  do i= 1,3