This paper presents a complete set of coupled equations that govern the bubble transport in three-dimensional gravity driven flow. The model accounts for bubbles growth or shrinkage due to pressure and temperature changes as well as for multiple gas diffusion in and out of the bubbles but neglects bubble coalescence, break-up, and nucleation. The model applies to glass melting furnaces but it could be extended to other two-phase flow applications such as metal and polymer processings, passive cooling systems, and two-phase flow around naval surface ship. Governing equations are given for the key variables which are, in the present case, (1) the refining agent concentration, (2) the gas species dissolved in the liquid phase, and (3) the bubble radius, gas molar fraction, and density function. The method of solution is discussed and based on the modified method of characteristics.