Regarding excited state dynamics, its non-adiabatic dynamical nature of IC and ISC requires intensive computations of excited states' energies, wave functions, and overlaps of the wave functions at adjacent time steps. Proper descriptions of various excited states are also challenging due to their distinctively different characters. Considering these, it is necessary to perform theoretical simulations using a methodology that: (i) provides a balanced and accurate treatment of both dynamic correlation and non-dynamic correlation (e.g. including the important double excitations) of the electronic states; (ii) is capable of correctly describing the conical intersections (CI) between ground and excited states as well as between different excited states; (iii) can provide efficient computations of both non-adiabatic and spin-orbit couplings; (iv) is efficient enough to perform statistical sampling for several picoseconds of large molecular systems. Perhaps, one of the most popular and practical electronic structure theories for ab initio NAMD is TDDFT, 7 although it produces topologically wrong CI 10 . Throughout this work, CI IJ is used to denote the conical intersection between adiabatic states I and J. Also, it was recently shown that both CASSCF and TDDFT fail to describe the S 2 -to-S 1 internal conversions (ICs) of s-trans-butadiene and s-trans-hexatriene via CI 21 , 9 which is due to the missing of either dynamic correlation or the doubly excited configurations.