Scopus (CiteScore 2022 =3.0, Q3) , ISC

Document Type : Original Research Article


1 School of Chemical Engineering, Iran University of Science and Technology, Tehran 16765-163, Iran

2 Materials and Nuclear Fuel Research School, Atomic Energy Organization, Tehran 11365-8486, Iran


In this study, a dynamic and equilibrium modeling of water distillation process for 18O isotope separation with multi-component distillation approach was studied. The mathematical model of this process was derived in order to obtain the equilibrium concentrations using isotopic exchange reactions through statistical methods based on stochastic exchanges and probability of molecular collisions. In this modeling, to calculate the liquid-vapor equilibriums, the separation factor (α) or volatility was considered as a function of column temperature. Then equations of this model were solved by numerical techniques in Matlab programming and distributions of isotopic water components were obtained in the distillation column under total reflux conditions. The simulation results showed that the separation factor has a significant effect on the 18O isotope separation as the H-labeled light components increased in the condenser above the column; the HD and D-labeled heavy components increased in the boiler due to the difference of the separation factor. Then experimental results of 18O in the condenser above the column and the reboiler from a semi-industrial plant were compared with the simulation results, which it showed a good conformity between the modeling and pilot results.

Graphical Abstract

Process modeling of 18o isotope separation using stochastic exchanges of isotopic components in a packed distillation column


Main Subjects


Oxygen naturally consists of three isotopes of 16O, 17O, and 18O. The isotopes are widely used in natural and physical sciences such as agricultural, environmental, biochemical researches, medical diagnosis, etc. Some typical methods to enrich the stable isotopes like 18O are chemical exchange processes, diffusion-based separation, laser separation methods and distillation [1-4]. Due to relative difference vapor pressure between the different isotopes of lights and heavy elements such as H216O and H218O, a practical method of isotope separation is distillation [5-9]. The abundance of 18O isotope in natural water is very small (about 0.2%) and more than 2000 equilibrium stages are required to separate it. So, to produce a high enrichment of 18O isotope, the distillation cascade columns are designed. By using this cascade columns, the vapor and liquid products with very high concentrations of (or even pure) more volatile and less volatile components, respectively, can be obtained (Figure 1).

Since process takes a lot of time (more than one year) to achieve steady state condition, the real experiments for process assessments are very difficult, so there are only a few experimental works on 18O isotope separation. The most studies on the 18O production were done by Dostrovsky and Wolf et al. [10-11]. Moreover, a few modeling studies were presented by some researchers. In 1963, Fraenkel and Raviv developed a modeling based on plate by plate calculations with considering three isotopic elements of D (Deuterium), 17O and 18O under steady state conditions [12]. In recent years, other researchers have developed the modeling with separation columns in series by ASPEN simulation software [13-14]. In 2016, Gao et al. studied static and dynamic modeling of 18O isotope separation by considering the effect of deuterium isotope on the 18O isotope separation [15]. In these simulations, three relative isotopes of D, 18O or 17O were considered using the constant values of 1.054, 1.006, and 1.003 for separation factor in the modeling.

The literature available indicated changing of process condition such as pressure or temperature in the distillation process which can change separation factor [16]. Furthermore, according to several kinds of isotopic molecules that exist in water, it should be considered as a multicomponent system. There are 18 isotopic molecules for the water by combination of H, D, T, 16O, 17O and 18O which T (tritium) isotope can be ignored due to its low concentration in the water. Due to this issue, 9 compounds have been studied in this study which their names and physical properties are listed in Table 1 [17-18].

Interactions and collisions between isotopic components in water lead to a complex isotope exchange reaction network [19]. These isotopic exchanges can change the concentration profile of the components along the column by creating a new equilibrium between the components and cause increase or decrease of the separation efficiency.

The objective of this paper is to present a distillation modeling with multicomponent approach of water distillation for 18O isotope separation. New work in this research includes: Estimation of isotopic equilibrium between liquid-vapor phases using statistical methods based on molecular collision, considering the separation coefficient as a function of column temperature (In the models presented by previous researchers, a constant separation factor is considered).

Due to the limitations of constructing a cascade of distillation columns, the simulation results reported in the literature have not been confirmed by experimental data. However in this study, the model was evaluated using experimental data from a pilot column operating under total reflux conditions.

Material and method

