Open access

Continuous Chromatography Modelling with 2D and 3D Networks and Stochastic Methods – Effects of Porous Structure and Solute Population

Written By

Leôncio Diógenes T. Câmara, Jader Lugon Junior, Flávio de Matos Silva, Guilherme Pereira de Oliveira, Lídice Camps Echevarria, Orestes Llanes Santiago and Antônio J. Silva Neto

Submitted: 08 February 2012 Published: 24 July 2013

DOI: 10.5772/51957

From the Edited Volume

Mass Transfer - Advances in Sustainable Energy and Environment Oriented Numerical Modeling

Edited by Hironori Nakajima

Chapter metrics overview

2,377 Chapter Downloads

View Full Metrics

1. Introduction

The modelling of separation chromatographic processes reported in the literature is, in general, related to macroscopic approaches for the phenomenological representation of the mass transfer mechanisms involved. In such models the microscopic aspects of the porous medium structure, related to the separation mechanisms, are incorporated implicitly, limiting the quality of the representation of the separation systems, which are strongly influenced by the micro porous adsorbent.

The modelling of fluid flow in porous media through the application of interconnected networks, which considers the global result from a system of interconnected microscopic elements, is related to the concepts of percolation theory.

The classical macroscopic models of chromatography have limitations in representing the structural parameters of the solid adsorbents, such as topology and morphology, as well as population effects of the molecules in the liquid phase, i.e. multi-molecules movement. Such important microscopic properties can be represented applying interconnected network models which can lead to a better understanding of the phenomenological aspects that contribute to the separation mechanisms in micro-porous media.

In the simulation of the molecules flow through the column porous medium, a stochastic approach is utilized to represent the adsorption, diffusion and convection phenomena.

The molecules can move freely in the network structure, from one neighbor site to another, being respected the requirement that the final position is not occupied by another molecule. Therefore, two molecules cannot occupy the same network site.

This chapter is dedicated to the modelling of continuous chromatography with a network approach combined with Monte Carlo like random walk stochastic methods. The porous structure of the solid adsorbent phase of the chromatographic column is represented by two and three dimensional networks, respectively square and cubic lattices (Oliveira et al., 2008; Biasse et al. 2010), in which population effects are taken into account, being related to the movement of multi-molecules modeled by stochastic phenomena of adsorption, desorption, diffusion and advection.

The use of network models to study chromatographic separation processes has been observed in the literature with different techniques and applications (Loh & Wang, 1995; Kier et al., 2000; Loh & Geng, 2003; Geng & Loh, 2004; Bryntensson, 2002; Oliveira et al., 2008, Biasse et al., 2010). In the work of Kier et al. (2000), a square network model was applied in the representation of the chromatographic column utilizing a cellular automata approach. The authors assumed arbitrarily pre-defined probabilities for the particles motion and interactions among them. Loh & Geng (2003) applied a cubic network model of interconnected cylindrical pores in the study of chromatographic systems of perfusion. Topological and morphological aspects, such as connectivity and pore size distribution, were analyzed, being observed a great influence of such porous adsorbent media characteristics on the mass transfer of the phenomena studied. In Geng & Loh (2004), the porous structure of the column was modeled considering three different Gaussian distributions of pore size, in order to represent the macro-pores, the micro-pores and the interstitial pores.

In this chapter, the fluid advection is assumed to be the main factor contributing to the molecules movement in the network structure. Such assumption is reasonable since the fluid movement in the chromatographic column comes from the external driving force provided by the pumping system. Two pore dimensions are assumed in the porous structure of adsorbent, the small and large cavities, which leads to steric and non-steric effects, respectively, due to the required space for the molecules movement.

The application of interconnected network models combined with stochastic phenomena of adsorption, diffusion and advection represents the main dynamical behaviors of the chromatographic processes of separation. The multi-molecules population effects allow the study of the dynamics of percolation through the chromatographic column, being therefore possible to evaluate the influence of molecules concentration on the mass transport phenomena along the chromatographic column.


2. Phenomenological and structural modelling

