For more than 150 years, scientists have advanced aberration theory to describe, analyze and eliminate imperfections that disturb the imaging quality of optical components and systems. Simultaneously, they have developed optical design methods for and manufacturing techniques of imaging systems with ever-increasing complexity and performance up to the point where they are now including optical elements that are unrestricted in their surface shape. These so-called optical freeform elements offer degrees of freedom that can greatly extend the functionalities and further boost the specifications of state-of-the-art imaging systems. However, the drastically increased number of surface coefficients of these freeform surfaces poses severe challenges for the optical design process, such that the deployment of freeform optics remained limited until today. In this paper, we present a deterministic direct optical design method for freeform imaging systems based on differential equations derived from Fermat's principle and solved using power series. The method allows calculating the optical surface coefficients that ensure minimal image blurring for each individual order of aberrations. We demonstrate the systematic, deterministic, scalable, and holistic character of our method with catoptric and catadioptric design examples. As such we introduce a disruptive methodology to design optical imaging systems from scratch, we largely reduce the"trial-and-error"approach in present-day optical design, and we pave the way to a fast-track uptake of freeform elements to create the next-generation high-end optics. We include a user application that allows users to experience this unique design method hands-on.