The distillation column used in this study had a diameter of 0.15 m and a height of 24 m and was filled with Sulzer-CY packing. The column was operated under total reflux condition and vacuum operating pressure for eight days to reach steady state condition. In this study, the distillation column was charged with 68 kg of the water, containing 12.7% of 18O and 10.2% of D isotopes that was supplied by a semi-industrial plant. The unit was operated for 8 days (192 hr) under total reflux condition and the samples were taken from the condenser and the reboiler of the column every 24 h and were analyzed with a LGR (isotope-ratio analyzer). Since the column was working at total reflux condition and also the fact that both the physical property and latent heat of all the isotope components are close to each other, the values of internal liquid and vapor flows were considered to be the same. The values of these internal flows were obtained from the reboiler and condenser duties so that maximum internal flow in the column was calculated 35 kg/hr, which was below the flooding point. The vacuum pressure of 30kPa was supplied by means of a vapor ejector at the top of the column. The pressure difference was considered linearly for each stage of column and related temperature was estimated from vapor pressure correlation such as Antoine equation. The holdup values of the column packing were determined through measuring the reboiler and condenser inventory and calculating the differences between them with the charged values.

Mathematical equations

The distillation process was modeled by deriving a series of algebraic and ordinary differential equations (ODEs), which were obtained by applying the mass balances between internal flows of the column. The following assumptions were used in dynamic simulation of the process:

  • The liquid and vapor of each stage are in equilibrium.
  • The operation is adiabatic and in this modeling was ignored ignores the heat of mixing and the difference between enthalpy of the two phases.
  • Because of the physical property similarities of the isotopic water component, heat and energy balance was not considered in this modeling.
  • Due to the operation of the column at vacuum pressure, the equilibrium calculations of the phases were performed using the rules of ideal systems.

Basic equations

The model equations were written with considering the above assumptions as follows:

Where L and V are the liquid and vapor flow rate in the column, M is the column inventory, x and y are the isotopes compositions, k is the equilibrium constant, i and j are the number of components and equilibrium stages.

The values of the isotope component concentrations in the vapor phase can be calculated from the concentrations in the liquid phase and separation factor values of each component from an equilibrium relation of a multicomponent system as [20]:

N is the equilibrium stage that was estimated about 150 by Fenski-Underwood method based on steady state concentration of D and 18O at the condenser and the reboiler of column [21]. In 1968, Van Hook presented a correlation for calculating vapor pressures of all the water isotopic components [22]:

where P is the vapor pressure. Through the equation 9 and the values for coefficients A, B, and C, the vapor pressures of each water isotopic components can be calculated. The coefficients are shown in Table 2.

Given that the left side of Equation 6, i.e. the ratio of vapor pressure of each component to the lightest component of water, is equal to the separation coefficient of water isotopes, it can be used to calculate the separation coefficient in modeling. The separation factor means the volatility differences between the components or ratio of their vapor pressures which was written as follows:      

Isotopic exchanges and equilibrium

In the water distillation process, the isotope exchanges of components are inevitable. This interaction is the same with reaction and can create lighter or heavier components at each equilibrium stage. A typical exchange that occurs among three isotopes is presented as follows:

The equilibrium constants of reaction for water isotopes and their kinetics are not clearly specified. So here the statistical procedure for calculating the equilibrium concentration of isotopic reactions were used [23]. Some of these equations are as follows: 

Equations 9-11 are obtained based on stochastic exchanges and probability of molecular collisions in the liquid phase. For example, HD18O is calculated as follows:

P here is the probability of collisions between molecules that forms a new component. In the above equations, XD, X17 and X18O have been defined as follows:

The relative concentrations of D and 18O or 17O measured by LGR were used in order to compare between the model results and the experimental data. The relative concentrations in the modeling were calculated from equations 14-16 and reversely by equation 9-11. The isotopic concentrations for feed were calculated by the same method from isotope components of D and 18O and 17O as shown in Table 3.

The equations from 1 to 8 express that the water distillation process is modelled by classical method as a conventional vapor/liquid equilibrium stage without any reaction or isotopic exchange. The outgoing liquid stream passes to a reactor where isotopic exchange is established (Equation 9-11).

Then, the stream leaving this reactor passes on to the next equilibrium stage. This modeling approach was used by Davis et al. in a distillation column in which formaldehyde reacts with water and methanol (Figure 2) [24].

Numerical solution

