program unresag implicit none include 'GA.inc' include 'common.inc' include 'io.inc' real*8 :: WEIGHTS(18) integer*4 :: i,j,m,n,iloc,pos,iter,maxiter,numarg,tmplen integer*4 :: WybierzOsobnika,Najlepszy,ZnajdzPodobnego,FindWorst integer*4 :: status integer*4 :: mnoznik = 10 character*16 :: mode,argn,argi,dir_nam character*15 :: filename real*8 :: loc,fc,sumfitness,sump,avrp,maxfc real :: rand, r real*8 :: maxfcever real*8 :: FunkcjaCeluB, FunkcjaCeluB2, FunkcjaCeluR character(80),dimension(5) :: wagi logical :: searching = .false. logical :: logdata = .true. logical :: fileex = .false. logical :: retry = .false. logical :: dir_ex = .false. logical :: debug = .true. character*200 :: text, tmptext,tmptext2 character(len=*), parameter :: FMT = "(19F10.5,E15.7,F10.5)" c ======================================================================= c Main program c ======================================================================= call ReadInput(status) if ((maxminstep.gt.0).and.(generation.eq.0)) then c Do first score do_fs=.true. endif c Is everything OK? if (status.gt.0) then write(*,*) "There ware input file errors - check log." goto 2010 endif call itime(timeArray) r = rand(timeArray(1)*10000+timeArray(2)*100+timeArray(3)) c Fair enough? c r = 527865426 c------------------------------------------------------- c bank & populacja arrays store this c 1-19 - weights c 20 - energy obtained from zscore c 21 - fittness c------------------------------------------------------- allocate (bank(banksize,21)) allocate (populacja(BANK_MULTIPLIER*banksize,21)) allocate (pairs(BANK_MULTIPLIER*banksize)) c c Zero the arrays c do i=1,BANK_MULTIPLIER*banksize do j=1,21 populacja(i,j)=0 enddo pairs(i)=0 enddo do i=1,banksize do j=1,21 bank(i,j)=0 enddo enddo c First time here? inquire(FILE=ostatefn,EXIST=fileex) c No. Prepere next generation if (fileex) then open(ostate,file=ostatefn) read(ostate,'(I4)') generation read(ostate,'(L2)') do_optima read(ostate,'(L2)') do_ga read(ostate,'(F7.5)') avrd close(ostate) if (do_ga) then call write2log("Doing GA in this step") endif if (.not.do_optima) then maxminstep=0 call write2log("ZSCORE weights minimalization disabled for now") generation=generation+1 endif write(tmptext,'(I4)') generation call write2log("This is genaration "//tmptext) call ReadOptimaW(BANK_MULTIPLIER*banksize,populacja) 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 c End of "first time here?" code endif c c Do we actual use a Genetic Algorithm? c if (do_ga.eqv..true.) then c Yes. This is only done when doing second pass with ZSCORE (without weight minimalizaton) c or with weight minimalization disabled (maxmin=0) just to obtain final score) c \\---//================================================================c c \\-// c c "// Genetic Algorithm code starts here c c //-\\ c c =//---\\===============================================================c call ReadZEnergy(BANK_MULTIPLIER*banksize,populacja) do i=1,BANK_MULTIPLIER*banksize c populacja(i,20)=FunkcjaCeluR(populacja(i,20)) write(*,*) "Ind ",i,populacja(i,20) enddo call write2log("==[ Doing GA now ]=============================") c c Populate bank depending on algorithm c select case(alg) case('simple') call GetNBest(populacja,bank,banksize) case('better') call GetNBest(populacja,bank,banksize) case('betterv2') call GetNBest(populacja,bank,banksize) case('csa') c --- debug begin --- if (debug) then write(*,*) "generation: ",generation write(*,*) "maxminstep: ",maxminstep write(*,*) "do_ga: ",do_ga write(*,*) "do_optima: ",do_optima write(*,*) "do_fs: ",do_fs c write(*,*) "cicutoff: ", cicutoff endif c --- debug end --- c csacutoff=cicutoff c c Fill the bank just after the first time we get the score c if (((generation.eq.1).and.(maxminstep.eq.0)).or.((generation.eq& &.1).and.(maxminstep.gt.0))) then c --- debug begin --- if (debug) then write(*,*) "Calling GetNBest..." endif c --- debug end --- call GetNBest(populacja,bank,banksize) call CalcAvgDist(bank,avrd) write(tmptext,'(F7.5)') avrd call write2log("Average distance between individuals in initial& & bank is "//trim(tmptext)) c c Cutoff c csacutoff=(maxco*avrd)-generation*avrd*(maxco-minco)/maxgen write(tmptext,'(F7.5)') csacutoff call write2log("CSA cutoff is now set to "//trim(tmptext)) c csacutoff=maxco*avrd c --- debug begin --- if (debug) then do i=1,banksize write(*,FMT) bank(i,:) enddo endif c --- debug end --- call CalcFitness(bank,sumfitness) call WriteBank(bank) c c Do normal CSA c else c --- debug begin --- if (debug) then write(*,*) "Szukam podobnego" endif c --- debug end --- call ReadBank(bank) write(tmptext,'(F7.5)') avrd call write2log("Average distance in bank is "//trim(tmptext)) csacutoff=maxco*avrd-generation*avrd*(maxco-minco)/maxgen write(tmptext,'(F7.5)') csacutoff call write2log("CSA cutoff is now set to "//trim(tmptext)) do i=1,BANK_MULTIPLIER*banksize write(tmptext,'(I4)') i call write2log("Checking ind "//trim(tmptext)) j=ZnajdzPodobnego(banksize,bank,populacja(i,:),csacutoff) if (j.gt.0) then ! W banku jest podobny if (populacja(i,20).lt.bank(j,20)) then write(tmptext,'(I4)') j write(tmptext2,'(I4)') i call write2log(" Swaping ind"//trim(tmptext)//" from bank t& &o ind "//trim(tmptext2)//" from population") write(tmptext2,'(19F8.5,E15.7,F8.5)') bank(j,:) call write2log(" BANK"//trim(tmptext)//":"//trim(tmptext2)) write(tmptext,'(I4)') i write(tmptext2,'(19F8.5,E15.7,F8.5)') populacja(i,:) call write2log(" POP "//trim(tmptext)//":"//trim(tmptext2)) bank(j,:)=populacja(i,:) else call write2log(" Found simialar but not better") endif else ! W banku nie ma podobnego j=FindWorst(banksize,bank) write(tmptext2,'(I4)') j if (populacja(i,20).lt.bank(j,20)) then call write2log(" Worst in bank is "//trim(tmptext2)) write(tmptext,'(I4)') i call write2log(" Swaping worst ind in bank to "//trim(tmpte& &xt)) write(tmptext,'(I4)') j write(tmptext2,'(19F8.5,E15.7,F8.5)') bank(j,:) call write2log(" BANK"//trim(tmptext)//":"//trim(tmptext2)) write(tmptext,'(I4)') i write(tmptext2,'(19F8.5,E15.7,F8.5)') populacja(i,:) call write2log(" POP "//trim(tmptext)//":"//trim(tmptext2)) bank(j,:)=populacja(i,:) else call write2log(" The worst in bank is better then this Ind"& &) endif endif enddo c ----------------------------------------------------------------------- c Calculate fitness for reproduction c ----------------------------------------------------------------------- call CalcFitness(bank,sumfitness) c --- debug begin --- if (debug) then write (*,*) "Dumping bank: " do i=1,banksize write(*,FMT) bank(i,:) enddo write (*,*) "Sumfitness: ", sumfitness endif c --- debug end --- call WriteBank(bank) endif case('cluster') write(*,*) "Well this is not implemented yet" goto 2010 end select c call WritePopSum() c Have we done last generation ? if (generation.eq.maxgen) then call write2log("Maximum number of genarations reached. QUITING.") goto 2010 endif c Prapering for next generation so increase counter iter=iter+1 c ----------------------------------------------------------------------- c Reproduction c ----------------------------------------------------------------------- call write2log("==[ Reproducting ]============================") select case(alg) case('simple','csa') do i=1,BANK_MULTIPLIER*banksize m=WybierzOsobnika(banksize,bank) write(tmptext,'(I4)') m write(tmptext2,'(I4)') i call write2log("Reproducting individual "//trim(tmptext)//" from& & bank as new individual "//trim(tmptext2)) populacja(i,:)=bank(m,:) end do case('better','betterv2') populacja(1,:)=bank(1,:) call write2log("Reproducting first (best) individual") do i=2,BANK_MULTIPLIER*banksize m=WybierzOsobnika(banksize,bank) write(tmptext,'(I4)') m write(tmptext2,'(I4)') i call write2log("Reproducting individual "//trim(tmptext)//" from& & bank as new individual "//trim(tmptext2)) populacja(i,:)=bank(m,:) end do end select c ----------------------------------------------------------------------- c Pairing c ----------------------------------------------------------------------- do i=1,BANK_MULTIPLIER*banksize pairs(i)=0 end do do i=1,BANK_MULTIPLIER*banksize if (pairs(i).eq.0) then 210 pos=1+int(rand(0)/(1/real(BANK_MULTIPLIER*banksize))) if (pos.le.i) then goto 210 endif do j=1,(i-1) if (pairs(j).eq.pos) then goto 210 endif end do pairs(i) = pos pairs(pos) = i endif end do call write2log("==[ Pairing ]==============================") do i=1,BANK_MULTIPLIER*banksize write(tmptext,'(A11,I3.3,A5,I3.3,A21)') "Individual ",i," and ",p& &airs(i)," will be crossed over" call write2log(trim(tmptext)) end do c ----------------------------------------------------------------------- c Crossover c ----------------------------------------------------------------------- call write2log("==[ Crossing ]==============================") do i=1,BANK_MULTIPLIER*banksize if (pairs(i) .gt. i) then call CrossoverS(populacja(i,:),populacja(pairs(i),:)) endif end do call write2log(" WLONG WSCP WELEC WBOND WANG & &WSCLOC WTOR WTORD WCORRH WCORR4 WCORR5 WCORR6 WEL_LOC WT& &URN3 WTURN4 WTURN6 WVDWPP WHPB WSCCOR") do i=1,BANK_MULTIPLIER*banksize tmplen = 9 text = "" write(text,'(A4,I3.3,A2)') "Ind",i,": " do j=1,19 write(tmptext,'(F7.5)') populacja(i,j) text = text(1:tmplen) // tmptext tmplen = tmplen + 8 end do call write2log(trim(text)) end do c ----------------------------------------------------------------------- c Mutation c ----------------------------------------------------------------------- call write2log("==[ Mutating ]==============================") do i = 1,19 if (IGNORE_WEIGHT(i).eq.1) then write(tmptext,'(I2)') i call write2log("Parameter "//trim(WNAME(i))//" will be ignored") endif end do do i=1,BANK_MULTIPLIER*banksize call MutationV(populacja(i,:),i,MUTRATIO) end do c for betterv2 and csa restore best individual select case(alg) case('betterv2','csa') populacja(1,:)=bank(1,:) end select c c Genetic Algorithm blokcode stops here c endif c --------------------------------- c Do some final stuff c --------------------------------- c c Set the variables before writing to state file c if ((.not.do_ga).and.(.not.do_optima)) then if (do_fs) then do_ga=.false. do_optima=.true. else do_ga=.true. do_optima=.false. endif else if (do_fs) then do_optima=.not.do_optima if (generation.eq.0) then do_ga=.false. else do_ga=.not.do_ga endif else do_ga=.true. do_optima=.false. endif endif c c Write te current state and population c call WriteState() c call WritePopSum() c c Create the inputs c write(tmptext,'(I3)') generation+1 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. c goto 2010 c ====================================================================== c c MAIN CODE ENDS HERE c c ====================================================================== 2000 deallocate(populacja) deallocate(pairs) 2010 end program c ====================================================================== c Clustering subroutines c ====================================================================== include 'cluster.inc' c ====================================================================== c Write2log subroutine c ====================================================================== c Creates logs c ---------------------------------------------------------------------- subroutine Write2log(text) include 'common.inc' include 'io.inc' character*200 :: tmptext character(len=*) :: text logical :: ex call itime(timeArray) inquire(file=ologfn,exist=ex) if (.not.ex) then open(olog, file=ologfn, access="append") write(olog,"(A)") "==============================================& &==========" tmptext = "= UNRES weights optimizer version "//version write(olog,"(A)") tmptext write(olog,"(A)") info write(olog,"(A)") "==============================================& &==========" close(olog) endif open(olog, file = ologfn , access = "append") write(olog,"(I2.2,A1,I2.2,A1,I2.2,A2,A)") timeArray(1),":",timeAr& &ray(2),":",timeArray(3),": ",text close(olog) end subroutine Write2log c ====================================================================== c GetNBest subroutine c ====================================================================== c Fills up the bank population up to banksize with individuals from c inputpop having the best(lowest) score c ---------------------------------------------------------------------- subroutine GetNBest(inputpop,bank,banksize) integer*4 :: banksize, i,j,idx real*8 :: inputpop, bank, best, last include 'GA.inc' dimension inputpop(banksize*BANK_MULTIPLIER,21) dimension bank(banksize,21) last=0.0 idx=1 do j=1,banksize best=huge(0.0d0) do i=1,banksize*BANK_MULTIPLIER if ((inputpop(i,20).lt.best).and.(inputpop(i,20).gt.last)) then best=inputpop(i,20) idx=i endif end do bank(j,:)=inputpop(idx,:) last=best end do end subroutine GetNBest c ====================================================================== c FindWorst subroutine c ====================================================================== c Finds the worst (highest score) individual in population c ---------------------------------------------------------------------- integer*4 function FindWorst(rozmiar,pop) integer*4 :: rozmiar, i, idx real*8 :: pop,last dimension pop(rozmiar,21) last=0.0 idx = 0 do i=1,rozmiar if (pop(i,20).gt.last) then idx = i last=pop(i,20) endif end do FindWorst=idx return end function c ====================================================================== c GenPopulation subroutine c ====================================================================== c Generates a population (pop) of nind individuals with 19 random c wights from a set of ranges (see GA.inc). c ---------------------------------------------------------------------- subroutine GenPopulation(nind,pop) integer*4 :: i,j,nind real*8 :: pop dimension pop(nind,21) include 'GA.inc' do j=1,nind do i=1,19 select case(i) case (1) pop(j,1)=WLONG_LOW+(rand(0)*(WLONG_HI-WLONG_LOW)) case (2) pop(j,2)=WSCP_LOW+(rand(0)*(WSCP_HI-WSCP_LOW)) case (3) pop(j,3)=WELEC_LOW+(rand(0)*(WELEC_HI-WELEC_LOW)) case (4) pop(j,4)=WBOND_LOW+(rand(0)*(WBOND_HI-WBOND_LOW)) case (5) pop(j,5)=WANG_LOW+(rand(0)*(WANG_HI-WANG_LOW)) case (6) pop(j,6)=WSCLOC_LOW+(rand(0)*(WSCLOC_HI-WSCLOC_LOW)) case (7) pop(j,7)=WTOR_LOW+(rand(0)*(WTOR_HI-WTOR_LOW)) case (8) pop(j,8)=WTORD_LOW+(rand(0)*(WTORD_HI-WTORD_LOW)) case (9) pop(j,9)=WCORRH_LOW+(rand(0)*(WCORRH_HI-WCORRH_LOW)) case (10) pop(j,10)=WCORR4_LOW+(rand(0)*(WCORR4_HI-WCORR4_LOW)) case (11) pop(j,11)=WCORR5_LOW+(rand(0)*(WCORR5_HI-WCORR5_LOW)) case (12) pop(j,12)=WCORR6_LOW+(rand(0)*(WCORR6_HI-WCORR6_LOW)) case (13) pop(j,13)=WEL_LOC_LOW+(rand(0)*(WEL_LOC_HI-WEL_LOC_LOW)) case (14) pop(j,14)=WTURN3_LOW+(rand(0)*(WTURN3_HI-WTURN3_LOW)) case (15) pop(j,15)=WTURN4_LOW+(rand(0)*(WTURN4_HI-WTURN4_LOW)) case (16) pop(j,16)=WTURN6_LOW+(rand(0)*(WTURN6_HI-WTURN6_LOW)) case (17) pop(j,17)=WVDWPP_LOW+(rand(0)*(WVDWPP_HI-WVDWPP_LOW)) case (18) pop(j,18)=WHPB_LOW+(rand(0)*(WHPB_HI-WHPB_LOW)) case (19) pop(j,19)=WSCCOR_LOW+(rand(0)*(WSCCOR_HI-WSCCOR_LOW)) end select end do end do end subroutine GenPopulation c ====================================================================== c CrossoverS subroutine c ====================================================================== c Does a simple value crossing-over on a pair c of two given individuals, giving two new individuals. c ---------------------------------------------------------------------- subroutine CrossoverS(ind1,ind2) real*8 :: ind1(19),ind2(19),temp(19) integer*4 :: loc loc = 1 + int(rand(0)/(1/18.0)) temp = ind2 do i=(loc),19 ind2(i)=ind1(i) ind1(i)=temp(i) end do end subroutine CrossoverS c ====================================================================== c MutationV subroutine c ====================================================================== c Does a value mutation for the given individual (ind) with the c given mutation ratio (mratio). c ---------------------------------------------------------------------- subroutine MutationV(ind,id,mratio) include 'GA.inc' c real*8 :: rand real :: mratio real*8 :: ind(19),ti integer*4 :: i,id character*100 :: tmptext character*150 :: logtext do i = 1,19 if (IGNORE_WEIGHT(i).eq.1) then c write(tmptext,'(I2)') i c call write2log("Skipping weight "//trim(tmptext)) goto 2300 endif if (rand(0)