A mixed quantum-classical method aimed at the study of nonadiabatic dynamics in the presence of external electromagnetic fields is developed within the framework of time-dependent density functional theory. To this end, we use a trajectory-based description of the quantum nature of the nuclear degrees of freedom according to Tully's fewest switches trajectories surface hopping, where both the nonadiabatic coupling elements between the different potential energy surfaces, and the coupling with the external field are given as functionals of the ground-state electron density or, equivalently, of the corresponding Kohn-Sham orbitals. The method is applied to the study of the photodissociation dynamics of some simple molecules in gas phase.TAVERNELLI, CURCHOD, AND ROTHLISBERGER PHYSICAL REVIEW A 81, 052508 (2010) expand the system wave function at a given time t in a linear combination of static electronic wave functions, (r,R,t) = ∞ J J (r; R) J (R,t), (2.4) with time-dependent coefficients J (R,t), which in the Born-Oppenheimer limit correspond to nuclear wave functions. { J (r; R)} describes a complete set of electronic basis functions that are solutions of the time-independent Schrödinger equation,Ĥ el (r; R) J (r; R) = E el J (R) J (r; R), (2.5) and depend parametrically on the nuclear coordinates R. E el J (R) = H J J (R) is called the J th potential energy surface (PES), which is a function of the nuclear coordinates R. Inserting Eq. (2.4) in the time-dependent Schrödinger Eq. (2.1) and multiplying from the left by * J (r; R) we get, after integration over r, ih ∂ J (R,t) ∂t = − γh 2 2M γ ∇ 2 γ J (R,t) + I H J I (R) I (R,t) + γ Ih 2 2M γ D γ J I (R) I (R,t) − γ,I =Jh 2 M γ d γ J I (R)∇ γ I (R,t), (2.6) 052508-2 MIXED QUANTUM-CLASSICAL DYNAMICS WITH TIME-. . . PHYSICAL REVIEW A 81, 052508 (2010) where H J I (R) = * J (r; R)Ĥ el I (r; R)d r. (2.7) d γ J I (R), the first-order coupling vectors, are defined as d γ J I (R) = { * J (r; R)[∇ γ I (r; R)]}d r, (2.8) and D γ J I (R), the second-order coupling elements, are given by D γ J I (R) = * J (r; R) ∇ 2 γ I (r; R) d r.(2.9)