ifdef poporawa
[unres4.git] / source / unres / io_config.f90
index 0dce2e1..0c2d645 100644 (file)
@@ -30,7 +30,7 @@
 !
 !-----------------------------------------------------------------------------
       contains
-#if .not. defined(WHAM_RUN) && .not. defined(CLUSTER)
+#if !defined(WHAM_RUN) && !defined(CLUSTER)
 !-----------------------------------------------------------------------------
 ! bank.F    io_csa
 !-----------------------------------------------------------------------------
       character(len=1) :: toronelet(-2:2) = (/"p","a","G","A","P"/)
       logical :: lprint,LaTeX
       real(kind=8),dimension(3,3,maxlob) :: blower     !(3,3,maxlob)
-      real(kind=8),dimension(13) :: b
+      real(kind=8),dimension(13) :: bN
       character(len=3) :: lancuch      !,ucase
 !el  local variables
       integer :: m,n,l,i,j,k,iblock,lll,llll,ll,nlobi,mm
       real(kind=8) :: v0ijsccor,v0ijsccor1,v0ijsccor2,v0ijsccor3,si,&
                 dwa16,rjunk,akl,v0ij,rri,epsij,rrij,sigeps,sigt1sq,&
                 sigt2sq,sigii1,sigii2,ratsig1,ratsig2,rsum_max,r_augm,&
-                res1,epsijlip
+                res1,epsijlip,epspeptube,epssctube,sigmapeptube,      &
+                sigmasctube
       integer :: ichir1,ichir2
 !      real(kind=8),dimension(maxterm,-maxtor:maxtor,-maxtor:maxtor,2) :: v1_el,v2_el !(maxterm,-maxtor:maxtor,-maxtor:maxtor,2)
 !el      allocate(v1_el(maxterm,-maxtor:maxtor,-maxtor:maxtor,2))
 
       do i=0,nloctyp-1
         read (ifourier,*,end=115,err=115)
-        read (ifourier,*,end=115,err=115) (b(ii),ii=1,13)
+        read (ifourier,*,end=115,err=115) (bN(ii),ii=1,13)
         if (lprint) then
           write (iout,*) 'Type',i
-          write (iout,'(a,i2,a,f10.5)') ('b(',ii,')=',b(ii),ii=1,13)
+          write (iout,'(a,i2,a,f10.5)') ('bN(',ii,')=',bN(ii),ii=1,13)
         endif
-        B1(1,i)  = b(3)
-        B1(2,i)  = b(5)
-        B1(1,-i) = b(3)
-        B1(2,-i) = -b(5)
+        B1(1,i)  = bN(3)
+        B1(2,i)  = bN(5)
+        B1(1,-i) = bN(3)
+        B1(2,-i) = -bN(5)
 !        b1(1,i)=0.0d0
 !        b1(2,i)=0.0d0
-        B1tilde(1,i) = b(3)
-        B1tilde(2,i) =-b(5)
-        B1tilde(1,-i) =-b(3)
-        B1tilde(2,-i) =b(5)
+        B1tilde(1,i) = bN(3)
+        B1tilde(2,i) =-bN(5)
+        B1tilde(1,-i) =-bN(3)
+        B1tilde(2,-i) =bN(5)
 !        b1tilde(1,i)=0.0d0
 !        b1tilde(2,i)=0.0d0
-        B2(1,i)  = b(2)
-        B2(2,i)  = b(4)
-        B2(1,-i)  =b(2)
-        B2(2,-i)  =-b(4)
+        B2(1,i)  = bN(2)
+        B2(2,i)  = bN(4)
+        B2(1,-i)  =bN(2)
+        B2(2,-i)  =-bN(4)
 
 !        b2(1,i)=0.0d0
 !        b2(2,i)=0.0d0
