Dear Cyrille. The main reason for this limitation is probably that it's not totally obvious that the fitted DFT-D parameters also work in the DFT+U case. But perhaps they will, and you can easily include a DFT-D2 or DFT-D3 correction by adding it to the LCAOCalcultor object by hand.
Simply combine the calculators from 2 scripts, one that uses DFT+U, and one that uses DFT-D. See attached script, which is straight from the VNL Scripter, except that I have manually added the DFT-D2 correction for PBE to the calculator:
correction_extension = GrimmeDFTD2(
global_scale_factor=0.75,
damping_factor=20.0,
maximum_neighbour_distance=30.0*Ang,
element_parameters={
Carbon: [1.452*Ang, 1.75*J*nm**6/mol],
},
)
calculator = LCAOCalculator(
basis_set=basis_set,
exchange_correlation=exchange_correlation,
numerical_accuracy_parameters=numerical_accuracy_parameters,
correction_extension=correction_extension,
)