From: Adam Sieradzan Date: Wed, 31 Jul 2024 14:17:13 +0000 (+0200) Subject: major change X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?a=commitdiff_plain;h=5c7bdad38e52abba274de9ec07f36ddb6c2dff5d;p=unres4.git major change --- diff --git a/source/cluster/cluster.F90 b/source/cluster/cluster.F90 index a47f418..3f5f4fd 100644 --- a/source/cluster/cluster.F90 +++ b/source/cluster/cluster.F90 @@ -9,7 +9,7 @@ use io_clust !#define CLUSTER use io_units - use io_base, only: permut + use io_base, only: permut, read_dist_constr,gen_dist_constr2 use geometry_data, only: nres,theta,phi,alph,omeg,& c,cref use energy_data, only: nnt,nct @@ -111,7 +111,7 @@ ! write (iout,*) "Before permut" ! write (iout,*) "symetr", symetr ! call flush(iout) - call permut(symetr) +! call permut(symetr) ! write (iout,*) "after permut" ! call flush(iout) print *,'MAIN: nnt=',nnt,' nct=',nct @@ -212,6 +212,7 @@ write(iout,*) "after daread_ccoords" ind1=0 DO I=1,NCON_work-1 +! print *,"Tu1",i if (mod(i,100).eq.0) print *,'Calculating RMS i=',i do k=1,2*nres do l=1,3 @@ -225,6 +226,7 @@ enddo enddo DO J=I+1,NCON_work +! print *,"tu2",j IND=IOFFSET(NCON_work,I,J) #ifdef MPI if (ind.ge.indstart(me) .and. ind.le.indend(me)) then @@ -490,17 +492,20 @@ real(kind=8) :: przes(3),obrot(3,3) real(kind=8) :: xx(3,2*nres),yy(3,2*nres) !(3,maxres2) integer :: i,ii,j,icon,jcon,kkk,nperm,chalen,zzz - integer :: iaperm,ibezperm,run + integer :: iaperm,ibezperm,run,nequiv(50) real(kind=8) :: rms,rmsmina -! write (iout,*) "tu dochodze" +! print *, "tu dochodze",symetr rmsmina=10d10 nperm=1 do i=1,symetr nperm=i*nperm enddo ! write (iout,*) "nperm",nperm - call permut(symetr) + if (symetr.ne.1) call permut(nequiv(i),nperm,tabperm) +! call permut(symetr) +! print *,"after permut" ! write (iout,*) "tabperm", tabperm(1,1) + if (symetr.ne.1) then do kkk=1,nperm if (lside) then ii=0 @@ -559,6 +564,20 @@ if (non_conv) print *,non_conv,icon,jcon if (rmsmina.gt.rms) rmsmina=rms enddo + else +! print *,"in else difconf" + nstart=2 + kkk=1 + do i=1,2*nres + do j=1,3 + c(j,i)=allcart(j,i,jcon) + enddo + enddo + call fitsq(rms,c(1,nstart),cref(1,nstart,kkk),nend-nstart+1,& + przes,& + obrot,non_conv) + rmsmina=rms + endif difconf=dsqrt(rmsmina) return end function difconf diff --git a/source/cluster/io_clust.F90 b/source/cluster/io_clust.F90 index 7f3e42d..a37bead 100644 --- a/source/cluster/io_clust.F90 +++ b/source/cluster/io_clust.F90 @@ -1405,6 +1405,7 @@ !c if (constr_homology) tole=dmax1(tole,1.5d0) write (iout,*) "with_homology_constr ",with_dihed_constr,& " CONSTR_HOMOLOGY",constr_homology + write(iout,*) "CONSTR_DIST=",constr_dist read_homol_frag = index(controlcard,"READ_HOMOL_FRAG").gt.0 out_template_coord = index(controlcard,"OUT_TEMPLATE_COORD").gt.0 out_template_restr = index(controlcard,"OUT_TEMPLATE_RESTR").gt.0 @@ -1929,7 +1930,7 @@ indpdb,constr_dist,raw_psipred, with_theta_constr use geometry, only: chainbuild,alloc_geo_arrays use energy, only: alloc_ener_arrays - use io_base, only: rescode + use io_base, only: rescode,read_dist_constr,gen_dist_constr2 use control, only: setup_var,init_int_table use conform_compar, only: contact ! implicit none @@ -2374,6 +2375,17 @@ !C enddo endif ! ntheta_constr.gt.0 endif! with_theta_constr + if (constr_dist.gt.0) then + call read_dist_constr + + link_start=1 + link_end=nhpb + clink_start=1 + clink_end=cnhpb + endif + + + if (constr_homology.gt.0) then !c write (iout,*) "About to call read_constr_homology" !c call flush(iout) @@ -3126,5 +3138,243 @@ 10 stop "Error in fragment file" end subroutine read_klapaucjusz -!----------------------------------------------------------------------------- + + +! subroutine read_dist_constr +! use MPI_data +!! use control +! use geometry, only: dist +! use geometry_data +! use control_data +! use energy_data +!! implicit real*8 (a-h,o-z) +!! include 'DIMENSIONS' +!#ifdef MPI +! include 'mpif.h' +!#endif +!! include 'COMMON.SETUP' +!! include 'COMMON.CONTROL' +!! include 'COMMON.CHAIN' +!! include 'COMMON.IOUNITS' +!! include 'COMMON.SBRIDGE' +! integer,dimension(2,100) :: ifrag_,ipair_ +! real(kind=8),dimension(100) :: wfrag_,wpair_ +! character(len=640) :: controlcard +! +!!el local variables +! integer :: i,k,j,ii,jj,itemp,mnumkk,mnumjj,k1,j1 +! integer :: nfrag_,npair_,ndist_ +! real(kind=8) :: dist_cut,ddjk +! +!! write (iout,*) "Calling read_dist_constr" +!! write (iout,*) "nres",nres," nstart_sup",nstart_sup," nsup",nsup +!! call flush(iout) +! if(.not.allocated(cihpb)) allocate(cihpb(maxdim)) +! if(.not.allocated(cjhpb)) allocate(cjhpb(maxdim)) +! if(.not.allocated(cdhpb)) allocate(cdhpb(maxdim)) +! if(.not.allocated(cdhpb1)) allocate(cdhpb1(maxdim)) +! if(.not.allocated(cforcon)) allocate(cforcon(maxdim)) +! if(.not.allocated(cfordepth)) allocate(cfordepth(maxdim)) +! if(.not.allocated(cibecarb)) allocate(cibecarb(maxdim)) +! if ((genconstr.gt.0).and.(constr_dist.eq.11)) then +! call gen_dist_constr2 +! go to 1712 +! endif +! cnhpb=0 +! call card_concat(controlcard,.true.) +! call readi(controlcard,"NFRAG",nfrag_,0) +! call readi(controlcard,"NPAIR",npair_,0) +! call readi(controlcard,"NDIST",ndist_,0) +! call reada(controlcard,'DIST_CUT',dist_cut,5.0d0) +! call multreadi(controlcard,"IFRAG",ifrag_(1,1),2*nfrag_,0) +! call multreadi(controlcard,"IPAIR",ipair_(1,1),2*npair_,0) +! call multreada(controlcard,"WFRAG",wfrag_(1),nfrag_,0.0d0) +! call multreada(controlcard,"WPAIR",wpair_(1),npair_,0.0d0) +! write (iout,*) "NFRAG",nfrag_," NPAIR",npair_," NDIST",ndist_ +!! write (iout,*) "IFRAG" +! do i=1,nfrag_ +! write (iout,*) i,ifrag_(1,i),ifrag_(2,i),wfrag_(i) +! enddo +!! write (iout,*) "IPAIR" +!! do i=1,npair_ +!! write (iout,*) i,ipair_(1,i),ipair_(2,i),wpair_(i) +!! enddo +!! if(.not.allocated(ihpb)) allocate(ihpb(maxdim)) +!! if(.not.allocated(jhpb)) allocate(jhpb(maxdim)) +!! if(.not.allocated(dhpb)) allocate(dhpb(maxdim)) +!! if(.not.allocated(forcon)) allocate(forcon(maxdim)) +! +! call flush(iout) +! do i=1,nfrag_ +! call flush(iout) +! if (wfrag_(i).gt.0.0d0) then +! do j=ifrag_(1,i),ifrag_(2,i)-1 +! do k=j+1,ifrag_(2,i) +! write (iout,*) "j",j," k",k,nres +! j1=j +! k1=k +! if (j.gt.nres) j1=j-nres +! if (k.gt.nres) k1=k-nres +! mnumkk=molnum(k1) +! mnumjj=molnum(j1) +! +! if ((itype(k1,mnumkk).gt.ntyp_molec(mnumkk)).or.& +! (itype(j1,mnumjj).gt.ntyp_molec(mnumjj))) cycle +! write (iout,*) "j",j," k",k,itype(k1,mnumkk),itype(j1,mnumjj) +! ddjk=dist(j,k) +! if (constr_dist.eq.1) then +! cnhpb=cnhpb+1 +! cihpb(cnhpb)=j +! cjhpb(cnhpb)=k +! cdhpb(cnhpb)=ddjk +! cforcon(cnhpb)=wfrag_(i) +! else if (constr_dist.eq.2) then +! if (ddjk.le.dist_cut) then +! print *,"tu",nhpb,i,j,k,maxdim +! cnhpb=cnhpb+1 +! cihpb(cnhpb)=j +! cjhpb(cnhpb)=k +! cdhpb(cnhpb)=ddjk +! cforcon(cnhpb)=wfrag_(i) +! endif +! else +! cnhpb=cnhpb+1 +! cihpb(cnhpb)=j +! cjhpb(cnhpb)=k +! cdhpb(cnhpb)=ddjk +! cforcon(cnhpb)=wfrag_(i)*dexp(-0.5d0*(ddjk/dist_cut)**2) +! endif +!#ifdef MPI +! if (.not.out1file .or. me.eq.king) & +! write (iout,'(a,3i5,f8.2,1pe12.2)') "+dist.constr ",& +! cnhpb,cihpb(cnhpb),cjhpb(cnhpb),cdhpb(cnhpb),cforcon(cnhpb) +!#else +! write (iout,'(a,3i5,f8.2,1pe12.2)') "+dist.constr ",& +! cnhpb,cihpb(cnhpb),cjhpb(cnhpb),cdhpb(cnhpb),cforcon(cnhpb) +!#endif +! enddo +! enddo +! endif +! enddo +! do i=1,npair_ +! if (wpair_(i).gt.0.0d0) then +! ii = ipair_(1,i) +! jj = ipair_(2,i) +! if (ii.gt.jj) then +! itemp=ii +! ii=jj +! jj=itemp +! endif +! do j=ifrag_(1,ii),ifrag_(2,ii) +! do k=ifrag_(1,jj),ifrag_(2,jj) +! cnhpb=cnhpb+1 +! cihpb(cnhpb)=j +! cjhpb(cnhpb)=k +! cforcon(cnhpb)=wpair_(i) +! cdhpb(cnhpb)=dist(j,k) +!#ifdef MPI +! if (.not.out1file .or. me.eq.king) & +! write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",& +! nhpb,cihpb(nhpb),cjhpb(nhpb),cdhpb(nhpb),cforcon(nhpb) +!#else +! write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",& +! nhpb,cihpb(nhpb),cjhpb(nhpb),cdhpb(nhpb),cforcon(nhpb) +!#endif +! enddo +! enddo +! endif +! enddo +! do i=1,ndist_ +!! read (inp,*) ihpb(nhpb+1),jhpb(nhpb+1),forcon(nhpb+1) +!! if (forcon(nhpb+1).gt.0.0d0) then +!! nhpb=nhpb+1 +!! dhpb(nhpb)=dist(ihpb(nhpb),jhpb(nhpb)) +! if (constr_dist.eq.11) then +! read (inp,*) cihpb(cnhpb+1),cjhpb(cnhpb+1),cdhpb(i),cdhpb1(i), & +! cibecarb(i),cforcon(cnhpb+1),cfordepth(cnhpb+1) +!!EL fordepth(nhpb+1)=fordepth(nhpb+1)/forcon(nhpb+1) +! cfordepth(cnhpb+1)=cfordepth(cnhpb+1)**(0.25d0) +! cforcon(cnhpb+1)=cforcon(cnhpb+1)**(0.25d0) +! else +! print *,"in else" +! read (inp,*) cihpb(cnhpb+1),cjhpb(cnhpb+1),cdhpb(i),cdhpb1(i), & +! cibecarb(i),cforcon(cnhpb+1) +! print *,"afterread",cihpb(cnhpb+1) +! endif +! if (cforcon(cnhpb+1).gt.0.0d0) then +! cnhpb=cnhpb+1 +! if (cibecarb(i).gt.0) then +! cihpb(cnhpb)=cihpb(cnhpb)+nres +! cjhpb(cnhpb)=cjhpb(cnhpb)+nres +! endif +! if (cdhpb(cnhpb).eq.0.0d0) & +! cdhpb(cnhpb)=dist(ihpb(cnhpb),jhpb(cnhpb)) +! endif +! +!#ifdef MPI +! if (.not.out1file .or. me.eq.king) & +! write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",& +! cnhpb,cihpb(cnhpb),cjhpb(cnhpb),cdhpb(cnhpb),cforcon(cnhpb) +!#else +! write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",& +! cnhpb,cihpb(cnhpb),cjhpb(cnhpb),cdhpb(cnhpb),cforcon(cnhpb) +!#endif +! enddo +! 1712 continue +! call flush(iout) +! return +! end subroutine read_dist_constr +! +!! if (ifrag_(1,i).lt.nstart_sup) ifrag_(1,i)=nstart_sup +!! if (ifrag_(2,i).gt.nstart_sup+nsup-1) & +!! ifrag_(2,i)=nstart_sup+nsup-1 +!! write (iout,*) i,ifrag_(1,i),ifrag_(2,i),wfrag_(i) +! subroutine gen_dist_constr2 +! use MPI_data +!! use control +! use geometry, only: dist +! use geometry_data +! use control_data +! use energy_data +! integer :: i,j +! real(kind=8) :: distance +! if (constr_dist.eq.11) then +! do i=nstart_sup,nstart_sup+nsup-1 +! do j=i+2,nstart_sup+nsup-1 +! distance=dist(i,j) +! if (distance.le.15.0) then +! cnhpb=cnhpb+1 +! cihpb(nhpb)=i+nstart_seq-nstart_sup +! cjhpb(nhpb)=j+nstart_seq-nstart_sup +! cforcon(nhpb)=sqrt(0.04*distance) +! cfordepth(nhpb)=sqrt(40.0/distance) +! cdhpb(nhpb)=distance-0.1d0 +! cdhpb1(nhpb)=distance+0.1d0 +! +!#ifdef MPI +! if (.not.out1file .or. me.eq.king) & +! write (iout,'(a,3i5,f8.2,f10.2)') "+dist.constr ", & +! cnhpb,cihpb(nhpb),cjhpb(nhpb),cdhpb(nhpb),cforcon(nhpb) +!#else +! write (iout,'(a,3i5,f8.2,f10.2)') "+dist.constr ", & +! cnhpb,cihpb(nhpb),cjhpb(nhpb),cdhpb(nhpb),cforcon(nhpb) +!#endif +! endif +! enddo +! enddo +! else +! do i=nstart_sup,nstart_sup+nsup-1 +! do j=i+2,nstart_sup+nsup-1 +! cnhpb=cnhpb+1 +! cihpb(nhpb)=i+nstart_seq-nstart_sup +! cjhpb(nhpb)=j+nstart_seq-nstart_sup +! cforcon(nhpb)=weidis +! cdhpb(nhpb)=dist(i,j) +! enddo +! enddo +! endif +! return +! end subroutine gen_dist_constr2 +! +!!----------------------------------------------------------------------------- end module io_clust diff --git a/source/cluster/main_clust.F b/source/cluster/main_clust.F index 650250e..be80144 100644 --- a/source/cluster/main_clust.F +++ b/source/cluster/main_clust.F @@ -371,7 +371,8 @@ c write (iout,*) "tu dochodze" nperm=i*nperm enddo c write (iout,*) "nperm",nperm - call permut(symetr) +! call permut(symetr) + call permut(nequiv(i),nperm,tabperm) c write (iout,*) "tabperm", tabperm(1,1) do kkk=1,nperm if (lside) then diff --git a/source/cluster/probabl.F90 b/source/cluster/probabl.F90 index 2f9c7d3..a732c05 100644 --- a/source/cluster/probabl.F90 +++ b/source/cluster/probabl.F90 @@ -166,7 +166,7 @@ enddo call int_from_cart1(.false.) ! call etotal(energia(0),fT) - call etotal(energia) + call etotal(energia(0)) totfree(i)=energia(0) ! write (iout,'(8f10.5)') ((c(l,k),l=1,3),k=1,nres) ! write (iout,'(8f10.5)') ((c(l,k+nres),l=1,3),k=nnt,nct) @@ -211,6 +211,7 @@ etors=enetb(13,i) etors_d=enetb(14,i) ehpb=enetb(15,i) +! print *,"ehpb",ehpb ! estr=enetb(18,i) estr=enetb(17,i) ! esccor=enetb(19,i) diff --git a/source/unres/CMakeLists.txt b/source/unres/CMakeLists.txt index 41d35fa..d9be5ea 100644 --- a/source/unres/CMakeLists.txt +++ b/source/unres/CMakeLists.txt @@ -82,7 +82,7 @@ if (Fortran_COMPILER_NAME STREQUAL "ifort") set (CMAKE_Fortran_FLAGS_RELEASE " ") set (CMAKE_Fortran_FLAGS_DEBUG "-O0 -g -traceback") # set(FFLAGS0 "-fpp -c -CB -g -ip " ) - set(FFLAGS0 "-O3 -ip -fpp -mcmodel=large" ) + set(FFLAGS0 "-CB -g -ip -fpp -mcmodel=large" ) # set(FFLAGS0 "-O0 -CB -CA -g" ) set(FFLAGS1 "-fpp -c -O " ) set(FFLAGS2 "-fpp -c -g -CA -CB ") diff --git a/source/unres/MD.F90 b/source/unres/MD.F90 index 52f96af..a676862 100644 --- a/source/unres/MD.F90 +++ b/source/unres/MD.F90 @@ -789,13 +789,13 @@ integer :: i,j,k,iti,mnum,term real(kind=8) :: KEt_p,KEt_sc,KEr_p,KEr_sc,incr(3),& mag1,mag2,v(3) -!#ifdef DEBUG +#ifdef DEBUG write (iout,*) "Velocities, kietic" do i=0,nres write (iout,'(i3,3f10.5,3x,3f10.5)') i,(d_t(j,i),j=1,3),& (d_t(j,i+nres),j=1,3) enddo -!#endif +#endif KEt_p=0.0d0 KEt_sc=0.0d0 ! write (iout,*) "ISC",(isc(itype(i,1)),i=1,nres) @@ -808,7 +808,7 @@ do i=nnt,term mnum=molnum(i) if (mnum.ge.5) mp(mnum)=msc(itype(i,mnum),mnum) - write (iout,*) "Kinetic trp:",i,(incr(j),j=1,3),mp(mnum) +! write (iout,*) "Kinetic trp:",i,(incr(j),j=1,3),mp(mnum) if (mnum.gt.4) then do j=1,3 v(j)=incr(j)+0.5d0*d_t(j,i) @@ -824,7 +824,7 @@ incr(j)=incr(j)+d_t(j,i) enddo enddo - write(iout,*) 'KEt_p', KEt_p +! write(iout,*) 'KEt_p', KEt_p ! The translational part for the side chain virtual bond ! Only now we can initialize incr with zeros. It must be equal ! to the velocities of the first Calpha. @@ -850,7 +850,7 @@ v(j)=incr(j)+d_t(j,nres+i) enddo endif - write (iout,*) "Kinetic trsc:",i,(incr(j),j=1,3) +! write (iout,*) "Kinetic trsc:",i,(incr(j),j=1,3) ! write (iout,*) "i",i," msc",msc(iti,mnum)," v",(v(j),j=1,3) KEt_sc=KEt_sc+mscab*(v(1)*v(1)+v(2)*v(2)+v(3)*v(3)) vtot(i+nres)=v(1)*v(1)+v(2)*v(2)+v(3)*v(3) @@ -859,7 +859,7 @@ enddo enddo ! goto 111 - write(iout,*) 'KEt_sc', KEt_sc +! write(iout,*) 'KEt_sc', KEt_sc ! The part due to stretching and rotation of the peptide groups KEr_p=0.0D0 do i=nnt,nct-1 @@ -869,12 +869,12 @@ do j=1,3 incr(j)=d_t(j,i) enddo - write (iout,*) "Kinetic rotp:",i,(incr(j),j=1,3),KEr_p,Ip(mnum) +! write (iout,*) "Kinetic rotp:",i,(incr(j),j=1,3),KEr_p,Ip(mnum) KEr_p=KEr_p+Ip(mnum)*(incr(1)*incr(1)+incr(2)*incr(2) & +incr(3)*incr(3)) enddo ! goto 111 - write(iout,*) 'KEr_p', KEr_p +! write(iout,*) 'KEr_p', KEr_p ! The rotational part of the side chain virtual bond KEr_sc=0.0D0 do i=nnt,nct @@ -886,16 +886,16 @@ do j=1,3 incr(j)=d_t(j,nres+i) enddo - write (iout,*) "Kinetic rotsc:",i,(incr(j),j=1,3) +! write (iout,*) "Kinetic rotsc:",i,(incr(j),j=1,3) KEr_sc=KEr_sc+Isc(iti,mnum)*(incr(1)*incr(1)+incr(2)*incr(2)+ & incr(3)*incr(3)) endif enddo ! The total kinetic energy 111 continue - write(iout,*) 'KEr_sc', KEr_sc +! write(iout,*) 'KEr_sc', KEr_sc KE_total=0.5d0*(KEt_p+KEt_sc+0.25d0*KEr_p+KEr_sc) - write (iout,*) "KE_total",KE_total +! write (iout,*) "KE_total",KE_total return end subroutine kinetic !----------------------------------------------------------------------------- @@ -3358,6 +3358,7 @@ innt=iposd_chain(ichain) ! innt_org= innt_org=chain_border(1,ichain) +! print *,"befor molnum0",innt_org if ((molnum(innt_org).eq.5).or.(molnum(innt_org).eq.4)) go to 137 if(.not.allocated(ghalf)) print *,"COCO" if(.not.allocated(Ghalf)) allocate(Ghalf(1300*(1300+1)/2)) @@ -3415,7 +3416,8 @@ write(iout,*) i,ghalf(i) enddo #endif - mnum=molnum(innt+1) +! print *,"befor molnum0.5",innt_org+1 + mnum=molnum(innt_org+1) call gldiag(1300*2,n,n,Ghalf,work,Geigen,Gvec,ierr,iwork) write(iout,*) "after internal",ierr,iwork #ifdef DEBUG @@ -3456,6 +3458,7 @@ do i=1,n do k=1,3 ii=ii+1 +! print *, "before molnum",innt_org+1 mnum=molnum(innt_org+1) if (molnum(innt_org+1).ge.4) geigen(i)=3.0/msc(itype(innt_org+i-1,mnum),mnum) !should not it be other way around @@ -3522,16 +3525,19 @@ do j=1,3 ind=ind+1 d_t(j,i)=d_t_work(ind) - enddo + enddo !j +! print *, "before molnum2",i mnum=molnum(i) + print *, "after molnum2",i if (itype(i,1).ne.10 .and. itype(i,mnum).ne.ntyp1_molec(mnum).and.mnum.le.2) then do j=1,3 ind=ind+1 d_t(j,i+nres)=d_t_work(ind) - enddo + enddo !j endif - enddo + enddo !i enddo + print *,"before large1" if (large) then write (iout,*) write (iout,*) "Random velocities in the Calpha,SC space" @@ -3541,7 +3547,10 @@ restyp(itype(i,mnum),mnum),i,(d_t(j,i),j=1,3),(d_t(j,i+nres),j=1,3) enddo endif + print *,"before kineticCASC" call kinetic_CASC(Ek1) + print *,"after kineticCASC" + ! ! Transform the velocities to virtual-bond space ! @@ -3636,6 +3645,7 @@ enddo ibeg=inct+1 do i=innt,inct +! print *,"before molnum3",i mnum=molnum(i) if (iabs(itype(i,1).eq.10).or.mnum.ge.3) then !c write (iout,*) "i",i,(d_t(j,i),j=1,3),(d_t(j,i+1),j=1,3) @@ -4491,7 +4501,7 @@ Im(2,2)=Im(2,2)+msc(iabs(iti),mnum)*(pr(3)*pr(3)+pr(1)*pr(1)) Im(3,3)=Im(3,3)+msc(iabs(iti),mnum)*(pr(1)*pr(1)+pr(2)*pr(2)) enddo - print *,"after msc",Im +! print *,"after msc",Im do i=nnt,nct-1 mnum=molnum(i) mnum1=molnum(i+1) @@ -4510,7 +4520,7 @@ Im(3,3)=Im(3,3)+Ip(mnum)*(1-dc_norm(3,i)*dc_norm(3,i))*& vbld(i+1)*vbld(i+1)*0.25d0 enddo - print *,"after Ip",Im +! print *,"after Ip",Im do i=nnt,nct mnum=molnum(i) mnum1=molnum(i+1) @@ -4529,12 +4539,12 @@ dc_norm(2,inres))*vbld(inres)*vbld(inres) Im(3,3)=Im(3,3)+Isc(iti,mnum)*(1-dc_norm(3,inres)*& dc_norm(3,inres))*vbld(inres)*vbld(inres) - print *,i,inres,vbld(inres),iti,dc_norm(1,inres),Im(1,1) +! print *,i,inres,vbld(inres),iti,dc_norm(1,inres),Im(1,1) endif enddo - print *,"after ISC",Im - print *,"before angnom",L - print *,"before angnom",Im +! print *,"after ISC",Im +! print *,"before angnom",L +! print *,"before angnom",Im call angmom(cm,L) Im(2,1)=Im(1,2) Im(3,1)=Im(1,3) @@ -4548,9 +4558,9 @@ enddo enddo !c Finding the eigenvectors and eignvalues of the inertia tensor - print *,"before djacob",Imcp,eigvec,eigval +! print *,"before djacob",Imcp,eigvec,eigval call djacob(3,3,10000,1.0d-10,Imcp,eigvec,eigval) - print *, "after djacob" +! print *, "after djacob" do i=1,3 if (dabs(eigval(i)).gt.1.0d-15) then Id(i,i)=1.0d0/eigval(i) diff --git a/source/unres/data/names.F90 b/source/unres/data/names.F90 index 32d7336..988be12 100644 --- a/source/unres/data/names.F90 +++ b/source/unres/data/names.F90 @@ -106,7 +106,7 @@ "WLT "," "," ","WTUBE " ,& "WVDWPPNUCL","WELPP ","WVDWPSB ","WELPSB ","WVDWSB ",& "WELSB ","WBOND_NUCL","WANG_NUCL ","WSBLOC ","WTOR_NUCL ",& - "WTORD_NUCL","WCORR_NUCL","WCORR3_NUC","WNULL ","WNULL ",& + "WTORD_NUCL","WCORR_NUCL","WCORR3NUCL","WNULL ","WNULL ",& "WCATPROT ","WCATCAT ","WNULL ","WNULL ","WNULL ",& "WSCBASE ","WPEPBASE ","WSCPHO ","WPEPPHO ","WCATNUCL ",& "H_CONS ","WNULL ","WNULL ","WNULL ","WNULL ",& diff --git a/source/unres/energy.F90 b/source/unres/energy.F90 index c14e4ad..ef82086 100644 --- a/source/unres/energy.F90 +++ b/source/unres/energy.F90 @@ -929,6 +929,7 @@ ! print *,"before",ees,evdw1,ecorr ! write(iout,*) ecorr_nucl,"ecorr_nucl",nres_molec(2) if (nres_molec(2).gt.0) then +! write(iout,*) "welpsb2,",welpsb call ebond_nucl(estr_nucl) call ebend_nucl(ebe_nucl) call etor_nucl(etors_nucl) @@ -1400,7 +1401,7 @@ wtor=weights(13)*fact(1) wtor_d=weights(14)*fact(2) wsccor=weights(21)*fact(1) - welpsb=weights(28)*fact(1) + welpsb=weights(29)*fact(1) wcorr_nucl= weights(37)*fact(1) wcorr3_nucl=weights(38)*fact(2) wtor_nucl= weights(35)*fact(1) @@ -5803,16 +5804,16 @@ !c Restraints from contact prediction dd=dist(ii,jj) if (constr_dist.eq.11) then - ehpb=ehpb+fordepth(i)**4.0d0 & + ehpb=ehpb+cfordepth(i)**4.0d0 & *rlornmr1(dd,cdhpb(i),cdhpb1(i),cforcon(i)) - fac=fordepth(i)**4.0d0 & + fac=cfordepth(i)**4.0d0 & *rlornmr1prim(dd,cdhpb(i),cdhpb1(i),cforcon(i))/dd if (energy_dec) write (iout,'(a6,2i5,3f8.3)') "edisl",ii,jj, & ehpb,cfordepth(i),dd else if (cdhpb1(i).gt.0.0d0) then ehpb=ehpb+2*cforcon(i)*gnmr1(dd,cdhpb(i),cdhpb1(i)) - fac=forcon(i)*gnmr1prim(dd,cdhpb(i),cdhpb1(i))/dd + fac=cforcon(i)*gnmr1prim(dd,cdhpb(i),cdhpb1(i))/dd !c write (iout,*) "beta nmr", !c & dd,2*forcon(i)*gnmr1(dd,dhpb(i),dhpb1(i)) else @@ -5845,15 +5846,15 @@ if (constr_dist.eq.11) then ehpb=ehpb+cfordepth(i)**4.0d0 & - *rlornmr1(dd,dhpb(i),cdhpb1(i),cforcon(i)) - fac=fordepth(i)**4.0d0 & + *rlornmr1(dd,cdhpb(i),cdhpb1(i),cforcon(i)) + fac=cfordepth(i)**4.0d0 & *rlornmr1prim(dd,cdhpb(i),cdhpb1(i),cforcon(i))/dd if (energy_dec) write (iout,'(a6,2i5,3f8.3)') "edisl",ii,jj, & - ehpb,fordepth(i),dd + ehpb,cfordepth(i),dd else if (cdhpb1(i).gt.0.0d0) then ehpb=ehpb+2*cforcon(i)*gnmr1(dd,cdhpb(i),cdhpb1(i)) - fac=forcon(i)*gnmr1prim(dd,cdhpb(i),cdhpb1(i))/dd + fac=cforcon(i)*gnmr1prim(dd,cdhpb(i),cdhpb1(i))/dd !c write (iout,*) "alph nmr", !c & dd,2*forcon(i)*gnmr1(dd,dhpb(i),dhpb1(i)) else @@ -5877,9 +5878,9 @@ waga=cforcon(i) !C Calculate the contribution to energy. ehpb=ehpb+waga*rdis*rdis - ! if (energy_dec) -! write (iout,'(a6,2i5,5f10.3)') "edisip",ii,jj, & -! ehpb,dd,cdhpb(i),waga,rdis + if (energy_dec) & + write (iout,'(a6,2i5,5f10.3)') "edisip",ii,jj, & + ehpb,dd,cdhpb(i),waga,rdis !c write (iout,*) "alpha reg",dd,waga*rdis*rdis !C @@ -6961,6 +6962,9 @@ !c AL 3/30/2022 handle the cases of an isolated-residue chain if (i.eq.nnt .and. itype(i+1,1).eq.ntyp1) cycle if (i.eq.nct .and. itype(i-1,1).eq.ntyp1) cycle + if (itype(i-1,1).eq.ntyp1 .and. itype(i+1,1).eq.ntyp1) cycle +! if (i.eq.nct .and. itype(i-1,1).eq.ntyp1) cycle + ! costtab(i+1) =dcos(theta(i+1)) if (it.eq.10) goto 1 #ifdef SC_END @@ -6977,10 +6981,10 @@ *cossc1)*gradene gsclocx(:,i)=gsclocx(:,i)+vbld_inv(i+nres)*& (dC_norm(:,i-1)-dC_norm(:,i+nres)*cossc1)*gradene -#ifdef ENERGY_DEC +!#ifdef ENERGY_DEC if (energy_dec) write (2,'(2hC ,a3,i6,2(a,f10.5))')& restyp(iti,1),i," angle",rad2deg*dacos(cossc1)," escloc",sumene -#endif +!#endif else if (i.eq.nnt .or. itype(i-1,1).eq.ntyp1) then !c AL 3/30/2022 handle a sidechain of a loose N-end cossc=scalar(dc_norm(1,i),dc_norm(1,i+nres)) @@ -6994,10 +6998,10 @@ *cossc)*gradene gsclocx(:,i)=gsclocx(:,i)+vbld_inv(i+nres)*& (dC_norm(:,i)-dC_norm(:,i+nres)*cossc)*gradene -#ifdef ENERGY_DEC - if (energy_dec) write (2,'(2hN ,a3,i6,2(a,f10.5))') - & restyp(iti),i," angle",rad2deg*dacos(cossc)," escloc",sumene -#endif +!#ifdef ENERGY_DEC + if (energy_dec) write (2,'(2hN ,a3,i6,2(a,f10.5))')& + restyp(iti,1),i," angle",rad2deg*dacos(cossc)," escloc",sumene +!#endif else #endif ! @@ -21652,6 +21656,7 @@ C !!!!!!!! NIE CZYTANE !!!!!!!!!!! !(maxres) allocate(jcont_hb(maxconts,nres)) #endif + allocate(jcont_hb(maxconts,nres)) allocate(num_cont_hb(nres)) !(maxconts,maxres) ! common /rotat/ @@ -22678,7 +22683,7 @@ C !!!!!!!! NIE CZYTANE !!!!!!!!!!! ! print *,"tu?",i,istart_nucl(i,iint),iend_nucl(i,iint) do j=istart_nucl(i,iint),iend_nucl(i,iint) ind=ind+1 -! print *,"JESTEM" +! print *,"JESTEM",i,j itypj=itype(j,2) if (itypj.eq.ntyp1_molec(2)) cycle dscj_inv=vbld_inv(j+nres) @@ -22770,9 +22775,9 @@ C !!!!!!!! NIE CZYTANE !!!!!!!!!!! !C Calculate angular part of the gradient. call sc_grad_nucl call eelsbij(eelij,num_conti2) - if (energy_dec .and. & - (j.eq.i+1.or.j.eq.nres-i+1.or.j.eq.nres-i.or.j.eq.nres-i+2)) & - write (istat,'(e14.5)') evdwij +! if (energy_dec .and. & +! (j.eq.i+1.or.j.eq.nres-i+1.or.j.eq.nres-i.or.j.eq.nres-i+2)) & +! write (istat,'(e14.5)') evdwij eelsb=eelsb+eelij enddo ! j enddo ! iint @@ -22855,14 +22860,14 @@ C !!!!!!!! NIE CZYTANE !!!!!!!!!!! !C 12/26/95 - for the evaluation of multi-body H-bonding interactions ees0ij=4.0D0+facfac-fac1 - if (energy_dec) then - if(j.eq.i+1.or.j.eq.nres-i+1.or.j.eq.nres-i.or.j.eq.nres-i+2) & - write (istat,'(2a1,i4,1x,2a1,i4,4f10.5,3e12.5,$)') & - sugartyp(istype(i)),restyp(itypi,2),i,sugartyp(istype(j)),& - restyp(itypj,2),j,1.0d0/rij,cosa,cosb,cosg,fac*r3ij, & - (4.0D0+facfac-fac1)*r6ij,(2.0d0-2.0d0*facfac+fac1)*r6ij - write (iout,'(a6,2i5,e15.3)') 'ees',i,j,eesij - endif +! if (energy_dec) then +! if(j.eq.i+1.or.j.eq.nres-i+1.or.j.eq.nres-i.or.j.eq.nres-i+2) & +! write (istat,'(2a1,i4,1x,2a1,i4,4f10.5,3e12.5,$)') & +! sugartyp(istype(i)),restyp(itypi,2),i,sugartyp(istype(j)),& +! restyp(itypj,2),j,1.0d0/rij,cosa,cosb,cosg,fac*r3ij, & +! (4.0D0+facfac-fac1)*r6ij,(2.0d0-2.0d0*facfac+fac1)*r6ij +! write (iout,'(a6,2i5,e15.3)') 'ees',i,j,eesij +! endif !C !C Calculate contributions to the Cartesian gradient. @@ -22913,9 +22918,8 @@ C !!!!!!!! NIE CZYTANE !!!!!!!!!!! gelsbc(k,j)=gelsbc(k,j)+ggg(k) gelsbc(k,i)=gelsbc(k,i)-ggg(k) enddo -! IF ( (wcorr_nucl.gt.0.0d0.or.wcorr3_nucl.gt.0.0d0) .and. - IF ( j.gt.i+1 .and.& - num_conti.le.maxcont) THEN + IF ( (wcorr_nucl.gt.0.0d0.or.wcorr3_nucl.gt.0.0d0) .and. & + ( j.gt.i+1 .and. num_conti2.le.maxcont)) THEN !C !C Calculate the contact function. The ith column of the array JCONT will !C contain the numbers of atoms that make contacts with the atom I (of numbers @@ -22929,11 +22933,11 @@ C !!!!!!!! NIE CZYTANE !!!!!!!!!!! num_conti=num_conti+1 num_conti2=num_conti2+1 - if (num_conti.gt.maxconts) then + if (num_conti2.gt.maxconts) then write (iout,*) 'WARNING - max. # of contacts exceeded;',& ' will skip next contacts for this conf.',maxconts else - jcont_hb(num_conti,i)=j + jcont_hb(num_conti2,i)=j !c write (iout,*) "num_conti",num_conti, !c & " jcont_hb",jcont_hb(num_conti,i) !C Calculate contact energies @@ -22977,7 +22981,7 @@ C !!!!!!!! NIE CZYTANE !!!!!!!!!!! ecosbm=ecosb1-ecosb2 ecosgm=ecosg1-ecosg2 !C End diagnostics - facont_hb(num_conti,i)=fcont + facont_hb(num_conti2,i)=fcont fprimcont=fprimcont/rij do k=1,3 gggp(k)=ecosbp*dcosb(k)+ecosgp*dcosg(k) @@ -22990,27 +22994,27 @@ C !!!!!!!! NIE CZYTANE !!!!!!!!!!! gggm(2)=gggm(2)+ees0mijp*yj gggm(3)=gggm(3)+ees0mijp*zj !C Derivatives due to the contact function - gacont_hbr(1,num_conti,i)=fprimcont*xj - gacont_hbr(2,num_conti,i)=fprimcont*yj - gacont_hbr(3,num_conti,i)=fprimcont*zj + gacont_hbr(1,num_conti2,i)=fprimcont*xj + gacont_hbr(2,num_conti2,i)=fprimcont*yj + gacont_hbr(3,num_conti2,i)=fprimcont*zj do k=1,3 !c !c Gradient of the correlation terms !c - gacontp_hb1(k,num_conti,i)= & + gacontp_hb1(k,num_conti2,i)= & (ecosap*(dc_norm(k,j+nres)-cosa*dc_norm(k,i+nres)) & + ecosbp*(erij(k)-cosb*dc_norm(k,i+nres)))*vbld_inv(i+nres) - gacontp_hb2(k,num_conti,i)= & + gacontp_hb2(k,num_conti2,i)= & (ecosap*(dc_norm(k,i+nres)-cosa*dc_norm(k,j+nres)) & + ecosgp*(erij(k)-cosg*dc_norm(k,j+nres)))*vbld_inv(j+nres) - gacontp_hb3(k,num_conti,i)=gggp(k) - gacontm_hb1(k,num_conti,i)= & + gacontp_hb3(k,num_conti2,i)=gggp(k) + gacontm_hb1(k,num_conti2,i)= & (ecosam*(dc_norm(k,j+nres)-cosa*dc_norm(k,i+nres)) & + ecosbm*(erij(k)-cosb*dc_norm(k,i+nres)))*vbld_inv(i+nres) - gacontm_hb2(k,num_conti,i)= & + gacontm_hb2(k,num_conti2,i)= & (ecosam*(dc_norm(k,i+nres)-cosa*dc_norm(k,j+nres))& + ecosgm*(erij(k)-cosg*dc_norm(k,j+nres)))*vbld_inv(j+nres) - gacontm_hb3(k,num_conti,i)=gggm(k) + gacontm_hb3(k,num_conti2,i)=gggm(k) enddo endif endif diff --git a/source/unres/io.F90 b/source/unres/io.F90 index 64b4ede..548802c 100644 --- a/source/unres/io.F90 +++ b/source/unres/io.F90 @@ -743,7 +743,7 @@ ! double precision x(maxvar) character(len=256) :: pdbfile character(len=800) :: weightcard - character(len=80) :: weightcard_t!,ucase + character(len=80) :: weightcard_t,chartest!,ucase ! integer,dimension(:),allocatable :: itype_pdb !(maxres) ! common /pizda/ itype_pdb logical :: fail !seq_comp, @@ -915,6 +915,7 @@ 'Scaling factor of 1,4 SC-p interactions:',scal14 write (iout,'(a,f8.3)') & 'General scaling factor of SC-p interactions:',scalscp + write(iout,*) "welpsb",welpsb endif r0_corr=cutoff_corr-delt_corr do i=1,ntyp @@ -1050,7 +1051,7 @@ enddo endif endif - print *,"CZY TU DOCHODZE" + print *,"CZY TU DOCHODZE",nucleic if (indpdb.eq.0) then nres_molec(:)=0 allocate(sequence(maxres,5)) @@ -1083,6 +1084,7 @@ if (nucleic) then ! Read sequence if not taken from the pdb file. molec=2 + print *,"jest before read nres" read (inp,*) nres_molec(molec) print *,'nres=',nres_molec(molec) ! allocate(sequence(maxres,5)) @@ -1654,11 +1656,11 @@ enddo do i=2,nres-1 alph(i)=110d0*deg2rad - if (molnum(i).eq.2) alph(i)=30d0*deg2rad + if (molnum(i).eq.2) alph(i)=60d0*deg2rad enddo do i=2,nres-1 omeg(i)=-120d0*deg2rad - if (molnum(i).eq.2) omeg(i)=60d0*deg2rad + if (molnum(i).eq.2) omeg(i)=30d0*deg2rad if (itype(i,1).le.0) omeg(i)=-omeg(i) enddo call chainbuild diff --git a/source/unres/io_base.F90 b/source/unres/io_base.F90 index a77c487..36799c9 100644 --- a/source/unres/io_base.F90 +++ b/source/unres/io_base.F90 @@ -658,20 +658,21 @@ read (inp,*) cihpb(cnhpb+1),cjhpb(cnhpb+1),cdhpb(i),cdhpb1(i), & cibecarb(i),cforcon(cnhpb+1),cfordepth(cnhpb+1) !EL fordepth(nhpb+1)=fordepth(nhpb+1)/forcon(nhpb+1) - fordepth(nhpb+1)=cfordepth(cnhpb+1)**(0.25d0) - forcon(nhpb+1)=cforcon(cnhpb+1)**(0.25d0) + cfordepth(cnhpb+1)=cfordepth(cnhpb+1)**(0.25d0) + cforcon(cnhpb+1)=cforcon(cnhpb+1)**(0.25d0) else -!C print *,"in else" + print *,"in else" read (inp,*) cihpb(cnhpb+1),cjhpb(cnhpb+1),cdhpb(i),cdhpb1(i), & cibecarb(i),cforcon(cnhpb+1) + print *,"afterread",cihpb(cnhpb+1) endif - if (cforcon(nhpb+1).gt.0.0d0) then + if (cforcon(cnhpb+1).gt.0.0d0) then cnhpb=cnhpb+1 - if (ibecarb(i).gt.0) then - cihpb(i)=ihpb(i)+nres - cjhpb(i)=jhpb(i)+nres + if (cibecarb(i).gt.0) then + cihpb(cnhpb)=cihpb(cnhpb)+nres + cjhpb(cnhpb)=cjhpb(cnhpb)+nres endif - if (cdhpb(nhpb).eq.0.0d0) & + if (cdhpb(cnhpb).eq.0.0d0) & cdhpb(cnhpb)=dist(ihpb(cnhpb),jhpb(cnhpb)) endif @@ -1145,7 +1146,7 @@ do i=nnt,nct ! write(iout,*), "coord",c(1,i),c(2,i),c(3,i) iti=itype(i,molnum(i)) - print *,i,molnum(i) +! print *,i,molnum(i) if (molnum(i+1).eq.0) then iti1=ntyp1_molec(molnum(i)) else diff --git a/source/unres/math.F90 b/source/unres/math.F90 index f6a800a..f78856d 100644 --- a/source/unres/math.F90 +++ b/source/unres/math.F90 @@ -28,7 +28,7 @@ COS45 = .70710678 S45SQ = 0.50 C45SQ = 0.50 - print *,"in math",N,NMAX +! print *,"in math",N,NMAX ! UNIT EIGENVECTOR MATRIX DO 70 I = 1,N DO 7 J = I,N @@ -102,7 +102,7 @@ 1000 FORMAT (/1X,'NONCONVERGENT JACOBI. LARGEST OFF-DIAGONAL ELE',& 'MENT = ',1PE14.7) ! ARRANGEMENT OF EIGENVALUES IN ASCENDING ORDER - print *,"before first AJJ?" +! print *,"before first AJJ?" 300 DO 14 I=1,N 14 AJJ(I)=A(I,I) LT=N+1 @@ -114,7 +114,7 @@ 17 AIIMIN=AJJ(I) IT=I 16 CONTINUE - print *,IT,"IT" +! print *,IT,"IT" IN=L AII(IN)=AIIMIN AJJ(IT)=1.0E+30 diff --git a/source/unres/minim.F90 b/source/unres/minim.F90 index c661960..9016aa7 100644 --- a/source/unres/minim.F90 +++ b/source/unres/minim.F90 @@ -6854,6 +6854,7 @@ !c set default parameters for the line search !c if (stpmax .eq. 0.0d0) stpmax = 5.0d0 + if (stpmax .eq. 0.0d0) stpmax = 2.0d0 stpmin = 1.0d-16 cappa = 0.9d0 slpmax = 10000.0d0 diff --git a/source/unres/search.F90 b/source/unres/search.F90 index 4439b32..89ee892 100644 --- a/source/unres/search.F90 +++ b/source/unres/search.F90 @@ -34,6 +34,7 @@ blank = ' ' if (stpmin .eq. 0.0d0) stpmin = 1.0d-16 if (stpmax .eq. 0.0d0) stpmax = 2.0d0 + if (stpmax .eq. 0.0d0) stpmax = 1.5d0 if (cappa .eq. 0.0d0) cappa = 0.1d0 if (slpmax .eq. 0.0d0) slpmax = 10000.0d0 if (angmax .eq. 0.0d0) angmax = 180.0d0 diff --git a/source/wham/enecalc.F90 b/source/wham/enecalc.F90 index 6ba59a0..0412693 100644 --- a/source/wham/enecalc.F90 +++ b/source/wham/enecalc.F90 @@ -136,7 +136,7 @@ integer :: errmsg_count,maxerrmsg_count=100000 !el real(kind=8) :: rmsnat,gyrate !el external rmsnat,gyrate - real(kind=8) :: tole=0.0d0 + real(kind=8) :: tole=0.5d0 integer i,itj,ii,iii,j,k,l,licz integer ir,ib,ipar,iparm integer iscor,islice @@ -301,6 +301,8 @@ write(iout,*)"enecalc_ i ntot",i,ntot wtor_d,wsccor,wbond #endif call restore_parm(iparm) + print *,"after restore parm" + print *,beta_h(ib,ipar) call rescale_weights(1.0d0/(beta_h(ib,ipar)*1.987D-3)) #ifdef DEBUG write (iout,*) "before etot w=",1.0d0/(beta_h(ib,ipar)*1.987D-3) @@ -309,7 +311,9 @@ write(iout,*)"enecalc_ i ntot",i,ntot wtor_d,wsccor,wbond #endif ! call etotal(energia(0),fT) + print *,"before etotal" call etotal(energia(0)) + print *,"after etotal" !write(iout,*)"check c and dc after etotal",1.0d0/(0.001987*beta_h(ib,ipar)) !do k=1,2*nres+2 !write(iout,*)k,"c=",(c(l,k),l=1,3) diff --git a/source/wham/io_wham.F90 b/source/wham/io_wham.F90 index d8cb743..e15f2fd 100644 --- a/source/wham/io_wham.F90 +++ b/source/wham/io_wham.F90 @@ -167,7 +167,7 @@ ! use energy_data use geometry_data, only:nres,deg2rad,c,dc,nres_molec,crefjlee,cref - use control_data, only:iscode,dyn_ss,pdbref,indpdb + use control_data, only:iscode,pdbref,indpdb use io_base, only:rescode use control, only:setup_var,init_int_table,hpb_partition use geometry, only:alloc_geo_arrays @@ -496,7 +496,9 @@ allocate(ww(max_eneW)) do i=1,n_eneW key = wname(i)(:ilen(wname(i))) + call reada(controlcard,key(:ilen(key)),ww(i),1.0d0) + write(iout,*) "NAMES ",key(:ilen(key)),ww(i),i enddo write (iout,*) "iparm",iparm," myparm",myparm @@ -528,6 +530,7 @@ allocate(ww(max_eneW)) wcatcat=ww(42) wcatprot=ww(41) wcorr3_nucl=ww(38) + write(iout,*) "CZY TU BLAD", ww(38) wcorr_nucl=ww(37) wtor_d_nucl=ww(36) wtor_nucl=ww(35) @@ -3137,7 +3140,7 @@ allocate(ww(max_eneW)) use control_data, only:indpdb,symetr,r_cut_ele,rlamb_ele,ions,& scelemode,TUBEmode,tor_mode,energy_dec,r_cut_ang,r_cut_mart,& - rlamb_mart + rlamb_mart,constr_dist use energy_data, only:distchainmax,tubeR0,tubecenter,dyn_ss,constr_homology use geometry_data, only:boxxsize,boxysize,boxzsize,bordtubetop,& diff --git a/source/wham/wham.F90 b/source/wham/wham.F90 index 78304d8..be3f2bf 100644 --- a/source/wham/wham.F90 +++ b/source/wham/wham.F90 @@ -10,7 +10,7 @@ use work_part ! use io_units - use control_data, only:indpdb + use control_data, only:indpdb,constr_dist #ifdef MPI use mpi_data ! use mpi_ diff --git a/source/wham/wham_data.F90 b/source/wham/wham_data.F90 index dc0a0cc..6192eb0 100644 --- a/source/wham/wham_data.F90 +++ b/source/wham/wham_data.F90 @@ -109,7 +109,7 @@ logical :: punch_dist,print_rms,caonly,verbose,merge_helices,& bxfile,cxfile,histfile,entfile,zscfile,rmsrgymap,& with_dihed_constr,check_conf,histout - integer :: icomparfunc,pdbint,ensembles,constr_dist,oldion,shield_mode + integer :: icomparfunc,pdbint,ensembles,oldion,shield_mode !--------------------------------------------------------------------------- ! COMMON.OBCINKA ! common /obcinka/