!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! Simulation - Competition for leadership !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! program Wave IMPLICIT REAL*8 (A-H,O-Z), INTEGER(I-N) Real(8), ALLOCATABLE ::y1(:),y2(:),y3(:),y4(:),y5(:) CHARACTER(15) :: name,w,name1 common /c1/ ncomb,nwaves,nsave,npe common /c2/ name1 !!!! Parameters !!!! npe=389 ! Period normalized (one wave) nesp=32 ! Number of species ndados=622400 ! npe*10*5*32 nrep=10000 ! Simulation repetitions nwaves=20 ! Waves by simulation nsave =0 !!!! Species !!!! !names: "lepto", !"burgersi", !"mordax", !"rapax", !"deichmanni", !"hirsutimanus", !"intermedia", !"beebei", !"thayeri", !"uruguayensis", !"heteropleura", !"seismella", !"polita", !"elegans", !"galapagensis", !"perplexa", !"terpsichores", !"capricornis", !"dussumieri", !"stenodactylus", !"flammula", !"signata", !"vomeris", !"mjoebergi", !"coarctata", !"longidigitum", !"herradurensis", !"batuenta", !"inaequalis", !"oerstedi", !"saltitanta", !"panamensis", name1="lepto" ! Change the species above OPEN(UNIT=53,FILE='apice_with_effect_'//Trim(adjustl(name1))//'.txt',STATUS='UNKNOWN') OPEN(UNIT=54,FILE='apice_without_effect_'//Trim(adjustl(name1))//'.txt',STATUS='UNKNOWN') n1=npe*10 n2=npe*10 n3=npe*10 n4=npe*10 n5=npe*10 Allocate (y1(n1),y2(n2),y3(n3),y4(n4),y5(n5)) OPEN(UNIT=61,FILE='normalized.txt',STATUS='old') j1=0 j2=0 j3=0 j4=0 j5=0 Do i=1,ndados read(61,*) w, name , nx, nxx , w0, yt if(name.eq.name1) then if(nx.eq.1) then j1=j1+1 y1(j1)=yt endif if(nx.eq.2)then j2=j2+1 y2(j2)=yt endif if(nx.eq.3)then j3=j3+1 y3(j3)=yt endif if(nx.eq.4)then j4=j4+1 y4(j4)=yt endif if(nx.eq.5)then j5=j5+1 y5(j5)=yt endif endif Enddo Close(61) ncomb=0 Do l=1, nrep Do i=1,5 Do ii=i,5 ncomb=ncomb+1 if (i.eq.1.and.ii.eq.2) call dinamic (n1,n2,y1,y2,i,ii) if (i.eq.1.and.ii.eq.3) call dinamic (n1,n3,y1,y3,i,ii) if (i.eq.1.and.ii.eq.4) call dinamic (n1,n4,y1,y4,i,ii) if (i.eq.1.and.ii.eq.5) call dinamic (n1,n5,y1,y5,i,ii) if (i.eq.2.and.ii.eq.3) call dinamic (n2,n3,y2,y3,i,ii) if (i.eq.2.and.ii.eq.4) call dinamic (n2,n4,y2,y4,i,ii) if (i.eq.2.and.ii.eq.5) call dinamic (n2,n5,y2,y5,i,ii) if (i.eq.3.and.ii.eq.4) call dinamic (n3,n4,y3,y4,i,ii) if (i.eq.3.and.ii.eq.5) call dinamic (n3,n5,y3,y5,i,ii) if (i.eq.4.and.ii.eq.5) call dinamic (n4,n5,y4,y5,i,ii) enddo enddo enddo close(53) close(54) 10 FORMAT(11(1X,I2)) END PROGRAM Wave !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! subroutine dinamic(n1,n2,y1,y2,i,ii) IMPLICIT REAL*8 (A-H,O-Z), INTEGER(I-N) CHARACTER(15) :: name1 common /c1/ ncomb,nwaves,nsave,npe common /c2/ name1 Dimension::y1(n1),y2(n2) call random_number(r1) call random_number(r2) n01=1+Int(r1*n1) n02=1+Int(r2*n2) m01=n01 m02=n02 nyant1=-1 nyant2=-1 nyant11=-1 nyant22=-1 napice1=0 napice2=0 napice=0 n=0 Do n=n+1 ncode=5 n01f=n01+1 n01ff=n01+2 if (mod(n01,npe).eq.0) then call random_number(r1) n01f=1+npe*(int(r1*10.0)) n01ff=n01f+1 endif if (mod(n01,(npe)).eq.(npe-1)) then call random_number(r1) n01ff=1+npe*(int(r1*10.0)) endif n02f=n02+1 n02ff=n02+2 if (mod(n02,npe).eq.0) then call random_number(r1) n02f=1+npe*(int(r1*10.0)) n02ff=n02f+1 endif if (mod(n02,(npe)).eq.(npe-1)) then call random_number(r1) n02ff=1+npe*(int(r1*10.0)) endif delta1=y1(n01f)- y1(n01) delta2=y2(n02f)- y2(n02) y10=y1(n01) y20=y2(n02) var=y1(n01)-y2(n02) if(delta1.gt.0.0.and.delta2.gt.0.0.and.var.ne.0.0) then if(var.gt.0.0) then y1int=y1(n01)+delta1/2 y2int=y2(n02f) n01=n01f n02=n02ff ncode=2 else y1int=y1(n01f) y2int=y2(n02)+delta2/2 n01=n01ff n02=n02f ncode=1 endif else if(delta1.lt.0.0.and.delta2.lt.0.0.and.var.ne.0)then if(var.lt.0.0) then y1int=y1(n01)+delta1/2 y2int=y2(n02f) n01=n01f n02=n02ff ncode=-2 else y1int=y1(n01f) y2int=y2(n02)+delta2/2 n01=n01ff n02=n02f ncode=-1 endif else y1int=y1(n01)+delta1/2 y2int=y2(n02)+delta2/2 n01=n01f n02=n02f ncode=0 endif endif if(int(y10).eq.-1.or.int(y1int).eq.-1) nyant1=-1 if((int(y10).eq.1.or.int(y1int).eq.1).and.nyant1.eq.-1) then if((napice.gt.nsave).or.(napice.gt.nsave)) then write(53,101) name1, ncomb, 1, n endif nyant1=1 napice=napice+1 if(napice.gt.nwaves) exit endif if(int(y20).eq.-1.or.int(y2int).eq.-1) nyant2=-1 if((int(y20).eq.1.or.int(y2int).eq.1).and.nyant2.eq.-1) then if((napice.gt.nsave).or.(napice.gt.nsave)) then write(53,101) name1, ncomb, 2, n endif nyant2=1 napice=napice+1 if(napice.gt.nwaves) exit endif if(int(y1(m01)).eq.-1) nyant11=-1 if(int(y1(m01)).eq.1.and.nyant11.eq.-1) then if((napice.gt.nsave).or.(napice.gt.nsave)) then write(54,101) name1, ncomb, 1, n endif nyant11=1 if(napice.gt.nwaves) exit endif if(int(y2(m02)).eq.-1) nyant22=-1 if(int(y2(m02)).eq.1.and.nyant22.eq.-1) then if((napice.gt.nsave).or.(napice.gt.nsave)) then write(54,101) name1, ncomb, 2, n endif nyant22=1 if(napice.gt.nwaves) exit endif if (mod(m01,npe).eq.0) then m01=1+npe*int(real(n01-1)/real(npe)) else m01=m01+1 endif if (mod(m02,npe).eq.0) then m02=1+npe*int(real(n02-1)/real(npe)) else m02=m02+1 endif enddo 101 FORMAT(A15,1x,I7,1x,I2,1x,I5) endsubroutine dinamic SUBROUTINE init_random_seed() INTEGER :: l, k, clock INTEGER, DIMENSION(:), ALLOCATABLE :: seed CALL RANDOM_SEED(size = k) ALLOCATE(seed(k)) CALL SYSTEM_CLOCK(COUNT=clock) seed = clock + 37 * (/ (l - 1, l = 1, k) /) CALL RANDOM_SEED(PUT = seed) DEALLOCATE(seed) END SUBROUTINE init_random_seed !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!