The Euler’s method in the Matlab programming was applied to solve the equations on a system of dual core CPU i5 and speed clock of 1.8GH and 6GB RAM. To increase precision, high discretization division applied to the time step (i.e. 0.00005 S). The simulation run time of the proposed model for 150 equilibrium stages was approximately 12 hr.

Results and discussion

Isotopic water components distribution

The concentration distributions of isotopic components in each equilibrium stage for time period of 192 hr are presented in Figures 3, 4, and 5. As can be seen from the figures, during the unsteady period, significant changes in isotopic components are observed which the changes are mainly influenced by separation factor of components. Among the water isotopic components, H-labeled components such as H216O, H217O and H218O are the lightest components. These components have a higher abundance and separation factor than other components. As shown in Figure 3a, the H216O increases at the condenser of the column to above 90% and decreases at the reboiler of the column below 50% respectively. The separation factors of the H217O and H218O come between the light (H216O) and heavy (Deuterium-labeled component) components of water. At the beginning of startup, these components are increased at the condenser of the column and heavy components increased at the reboiler of the column. With the progress of the process, the competition between H217O and H218O with H216O causes decreasing them at the condenser of the column and relatively decreasing at the reboiler of the column. As it can be seen in the Figure 3b-c, H217O and H218O concentrate in the middle of column.

Figure 4a-c shows the HD-components distributions in the column. The HD components such as HD16O, HD17O and HD18O have a low separation factor and tend to concentrate at the reboiler of the column.

Figure 5a-c shows the distributions of deuterium-labeled components in the column. The deuterium-labeled components such as D216O, D217O, and D218O have a lower separation factor than other components and tend to concentrate at the bottom of the column. As observed in these figures, higher concentrations of deuterium-labeled compounds such as HD18O and D218O in water, due to the low separation factor, help to increase the isotope of 18O at the bottom of the column.

The results obtained from the distribution of isotopic compounds show that the separation factor of components, which is a function of temperature, has an important role in concentration changes of isotopic components in  a water distillation column. The compounds with low separation factor are concentrated at the bottom of the column and compounds with high separation factor are concentrated at the top of the column. Also, the interaction between these components can change the concentration values of each component in the column and increase or decrease the separation efficiency. To investigate these issues on the modeling and evaluate the model, the model results have been compared with experimental data obtained from a semi-industrial pilot, which is mentioned in the next section.

Model validation

D and 18O isotopes were analyzed by LGR at the top and bottom of the column as the zones for controlling component concentrations in the distillation column and compared with the results obtained from the concentrations of D and 18O isotopes at this locations in unsteady state periods. According to Figures 6a-d, the concentration of D and 18O isotopes is increased at the reboiler of the column and decreased at the condenser of the column. As can be seen in these figures, there is a well adaptation among the modeling results with the pilot data.


In this paper, modeling of the 18O isotope separation based on multicomponent approach of water distillation process was studied. The effect of temperature was considered in calculating the separation factor and the isotopic reactions of water isotopic components were estimated by statistical method based on stochastic exchanges of isotopic elements. The modeling equations were solved using the numerical procedure in the Matlab programming for the distribution of isotopic components in the distillation column. The simulation results indicated that increasing the deuterium-labeled components in the water have a significant effects on the 18O isotope separation. The results of concentrations for 18O, and deuterium isotopes at the condenser and the reboiler of the column were compared with the corresponding pilot data which obtained from a semi-industrial packed column. It resulted in a well adaptation between the modeling results and the pilot data. The proposed model for this case study can be used for simulation and optimization of the 18O isotope separation in the real systems.



The authors are thankful to all the professors and colleagues in the field of academia and industry who participated in this study.


k                           k-value

L                           Liquid Rate (mol/s)

M                          Hold Up (mole)

16O, 17O, 18O         Oxygen Isotopes

P                           Pressure (kPa)

p                           Probability

T                          Temperature (K)

V                          Vapor Rate (mol/s)

M                          Hold Up (mole)

x                            Liquid Mole fraction

y                           Vapor Mole fraction

L                           Liquid Flow Rate (mol/s)

α                           Separation Factor Subscripts

i                           Component

j                           Stage



H. Mallah:


How to cite this article: R. Lashkarboluki, M. H. Mallah*. Process modeling of 18o isotope separation using stochastic exchanges of isotopic components in a packed distillation column. Eurasian Chemical Communications, 2022, 4(2), 137-146. Link:


