2 !-----------------------------------------------------------------------------
5 !-----------------------------------------------------------------------------
8 !-----------------------------------------------------------------------------
10 !----------------------------------------------------------------------------
12 !----------------------------------------------------------------------------
13 subroutine probabl(ib,nlist,ncon,*)
14 ! construct the conformational ensembles at REMD temperatures
16 ! include "DIMENSIONS"
17 ! include "sizesclu.dat"
19 use io_clust, only: daread_ccoords
20 use geometry_data, only: nres,c
21 use energy_data, only: nss,ihpb,jhpb,max_ene
22 use geometry, only: int_from_cart1
23 use energy, only: etotal,rescale_weights
24 use energy, only: rescale_mode,enerprint,weights
28 ! include "COMMON.MPI"
29 integer :: ierror,errcode !,status(MPI_STATUS_SIZE)
31 ! include "COMMON.IOUNITS"
32 ! include "COMMON.FREE"
33 ! include "COMMON.FFIELD"
34 ! include "COMMON.INTERACT"
35 ! include "COMMON.SBRIDGE"
36 ! include "COMMON.CHAIN"
37 ! include "COMMON.CLUSTER"
38 real(kind=4) :: csingle(3,2*nres)
39 real(kind=8) :: fT(6),fTprim(6),fTbis(6),quot,quotl1,quotl,kfacl,&
40 eprim,ebis,temper,kfac=2.4d0,T0=300.0d0
41 real(kind=8) :: etot,evdw,evdw2,ees,evdw1,ebe,etors,escloc,&
42 ehpb,ecorr,ecorr5,ecorr6,eello_turn4,eello_turn3,&
43 eturn6,eel_loc,edihcnstr,etors_d,estr,evdw2_14,esccor,&
45 integer :: i,ii,ik,iproc,iscor,j,k,l,ib,nlist,ncon
46 real(kind=8) :: qfree,sumprob,eini,efree,rmsdev
47 character(len=80) :: bxname
48 character(len=2) :: licz1
49 character(len=5) :: ctemper
52 real(kind=4) :: Fdimless(maxconf),Fdimless_(maxconf)
53 real(kind=8) :: totfree_(0:maxconf),entfac_(maxconf)
54 real(kind=8) :: energia(0:max_ene)
55 integer,dimension(0:nprocs) :: scount_
61 ! write (iout,*) i,list_conf(i)
64 write (iout,*) me," indstart",indstart(me)," indend",indend(me)
65 call daread_ccoords(indstart(me),indend(me))
67 ! write (iout,*) "ncon",ncon
68 temper=1.0d0/(beta_h(ib)*1.987D-3)
69 !elwrite(iout,*)"rescale_mode", rescale_mode
70 ! write (iout,*) "ib",ib," beta_h",beta_h(ib)," temper",temper
71 ! quot=1.0d0/(T0*beta_h(ib)*1.987D-3)
78 ! fT(l)=kfacl/(kfacl-1.0d0+quotl)
80 ! EL start old rescale-------------------------------
81 ! if (rescale_mode.eq.1) then
82 ! quot=1.0d0/(T0*beta_h(ib)*1.987D-3)
89 ! fT(l)=kfacl/(kfacl-1.0d0+quotl)
92 ! ft(6)=(320.0d0+80.0d0*dtanh((betaT-320.0d0)/80.0d0))/ &
94 ! ftprim(6)=1.0d0/(320.0d0*dcosh((betaT-320.0d0)/80.0d0)**2)
95 ! ftbis(6)=-2.0d0*dtanh((betaT-320.0d0)/80.0d0) &
96 ! /(320.0d0*80.0d0*dcosh((betaT-320.0d0)/80.0d0)**3)
107 ! else if (rescale_mode.eq.2) then
108 ! quot=1.0d0/(T0*beta_h(ib)*1.987D-3)
112 ! fT(l)=1.12692801104297249644d0/ &
113 ! dlog(dexp(quotl)+dexp(-quotl))
115 !el write (iout,*) 1.0d0/(beta_h(ib)*1.987D-3)
118 ! ft(6)=(320.0d0+80.0d0*dtanh((betaT-320.0d0)/80.0d0))/ &
120 ! ftprim(6)=1.0d0/(320.0d0*dcosh((betaT-320.0d0)/80.0d0)**2)
121 ! ftbis(6)=-2.0d0*dtanh((betaT-320.0d0)/80.0d0) &
122 ! /(320.0d0*80.0d0*dcosh((betaT-320.0d0)/80.0d0)**3)
133 !EL end old rescele----------------------
135 call rescale_weights(1.0d0/(beta_h(ib)*1.987D-3))
144 ! write (iout,*) "i",i," ii",ii
149 c(k,j)=allcart(k,j,i)
150 c(k,j+nres)=allcart(k,j+nres,i)
155 c(k,nres+nres)=c(k,nres)
159 !write(iout,*)"c, k, j",k,j,c(k,j)
164 ihpb(j)=ihpb_all(j,i)
165 jhpb(j)=jhpb_all(j,i)
167 call int_from_cart1(.false.)
168 ! call etotal(energia(0),fT)
170 totfree(i)=energia(0)
171 ! write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres)
172 ! write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct)
173 ! call enerprint(energia)
174 ! call pdbout(totfree(i),16,i)
176 write (iout,*) i," energia",(energia(j),j=0,21)
177 write (iout,*) "etot", etot
178 ! write (iout,*) "ft(6)", ft(6)
181 enetb(k,i)=energia(k)
185 ! write (iout,*) evdw
188 !el evdw2_14=enetb(17,i)
190 evdw2=enetb(2,i)+evdw2_14
206 eello_turn3=enetb(8,i)
207 eello_turn4=enetb(9,i)
218 ! edihcnstr=enetb(20,i)
219 edihcnstr=enetb(19,i)
221 ! etot=wsc*(evdw+ft(6)*evdw_t)+wscp*evdw2+
222 ! &ft(1)*welec*ees+wvdwpp*evdw1
223 ! & +wang*ebe+ft(1)*wtor*etors+wscloc*escloc
224 ! & +wstrain*ehpb+nss*ebr+ft(3)*wcorr*ecorr+ft(4)*wcorr5*ecorr5
225 ! & +ft(5)*wcorr6*ecorr6+ft(3)*wturn4*eello_turn4
226 ! & +ft(2)*wturn3*eello_turn3
227 ! & +ft(5)*wturn6*eturn6+ft(2)*wel_loc*eel_loc
228 ! & +edihcnstr+ft(2)*wtor_d*etors_d+ft(1)*wsccor*esccor
231 ! etot=wsc*(evdw+ft(6)*evdw_t)+wscp*evdw2+ft(1)*welec*(ees+evdw1)
232 ! & +wang*ebe+ft(1)*wtor*etors+wscloc*escloc
233 ! & +wstrain*ehpb+nss*ebr+ft(3)*wcorr*ecorr+ft(4)*wcorr5*ecorr5
234 ! & +ft(5)*wcorr6*ecorr6+ft(3)*wturn4*eello_turn4
235 ! & +ft(2)*wturn3*eello_turn3
236 ! & +ft(5)*wturn6*eturn6+ft(2)*wel_loc*eel_loc+edihcnstr
237 ! & +ft(2)*wtor_d*etors_d+ft(1)*wsccor*esccor
241 ! write (iout,*) "etot2", etot
242 ! write (iout,*) "evdw", wsc, evdw,evdw_t
243 ! write (iout,*) "evdw2", wscp, evdw2
244 ! write (iout,*) "welec", ft(1),welec,ees
245 ! write (iout,*) "evdw1", wvdwpp,evdw1
246 ! write (iout,*) "ebe",ebe,wang
248 write(iout,*) "nss",nss
249 Fdimless(i)=beta_h(ib)*etot+entfac(ii)
252 write (iout,*) "fdim calc", i,ii,ib,&
253 1.0d0/(1.987d-3*beta_h(ib)),totfree(i),&
254 entfac(ii),Fdimless(i)
256 Fdimless_(i)=Fdimless(i)
257 totfree_(i)=totfree(i)
258 call enerprint(energia(0)) !el
269 call MPI_Gatherv(Fdimless_(1),scount_(me),&
270 MPI_REAL,Fdimless(1),&
271 scount(0),idispl(0),MPI_REAL,Master,&
272 MPI_COMM_WORLD, IERROR)
273 call MPI_Gatherv(totfree_(1),scount(me),&
274 MPI_DOUBLE_PRECISION,totfree(1),&
275 scount_(0),idispl(0),MPI_DOUBLE_PRECISION,Master,&
276 MPI_COMM_WORLD, IERROR)
277 call MPI_Gatherv(entfac_(indstart(me)+1),scount_(me),&
278 MPI_DOUBLE_PRECISION,entfac(1),&
279 scount(0),idispl(0),MPI_DOUBLE_PRECISION,Master,&
280 MPI_COMM_WORLD, IERROR)
282 if (me.eq.Master) then
285 write (iout,*) "The FDIMLESS array before sorting"
287 write (iout,*) i,fdimless(i)
290 call mysort1(ncon,Fdimless,list_conf)
292 write (iout,*) "The FDIMLESS array after sorting"
294 write (iout,*) i,list_conf(i),fdimless(i)
298 totfree(i)=fdimless(i)
302 qfree=qfree+exp(-fdimless(i)+fdimless(1))
303 ! write (iout,*) "fdimless", fdimless(i)
305 write (iout,*) "qfree",qfree !d
308 write (iout,*) "ncon", ncon,maxstr_proc
309 do i=1,min0(ncon,maxstr_proc)-1
310 sumprob=sumprob+exp(-fdimless(i)+fdimless(1))/qfree
312 write (iout,*) i,ib,beta_h(ib),&
313 1.0d0/(1.987d-3*beta_h(ib)),list_conf(i),&
314 totfree(list_conf(i)),&
315 -entfac(list_conf(i)),fdimless(i),sumprob
317 if (sumprob.gt.prob_limit) goto 122
318 ! if (sumprob.gt.1.00d0) goto 122
324 call MPI_Bcast(nlist, 1, MPI_INTEGER, Master, MPI_COMM_WORLD, &
326 call MPI_Bcast(list_conf,nlist,MPI_INTEGER,Master,MPI_COMM_WORLD,&
329 ! write (iout,*) "iproc",iproc," indstart",indstart(iproc),
330 ! & " indend",indend(iproc)
333 !write(iout,*)"koniec probabl"
335 end subroutine probabl
336 !--------------------------------------------------
337 subroutine mysort1(n, x, ipermut)
339 integer :: i,j,imax,ipm,n
341 integer :: ipermut(n)
342 real(kind=4) :: xtemp
347 if (x(j).lt.xtemp) then
355 ipermut(imax)=ipermut(i)
359 end subroutine mysort1
360 !-----------------------------------------------------------------------------
361 !-----------------------------------------------------------------------------
362 end module probability