What is the covariance of $\min(X,Y)$ and $\min(X,Z)$ when $X$, $Y$, $Z$ are i.i.d. standard normal? Ecological growth can be thought of as being determined by the most limiting factor, and I'd like to have a solution for the covariance of growth records subject to common and independent limitations.
The exact question is the $\operatorname{cov}(\min(X,Y), \min(X,Z))$ where $X,Y,Z$ are each iid and $N(0,1)$. 
This is similar to an earlier question, but now with a double min: Compute $\operatorname{cov}(X, \max(X,Y))$ and $\operatorname{cov}(X, \min(X,Y))$ when $X$ and $Y$ are standard normal
I had thought to parallel the answer posted there by Did, letting $A=\min(X,Y)$ and $B=\min(X,Z)$, such that the covariance of $A$ and $B$ is, 
$$\mathbb E(AB)-\mathbb E(A)\mathbb E(B).$$  
Using $A=X\mathbf 1_{X\gt Y}+Y\mathbf 1_{X\lt Y}$ and $B=X\mathbf 1_{X\gt Z}+Z\mathbf 1_{X\lt Z}$ permits for partitioning expectations, but I'm stuck because many of the terms no longer go to zero.
 A: The result is: 

$$\color{red}{\mathrm{cov}(\min(X,Y),\min(X,Z))=\frac13-\frac{2-\sqrt3}{2\pi}\approx0.290687895}$$


To show this, note that the identities $$\min(X,Y)=X\mathbf 1_{X<Y}+Y\mathbf 1_{Y<X}\qquad\min(X,Z)=X\mathbf 1_{X<Z}+Z\mathbf 1_{Z<X}$$ yield $$E(\min(X,Y)\min(X,Z))=u+v+w$$ with $$u=E(X^2;X<Y,X<Z)\qquad v=E(YZ;Y<X,Z<X)$$ and $$w=E(XY;Y<X<Z)+E(XZ;Z<X<Y)$$ Considering the common CDF $F$ of $X$, $Y$ and $Z$, one sees that $$u=E(T^2)\qquad\text{with}\ T=X(1-F(X))$$ Using the symmetry $(X,Y,Z)\to(X,Z,Y)$, $$w=2E(XY;Y<X<Z)=E(XY;X<Z,Y<Z)=v$$ 
For every $z$, $$E(XY;X<z,Y<z)=G(z)^2\qquad\text{with}\ G(z)=E(X;X<z)$$ hence $$w=v=E(G(X)^2)$$ Likewise, $$E(\min(X,Y))=E(X;X<Y)+E(Y;Y<X)=2E(T)$$ hence 

$$\mathrm{cov}(\min(X,Y),\min(X,Z))=E(X^2(1-F(X))^2)+2E(G(X)^2)-4E(X(1-F(X)))^2$$ 

where we recall that $$G(x)=E(X;X<x)$$

The preceding steps hold for any i.i.d. triplet. From now on, we use that $X$ is standard normal, with PDF $\varphi$ and CDF $\Phi$, hence $E(X)=0$ and $$E(T)=-E(X\Phi(X))=-\int_\mathbb R x\Phi(x)\varphi(x)dx$$ Since $\Phi'(x)=\varphi(x)$ and $\varphi'(x)=-x\varphi(x)$, an integration by parts yields $$E(T)=\left.\varphi(x)\Phi(x)\right|_{-\infty}^\infty-\int_\mathbb R\varphi(x)^2dx=-\int_\mathbb R\frac1{\sqrt{2\pi}}\varphi(\sqrt2x)dx=-\frac1{2\sqrt{\pi}}$$ Likewise, $X$ is symmetric hence $1-\Phi(X)=\Phi(-X)$ and $$E(T^2)=E(X^2\Phi(X)^2)=\int_\mathbb R x^2\varphi(x)\Phi(x)^2dx$$ Using $\varphi''(x)=(x^2-1)\varphi(x)$, one gets $E(T^2)=s+t$ with $$s=\int_\mathbb R \varphi''(x)\Phi(x)^2dx=\left.\varphi'(x)\Phi^2(x)\right|_{-\infty}^\infty-\int_\mathbb R\varphi(x)\cdot2\varphi'(x)\Phi(x)dx=-\int_\mathbb R(\varphi(x)^2)'\Phi(x)dx$$ hence $$s=\left.-\varphi(x)^2\Phi(x)\right|_{-\infty}^\infty+\int_\mathbb R\varphi(x)^3dx=\int_\mathbb R\frac1{2\pi}\varphi(\sqrt3x)dx=\frac1{2\pi\sqrt3}$$
and $$t=\int_\mathbb R \varphi(x)\Phi(x)^2dx=\left.\frac13\Phi(x)^3\right|_{-\infty}^\infty=\frac13$$ Finally, $$G(x)=\int_{-\infty}^xz\varphi(z)dz=-\varphi(x)$$ hence $$E(G(X)^2)=\int_\mathbb R\varphi(x)^3dx=\int_\mathbb R\frac1{2\pi}\varphi(\sqrt3x)dx=\frac1{2\pi\sqrt3}$$
Summing these contributions, one gets the numerical value at the beginning of this post.
