The outline of this work is as follows. In Section II, the setup of the system under consideration is briefly described. The canonical rMC simulations with multivalent anions are reviewed in Section III. Next, the theoretical approach of canonical CR is discussed in Section IV. Results are shown and discussed in Section V, followed by the closing remarks and perspective of new developments in Section VI. Finally, the detailed calculations of the electrostatic correlations contributions to the free energy, which are included in the modified PB equation, are provided in the Appendix A.