In the adsorption phenomenon, the molecule arriving at a new site can be adsorbed according to a probability of adsorption (pads). In this case, the system generates a random number R, from a uniform distribution between 0 and 1, being it compared to a pre-defined adsorption probability (pads). The molecule is then adsorbed at the new site if the random number R is smaller than the adsorption probability (R<pads). The same procedure is applied in the desorption of a previously adsorbed molecule, with a desorption probability (pdes). The relation between the adsorption and desorption probabilities considered in the simulations is

pads+pdes=1,  (0pads1)E1

The adsorption and desorption phenomena are represented by

A+s k1k2 A.sE2

in which a molecule of solute A can be adsorbed at or desorbed from a site s according to a kinetic constant of adsorption (k1) or desorption (k2), respectively. From Eq. (2) it can be observed that the adsorption rate is proportional to the concentration of both solute A and vacant sites s.

The simulation of the dynamic process of chromatographic separation in adsorption columns was performed combining the network modelling of the column porous structure with the stochastic modelling of molecules movement and interactions percolating the system.

In the schematic representation of the adsorption column shown in Fig. 1, the symbols C0 and C represent the solute concentrations at the entrance and at the exit of the column, respectively.

A square network was used to model surface phenomenon (studying equilibrium isotherms) and a cubic one to model the diffusion phenomenon. The porous medium of the chromatography column was represented by a two or three-dimensional cubic network model of interconnected sites (Vide Fig. 2). In such structure the connectivity, i.e. the number of neighbors connected to each site, is equal to four (2D) or six (3D). Each intersection node or site in the network corresponds to a potential adsorption location, in which the solute molecule may be adsorbed, being permitted only one adsorbed molecule per intersection.

In Fig. 2 one can see the graphical representation of both models used with the percolation threshold values (pcs) considered over the site approach. The percolation threshold is the minimum probability of site occupation that represents the percolation through the whole system, becoming a cluster that goes from one extreme of the network to the other.

Figure 1.

Chromatographic column represented by a three-dimensional cubic network.

Figure 2.

Graphical representation of two-dimension square network (A) and cubic network (B) with respective percolation threshold values (pcs)

2.1. Adsorption isotherm model

In this stage, molecules adsorption in stirred tanks was modeled, so that the adsorption surface was represented as a two-dimensional network. In Fig. 3 the process is represented, the liquid phase molecule can be adsorbed at the adsorbent surface material, being the latter represented as the network at the tank bottom.

Figure 3.

Stirred tank represented by a two-dimensional square network.

In this model it is assumed that the liquid phase molecules can be adsorbed at the network according to an adsorption probability (pads), defined as


where Cit and Ci0 represent the liquid phase molecules concentration at adsorption times previously prescribed as t and t=0 respectively. It is known that at time t=0, the concentration is Cit=Ci0. The network used was a 100 x 100 sites, so that it takes to a number of 10.000 adsorbent sites that can be occupied, or not, by the molecules in the system. For each time step (t) a liquid phase molecule can adsorb to the surface according to a previously prescribed adsorption probability (pads). The adsorption takes place if the random number (R), taken by a uniform random number generator, is lower than (pads), that is, R<pads. Another condition is that the random chosen position is not occupied by another molecule. The adsorbed molecules are subtracted from the liquid phase. It is considered a final time equal to 106, so that a infinite adsorption time, or equilibrium, is represented. It is not considered the interaction between adsorbed molecules, i.e. the superficial diffusion processes are not taken into account.

2.2. Diffusion phenomenon modelling

The diffusion phenomenon model was preformed through the "random walks" technique both for the two-dimensional and three-dimensional percolation at the square network and at the cubic network, respectively. Using this procedure one can use diffusion parameters for chromatography porous media and topological properties from network models, and obtain a fundamental universal correlation for these complex phenomena. Two different diffusion phenomena were studied at this stage: the first assumes that a certain molecule may follow any direction at the network; and the second that assumes that a molecule is only able to diffuse in the axial direction. The second case establishes diffusion relations for the axial dispersion for a chromatographic column, which is a important phenomenon that draws a lot of attention in separation chromatographic processes. In general, on the macroscopic models, only dispersion at the axial direction is taken into account, neglecting the transversal dispersion, that makes the numerical solutions to be simpler. In the present modelling, both axial and transversal dispersion are implicit on the stochastic "random walk" model, since the molecule is able to diffuse on both directions.

