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