This study introduces a unified thermodynamic methodology for modeling chemical reactions in heterogeneous equilibria. The approach addresses complex multicomponent systems where both soluble and insoluble species coexist, enabling accurate prediction of mineral phase formation, particularly struvite and vivianite, under variable environmental conditions. A central innovation of the method lies in its ability to account for the simultaneous precipitation of multiple insoluble metal-containing species, often governed by competitive equilibria. Additionally, the methodology introduces a generalized chemical equation that captures the overall process, including metal ion hydrolysis, complex formation, ligand protonation, and associated reactions. By integrating derived thermodynamic functions with originally developed mass balance equations for solid phases, the model reflects the coupled and dynamic nature of chemical interactions in heterogeneous systems. Special emphasis is placed on the influence of pH and the presence of competing metal ions such as calcium, sodium, and potassium, which strongly affect mineral solubility and system buffering. The methodology enables the prediction of mineral formation as a function of initial wastewater composition and environmental parameters, offering valuable insight into the conditions favorable for contaminant removal through controlled precipitation. Application of the model to real wastewater systems confirms its effectiveness in identifying optimal pH values and ionic conditions for maximizing mineral recovery. This thermodynamic framework provides a robust predictive tool for optimizing wastewater treatment processes. It reduces the need for extensive laboratory experimentation and supports the development of sustainable treatment strategies that promote resource recovery and ensure chemical stability in treated effluents.