From: Dawid Jagiela Date: Mon, 12 Mar 2012 12:12:06 +0000 (+0100) Subject: csacutoff is now calculated from avrage distance (avrd) between individuals in first... X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?a=commitdiff_plain;h=b38671b328206bf6c62f0935c8dba18f3c906250;p=unres.git csacutoff is now calculated from avrage distance (avrd) between individuals in first bank. Border values are taken as a factor off this value. Factors are read from main input (minco,maxco). --- diff --git a/source/ga/GA.f b/source/ga/GA.f index 8d5b1aa..cfb58a3 100644 --- a/source/ga/GA.f +++ b/source/ga/GA.f @@ -11,7 +11,7 @@ integer*4 :: mnoznik = 10 character*16 :: mode,argn,argi,dir_nam character*15 :: filename - real*8 :: loc,fc,sumfitness,sump,avrp,maxfc,avrd + real*8 :: loc,fc,sumfitness,sump,avrp,maxfc real :: rand, r real*8 :: maxfcever real*8 :: FunkcjaCeluB, FunkcjaCeluB2, FunkcjaCeluR @@ -31,7 +31,7 @@ c ======================================================================= call ReadInput(status) - if ((maxminstep.gt.0).and.(generation.eq.1)) then + if ((maxminstep.gt.0).and.(generation.eq.0)) then c Do first score do_fs=.true. endif @@ -78,8 +78,11 @@ c No. Prepere next generation open(ostate,file=ostatefn) read(ostate,'(I4)') generation read(ostate,'(L2)') do_optima - read(ostate,'(L2)') do_ga + read(ostate,'(L2)') do_ga + read(ostate,'(F7.5)') avrd close(ostate) + write(tmptext,'(I4)') generation + call write2log("This is genaration "//tmptext) if (do_ga) then call write2log("Doing GA in this step") endif @@ -94,6 +97,8 @@ c No. Prepere next generation c Yes. Generate random zero-population else call GenPopulation(BANK_MULTIPLIER*banksize,populacja) + write(tmptext,'(I4)') generation + call write2log("This is genaration "//tmptext) if (do_optima) then do_ga=.false. endif @@ -142,17 +147,17 @@ c --- debug begin --- write(*,*) "do_ga: ",do_ga write(*,*) "do_optima: ",do_optima write(*,*) "do_fs: ",do_fs - write(*,*) "cicutoff: ", cicutoff +c write(*,*) "cicutoff: ", cicutoff endif c --- debug end --- - csacutoff=cicutoff +c csacutoff=cicutoff c c Fill the bank just after the first time we get the score c - if (((generation.eq.2).and.(maxminstep.eq.0)).or.((generation.eq& - &.2).and.(maxminstep.gt.0))) then + if (((generation.eq.1).and.(maxminstep.eq.0)).or.((generation.eq& + &.1).and.(maxminstep.gt.0))) then c --- debug begin --- if (debug) then @@ -160,6 +165,14 @@ c --- debug begin --- endif c --- debug end --- call GetNBest(populacja,bank,banksize) + call CalcAvgDist(bank,avrd) + write(tmptext,'(F7.5)') csacutoff + call write2log("CSA cutoff is now set to "//tmptext) + csacutoff=(maxco*avrd)-generation*avrd*(maxco-minco)/maxgen +c csacutoff=maxco*avrd + + + c --- debug begin --- if (debug) then do i=1,banksize @@ -183,7 +196,6 @@ c --- debug begin --- c --- debug end --- call ReadBank(bank) - call CalcAvgDist(bank,avrd) write(tmptext,'(F7.5)') avrd call write2log("Average distance in bank is "//trim(tmptext)) @@ -228,11 +240,12 @@ c --- debug begin --- c --- debug end --- call WriteBank(bank) - csacutoff=csacutoff-(generation*cicutoff/maxgen) + + csacutoff=maxco*avrd-generation*avrd*(maxco-minco)/maxgen +c csacutoff=csacutoff-(generation*cicutoff/maxgen) c csacutoff=cicutoff*(0.8**(iter-1)) write(tmptext,'(F7.5)') csacutoff - call write2log("CSA cutoff is now set to "//tmptext) endif case('cluster') @@ -396,8 +409,8 @@ c c Create the inputs c write(tmptext,'(I)') generation - call write2log("Preparing inputs for generation "//trim(adjustl(tm& - &ptext))) + call write2log("Preparing inputs for next generation ("//trim(adju& + &stl(tmptext))//")") call CreateInputs(BANK_MULTIPLIER*banksize,populacja) c c All done? Then let's quit. @@ -818,29 +831,31 @@ c ====================================================================== do read(inp, '(A)', iostat=stat) wiersz if (stat /= 0) exit - if ((wiersz(1:4).eq.'PDB=').or.(wiersz(1:4).eq.'pdb=')) then - npdb=1 - tmp = wiersz(5:len_trim(wiersz)) - do i=1,len_trim(tmp) - if (tmp(i:i).eq.' ') then - npdb=npdb+1 - endif - end do - if (npdb.gt.maxnpdb) then - call write2log("Number of input PDB exceeds maxnpdb!") - status = 1 - exit - endif - do i=1,npdb - if (index(trim(tmp),' ').gt.0) then - pdbfiles(i)=tmp(1:index(trim(tmp),' ')) - else - pdbfiles(i)=tmp(1:len_trim(tmp)) - endif - tmp=tmp(index(trim(tmp),' ')+1:len_trim(tmp)) +c PDB= + if ((wiersz(1:4).eq.'PDB=').or.(wiersz(1:4).eq.'pdb=')) then + npdb=1 + tmp = wiersz(5:len_trim(wiersz)) + do i=1,len_trim(tmp) + if (tmp(i:i).eq.' ') then + npdb=npdb+1 + endif + end do + if (npdb.gt.maxnpdb) then + call write2log("Number of input PDB exceeds maxnpdb!") + status = 1 + exit + endif + do i=1,npdb + if (index(trim(tmp),' ').gt.0) then + pdbfiles(i)=tmp(1:index(trim(tmp),' ')) + else + pdbfiles(i)=tmp(1:len_trim(tmp)) + endif + tmp=tmp(index(trim(tmp),' ')+1:len_trim(tmp)) end do endif ! Koniec czytania "PDB=" - + +c ALG= if ((wiersz(1:4).eq.'ALG=').or.(wiersz(1:4).eq.'alg=')) then alg=wiersz(5:len_trim(wiersz)) select case(alg) @@ -858,22 +873,40 @@ c ====================================================================== alg="csa" call write2log ("Unknown algorithm. Using 'csa' as default") end select - endif ! Koniec czytania "ALG=" + endif +c GENERATIONS= select case (wiersz(1:12)) case('GENERATIONS=','generations=') tmp = wiersz(13:len_trim(wiersz)) read(tmp,'(I)') maxgen - end select ! Koniec czytania "GENERATIONS=" + end select c CICUTOFF= - select case(wiersz(1:9)) - case('CICUTOFF=','cicutoff=') - tmp = wiersz(10:len_trim(wiersz)) - read(tmp(1:len_trim(tmp)),'(F7.5)') cicutoff - call write2log("Initial CSA cutoff is set to "//tmp) +c select case(wiersz(1:9)) +c case('CICUTOFF=','cicutoff=') +c tmp = wiersz(10:len_trim(wiersz)) +c read(tmp(1:len_trim(tmp)),'(F7.5)') cicutoff +c call write2log("Initial CSA cutoff is set to "//tmp) +c end select + +c MINCO= + select case(wiersz(1:6)) + case('MINCO=','minco=') + tmp = wiersz(7:len_trim(wiersz)) + read(tmp(1:len_trim(tmp)),'(F7.5)') minco + call write2log("Minimal CSA cutoff factor is set to "//tmp) end select +c MAXCO= + select case(wiersz(1:6)) + case('MAXCO=','maxco=') + tmp = wiersz(7:len_trim(wiersz)) + read(tmp(1:len_trim(tmp)),'(F7.5)') maxco + call write2log("Maximal CSA cutoff factor is set to "//tmp) + end select + + c POPULATION= select case(wiersz(1:11)) case('POPULATION=','population=') @@ -1338,6 +1371,7 @@ c ---------------------------------------------------------------------- write(ostate,'(I4)') generation write(ostate,'(L2)') do_optima write(ostate,'(L2)') do_ga + write(ostate,'(F7.5)') avrd close(ostate) end subroutine WriteState diff --git a/source/ga/common.inc b/source/ga/common.inc index c484ece..b392108 100644 --- a/source/ga/common.inc +++ b/source/ga/common.inc @@ -5,12 +5,14 @@ integer*4 :: nscripts = 0 ! number of shell scripts to copy integer*4 :: maxgen = 0 ! maximum number of generations integer*4 :: banksize = 0 ! size of bank - integer*4 :: generation = 1 ! current generation + integer*4 :: generation = 0 ! current generation integer*4 :: maxminstep = 0 ! max minimalization to be done by zscore integer*4,parameter :: maxnpdb = 10 ! hard limit for maximum number of proteins integer*4,parameter :: maxscripts = 10 - real*8 :: cicutoff = 5.0 ! CSA initial cutoff real*8 :: csacutoff = 0.0 ! CSA cutoff + real*8 :: minco = 0.0 ! minimal CSA cutoff factor + real*8 :: maxco = 0.0 ! maximal CSA cutoff factor + real*8 :: avrd = 0.0 ! average distance between ind in first bank logical :: do_optima = .false. logical :: do_ga = .false. logical :: do_fs = .false. @@ -20,9 +22,9 @@ character*32 :: mremdtemplate(maxnpdb) character*16 :: alg - common /inputy/ npdb,nscripts,ntwham,ntmremd,maxgen,banksize,cicut& - &off,csacutoff,alg,pdbfiles,scripts,whamtemplate,mremdtemplate,gene& - &ration,maxminstep,do_optima, do_ga + common /inputy/ npdb,nscripts,ntwham,ntmremd,maxgen,banksize,csacu& + &toff,avrd,minco,maxco,alg,pdbfiles,scripts,whamtemplate,mremdtempl& + &ate,generation,maxminstep,do_optima, do_ga character*7 :: version = "1.1.1" character*50 :: info = "= Last modified by Lightnir 09/03/2012" real*8,allocatable :: bank(:,:),populacja(:,:),temppopulacja(:,:)