The catenary equation can be derived by considering the minimization of the functionals,
\begin{align}
f[y]&=\int_0^{l_2}\mu gy\sqrt{1+\left(\frac{\mathrm dy}{\mathrm dx}\right)^2}\,\mathrm dx \tag{pot. e.}\\
g[y]&=\int_0^{l_2}\sqrt{1+\left(\frac{\mathrm dy}{\mathrm dx}\right)^2}\,\mathrm dx\tag{length}
\end{align}
where $l_2$ is the point we're hanging a chain of length $l_1$.
Next, we can define the functional,
$$
S[y]=f[y]+\lambda g[y]=\int_0^{l_2}\mathcal{L}(y,\,y',\,x)\,\mathrm{d}x
$$
where $\lambda$ is a parameter that helps us determine the correct chain length and $y'=\mathrm{d}y/\mathrm{d}x$. From this definition, we have that our Lagrangian is,
$$\mathcal{L}\triangleq\left(\mu gy+\lambda\right)\sqrt{1+(y')^2}.$$
So by applying the Legendre transform, we get our Hamiltonian,
\begin{align}
\mathcal{H}&=y'\frac{\partial\mathcal{L}}{\partial y'}-\mathcal{L} \\
&=y'\left(\mu gy+\lambda\right)\cdot\frac{y'}{\sqrt{1+(y')^2}}-\left(\mu gy+\lambda\right)\sqrt{1+(y')^2}.
\end{align}
If we then let the Hamiltonian be the constant total energy (density), $\varepsilon$, then we end up with,
$$\left(\frac{\mu gy+\lambda}{\varepsilon}\right)^2=1+(y')^2$$
Under the constant mass density case, we can trivially solve this to get that common $\propto\cosh(x)$ function.
In the non-constant density, the ODE is a little different,
$$\frac{\mathrm dy}{\mathrm dx}=\left(\left(\frac{g\mu(x)y(x)+\lambda}{\varepsilon}\right)^2-1\right)^{1/2}\tag{1}$$
and may be a bit more complicated to solve, possibly requiring some numerical solution to fit the two boundary conditions ($y(0)=y_0$ and $y(l_2)=y_2$) along with the condition that $g[y]=l_1$ (i.e., length is fixed).
Since you are interested in an effective mass density when given two different-but-uniform ropes, your mass density is,
$$
\mu(x)=\begin{cases} \mu_1 & x<x_p \\ \mu_2 & \text{otherwise} \end{cases}
$$
where $x_p$ is the joint. This seems to point you towards having two separate linear cases and an additional constraint that $y(x_p^-)=y(x_p^+)$ (i.e., the heights at the joint as measured from above and from below must match). This likely will require some numerical minimization methods to satisfy.
Once you have that solution, it may be possible to add an additional parameter of the equivalent mass density, $\mu_\text{equiv}$, by fitting the constant-density case with the additional condition that $y(x_p)=y_p$ (i.e., the curve must also pass through the joint), though I don't know if this would result in a solution or not.