Equation 4 represents the diffusion mechanism, in which the molecule can diffuse in (n) directions with equal diffusion probability (pi) for each direction. In general, for the square and cubic network, with axial dispersion, (n) corresponds respectively to 4 and 6, but when axial dispersion is considered (n) corresponds to 3 and 5, respectively, being neglected the backward diffusion along the axial direction. In all cases, the moving probability for each direction (i) is considered to be equal.


At first, a probability of occupation for network elements (poc) in each direction i is considered to be equal. Simulations for poc<1 were also done, in such a way to represent the porosity (ε) of the porous adsorbent media. A porosity represents the empty space fraction in the chromatographic column, being proportional to the occupation probability (pocε). In the case that poc=1, it is considered that the column is completely empty, and therefore (ε=1). Networks with poc<1 were built randomly, being the existence of a certain site conditioned to the generation of a random number R, and the satisfaction of R<poc condition.

Simulations for different time ranges were made (10, 25, 50 and 75), being monitored the distance (dp) by the "random walk" and also the area (Ap), the latter counted in site units. In the case of 3D, the latter has a relationship with the volume (Vd), also counted in site units. For each time-step, a different molecule direction was randomly chosen, in a such way that the destination site was forced to exist; A random number of events N=106 was performed in order to obtain a distance (dp) by the "random walk" with smaller dispersion.

In the case of smaller molecular diffusion (Dm), it was observed a certain proportion Dmdp with the distance, obeying to the power law (Stauffer and Aharony, 1992, Biasse et al., 2010).


where c and k are constants and t represents time.

2.3. Movement rules and steric effect

In the simulation of the molecules movements, four rules (MR – movement rules) are considered to be representative of the diffusion and advection mechanisms. Such rules are schematically represented in Fig. 4.

From Fig. 4 it can be observed that the movement rules are determined by the directions considered in each situation. In the movement rules (I) and (II), the molecules can move in all directions of the structure while in the movement rules (III) and (IV) the molecules cannot move in the direction 6, which corresponds to the movement against the longitudinal main stream flow direction, i.e. the upstream movement is not allowed. The MR-(II) and (IV) are similar to MR-(I) and (III), respectively, with a higher chance of movement in the direction 5 (longitudinal or axial direction).

The MR-(I) is considered to be representative of the diffusion mechanism of solute molecules through the porous structure.

This assumption is coherent as there is not a driving mechanism forcing the flow in any particular direction of the network. In this case, according to the first law of Fick (Bird et al., 2002), the solute flow is determined by the concentration gradient of molecules, without significant effects of external forces. The diffusion mechanism of MR-(I) is governed by

pi=16,  i=1,2,,6E6

in which Pi represents the probability of the solute molecule to move in the i-th direction of the network. According to Eq. (6), the chance of the molecules to follow any direction (six directions) in the network is the same.

The MR-(II), (III) and (IV) are considered to be representative of the advective mechanism as these configurations favor the flow in the axial downstream direction of the column.

The advective mechanism MR-(II) is governed by

pi=1-p55, i=1,2,,6ii5, p5>16 E7

in which p5 represents the probability of the solute molecule to move in the axial downward direction of the column. In this case, the chance of the solute molecule to move in the axial downward direction is greater than the other directions, and the probabilities related to the movements in the remaining directions are the same.

Figure 4.

Movement rules (MR) for the solute molecules with the corresponding directions.

The advective mechanisms MR-(III) and (IV) are determined, respectively, by

pi=15,  i=1,2,,5E8
pi=1-p54, i=1,2,,4,   p5>15E9

From Eq. (8) it can be observed that the solute movement in the direction 6 is prohibited, having the same probability for the moves in the other directions. In Eq. (9) as in Eq. (7), the probability of moving in the direction 5 is greater than in the other directions.

