#====================================================================== #============ Simulation parameters ================= preydistribution<- 0 #choose the prey distribution to simulate # 0 = random prey distribution # 1 = aggregated marinelandscape<- 1 #choose the seascape to simulate # 0 = no devices # 1 = devices Nprey<-700 # insert number of prey Npred<- 500 # insernt number of predators # Prey aggregation ncluster<- 30 # insert number of clusters in case of aggregated # prey distribution and no devices in the landscape boxS<-37 # size of the area for the aggregation aggWidth <- 19 # side distance of aggregation area around devices aggHeight <- 10 # top/bottom distance of aggreagation # area around devices #Marine landscape size x1<-100 # depth of the Marine landscape y1<-700 # width of the Marine landscape #Energetic costs zones def.prm<- -8 surface<-95 med<-80 bottom<-x1 #Energetic costs values h.cost<-3 m.cost<-2 low.cost<-1 #Prey energetic value preyE<-1000 #Device value devE<- -8 #Device size n.devices<-10 #number of devices sizedevices<- 15 # size of devices (in paper small devices== 15, # large devices== 50) firstdistance<-20 #distance of the first device that will be placed #from the left edge bottomdistance<- 45 #distance of the device from the bottom edge #(in paper small devices== 45, # large devices== 20) totdevices<-sizedevices * n.devices distance<-(y1 - totdevices)/n.devices rightEdge<-firstdistance + n.devices * sizedevices + (n.devices-1) * distance if(rightEdge>y1){ print("the devices don't fit in the landscape") } #Dive time des<-23 #descending phase sea<-seq(10,300,10) #searching phase tot<-seq(210,500,10) #tot Maximum Dive Time nDiveTime<-length(sea) # number of different dive times maxFish<- 3 #assign maximum number of prey that the predator can catch # -1 <- unlimited (PN in paper) # 3 <- P3 in paper # 1 <- P1 in paper #Example of the matrix with the position of the predator (NA) #and the position of each step around the predator within its percepual range # [,1] [,2] [,3] [,4] [,5] # [1,] 20 21 22 23 24 # [2,] 15 16 17 18 19 # [3,] 11 12 NA 13 14 # [4,] 6 7 8 9 10 # [5,] 1 2 3 4 5 #====== functions for the prey distribution build.marinelandscape<-function() { #------------ Random distribution no devices if (preydistribution==0 && marinelandscape== 0){ for(y in 1:Nprey){ fr[y]<-sample(c(1:y1),1, replace=T) #random sample of X and Y coordinates fci[y]<-sample(c(1:(med-1)),1, replace=T) if ((landscape[fr[y],fci[y]])== 1){ landscape[fr[y],fci[y]]<-preyE} else if ((landscape[fr[y],fci[y]]) == preyE){ #if in the cell there is already fr[y]<-sample(c(1:y1),1, replace=T) #another prey, sample fci[y]<-sample(c(1:(med-1)),1, replace=T) #again landscape[fr[y],fci[y]]<-preyE } } } #------------- Random distribution and devices else if (preydistribution==0 && marinelandscape>= 1) { for (i in 1:n.devices){ #draw the devices leftSide<-firstdistance + ((i-1) * (sizedevices + distance)) landscape[leftSide:(leftSide+sizedevices), bottomdistance:(bottomdistance+sizedevices)]<- devE } #prey distribution for(y in 1:Nprey){ fr[y]<-sample(c(1:y1),1, replace=T) #random sample of X and Y coordinates fci[y]<-sample(c(1:(med-1)),1, replace=T) if ((landscape[fr[y],fci[y]])== 1){ landscape[fr[y],fci[y]]<-preyE} else { if ((landscape[fr[y],fci[y]]) == devE){ #if in the cell there is the fr[y]<-sample(c(1:y1),1, replace=T) # "no through" area sample fci[y]<-sample(c(1:(med-1)),1, replace=T) #another position in while((landscape[fr[y],fci[y]]) == devE){ #the matrix fr[y]<-sample(c(1:y1),1, replace=T) fci[y]<-sample(c(1:(med-1)),1, replace=T) } landscape[fr[y],fci[y]]<-preyE} else { if ((landscape[fr[y],fci[y]]) == preyE){ # if in the cell there is already fr[y]<-sample(c(1:y1),1, replace=T) #another prey, sample fci[y]<-sample(c(1:(med-1)),1, replace=T) #again landscape[fr[y],fci[y]]<-preyE }}} } } #---------------- Aggregated Distribution no devices else if (preydistribution==1 && marinelandscape==0) { numbercl<-sample(seq(1,50,by=1),ncluster,replace=T) for(i in 1) { if (sum(numbercl)==Nprey) { #print(sum(numbercl)) break } repeat { #code if (sum(numbercl)!=Nprey) { numbercl<-sample(seq(1,50,by=1),ncluster,replace=T) if (sum(numbercl)==Nprey) { #print(sum(numbercl)) break } } } } drawbo<-data.frame(x=NA,y=NA) swormc<- 1 # draw the boxS while (swormc<= ncluster) { tempx<-sample(seq(1,(y1-boxS),by=1),1,replace=T) tempy<-sample(seq(1,((med-1)-boxS),by=1),1,replace=T) if(swormc==1){ drawbo[swormc,1]<-tempx drawbo[swormc,2]<-tempy swormc<- swormc+1 } else { for(i in 1:(swormc-1)) { if (tempx > (drawbo[i,1]-boxS) && tempx < (drawbo[i,1]+boxS) && tempy > (drawbo[i,2]-boxS) && tempy < (drawbo[i,2]+boxS)) { break } } drawbo[swormc,1]<-tempx drawbo[swormc,2]<-tempy swormc<- swormc+1 } } for (m in 1:ncluster) { fr<-rep(NA,numbercl[m]) fci<-rep(NA,numbercl[m]) y<-1 while (y <= numbercl[m]) { tempfishx<-sample(c(drawbo[m,1]:(drawbo[m,1]+ boxS)),1, replace=T) #random sample of X and Y coordinates tempfishy<-sample(c(drawbo[m,2]:(drawbo[m,2]+ boxS)),1, replace=T) if((landscape[tempfishx,tempfishy])== low.cost) { fr[y]<-tempfishx fci[y]<-tempfishy landscape[fr[y],fci[y]]<- preyE y<- y + 1 } } } } #----------- Aggregated distribution and devices else if (preydistribution==1 && marinelandscape==1) { prey.devices <- Nprey / n.devices for (i in 1:n.devices){ #draw the devices leftSide<-firstdistance + ((i-1) * (sizedevices + distance)) landscape[leftSide:(leftSide+sizedevices), bottomdistance:(bottomdistance+sizedevices)]<- devE numFishFound <- 0 while (numFishFound < prey.devices) { #newFish <- rnorm2d(1, rho = 0) newFish<-data.frame(NA,NA) colnames (newFish)[1]<- "x" colnames (newFish)[2]<- "y" newFish[1,1] <- runif(1, min = (leftSide - aggWidth), max = (leftSide + sizedevices+aggWidth)) newFish[1,2] <- runif(1, min = (bottomdistance - aggHeight), max = (bottomdistance + sizedevices + aggHeight)) #newFish[1,1] <- newFish[1,1] * ((sizedevices + 2* aggWidth)/2) + (leftSide + (sizedevices/2)) #newFish[1,2] <- newFish[1,2] * ((sizedevices + 2* aggHeight)/2) + (bottomdistance + (sizedevices/2)) newFish[1,1]<-round(newFish[1,1]) newFish[1,2]<-round(newFish[1,2]) #print(newFish) if(newFish[1,1]<=y1 && newFish[1,1]>0 && newFish[1,2]<=x1 && newFish[1,2]>0){ if (landscape[newFish[1,1],newFish[1,2]] == low.cost) { numFishFound <- numFishFound + 1 landscape[newFish[1,1],newFish[1,2]] <- preyE #print(numFishFound) } } } } } return(landscape) } #============================================================================ #============= Starting simulation ========================================== nam<-as.character(seq(1,length(sea),1)) for(ss in 1:nDiveTime){ nft <- array(0,c(Npred,1)) food1 <-array(0,c(Npred,1)) time.dis <- array(0,c(Npred,1)) time.search <- array(0,c(Npred,1)) time.asc <- array(0,c(Npred,1)) for(w in 1:Npred){ fr<-rep(NA,Nprey) fci<-rep(NA,Nprey) prey.cord <-data.frame(x=NA,y=NA) xy.cord <- data.frame(x = NA, y = NA) landscape<- matrix(0, ncol=x1, nrow=y1) landscape[1:y1,surface:x1]<- h.cost landscape[1:y1,med:(surface-1)]<- m.cost landscape[1:y1,1:(med-1)]<- low.cost landscape <- build.marinelandscape() # choose the marinelandscape to simulate # calculate total number of prey for (a in 1:y1){ for(b in 1:(med-1)){ if((landscape[a,b])==preyE){ nft[w,]<- nft[w,]+1 prey.cord[nft[w,],1] <- a prey.cord[nft[w,],2] <- b }}} # dev.co<-data.frame(NA,NA) deviceFieldCounter <- 1 for (a in 1:y1){ for(b in 1:(med-1)){ if((landscape[a,b])==devE){ dev.co[deviceFieldCounter,1] <- a dev.co[deviceFieldCounter,2] <- b deviceFieldCounter <- deviceFieldCounter + 1 }}} ########Move the predator--> FIRST STEP r<-sample(c(1:y1),1, replace=F) ci<-100 landscape[r,ci]<-NA ##DESCENDING PHASE for(gen in 1:des){ xy.cord[gen, ] <- cbind(r, ci) for (r in r){ for(ci in ci){ #create a vector called "value" of "9" in which to store # the values of the perception matrix of the predator value <- array(def.prm,24) # as above to store rows and columns position of each value rows <- array(def.prm,24) cols <- array(def.prm,24) i <- 1 # assign the perceptual range of 2 cells for(r1 in -2:2){ for(c1 in -2:2){ if(r1!=0 || c1!=0){ c2<-ci+c1 r2<-r+r1 if(r2>0 && r2<=y1 && c2>0 && c2<=x1){ value[i]<-landscape[r2,c2] rows[i]<-r2 cols[i]<-c2 } i<-i+1 } } } time.dis[w,]<-time.dis[w,]+1 #calculate time in the descending phase if((value[8])>= h.cost){ landscape[r,ci]<- h.cost #erase predator rr <-12 #if the position choosen is outside the matrix (value 9) while( value[rr] == def.prm){ rr <-sample(c(12,8,17),1, replace=T) } r <- rows[rr] ci <- cols[rr] landscape[r,ci]<-NA } else { if((value[8])== m.cost){ landscape[r,ci]<- m.cost rr <-12 #if the position choosen is outside the matrix (value 9) while( value[rr] == def.prm){ rr <-sample(c(12,8,17),1, replace=T) } r <- rows[rr] ci <- cols[rr] landscape[r,ci]<-NA } else { if(min(value)== low.cost){ landscape[r,ci]<- low.cost rr <-12 while( value[rr] == def.prm){ rr <-sample(c(12,8,17),1, replace=T) } r <- rows[rr] ci <- cols[rr] landscape[r,ci]<-NA } } } } } # jpeg(file = paste("C://Users//Utente//Desktop//Video//",gen,"myplot.jpg")) # plot(prey.cord, xlim=c(1,y1),ylim=c(1,x1), pch=19, col="blue", cex=0.2, xlab="x", ylab="y") # lines(xy.cord, type="l", col="red", lwd=2) # points(dev.co, col="blue",pch=15) # dev.off() #image(landscape, col=c(topo.colors(20),heat.colors(8),"black")) } ####SEARCHING PHASE steps <- 0 pow <-0 xy.moo <- c() for(gen in (des+1):((des+1)+sea[ss])) { xy.cord[gen, ] <- cbind(r, ci) f <- xy.cord[gen,]-xy.cord[gen-1,] for (r in r){ for(ci in ci){ value <- array(def.prm,24) rows <- array(def.prm,24) cols <- array(def.prm,24) i <- 1 for(r1 in -2:2){ for(c1 in -2:2){ if(r1!=0 || c1!=0){ c2<-ci+c1 r2<-r+r1 if(r2>0 && r2<=y1 && c2>0 && c2<=x1){ value[i]<-landscape[r2,c2] rows[i]<-r2 cols[i]<-c2 } i<-i+1 } } } #count time in the searching phase time.search[w,]<-time.search[w,] +1 landscape[r,ci]<- low.cost #erase the predator # calculates means for the 'zones' top.left<- mean(c(value[1],value[2],value[6],value[7]))*sqrt(2) top<- mean(c(value[11],value[12],value[6],value[7],value[16],value[15])) top.right<- mean(c(value[16],value[15],value[21],value[20]))*sqrt(2) middle.left<- mean(c(value[3],value[2],value[8],value[7],value[4],value[9])) bottom.left<- mean(c(value[9],value[10],value[4],value[5]))*sqrt(2) bottom<- mean(c(value[9],value[13],value[19],value[14],value[10],value[18])) bottom.right<- mean(c(value[18],value[24],value[19],value[23]))*sqrt(2) middle.right<- mean(c(value[17],value[18],value[16],value[22],value[23],value[21])) means2 <- c(top.left, middle.left, bottom.left, top, bottom, top.right, middle.right, bottom.right) if(min(value)==devE && (food1[w,])!=maxFish){ rr <-sample(c(7,8,9,12,13,16,17,18),1, replace=T) while( value[rr] == devE){ rr <-sample(c(7,8,9,12,13,16,17,18),1, replace=T) } r <- rows[rr] ci <- cols[rr] landscape[r,ci]<-NA } else { if((min(value)>= low.cost) && (max(value)=devE) && (max(value)==preyE) && (food1[w,])< maxFish)) { j.v<- which(means2==max(means2,na.rm=TRUE),arr.ind=TRUE) ifelse(length(j.v)>1, j<- sample(j.v,1,replace=F), j<-j.v) # print(value) # print(means2) # print(j) if(j==1){ r<-rows[7] ci<-cols[7] steps <- 0 if((landscape[r,ci])==preyE) { food1[w,]<-food1[w,]+1 } landscape[r,ci]<-NA } else { if(j==2){ r<-rows[8] ci<-cols[8] steps <- 0 if((landscape[r,ci])==preyE) { food1[w,]<-food1[w,]+1 } landscape[r,ci]<-NA } else { if(j==3){ r<-rows[9] ci<-cols[9] steps <- 0 if((landscape[r,ci])==preyE) { food1[w,]<-food1[w,]+1 } landscape[r,ci]<-NA } else { if(j==4){ r<-rows[12] ci<-cols[12] steps <- 0 if((landscape[r,ci])==preyE) { food1[w,]<-food1[w,]+1 } landscape[r,ci]<-NA } else { if(j==5){ r<-rows[13] ci<-cols[13] steps <- 0 if((landscape[r,ci])==preyE) { food1[w,]<-food1[w,]+1 } landscape[r,ci]<-NA } else { if(j==6){ r<-rows[16] ci<-cols[16] steps <- 0 if((landscape[r,ci])==preyE) { food1[w,]<-food1[w,]+1 } landscape[r,ci]<-NA } else { if(j==7){ r<-rows[17] ci<-cols[17] steps <- 0 if((landscape[r,ci])==preyE) { food1[w,]<-food1[w,]+1 } landscape[r,ci]<-NA } else { if(j==8){ r<-rows[18] ci<-cols[18] steps <- 0 if((landscape[r,ci])==preyE) { food1[w,]<-food1[w,]+1 } landscape[r,ci]<-NA } } } } } } } } } else { if(maxFish < 0 | ((min(value)>=low.cost) && (max(value)==preyE) && (food1[w,])< maxFish)) { j.v<- which(means2==max(means2,na.rm=TRUE),arr.ind=TRUE) ifelse(length(j.v)>1, j<- sample(j.v,1,replace=F), j<-j.v) # print(value) # print(means2) # print(j) if(j==1){ r<-rows[7] ci<-cols[7] steps <- 0 if((landscape[r,ci])==preyE) { food1[w,]<-food1[w,]+1 } landscape[r,ci]<-NA } else { if(j==2){ r<-rows[8] ci<-cols[8] steps <- 0 if((landscape[r,ci])==preyE) { food1[w,]<-food1[w,]+1 } landscape[r,ci]<-NA } else { if(j==3){ r<-rows[9] ci<-cols[9] steps <- 0 if((landscape[r,ci])==preyE) { food1[w,]<-food1[w,]+1 } landscape[r,ci]<-NA } else { if(j==4){ r<-rows[12] ci<-cols[12] steps <- 0 if((landscape[r,ci])==preyE) { food1[w,]<-food1[w,]+1 } landscape[r,ci]<-NA } else { if(j==5){ r<-rows[13] ci<-cols[13] steps <- 0 if((landscape[r,ci])==preyE) { food1[w,]<-food1[w,]+1 } landscape[r,ci]<-NA } else { if(j==6){ r<-rows[16] ci<-cols[16] steps <- 0 if((landscape[r,ci])==preyE) { food1[w,]<-food1[w,]+1 } landscape[r,ci]<-NA } else { if(j==7){ r<-rows[17] ci<-cols[17] steps <- 0 if((landscape[r,ci])==preyE) { food1[w,]<-food1[w,]+1 } landscape[r,ci]<-NA } else { if(j==8){ r<-rows[18] ci<-cols[18] steps <- 0 if((landscape[r,ci])==preyE) { food1[w,]<-food1[w,]+1 } landscape[r,ci]<-NA } } } } } } } } } else { if ((food1[w,])==maxFish) { time.search[w,]<-time.search[w,] - 1 if ((value[13])== low.cost | (value[13])== preyE){ landscape[r,ci]<- low.cost time.asc[w,]<-time.asc[w,] + 1 rr <-sample(c(8,9,13,17,18),1, replace=T) while( value[rr] == def.prm){ rr <-sample(c(8,9,13,17,18),1, replace=T) } r <- rows[rr] ci <- cols[rr] landscape[r,ci]<-NA } else { if ((max(value)== low.cost &&(value[13])== devE) | (max(value)== preyE &&(value[13])== devE) ){ landscape[r,ci]<- low.cost time.asc[w,]<-time.asc[w,] + 1 rr <-sample(c(7,8,17,16),1, replace=T) while( value[rr] == def.prm){ rr <-sample(c(7,8,17,16),1, replace=T) } r <- rows[rr] ci <- cols[rr] landscape[r,ci]<-NA } else { if ((value[13])== m.cost){ landscape[r,ci]<- value[8] time.asc[w,]<-time.asc[w,] + 1 rr <-13 r <- rows[rr] ci <- cols[rr] landscape[r,ci]<-NA } else { if ((value[13])== h.cost){ landscape[r,ci]<- h.cost time.asc[w,]<-time.asc[w,] + 1 rr <- 13 while( value[rr] == def.prm){ rr <-sample(c(8,17),1, replace=T) } r <- rows[rr] ci <- cols[rr] landscape[r,ci]<-NA } else { if ((value[13])== def.prm && (max(value)== h.cost)){ landscape[r,ci]<- h.cost rr <-sample(c(8,17),1, replace=T) while( value[rr] == def.prm){ rr <-sample(c(8,17),1, replace=T) } r <- rows[rr] ci <- cols[rr] landscape[r,ci]<-NA } } } } } } } } } } } } #image(landscape, col=c(topo.colors(20),heat.colors(8),"black")) # jpeg(file = paste("C://Users//Utente//Desktop//Video//",gen,"myplot.jpg")) #plot(prey.cord, xlim=c(1,y1),ylim=c(1,x1), pch=19, col="blue", cex=0.2, xlab="x", ylab="y") # lines(xy.cord, type="l", col="red", lwd=2) # points(dev.co, col="blue",pch=15) # dev.off() } #} ##ASCHENDING PHASE tse1<-((des+1)+sea[ss])+1 tse2<-(tot[ss]) for(gen in tse1:tse2){ xy.cord[gen, ] <- cbind(r, ci) for (r in r){ for(ci in ci){ value <- array(def.prm,24) rows <- array(def.prm,24) cols <- array(def.prm,24) i <- 1 for(r1 in -2:2){ for(c1 in -2:2){ if(r1!=0 || c1!=0){ c2<-ci+c1 r2<-r+r1 if(r2>0 && r2<=y1 && c2>0 && c2<=x1){ value[i]<-landscape[r2,c2] rows[i]<-r2 cols[i]<-c2 } i<-i+1 } } } if ((value[13])== low.cost | (value[13])== preyE){ landscape[r,ci]<- low.cost time.asc[w,]<-time.asc[w,] + 1 rr <-sample(c(8,9,13,17,18),1, replace=T) while( value[rr] == def.prm){ rr <-sample(c(8,9,13,17,18),1, replace=T) } r <- rows[rr] ci <- cols[rr] landscape[r,ci]<-NA } else { if ((max(value)== low.cost &&(value[13])== devE) | (max(value)== preyE &&(value[13])== devE) ){ landscape[r,ci]<- low.cost time.asc[w,]<-time.asc[w,] + 1 rr <-sample(c(7,8,17,16),1, replace=T) while( value[rr] == def.prm){ rr <-sample(c(7,8,17,16),1, replace=T) } r <- rows[rr] ci <- cols[rr] landscape[r,ci]<-NA } else { if ((value[13])== m.cost){ landscape[r,ci]<- value[8] time.asc[w,]<-time.asc[w,] + 1 rr <-13 r <- rows[rr] ci <- cols[rr] landscape[r,ci]<-NA } else { if ((value[13])== h.cost){ landscape[r,ci]<- h.cost time.asc[w,]<-time.asc[w,] + 1 rr <- 13 while( value[rr] == def.prm){ rr <-sample(c(8,17),1, replace=T) } r <- rows[rr] ci <- cols[rr] landscape[r,ci]<-NA } else { if ((value[13])== def.prm && (max(value)== h.cost)){ landscape[r,ci]<- h.cost rr <-sample(c(8,17),1, replace=T) while( value[rr] == def.prm){ rr <-sample(c(8,17),1, replace=T) } r <- rows[rr] ci <- cols[rr] landscape[r,ci]<-NA } } } } } } } #image(landscape, col=c(topo.colors(20),heat.colors(8),"black")) # jpeg(file = paste("C://Users//Utente//Desktop//Video//",gen,"myplot.jpg")) # plot(prey.cord, xlim=c(1,y1),ylim=c(1,x1), pch=19, col="black", cex=0.2, xlab="x", ylab="y") # lines(xy.cord, type="l", col="red", lwd=2) # points(dev.co, col="blue",pch=15) #dev.off() } #} sim<-data.frame(time.dis,time.search,time.asc,food1) write.table(sim, paste("/home/severin/MariannaR/shortHDRanDistSmallDev.TimeDive.",nam[ss],".P3" , ".csv", collapse=NULL, sep=""), sep=",", col.names=T,row.names=F) } }