Thank you for your help. But after I re-modified the code, I found some problems. How should nested structures be parallelized? And how should global variables be passed in foreach? In the code below I simulate the data I am using and my updated code. But I'm still confused about how to convert to foreach for parallel computation. Can you help me?
#######import package####
library(ff)
library(raster)
####Here is the simulated data#########
TDD_MAT<-matrix(1, nrow = 50, ncol = 50)
FDD_MAT<-matrix(1, nrow = 50, ncol = 50)
for(k in 1:24){
matrix<-matrix(runif(250,min = -15+k,max = -14+k), nrow = 50, ncol = 50)
raster<-raster(matrix)
extent(raster)<-extent(1,50,1,50)
raster_stack<-addLayer(raster_stack,raster,k)
}
MAT_2015 <- ff(vmode="double",dim=c(ncell(raster_stack),nlayers(raster_stack)),filename=paste0(getwd(),"/stack.ffdata"))
for(i in 1:nlayers(raster_stack)){
MAT_2015[,i] <- raster_stack[[i]][]
}
doy<-c(2,4,6,35,53,65,72,105,125,132,145,165,178,184,195,243,255,265,288,293,302,333,355,362)
#########################The above is a data simulation#############
##########################The following is a for loop that wants to calculate in parallel###
for (i in 1:10000) {
if(any(is.na(MAT_2015[i,]))){
FDD_MAT[i]<-NA ##NA cases do not appear in simulated data
TDD_MAT[i]<-NA
}else{
day<-doy ##doy is a global variable I defined,just like c(2,65,160,243,325)
x<-c(day*6.28/365)
y<-MAT_2015[i,]
fm <- nls(y ~ cbind(a = 1, b = sin(x + c)), start = list(c=30), alg = "plinear"); ##perform regression in the specified format,to get c(a,b,c)
yy<-fm$m$getAllPars() ##get the coefficients
TDD<-0
FDD<-0
for(j in 1:365){
y<-yy[2]+yy[3]*sin((j*6.28/365)+yy[1])
if(y>0){
TDD<-TDD+y}
else{
FDD<-FDD+y}
}
FDD_MAT[i]<-FDD
TDD_MAT[i]<-TDD
}
}