Abstract: We present direct methods, algorithms, and symbolic software for the computation of conservation laws of nonlinear partial differential equations (PDEs) and differential-difference equations (DDEs). Our method for PDEs is based on calculus, linear algebra, and variational calculus. First, we compute the dilation symmetries of the given nonlinear system. Next, we build a candidate density as a linear combination with undetermined coefficients of terms that are scaling invariant. The variational derivative (Euler operator) is used to derive a linear system for the undetermined coefficients. This system is then analyzed and solved. Finally, we compute the flux with the homotopy operator. The method is applied to nonlinear PDEs in (1+1) dimensions with polynomial nonlinearities which include the Korteweg-de Vries (KdV), Boussinesq, and Drinfel'd-Sokolov-Wilson equations. An adaptation of the method is applied to PDEs with transcendental nonlinearities. Examples include the sine-Gordon, sinh-Gordon, and Liouville equations. For equations in laboratory coordinates, the coefficients of the candidate density are undetermined functions which must satisfy a mixed linear system of algebraic and ordinary differential equations. For the computation of conservation laws of nonlinear DDEs we use a splitting of the identity operator. This method is more efficient that an approach based on the discrete Euler and homotopy operators. We apply the method of undetermined coefficients to the Kac-van Moerbeke, Toda, and Ablowitz-Ladik lattices. To overcome the shortcomings of the undetermined coefficient technique, we designed a new method that first calculates the leading order term and then the required terms of lower order. That method, which is no longer restricted to polynomial conservation laws, is applied to discretizations of the KdV and modified KdV equations, and a combination thereof. Additional examples include lattices due to Bogoyavlenskii, Belov-Chaltikian, and Blaszak-Marciniak. The undetermined coefficient methods for PDEs and DDEs have been implemented in Mathematica. The code TransPDEDensityFlux.m computes densities and fluxes of systems of PDEs with or without transcendental nonlinearities. The code {\tt DDEDensity\\Flux.m} does the same for polynomial nonlinear DDEs. Starting from the leading order terms, the new {\tt Maple} library {\tt discrete} computes densities and fluxes of nonlinear DDEs. The software can be used to answer integrability questions and to gain insight in the physical and mathematical properties of nonlinear models. When applied to nonlinear systems with parameters, the software computes the conditions on the parameters for conservation laws to exist. The existence of a hierarchy of conservation laws is a predictor for complete integrability of the system and its solvability with the Inverse Scattering Transform.