-        CC(1,1,i)= b(7)
-        CC(2,2,i)=-b(7)
-        CC(2,1,i)= b(9)
-        CC(1,2,i)= b(9)
-        CC(1,1,-i)= b(7)
-        CC(2,2,-i)=-b(7)
-        CC(2,1,-i)=-b(9)
-        CC(1,2,-i)=-b(9)
+        CC(1,1,i)= bN(7)
+        CC(2,2,i)=-bN(7)
+        CC(2,1,i)= bN(9)
+        CC(1,2,i)= bN(9)
+        CC(1,1,-i)= bN(7)
+        CC(2,2,-i)=-bN(7)
+        CC(2,1,-i)=-bN(9)
+        CC(1,2,-i)=-bN(9)
 !        CC(1,1,i)=0.0d0
 !        CC(2,2,i)=0.0d0
 !        CC(2,1,i)=0.0d0
 !        CC(1,2,i)=0.0d0
-        Ctilde(1,1,i)=b(7)
-        Ctilde(1,2,i)=b(9)
-        Ctilde(2,1,i)=-b(9)
-        Ctilde(2,2,i)=b(7)
-        Ctilde(1,1,-i)=b(7)
-        Ctilde(1,2,-i)=-b(9)
-        Ctilde(2,1,-i)=b(9)
-        Ctilde(2,2,-i)=b(7)
+        Ctilde(1,1,i)=bN(7)
+        Ctilde(1,2,i)=bN(9)
+        Ctilde(2,1,i)=-bN(9)
+        Ctilde(2,2,i)=bN(7)
+        Ctilde(1,1,-i)=bN(7)
+        Ctilde(1,2,-i)=-bN(9)
+        Ctilde(2,1,-i)=bN(9)
+        Ctilde(2,2,-i)=bN(7)
 
 !        Ctilde(1,1,i)=0.0d0
 !        Ctilde(1,2,i)=0.0d0
 !        Ctilde(2,1,i)=0.0d0
 !        Ctilde(2,2,i)=0.0d0
-        DD(1,1,i)= b(6)
-        DD(2,2,i)=-b(6)
-        DD(2,1,i)= b(8)
-        DD(1,2,i)= b(8)
-        DD(1,1,-i)= b(6)
-        DD(2,2,-i)=-b(6)
-        DD(2,1,-i)=-b(8)
-        DD(1,2,-i)=-b(8)
+        DD(1,1,i)= bN(6)
+        DD(2,2,i)=-bN(6)
+        DD(2,1,i)= bN(8)
+        DD(1,2,i)= bN(8)
+        DD(1,1,-i)= bN(6)
+        DD(2,2,-i)=-bN(6)
+        DD(2,1,-i)=-bN(8)
+        DD(1,2,-i)=-bN(8)
 !        DD(1,1,i)=0.0d0
 !        DD(2,2,i)=0.0d0
 !        DD(2,1,i)=0.0d0
 !        DD(1,2,i)=0.0d0
-        Dtilde(1,1,i)=b(6)
-        Dtilde(1,2,i)=b(8)
-        Dtilde(2,1,i)=-b(8)
-        Dtilde(2,2,i)=b(6)
-        Dtilde(1,1,-i)=b(6)
-        Dtilde(1,2,-i)=-b(8)
-        Dtilde(2,1,-i)=b(8)
-        Dtilde(2,2,-i)=b(6)
+        Dtilde(1,1,i)=bN(6)
+        Dtilde(1,2,i)=bN(8)
+        Dtilde(2,1,i)=-bN(8)
+        Dtilde(2,2,i)=bN(6)
+        Dtilde(1,1,-i)=bN(6)
+        Dtilde(1,2,-i)=-bN(8)
+        Dtilde(2,1,-i)=bN(8)
+        Dtilde(2,2,-i)=bN(6)
 
 !        Dtilde(1,1,i)=0.0d0
 !        Dtilde(1,2,i)=0.0d0
 !        Dtilde(2,1,i)=0.0d0
 !        Dtilde(2,2,i)=0.0d0
