In this paper we study strongly singular problems with Dirichlet boundary condition on bounded domains given by \(\begin{aligned} -\operatorname {div} \left( |\nabla u|^{p-2}\nabla u+\mu (x)|\nabla u|^{q-2}\nabla u \right) = \frac{h(x)}{{\left( u^+\right) }^r} \quad \text {in } \Omega , \end{aligned}\) where \(1<p<N\) , \(p<q<p^*=\frac{Np}{N-p}\) , \(0 \le \mu (\cdot ) \in L^\infty (\Omega )\) , \(1<r\) and \(h\in L^1(\Omega )\) with \(h(x)>0\) for a.a. \(x\in \Omega \) . Since the exponent r is larger than one, the corresponding energy functional is not continuous anymore and so the related Nehari manifold \(\begin{aligned} \mathcal {N} = \left\{ u \in W^{1,\mathcal {H}}_0(\Omega ):\Vert \nabla u\Vert _p^p+\Vert \nabla u\Vert _{q,\mu }^q- \int _\Omega h(x){\left( u^+\right) }^{1-r} \,\textrm{d}x = 0\right\} \end{aligned}\) is not closed in the Musielak-Orlicz Sobolev space \(W^{1,\mathcal {H}}_0(\Omega )\) . Instead we are minimizing the energy functional over the constraint set \(\begin{aligned} \mathcal {M} = \left\{ u \in W^{1,\mathcal {H}}_0(\Omega ):\Vert \nabla u\Vert _p^p+\Vert \nabla u\Vert _{q,\mu }^q- \int _\Omega h(x){\left( u^+\right) }^{1-r} \,\textrm{d}x \ge 0\right\} , \end{aligned}\) which turns out to be closed in \(W^{1,\mathcal {H}}_0(\Omega )\) and prove the existence of at least one weak solution. Our result is even new in the case when the weight function \(\mu \) is away from zero.