  The use of orthonormal polynomial bases has been found to be efficient in
preventing ill-conditioning of the system matrix in the primal formulation of
Virtual Element Methods (VEM) for high values of polynomial degree and in
presence of badly-shaped polygons. However, we show that using the natural
extension of a orthogonal polynomial basis built for the primal formulation is
not sufficient to cure ill-conditioning in the mixed case. Thus, in the present
work, we introduce an orthogonal vector-polynomial basis which is built ad hoc
for being used in the mixed formulation of VEM and which leads to very
high-quality solution in each tested case. Furthermore, a numerical experiment
related to simulations in Discrete Fracture Networks (DFN), which are often
characterised by very badly-shaped elements, is proposed to validate our
procedures.
