1 subroutine probabl(ib,nlist,ncon,*)
2 ! construct the conformational ensembles at REMD temperatures
9 integer ierror,errcode,status(MPI_STATUS_SIZE)
11 include "COMMON.CONTROL"
12 include "COMMON.IOUNITS"
14 include "COMMON.FFIELD"
15 include "COMMON.INTERACT"
16 include "COMMON.SBRIDGE"
17 include "COMMON.CHAIN"
18 include "COMMON.CLUSTER"
19 real*4 csingle(3,maxres2)
20 double precision fT(6),fTprim(6),fTbis(6),quot,quotl1,quotl,kfacl,
21 & eprim,ebis,temper,kfac/2.4d0/,T0/300.0d0/
22 double precision etot,evdw,evdw2,ees,evdw1,ebe,etors,escloc,
23 & ehpb,ecorr,ecorr5,ecorr6,eello_turn4,eello_turn3,
24 & eturn6,eel_loc,edihcnstr,etors_d,estr,evdw2_14,esccor,
26 integer i,ii,ik,iproc,iscor,j,k,l,ib,nlist,ncon
27 double precision qfree,sumprob,eini,efree,rmsdev
33 character*80 structure/'Structure'/
34 real*4 Fdimless(maxconf), Fdimless_buf(maxconf)
35 double precision energia(0:max_ene), totfree_buf(0:maxconf),
37 double precision buffer(maxconf)
42 c write (iout,*) i,list_conf(i)
45 write (iout,*) me," indstart",indstart(me)," indend",indend(me)
46 call daread_ccoords(indstart(me),indend(me))
48 C write (iout,*) "ncon",ncon
50 temper=1.0d0/(beta_h(ib)*1.987D-3)
51 c write (iout,*) "ib",ib," beta_h",beta_h(ib)," temper",temper
52 c quot=1.0d0/(T0*beta_h(ib)*1.987D-3)
59 c fT(l)=kfacl/(kfacl-1.0d0+quotl)
62 if (rescale_mode.eq.1) then
63 quot=1.0d0/(T0*beta_h(ib)*1.987D-3)
70 fT(l)=kfacl/(kfacl-1.0d0+quotl)
73 ft(6)=(320.0d0+80.0d0*dtanh((betaT-320.0d0)/80.0d0))/
75 ftprim(6)=1.0d0/(320.0d0*dcosh((betaT-320.0d0)/80.0d0)**2)
76 ftbis(6)=-2.0d0*dtanh((betaT-320.0d0)/80.0d0)
77 & /(320.0d0*80.0d0*dcosh((betaT-320.0d0)/80.0d0)**3)
88 else if (rescale_mode.eq.2) then
89 quot=1.0d0/(T0*beta_h(ib)*1.987D-3)
93 fT(l)=1.12692801104297249644d0/
94 & dlog(dexp(quotl)+dexp(-quotl))
96 c write (iout,*) 1.0d0/(beta_h(ib)*1.987D-3),ft
99 ft(6)=(320.0d0+80.0d0*dtanh((betaT-320.0d0)/80.0d0))/
101 ftprim(6)=1.0d0/(320.0d0*dcosh((betaT-320.0d0)/80.0d0)**2)
102 ftbis(6)=-2.0d0*dtanh((betaT-320.0d0)/80.0d0)
103 & /(320.0d0*80.0d0*dcosh((betaT-320.0d0)/80.0d0)**3)
122 C write (iout,*) "i",i," ii",ii,"ib",ib,scount(me)
127 c(k,j)=allcart(k,j,i)
128 c(k,j+nres)=allcart(k,j+nres,i)
129 C write(iout,*) "coord",i,j,k,allcart(k,j,i),c(k,j),
130 C & c(k,j+nres),allcart(k,j+nres,i)
133 C write(iout,*) "out of j loop"
137 c(k,nres+nres)=c(k,nres)
139 C write(iout,*) "after nres+nres",nss_all(i)
143 ihpb(j)=ihpb_all(j,i)
144 jhpb(j)=jhpb_all(j,i)
146 call int_from_cart1(.false.)
147 call etotal(energia(0),fT)
149 c write (structure(9:),'(bz,i6.6)') i
150 call TMscore_sub(rmsdev,gdt_ts_tb(i),
151 & gdt_ha_tb(i),tmscore_tb(i),Structure,.false.)
153 write (iout,*) i,rmsdev,gdt_ts_tb(i),gdt_ha_tb(i),
157 totfree(i)=energia(0)
158 totfree_buf(i)=totfree(i)
159 c write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
160 c write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
161 c call pdbout(totfree(i),16,i)
165 write (iout,*) "conformation", i
166 call enerprint(energia(0),fT)
170 enetb(k,i)=energia(k)
174 c write (iout,*) evdw
178 evdw2=enetb(2,i)+evdw2_14
194 eello_turn3=enetb(8,i)
195 eello_turn4=enetb(9,i)
204 edihcnstr=enetb(20,i)
206 c etot=wsc*(evdw+ft(6)*evdw_t)+wscp*evdw2+
207 c &ft(1)*welec*ees+wvdwpp*evdw1
208 c & +wang*ebe+ft(1)*wtor*etors+wscloc*escloc
209 c & +wstrain*ehpb+nss*ebr+ft(3)*wcorr*ecorr+ft(4)*wcorr5*ecorr5
210 c & +ft(5)*wcorr6*ecorr6+ft(3)*wturn4*eello_turn4
211 c & +ft(2)*wturn3*eello_turn3
212 c & +ft(5)*wturn6*eturn6+ft(2)*wel_loc*eel_loc
213 c & +edihcnstr+ft(2)*wtor_d*etors_d+ft(1)*wsccor*esccor
216 c etot=wsc*(evdw+ft(6)*evdw_t)+wscp*evdw2+ft(1)*welec*(ees+evdw1)
217 c & +wang*ebe+ft(1)*wtor*etors+wscloc*escloc
218 c & +wstrain*ehpb+nss*ebr+ft(3)*wcorr*ecorr+ft(4)*wcorr5*ecorr5
219 c & +ft(5)*wcorr6*ecorr6+ft(3)*wturn4*eello_turn4
220 c & +ft(2)*wturn3*eello_turn3
221 c & +ft(5)*wturn6*eturn6+ft(2)*wel_loc*eel_loc+edihcnstr
222 c & +ft(2)*wtor_d*etors_d+ft(1)*wsccor*esccor
226 write (iout,*) "etot2", etot
227 write (iout,*) "evdw", wsc, evdw,evdw_t
228 write (iout,*) "evdw2", wscp, evdw2
229 write (iout,*) "welec", ft(1),welec,ees
230 write (iout,*) "evdw1", wvdwpp,evdw1
231 write (iout,*) "ebe", ebe,wang
233 Fdimless(i)=beta_h(ib)*etot+entfac(ii)
234 Fdimless_buf(i)=Fdimless(i)
236 totfree_buf(i)=totfree(i)
238 write (iout,*) "fdim calc", i,ii,ib,
239 & 1.0d0/(1.987d-3*beta_h(ib)),totfree(i),
240 & entfac(ii),Fdimless(i)
245 entfac_buf(ijk)=entfac(ijk)
246 Fdimless_buf(ijk)=Fdimless(ijk)
249 totfree_buf(ijk)=totfree(ijk)
253 c scount_buf=scount(me)
254 c scount_buf2=scount(0)
256 c entfac_buf(indstart(me)+1)=entfac(indstart(me)+1)
259 c WRITE (iout,*) "Wchodze do call MPI_Gatherv1 (Propabl)"
260 call MPI_Gatherv(Fdimless_buf(1),scount(me),
261 & MPI_REAL,Fdimless(1),
262 & scount(0),idispl(0),MPI_REAL,Master,
263 & MPI_COMM_WORLD, IERROR)
264 c WRITE (iout,*) "Wchodze do call MPI_Gatherv2 (Propabl)"
265 call MPI_Gatherv(totfree_buf(1),scount(me),
266 & MPI_DOUBLE_PRECISION,totfree(1),
267 & scount(0),idispl(0),MPI_DOUBLE_PRECISION,Master,
268 & MPI_COMM_WORLD, IERROR)
269 c WRITE (iout,*) "Wchodze do call MPI_Gatherv3 (Propabl)"
270 call MPI_Gatherv(entfac_buf(indstart(me)+1),scount(me),
271 & MPI_DOUBLE_PRECISION,entfac(1),
272 & scount(0),idispl(0),MPI_DOUBLE_PRECISION,Master,
273 & MPI_COMM_WORLD, IERROR)
274 c WRITE (iout,*) "Wychodze z call MPI_Gatherv (Propabl)"
277 buffer(i)=gdt_ts_tb(i)
279 call MPI_Gatherv(buffer(1),scount(me),MPI_DOUBLE_PRECISION,
280 & gdt_ts_tb(1),scount(0),idispl(0),MPI_DOUBLE_PRECISION,Master,
281 & MPI_COMM_WORLD,IERROR)
283 buffer(i)=gdt_ha_tb(i)
285 call MPI_Gatherv(buffer(1),scount(me),MPI_DOUBLE_PRECISION,
286 & gdt_ha_tb(1),scount(0),idispl(0),MPI_DOUBLE_PRECISION,Master,
287 & MPI_COMM_WORLD,IERROR)
289 buffer(i)=tmscore_tb(i)
291 call MPI_Gatherv(buffer(1),scount(me),MPI_DOUBLE_PRECISION,
292 & tmscore_tb(1),scount(0),idispl(0),MPI_DOUBLE_PRECISION,Master,
293 & MPI_COMM_WORLD,IERROR)
295 if (me.eq.Master) then
296 c WRITE (iout,*) "me.eq.Master"
299 write (iout,*) "The FDIMLESS array before sorting"
301 write (iout,'(2i5,4f10.5)') i,list_conf(i),fdimless(i),
302 & gdt_ts_tb(i),gdt_ha_tb(i),tmscore_tb(i)
305 c WRITE (iout,*) "Wchodze do call mysort1"
306 call mysort1(ncon,Fdimless,list_conf)
307 c WRITE (iout,*) "Wychodze z call mysort1"
309 write (iout,*) "The FDIMLESS array after sorting"
311 write (iout,'(2i5,4f10.5)') i,list_conf(i),fdimless(i),
312 & gdt_ts_tb(i),gdt_ha_tb(i),tmscore_tb(i)
315 c WRITE (iout,*) "Wchodze do petli i=1,ncon totfree(i)=fdimless(i)"
317 totfree(i)=fdimless(i)
321 qfree=qfree+exp(-fdimless(i)+fdimless(1))
322 c write (iout,*) "fdimless", fdimless(i)
324 c write (iout,*) "qfree",qfree
327 write (iout,*) "ncon", ncon,maxstr_proc
328 do i=1,min0(ncon,maxstr_proc)-1
329 sumprob=sumprob+exp(-fdimless(i)+fdimless(1))/qfree
331 write (iout,*) i,ib,beta_h(ib),
332 & 1.0d0/(1.987d-3*beta_h(ib)),list_conf(i),
333 & totfree(list_conf(i)),
334 & -entfac(list_conf(i)),fdimless(i),sumprob
336 if (sumprob.gt.prob_limit) goto 122
337 c if (sumprob.gt.1.00d0) goto 122
343 call MPI_Bcast(nlist, 1, MPI_INTEGER, Master, MPI_COMM_WORLD,
345 call MPI_Bcast(list_conf,nlist,MPI_INTEGER,Master,MPI_COMM_WORLD,
348 c write (iout,*) "iproc",iproc," indstart",indstart(iproc),
349 c & " indend",indend(iproc)
351 write (iout,*) "nlist",nlist
355 !--------------------------------------------------
356 subroutine mysort1(n, x, ipermut)
358 integer i,j,imax,ipm,n
366 if (x(j).lt.xtemp) then
374 ipermut(imax)=ipermut(i)