-        EE(1,1,i)= b(10)+b(11)
-        EE(2,2,i)=-b(10)+b(11)
-        EE(2,1,i)= b(12)-b(13)
-        EE(1,2,i)= b(12)+b(13)
-        EE(1,1,-i)= b(10)+b(11)
-        EE(2,2,-i)=-b(10)+b(11)
-        EE(2,1,-i)=-b(12)+b(13)
-        EE(1,2,-i)=-b(12)-b(13)
+        EE(1,1,i)= bN(10)+bN(11)
+        EE(2,2,i)=-bN(10)+bN(11)
+        EE(2,1,i)= bN(12)-bN(13)
+        EE(1,2,i)= bN(12)+bN(13)
+        EE(1,1,-i)= bN(10)+bN(11)
+        EE(2,2,-i)=-bN(10)+bN(11)
+        EE(2,1,-i)=-bN(12)+bN(13)
+        EE(1,2,-i)=-bN(12)-bN(13)
 
 !        ee(1,1,i)=1.0d0
 !        ee(2,2,i)=1.0d0
       allocate(aa_aq(ntyp1,ntyp1),bb_aq(ntyp1,ntyp1),chi(ntyp1,ntyp1)) !(ntyp,ntyp)
       allocate(aa_lip(ntyp1,ntyp1),bb_lip(ntyp1,ntyp1)) !(ntyp,ntyp)
       allocate(sigma(0:ntyp1,0:ntyp1),r0(ntyp1,ntyp1)) !(0:ntyp1,0:ntyp1)
+      allocate(acavtub(ntyp1),bcavtub(ntyp1),ccavtub(ntyp1),&
+        dcavtub(ntyp1))
+      allocate(sc_aa_tube_par(ntyp1),sc_bb_tube_par(ntyp1),&
+        tubetranene(ntyp1))
       aa_aq(:,:)=0.0D0
       bb_aq(:,:)=0.0D0
       aa_lip(:,:)=0.0D0
       chi(:,:)=0.0D0
       sigma(:,:)=0.0D0
       r0(:,:)=0.0D0
+      acavtub(:)=0.0d0
+      bcavtub(:)=0.0d0
+      ccavtub(:)=0.0d0
+      dcavtub(:)=0.0d0
+      sc_aa_tube_par(:)=0.0d0
+      sc_bb_tube_par(:)=0.0d0
+
 !--------------------------------
 
       do i=2,ntyp
          endif
         enddo
       enddo
+      write(iout,*) "tube param"
+      read(itube,*) epspeptube,sigmapeptube,acavtubpep,bcavtubpep, &
+      ccavtubpep,dcavtubpep,tubetranenepep
+      sigmapeptube=sigmapeptube**6
+      sigeps=dsign(1.0D0,epspeptube)
+      epspeptube=dabs(epspeptube)
+      pep_aa_tube=4.0d0*epspeptube*sigmapeptube**2
+      pep_bb_tube=-sigeps*4.0d0*epspeptube*sigmapeptube
+      write(iout,*) pep_aa_tube,pep_bb_tube,tubetranenepep
+      do i=1,ntyp
+       read(itube,*) epssctube,sigmasctube,acavtub(i),bcavtub(i), &
+      ccavtub(i),dcavtub(i),tubetranene(i)
+       sigmasctube=sigmasctube**6
+       sigeps=dsign(1.0D0,epssctube)
+       epssctube=dabs(epssctube)
+       sc_aa_tube_par(i)=4.0d0*epssctube*sigmasctube**2
+       sc_bb_tube_par(i)=-sigeps*4.0d0*epssctube*sigmasctube
+      write(iout,*) sc_aa_tube_par(i), sc_bb_tube_par(i),tubetranene(i)
+      enddo
 
       allocate(aad(ntyp,2),bad(ntyp,2)) !(ntyp,2)
       bad(:,:)=0.0D0
 
       return
       end subroutine readpdb
