  MHD simulations show that the magnetic field can drive molecular outflows
during the formation of massive protostars. The best probe to observationally
measure both the morphology and the strength of this magnetic field at scales
of 10-100 au is maser polarization. We measure the direction of magnetic fields
at mas resolution around a sample of massive star-forming regions to determine
whether there is a relation between the orientation of the magnetic field and
of the outflows. In addition, by estimating the magnetic field strength via the
Zeeman splitting measurements, the role of magnetic field in the dynamics of
the massive star-forming region is investigated. We selected a flux-limited
sample of 31 massive star-forming regions to perform a statistical analysis of
the magnetic field properties with respect to the molecular outflows
characteristics. We report the linearly and circularly polarized emission of
6.7 GHz CH3OH masers towards seven massive star-forming regions of the total
sample with the EVN. The sources are: G23.44-0.18, G25.83-0.18, G25.71-0.04,
G28.31-0.39, G28.83-0.25, G29.96-0.02, and G43.80-0.13. We identified a total
of 219 CH3OH maser features, 47 and 2 of which showed linearly and circularly
polarized emission, respectively. We measured well-ordered linear polarization
vectors around all the massive YSOs and Zeeman splitting towards G25.71-0.04
and G28.83-0.25. Thanks to recent theoretical results, we were able to provide
lower limits to the magnetic field strength from our Zeeman splitting
measurements. We further confirm (based on ~80% of the total flux-limited
sample) that the magnetic field on scales of 10-100 au is preferentially
oriented along the outflow axes. The estimated magnetic field strength of
|B_{||}|>61 mG and >21 mG towards G25.71-0.04 and G28.83-0.2, respectively,
indicates that it dominates the dynamics of the gas in both regions.
