I have a 3 by 3 matrix of random variables and I would like to generate the covariance matrix of these 9 random variables. The restrictions are as follows:

\begin{aligned} \operatorname{Corr}\left(b_{u_{1}, d_{1}}, b_{u_{2}, d_{2}}\right) &=\rho_{1}, \text { if } u_{1} \neq u_{2}, d_{1} \neq d_{2} \\ \operatorname{Corr}\left(b_{u_{1}, d_{1}}, b_{u_{2}, d_{1}}\right) &=\rho_{2}, \text { if } u_{1} \neq u_{2} \\ \operatorname{Corr}\left(b_{u_{1}, d_{1}}, b_{u_{1}, d_{2}}\right) &=\rho_{3}, \text { if } d_{1} \neq d_{2} \\ \operatorname{SD}\left(b_{u_{1}, d_{1}}\right) &=\sigma \end{aligned}

where \rho and \sigma are parameters to be estimated. I write a function to check if two variables are in the same row or column:

```
# Assign covariances and variances
function fun_cov(p1::Tuple{Int64,Int64}, p2::Tuple{Int64,Int64}, rho1_est::Float64, rho2_est::Float64, rho3_est::Float64, sigma_est::Float64)
agreements = sum(p1 .== p2)
ifelse(agreements==2, sigma_est^2, ifelse(agreements==0, rho1_est * sigma_est^2, ifelse(p1[2]==p2[2], rho2_est * sigma_est^2, rho3_est * sigma_est^2)))
end
```

Then, I generate the covariance matrix:

```
ind = vec(collect(Iterators.product(1:N, 1:N))) # initialize the index
cov =[fun_cov(p1, p2, rho1_est, rho2_est, rho3_est, sigma_est) for p1 in ind, p2 in ind]
```

Here is the problem, I find out that some values of \rho and \sigma do not generate a covariance matrix. For example, when

```
rho1_est=0.6; rho2_est=0.4; rho3_est=0.2; sigma_est=2.0; N=3;
```

the `cov`

I get has negative eigenvalues. I am not sure what is wrong in my code and I do not see why this setting of correlations is not possible. In my code, I need to guess the value of these parameters to generate the covariance matrix, which is then used to compute the objective function. I wonder how should I make sure `cov`

is indeed a covariance matrix while satisfying the restrictions on correlations. Thanks in advance.