We have developed robust and efficient numerical methods for solving integral equation theories for electrolyte solutions. These methods are hybrids of Newton-Raphson and Picard iterations and have been obtained as extended versions of the previous methods for pure solvents by solving nontrivial problems posed by the inclusion of ions. Bulk electrolytes and electrolytes near both inert and metallic surfaces are considered. The basic equations previously derived for a one-component fluid near a planar wall are extended to a multicomponent fluid. Analytical expressions for elements of the Jacobian matrices are arranged in compact form. A striking feature of the method for surface problems is that the Jacobian is determined only from bulk properties. A discussion of some special treatments that need to be considered for asymmetric anions and cations is included. These methods have been demonstrated using the full reference hypernetted-chain theory for various sizes of ions in a wide range of ionic concentrations. (C) 1996 Academic Press, Inc.