  We propose a new loss function for supervised and physics-informed training
of neural networks and operators that incorporates a posteriori error estimate.
More specifically, during the training stage, the neural network learns
additional physical fields that lead to rigorous error majorants after a
computationally cheap postprocessing stage. Theoretical results are based upon
the theory of functional a posteriori error estimates, which allows for the
systematic construction of such loss functions for a diverse class of
practically relevant partial differential equations. From the numerical side,
we demonstrate on a series of elliptic problems that for a variety of
architectures and approaches (physics-informed neural networks,
physics-informed neural operators, neural operators, and classical
architectures in the regression and physics-informed settings), we can reach
better or comparable accuracy and in addition to that cheaply recover
high-quality upper bounds on the error after training.
