source/unres/src_MD-M/checkder_p.F

   1       subroutine check_cartgrad
   2 C Check the gradient of Cartesian coordinates in internal coordinates.
   3       implicit real*8 (a-h,o-z)
   4       include 'DIMENSIONS'
   5       include 'COMMON.IOUNITS'
   6       include 'COMMON.VAR'
   7       include 'COMMON.CHAIN'
   8       include 'COMMON.GEO'
   9       include 'COMMON.LOCAL'
  10       include 'COMMON.DERIV'
  11       dimension temp(6,maxres),xx(3),gg(3)
  12       indmat(i,j)=((2*(nres-2)-i)*(i-1))/2+j-1
  13 *
  14 * Check the gradient of the virtual-bond and SC vectors in the internal
  15 * coordinates.
  16 *
  17       aincr=1.5d-7
  18       aincr2=2.5d-8
  19       call cartder
  20       write (iout,'(a)') '**************** dx/dalpha'
  21       write (iout,'(a)')
  22       do i=2,nres-1
  23         alphi=alph(i)
  24         alph(i)=alph(i)+aincr
  25         do k=1,3
  26           temp(k,i)=dc(k,nres+i)
  27         enddo
  28         call chainbuild
  29         do k=1,3
  30           gg(k)=(dc(k,nres+i)-temp(k,i))/aincr
  31           xx(k)=dabs((gg(k)-dxds(k,i))/(aincr*dabs(dxds(k,i))+aincr))
  32         enddo
  33         write (iout,'(i4,3e15.6/4x,3e15.6,3f9.3)')
  34      &  i,(gg(k),k=1,3),(dxds(k,i),k=1,3),(xx(k),k=1,3)
  35         write (iout,'(a)')
  36         alph(i)=alphi
  37         call chainbuild
  38       enddo
  39       write (iout,'(a)')
  40       write (iout,'(a)') '**************** dx/domega'
  41       write (iout,'(a)')
  42       do i=2,nres-1
  43         omegi=omeg(i)
  44         omeg(i)=omeg(i)+aincr
  45         do k=1,3
  46           temp(k,i)=dc(k,nres+i)
  47         enddo
  48         call chainbuild
  49         do k=1,3
  50           gg(k)=(dc(k,nres+i)-temp(k,i))/aincr
  51           xx(k)=dabs((gg(k)-dxds(k+3,i))/
  52      &          (aincr*dabs(dxds(k+3,i))+aincr))
  53         enddo
  54         write (iout,'(i4,3e15.6/4x,3e15.6,3f9.3)')
  55      &      i,(gg(k),k=1,3),(dxds(k+3,i),k=1,3),(xx(k),k=1,3)
  56         write (iout,'(a)')
  57         omeg(i)=omegi
  58         call chainbuild
  59       enddo
  60       write (iout,'(a)')
  61       write (iout,'(a)') '**************** dx/dtheta'
  62       write (iout,'(a)')
  63       do i=3,nres
  64         theti=theta(i)
  65         theta(i)=theta(i)+aincr
  66         do j=i-1,nres-1
  67           do k=1,3
  68             temp(k,j)=dc(k,nres+j)
  69           enddo
  70         enddo
  71         call chainbuild
  72         do j=i-1,nres-1
  73           ii = indmat(i-2,j)
  74 c         print *,'i=',i-2,' j=',j-1,' ii=',ii
  75           do k=1,3
  76             gg(k)=(dc(k,nres+j)-temp(k,j))/aincr
  77             xx(k)=dabs((gg(k)-dxdv(k,ii))/
  78      &            (aincr*dabs(dxdv(k,ii))+aincr))
  79           enddo
  80           write (iout,'(2i4,3e14.6/8x,3e14.6,3f9.3)')
  81      &        i,j,(gg(k),k=1,3),(dxdv(k,ii),k=1,3),(xx(k),k=1,3)
  82           write(iout,'(a)')
  83         enddo
  84         write (iout,'(a)')
  85         theta(i)=theti
  86         call chainbuild
  87       enddo
  88       write (iout,'(a)') '***************** dx/dphi'
  89       write (iout,'(a)')
  90       do i=4,nres
  91         phi(i)=phi(i)+aincr
  92         do j=i-1,nres-1
  93           do k=1,3
  94             temp(k,j)=dc(k,nres+j)
  95           enddo
  96         enddo
  97         call chainbuild
  98         do j=i-1,nres-1
  99           ii = indmat(i-2,j)
 100 c         print *,'ii=',ii
 101           do k=1,3
 102             gg(k)=(dc(k,nres+j)-temp(k,j))/aincr
 103             xx(k)=dabs((gg(k)-dxdv(k+3,ii))/
 104      &            (aincr*dabs(dxdv(k+3,ii))+aincr))
 105           enddo
 106           write (iout,'(2i4,3e14.6/8x,3e14.6,3f9.3)')
 107      &        i,j,(gg(k),k=1,3),(dxdv(k+3,ii),k=1,3),(xx(k),k=1,3)
 108           write(iout,'(a)')
 109         enddo
 110         phi(i)=phi(i)-aincr
 111         call chainbuild
 112       enddo
 113       write (iout,'(a)') '****************** ddc/dtheta'
 114       do i=1,nres-2
 115         thet=theta(i+2)
 116         theta(i+2)=thet+aincr
 117         do j=i,nres
 118           do k=1,3
 119             temp(k,j)=dc(k,j)
 120           enddo
 121         enddo
 122         call chainbuild
 123         do j=i+1,nres-1
 124           ii = indmat(i,j)
 125 c         print *,'ii=',ii
 126           do k=1,3
 127             gg(k)=(dc(k,j)-temp(k,j))/aincr
 128             xx(k)=dabs((gg(k)-dcdv(k,ii))/
 129      &           (aincr*dabs(dcdv(k,ii))+aincr))
 130           enddo
 131           write (iout,'(2i4,3e14.6/8x,3e14.6,3f9.3)')
 132      &           i,j,(gg(k),k=1,3),(dcdv(k,ii),k=1,3),(xx(k),k=1,3)
 133           write (iout,'(a)')
 134         enddo
 135         do j=1,nres
 136           do k=1,3
 137             dc(k,j)=temp(k,j)
 138           enddo
 139         enddo
 140         theta(i+2)=thet
 141       enddo
 142       write (iout,'(a)') '******************* ddc/dphi'
 143       do i=1,nres-3
 144         phii=phi(i+3)
 145         phi(i+3)=phii+aincr
 146         do j=1,nres
 147           do k=1,3
 148             temp(k,j)=dc(k,j)
 149           enddo
 150         enddo
 151         call chainbuild
 152         do j=i+2,nres-1
 153           ii = indmat(i+1,j)
 154 c         print *,'ii=',ii
 155           do k=1,3
 156             gg(k)=(dc(k,j)-temp(k,j))/aincr
 157             xx(k)=dabs((gg(k)-dcdv(k+3,ii))/
 158      &           (aincr*dabs(dcdv(k+3,ii))+aincr))
 159           enddo
 160           write (iout,'(2i4,3e14.6/8x,3e14.6,3f9.3)')
 161      &         i,j,(gg(k),k=1,3),(dcdv(k+3,ii),k=1,3),(xx(k),k=1,3)
 162           write (iout,'(a)')
 163         enddo
 164         do j=1,nres
 165           do k=1,3
 166             dc(k,j)=temp(k,j)
 167           enddo
 168         enddo
 169         phi(i+3)=phii
 170       enddo
 171       return
 172       end
 173 C----------------------------------------------------------------------------
 174       subroutine check_ecart
 175 C Check the gradient of the energy in Cartesian coordinates.
 176       implicit real*8 (a-h,o-z)
 177       include 'DIMENSIONS'
 178       include 'COMMON.CHAIN'
 179       include 'COMMON.DERIV'
 180       include 'COMMON.IOUNITS'
 181       include 'COMMON.VAR'
 182       include 'COMMON.CONTACTS'
 183       common /srutu/ icall
 184       dimension ggg(6),cc(3),xx(3),ddc(3),ddx(3),x(maxvar),g(maxvar)
 185       dimension grad_s(6,maxres)
 186       double precision energia(0:n_ene),energia1(0:n_ene)
 187       integer uiparm(1)
 188       double precision urparm(1)
 189       external fdum
 190       icg=1
 191       nf=0
 192       nfl=0
 193       call zerograd
 194       aincr=1.0D-7
 195       print '(a)','CG processor',me,' calling CHECK_CART.'
 196       nf=0
 197       icall=0
 198       call geom_to_var(nvar,x)
 199       call etotal(energia(0))
 200       etot=energia(0)
 201       call enerprint(energia(0))
 202       call gradient(nvar,x,nf,g,uiparm,urparm,fdum)
 203       icall =1
 204       do i=1,nres
 205         write (iout,'(i5,3f10.5)') i,(gradxorr(j,i),j=1,3)
 206       enddo
 207       do i=1,nres
 208         do j=1,3
 209           grad_s(j,i)=gradc(j,i,icg)
 210           grad_s(j+3,i)=gradx(j,i,icg)
 211         enddo
 212       enddo
 213       call flush(iout)
 214       write (iout,'(/a/)') 'Gradient in virtual-bond and SC vectors'
 215       do i=1,nres
 216         do j=1,3
 217           xx(j)=c(j,i+nres)
 218           ddc(j)=dc(j,i)
 219           ddx(j)=dc(j,i+nres)
 220         enddo
 221         do j=1,3
 222           dc(j,i)=dc(j,i)+aincr
 223           do k=i+1,nres
 224             c(j,k)=c(j,k)+aincr
 225             c(j,k+nres)=c(j,k+nres)+aincr
 226           enddo
 227           call etotal(energia1(0))
 228           etot1=energia1(0)
 229           ggg(j)=(etot1-etot)/aincr
 230           dc(j,i)=ddc(j)
 231           do k=i+1,nres
 232             c(j,k)=c(j,k)-aincr
 233             c(j,k+nres)=c(j,k+nres)-aincr
 234           enddo
 235         enddo
 236         do j=1,3
 237           c(j,i+nres)=c(j,i+nres)+aincr
 238           dc(j,i+nres)=dc(j,i+nres)+aincr
 239           call etotal(energia1(0))
 240           etot1=energia1(0)
 241           ggg(j+3)=(etot1-etot)/aincr
 242           c(j,i+nres)=xx(j)
 243           dc(j,i+nres)=ddx(j)
 244         enddo
 245         write (iout,'(i3,6(1pe12.5)/3x,6(1pe12.5)/)')
 246      &   i,(ggg(k),k=1,6),(grad_s(k,i),k=1,6)
 247       enddo
 248       return
 249       end
 250 c----------------------------------------------------------------------------
 251       subroutine check_ecartint
 252 C Check the gradient of the energy in Cartesian coordinates.
 253       implicit real*8 (a-h,o-z)
 254       include 'DIMENSIONS'
 255       include 'COMMON.CONTROL'
 256       include 'COMMON.CHAIN'
 257       include 'COMMON.DERIV'
 258       include 'COMMON.IOUNITS'
 259       include 'COMMON.VAR'
 260       include 'COMMON.CONTACTS'
 261       include 'COMMON.MD'
 262       include 'COMMON.LOCAL'
 263       include 'COMMON.SPLITELE'
 264       common /srutu/ icall
 265       dimension ggg(6),ggg1(6),cc(3),xx(3),ddc(3),ddx(3),x(maxvar),
 266      &  g(maxvar)
 267       dimension dcnorm_safe(3),dxnorm_safe(3)
 268       dimension grad_s(6,0:maxres),grad_s1(6,0:maxres)
 269       double precision phi_temp(maxres),theta_temp(maxres),
 270      &  alph_temp(maxres),omeg_temp(maxres)
 271       double precision energia(0:n_ene),energia1(0:n_ene)
 272       integer uiparm(1)
 273       double precision urparm(1)
 274       external fdum
 275 c      r_cut=2.0d0
 276 c      rlambd=0.3d0
 277       icg=1
 278       nf=0
 279       nfl=0
 280       print *,"ATU 3"
 281       call intout
 282 c      call intcartderiv
 283 c      call checkintcartgrad
 284       call zerograd
 285       aincr=8.0D-7
 286       write(iout,*) 'Calling CHECK_ECARTINT.'
 287       nf=0
 288       icall=0
 289       call geom_to_var(nvar,x)
 290       if (.not.split_ene) then
 291         call etotal(energia(0))
 292         etot=energia(0)
 293         call enerprint(energia(0))
 294         call flush(iout)
 295         write (iout,*) "enter cartgrad"
 296         call flush(iout)
 297         call cartgrad
 298         write (iout,*) "exit cartgrad"
 299         call flush(iout)
 300         icall =1
 301         do i=1,nres
 302           write (iout,'(i5,3f10.5)') i,(gradxorr(j,i),j=1,3)
 303         enddo
 304         do j=1,3
 305           grad_s(j,0)=gcart(j,0)
 306         enddo
 307         do i=1,nres
 308           do j=1,3
 309             grad_s(j,i)=gcart(j,i)
 310             grad_s(j+3,i)=gxcart(j,i)
 311           enddo
 312         enddo
 313       else
 314 !- split gradient check
 315         call zerograd
 316         call etotal_long(energia(0))
 317         call enerprint(energia(0))
 318         call flush(iout)
 319         write (iout,*) "enter cartgrad"
 320         call flush(iout)
 321         call cartgrad
 322         write (iout,*) "exit cartgrad"
 323         call flush(iout)
 324         icall =1
 325         write (iout,*) "longrange grad"
 326         do i=1,nres
 327           write (iout,'(i5,3f10.5,5x,3f10.5)') i,(gcart(j,i),j=1,3),
 328      &    (gxcart(j,i),j=1,3)
 329         enddo
 330         do j=1,3
 331           grad_s(j,0)=gcart(j,0)
 332         enddo
 333         do i=1,nres
 334           do j=1,3
 335             grad_s(j,i)=gcart(j,i)
 336             grad_s(j+3,i)=gxcart(j,i)
 337           enddo
 338         enddo
 339         call zerograd
 340         call etotal_short(energia(0))
 341         call enerprint(energia(0))
 342         call flush(iout)
 343         write (iout,*) "enter cartgrad"
 344         call flush(iout)
 345         call cartgrad
 346         write (iout,*) "exit cartgrad"
 347         call flush(iout)
 348         icall =1
 349         write (iout,*) "shortrange grad"
 350         do i=1,nres
 351           write (iout,'(i5,3f10.5,5x,3f10.5)') i,(gcart(j,i),j=1,3),
 352      &    (gxcart(j,i),j=1,3)
 353         enddo
 354         do j=1,3
 355           grad_s1(j,0)=gcart(j,0)
 356         enddo
 357         do i=1,nres
 358           do j=1,3
 359             grad_s1(j,i)=gcart(j,i)
 360             grad_s1(j+3,i)=gxcart(j,i)
 361           enddo
 362         enddo
 363       endif
 364       write (iout,'(/a/)') 'Gradient in virtual-bond and SC vectors'
 365       do i=0,nres
 366         print *,i
 367         do j=1,3
 368           xx(j)=c(j,i+nres)
 369           ddc(j)=dc(j,i)
 370           ddx(j)=dc(j,i+nres)
 371           do k=1,3
 372             dcnorm_safe(k)=dc_norm(k,i)
 373             dxnorm_safe(k)=dc_norm(k,i+nres)
 374           enddo
 375         enddo
 376         do j=1,3
 377           dc(j,i)=ddc(j)+aincr
 378           call chainbuild_cart
 379 #ifdef MPI
 380 c Broadcast the order to compute internal coordinates to the slaves.
 381 c          if (nfgtasks.gt.1)
 382 c     &      call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR)
 383 #endif
 384 c          call int_from_cart1(.false.)
 385           if (.not.split_ene) then
 386             call etotal(energia1(0))
 387             etot1=energia1(0)
 388           else
 389 !- split gradient
 390             call etotal_long(energia1(0))
 391             etot11=energia1(0)
 392             call etotal_short(energia1(0))
 393             etot12=energia1(0)
 394 c            write (iout,*) "etot11",etot11," etot12",etot12
 395           endif
 396 !- end split gradient
 397 c          write(iout,'(2i5,2(a,f15.10))')i,j," etot",etot," etot1",etot1
 398           dc(j,i)=ddc(j)-aincr
 399           call chainbuild_cart
 400 C          print *,c(j,i)
 401 c          call int_from_cart1(.false.)
 402           if (.not.split_ene) then
 403             call etotal(energia1(0))
 404             etot2=energia1(0)
 405             ggg(j)=(etot1-etot2)/(2*aincr)
 406           else
 407 !- split gradient
 408             call etotal_long(energia1(0))
 409             etot21=energia1(0)
 410             ggg(j)=(etot11-etot21)/(2*aincr)
 411             call etotal_short(energia1(0))
 412             etot22=energia1(0)
 413             ggg1(j)=(etot12-etot22)/(2*aincr)
 414 !- end split gradient
 415 c            write (iout,*) "etot21",etot21," etot22",etot22
 416           endif
 417 c          write(iout,'(2i5,2(a,f15.10))')i,j," etot",etot," etot2",etot2
 418           dc(j,i)=ddc(j)
 419           call chainbuild_cart
 420         enddo
 421         do j=1,3
 422           dc(j,i+nres)=ddx(j)+aincr
 423           call chainbuild_cart
 424 c          write (iout,*) "i",i," j",j," dxnorm+ and dxnorm"
 425 c          write (iout,'(3f15.10)') (dc_norm(k,i+nres),k=1,3)
 426 c          write (iout,'(3f15.10)') (dxnorm_safe(k),k=1,3)
 427 c          write (iout,*) "dxnormnorm",dsqrt(
 428 c     &  dc_norm(1,i+nres)**2+dc_norm(2,i+nres)**2+dc_norm(3,i+nres)**2)
 429 c          write (iout,*) "dxnormnormsafe",dsqrt(
 430 c     &      dxnorm_safe(1)**2+dxnorm_safe(2)**2+dxnorm_safe(3)**2)
 431 c          write (iout,*)
 432           if (.not.split_ene) then
 433             call etotal(energia1(0))
 434             etot1=energia1(0)
 435           else
 436 !- split gradient
 437             call etotal_long(energia1(0))
 438             etot11=energia1(0)
 439             call etotal_short(energia1(0))
 440             etot12=energia1(0)
 441           endif
 442 !- end split gradient
 443 c          write(iout,'(2i5,2(a,f15.10))')i,j," etot",etot," etot1",etot1
 444           dc(j,i+nres)=ddx(j)-aincr
 445           call chainbuild_cart
 446 c          write (iout,*) "i",i," j",j," dxnorm- and dxnorm"
 447 c          write (iout,'(3f15.10)') (dc_norm(k,i+nres),k=1,3)
 448 c          write (iout,'(3f15.10)') (dxnorm_safe(k),k=1,3)
 449 c          write (iout,*)
 450 c          write (iout,*) "dxnormnorm",dsqrt(
 451 c     &  dc_norm(1,i+nres)**2+dc_norm(2,i+nres)**2+dc_norm(3,i+nres)**2)
 452 c          write (iout,*) "dxnormnormsafe",dsqrt(
 453 c     &      dxnorm_safe(1)**2+dxnorm_safe(2)**2+dxnorm_safe(3)**2)
 454           if (.not.split_ene) then
 455             call etotal(energia1(0))
 456             etot2=energia1(0)
 457             ggg(j+3)=(etot1-etot2)/(2*aincr)
 458           else
 459 !- split gradient
 460             call etotal_long(energia1(0))
 461             etot21=energia1(0)
 462             ggg(j+3)=(etot11-etot21)/(2*aincr)
 463             call etotal_short(energia1(0))
 464             etot22=energia1(0)
 465             ggg1(j+3)=(etot12-etot22)/(2*aincr)
 466 !- end split gradient
 467           endif
 468 c          write(iout,'(2i5,2(a,f15.10))')i,j," etot",etot," etot2",etot2
 469           dc(j,i+nres)=ddx(j)
 470           call chainbuild_cart
 471         enddo
 472         write (iout,'(i3,6(1pe12.5)/3x,6(1pe12.5)/3x,6(1pe12.5)/)')
 473      &   i,(ggg(k),k=1,6),(grad_s(k,i),k=1,6),(ggg(k)/grad_s(k,i),k=1,6)
 474         if (split_ene) then
 475           write (iout,'(i3,6(1pe12.5)/3x,6(1pe12.5)/3x,6(1pe12.5)/)')
 476      &   i,(ggg1(k),k=1,6),(grad_s1(k,i),k=1,6),(ggg1(k)/grad_s1(k,i),
 477      &   k=1,6)
 478          write (iout,'(i3,6(1pe12.5)/3x,6(1pe12.5)/3x,6(1pe12.5)/)')
 479      &   i,(ggg(k)+ggg1(k),k=1,6),(grad_s(k,i)+grad_s1(k,i),k=1,6),
 480      &   ((ggg(k)+ggg1(k))/(grad_s(k,i)+grad_s1(k,i)),k=1,6)
 481         endif
 482       enddo
 483       return
 484       end
 485 c-------------------------------------------------------------------------
 486       subroutine int_from_cart1(lprn)
 487       implicit real*8 (a-h,o-z)
 488       include 'DIMENSIONS'
 489 #ifdef MPI
 490       include 'mpif.h'
 491       integer ierror
 492 #endif
 493       include 'COMMON.IOUNITS'
 494       include 'COMMON.VAR'
 495       include 'COMMON.CHAIN'
 496       include 'COMMON.GEO'
 497       include 'COMMON.INTERACT'
 498       include 'COMMON.LOCAL'
 499       include 'COMMON.NAMES'
 500       include 'COMMON.SETUP'
 501       include 'COMMON.TIME1'
 502       logical lprn
 503       if (lprn) write (iout,'(/a)') 'Recalculated internal coordinates'
 504 #ifdef TIMING
 505       time01=MPI_Wtime()
 506 #endif
 507 #if defined(PARINT) && defined(MPI)
 508       do i=iint_start,iint_end
 509 #else
 510       do i=2,nres
 511 #endif
 512 C        print *,i
 513         dnorm1=dist(i-1,i)
 514         dnorm2=dist(i,i+1)
 515 C        print *,i,dnorm1,dnorm2
 516         do j=1,3
 517           c(j,maxres2)=0.5D0*(2*c(j,i)+(c(j,i-1)-c(j,i))/dnorm1
 518      &     +(c(j,i+1)-c(j,i))/dnorm2)
 519         enddo
 520         be=0.0D0
 521         if (i.gt.2) then
 522         if (i.le.nres) phi(i+1)=beta(i-2,i-1,i,i+1)
 523         if ((itype(i).ne.10).and.(itype(i-1).ne.10)) then
 524          tauangle(3,i+1)=beta(i+nres-1,i-1,i,i+nres)
 525         endif
 526         if (itype(i-1).ne.10) then
 527          tauangle(1,i+1)=beta(i-1+nres,i-1,i,i+1)
 528          omicron(1,i)=alpha(i-2,i-1,i-1+nres)
 529          omicron(2,i)=alpha(i-1+nres,i-1,i)
 530         endif
 531         if (itype(i).ne.10) then
 532          tauangle(2,i+1)=beta(i-2,i-1,i,i+nres)
 533         endif
 534         endif
 535         omeg(i)=beta(nres+i,i,maxres2,i+1)
 536 C        print *,omeg(i)
 537         alph(i)=alpha(nres+i,i,maxres2)
 538 C        print *,alph(i)
 539         theta(i+1)=alpha(i-1,i,i+1)
 540         vbld(i)=dist(i-1,i)
 541 C        print *,vbld(i)
 542         vbld_inv(i)=1.0d0/vbld(i)
 543         vbld(nres+i)=dist(nres+i,i)
 544 C            print *,vbld(i+nres)
 545
 546         if (itype(i).ne.10) then
 547           vbld_inv(nres+i)=1.0d0/vbld(nres+i)
 548         else
 549           vbld_inv(nres+i)=0.0d0
 550         endif
 551       enddo
 552 #if defined(PARINT) && defined(MPI)
 553        if (nfgtasks1.gt.1) then
 554 cd       write(iout,*) "iint_start",iint_start," iint_count",
 555 cd     &   (iint_count(i),i=0,nfgtasks-1)," iint_displ",
 556 cd     &   (iint_displ(i),i=0,nfgtasks-1)
 557 cd       write (iout,*) "Gather vbld backbone"
 558 cd       call flush(iout)
 559        time00=MPI_Wtime()
 560        call MPI_Allgatherv(vbld(iint_start),iint_count(fg_rank1),
 561      &   MPI_DOUBLE_PRECISION,vbld(1),iint_count(0),iint_displ(0),
 562      &   MPI_DOUBLE_PRECISION,FG_COMM1,IERR)
 563 cd       write (iout,*) "Gather vbld_inv"
 564 cd       call flush(iout)
 565        call MPI_Allgatherv(vbld_inv(iint_start),iint_count(fg_rank1),
 566      &   MPI_DOUBLE_PRECISION,vbld_inv(1),iint_count(0),iint_displ(0),
 567      &   MPI_DOUBLE_PRECISION,FG_COMM1,IERR)
 568 cd       write (iout,*) "Gather vbld side chain"
 569 cd       call flush(iout)
 570        call MPI_Allgatherv(vbld(iint_start+nres),iint_count(fg_rank1),
 571      &  MPI_DOUBLE_PRECISION,vbld(nres+1),iint_count(0),iint_displ(0),
 572      &   MPI_DOUBLE_PRECISION,FG_COMM1,IERR)
 573 cd       write (iout,*) "Gather vbld_inv side chain"
 574 cd       call flush(iout)
 575        call MPI_Allgatherv(vbld_inv(iint_start+nres),
 576      &   iint_count(fg_rank1),MPI_DOUBLE_PRECISION,vbld_inv(nres+1),
 577      &   iint_count(0),iint_displ(0),MPI_DOUBLE_PRECISION,FG_COMM1,IERR)
 578 cd       write (iout,*) "Gather theta"
 579 cd       call flush(iout)
 580        call MPI_Allgatherv(theta(iint_start+1),iint_count(fg_rank1),
 581      &   MPI_DOUBLE_PRECISION,theta(2),iint_count(0),iint_displ(0),
 582      &   MPI_DOUBLE_PRECISION,FG_COMM1,IERR)
 583 cd       write (iout,*) "Gather phi"
 584 cd       call flush(iout)
 585        call MPI_Allgatherv(phi(iint_start+1),iint_count(fg_rank1),
 586      &   MPI_DOUBLE_PRECISION,phi(2),iint_count(0),iint_displ(0),
 587      &   MPI_DOUBLE_PRECISION,FG_COMM1,IERR)
 588 #ifdef CRYST_SC
 589 cd       write (iout,*) "Gather alph"
 590 cd       call flush(iout)
 591        call MPI_Allgatherv(alph(iint_start),iint_count(fg_rank1),
 592      &   MPI_DOUBLE_PRECISION,alph(1),iint_count(0),iint_displ(0),
 593      &   MPI_DOUBLE_PRECISION,FG_COMM1,IERR)
 594 cd       write (iout,*) "Gather omeg"
 595 cd       call flush(iout)
 596        call MPI_Allgatherv(omeg(iint_start),iint_count(fg_rank1),
 597      &   MPI_DOUBLE_PRECISION,omeg(1),iint_count(0),iint_displ(0),
 598      &   MPI_DOUBLE_PRECISION,FG_COMM1,IERR)
 599 #endif
 600        time_gather=time_gather+MPI_Wtime()-time00
 601       endif
 602 #endif
 603       do i=1,nres-1
 604         do j=1,3
 605           dc_norm(j,i)=dc(j,i)*vbld_inv(i+1)
 606         enddo
 607       enddo
 608       do i=2,nres-1
 609         do j=1,3
 610           dc_norm(j,i+nres)=dc(j,i+nres)*vbld_inv(i+nres)
 611         enddo
 612       enddo
 613       if (lprn) then
 614       do i=2,nres
 615        write (iout,1212) restyp(itype(i)),i,vbld(i),
 616      &rad2deg*theta(i),rad2deg*phi(i),vbld(nres+i),
 617      &rad2deg*alph(i),rad2deg*omeg(i)
 618       enddo
 619       endif
 620  1212 format (a3,'(',i3,')',2(f15.10,2f10.2))
 621 #ifdef TIMING
 622       time_intfcart=time_intfcart+MPI_Wtime()-time01
 623 #endif
 624       return
 625       end
 626 c----------------------------------------------------------------------------
 627       subroutine check_eint
 628 C Check the gradient of energy in internal coordinates.
 629       implicit real*8 (a-h,o-z)
 630       include 'DIMENSIONS'
 631       include 'COMMON.CHAIN'
 632       include 'COMMON.DERIV'
 633       include 'COMMON.IOUNITS'
 634       include 'COMMON.VAR'
 635       include 'COMMON.GEO'
 636       common /srutu/ icall
 637       dimension x(maxvar),gana(maxvar),gg(maxvar)
 638       integer uiparm(1)
 639       double precision urparm(1)
 640       double precision energia(0:n_ene),energia1(0:n_ene),
 641      &  energia2(0:n_ene)
 642       character*6 key
 643       external fdum
 644       call zerograd
 645       aincr=1.0D-7
 646       print '(a)','Calling CHECK_INT.'
 647       nf=0
 648       nfl=0
 649       icg=1
 650       call geom_to_var(nvar,x)
 651       call var_to_geom(nvar,x)
 652       call chainbuild
 653       icall=1
 654       print *,'ICG=',ICG
 655       call etotal(energia(0))
 656       etot = energia(0)
 657       call enerprint(energia(0))
 658       print *,'ICG=',ICG
 659 #ifdef MPL
 660       if (MyID.ne.BossID) then
 661         call mp_bcast(x(1),8*(nvar+3),BossID,fgGroupID)
 662         nf=x(nvar+1)
 663         nfl=x(nvar+2)
 664         icg=x(nvar+3)
 665       endif
 666 #endif
 667       nf=1
 668       nfl=3
 669 cd    write (iout,'(10f8.3)') (rad2deg*x(i),i=1,nvar)
 670       call gradient(nvar,x,nf,gana,uiparm,urparm,fdum)
 671 cd    write (iout,'(i3,1pe14.4)') (i,gana(i),i=1,nvar)
 672       icall=1
 673       do i=1,nvar
 674         xi=x(i)
 675         x(i)=xi-0.5D0*aincr
 676         call var_to_geom(nvar,x)
 677         call chainbuild
 678         call etotal(energia1(0))
 679         etot1=energia1(0)
 680         x(i)=xi+0.5D0*aincr
 681         call var_to_geom(nvar,x)
 682         call chainbuild
 683         call etotal(energia2(0))
 684         etot2=energia2(0)
 685         gg(i)=(etot2-etot1)/aincr
 686         write (iout,*) i,etot1,etot2
 687         x(i)=xi
 688       enddo
 689       write (iout,'(/2a)')' Variable        Numerical       Analytical',
 690      &    '     RelDiff*100% '
 691       do i=1,nvar
 692         if (i.le.nphi) then
 693           ii=i
 694           key = ' phi'
 695         else if (i.le.nphi+ntheta) then
 696           ii=i-nphi
 697           key=' theta'
 698         else if (i.le.nphi+ntheta+nside) then
 699            ii=i-(nphi+ntheta)
 700            key=' alpha'
 701         else
 702            ii=i-(nphi+ntheta+nside)
 703            key=' omega'
 704         endif
 705         write (iout,'(i3,a,i3,3(1pd16.6))')
 706      & i,key,ii,gg(i),gana(i),
 707      & 100.0D0*dabs(gg(i)-gana(i))/(dabs(gana(i))+aincr)
 708       enddo
 709       return
 710       end