In conclusion, we integrated CNEO with the QM/MM electrostatic embedding scheme and achieved the accurate and efficient incoporation of nuclear quantum effects, particularly quantum delocalization effects, in the QM region of QM/MM simulations. We applied the resulting CNEO QM/MM theory to the calculation of a phenol-water complex and a glutamic acid-glutamate complex in both gas phase and aqueous phase. We investigated the impact of nuclear quantum effects on both optimized geometries and molecular dynamics. For the neutral phenol-water complex, solvation reduces the hydrogen bond distance and the incorporation of nuclear quantum effects through CNEO-DFT leads to a shorter hydrogen bond than that of DFT simulations. In contrast, for the negatively charged glutamic acid-glutamate complex, solvation leads to a slight increase for the hydrogen bond distance as predicted by DFT and DFT QM/MM but a negligible change as predicted by CNEO-DFT and CNEO-DFT QM/MM. Through dynamics simulations, we observed much more frequent proton transfer in CNEO-DFT QM/MM simulations than in DFT QM/MM simulations for the glutamic acidglutamate complex in both gas phase and aqueous phase due to the incorporation of quantum nuclear delocalization effects. Additionally, the location of the shared proton predicted by CNEO QM/MM is in great agreement with the experimental observations. All of these results demonstrate the significant impact of the solvation environment and quantum nuclear delocalization effects, both of which are key features of the CNEO QM/MM approach. As an accurate and efficient method, CNEO QM/MM holds great promise for future investigation of hydrogenrelated processes in complex chemical and biological environments.