+++ /dev/null
- 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
-c 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
-cmd 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
-cmd 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
-cmd 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
-cmd call fricmat_mult(z,d_a_tmp)
- else if (iorder.eq.10) then
-cmd 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)
-c 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
projects / unres.git / blobdiff