correction in readpdb (bugfix)
[unres.git] / source / unres / src_MIN / unres_min.F
1 CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
2 C                                                                              C
3 C                                U N R E S                                     C
4 C                                                                              C
5 C Program to carry out conformational search of proteins in an united-residue  C
6 C approximation.                                                               C
7 C                                                                              C
8 CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
9       implicit real*8 (a-h,o-z)
10       include 'DIMENSIONS'
11
12
13 #ifdef MPI
14       include 'mpif.h'
15       include 'COMMON.SETUP'
16 #endif
17       include 'COMMON.TIME1'
18       include 'COMMON.INTERACT'
19       include 'COMMON.NAMES'
20       include 'COMMON.GEO'
21       include 'COMMON.HEADER'
22       include 'COMMON.CONTROL'
23       include 'COMMON.CONTACTS'
24       include 'COMMON.CHAIN'
25       include 'COMMON.VAR'
26       include 'COMMON.IOUNITS'
27       include 'COMMON.FFIELD'
28 c      include 'COMMON.REMD'
29 c      include 'COMMON.MD'
30       include 'COMMON.SBRIDGE'
31       double precision hrtime,mintime,sectime
32       character*64 text_mode_calc(-2:14) /'test',
33      & 'SC rotamer distribution',
34      & 'Energy evaluation or minimization',
35      & 'Regularization of PDB structure',
36      & 'Threading of a sequence on PDB structures',
37      & 'Monte Carlo (with minimization) ',
38      & 'Energy minimization of multiple conformations',
39      & 'Checking energy gradient',
40      & 'Entropic sampling Monte Carlo (with minimization)',
41      & 'Energy map',
42      & 'CSA calculations',
43      & 'Not used 9',
44      & 'Not used 10',
45      & 'Soft regularization of PDB structure',
46      & 'Mesoscopic molecular dynamics (MD) ',
47      & 'Not used 13',
48      & 'Replica exchange molecular dynamics (REMD)'/
49       external ilen
50
51 c      call memmon_print_usage()
52
53       call init_task
54       if (me.eq.king)
55      & write(iout,*)'### LAST MODIFIED  11/03/09 1:19PM by czarek'  
56       if (me.eq.king) call cinfo
57 C Read force field parameters and job setup data
58       call readrtns
59       call flush(iout)
60 C
61       if (me.eq.king .or. .not. out1file) then
62        write (iout,'(2a/)') 
63      & text_mode_calc(modecalc)(:ilen(text_mode_calc(modecalc))),
64      & ' calculation.' 
65        if (minim) write (iout,'(a)') 
66      &  'Conformations will be energy-minimized.'
67        write (iout,'(80(1h*)/)') 
68       endif
69       call flush(iout)
70 #ifdef MPI
71       if (fg_rank.gt.0) then
72 C Fine-grain slaves just do energy and gradient components.
73         call ergastulum ! slave workhouse in Latin
74       else
75 #endif
76       if (modecalc.eq.0) then
77         call exec_eeval_or_minim
78       else if (modecalc.eq.5) then
79          call exec_checkgrad
80       else
81         write (iout,'(a)') 'This calculation type is not supported',
82      &   ModeCalc
83       endif
84 #ifdef MPI
85       endif
86 C Finish task.
87       if (fg_rank.eq.0) call finish_task
88 c      call memmon_print_usage()
89 #ifdef TIMING
90        call print_detailed_timing
91 #endif
92       call MPI_Finalize(ierr)
93       stop 'Bye Bye...'
94 #else
95       call dajczas(tcpu(),hrtime,mintime,sectime)
96       stop '********** Program terminated normally.'
97 #endif
98       end
99 c---------------------------------------------------------------------------
100       subroutine exec_eeval_or_minim
101       implicit real*8 (a-h,o-z)
102       include 'DIMENSIONS'
103 #ifdef MPI
104       include 'mpif.h'
105 #endif
106       include 'COMMON.SETUP'
107       include 'COMMON.TIME1'
108       include 'COMMON.INTERACT'
109       include 'COMMON.NAMES'
110       include 'COMMON.GEO'
111       include 'COMMON.HEADER'
112       include 'COMMON.CONTROL'
113       include 'COMMON.CONTACTS'
114       include 'COMMON.CHAIN'
115       include 'COMMON.VAR'
116       include 'COMMON.IOUNITS'
117       include 'COMMON.FFIELD'
118       include 'COMMON.SBRIDGE'
119       common /srutu/ icall
120       double precision energy(0:n_ene),varia(maxvar)
121       double precision energy_long(0:n_ene),energy_short(0:n_ene)
122       if (indpdb.eq.0) call chainbuild
123 #ifdef MPI
124       time00=MPI_Wtime()
125 #else
126       time00=tcpu()
127 #endif
128       call chainbuild_cart
129       call etotal(energy(0))
130 #ifdef MPI
131       time_ene=MPI_Wtime()-time00
132 #else
133       time_ene=tcpu()
134 #endif
135       write (iout,*) "Time for energy evaluation",time_ene
136       print *,"after etotal"
137       etota = energy(0)
138       etot =etota
139       call enerprint(energy(0))
140 c      call hairpin(.true.,nharp,iharp)
141 c      call secondary2(.true.)
142       if (minim) then
143
144 crc overlap test
145         if (overlapsc) then
146           print *, 'Calling OVERLAP_SC'
147           call overlap_sc(fail)
148         endif
149
150         if (searchsc) then
151           call sc_move(2,nres-1,10,1d10,nft_sc,etot)
152           print *,'SC_move',nft_sc,etot
153           write(iout,*) 'SC_move',nft_sc,etot
154         endif
155
156         if (dccart) then
157           print *, 'Calling MINIM_DC'
158 #ifdef MPI
159           time1=MPI_WTIME()
160 #else
161           time1=tcpu()
162 #endif
163           call minim_dc(etot,iretcode,nfun)
164         else
165           if (indpdb.ne.0) then 
166             call bond_regular
167             call chainbuild
168           endif
169           call geom_to_var(nvar,varia)
170           print *,'Calling MINIMIZE.'
171 #ifdef MPI
172           time1=MPI_WTIME()
173 #else
174           time1=tcpu()
175 #endif
176           call minimize(etot,varia,iretcode,nfun)
177         endif
178         print *,'SUMSL return code is',iretcode,' eval ',nfun
179 #ifdef MPI
180         evals=nfun/(MPI_WTIME()-time1)
181 #else
182         evals=nfun/(tcpu()-time1)
183 #endif
184         print *,'# eval/s',evals
185         print *,'refstr=',refstr
186 c        call hairpin(.true.,nharp,iharp)
187 c        call secondary2(.true.)
188         call etotal(energy(0))
189         etot = energy(0)
190         call enerprint(energy(0))
191
192         call intout
193         call briefout(0,etot)
194 c        if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
195           write (iout,'(a,i3)') 'SUMSL return code:',iretcode
196           write (iout,'(a,i20)') '# of energy evaluations:',nfun+1
197           write (iout,'(a,f16.3)')'# of energy evaluations/sec:',evals
198 c      else
199 c        print *,'refstr=',refstr
200 c        if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
201 c        call briefout(0,etot)
202       endif
203       if (outpdb) call pdbout(etot,titel(:32),ipdb)
204       if (outmol2) call mol2out(etot,titel(:32))
205       return
206       end
207
208       subroutine exec_checkgrad
209       implicit real*8 (a-h,o-z)
210       include 'DIMENSIONS'
211 #ifdef MPI
212       include 'mpif.h'
213 #endif
214       include 'COMMON.SETUP'
215       include 'COMMON.TIME1'
216       include 'COMMON.INTERACT'
217       include 'COMMON.NAMES'
218       include 'COMMON.GEO'
219       include 'COMMON.HEADER'
220       include 'COMMON.CONTROL'
221       include 'COMMON.CONTACTS'
222       include 'COMMON.CHAIN'
223       include 'COMMON.VAR'
224       include 'COMMON.IOUNITS'
225       include 'COMMON.FFIELD'
226 c      include 'COMMON.REMD'
227       include 'COMMON.MD_'
228       include 'COMMON.SBRIDGE'
229       common /srutu/ icall
230       double precision energy(0:max_ene)
231 c      do i=2,nres
232 c        vbld(i)=vbld(i)+ran_number(-0.1d0,0.1d0)
233 c        if (itype(i).ne.10) 
234 c     &      vbld(i+nres)=vbld(i+nres)+ran_number(-0.001d0,0.001d0)
235 c      enddo
236       if (indpdb.eq.0) call chainbuild
237 c      do i=0,nres
238 c        do j=1,3
239 c          dc(j,i)=dc(j,i)+ran_number(-0.2d0,0.2d0)
240 c        enddo
241 c      enddo
242 c      do i=1,nres-1
243 c        if (itype(i).ne.10) then
244 c          do j=1,3
245 c            dc(j,i+nres)=dc(j,i+nres)+ran_number(-0.2d0,0.2d0)
246 c          enddo
247 c        endif
248 c      enddo
249 c      do j=1,3
250 c        dc(j,0)=ran_number(-0.2d0,0.2d0)
251 c      enddo
252       usampl=.true.
253       totT=1.d0
254       eq_time=0.0d0
255 c      call read_fragments
256       call chainbuild_cart
257       call cartprint
258       call intout
259       icall=1
260       call etotal(energy(0))
261       etot = energy(0)
262       call enerprint(energy(0))
263       write (iout,*) "Uconst",Uconst," Uconst_back",uconst_back
264       print *,'icheckgrad=',icheckgrad
265       goto (10,20,30) icheckgrad
266   10  call check_ecartint
267       return
268   20  call check_cartgrad
269       return
270   30  call check_eint
271       return
272       end