Multiple Hermite polynomials are an extension of the classical Hermite polynomials for which orthogonality conditions are imposed with respect to r > 1 normal (Gaussian) weights w j (x) = e −x 2 +c j x with different means c j /2, 1 ≤ j ≤ r. These polynomials have a number of properties, such as a Rodrigues formula, recurrence relations (connecting polynomials with nearest neighbor multi-indices), a differential equation, etc. The asymptotic distribution of the (scaled) zeros is investigated and an interesting new feature happens: depending on the distance between the c j , 1 ≤ j ≤ r, the zeros may accumulate on s disjoint intervals, where 1 ≤ s ≤ r. We will use the zeros of these multiple Hermite polynomials to approximate integrals of the form ∞ −∞ f (x) exp(−x 2 + c j x) dx simultaneously for 1 ≤ j ≤ r for the case r = 3 and the situation when the zeros accumulate on three disjoint intervals. We also give some properties of the corresponding quadrature weights. * Work supported EOS project 30889451 and FWO project G.0864.16N