The bottom panels of Figure (a)-(h) compare the experimental and computed XAS spectra of our benchmark systems. The computed XAS spectra are shifted to align with the respective first experimental peak; the shifts are given in the figure legends. These shifts are <1.1 eV and significantly smaller than the corresponding shifts computed at the lower TDDFT/CAM-B3LYP 100% level of theory in Ref. 63. The XAS spectra yielded by (fc-CVS-)EOM-EE-CCSD show only small discrepancies with respect to experimental spectra in terms of peak positions and intensities of dominant peaks, compared to similar discrepancies with the TDDFT approach in Ref. 63, further highlighting its robustness relative to TDDFT in modeling XAS. Below, all computed energy values are shifted values. The SI provide the non-shifted peak po- The top two panels of Figure (a)-(h) provide the X2PA spectra computed, respectively, with the CVS-0 and CVS-uS resolvents within the (fc-CVS-)EOM-EE-CCSD framework for the benchmark systems (shifted similar to XAS). Compared to CVS-0, the intensities of the spectra computed with CVS-uS are one order of magnitude smaller across all systems. This is a consequence of cancellation of contributions from the valence states as intermediates and the core states as intermediates. Further, the shapes of these spectra with the two approaches show significant differences. We also note that the magnitudes of the X2PA cross sections (∼10 -6 ) in a.u. are orders of magnitude smaller than typical valence 2PA cross sections for these systems, a consequence of small overlaps of core orbitals with unoccupied orbitals. However, the X2PA moments used for computing the cross sections are greater than the precision thresholds set in our calculations and, therefore, these moments are numerically reliable.