Copyright © 2022 by SPC (Sami Publishing Company) + is an open access article distributed under the Creative Commons Attribution License(CC BY)  license  (, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.



[1] A.P. De Groot, Elsevier, Handbook of Stable Isotope Analytical Techniques Vol II Amsterdam, 2008, 1025–1032. [Pdf], [Google Scholar], [Publisher]
[2] S. Rotering, K. Franke, J. Zessin, P. Brust, F. Füchtner, S. Fischer, J. Steinbach, Appl. Radiat. Isot., 2015, 101, 44–52. [crossref], [Google Scholar], [Publisher]
[3] E. Schwartz, M. Hayer, B.A. Hungate, B.J. Koch, T.A. McHugh, W. Mercurio, E.M. Morrissey, K. Soldanova, Curr. Opin. Biotechnol., 2016, 41, 14–18. [crossref], [Google Scholar], [Publisher]
[4] W.A. Van Hook, Isotope separation, Kluwer Academic Publishers, Dordrecht, 2003, 177–211. [crossref], [Publisher]
[5] J. Horita, K. Rozanski, S. Cohen, Isotopes Environ. Health Stud., 2008, 44, 23–49. [crossref], [Google Scholar], [Publisher]   
[6] H. London, Separation of Isotopes, George Newnes Limited, London, 1961.
[7] H.G. Thode, S.R. Smith, F.O. Walking, Canadian Journal of Research, 1944, 22b, 127–136. [crossref], [Google Scholar], [Publisher]  
[8] S. Szapiro, F. Steckel, Trans. Faraday Soc., 1967, 63, 883-894. [crossref], [Google Scholar], [Publisher]  
[9] N. Matsunaga, A. Nagashima, I.J.T., 1987, 8, 681–694. [crossref], [Google Scholar], [Publisher]  
[10] I. Dostrovsky, E.D. Hughes, D.R. Llewellyn, Nature, 1948, 161, 858-859. [crossref], [Google Scholar], [Publisher]  
[11] D. Wolf, H. Cohen, Can. J. Chem. Eng., 1972, 50, 621–627. [crossref], [Google Scholar], [Publisher]  
[12] Z. Fraenkel, A. Raviv, W. Klein, Chem. Eng. Sci., 1963, 18, 697-709. [crossref], [Google Scholar], [Publisher]  
[13] Y.H. Gao, Z.H. Xu, Z.J. Yu, W.Y. Fei, CIESC J., 2014, 65, 213–219. [Google Scholar]
[14] Y.Y. Jiang, Y.Y. Chen, C.J. Qin, Y. Liu, H.S. Gu, J. Isot., 2011, 102–105. [crossref]
[15] Y.H. Gao, Z.H. Xu, K. Wu, Z. Yu, W.Y. Fei, Chinese J. Chem. Eng., 2016, 24, 979–988. [crossref], [Google Scholar], [Publisher]
[16] B.M. Andreev, E.P. Magomedbekov, A.A. Raitman, A.Y. Pozenkevich, Sakharovsky, A.V. Khoroshilov, Elsevier, Amsterdam, 2007. [Google Scholar], [Publisher]
[17] G.S. Kell, J. Phys. Chem. Ref. Data., 1977, 6, 1109-1131. [crossref], [Google Scholar], [Publisher]
[18] D. Palmer, A.R. Fernandez-Prini, A.H. Harvey, steam and hydrothermal solutions, 2004. [Google Scholar], [Publisher]
[19] W.A. Van Hook, J. Phys. Chem., 1968, 72, 1234-1243. [crossref], [Google Scholar], [Publisher]
[20] R. TaylorR. Krishna, Multicomponent Mass Transfer, Wiley, New York, 1993, 395-417. [Pdf], [Google Scholar], [Publisher]
[21] R.W. Rousseau, Handbook of Separation Process Technology, New York: J. Wiley, 1987, 70-90. [Pdf], [Google Scholar], [Publisher]
[22] G. Jancso, W.A. Van Hook, Condensed phase isotope effects, Chemical Reviews, 1974, 689–750. [crossref], [Google Scholar], [Publisher]
[23] R. Hanson, S. Green, Introduction to Molecular Thermodynamics, University Science Books, 2008, 1-17. [Pdf], [Google Scholar], [Publisher]
[24] S. Dilfanian, Ph.D. Thesis, University of Aston in Birmingham, 1978, 56-65. [Google Scholar], [Publisher]