The simulation of the flow of solute molecules was carried out according to the procedure described next. The molecules were introduced randomly at the column entrance (at the nodes with k = 1), maintaining the concentration constant at this section (C0). The molecules concentration (C) at the column exit (k = 30) was calculated as the ratio of the number of molecules occupying the network intersections (sites) and the total number of intersections. In each step of the simulation, corresponding to a discrete time (t*), all the network sites occupied by molecules were checked, being assumed for each molecule a direction of movement according to the MR to be followed. A molecule was kept at its original location if the new site, it was supposed to be move to, was occupied by another molecule or was located outside the lateral limits of the network domain. The molecules at the entrance could move only in the axial downward direction, being kept for every discrete time t*, a constant value for the concentration at this section (C0, k = 1). At the column exit (k = 30) the concentration (C) was monitored as a function of the discrete time (t*).

Figure 5.

Representation of molecule movement along the column: (A) the non-steric model and (B) with steric model. Obs.: The arrows represent the allowed movement that the molecule can do, and crossed arrows represent movement that can´t be done.

An important parameter to be taken into account in the stochastic modelling is the number of simulations (N), which indicates the number of times that the same procedure of the simulation is performed using the same control parameters. The increase in the number of simulations (N) leads to a decrease in the dispersion of the calculated value.

The steric effects were also investigated, that is, while the solid phase is able to adsorb one molecule at each site, one possibility is that the liquid phase site is able to contain one molecule and another possibility is that it is able to contain a unlimited number of molecules. In Fig. 5 are represented those two possibilities, being (A) for the case without the steric limitation and (B) representing the steric restriction.


3. Results and discussions

3.1. Langmuir isotherm model for surface adsorption

In Fig. 6 are represented the adsorption isotherms: in (A) are the results obtained for the simulations using a square network and in (B) are the experimental and deterministic model adapted from Silva (2004). One can observe that the results obtained using the stochastic adsorption method are representative for the studied phenomenon with the classic Langmuir isotherm.

These results show that the stochastic phenomenology is determinant to the behavior of equilibrium systems with multimolecules, and the overall result is governed by the individual actions of each component.

Figure 6.

Results obtained with the stochastic simulation (A) and the experimental and classical modelling (B) by Silva (2004).

3.2. Results for the diffusion phenomenon (2D and 3D modelling)

In Figs. 7, 8 and 9 are represented the percolation evolution both in the square and cubic networks using different moving mechanisms. In Fig. 7 are represented the percolation using a square network of 50x50 nodes, for a time equal to 2000 steps, for four directions, and occupation probability poc=1 (Fig. 7 A) andpoc=0.7 (Fig. 7 B).

In Fig. 8 are represented results for the same conditions considered before, but allowing only 3 directions for the dispersion mechanism. One can observe a greater tendency towards the axial diffusion.

A evolution percolation method for 3D is presented in Fig. 9 for two movement mechanisms.

Figure 7.

Diffusion through the square network for four directions with poc=1 (A) and poc=0.7 (B).

Figure 8.

Diffusion through the square network for three directions with poc=1 (A) and poc=0.7 (B).

Figure 9.

Representation for cubic percolation considering six directions (A) and five directions (B).

In Figs. 10 and 11 are presented the distance (dp) with time in log scale. Using this scale, the power law (Biasse et al., 2010) represented by Eq. (5) becomes


We can obtain the expoent k, that corresponds to the angular coefficient of the expression above. The k exponent is important because of its relation to the molecular diffusivity and time, being this a fundamental parameter for this relation and sensitive to those values.

From the analysis of Figs. 10 and 11 we are able to observe that diffusivity, which presents a relation with dp, has a direct link to the power law established by Eq. (5), once the model leads to a good correlation with the simulated data for the percolation both for the square (2D) and cubic (3D) networks. The results presented high correlation coefficients, which indicates a good agreement to the power law. Both situations for the regular network with poc=1, for the 2D (Fig. 10A) and 3D (Fig. 11A) networks, presented ideal percolation with universal exponents equal to 0.5014 and 0.4987, respectively. In those cases, because both are regular, with poc=1 and movement in all directions, the correlation coefficient was very closed to 1 (R2=1).

