The nonideal mixing of phosphatidylserine (PS) and phosphatidylcholine (PC) binary lipid mixtures was studied by computer simulation based on a model wherein the excess energy of mixing is divided between an electrostatic term and one adjustable term delta Em that includes all other nonideal interactions. The lateral distribution of the lipids and the energy of the mixtures were obtained by using Kawasaki relaxation in a canonical ensemble. The Gibbs free energies were calculated by Kirkwood's coupling parameter method. The simulation results are strongly dependent on simulation size for sizes smaller than about 1000 lipids. Nonideal interaction between lipids can result in large scale separation of lipid phases of different composition at reasonable delta Em values as well as clustering of like lipids. In plots of total Gibbs free energy of mixing versus PS mole fraction in PS/PC, the boundaries of the two phase region could be accurately determined. The electrostatic interaction influences cluster size and shape, and also the composition of phases in the two-phase region.