the functional defined on $\M({\overline}{\Omega},{\mathbf{R}}^2)\to{\mathbf{R}}_+\cup\{+\infty\}$ as $${\mathcal{E}}_f(\s):=
\begin{dcases}
\int_Mf(\t)\dH^1&\text{if }\s=U(\t,\xi,M),\\
\quad +{\infty}&\text{in the other cases.}
\end{dcases}$$ Given a sequence of $N+1$ distinct points $S=(x_0,\dots,x_N)\in {\Omega}^{N+1}$, we consider the minimization of ${\mathcal{E}}_f(\s)$ for $\s\in \M({\overline}{\Omega},{\mathbf{R}}^2)$ satisfying the constraint $$\label{eq:divconstraint1}
{\nabla}\cdot\s=N\d_{x_0}-\sum_{i=1}^N\d_{x_i}\quad\mbox{ in }{\mathcal{D}}'({\mathbf{R}}^2).$$ The distributional support of such $\sigma$ connects the source in $x_0$ to the sinks in $x_1,\cdots,x_N$. In general, a model for branched transport connecting a set of sources to a set of sinks (represented by two discrete probabilistic measures supported on a set of points in ${\Omega}$) is obtained by choosing $f(\t)=|\t|^\a$ with $0<\a<1$ and minimizing the associated functional under a divergence constraint similar to equation . The direct numerical optimization of the functional ${\mathcal{E}}_f$ is not easy because we do not know [*a priori*]{} the topological properties of the tree $M$. For this reason it is interesting to optimize an “approximate" functional defined on more flexible objects such as functions. Such approximate model has been introduced in [@OS2011] where the authors study the $\Gamma$-convergence (see [@Bra2]) of a family of functionals inspired by the well known work of Modica and Mortola [@Mod_Mort]. Another effort in this direction can be found in the work [@Bon_Lem_San] where is studied an 