Not really sure that this is what you are asking for. Give it a check.
x<-rnorm(100)
y<-rnorm(100)
z<-rnorm(100)
xz<-x*z
yz<-y*z
co<-data.frame(nrow=100, ncol=2)
for (i in 1:100){
z[i]<-rnorm(1)
xz[i]<-z[i]*x[i]
yz[i]<-z[i]*y[i]
co[i,]<-c(cor(z,x), cor(z,y))
}