subroutine minimize(etot,x,iretcode,nfun) implicit real*8 (a-h,o-z) include 'DIMENSIONS' parameter (liv=60,lv=(77+maxvar*(maxvar+17)/2)) ********************************************************************* * OPTIMIZE sets up SUMSL or DFP and provides a simple interface for * * the calling subprogram. * * when d(i)=1.0, then v(35) is the length of the initial step, * * calculated in the usual pythagorean way. * * absolute convergence occurs when the function is within v(31) of * * zero. unless you know the minimum value in advance, abs convg * * is probably not useful. * * relative convergence is when the model predicts that the function * * will decrease by less than v(32)*abs(fun). * ********************************************************************* include 'COMMON.IOUNITS' include 'COMMON.VAR' include 'COMMON.GEO' include 'COMMON.MINIM' common /srutu/ icall dimension iv(liv) double precision minval,x(maxvar),d(maxvar),v(1:lv),xx(maxvar) double precision energia(0:n_ene) external func,gradient,fdum external func_restr,grad_restr logical not_done,change,reduce common /przechowalnia/ v icall = 1 NOT_DONE=.TRUE. c DO WHILE (NOT_DONE) call deflt(2,iv,liv,lv,v) * 12 means fresh start, dont call deflt iv(1)=12 * max num of fun calls if (maxfun.eq.0) maxfun=500 iv(17)=maxfun * max num of iterations if (maxmin.eq.0) maxmin=1000 iv(18)=maxmin * controls output iv(19)=2 * selects output unit iv(21)=0 if (print_min_ini+print_min_stat+print_min_res.gt.0) iv(21)=iout * 1 means to print out result iv(22)=print_min_res * 1 means to print out summary stats iv(23)=print_min_stat * 1 means to print initial x and d iv(24)=print_min_ini * min val for v(radfac) default is 0.1 v(24)=0.1D0 * max val for v(radfac) default is 4.0 v(25)=2.0D0 c v(25)=4.0D0 * check false conv if (act fnctn decrease) .lt. v(26)*(exp decrease) * the sumsl default is 0.1 v(26)=0.1D0 * false conv if (act fnctn decrease) .lt. v(34) * the sumsl default is 100*machep v(34)=v(34)/100.0D0 * absolute convergence if (tolf.eq.0.0D0) tolf=1.0D-4 v(31)=tolf * relative convergence if (rtolf.eq.0.0D0) rtolf=1.0D-4 v(32)=rtolf * controls initial step size v(35)=1.0D-1 * large vals of d correspond to small components of step do i=1,nphi d(i)=1.0D-1 enddo do i=nphi+1,nvar d(i)=1.0D-1 enddo cd print *,'Calling SUMSL' c call var_to_geom(nvar,x) c call chainbuild c call etotal(energia(0)) c etot = energia(0) IF (mask_r) THEN call x2xx(x,xx,nvar_restr) call sumsl(nvar_restr,d,xx,func_restr,grad_restr, & iv,liv,lv,v,idum,rdum,fdum) call xx2x(x,xx) ELSE call sumsl(nvar,d,x,func,gradient,iv,liv,lv,v,idum,rdum,fdum) ENDIF etot=v(10) iretcode=iv(1) cd print *,'Exit SUMSL; return code:',iretcode,' energy:',etot cd write (iout,'(/a,i4/)') 'SUMSL return code:',iv(1) c call intout c change=reduce(x) call var_to_geom(nvar,x) c if (change) then c write (iout,'(a)') 'Reduction worked, minimizing again...' c else c not_done=.false. c endif call chainbuild c call etotal(energia(0)) c etot=energia(0) c call enerprint(energia(0)) nfun=iv(6) c write (*,*) 'Processor',MyID,' leaves MINIMIZE.' c ENDDO ! NOT_DONE return end #ifdef MPI c---------------------------------------------------------------------------- subroutine ergastulum implicit real*8 (a-h,o-z) include 'DIMENSIONS' #ifdef MPI include "mpif.h" #endif include 'COMMON.SETUP' include 'COMMON.DERIV' include 'COMMON.VAR' include 'COMMON.IOUNITS' include 'COMMON.FFIELD' include 'COMMON.INTERACT' include 'COMMON.MD' include 'COMMON.TIME1' double precision z(maxres6),d_a_tmp(maxres6) double precision edum(0:n_ene),time_order(0:10) double precision Gcopy(maxres2,maxres2) common /przechowalnia/ Gcopy integer icall /0/ C Workers wait for variables and NF, and NFL from the boss iorder=0 do while (iorder.ge.0) c write (*,*) 'Processor',fg_rank,' CG group',kolor, c & ' receives order from Master' time00=MPI_Wtime() call MPI_Bcast(iorder,1,MPI_INTEGER,king,FG_COMM,IERR) time_Bcast=time_Bcast+MPI_Wtime()-time00 if (icall.gt.4 .and. iorder.ge.0) & time_order(iorder)=time_order(iorder)+MPI_Wtime()-time00 icall=icall+1 c write (*,*) c & 'Processor',fg_rank,' completed receive MPI_BCAST order',iorder if (iorder.eq.0) then call zerograd call etotal(edum) c write (2,*) "After etotal" c write (2,*) "dimen",dimen," dimen3",dimen3 c call flush(2) else if (iorder.eq.2) then call zerograd call etotal_short(edum) c write (2,*) "After etotal_short" c write (2,*) "dimen",dimen," dimen3",dimen3 c call flush(2) else if (iorder.eq.3) then call zerograd call etotal_long(edum) c write (2,*) "After etotal_long" c write (2,*) "dimen",dimen," dimen3",dimen3 c call flush(2) else if (iorder.eq.1) then call sum_gradient c write (2,*) "After sum_gradient" c write (2,*) "dimen",dimen," dimen3",dimen3 c call flush(2) else if (iorder.eq.4) then call ginv_mult(z,d_a_tmp) else if (iorder.eq.5) then c Setup MD things for a slave dimen=(nct-nnt+1)+nside dimen1=(nct-nnt)+(nct-nnt+1) dimen3=dimen*3 c write (2,*) "dimen",dimen," dimen3",dimen3 c call flush(2) call int_bounds(dimen,igmult_start,igmult_end) igmult_start=igmult_start-1 call MPI_Allgather(3*igmult_start,1,MPI_INTEGER, & ng_start(0),1,MPI_INTEGER,FG_COMM,IERROR) my_ng_count=igmult_end-igmult_start call MPI_Allgather(3*my_ng_count,1,MPI_INTEGER,ng_counts(0),1, & MPI_INTEGER,FG_COMM,IERROR) c write (2,*) "ng_start",(ng_start(i),i=0,nfgtasks-1) c write (2,*) "ng_counts",(ng_counts(i),i=0,nfgtasks-1) myginv_ng_count=maxres2*my_ng_count c write (2,*) "igmult_start",igmult_start," igmult_end", c & igmult_end," my_ng_count",my_ng_count c call flush(2) call MPI_Allgather(maxres2*igmult_start,1,MPI_INTEGER, & nginv_start(0),1,MPI_INTEGER,FG_COMM,IERROR) call MPI_Allgather(myginv_ng_count,1,MPI_INTEGER, & nginv_counts(0),1,MPI_INTEGER,FG_COMM,IERROR) c write (2,*) "nginv_start",(nginv_start(i),i=0,nfgtasks-1) c write (2,*) "nginv_counts",(nginv_counts(i),i=0,nfgtasks-1) c call flush(2) c call MPI_Barrier(FG_COMM,IERROR) time00=MPI_Wtime() call MPI_Scatterv(ginv(1,1),nginv_counts(0), & nginv_start(0),MPI_DOUBLE_PRECISION,gcopy(1,1), & myginv_ng_count,MPI_DOUBLE_PRECISION,king,FG_COMM,IERR) #ifdef TIMING time_scatter_ginv=time_scatter_ginv+MPI_Wtime()-time00 #endif do i=1,dimen do j=1,2*my_ng_count ginv(j,i)=gcopy(i,j) enddo enddo c write (2,*) "dimen",dimen," dimen3",dimen3 c write (2,*) "End MD setup" c call flush(2) c write (iout,*) "My chunk of ginv_block" c call MATOUT2(my_ng_count,dimen3,maxres2,maxers2,ginv_block) else if (iorder.eq.6) then call int_from_cart1(.false.) else if (iorder.eq.7) then call chainbuild_cart else if (iorder.eq.8) then call intcartderiv else if (iorder.eq.9) then call fricmat_mult(z,d_a_tmp) else if (iorder.eq.10) then call setup_fricmat endif enddo write (*,*) 'Processor',fg_rank,' CG group',kolor, & ' absolute rank',myrank,' leves ERGASTULUM.' write(*,*)'Processor',fg_rank,' wait times for respective orders', & (' order[',i,']',time_order(i),i=0,10) return end #endif ************************************************************************ subroutine func(n,x,nf,f,uiparm,urparm,ufparm) implicit real*8 (a-h,o-z) include 'DIMENSIONS' include 'COMMON.DERIV' include 'COMMON.IOUNITS' include 'COMMON.GEO' common /chuju/ jjj double precision energia(0:n_ene) integer jjj double precision ufparm external ufparm integer uiparm(1) real*8 urparm(1) dimension x(maxvar) c if (jjj.gt.0) then c write (iout,'(10f8.3)') (rad2deg*x(i),i=1,n) c endif nfl=nf icg=mod(nf,2)+1 cd print *,'func',nf,nfl,icg call var_to_geom(n,x) call zerograd call chainbuild cd write (iout,*) 'ETOTAL called from FUNC' call etotal(energia(0)) call sum_gradient f=energia(0) c if (jjj.gt.0) then c write (iout,'(10f8.3)') (rad2deg*x(i),i=1,n) c write (iout,*) 'f=',etot c jjj=0 c endif return end ************************************************************************ subroutine func_restr(n,x,nf,f,uiparm,urparm,ufparm) implicit real*8 (a-h,o-z) include 'DIMENSIONS' include 'COMMON.DERIV' include 'COMMON.IOUNITS' include 'COMMON.GEO' common /chuju/ jjj double precision energia(0:n_ene) integer jjj double precision ufparm external ufparm integer uiparm(1) real*8 urparm(1) dimension x(maxvar) c if (jjj.gt.0) then c write (iout,'(10f8.3)') (rad2deg*x(i),i=1,n) c endif nfl=nf icg=mod(nf,2)+1 call var_to_geom_restr(n,x) call zerograd call chainbuild cd write (iout,*) 'ETOTAL called from FUNC' call etotal(energia(0)) call sum_gradient f=energia(0) c if (jjj.gt.0) then c write (iout,'(10f8.3)') (rad2deg*x(i),i=1,n) c write (iout,*) 'f=',etot c jjj=0 c endif return end c------------------------------------------------------- subroutine x2xx(x,xx,n) implicit real*8 (a-h,o-z) include 'DIMENSIONS' include 'COMMON.VAR' include 'COMMON.CHAIN' include 'COMMON.INTERACT' double precision xx(maxvar),x(maxvar) do i=1,nvar varall(i)=x(i) enddo ig=0 igall=0 do i=4,nres igall=igall+1 if (mask_phi(i).eq.1) then ig=ig+1 xx(ig)=x(igall) endif enddo do i=3,nres igall=igall+1 if (mask_theta(i).eq.1) then ig=ig+1 xx(ig)=x(igall) endif enddo do ij=1,2 do i=2,nres-1 if (itype(i).ne.10) then igall=igall+1 if (mask_side(i).eq.1) then ig=ig+1 xx(ig)=x(igall) endif endif enddo enddo n=ig return end subroutine xx2x(x,xx) implicit real*8 (a-h,o-z) include 'DIMENSIONS' include 'COMMON.VAR' include 'COMMON.CHAIN' include 'COMMON.INTERACT' double precision xx(maxvar),x(maxvar) do i=1,nvar x(i)=varall(i) enddo ig=0 igall=0 do i=4,nres igall=igall+1 if (mask_phi(i).eq.1) then ig=ig+1 x(igall)=xx(ig) endif enddo do i=3,nres igall=igall+1 if (mask_theta(i).eq.1) then ig=ig+1 x(igall)=xx(ig) endif enddo do ij=1,2 do i=2,nres-1 if (itype(i).ne.10) then igall=igall+1 if (mask_side(i).eq.1) then ig=ig+1 x(igall)=xx(ig) endif endif enddo enddo return end c---------------------------------------------------------- subroutine minim_dc(etot,iretcode,nfun) implicit real*8 (a-h,o-z) include 'DIMENSIONS' parameter (liv=60,lv=(77+maxvar*(maxvar+17)/2)) #ifdef MPI include 'mpif.h' #endif include 'COMMON.SETUP' include 'COMMON.IOUNITS' include 'COMMON.VAR' include 'COMMON.GEO' include 'COMMON.MINIM' include 'COMMON.CHAIN' dimension iv(liv) double precision minval,x(maxvar),d(maxvar),v(1:lv),xx(maxvar) common /przechowalnia/ v double precision energia(0:n_ene) external func_dc,grad_dc,fdum logical not_done,change,reduce double precision g(maxvar),f1 call deflt(2,iv,liv,lv,v) * 12 means fresh start, dont call deflt iv(1)=12 * max num of fun calls if (maxfun.eq.0) maxfun=500 iv(17)=maxfun * max num of iterations if (maxmin.eq.0) maxmin=1000 iv(18)=maxmin * controls output iv(19)=2 * selects output unit iv(21)=0 if (print_min_ini+print_min_stat+print_min_res.gt.0) iv(21)=iout * 1 means to print out result iv(22)=print_min_res * 1 means to print out summary stats iv(23)=print_min_stat * 1 means to print initial x and d iv(24)=print_min_ini * min val for v(radfac) default is 0.1 v(24)=0.1D0 * max val for v(radfac) default is 4.0 v(25)=2.0D0 c v(25)=4.0D0 * check false conv if (act fnctn decrease) .lt. v(26)*(exp decrease) * the sumsl default is 0.1 v(26)=0.1D0 * false conv if (act fnctn decrease) .lt. v(34) * the sumsl default is 100*machep v(34)=v(34)/100.0D0 * absolute convergence if (tolf.eq.0.0D0) tolf=1.0D-4 v(31)=tolf * relative convergence if (rtolf.eq.0.0D0) rtolf=1.0D-4 v(32)=rtolf * controls initial step size v(35)=1.0D-1 * large vals of d correspond to small components of step do i=1,6*nres d(i)=1.0D-1 enddo k=0 do i=1,nres-1 do j=1,3 k=k+1 x(k)=dc(j,i) enddo enddo do i=2,nres-1 if (ialph(i,1).gt.0) then do j=1,3 k=k+1 x(k)=dc(j,i+nres) enddo endif enddo call sumsl(k,d,x,func_dc,grad_dc,iv,liv,lv,v,idum,rdum,fdum) k=0 do i=1,nres-1 do j=1,3 k=k+1 dc(j,i)=x(k) enddo enddo do i=2,nres-1 if (ialph(i,1).gt.0) then do j=1,3 k=k+1 dc(j,i+nres)=x(k) enddo endif enddo call chainbuild_cart cd call zerograd cd nf=0 cd call func_dc(k,x,nf,f,idum,rdum,fdum) cd call grad_dc(k,x,nf,g,idum,rdum,fdum) cd cd do i=1,k cd x(i)=x(i)+1.0D-5 cd call func_dc(k,x,nf,f1,idum,rdum,fdum) cd x(i)=x(i)-1.0D-5 cd print '(i5,2f15.5)',i,g(i),(f1-f)/1.0D-5 cd enddo etot=v(10) iretcode=iv(1) nfun=iv(6) return end subroutine func_dc(n,x,nf,f,uiparm,urparm,ufparm) implicit real*8 (a-h,o-z) include 'DIMENSIONS' #ifdef MPI include 'mpif.h' #endif include 'COMMON.SETUP' include 'COMMON.DERIV' include 'COMMON.IOUNITS' include 'COMMON.GEO' include 'COMMON.CHAIN' include 'COMMON.VAR' double precision energia(0:n_ene) double precision ufparm external ufparm integer uiparm(1) real*8 urparm(1) dimension x(maxvar) nfl=nf cbad icg=mod(nf,2)+1 icg=1 k=0 do i=1,nres-1 do j=1,3 k=k+1 dc(j,i)=x(k) enddo enddo do i=2,nres-1 if (ialph(i,1).gt.0) then do j=1,3 k=k+1 dc(j,i+nres)=x(k) enddo endif enddo call chainbuild_cart call zerograd call etotal(energia(0)) f=energia(0) cd print *,'func_dc ',nf,nfl,f return end subroutine grad_dc(n,x,nf,g,uiparm,urparm,ufparm) implicit real*8 (a-h,o-z) include 'DIMENSIONS' #ifdef MPI include 'mpif.h' #endif include 'COMMON.SETUP' include 'COMMON.CHAIN' include 'COMMON.DERIV' include 'COMMON.VAR' include 'COMMON.INTERACT' include 'COMMON.FFIELD' include 'COMMON.MD' include 'COMMON.IOUNITS' external ufparm integer uiparm(1),k double precision urparm(1) dimension x(maxvar),g(maxvar) c c c cbad icg=mod(nf,2)+1 icg=1 cd print *,'grad_dc ',nf,nfl,nf-nfl+1,icg if (nf-nfl+1) 20,30,40 20 call func_dc(n,x,nf,f,uiparm,urparm,ufparm) cd print *,20 if (nf.eq.0) return goto 40 30 continue cd print *,30 k=0 do i=1,nres-1 do j=1,3 k=k+1 dc(j,i)=x(k) enddo enddo do i=2,nres-1 if (ialph(i,1).gt.0) then do j=1,3 k=k+1 dc(j,i+nres)=x(k) enddo endif enddo call chainbuild_cart C C Evaluate the derivatives of virtual bond lengths and SC vectors in variables. C 40 call cartgrad cd print *,40 k=0 do i=1,nres-1 do j=1,3 k=k+1 g(k)=gcart(j,i) enddo enddo do i=2,nres-1 if (ialph(i,1).gt.0) then do j=1,3 k=k+1 g(k)=gxcart(j,i) enddo endif enddo return end