MPA in 2D systems: BZ sampling?

[csk]

Dear all!

When trying to use the MPA with 2D systems (MoS2 monolayer in this case), I encounter the following behavior (gap and VBM/CBM at the Gamma point!):

Convergence w.r.t. the BZ sampling within the Godby-Needs PPM:
kpt bnd_x bnd_g cutoff (Ry) gap vbm cbm
6x6x1 200 200 4 3.91217 -0.87991 3.03226
12x12x1 200 200 4 3.88638 -0.83605 3.05033

Using the RIM and RIM-W parameters for 2D systems, the QP energies are already converged (to within 50 meV here) with a 12x12x1 k-mesh, as it should be. Now the same thing for MPA with 'ETStpsXm = 8' (converged vs. 16 poles):

kpt bnd_x bnd_g cutoff (Ry) gap vbm cbm
6x6x1 600 600 16 5.12461 -1.82519 3.29942
6x6x1 200 200 4 5.20877 -1.37231 3.83646
12x12x1 200 200 4 4.37506 -0.96955 3.40551
24x24x1 200 200 4 3.98923 -0.79276 3.19647
48x48x1 200 200 4 3.82818 -0.66470 3.16348

The first line, in comparison to the second line, shows the convergence w.r.t. the bands (both for X and G) and the cutoff, both using a 6x6x1 k-mesh.

What am I doing wrong here? Why is the convergence w.r.t. the k-mesh much slower in this case? Is there something else to consider for the MPA in 2D systems? I attach the QE and Yambo input files, as well as the respective report files to this post.

Thanks for your explanations!
Christian