A differential evolution based methodology is introduced for the solution of elliptic partial differential equations (PDEs) with Dirichlet and/or Neumann boundary conditions. The solutions evolve over bounded domains throughout the interior nodes by minimization of nodal deviations among the population. The elliptic PDEs are replaced by the corresponding system of finite difference approximation, yielding an expression for nodal residues. The global residue is declared as the root-mean-square value of the nodal residues and taken as the cost function. The standard differential evolution is then used for the solution of elliptic PDEs by conversion to a minimization problem of the global residue. A set of benchmark problems consisting of both linear and nonlinear elliptic PDEs has been considered for validation, proving the effectiveness of the proposed algorithm. To demonstrate its robustness, sensitivity analysis has been carried out for various differential evolution operators and parameters. Comparison of the differential evolution based computed nodal values with the corresponding data obtained using the exact analytical expressions shows the accuracy and convergence of the proposed methodology.