How can I run this parallel code in R/RStan?

Question

i want to run a code, but i have an error, (both descrives below). The main idea is to run each cycle independently, because doing it sequentially takes a long time, plus I'm running the code on a PC with win10 and 24 cores.

the code:

library(readxl)
library(rstan)
library(ggplot2)
library(shinystan)
library(parallel)
library(foreach)

options(mc.cores = parallel::detectCores())


setwd('### dir of data')

# data

U<- read_excel(path = 'My_data_cell.xlsx',sheet = NULL,range = NULL)
Um=as.matrix(U)
Um_datos=Um[,-1]
Metadata<-read_excel(path='My_other_data_cell.xlsx',sheet = NULL,range = NULL)
MData=as.matrix(Metadata)
MData_datos=MData[,-1]
freq=Um[,1]
LC<-length(freq)
num=as.character(1:LC)

# stan Code

C_model <- stan_model("My_code.stan",
                      model_name='model_1',
                      auto_write = TRUE,
                      verbose    = TRUE)

# trying to parallel 

cl <- parallel::makeCluster(20)
doParallel::registerDoParallel(cl)

time_in1<-Sys.time()
#:LC
foreach(i =1:LC) %dopar% {
  
  
  #usefull data
  
  data_Y<-log(Um_datos[i,])
  nonNAdata=!is.na(data_Y)
  data_y=data_Y[nonNAdata]
    
  Mag<-as.numeric(MData_datos[1,])
  Mag=Mag[nonNAdata]
  
  Vs<-as.numeric(MData_datos[2,])
  Vs=Vs[nonNAdata]
  
  f0<-as.numeric(MData_datos[3,])
  f0=Vs[nonNAdata]
  
  Rrup<-as.numeric(MData_datos[4,])
  Rrup=Rrup[nonNAdata]
  
  F_event=as.numeric(MData_datos[5,])
  F_event=F_event[nonNAdata]
  
  #ID EST and EVN #####
  EVNId<-as.numeric(MData_datos[6,])  #ID_EQ
  EVNID=EVNId[nonNAdata]
  
  ESTId<-as.numeric(MData_datos[7,]) #ID_EST
  ESTID=ESTId[nonNAdata]
 
  
  ESTID_nr<-ESTID[!duplicated(ESTID)]
  EVNID_nr<-EVNID[!duplicated(EVNID)]
  
  eqid=EVNID
  statid=ESTID
  
  neq=max(as.numeric(EVNID_nr))
  nstat=max(as.numeric(ESTID_nr))
  
  
  data_list <- list(N = length(Mag),
                    NEQ = neq,
                    NSTAT = nstat,
                    M = Mag,
                    R = Rrup,
                    f0 = f0,
                    VS = Vs,
                    Fevent = F_event,
                    Y = data_y,
                    freq = freq[i],
                    idx_eq = eqid,
                    idx_st = statid)
    
  fit <- sampling(C_model,
                  data = data_list,
                  cores = 20,
                  chains = 1,
                  iter = 100,
                  warmup= 200,
                  init_r = 0.5, 
                  control = list(adapt_delta = 0.85, max_treedepth = 16),
                  seed=123)
  
  nombre=paste('fit_NUM_',num[i],'.rds',sep="")
  saveRDS(fit,nombre)

}

time_end1<-Sys.time()
deltatime1=time_in1-time_end1

the error:

Error in { : task 1 failed - "could not find function 'sampling'"

What is the correct way to run a STAN code in several loops but parallelizing each one?