The present paper describes an evaluation method of the shape derivatives, which is defined as the Frechet derivative with respect to the domain variation, of cost functions in a shape optimization problem of a domain in which a boundary value problem of partial differential equation is defined. Instead of the known method using the equation of boundary integral type, an equation of domain integral type is proposed to evaluate the shape derivative. This method is applied to the mean compliance minimization problem of linear elastic body. Numerical example shows the validity of the method.