We consider PDE constrained nonparametric regression problems in which the parameter f is the unknown coefficient function of a second order elliptic partial differential operator L-f, and the unique solution u(f) of the boundary value problem L(f)u = g(1) on O, u = g(2) on partial derivative O, is observed corrupted by additive Gaussian white noise. Here O is a bounded domain in R-d with smooth boundary partial derivative O, and g(1), g(2) are given functions defined on O, partial derivative O, respectively. Concrete examples include L-f u = Delta u-2fu (Schrodinger equation with attenuation potential f) and L(f)u = div(f del u) (divergence form equation with conductivity f). In both cases, the parameter space F = {f is an element of H-alpha(O)vertical bar f > 0}, alpha > 0, where H-alpha (O) is the usual order alpha Sobolev space, induces a set of nonlinearly constrained regression functions {u(f): f is an element of F}. We study Tikhonov-type penalized least squares estimators (f) over cap for f. The penalty functionals are of squared Sobolev-norm type and thus (f) over cap can also be interpreted as a Bayesian "maximum a posteriori" estimator corresponding to some Gaussian process prior. We derive rates of convergence of (f) over cap and of u((f) over cap), to f, u(f), respectively. We prove that the rates obtained are minimax-optimal in prediction loss. Our bounds are derived from a general convergence rate result for nonlinear inverse problems whose forward map satisfies a modulus of continuity condition, a result of independent interest that is applicable also to linear inverse problems, illustrated in an example with the Radon transform.