existence of feasible solutions for the outputs of interest is only given in two dimensions. The computed bound gaps are found to converge optimally, that is, at the same rate as the finite element approximation. An attractive feature of the proposed approach is that it allows for a data set to be generated that can be used to certify and document the computed bounds. Using this data set and a simple algorithm, the correctness of the computed bounds can be established without recourse to the original code used to compute them. In the present paper, only computational domains whose boundary is made up of straight sided segments and polynomially varying loads are considered. Two examples are given to illustrate the proposed methodology.