Recent polariton quantum dynamics simulations are focusing on obtaining these expensive derivatives using semi-empirical electronic structure methods. For example, Zhang et al. derived the derivative of transition dipoles at the AM1/CIS level of theory 22 and applied it to the TSH simulations of stilbene molecule coupled to the cavity based on the Jaynes-Cummings (JC) and Tavis-Cummings (TC) models. Fregoni et al. derived the derivative of transition dipoles at the AM1/FOMO-CI level and performed the TSH simulations for azobenzene molecule based on the Rabi-type model (that excludes the permanent dipole moment and dipole selfenergy). Groenhof et al. approximated transition dipoles by a first-order Taylor expansion around the ground state equilibrium and obtained the derivatives through least-squares fitting with a large number geometrical conformations. The lack of accurate derivatives on dipoles (at the correlated wavefunction level of theory) had become the major theoretical bottleneck in performing ab initio on-the-fly simulations of polariton chemistry. To the best of our knowledge, there is no previous work on using correlated wavefunction level of theory to perform ab initio on-the-fly quantum dynamics simulations of polariton chemistry with the rigorous nuclear gradient. In this work, we construct the dipoles (including permanent and transition dipoles) and their derivatives for realistic molecules using machine-learning techniques. We apply the machine-learned dipoles and derivatives to ab initio on-the-fly polariton dynamics simulations of a realistic polyatomic molecular system, azomethane, coupled to an optical cavity. We use the accurate nuclear gradient expression derived from the Pauli-Fierz (PF) Quantum Electrodynamics Hamiltonian without making further approximations. The energies, gradients, and derivative couplings of the molecular systems are obtained from the on-the-fly simulations at the level of complete active space self-consistent field (CASSCF), which are used to compute the polariton energies and nuclear gradients. The derivatives of dipoles are obtained from the machine learning model, which is trained with the data obtained from the CASSCF level of electronic structure calculations. The cavity loss process is modeled with the Lindblad jump superoperator on the reduced density of the electronicphotonic quantum subsystem. Our results show the accuracy of the machine-learned dipoles and their usage in simulating polariton dynamics.