1 subroutine make_list(lprn,first_call,nvarr,x)
4 include "DIMENSIONS.ZSCOPT"
7 integer IERROR,ErrCode,Status(MPI_STATUS_SIZE,10)
8 integer req(10),msg_in(5),msg_out(5),address,size
9 character*1 buffer(8*(2*maxstr_proc*nntyp+8000))
12 include "COMMON.WEIGHTS"
13 include "COMMON.WEIGHTDER"
14 include "COMMON.ENERGIES"
15 include "COMMON.IOUNITS"
16 include "COMMON.VMCPAR"
17 include "COMMON.NAMES"
18 include "COMMON.INTERACT"
19 include "COMMON.TIME1"
20 include "COMMON.CHAIN"
21 include "COMMON.PROTFILES"
24 include "COMMON.OPTIM"
25 include "COMMON.CLASSES"
26 include "COMMON.COMPAR"
27 include "COMMON.TORSION"
28 C Define local variables
30 integer i,ii,iii,ibatch,kkk,jj,j,k,kk,l,iprot,ib,nn
31 integer ipass_conf,istart_conf,iend_conf,Previous,Next
32 double precision energia(0:max_ene)
33 double precision etoti,elowesti,dene
34 double precision tcpu_ini,tcpu_fin,tcpu
35 double precision elowest_t(2,maxT),etot_aux,enepsjk,
36 & emini,elowest_aux(2,maxT)
39 logical lprn,first_call
40 logical*1 lflag(maxstr)
41 integer ipermut(maxstr)
42 integer list_conf_(maxstr,maxprot)
43 double precision x(nvarr)
44 double precision ftune_eps
47 double precision e_total_(maxstr_proc)
49 if (me.eq.Master) then
50 write (iout,*) "Making the worklist of conformations."
51 write (iout,*) "enecut",(enecut(i),i=1,nprot)
52 write (iout,*) "first_call ",first_call," nvarr",nvarr
59 ntot_work(iprot)=ntot(iprot)
63 c Divide the whole database between processors
67 if (Previous.lt.0) Previous = Nprocs-1
68 if (Next.ge.Nprocs) Next = 0
69 call work_partition(lprn)
76 call restore_molinfo(iprot)
77 c write (iout,*) "Processor",me," iprot",iprot,
78 c & " indstart",indstart(me1,iprot)," indend",indend(me1,iprot),
79 c & " init_ene",init_ene," mod_fourier",mod_fourier(0),
80 c & " mod_elec",mod_elec," mod_scp",mod_scp
81 call restore_molinfo(iprot)
83 c Loop over the conformations of protein IPROT assigned to the current processor.
84 c The conformations are read off a DA scratchfile and processed in passes, each
85 c of which requires no more than MAXSTR_PROC conformations to be stored in memory
88 if (.not.init_ene .and. mod_fourier(nloctyp).eq.0
89 & .and. .not. mod_elec .and. .not. mod_scp) then
90 open (ientin,file=benefiles(iprot),status="old",
91 & form="unformatted",access="direct",recl=lenrec_ene(iprot))
93 open (icbase,file=bprotfiles(iprot),status="old",
94 & form="unformatted",access="direct",recl=lenrec(iprot))
95 open (ientout,file=benefiles(iprot),status="old",
96 & form="unformatted",access="direct",recl=lenrec_ene(iprot))
100 nchunk_ene(iprot) = iroof(scount(me1,iprot),maxstr)
101 nchunk_conf(iprot) = iroof(scount(me1,iprot),maxstr_proc)
103 & write (iout,*)"Protein",iprot," energy evaluation in",
104 & nchunk_conf(iprot)," passes."
107 do i=indstart(me1,iprot),indend(me1,iprot),maxstr_proc
108 ipass_conf=ipass_conf+1
109 if (me.eq.Master) write (iout,*) "MAKE_LIST: Pass",ipass_conf
111 iend_conf=min0(i+maxstr_proc-1,indend(me1,iprot))
113 nchunk_ene(iprot) = iroof(ntot(iprot),maxstr)
114 nchunk_conf(iprot) = iroof(ntot(iprot),maxstr_proc)
116 & write (iout,*)"Protein",iprot," energy evaluation in",
117 & nchunk_conf," passes."
119 do i=1,ntot(iprot),maxstr_proc
120 ipass_conf=ipass_conf+1
121 if (me.eq.Master) write (iout,*) "MAKE_LIST: Pass",ipass_conf
123 iend_conf=min0(i+maxstr_proc-1,ntot(iprot))
127 c Read the chunk of conformations off a DA scratchfile.
129 if (.not.init_ene .and. mod_fourier(nloctyp).eq.0
130 & .and. .not. mod_elec .and. .not. mod_scp) then
132 c If energy components have been pre-computed read them off a DA file.
134 call daread_ene(iprot,istart_conf,iend_conf)
135 do iii=istart_conf,iend_conf
139 enetb(ii,1,iprot)=0.0d0
142 enetb(ii,1,iprot)=enetb(ii,1,iprot)+ftune_eps(eps(j,k))*
143 & eneps(1,icant(j,k),ii,iprot)+
144 & eps(j,k)*eneps(2,icant(j,k),ii,iprot)
150 etot_aux=etot_aux+ww(j)*enetb(ii,j,iprot)
152 e_total(iii,iprot)=etot_aux
154 write (iout,'(i5,16(1pe12.4))') iii,
155 & (enetb(ii,j,iprot),j=1,n_ene),e_total(iii,iprot)
159 if (first_call .and. mod_side) then
160 c write (iout,*) "Callig x2w"
162 call x2w(nvarr,x_orig)
163 c write (iout,*) "After x2w"
166 do iii=istart_conf,iend_conf
168 enetb_oorig(ii,1,iprot)=0.0d0
171 enetb_oorig(ii,1,iprot)=enetb_oorig(ii,1,iprot)+
172 & ftune_eps(eps(j,k))*eneps(1,icant(j,k),ii,iprot)+
173 & eps(j,k)*eneps(2,icant(j,k),ii,iprot)
177 enetb_oorig(ii,j,iprot)=enetb(ii,j,iprot)
183 call daread_ccoords(iprot,istart_conf,iend_conf)
185 c Compute the energies of the conformations currently in memory and compute
186 c the lowest energies.
188 do iii=istart_conf,iend_conf
191 call restore_ccoords(iprot,ii)
192 call int_from_cart1(.false.)
194 write (iout,*) "Before etotal",iii,i
197 call etotal(energia(0))
199 write (iout,*) "After etotal",i
201 call enerprint(energia(0))
202 c write (iout,'(i5,16(1pe12.4))') i,
203 c & (energia(j),j=1,n_ene),energia(0)
207 write (iout,*) "Conformation:",i
208 write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
209 write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
210 write (iout,'(8f10.4)') (vbld(k),k=2,nres)
211 write (iout,'(8f10.4)') (vbld(k+nres),k=nnt,nct)
212 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
213 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
214 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
215 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
216 call enerprint(energia(0))
218 e_total(iii,iprot)=energia(0)
220 enetb(ii,j,iprot)=energia(j)
224 eneps(1,icant(j,k),ii,iprot)=eneps_temp(1,icant(j,k))
225 eneps(2,icant(j,k),ii,iprot)=eneps_temp(2,icant(j,k))
229 write (iout,'(i5,20(1pe12.4))') iii,
230 & (energia(j),j=1,n_ene),energia(0),eini(iii,iprot),
234 if (energia(0).ge.1.0d99) then
235 if (me.eq.Master) then
237 & "MAKE_LIST:CHUJ NASTAPIL in energy evaluation for",
238 & " point",i,". Probably NaNs in some of the energy components."
239 write (iout,*) "The components of the energy are:"
240 call enerprint(energia(0))
241 write (iout,*) "Conformation:",i
242 write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
243 write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
244 write (iout,'(8f10.4)') (vbld(k),k=2,nres)
245 write (iout,'(8f10.4)') (vbld(k+nres),k=nnt,nct)
246 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
247 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
248 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
249 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
250 write (iout,*) "Calculation terminated at this point.",
251 & " Check the database of conformations"
255 call MPI_Abort(ALL_COMM,ErrCode,IERROR)
257 stop "SEVERE error in energy calculation"
261 do ii=1,iend_conf-istart_conf+1
263 enetb_orig(ii,j,iprot)=enetb(ii,j,iprot)
267 call x2w(nvarr,x_orig)
269 write (iout,*) "x,xorig"
271 write (iout,'(i5,2f10.5)') k,x(k),x_orig(k)
274 do iii=istart_conf,iend_conf
276 call restore_ccoords(iprot,ii)
277 call int_from_cart1(.false.)
278 call etotal(energia(0))
280 write (iout,*) "Conformation:",iii,ii
281 write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
282 write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
283 write (iout,'(8f10.4)') (vbld(k),k=2,nres)
284 write (iout,'(8f10.4)') (vbld(k+nres),k=nnt,nct)
285 write (iout,'(8f10.4)') (rad2deg*theta(k),k=3,nres)
286 write (iout,'(8f10.4)') (rad2deg*phi(k),k=4,nres)
287 write (iout,'(8f10.4)') (rad2deg*alph(k),k=2,nres-1)
288 write (iout,'(8f10.4)') (rad2deg*omeg(k),k=2,nres-1)
289 call enerprint(energia(0))
292 enetb_oorig(ii,j,iprot)=energia(j)
295 write (iout,'(2i5,20(1pe12.4))') iii,ii,
296 & (energia(j),j=1,n_ene),energia(0)
297 write (iout,'(2i5,20(1pe12.4))') iii,ii,
298 & (enetb(ii,j,iprot),j=1,n_ene)
301 if (energia(0).ge.1.0d99) then
302 if (me.eq.Master) then
303 write (iout,*) "CHUJ NASTAPIL in energy evaluation for",
304 & " point",i,". Probably NaNs in some of the energy components."
305 write (iout,*) "The components of the energy are:"
306 call enerprint(energia(0))
307 write (iout,*) "Calculation terminated at this point.",
308 & " Check the database of conformations"
312 call MPI_Abort(ALL_COMM,ErrCode,IERROR)
314 stop "SEVERE error in energy calculation"
317 c write (iout,*) "MAKE_LIST Callig x2w"
320 c write (iout,*) "After x2w"
323 c write (iout,*) "make_list: calling dawrite_ene"
324 c write (iout,*) "istart_conf",istart_conf,
325 c & " iend_conf",iend_conf
326 call dawrite_ene(iprot,istart_conf,iend_conf,ientout)
328 c Distribute energy components through ring
329 call MPI_Barrier(WHAM_COMM,IERROR)
330 c write (iout,*) "Processes synchronized in make_list"
332 msg_out(1)= 5*me1+1000*ipass_conf+1
333 msg_out(2)= 5*me1+1000*ipass_conf+2
334 msg_out(3)= 5*me1+1000*ipass_conf+3
335 msg_out(4)= 5*me1+1000*ipass_conf+4
336 msg_out(5)= 5*me1+1000*ipass_conf+5
338 c Send the current energy tables to the right neighbor
339 c Receive the energy tables produced by processor kkk from the left neighbor
340 kkk = mod(me1-iii+NProcs,Nprocs)
341 msg_in(1)= 5*kkk+1000*ipass_conf+1
342 msg_in(2)= 5*kkk+1000*ipass_conf+2
343 msg_in(3)= 5*kkk+1000*ipass_conf+3
344 msg_in(4)= 5*kkk+1000*ipass_conf+4
345 msg_in(5)= 5*kkk+1000*ipass_conf+5
346 if (me.eq.Master) then
347 write (iout,*) "me1",me1," iii",iii," Previous",Previous,
348 & " Next",Next," kkk",kkk
349 write (iout,*) "msg_in",msg_in
350 write (iout,*) "msg_out",msg_out
352 write (iout,*) "Processor",me1," Start Send and receive"
355 call MPI_Send(istart_conf,1,MPI_INTEGER,Next,msg_out(1),
357 if (me.eq.Master)write (iout,*) "Send",msg_out(1)," complete"
358 if (me.eq.Master)call flush(iout)
359 call MPI_Recv(istart_conf,1,MPI_INTEGER,Previous,
360 & msg_in(1),WHAM_COMM,STATUS(1,6),IERROR)
361 if (me.eq.Master)write (iout,*) "Recv",msg_in(1)," complete"
362 if (me.eq.Master)call flush(iout)
363 call MPI_Send(iend_conf,1,MPI_INTEGER,Next,msg_out(2),
365 if (me.eq.Master)write (iout,*) "Send",msg_out(2)," complete"
366 if (me.eq.Master)call flush(iout)
367 call MPI_Recv(iend_conf,1,MPI_INTEGER,Previous,
368 & msg_in(2),WHAM_COMM,STATUS(1,7),IERROR)
369 if (me.eq.Master)write (iout,*) "Recv",msg_in(2)," complete"
370 if (me.eq.Master)call flush(iout)
371 call MPI_Buffer_Attach(buffer(1),8*(2*maxstr_proc*n_ene+800),
373 call MPI_BSend(enetb(1,1,iprot),maxstr_proc*n_ene,
374 & MPI_DOUBLE_PRECISION,Next,msg_out(3),
376 if (me.eq.Master) then
377 write (iout,*) "Send",msg_out(3)," complete (enetb)"
380 call MPI_Recv(enetb(1,1,iprot),maxstr_proc*n_ene,
381 & MPI_DOUBLE_PRECISION,Previous,msg_in(3),WHAM_COMM,
382 & STATUS(1,8),IERROR)
383 if (me.eq.Master) then
384 write (iout,*) "Recv",msg_in(3)," complete (enetb)"
387 call MPI_Buffer_Detach(address,size,IERROR)
388 c write (iout,*) "MPI_Buffer_Detach complete (enetb)"
390 call MPI_Buffer_Attach(buffer(1),8*(2*maxstr_proc*nntyp+800),
392 c write (iout,*) "MPI_Buffer_Attach complete (eneps)"
394 call MPI_BSend(eneps(1,1,1,iprot),2*maxstr_proc*nntyp,
395 & MPI_DOUBLE_PRECISION,Next,msg_out(4),
397 if (me.eq.Master) then
398 write (iout,*) "Send",msg_out(4)," complete (eneps)"
401 call MPI_Recv(eneps(1,1,1,iprot),2*maxstr_proc*nntyp,
402 & MPI_DOUBLE_PRECISION,Previous,msg_in(4),WHAM_COMM,
403 & STATUS(1,9),IERROR)
404 if (me.eq.Master) then
405 write (iout,*) "Recv",msg_in(4)," complete (eneps)"
408 call MPI_Buffer_Detach(address,size,IERROR)
409 call MPI_Buffer_Attach(buffer(1),8*(2*maxstr_proc*maxnatlike
411 call MPI_BSend(nu(1,1,1,iprot),
412 & maxstr_proc*maxnatlike*maxdimnat,
413 & MPI_DOUBLE_PRECISION,Next,msg_out(5),
415 if (me.eq.Master) then
416 write (iout,*) "Send",msg_out(5)," complete (nu)"
419 call MPI_Recv(nu(1,1,1,iprot),
420 & maxstr_proc*maxnatlike*maxdimnat,
421 & MPI_DOUBLE_PRECISION,Previous,msg_in(5),WHAM_COMM,
422 & STATUS(1,10),IERROR)
423 if (me.eq.Master) then
424 write (iout,*) "Recv",msg_in(5)," complete (nu)"
427 call MPI_Buffer_Detach(address,size,IERROR)
428 if (me.eq.Master) then
429 write (iout,*) "Send and receive complete"
431 write (iout,*) "Processor",me1," calling dawrite_ene",
432 & " istart_conf",istart_conf," iend_conf",iend_conf
436 c write (iout,*) "assignment of enetb_orig"
437 do ii=1,iend_conf-istart_conf+1
439 enetb_orig(ii,j,iprot)=enetb(ii,j,iprot)
441 c write (iout,'(i5,20f8.2)')
442 c & ii,(enetb_orig(ii,j,iprot),j=1,n_ene)
445 call dawrite_ene(iprot,istart_conf,iend_conf,ientout)
453 if (.not.init_ene .and. mod_fourier(nloctyp).eq.0
454 & .and. .not. mod_elec .and. .not. mod_scp) then
463 c Lowest free energies of structural classes
467 IF (NCHUNK_CONF(IPROT).EQ.1) THEN
470 do i=1,ntot_work(iprot)
474 do i=indstart(me1,iprot),indend(me1,iprot)
478 istart_conf=indstart(me1,iprot)
479 iend_conf=indend(me1,iprot)
481 do i=1,ntot_work(iprot)
485 iend_conf=ntot_work(iprot)
488 write (iout,*) "i2ii at make_list"
489 do i=1,ntot_work(iprot)
490 write (iout,*) "i",i," i2ii",i2ii(i,iprot)
494 open (ientin,file=benefiles(iprot),status="old",
495 & form="unformatted",access="direct",recl=lenrec_ene(iprot))
496 call daread_ene(iprot,istart_conf,iend_conf)
497 call emin_search(iprot)
501 open (ientin,file=benefiles(iprot),status="old",
502 & form="unformatted",access="direct",recl=lenrec_ene(iprot))
505 do istart_conf=indstart(me1,iprot),indend(me1,iprot),maxstr_proc
506 iend_conf=min0(istart_conf+maxstr_proc-1,indend(me1,iprot))
508 do istart_conf=1,ntot_work(iprot),maxstr_proc
509 iend_conf=min0(istart_conf+maxstr_proc-1,ntot_work(iprot))
512 c Read the chunk of energies and derivatives off a DA scratchfile.
514 ipass_conf=ipass_conf+1
515 do i=1,ntot_work(iprot)
519 do i=istart_conf,iend_conf
524 write (iout,*) "ipass_conf",ipass_conf,
525 & " istart_conf",istart_conf," iend_conf",iend_conf
526 do i=1,ntot_work(iprot)
527 write (iout,*) "i",i," i2ii",i2ii(i,iprot)
531 call daread_ene(iprot,istart_conf,iend_conf)
532 call emin_search(iprot)
540 c Complete the calculation of the lowest energies over all classes and
541 c distribute the values to all procs
543 do ibatch=1,natlike(iprot)+2
545 do ib=1,nbeta(ibatch,iprot)
546 write (iout,'(7hELOWEST,3i3,f15.3,i12)') iprot,ibatch,ib,
547 & elowest(ib,ibatch,iprot),ind_lowest(ib,ibatch,iprot)
550 do ib=1,nbeta(ibatch,iprot)
551 elowest_aux(1,ib)=elowest(ib,ibatch,iprot)
552 elowest_aux(2,ib)=ind_lowest(ib,ibatch,iprot)
554 call MPI_Allreduce(elowest_aux(1,1),elowest_t(1,1),
555 & nbeta(ibatch,iprot),
556 & MPI_2DOUBLE_PRECISION, MPI_MINLOC, Comm1, IERROR)
558 do ib=1,nbeta(ibatch,iprot)
559 write (iout,*) "beta=",betaT(ib,ibatch,iprot)
560 write (iout,'(9helowest_t,10f15.3)')
561 & elowest_t(1,ib),elowest_t(2,ib)
563 write (iout,*) "Processor",me,me1," finished MPI_Reduce: 2"
565 do ib=1,nbeta(ibatch,iprot)
566 elowest(ib,ibatch,iprot)=elowest_t(1,ib)
567 ind_lowest(ib,ibatch,iprot)=elowest_t(2,ib)
571 if (me.eq.Master) then
573 do ibatch=1,natlike(iprot)+2
574 do ib=1,nbeta(ibatch,iprot)
575 write (iout,*) "iprot",iprot," ibatch",ibatch," ib",ib,
576 & " elowest",elowest(ib,ibatch,iprot)
582 c Allgather to provide all energies to all processors
585 do i=1,scount(me1,iprot)
586 e_total_(i)=e_total(indstart(me1,iprot)+i-1,iprot)
588 call MPI_Allgatherv(e_total_(1),
589 & scount(me1,iprot),MPI_DOUBLE_PRECISION,e_total(1,iprot),
590 & scount(0,iprot),idispl(0,iprot),MPI_DOUBLE_PRECISION,
592 c call MPI_Allgatherv(e_total(indstart(me1,iprot),iprot),
593 c & scount(me1,iprot),MPI_DOUBLE_PRECISION,e_total(1,iprot),
594 c & scount(0,iprot),idispl(0,iprot),MPI_DOUBLE_PRECISION,
599 c Now determine which conformations will enter the database.
602 call restore_molinfo(iprot)
603 c Clear the list of conformations
604 do i=1,min0(ntot(iprot),maxstr)
610 c Make the list of conformations based on energy cut-off.
613 write (iout,*) "iprot",iprot," ibatch",ibatch," betmin",
614 & betmin(ibatch,iprot)
617 write (iout,*) "e_lowb",e_lowb(iprot)
618 write (iout,*) "t_lowb",t_lowb(iprot)
620 do ibatch=1,natlike(iprot)+2
624 c write (iout,*) "i",i," ii",ii," indstart",
625 c & indstart(me1,iprot)," indend",indend(me1,iprot)
626 if (i.ge.indstart(me1,iprot).and.i.le.indend(me1,iprot))
629 c write (iout,*) "i",i," ii",ii," kbatch",kbatch(i,iprot),
631 if (.not.lflag(i)) then
632 c if (eini(i,iprot).lt.e_lowb(iprot) .or.
633 c & entfac(i,iprot).lt.t_lowb(iprot)) then
635 c write (iout,*) "Conformation",i," eini",eini(i,iprot),
636 c & " entfac",entfac(i,iprot)," e_lowb",
638 c & " t_lowb",t_lowb(iprot)
643 do ib=1,nbeta(ibatch,iprot)
644 dene=betaT(ib,ibatch,iprot)*(e_total(i,iprot)
645 & -elowest(ib,ibatch,iprot))+entfac(i,iprot)
647 write (iout,*) "beta",betaT(ib,ibatch,iprot),
650 & " elowest",elowest(ib,ibatch,iprot)," dene",dene,
651 & " enecut",enecut(iprot)
653 if (dene.lt.enecut(iprot)) then
655 list_conf(nn,iprot)=i
668 if (nn.gt.maxstr) then
669 if (me.eq.Master) then
670 write (iout,*) "Error - after applying cutoff the number",
671 & " of conformations for protein ",i," exceeds MAXSTR:",
673 write (iout,*) "The calculation is terminating."
677 call MPI_Finalize(ierror)
681 call imysort(ntot_work(iprot),list_conf(1,iprot),ipermut)
683 & write (iout,*) "Protein",iprot,ntot_work(iprot),
684 & " conformations within scaled energy cut-off=",enecut(iprot),
685 & " found at processor",me
688 c All workers get the complete list of conformations.
690 call MPI_Allgather(ntot_work(iprot),1,MPI_INTEGER,
691 & scount(0,iprot),1,MPI_INTEGER,Comm1,IERROR)
694 idispl(i,iprot)=idispl(i-1,iprot)+scount(i-1,iprot)
697 write (iout,*) "Protein",iprot," Scount and Idispl"
699 write (iout,*) i,scount(i,iprot),idispl(i,iprot)
701 write (iout,*) "Protein",i,
702 & " local list of conformations of processor",me
703 do i=1,ntot_work(iprot)
704 write(iout,*) i,list_conf(i,iprot)
706 write (iout,*) "Before REDUCE: ntot_work",ntot_work(iprot)
709 call MPI_Allreduce(ntot_work(iprot),nn,1,MPI_INTEGER,MPI_SUM,
713 write (iout,*) "After REDUCE: ntot_work",ntot_work(iprot)
716 call MPI_Allgatherv(list_conf(1,iprot),
717 & scount(me1,iprot),MPI_INTEGER,list_conf_(1,iprot),
718 & scount(0,iprot),idispl(0,iprot),MPI_INTEGER,Comm1,IERROR)
719 do i=1,ntot_work(iprot)
720 list_conf(i,iprot)=list_conf_(i,iprot)
723 write (iout,*) "Protein",i,
724 & " global list of conformations of processor",me
725 do i=1,ntot_work(iprot)
726 write(iout,'(2i5,e15.5,33f7.3)')i,list_conf(i,iprot),
727 & e_total(list_conf(i,iprot),iprot)
728 c & ,(nu(k,i,iprot),k=1,natlike(iprot))
734 c Construct the mapping of the new list to the original numbers of
740 do i=1,ntot_work(iprot)
741 tsil(list_conf(i,iprot),iprot)=i
744 write (iout,*) "Protein",i," List-to-conformation mapping"
746 write(iout,*) i,tsil(i,iprot)
751 c Divide the work again based on the current lists of conformations
753 if (me.eq.Master) then
754 call work_partition(.true.)
756 call work_partition(.false.)
759 c If the conformations fit into memory, read them off a DA scratchfile.
762 call restore_molinfo(iprot)
764 nchunk_conf(iprot)=iroof(scount(me1,iprot),maxstr_proc)
766 nchunk_conf(iprot)=iroof(ntot_work(iprot),maxstr_proc)
768 if (nchunk_conf(iprot).eq.1) then
769 if (me.eq.Master) write (iout,*) "Protein",iprot,
770 & " in-memory storage of conformations."
771 if (init_ene .or. mod_fourier(nloctyp).gt.0
772 & .or. mod_elec .or. mod_scp) then
773 open (icbase,file=bprotfiles(iprot),status="old",
774 & form="unformatted",access="direct",recl=lenrec(iprot))
776 call daread_ccoords(iprot,indstart(me1,iprot),
779 call daread_ccoords(iprot,1,ntot_work(iprot))
783 open (ientin,file=benefiles(iprot),status="old",
784 & form="unformatted",access="direct",recl=lenrec_ene(iprot))
786 call daread_ene(iprot,indstart(me1,iprot),
789 call daread_ene(iprot,1,ntot_work(iprot))
792 open (icbase,file=bprotfiles(iprot),status="old",
793 & form="unformatted",access="direct",recl=lenrec(iprot))
795 call daread_ccoords(iprot,indstart(me1,iprot),
798 call daread_ccoords(iprot,1,ntot_work(iprot))
803 write (iout,*) "Protein",i,
804 & " global list of conformations of processor",me
805 do i=1,ntot_work(iprot)
806 write(iout,'(2i5,e15.5,33f7.3)')i,list_conf(i,iprot),
807 & e_total(list_conf(i,iprot),iprot)
811 if (me.eq.Master) write (iout,*) "Protein",iprot,
812 & " off-memory storage of conformations; ",
813 & "energy will be evaluated in",nchunk_conf(iprot)," passes."
821 write (iout,*) "E_TOTAL and ETOT_ORIG of protein",iprot
823 do i=1,ntot_work(iprot)
824 etot_orig(i,iprot)=e_total(list_conf(i,iprot),iprot)
826 write (iout,*) i,list_conf(i,iprot),
827 & e_total(list_conf(i,iprot),iprot),etot_orig(i,iprot)
830 do ibatch=1,natlike(iprot)+2
831 do ib=1,nbeta(ibatch,iprot)
832 emini = elowest(ib,ibatch,iprot)
833 if (elowest(ib,ibatch,iprot) .lt. emini)
834 & emini = elowest(ib,ibatch,iprot)
835 emin_orig(ib,ibatch,iprot)=emini
840 if (me1 .eq. Master) call write_conf_count
842 call write_conf_count
845 write (iout,*) "ELOWEST at the end of MAKE_LIST"
847 do ibatch=1,natlike(iprot)+2
848 do ib=1,nbeta(ibatch,iprot)
849 write (iout,*) "iprot",iprot," ibatch",ibatch," ib",ib,
850 & " elowest",elowest(ib,ibatch,iprot)
855 tcpu_fin = tcpu() - tcpu_ini
856 if (me.eq.Master) then
857 write (iout,*) "Time for creating list of conformations",tcpu_fin
860 t_func = t_func + tcpu_fin
863 c----------------------------------------------------------------------------
864 subroutine emin_search(iprot)
867 include "DIMENSIONS.ZSCOPT"
870 integer IERROR,ErrCode,Status(MPI_STATUS_SIZE,10)
871 integer req(10),msg_in(5),msg_out(5),address,size
872 character*1 buffer(8*(2*maxstr_proc*nntyp+8000))
875 include "COMMON.WEIGHTS"
876 include "COMMON.WEIGHTDER"
877 include "COMMON.COMPAR"
878 include "COMMON.ENERGIES"
879 include "COMMON.IOUNITS"
880 include "COMMON.VMCPAR"
881 include "COMMON.NAMES"
882 include "COMMON.INTERACT"
883 include "COMMON.TIME1"
884 include "COMMON.CHAIN"
885 include "COMMON.PROTFILES"
888 include "COMMON.CLASSES"
889 C Define local variables
890 integer i,ii,iii,kkk,jj,j,k,kk,l,iprot,ib,ibatch,nn
891 integer ipass_conf,istart_conf,iend_conf,Previous,Next
892 double precision energia(0:max_ene)
893 double precision etoti,elowesti,dene
894 double precision tcpu_ini,tcpu_fin,tcpu
895 double precision etot_aux,enepsjk,
896 & emini,elowest_aux(2,maxT)
901 do ibatch=1,natlike(iprot)+2
902 do ib=1,nbeta(ibatch,iprot)
903 elowest(ib,ibatch,iprot)=1.0d20
906 do ibatch=1,natlike(iprot)+2
912 do ib=1,nbeta(ibatch,iprot)
915 etoti=etoti+ww(kk)*enetb(jj,kk,iprot)
916 & *escal(kk,ib,ibatch,iprot)
918 if (ib.eq.1 .and. etoti.lt.elowest_ent(1,ibatch,iprot))
920 elowest_ent(1,ibatch,iprot)=etoti
921 elowest_ent(2,ibatch,iprot)=entfac(k,iprot)
924 c & etoti=etoti+entfac(k,iprot)/betaT(ib,ibatch,iprot)
925 etoti=etoti+entfac(k,iprot)/betaT(ib,ibatch,iprot)
926 if (etoti.lt.elowest(ib,ibatch,iprot)) then
927 elowest(ib,ibatch,iprot)=etoti
928 ind_lowest(ib,ibatch,iprot)=k
930 c write (iout,*) ib,betaT(ib,ibatch,iprot),etoti,
933 write (iout,'(2i5,2e15.5,f8.3)') k,jj,
934 & etoti,elowest(ib,ibatch,iprot),betaT(ib,ibatch,iprot)
938 write (iout,'(2i5,20f8.2)') j,jj,(enetb(jj,kk,iprot),