We analyze the flux conservation property of the finite element method. It is shown that the finite element solution does approximate the flux locally in the optimal order, i.e., the same order as that of the nodal interpolation operator. We propose two methods, post-processing the finite element solutions locally. The new solutions, remaining as optimal-order solutions, are flux-conserving elementwise. In one of our methods, the processed solution also satisfies the original finite element equations. While the high-order finite volume schemes are still under construction, our methods produce finite-volume-like finite element solution of any order. In particular, our methods avoid solving non-symmetric finite volume equations. Numerical tests in 2D and 3D verify our findings.