Predictive modelling in physical science and engineering is mostly based on solving certain partial differential equations where the complexity of solutions is dictated by the geometry of the domain. Motivated by the broad applications of explicit solutions for spherical and ellipsoidal domains, in particular, the Eshelby's solution in elasticity, we propose a generalization of ellipsoidal shapes called polynomial inclusions. A polynomial inclusion (or p-inclusion for brevity) of degree k is defined as a smooth, connected, and bounded body whose Newtonian potential is a polynomial of degree k inside the body. From this viewpoint, ellipsoids are identified as the only p-inclusions of degree two, and many fundamental problems in various physical settings admit simple closed-form solutions for general p-inclusions as for ellipsoids. Therefore, we anticipate that p-inclusions will be useful for applications including predictive materials models, optimal designs, and inverse problems. However, the existence of p-inclusions beyond degree two is not obvious, not to mention their direct algebraic parameterizations.
In this work, we explore alternative definitions and properties of p-inclusions in the context of potential theory. Based on the theory of variational inequalities, we show that p-inclusions do exist for certain polynomials, though a complete characterization remains open. We reformulate the determination of surfaces of p-inclusions as nonlocal geometric flows which are convenient for numerical simulations and studying geometric properties of p-inclusions. In two dimensions, by the method of conformal mapping we find an explicit algebraic parameterization of p-inclusions. We also propose a few open problems whose solution will deepen our understanding of relations between domain geometry, Newtonian potentials, and solutions to general partial differential equations. We conclude by presenting examples of applications of p-inclusions in the context of Eshelby inclusion problems and magnet designs.