Figure 10.

Results for 2D networks for four (A) and three (B) moving directions

It is observed, in axial and radial dispersion situations, represented by Figs. 10B and 11B, a significant increase in the value of k exponent when compared to the conditions with diffusion in all directions, demonstrating that diffusion was higher in those cases. For example, we can point out the difference between Figs. 10A and 10B (with poc=1), where the k exponent was increased by 79,7%. With respect to the porosity ε effect, that has a relation to the poc was observed that the reduction of poc from 1 to 0.7 had a more significant effect over the k exponent on 2D networks than on 3D ones. This effect can be observed by the greater slope in Fig. 10. Therefore the porosity reduction leads to a more significant slope on the diffusion for the 2D network represented by lower k exponents. The 3D network, has a lower percolation threshold when compared to the 2D network, presenting a structure with a larger number of possible percolation ways, that with the possibility for the "random walk" to go over larger distances, leading therefore to the larger values of k exponents when compared to 2D networks with the same value of poc<1.

Figure 11.

Results for 3D network for six (A) and five (B) movement directions.

The analysis of the power law Eq. 5 (Biasse et al., 2010), not only related to the distance (dp) but mainly in terms of (Ap) and (Vp), took to equivalent results to those previously reported, showing a good agreement of the experimental data with the power law model assumed. Table 1 presents the values of k exponents and the corresponding correlation coefficients. The values for Ap and Vp were obtained in percolated site units, being these units respectively related to 2D and 3D. Such as in the previous cases reported, related to dp, it was observed a reduction of exponent k for the axial and radial dispersion (three and five directions) situations. It was also observed a more significant reduction of k with the porosity reduction in 2D networks.

Table 1.

Exponents (k) from power law related to the area Ap (2D) and volume Vd (3D).

3.3. Results of 3D network column with advective phenomenon

All results presented in this section were obtained for 5000 simulations. In Figs. 12 and 13 the test results varying pads are presented. One can observe that for pads=0 there is a concentration step, because since no molecule is retained at each time, all of them are allowed to move to the next column section until they reach the exit section. Another observation is that when pads is increased there is a delay in the variation of the concentration at the column exit, because the chances for the molecule to be adsorbed are higher.

Figure 12.

Chromatography column exit for pads=0 and pads=1.

Figure 13.

Chromatography column exit for different values of pads.

In Fig. 14 it is presented a comparison of the simulation results, without considering axial dispersion, with the experimental data acquired by Cruz (1997). One can observe the good results obtained by the network simulations, despite the breakthrough curve too sharp when equilibrium is reached (CC0=1), in which the curve does not fit well, probably because of the effects of axial dispersion.

Figure 14.

Dimensionless concentration at the column exit as a function of time. A comparison between the experimental data of Cruz (1997) and network simulations without axial dispersion

Another important parameter to be considered here is the equilibrium fraction (keq), i.e. the fraction between adsorption and desorption. Another way to find the value for that ratio is comparing the amount of particles adsorbed and the amount of particles in the liquid phase (qece). Using this model one can confirm that there is a direct relationship between the concentration fraction and the padspdes. See Table 2.

Table 2.

Relationship between keq and padspdes

In Fig. 15 the same network, as the one used to obtain the results shown in Fig. 14 and Table 2, was considered assuming also the axial dispersion phenomenon.

Finally, in Fig. 16 are presented the simulation results (3D network model) for the discrete relative concentration (C/C0) as a function of the adsorption probability ( CC0) taking into account the steric and non-steric effects. The relative concentration indicates a measure of the porous medium resistance to the molecules percolating through the system.

The resistance to the percolation through the chromatographic column decreases as the relative concentration (pads) increases and vice-versa. Therefore, values close to one show a structure without any resistance to the molecules to percolate the column structure.

Figure 15.

Dimensionless Concentration at the column exit as a function of time. Comparison of experimental data acquisition by Cruz (1997) and network simulation considering axial dispersion.