-#if .not. defined(WHAM_RUN) && .not. defined(CLUSTER)
+#if !defined(WHAM_RUN) && !defined(CLUSTER)
 !-----------------------------------------------------------------------------
 ! readrtns_CSA.F
 !-----------------------------------------------------------------------------
 !C      if(me.eq.king .or. .not. out1file .and. fg_rank.eq.0) then
         write(iout,*) "shield_mode",shield_mode
 !C  Varibles set size of box
+      with_theta_constr = index(controlcard,"WITH_THETA_CONSTR").gt.0
+      write (iout,*) "with_theta_constr ",with_theta_constr
+      AFMlog=(index(controlcard,'AFM'))
+      selfguide=(index(controlcard,'SELFGUIDE'))
+      print *,'AFMlog',AFMlog,selfguide,"KUPA"
+      call readi(controlcard,'GENCONSTR',genconstr,0)
       call reada(controlcard,'BOXX',boxxsize,100.0d0)
       call reada(controlcard,'BOXY',boxysize,100.0d0)
       call reada(controlcard,'BOXZ',boxzsize,100.0d0)
+      call readi(controlcard,'TUBEMOD',tubemode,0)
+      write (iout,*) TUBEmode,"TUBEMODE"
+      if (TUBEmode.gt.0) then
+       call reada(controlcard,"XTUBE",tubecenter(1),0.0d0)
+       call reada(controlcard,"YTUBE",tubecenter(2),0.0d0)
+       call reada(controlcard,"ZTUBE",tubecenter(3),0.0d0)
+       call reada(controlcard,"RTUBE",tubeR0,0.0d0)
+       call reada(controlcard,"TUBETOP",bordtubetop,boxzsize)
+       call reada(controlcard,"TUBEBOT",bordtubebot,0.0d0)
+       call reada(controlcard,"TUBEBUF",tubebufthick,1.0d0)
+       buftubebot=bordtubebot+tubebufthick
+       buftubetop=bordtubetop-tubebufthick
+      endif
+
 ! CUTOFFF ON ELECTROSTATICS
       call reada(controlcard,"R_CUT_ELE",r_cut_ele,15.0d0)
       call reada(controlcard,"LAMBDA_ELE",rlamb_ele,0.3d0)
       open (isidep,file=sidename,status='old',action='read')
       call getenv_loc('LIPTRANPAR',liptranname)
       open (iliptranpar,file=liptranname,status='old',action='read')
+      call getenv_loc('TUBEPAR',tubename)
+      open (itube,file=tubename,status='old',action='read')
+
 !      print *,"Processor",myrank," opened file ISIDEP" 
 !      print *,"Processor",myrank," opened parameter files" 
 #elif (defined G77)
       open (isidep,file=sidename,status='old')
       call getenv_loc('LIPTRANPAR',liptranname)
       open (iliptranpar,file=liptranname,status='old')
+      call getenv_loc('TUBEPAR',tubename)
+      open (itube,file=tubename,status='old')
 #else
       open(1,file=pref_orig(:ilen(pref_orig))//'.inp',status='old',&
         readonly)
       open (isidep,file=sidename,status='old',readonly)
       call getenv_loc('LIPTRANPAR',liptranname)
       open (iliptranpar,file=liptranname,status='old',action='read')
+      call getenv_loc('TUBEPAR',tubename)
+      open (itube,file=tubename,status='old',action='read')
+
 #ifndef CRYST_SC
       call getenv_loc('ROTPARPDB',rotname_pdb)
       open (irotam_pdb,file=rotname_pdb,status='old',action='read')
 
       open(irest2,file=rest2name,status='unknown')
       read(irest2,*) totT,EK,potE,totE,t_bath
+      totTafm=totT
 !      do i=1,2*nres
 ! AL 4/17/17: Now reading d_t(0,:) too
       do i=0,2*nres