All variational parameters of the Jastrow factor were variationally optimized at the VMC level with the stochastic reconfiguration optimization method , while molecular orbital coefficients, basis set contraction coefficients, and exponents were kept fixed from the initial DFT calculation. The optimized VMC wave functions were taken as guiding functions in the FN-DMC calculation, in which we also employed the ccECP pseudopotentials to approximate core electrons for each atom, integrated with the determinant localization approximation (DLA) . Consequently, the FN-DMC calculations were computed at two time-steps of 0.050 and 0.025 (a.u.), using 12800 walkers divided into 300 blocks, each 100 steps long, for a total of 4 × 10 8 sampled configurations. For some systems, it was required to run an additional third calculation with a time-step of 0.015 a.u. to get statistical agreement within 1σ in the observed binding energies. In Fig. we displayed the time-step convergence against the CCSD(T) reference values. Both VMC wave function optimization and DMC calculations were performed with the QMeCha code .