The highest medium resistance corresponds to an adsorption probability near 0.68 which is the minimum point of the curve. At this point the number of molecules percolating the system is reduced to the lowest level.

Figure 16.

Relative concentration as a function of the adsorption probability considering steric and non-steric effects.

It must be noted that the medium resistance disappears at adsorption probability equal to 1 as at this point the adsorption is irreversible (no desorption is observed), and no molecule already adsorbed goes back to the liquid phase of the cavity, and therefore there is no interference in the movement of the other molecules.


4. Conclusions

The network stochastic surface multimolecular modelling was able to represent the behavior of Langmuir type adsorption isotherms. Such computational tool can be used to better understand the adsorption mechanisms to surfaces, resulting in the improvement of adsorbent materials.

The diffusion modelling through the 2D and 3D network "random walks" presented results obeying the power law with universal exponents well defined, being the latter related to the diffusive coefficients.

The effects of the power law were observed through the distance  CC0, the area dp and the volume Ap, being these values related to the diffusion phenomenon. Axial and radial dispersion mechanisms were also represented with a power law modelling, being that behavior related the variations of porosity Vd and probability of occupation ε.

A final important result is the observation of the existence of a direct relationship between the adsorption and desorption probabilities ratio with the equilibrium constant poc.



The authors acknowledge the financial support provided by CNPq, Conselho Nacional de Desenvolvimento Científico e Tecnológico, FAPERJ, Fundação Carlos Chagas Filho de Amparo à Pesquisa do Estado do Rio de Janeiro, and CAPES, Coordenação de Aperfeiçoamento de Pessoal de Nível Superior.


  1. 1. BiasseA. D.FonsecaE. X.SilvaNeto. A. J.CâmaraL. D. T.2010Power-Laws of Diffusion Phenomena In Chromatographic Columns, Advances and Applications in Mathematical Sciences 22313320
  2. 2. BirdR. B.StewartW. E.LightfootE. N.2002Transport Phenomena, second ed., Wiley.
  3. 3. BryntessonL. M.2002Pore Network Modelling of the Behaviour of a Solute in Chromatography Media: Transient and Steady-State Diffusion Properties,J. Chromatogr. A, 945, 103115
  4. 4. CruzJ. M.1997Insulin Adsorption in Ion Exchanger Resins using Fix and Fluidzed Beds, M.Sc. Dissertation, Chemical Engineering Faculty, UNICAMP. (In Portuguese).
  5. 5. GengA.LohK. C.2004Effects of Adsorption Kinetics and Surface Heterogeneity on Band Spreading in Perfusion Chromatography- A Network Model Analysis, Chem. Eng. Sci. 59, 24472456
  6. 6. KierL. B.ChengC. K.KarnesH. T.2000A Cellular Automata Model of Chromatography,Biomed. Chromatogr. 14, 530534
  7. 7. LohK. C.GengA.2003Hydrodynamic Dispersion in Perfusion Chromatography- A Network Model Analysis, Chem. Eng. Sci. 58, 34393451
  8. 8. LohK. C.WangD. I. C.1995Characterization of Pore Size Distribution of Packing Materials Used in Perfusion Chromatography Using Network Model,J. Chromatogr. A, 718, 239255
  9. 9. OliveiraG. P.CâmaraL. D. T.SilvaNeto. A. J.2008Modelling of the Separation Mechanisms in Chromatographic Columns through Adsorption and Advection in Two Dimensional Networks, XI Computational Modelling Meeting, Volta Redonda, Brazil. (In Portuguese)
  10. 10. Silva,2000Study of Insulin Adsorption in Ion Exclange Resin Columns: Experimental Parameters and Modelling, D.Sc. Thesis, Chemical Engineering Faculty, UNICAMP. (in Portuguese)
  11. 11. StaufferD.AharonyA.1992Introduction to Percolation TheorySecond Edition.

Written By

Leôncio Diógenes T. Câmara, Jader Lugon Junior, Flávio de Matos Silva, Guilherme Pereira de Oliveira, Lídice Camps Echevarria, Orestes Llanes Santiago and Antônio J. Silva Neto

Submitted: 08 February 2012 Published: 24 July 2013