I have two endogenous latent variables
Your syntax is correct, except for one detail. How you set the residual-covariance parameters is to calculate the covariance among their explained variances (i.e., expected values given the common-factor effects):
Even if these were the correct values, it would be sufficient to simply take the covariance matrix returned by cov()
cov(factors[ , c("endog1","endog2")])
But the residual covariance is the unexplained variances. You can just choose any values you want. If the explained variances are < 1 and you want the population slopes to be interpretable as standardized slopes, then you can use an identity matrix to calculate 1 minus the explain variance:
exoMat <- cov(factors[ , c("endog1","endog2")])
## subtract from 2 × 2 identity matrix
endoMat <- diag(2) - exoMat
## choose a residual correlation
resCor <- 0.3 # or whatever value you want
## convert to residual covariance by multiplying by the SDs
endoMat[1,2] <- endoMat[2,1] <- resCor * prod(sqrt(diag(endoMat)))
Then this line will correctly simulate the latent-variable residuals: