Using the description of hypermaps with matchings, Goulden and Jackson have given an expression of the generating series of rooted bipartite maps in terms of the zonal polynomials. We generalize this approach to the case of constellations on non-oriented surfaces that have recently been introduced by Chapuy and Dołęga. A key step in the proof is an encoding of constellations with tuples of matchings. We consider a one parameter deformation of the generating series of constellations using Jack polynomials and we introduce the coefficients $c^λ_{μ^0,...,μ^k}(b)$ obtained by the expansion of these functions in the power-sum basis. These coefficients are indexed by $k+2$ integer partitions and the deformation parameter $b$, and can be considered as a generalization for $k\geq1$ of the connection coefficients introduced by Goulden and Jackson. We prove that when we take some marginal sums, these coefficients enumerate $b$-weighted $k$-tuples of matchings. This can be seen as an "unrooted" version of a recent result of Chapuy and Dołęga for constellations. For $k=1$, this gives a partial answer to Goulden and Jackson Matching-Jack conjecture. Lassale has formulated a positivity conjecture for the coefficients $θ^{(α)}_μ(λ)$, defined as the coefficient of the Jack polynomial $J_λ^{(α)}$ in the power-sum basis. We use the second main result of this paper to give a proof of this conjecture in the case of partitions $λ$ with rectangular shape.