subroutine averages(iprot) implicit none include "DIMENSIONS" include "DIMENSIONS.ZSCOPT" #ifdef MPI include "mpif.h" integer IERROR,ErrCode include "COMMON.MPI" #endif integer n,nf,What include "COMMON.WEIGHTS" include "COMMON.WEIGHTDER" include "COMMON.ENERGIES" include "COMMON.IOUNITS" include "COMMON.VMCPAR" include "COMMON.NAMES" include "COMMON.INTERACT" include "COMMON.TIME1" include "COMMON.CHAIN" include "COMMON.PROTFILES" include "COMMON.COMPAR" include "COMMON.CLASSES" C Define local variables integer i,ii,iii,inat,jj,j,k,l,ik,jk,iprot,ib integer ipass_conf,istart_conf,iend_conf double precision energia(0:max_ene) double precision etoti,sumpart,esquarei,emeani,elowesti,enepsjk double precision aux,fac,ef1,ef2,em1,em2,var1,var2 double precision gnmr,tcpu_ini,tcpu logical lprn integer icant external icant lprn=.true. c c Zero out tables. c c Mazimum likeliood do ib=1,nbeta(1,iprot) do k=1,n_ene sumlikder(k,ib,iprot)=0.0d0 sumqder(k,ib,iprot)=0.0d0 enddo do k=1,nntyp sumlikeps(k,ib,iprot)=0.0d0 sumqeps(k,ib,iprot)=0.0d0 enddo do k=1,maxtor_var sumliktor(k,ib,iprot)=0.0d0 sumqtor(k,ib,iprot)=0.0d0 enddo sumlik(ib,iprot)=0.0d0 efree(ib,1,iprot)=0.0d0 enddo c Heat capacity do ib=1,nbeta(2,iprot) do k=1,n_ene eave_pftot(k,ib,iprot)=0.0d0 eave_pfprimtot(k,ib,iprot)=0.0d0 eave_pfbistot(k,ib,iprot)=0.0d0 emix_pftot(k,ib,iprot)=0.0d0 emix_pfprimtot(k,ib,iprot)=0.0d0 emix_pfbistot(k,ib,iprot)=0.0d0 emixsq_pftot(k,ib,iprot)=0.0d0 enddo do k=1,nntyp enepsave_ftot(k,ib,iprot)=0.0d0 eneps_mix_ftot(k,ib,iprot)=0.0d0 eneps_mix_fbistot(k,ib,iprot)=0.0d0 eneps_mixsq_ftot(k,ib,iprot)=0.0d0 enddo do k=1,maxtor_var entorave_ftot(k,ib,iprot)=0.0d0 entorave_fprimtot(k,ib,iprot)=0.0d0 entorave_fbistot(k,ib,iprot)=0.0d0 entor_mix_ftot(k,ib,iprot)=0.0d0 entor_mix_fprimtot(k,ib,iprot)=0.0d0 entor_mix_fbistot(k,ib,iprot)=0.0d0 entor_mixsq_ftot(k,ib,iprot)=0.0d0 enddo emean_ftot(ib,iprot)=0.0d0 ebis_ftot(ib,iprot)=0.0d0 esquare_ftot(ib,iprot)=0.0d0 enddo c Conformation-dependent averages do inat=1,natlike(iprot) do ib=1,nbeta(inat+2,iprot) do l=1,natdim(inat,iprot) do k=1,n_ene nu_pf(k,l,ib,inat,iprot)=0.0d0 enddo do k=1,nntyp nuepsave_f(k,l,ib,inat,iprot)=0.0d0 enddo do k=1,maxtor_var nutorave_f(k,l,ib,inat,iprot)=0.0d0 enddo nuave(l,ib,inat,iprot)=0.0d0 enddo enddo enddo c c Calculate the contributions to averages from each processor or all c contributions if calculations are run in uniprocessor mode. c The derivatives of energy in epsilons are dumped to disk, if needed. c #ifdef DEBUG write (iout,*) "Protein",iprot," nchunk_conf",nchunk_conf(iprot) #endif IF (NCHUNK_CONF(IPROT).EQ.1) THEN #ifdef MPI do i=1,ntot_work(iprot) i2ii(i,iprot)=0 enddo ii=0 do i=indstart(me1,iprot),indend(me1,iprot) ii=ii+1 i2ii(i,iprot)=ii enddo #else do i=1,ntot_work(iprot) i2ii(i,iprot)=i endif #endif #ifdef DEBUG do i=1,ntot_work(iprot) write (iout,*) "i",i," i2ii",i2ii(i,iprot) enddo call flush(iout) #endif call ave_eval(iprot) ELSE open (ientin,file=benefiles(iprot),status="old", & form="unformatted",access="direct",recl=lenrec_ene(iprot)) ipass_conf=0 #ifdef MPI do istart_conf=indstart(me1,iprot),indend(me1,iprot),maxstr_proc iend_conf=min0(istart_conf+maxstr_proc-1,indend(me1,iprot)) #else do istart_conf=1,ntot_work(iprot),maxstr_proc iend_conf=min0(istart_conf+maxstr_proc-1,ntot_work(iprot)) #endif c c Read the chunk of energies and derivatives off a DA scratchfile. c ipass_conf=ipass_conf+1 do i=1,ntot_work(iprot) i2ii(i,iprot)=0 enddo ii=0 do i=istart_conf,iend_conf ii=ii+1 i2ii(i,iprot)=ii enddo #ifdef DEBUG write (iout,*) "ipass_conf",ipass_conf, & " istart_conf",istart_conf," iend_conf",iend_conf do i=1,ntot_work(iprot) write (iout,*) "i",i," i2ii",i2ii(i,iprot) enddo call flush(iout) #endif call daread_ene(iprot,istart_conf,iend_conf) call ave_eval(iprot) enddo close(ientin) ENDIF c c Put tohether all contributions. c call ave_sum(iprot) return end c------------------------------------------------------------------------------- subroutine ave_eval(iprot) implicit none include "DIMENSIONS" include "DIMENSIONS.ZSCOPT" #ifdef MPI include "mpif.h" integer IERROR,ErrCode include "COMMON.MPI" #endif integer n,nf,What include "COMMON.WEIGHTS" include "COMMON.WEIGHTDER" include "COMMON.ENERGIES" include "COMMON.IOUNITS" include "COMMON.VMCPAR" include "COMMON.NAMES" include "COMMON.INTERACT" include "COMMON.TIME1" include "COMMON.CHAIN" include "COMMON.PROTFILES" include "COMMON.THERMAL" include "COMMON.VAR" include "COMMON.CLASSES" include "COMMON.COMPAR" C Define local variables integer i,ii,iii,jj,j,k,l,m,ik,jk,iprot,ib,inat, & ic integer ipass_conf,istart_conf,iend_conf double precision energia(0:max_ene) double precision etoti,sumpart,esquarei,emeani, & elowesti,enepsjk,eprim,ebis,eprimi,ebisi,etoti_orig, & entorjk,eave_pft(max_ene),emix_pft(max_ene), & esquare_ft,efree_t,emixsq_pft(max_ene), & eneps_mixt(nntyp),eneps_meant(nntyp), & enepsave_ft(nntyp),eneps_mix_ft(nntyp), & eneps_mixsq_ft(nntyp),emean_ft double precision aux,auxf,fac,facf,ef1,ef2,em1,em2,var1,var2, & deltei,deltei_orig,temper,elowest_all double precision gnmr,tcpu_ini,tcpu double precision ftune_epsprim external ftune_epsprim logical lprn integer icant external icant lprn=.false. C Compute the likelihood sum DO IB=1,NBETA(1,IPROT) elowest_all=elowest(ib,1,iprot) #ifdef DEBUG write(iout,*) "iprot",iprot," ib",ib, & " elowest",elowest_l(ib,iprot) #endif fac=betaT(ib,1,iprot) temper=1.0d0/(fac*Rgas) sumpart=0.0d0 sumlik(ib,iprot)=0.0d0 do k=1,n_ene sumlikder(k,ib,iprot)=0.0d0 sumqder(k,ib,iprot)=0.0d0 enddo do j=1,ntot_work(iprot) i = tsil(j,iprot) #ifdef DEBUG if (i.gt.0) write (iout,*) "i",i," iprot",iprot, & " indstart",indstart(me1,iprot), & " indend",indend(me1,iprot) #endif if (i.eq.0) goto 221 jj = i2ii(i,iprot) if (jj.gt.0) then c Temperature-dependent energy etoti=0.0d0 do k=1,n_ene etoti=etoti+ww(k)*escal(k,ib,1,iprot) & *enetb(jj,k,iprot) enddo deltei=etoti-elowest(ib,1,iprot) #ifdef DEBUG write (iout,*) "etoti",etoti," deltei",deltei write (iout,'(20f8.3)') (ww(k),k=1,n_ene) write (iout,'(20f8.3)') (enetb(jj,k,iprot), & k=1,n_ene) #endif aux=entfac(i,iprot)+fac*deltei #ifdef DEBUG write (iout,'(2i5,f5.2,7(a,e15.5))') & i,ib,fac," e_total",etoti, & " eini",eini(i,iprot)," entfac",entfac(i,iprot), & " e_lowest",elowest(ii,ib,iprot)," aux",aux #endif sumlik(ib,iprot)=sumlik(ib,iprot)+aux*Ptab(jj,ib,iprot) if (aux.le.150.0) then aux=dexp(-aux) sumpart=sumpart+aux do k=1,n_ene sumlikder(k,ib,iprot)=sumlikder(k,ib,iprot)+ & +enetb(jj,k,iprot)*escal(k,ib,1,iprot)*Ptab(jj,ib,iprot) sumqder(k,ib,iprot)=sumqder(k,ib,iprot) & +aux*enetb(jj,k,iprot)*escal(k,ib,1,iprot) enddo k=0 do ik=1,ntyp do jk=1,ik k=k+1 sumlikeps(k,ib,iprot)=sumlikeps(k,ib,iprot)+ & (ftune_epsprim(eps(ik,jk))* & eneps(1,k,jj,iprot)+eneps(2,k,jj,iprot)) & *Ptab(jj,ib,iprot) sumqeps(k,ib,iprot)=sumqeps(k,ib,iprot)+ & (ftune_epsprim(eps(ik,jk))* & eneps(1,k,jj,iprot)+eneps(2,k,jj,iprot))*aux enddo enddo k=0 do ik=1,ntyp do jk=1,ntyp if (mask_tor(jk,ik).gt.0) then k=k+1 entorjk=entor(k,jj,iprot) sumliktor(k,ib,iprot)=sumliktor(k,ib,iprot) & +Ptab(jj,ib,iprot)*entorjk*escal(13,ib,1,iprot) sumqtor(k,ib,iprot)= & sumqtor(k,ib,iprot)+aux*entorjk*escal(13,ib,1,iprot) endif enddo enddo endif endif 221 continue enddo ! j efree(ib,1,iprot)=sumpart #define DEBUG #ifdef DEBUG write (iout,*) " ib",ib," iprot",iprot, & " sumlik",sumlik(ib,iprot), & " sumq",efree(ib,1,iprot) #endif #undef DEBUG ENDDO ! ib C Calculation of heat capacity DO IB=1,NBETA(2,IPROT) elowest_all=elowest(ib,2,iprot) #define DEBUG #ifdef DEBUG write(iout,*) "iprot",iprot," ib",ib, & " elowest",elowest(ib,2,iprot) write(iout,*) "escalbis",(escalbis(k,ib,2,iprot),k=1,n_ene) #endif #undef DEBUG fac=betaT(ib,2,iprot) temper=1.0d0/(fac*Rgas) sumpart=0.0d0 emeani=0.0d0 esquarei=0.0d0 ebisi=0.0d0 do j=1,ntot_work(iprot) i = tsil(j,iprot) #ifdef DEBUG if (i.gt.0) write (iout,*) "i",i," iprot",iprot, & " indstart",indstart(me1,iprot), & " indend",indend(me1,iprot) #endif if (i.eq.0) cycle jj = i2ii(i,iprot) #ifdef DEBUG if (jj.gt.0) then write (iout,*) "ib",ib," j",j, write (iout,*) "nu",(nu(k,jj,iprot),k=1, & natconstr(iprot)) endif #endif if (jj.gt.0) then c Temperature-dependent energy etoti=0.0d0 etoti_orig=0.0d0 eprim=0.0d0 ebis=0.0d0 do k=1,n_ene etoti=etoti+ww(k)*escal(k,ib,2,iprot) & *enetb(jj,k,iprot) eprim=eprim+ww(k)*escalprim(k,ib,2,iprot) & *enetb(jj,k,iprot) ebis=ebis+ww(k)*escalbis(k,ib,2,iprot) & *enetb(jj,k,iprot) enddo deltei=etoti-elowest(ib,2,iprot) #ifdef DEBUG write (iout,*) "etoti",etoti," deltei",deltei write (iout,'(20f8.3)') (ww(k),k=1,n_ene) write (iout,'(20f8.3)') (enetb(jj,k,iprot), & k=1,n_ene) #endif aux=entfac(i,iprot)+fac*deltei #ifdef DEBUG write (iout,'(f5.2,7(a,e15.5))') & fac," e_total",etoti, & " eini",eini(i,iprot)," entfac",entfac(i,iprot), & " eprim",eprim," ebis",ebis, & " e_lowest",elowest(ib,2,iprot)," aux",aux #endif if (aux.le.150.0) then aux=dexp(-aux) sumpart=sumpart+aux etoti=etoti-temper*eprim emeani=emeani+etoti*aux esquarei=esquarei+aux*etoti**2 ebisi=ebisi+aux*ebis do k=1,n_ene eave_pftot(k,ib,iprot)=eave_pftot(k,ib,iprot) & +aux*enetb(jj,k,iprot)*escal(k,ib,2,iprot) c write (iout,*) "eave_pf",eave_pf(k,ii,ib,iprot) eave_pfprimtot(k,ib,iprot)=eave_pfprimtot(k,ib,iprot) & +aux*enetb(jj,k,iprot)*(escal(k,ib,2,iprot)- & temper*escalprim(k,ib,2,iprot)) eave_pfbistot(k,ib,iprot)=eave_pfbistot(k,ib,iprot) & +aux*enetb(jj,k,iprot)*escalbis(k,ib,2,iprot) emix_pftot(k,ib,iprot)=emix_pftot(k,ib,iprot) & +aux*enetb(jj,k,iprot)*etoti*escal(k,ib,2,iprot) emix_pfprimtot(k,ib,iprot)=emix_pfprimtot(k,ib,iprot) & +aux*enetb(jj,k,iprot)*etoti*(escal(k,ib,2,iprot) & -temper*escalprim(k,ib,2,iprot)) emix_pfbistot(k,ib,iprot)=emix_pfbistot(k,ib,iprot) & +aux*enetb(jj,k,iprot)*ebis*escal(k,ib,2,iprot) emixsq_pftot(k,ib,iprot)=emixsq_pftot(k,ib,iprot) & +aux*enetb(jj,k,iprot)*etoti**2*escal(k,ib,2,iprot) enddo k=0 do ik=1,ntyp do jk=1,ik k=k+1 enepsjk=ftune_epsprim(eps(ik,jk))* & eneps(1,k,jj,iprot)+eneps(2,k,jj,iprot) enepsave_ftot(k,ib,iprot)=enepsave_ftot(k,ib,iprot) & +aux*enepsjk eneps_mix_ftot(k,ib,iprot)= & eneps_mix_ftot(k,ib,iprot)+aux*enepsjk*etoti eneps_mix_fbistot(k,ib,iprot)= & eneps_mix_fbistot(k,ib,iprot)+aux*enepsjk*ebis eneps_mixsq_ftot(k,ib,iprot)= & eneps_mixsq_ftot(k,ib,iprot)+aux*enepsjk*etoti**2 enddo enddo k=0 do ik=1,ntyp do jk=1,ntyp if (mask_tor(jk,ik).gt.0) then k=k+1 entorjk=entor(k,jj,iprot) c write (iout,*) " k"," ik",ik," jk",jk, c & " entor",entorjk," contirb", c & aux*entorjk*escal(13,ib,2,iprot) entorave_ftot(k,ib,iprot)= & entorave_ftot(k,ib,iprot) & +aux*entorjk*escal(13,ib,2,iprot) c write (iout,*) "entorave_f", c & entorave_ftot(k,ib,iprot) entorave_fprimtot(k,ib,iprot)= & entorave_fprimtot(k,ib,iprot) & +aux*entorjk*(escal(13,ib,2,iprot)- & temper*escalprim(13,ib,2,iprot)) entorave_fbistot(k,ib,iprot)= & entorave_fbistot(k,ib,iprot) & +aux*entorjk*escalbis(13,ib,2,iprot) entor_mix_ftot(k,ib,iprot)= & entor_mix_ftot(k,ib,iprot) & +aux*entorjk*etoti*escal(13,ib,2,iprot) entor_mix_fprimtot(k,ib,iprot)= & entor_mix_fprimtot(k,ib,iprot) & +aux*entorjk*etoti*(escal(13,ib,2,iprot)- & temper*escalprim(13,ib,2,iprot)) entor_mix_fbistot(k,ib,iprot)= & entor_mix_fbistot(k,ib,iprot) & +aux*entorjk*ebis*escal(13,ib,2,iprot) entor_mixsq_ftot(k,ib,iprot)= & entor_mixsq_ftot(k,ib,iprot) & +aux*entorjk*etoti**2*escal(13,ib,2,iprot) endif enddo enddo endif endif efree(ib,2,iprot)=efree(ib,2,iprot)+sumpart emean_ftot(ib,iprot)=emean_ftot(ib,iprot)+emeani ebis_ftot(ib,iprot)=ebis_ftot(ib,iprot)+ebisi esquare_ftot(ib,iprot)=esquare_ftot(ib,iprot)+esquarei #ifdef DEBUG write (iout,*) "efree",efree(ib,2,iprot), & " emean",emean_ftot(ib,iprot), & " esquare",esquare_ftot(ib,iprot), & " ebis",ebis_ftot(ib,iprot), & " emix_pfprim",(emix_pfprimtot(k,ib,iprot),k=1,n_ene) #endif enddo ebis_ftot(ib,iprot)=ebis_ftot(ib,iprot)*temper do k=1,n_ene eave_pfbistot(k,ib,iprot)= & eave_pfbistot(k,ib,iprot)*temper emix_pfbistot(k,ib,iprot)= & emix_pfbistot(k,ib,iprot)*temper enddo do k=1,nntyp eneps_mix_fbistot(k,ib,iprot)= & eneps_mix_fbistot(k,ib,iprot)*temper enddo do k=1,ntor_var entorave_fbistot(k,ib,iprot)= & entorave_fbistot(k,ib,iprot)*temper entor_mix_fbistot(k,ib,iprot)= & entor_mix_fbistot(k,ib,iprot)*temper enddo #ifdef DEBUG write (iout,*) "eave_pf",(eave_pftot(k,ib,iprot), & k=1,n_ene) write (iout,*) "entorave_f",(entorave_ftot(k,ib,iprot), & k=1,ntor_var) #endif #define DEBUG #ifdef DEBUG write (iout,*) "ib",ib," temper",temper, & " ebis",ebis_ftot(ib,iprot) #endif #undef DEBUG ENDDO ! ib C Calculation of conformation-dependent averages DO INAT=1,NATLIKE(IPROT) DO IB=1,NBETA(I+2,IPROT) elowest_all=elowest(ib,inat+2,iprot) #ifdef DEBUG write(iout,*) "iprot",iprot," ib",ib, & " elowest",elowest(ib,iprot) #endif fac=betaT(ib,2,iprot) temper=1.0d0/(fac*Rgas) emeani=0.0d0 do j=1,ntot_work(iprot) i = tsil(j,iprot) #ifdef DEBUG if (i.gt.0) write (iout,*) "i",i," iprot",iprot, & " indstart",indstart(me1,iprot), & " indend",indend(me1,iprot) #endif if (i.eq.0) cycle jj = i2ii(i,iprot) #ifdef DEBUG if (jj.gt.0) then write (iout,*) "ib",ib," j",j, write (iout,*) "nu",(nu(k,jj,iprot),k=1, & natconstr(iprot)) endif #endif if (jj.gt.0) then c Temperature-dependent energy etoti=0.0d0 do k=1,n_ene etoti=etoti+ww(k)*escal(k,ib,inat+2,iprot) & *enetb(jj,k,iprot) enddo deltei=etoti-elowest(ib,inat+2,iprot) #ifdef DEBUG write (iout,*) "etoti",etoti," deltei",deltei write (iout,'(20f8.3)') (ww(k),k=1,n_ene) write (iout,'(20f8.3)') (enetb(jj,k,iprot), & k=1,n_ene) #endif aux=entfac(i,iprot)+fac*deltei #ifdef DEBUG write (iout,'(f5.2,7(a,e15.5))') & fac," e_total",etoti, & " eini",eini(i,iprot)," entfac",entfac(i,iprot), & " eprim",eprim," ebis",ebis, & " e_lowest",elowest(ib,inat+2,iprot)," aux",aux #endif if (aux.le.150.0) then aux=dexp(-aux) sumpart=sumpart+aux c 4/13/04 AL Components of the conformation-dependent averages do k=1,n_ene eave_nat_pftot(k,ib,inat,iprot)= & eave_nat_pftot(k,ib,inat,iprot) & +aux*enetb(jj,k,iprot)*escal(k,ib,inat+2,iprot) enddo k=0 do ik=1,ntyp do jk=1,ik k=k+1 enepsjk=ftune_epsprim(eps(ik,jk))* & eneps(1,k,jj,iprot)+eneps(2,k,jj,iprot) enepsave_nat_ftot(k,ib,inat,iprot)= & enepsave_nat_ftot(k,ib,inat,iprot) & +aux*enepsjk enddo enddo k=0 do ik=1,ntyp do jk=1,ntyp if (mask_tor(jk,ik).gt.0) then k=k+1 entorjk=entor(k,jj,iprot) c write (iout,*) " k"," ik",ik," jk",jk, c & " entor",entorjk," contirb", c & aux*entorjk*escal(13,ib,inat+2,iprot) entorave_nat_ftot(k,ib,inat,iprot)= & entorave_nat_ftot(k,ib,inat,iprot) & +aux*entorjk*escal(13,ib,inat+2,iprot) endif enddo enddo do l=1,natdim(inat,iprot) nuave(l,ib,inat,iprot)=nuave(l,ib,inat,iprot) & +aux*nu(l,inat,jj,iprot) do k=1,n_ene nu_pf(k,l,ib,inat,iprot)=nu_pf(k,l,ib,inat,iprot) & +aux*enetb(jj,k,iprot)*nu(k,inat,jj,iprot)* & escal(k,ib,inat+2,iprot) enddo k=0 do ik=1,ntyp do jk=1,ik k=k+1 nuepsave_f(k,l,ib,inat,iprot)= & nuepsave_f(k,l,ib,inat,iprot)+aux*enepsjk* & nu(l,inat,jj,iprot) enddo enddo k=0 do ik=1,ntyp do jk=1,ntyp if (mask_tor(jk,ik).gt.0) then k=k+1 nutorave_f(k,l,ib,inat,iprot)= & nutorave_f(k,l,ib,inat,iprot) & +aux*entorjk*escal(13,ib,inat+2,iprot)* & nu(l,inat,jj,iprot) endif enddo enddo enddo endif endif #ifdef DEBUG write (iout,*) "iprot",iprot," ib",ib, write (iout,*) "nuave" write (iout,'(20f10.5)') (nuave(k,ib,iprot), & k=1,natconstr(iprot)) #endif #ifdef DEBUG write (iout,*) "inat",inat, & " efree_nat",efree_nat(ib,inat,iprot) #endif enddo #ifdef DEBUG write (iout,*) "iprot",iprot," ib",ib write (iout,*) "nuave0" write (iout,'(20f10.5)') (nuave(k,ib,inat,iprot), & k=1,natdim(inat,iprot)) #endif ENDDO ! ib ENDDO ! INAT return end c------------------------------------------------------------------------------- subroutine ave_sum(iprot) implicit none include "DIMENSIONS" include "DIMENSIONS.ZSCOPT" #ifdef MPI include "mpif.h" integer IERROR,ErrCode include "COMMON.MPI" #endif integer n,nf,What include "COMMON.WEIGHTS" include "COMMON.WEIGHTDER" include "COMMON.ENERGIES" include "COMMON.IOUNITS" include "COMMON.VMCPAR" include "COMMON.NAMES" include "COMMON.INTERACT" include "COMMON.TIME1" include "COMMON.CHAIN" include "COMMON.PROTFILES" include "COMMON.VAR" include "COMMON.COMPAR" C Define local variables integer i,ii,iii,jj,j,k,l,ik,jk,iprot,ib,inat integer ipass_conf,istart_conf,iend_conf double precision energia(0:max_ene) double precision etoti,sumpart,esquarei,emeani,elowesti,enepsjk, & eave_pft(max_ene),emix_pft(max_ene),eave_pfprimt(max_ene), & eave_pfbist(max_ene),emix_pfprimt(max_ene),emix_pfbist(max_ene), & esquare_ft,efree_t, & emixsq_pft(max_ene),eneps_mixt(nntyp),eneps_meant(nntyp), & enepsave_ft(nntyp),eneps_mix_ft(nntyp),entorave_ft(maxtor_var), & entor_mix_ft(maxtor_var), & eneps_mix_fbist(nntyp),entorave_fbist(maxtor_var), & entorave_fprimt(maxtor_var),entor_mix_fprimt(maxtor_var), & entor_mix_fbist(maxtor_var),eneps_mixsq_ft(nntyp), & entor_mixsq_ft(maxtor_var),emean_ft,ebis_ft,nuave_t(maxdimnat), & nu_pft(max_ene,maxdimnat), & nuepsave_ft(nntyp,maxdimnat), & nutorave_ft(maxtor_var,maxdimnat), & sumlik_t,sumlikder_t(max_ene),sumlikeps_t(nntyp), & sumliktor_t(maxtor_var), & sumq_t,sumqder_t(max_ene),sumqeps_t(nntyp),sumqtor_t(maxtor_var) double precision aux,fac,ef1,ef2,em1,em2,var1,var2,efree_tot,facF, & elowest_all double precision gnmr,tcpu_ini,tcpu logical lprn integer icant external icant lprn=.true. c Maximum likelihood contribution DO IB=1,NBETA(1,IPROT) fac=betaT(ib,1,iprot) #ifdef MPI #define DEBUG #ifdef DEBUG write (iout,*) "Processor",me,me1," calling MPI_Reduce: 3" write (iout,*) "iprot",iprot," ib",ib write (iout,*) "sumlik",sumlik(ib,iprot), & " sumq",efree(ib,1,iprot) call flush(iout) #endif call MPI_Reduce( sumlik(ib,iprot), sumlik_t,1, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce( efree(ib,1,iprot), sumq_t,1, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) #ifdef DEBUG write (iout,*) "sumlik",sumlik(ib,iprot)," sumlik_t",sumlik_t write (iout,*) "sumq",efree(ib,1,iprot)," sumq_t",sumq_t call flush(iout) #endif call MPI_Reduce(sumlikder(1,ib,iprot),sumlikder_t(1), & n_ene, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(sumlikeps(1,ib,iprot), & sumlikeps_t(1), nntyp, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(sumliktor(1,ib,iprot), & sumliktor_t(1), ntor_var, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(sumqder(1,ib,iprot),sumqder_t(1),n_ene, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(sumqeps(1,ib,iprot),sumqeps_t(1), nntyp, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(sumqtor(1,ib,iprot),sumqtor_t(1), ntor_var, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) #ifdef DEBUG write (iout,*) "Processor",me,me1," finished MPI_Reduce: 3" call flush(iout) #endif if (me1.eq.master) then #ifdef DEBUG write (iout,*) "ib",ib, & " elowest",elowest(ib,1,iprot), & " sumlik",sumlik_t," qsum",sumq_t," fac",fac #endif do k=1,n_ene sumlikder(k,ib,iprot)=-sumlikder_t(k) & +sumqder_t(k)/sumq_t enddo do k=1,nntyp sumlikeps(k,ib,iprot)=-sumlikeps_t(k)/sumq_t enddo c write (iout,*) "eavetor",eave_pftot(13,ib,iprot), c & eave_pftot(19,ib,iprot) do k=1,ntor_var sumliktor(k,ib,iprot)=-sumliktor_t(k)/sumq_t enddo sumlik(ib,iprot)=-sumlik(ib,iprot)+dlog(sumq_t) & -elowest(ib,1,iprot)*fac efree(ib,1,iprot)=sumq_t #ifdef DEBUUG write (iout,*) "ib",ib," iprot",iprot," final sumlik", & sumlik(ib,iprot)," sumq",efree(ib,1,iprot) #endif endif #undef DEBUG #else do k=1,n_ene sumlikder(k,ib,iprot)=sumlikder_t(k,ib,iprot) & +sumqder(k,ib,iprot)/efree(ib,1,iprot) enddo do k=1,nntyp sumlikeps(k,ib,iprot)= & sumlikeps(k,ib,iprot)/efree(ib,1,iprot) enddo do k=1,ntor_var sumliktor(k,ib,iprot)= & sumliktor(k,ib,iprot)/efree(ib,1,iprot) enddo sumlik(ib,iprot)=sumlik(ib,iprot)+dlog(efree(ib,1,iprot) & -elowest(nbeta(iprot)+ib,iprot)*fac #endif ENDDO ! ib c Heat capacity and averages DO IB=1,NBETA(2,IPROT) fac=betaT(ib,2,iprot) #ifdef MPI #ifdef DEBUG write (iout,*) "Processor",me,me1," calling MPI_Reduce: 3" write (iout,*) "iprot",iprot," ib",ib call flush(iout) #endif call MPI_Reduce( efree(ib,2,iprot), efree_t,1, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) #ifdef DEBUG write (iout,*) "efree",efree(iprot write (iout,*) "efree_t",efree_t call flush(iout) #endif call MPI_Reduce(emean_ftot(ib,iprot),emean_ft,1, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(ebis_ftot(ib,iprot),ebis_ft,1, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(esquare_ftot(ib,iprot),esquare_ft,1, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(eave_pftot(1,ib,iprot),eave_pft, & n_ene, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(eave_pfprimtot(1,ib,iprot), & eave_pfprimt,n_ene, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(eave_pfbistot(1,ib,iprot), & eave_pfbist,n_ene, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) c write (iout,*) "eave_pf",(eave_pf(k,iprot), c & k=1,n_ene) c write (iout,*) "eave_pft",(eave_pft(k),k=1,n_ene) call MPI_Reduce(emix_pftot(1,ib,iprot),emix_pft, & n_ene, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(emix_pfprimtot(1,ib,iprot), & emix_pfprimt,n_ene, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(emix_pfbistot(1,ib,iprot), & emix_pfbist,n_ene, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(emixsq_pftot(1,ib,iprot), & emixsq_pft, n_ene, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(enepsave_ftot(1,ib,iprot), & enepsave_ft, nntyp, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(eneps_mix_ftot(1,ib,iprot), & eneps_mix_ft,nntyp, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(eneps_mix_fbistot(1,ib,iprot), & eneps_mix_fbist,nntyp, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(eneps_mixsq_ftot(1,ib,iprot), & eneps_mixsq_ft,nntyp, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) c write (iout,*) "enepsave_f",(enepsave_f(k,iprot), c & k=1,nntyp) c write (iout,*) "enepsave_ft",(enepsave_ft(k), c & k=1,nntyp) call MPI_Reduce(entorave_ftot(1,ib,iprot), & entorave_ft, ntor_var, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(entorave_fprimtot(1,ib,iprot), & entorave_fprimt, ntor_var, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(entorave_fbistot(1,ib,iprot), & entorave_fbist, ntor_var, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(entor_mix_ftot(1,ib,iprot), & entor_mix_ft,ntor_var, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(entor_mix_fprimtot(1,ib,iprot), & entor_mix_fprimt,ntor_var, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(entor_mix_fbistot(1,ib,iprot), & entor_mix_fbist,ntor_var, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(entor_mixsq_ftot(1,ib,iprot), & entor_mixsq_ft,ntor_var, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) c write (iout,*) "entorave_f",(entorave_f(k,iprot), c & k=1,ntor_var) c write (iout,*) "entorave_ft",(entorave_ft(k), c & k=1,ntor_var) #ifdef DEBUG write (iout,*) "Processor",me,me1," finished MPI_Reduce: 3" call flush(iout) #endif if (me1.eq.master) then elowest_all=elowest(ib,2,iprot) #ifdef DEBUG write (iout,*) "ib",ib, & " elowest",elowest(ib,iprot), & "efree",efree_t," fac",fac," facF",facF, & " efree_tot",efree_tot #endif do k=1,n_ene eave_pftot(k,ib,iprot)=eave_pft(k)/efree_t eave_pfprimtot(k,ib,iprot)=eave_pfprimt(k)/efree_t eave_pfbistot(k,ib,iprot)=eave_pfbist(k)/efree_t emix_pftot(k,ib,iprot)=emix_pft(k)/efree_t emix_pfprimtot(k,ib,iprot)=emix_pfprimt(k)/efree_t emix_pfbistot(k,ib,iprot)=emix_pfbist(k)/efree_t emixsq_pftot(k,ib,iprot)=emixsq_pft(k)/efree_t enddo do k=1,nntyp enepsave_ftot(k,ib,iprot)=enepsave_ft(k)/efree_t eneps_mix_ftot(k,ib,iprot)=eneps_mix_ft(k)/efree_t eneps_mix_fbistot(k,ib,iprot)=eneps_mix_fbist(k)/efree_t eneps_mixsq_ftot(k,ib,iprot)=eneps_mixsq_ft(k)/efree_t enddo c write (iout,*) "eavetor",eave_pftot(13,ib,iprot), c & eave_pftot(19,ib,iprot) do k=1,ntor_var entorave_ftot(k,ib,iprot)=entorave_ft(k)/ & efree_t c write (iout,*) "iprot",iprot," ib",ib," k",k, c & " entorave_ftot", entorave_ftot(k,ib,iprot) entorave_fprimtot(k,ib,iprot)=entorave_fprimt(k)/efree_t entorave_fbistot(k,ib,iprot)=entorave_fbist(k)/efree_t entor_mix_ftot(k,ib,iprot)=entor_mix_ft(k)/efree_t entor_mix_fprimtot(k,ib,iprot)=entor_mix_fprimt(k)/efree_t entor_mix_fbistot(k,ib,iprot)=entor_mix_fbist(k)/efree_t entor_mixsq_ftot(k,ib,iprot)=entor_mixsq_ft(k)/efree_t enddo emean_ftot(ib,iprot)=emean_ft/efree_t ebis_ftot(ib,iprot)=ebis_ft/efree_t esquare_ftot(ib,iprot)=esquare_ft/efree_t efree(ib,2,iprot)=-dlog(efree_t)/fac+elowest(ib,2,iprot) endif #else do k=1,n_ene eave_pftot(k,ib,iprot)=eave_pftot(k,ib,iprot) & /efree(ib,2,iprot) eave_pfprimtot(k,ib,iprot)=eave_pfprimtot(k,ib,iprot) & /efree(ib,2,iprot) eave_pfbistot(k,ib,iprot)=eave_pfbistot(k,ib,iprot) & /efree(ib,2,iprot) emix_pftot(k,ib,iprot)=emix_pftot(k,ib,iprot)/efree(ib,iprot) emix_pfprimtot(k,ib,iprot)=emix_pfprimtot(k,ib,iprot) & /efree(ib,2,iprot) emix_pfbistot(k,ib,iprot)=emix_pfbistot(k,ib,iprot) & /efree(ib,2,iprot) emixsq_pftot(k,ib,iprot)=emixsq_pftot(k,ib,iprot) & /efree(ib,2,iprot) enddo do k=1,nntyp enepsave_ftot(k,ib,iprot)=enepsave_ftot(k,ib,iprot) & /efree(ib,2,iprot) eneps_mix_ftot(k,ib,iprot)=eneps_mix_ftot(k,ib,iprot) & /efree(ib,2,iprot) eneps_mixsq_ftot(k,ib,iprot)=eneps_mixsq_ftot(k,ib,iprot) & /efree(ib,2,iprot) enddo do k=1,ntor_var entorave_ftot(k,ib,iprot)=entorave_f(k,ib,iprot) & /efree(ib,2,iprot) entor_mix_ftot(k,ib,iprot)=entor_mix_ftot(k,ib,iprot) & /efree(ib,2,iprot) entor_mixsq_ftot(k,ib,iprot)=entor_mixsq_ftot(k,ib,iprot) & /efree(ib,2,iprot) enddo emean_ftot(ib,iprot)=emean_ftot(ib,iprot)/efree(ib,2,iprot) ebis_ftot(ib,iprot)=ebis_ftot(ib,iprot)/efree(ib,2,iprot) esquare_ftot(ib,iprot)=esquare_ftot(ib,iprot)/efree(ib,2,iprot) #endif ENDDO ! IB c 4/13/04 AL Components of the correlation coefficients and their derivatives DO INAT=1,NATLIKE(IPROT) DO IB=1,NBETA(INAT+2,IPROT) #ifdef MPI call MPI_Reduce( efree(ib,inat+2,iprot), efree_t,1, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce( eave_nat_pftot(1, ib,inat,iprot), eave_pft, & n_ene, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce( enepsave_nat_ftot(1, ib,inat,iprot), & enepsave_ft, nntyp, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(entorave_nat_ftot(1,ib,inat,iprot), & entorave_ft, ntor_var, & MPI_DOUBLE_PRECISION, MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(nuave(1,ib,inat,iprot),nuave_t(1), & natdim(inat,iprot), MPI_DOUBLE_PRECISION, MPI_SUM, & Master, Comm1, IERROR) #ifdef DEBUG write (iout,*) "After REDUCE nuave",iprot,ib write (iout,'(20f10.5)') & (nuave(l,ib,iprot),l=1,natconstr(iprot)) #endif call MPI_Reduce(nu_pf(1,1,ib,inat,iprot),nu_pft(1,1), & max_ene*natdim(inat,iprot), MPI_DOUBLE_PRECISION, & MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(nuepsave_f(1,1,ib,inat,iprot), & nuepsave_ft(1,1), & nntyp*natdim(inat,iprot), MPI_DOUBLE_PRECISION, & MPI_SUM, Master, Comm1, IERROR) call MPI_Reduce(nutorave_f(1,1,ib,inat,iprot), & nutorave_ft(1,1), & maxtor_var*natdim(inat,iprot), MPI_DOUBLE_PRECISION, & MPI_SUM, Master, Comm1, IERROR) #ifdef DEBUG write (iout,*) "Processor",me,me1," finished MPI_Reduce: 3" call flush(iout) #endif if (me1.eq.master) then #ifdef DEBUG write (iout,*) "ib",ib, & " elowest",elowest(ib,iprot), & "efree",efree_t," fac",fac, & " efree_tot",efree_tot #endif do l=1,natdim(inat,iprot) do k=1,n_ene nu_pf(k,l,ib,inat,iprot)=nu_pft(k,l)/efree_t enddo do k=1,nntyp nuepsave_f(k,l,ib,inat,iprot)=nuepsave_ft(k,l)/efree_t enddo do k=1,ntor_var nutorave_f(k,l,ib,inat,iprot)=nutorave_ft(k,l)/efree_t enddo nuave(l,ib,inat,iprot)=nuave_t(l)/efree_t enddo do k=1,n_ene eave_nat_pftot(k,ib,inat,iprot)=eave_pft(k)/efree_t enddo do k=1,nntyp enepsave_nat_ftot(k,ib,inat,iprot)=enepsave_ft(k)/efree_t enddo do k=1,ntor_var entorave_nat_ftot(k,ib,inat,iprot)=entorave_ft(k)/efree_t enddo endif #else do l=1,natdim(inat,iprot) do k=1,n_ene nu_pf(k,l,ib,inat,iprot)=nu_pf(k,l,ib,inat,iprot) & /efree(ib,inat+2,iprot) enddo do k=1,nntyp nuepsave_f(k,l,ib,inat,iprot)=nuepsave_f(k,l,ib,inat,iprot) & /efree(ib,inat+2,iprot) enddo do k=1,ntor_var nutorave_ftot(k,l,ib,inat,iprot)= & nutorave_ftot(k,l,ib,inat,iprot) & /efree(ib,inat+2,iprot) enddo nuave(l,ib,inat,iprot)=nuave(l,ib,inat,iprot) & /efree(ib,inat+2,iprot) enddo ! l do k=1,n_ene eave_nat_pftot(k,ib,inat,iprot)= & eave_nat_pftot(k,ib,inat,iprot) & /efree(ib,inat+2,iprot) enddo do k=1,nntyp enepsave_nat_ftot(k,ib,inat,iprot)= & enepsave_nat_ftot(k,ib,inat,iprot)/efree(ib,inat+2,iprot) enddo do k=1,ntor_var enetorave_nat_ftot(k,ib,inat,iprot)= & enetorave_nat_ftot(k,ib,inat,iprot)/efree(ib,inat+2,iprot) enddo #endif ENDDO ! IB #ifdef DEBUG write (iout,*) "ib",ib," efree_tot",efree_tot #endif ENDDO ! INAT return end