In the present work, model equations describing the simultaneous removal of heavy metal ions, such as Zn(II), Cu(II), and Cd(II) ions from aqueous solutions by adsorption with phosphogypsum (PG) were calculated, employing response surface methodology (RSM), concerning operating conditions such as contact time (X1), temperature (X2), and initial pH (X3). The metal ions–PG interactions were characterised by thermal analysis techniques (TGA-DTG), scanning electron microscopy (SEM), Fourier transform infrared (FT-IR), and energy-dispersive X-ray spectroscopy (EDX). The results proved that the equations of the determined model agree well with the quadratic equations with R2(Adj) values that reach up to 0.876 (p < 0.05 with a 95% confidence level). It has been demonstrated that from aqueous solutions, phosphogypsum able to simultaneously remove more than 65.51% of Cd(II), 45.57% of Cu(II), and up to 45.44% of Zn(II) ions, and the optimums were obtained at the set of parameters of a contact time of 28.79 min, the temperature of 40.62 °C, and pH 9.17. Under these conditions, phosphogypsum was more selective towards cadmium (II) ions. Furthermore, the calculated models could provide appropriate data to strengthen the adsorption process.