An in-silico study on transient enzyme diffusion and adsorption within lignocellulosic biomass using a multi-scale model

1 Introduction

Hydrolysis of lignocellulosic (woody) biomass to yield monomer sugars is an essential step in the production of low-emission biofuels, such as cellulosic ethanol, hydrocarbon fuels from genetically engineered microorganisms, or high-value bioplastics and other bioproducts (Tyagi and Anand, 2024). Lignocellulosic or woody feedstocks include energy crops such as hybrid poplar, eucalyptus, willow, switchgrass, and residues such as wheat straw, rice straw, logging residues, and corn stover. The conversion pathway shown in Figure 1 from lignocellulosic biomass feedstock to biofuels involves several significant steps: deconstruction and fractionation, followed by synthesis and upgrading (Nanda et al., 2014). The first step is the deconstruction of biomass into its component sugars using catalysts and heat (pretreatment). The pretreated biomass is then hydrolyzed to produce fermentable sugars in a process called enzymatic hydrolysis (De et al., 2025; Zhang et al., 2021) using adsorbing enzymes (that work at the pore fluid-cellulose interface such as cellobiohydrolase and endoglucanase) and non-adsorbing enzymes (that work in the pore and bulk solution) such as β-glucosidase (Kumar and Murthy, 2017). After enzymatic hydrolysis, the produced sugars are fermented using microorganisms to produce biofuel. The conversion of lignocellulosic materials using these processes is expensive, and to commercialize the production of biofuels from lignocellulosic biomass, the efficiency of the process must be improved (Mustafa et al., 2024). Critical to enhancing enzymatic hydrolysis efficiency is understanding the effects of biomass particle morphology on enzyme accessibility and intraparticle diffusional limitations. The spatial stochastic modeling study on enzyme saccharification of biomass revealed sugar release is the product of crystallinity fraction of biomass particles and enzyme kinetics, inhibition, and steric, which lead to further studies on diffusion limitation of enzymes (De et al., 2025).

Diffusional limitations in biomass conversions have been studied experimentally and theoretically for many bioconversion processes, such as enzymatic hydrolysis, dilute acid hydrolysis, lignin methanolysis, and others (De et al., 2024). Sensitivity analysis of multiple parameters on enzymatic saccharification of pre-treated lignocellulose by De et al., demonstrated Cellobiohydrolase activity, cellulose crystallinity, and crystalline bond digestibility are the major factors influencing enzymatic saccharification (De et al., 2024). Kim and Lee studied diffusion of sulfuric acid into cane bagasse, corn stover, rice straw, and yellow poplar to measure effective diffusivity (Kim and Lee, 2002). They developed a transient diffusion model assuming a flat-plate geometry (of half thickness L) to predict critical diffusion time L ²/D _*e*, with effective diffusivities (D _*e*) extrapolated to elevated temperatures of dilute acid hydrolysis. A criterion for determining diffusional limitations on dilute acid hydrolysis was then derived by setting the Thiele Modules = 1, and solving for a critical particle size, L< (D _*e*/k)^½, where k is a first order intrinsic hemicellulose hydrolysis rate constant. Chen et al. developed an acid diffusion biphasic hydrolysis kinetic model for dilute acid pretreatment of corn stover assuming a cylindrical biomass particle, accounted for both fast and slow hydrolysis fractions of hemicellulose, an evolving porosity, and coupled the acid particle diffusion model to bulk liquid mass balance equations (Chen et al., 2015). Models were developed for production rates of xylose and furfural and compared them to time-dependent data to extract model parameter values. Separate model runs for particle sizes of 0.42 mm, 2 mm, and 4 mm with the same biomass loading confirmed acid diffusional limitations on rates of xylose production, which were most pronounced early during the fast hydrolysis phase, but during the hydrolysis of the slow fraction of hemicellulose rates of production were equal and not affected by diffusional limitations. The study concluded that for pretreatment processes required less than 3 min of reaction time, the bioconversion will be limited by acid diffusion for particle sizes >0.4 mm.

In a study of methanol-based lignin extraction of poplar, Thornburg et al. developed a mesoscale reaction and diffusion model to determine intrinsic reaction kinetics parameters from batch experiments (Thornburg et al., 2020). The model accounted for both forward and reverse reactions for lignin extraction/deposition and diffusion within cell wall and axially in lumen pores. A Thiele Modules analysis concluded that reaction and diffusion processes operate on similar time scales for lignin methanolysis, and that biomass particles of greater than 2 mm in length experience mass transfer limitations. Luterbacher et al. developed a reaction-diffusion model for enzymatic hydrolysis of lignocellulosic biomass and showed that enzyme mass transfer limits glucose yields for particle diameters above 0.5 mm due to enzyme diffusion into biomass particles that occurs prior to hydrolysis (Luterbacher et al., 2013). They also showed that diffusional limitations reduce glucose yield for even the smallest particle sizes during early stages of the hydrolysis process. Zhang et al. studied ultrafine grinding pretreatment of corn stover, achieving size reduction from 218 to 17 μm, which increased specific surface area, pore volume, and enzyme adsorption capacity and affinity constant (Zhang et al., 2016). Behle and Raguin Behle and Raguin developed an in silico stochastic model to predict the enzymatic saccharification of corn stover that incorporated the synergistic effect of enzyme cocktail and effect of corn stover crystallinity and composition (Behle and Raguin, 2021). Rohrbach and Luterbacher developed a multi-scale enzyme diffusion, adsorption, and hydrolysis kinetics model and predicted glucose yields in the presence of an evolving particle porosity (Rohrbach and Luterbacher, 2021). Using this multi-scale model, they studied in silico the effect of particle size, enzyme loading, and biomass loading. One of the key model predictions showed that ultimate glucose yields for biomass particle sizes below 500 μ m was not diffusion limited. This prediction was for the stoichiometric enzyme to biomass loadings and at 30% biomass loading in solution. However, further predictions showed that the particle size threshold for diffusion limitations decreased at lower biomass loadings with stoichiometric enzyme to biomass loading due to depletion of enzymes in the bulk solution.

While most studies involved stochastic modeling, Kumar et al. used Monte Carlo simulations to model the hydrolysis events to account for enzyme crowding effect for various enzymes which captures more realistic interactions between enzymes and cellulose during hydrolysis but does not account for adsorption kinetics of the enyzmes (Kumar and Murthy, 2017). In addition to the studies discussed above, several studies investigated the kinetics of enzyme adsorption within different biomass types and particle sizes. Table 1 presents a summary of these adsorption kinetic studies. In these studies, the biomass shape was reported as cylindrical, and particle radius was expressed as smaller than a particular screen size rather than one specific size. Based on the cellulase adsorption data reported in these studies, we approximated the characteristic time required to achieve equilibrium adsorption (by observation of the time course of the adsorption data), which occurs after diffusion of enzyme into the particle’s ceases, and the enzyme adsorption capacity was achieved.

TABLE 1. Effects of biomass and adsorption parameters on diffusion characteristic time.

(A) - PWS (Pretreated Wheat Straw), Qi et al. (2011); (B) - Machado et al. (2015).

The dataset shown in Table 1 shows a strong dependence of biomass particle size on the characteristic time (CT) to achieve equilibrium adsorption. For the largest particle size in the datasets of 0.025 cm in radius (500 μ m in diameter), the characteristic time to establish adsorption equilibrium was approximately 3–4 h. This time period is short compared to the expected hydrolysis time typical for cellulose hydrolysis of 24–72 h, and therefore a diffusion characteristic time of up to 6 h is not expected to significantly affect the observable kinetics of glucose production. This is consistent with observations and modeling from other studies in which diffusion limitations on hydrolysis kinetics are not expected for particle diameters less than approximately 500 μ m (Rohrbach and Luterbacher, 2021).

The numerous studies cited above have demonstrated the importance of particle diffusion characteristic time in enzyme adsorption and bioconversion processes. However, there is a lack of studies with a focus on understanding the effects of transport parameters on characteristic time of diffusion in biomass enzymatic conversion processes. Conversely, Nill et al. (Nill et al., 2018) conducted a modeling study of system parameters on the kinetics of biomass enzymatic hydrolysis, separate from diffusional effects, revealing insights within the enzyme kinetics domain. Here we report on a similar study with a focus on diffusion and adsorption phenomena only, while omitting enzyme hydrolysis kinetics. Therefore, the evolution of particle parameters, such as porosity due to enzyme reaction, is not considered in this study. The goal is to understand which parameters (particle size, porosity, adsorption strength, and affinity) most strongly affect diffusional characteristic time, individually and in combination. To our knowledge, this kind of in silico study has not been conducted before for biomass enzymatic hydrolysis.

2 Model development

2.1 Particle level model

The purpose of developing these particle-level equations is to determine the relative importance of various diffusion and adsorption parameters on the characteristic time for the diffusion of enzymes into non-shrinking porous cylindrical biomass particles. We also employ simplified boundary conditions by assuming that the enzyme concentration at the particle surface remains constant for this purpose.

2.1.1 Non-adsorbing enzyme: β-glucosidase (BG)

The diffusion of enzymes is influenced by biomass physical properties such as particle size and porosity (Zhang et al., 2021). The porosity ε of a biomass particle can be calculated using the equation adapted from (Zhang et al., 2016):where V P is the pore volume per unit mass of particle and ρ s is the solid density of the particle (density of the wood itself, not accounting for the voids), calculated using an equation adapted from (Zhang et al., 2016):where M c , M h , M l a n d M o are the mass fractions of cellulose (c), hemicellulose (h), lignin (L), and other components (o) on a dry weight basis, respectively. ρ c , ρ h , ρ l , a n d ρ o stand for the solid densities of cellulose, hemicellulose, lignin, and other components on a dry weight basis, respectively. The molecular diffusivity of an enzyme in water (D _*w*) is a fundamental transport property and is affected by temperature, size of the enzyme, and the properties of water. The effective diffusion coefficient of the enzyme within the biomass particle (D _*p*) is lower than D _*w* because of the tortuous path, dead-end path, and non-uniform pore widths that exist in biomass particles. D _*p* is calculated using Equation 3, a correlation derived from data by Whitaker (Whitaker, 2013) for diffusion of liquid solutes in an isotropic porous media:where D w is the diffusion coefficient of enzymes in pure water with no porous media present, and ϵ is the biomass internal porosity. The use of this correlation assumes isotropic diffusion in the biomass particle. Because there is no bulk flow of fluid through the small pores in the biomass, the mass transfer in this system occurs through diffusion only. In a porous media in which cellulose-degrading enzymes of typical size of less than 5 nm are only allowed into pores of sufficient size, there will be a fraction of the pore volume, φ, that is accessible to the enzymes, and therefore in which, diffusive transport can occur. Assuming that the volume average concentration of enzymes (E) within the biomass particle is based on only the accessible pore volume within the particle, then the radial mass diffusive flux of enzymes in an assumed cylindrical-shaped biomass particle is defined as

ϵ = V P V P + 1 ρ s (1)

ρ s = 1 M c ρ c + M h ρ h + M l ρ l + M o ρ 0 (2)

D p = D w ϵ 2 − ϵ (3)

j E = − φ D p ∂ E ∂ r (4)

Thus, in the special case of φ = 0 (no accessible pores to the enzyme), the diffusive flux is zero, and no enzyme can enter the biomass particles from the bulk solution during enzymatic hydrolysis of biomass. It seems correct to base E on the accessible pore volume in the particle because enzyme adsorption and hydrolysis kinetics depend on conditions in the accessible pores. Also, this will assure that gradients of E with respect to r will be the same regardless of the value of φ. In addition, the enzyme can only accumulate in accessible pores in the particle, and therefore the accumulation term in the governing equation of diffusive transport must also include φ. The governing diffusion equation for a non-absorbing enzyme, such as β-glucosidase (BG), in which ϵ and φ are dependent on time, and radial position, is given by Equation 4, an adaptation from Whitaker (Whitaker, 2013):where BG is the intrinsic volume average concentration of non-adsorbing enzymes in the accessible pore fluid of the particle at radial position r, φ is the ratio of pore volume accessible to the enzymes V P A to the total pore volume V P, r is a radial position inside the particle, and D _*p* is the pore effective diffusion coefficient (Equation 3). In the case of constant ϵ and φ, the case we are interested in here, the diffusion equation simplifies to

∂ φ ϵ B G ∂ t = 1 r ∂ ∂ r r φ D p ∂ B G ∂ r (5)

φ ϵ ∂ B G ∂ t = φ 1 r ∂ ∂ r r D p ∂ B G ∂ r (6)

And with further simplification to

ϵ ∂ B G ∂ t = 1 r ∂ ∂ r r D p ∂ B G ∂ r (7)

The initial condition for the model is given by:

At t = 0 ∀ 0 ≤ r ≤ R,

B G = 0 (8)

The enzymatic hydrolysis process and the assumed geometry of the biomass particle are shown in Figure 2. Boundary conditions for the model are given by Equations 9, 10:

For t > 0, At r = 0,

∂ B G ∂ r = 0 (9)

At the particle surface, r = R,where C B G 0 is the bulk fluid concentration of BG, which is assumed to be a constant for this derivation.

B G = C B G 0 (10)

2.1.2 Adsorbing enzymes: cellobiohydrolase and endoglucanase (CE)

The solid density does not account for the pore spaces; hence, to further incorporate the internal pore volume, we use apparent density, which is estimated using Equation 11 adapted from Zhang et al. (Zhang et al., 2016).where V p is the total pore volume, ρ s is the solid density. In this study, the Langmuir adsorption isotherm (Zhang et al., 2016), is used to model the adsorption of enzymes. The Langmuir adsorption isotherm is given by:where q is the equilibrium concentration of solid-bound enzyme inside the particle pores (mg enzyme/g solid), q m is the maximum solid-bound capacity when all pores are accessible to the enzyme (mg enzyme/g solid), K a is the affinity constant (cm³ pore fluid/g enzyme), and C E is the enzyme intrinsic volume average mass concentration inside the accessible pores of the particle (mg enzyme/cm³ pore fluid). The mass balance equation of the adsorbing enzyme with internal pore diffusion inside the cylindrical particles is given by:where D P C E is the effective diffusion coefficient of adsorbing enzyme. In this study, we assume that pore diffusion only, and that surface diffusion has no significant effect on the characteristic time of diffusion of hydrolytic adsorbing enzymes. Furthermore, it is assumed, as before, that ϵ and φ are constant, for which the mass diffusion equation for adsorbing enzymes is:

ρ a = ρ s 1 + V P ρ s (11)

q = φ q m K a C E 1 + K a C E (12)

∂ ϵ φ C E ∂ t + ρ a ∂ q ∂ t = 1 r ∂ ∂ r r φ D P C E ∂ C E ∂ r (13)

ϵ φ ∂ C E ∂ t + ρ a ∂ q ∂ t = φ 1 r ∂ ∂ r r D P C E ∂ C E ∂ r (14)

The model can be further simplified by differentiating Equation 12 with respect to C E, which yields:where, f ′ C E = φ q m K a 1 + K a C E 2 is the differential of the Langmuir isotherm with respect to the concentration of enzyme in the pore fluid (CE). In this model, instantaneous adsorption equilibrium is assumed for modeling the effects of adsorption on the diffusion characteristic time. With d q / d t = d q / d C E d C E / d t = f ′ C E d C E / d t, substituting Equation 15 into Equation 14 yields:

d q = f ′ C E d C E (15)

ϵ + ρ a q m C E K a 1 + K a C E 2 ∂ C E ∂ t = 1 r ∂ ∂ r r D P C E ∂ C E ∂ r (16)

The initial condition for the model is given by:

At t = 0 ∀ 0 ≤ r ≤ R,

C E = 0 (17)

Boundary conditions for the model are given by Equations 18, 19: t > 0 ,

At r = 0,

∂ C E ∂ r = 0 (18)

At r = R,where C C E 0 is the bulk fluid concentration of adsorbing enzyme, which is assumed constant.

C E = C C E 0 (19)

2.2 Multi-scale model for diffusion and adsorption

The multi-scale model couples the diffusion/adsorption mass transfer within the biomass particles with the changes in concentration of enzymes in the bulk solution containing the biomass particles. This model includes two spatial scales; the scale of the biomass particle in addition to the scale of the bulk solution, which spans many 10s or 100s the length scale of the biomass particle. The changing bulk concentration affects the particle surface boundary condition, which in turn affects the rate of diffusion into the particles. Thus, these are coupled transport phenomena. The purpose of developing these coupled model equations is to predict the effects of changing the initial bulk concentration of enzyme on characteristic time for diffusion into the particle, as well as to apply the multi-scale model to map the parameter space for which diffusive mass transport is expected to begin limiting observable enzymatic hydrolysis kinetics. The particle diffusion equations are the same as presented in Equations 7, 16 for non-absorbing and adsorbing enzymes, respectively.

The equation for the flux boundary condition for both non-adsorbing and adsorbing enzymes at r = R is identical, thus only the equation for non-adsorbing enzymes is given next:where k represents the external-film mass transfer coefficient, D p is the effective diffusion coefficient of the non-adsorbing enzyme inside the biomass particle, C B G is the concentration of non-adsorbing enzyme in the bulk solution. The change in concentration of enzymes in the bulk solution takes the same form for both non-adsorbing and adsorbing enzymes; therefore, we only show the equation for the non-adsorbing enzyme, B G, given by the mass balance equation:where V is the initial volume of the fluid inside the bulk solution and A is the total outer surface area of the biomass particles inside the bulk solution, given by:where D is the diameter of the cylindrical biomass particle and L is the length of the biomass particle, n p a r t i c l e is the total number of particles inside the flask given by:where, W b i o m a s s is the total dry weight of biomass particles inside the bulk solution, and W p a r t i c l e is the weight of each particle inside the flask given by:where, V p a r t i c l e is the volume of the particle and ρ b u l k is the bulk density of the biomass particles given by: and where, ρ s is the solid density of the biomass particle. It is noted that A is proportional to the biomass loading in the bulk solution.

φ D p ∂ B G ∂ r r = R = − k C B G − B G r = R (20)

V d C B G d t = − k A C B G − B G r = R (21)

A = n p a r t i c l e π D L (22)

n p a r t i c l e = W b i o m a s s W p a r t i c l e (23)

W p a r t i c l e = V p a r t i c l e × ρ b u l k (24)

V p a r t i c l e = π R 2 L (25)

ρ b u l k = 1 − ε ρ s (26)

3 Method of solution

The particle-level model for adsorbing and non-adsorbing enzymes is numerically solved using Equations 1–19 and Avicel data from Table 2. The characteristic time (i.e., the time needed for the enzyme concentration inside the particle at r = 0 to reach 99% of the enzyme concentration in the bulk solution) is calculated using a custom-written program in MATLAB (MathWorks, Natick, MA, United States). The numerical solution was successfully validated in two ways: (i) by comparing the characteristic times predicted using the model with characteristic times obtained from the analytical expression (Equation 27) by Crank (Crank, 1979) at the limiting case of no adsorption as shown in Figure 3. The numerical model is validated by comparing the concentration of enzymes at different points of time and different radial positions predicted by the particle level model with analytical solution provided by Crank (Crank, 1979) for a cylinder with constant surface concentration.

C = 2 a 2 ∑ n = 1 ∞ exp − D α n 2 t J 0 r α n J 1 2 a α n ∫ 0 a r f r J 0 r α n d r (27)

TABLE 2. Material and parameter data used in the model.

A - Machado et al. (2015), B - assumption, C - All particles are of the same size, D - calculated using Equation 11, E − calculated using Equation 2, F - Zhang et al. (Zhang et al., 2016).

G - Tsai et al. (Tsai et al., 2014).

H - calculated using Equation 3, I - calculated using Equation 1, J - Tantasucharit (Tantasucharit, 1995).

K - Yohana et al. (Chaerunisaa et al., 2019).

L - D w ranges from 10^–6 to 10^–7 cm²/s and assumption, M - Ankathi et al. (Ankathi et al., 2019).

N - Lowest particle size from a particle size distribution (28 mesh, Tyler), O - Chen et al. (Chen, 2015).

P - calculated using Equation 22, Q - assuming aspect ratio L = 4 D, R - calculated using Equation 23, S - calculated using Equation 24, T - calculated by Equation 25, U - calculated using Equation 26, V - (BG enzyme loading: 13.5 μL for 0.15 g of dry biomass in each flask), (CE enzyme loading: 37.5 μL for 0.15 g of dry biomass in each flask), the protein concentration of Accelerase 1,500 and Accelerase BG enzyme is 82 mg/mL.

The effect of porosity, pore size, maximum adsorption capacity, and affinity constant on the characteristic time for enzyme diffusion into a cylindrical Avicel particle is determined here by varying the porosity from 0.1 to 1 with intervals of 0.1, varying the particle radius from 0.0025 cm to 0.025 cm with intervals of 0.0025 cm, varying maximum adsorption capacity from 0 to 100 mg/g substrate with intervals of 20 mg/g substrate and varying the affinity constant from 0 to 20 mL/mg with intervals of 4 mL/mg, respectively. Each parameter is varied while holding other parameters constant, in order to isolate the effects of one parameter at a time. The accessibility factor is the ratio of accessible pore volume to total pore volume. In the case of Avicel biomass particle, all the pores are assumed to be accessible and φ value was calculated to be 1 based on this assumption. But for DAP poplar, Ankathi et al. (Ankathi et al., 2019) measured the accessibility factor to be 0.92.In the case for DAP poplar, enzymatic hydrolysis occurred only for pore volumes above 51 Å, hence anything above 51 Å was considered accessible pore volume and anything below was not accessible by the enzymes. The total pore volume of the particles was measured using NMR cryoporometry.

The explicit finite difference method was used to solve the coupled equations for the multi-scale model using Equations 1–26, using 21-min dilute acid (H₂SO₄) pretreatment (DAP) poplar data from Table 2. DAP poplar data was used because these experiments included a measure of porosity change over pretreatment severity, which established bounds on porosity for this study. Numerical finite difference expressions with an accuracy of order (Δx) (Nanda et al., 2014) were utilized for all derivative terms. The biomass particle was segmented into 20 equally spaced sections, totaling 21 nodes in the radial direction. Time steps were kept below the stability criteria as shown in Equation 28 to assure a stable and accurate solution.

0.5 = Δ t D P C E Δ r 2 ϵ + ρ a q m C E K a 1 + K a C E 2 − 1 (28)

The explicit finite difference numerical solution was implemented using Microsoft Excel to predict transient non-adsorbing and adsorbing enzyme concentrations at different radial locations at the particle level and in the bulk solution. The accuracy of the explicit finite difference solution was confirmed for the non-adsorbing enzyme in the dilute limit of biomass loading (bulk concentration of enzyme remaining constant) by comparing it to the analytical solution, as expressed above, where the agreement was within 1% for characteristic time (time required for enzyme concentration at r = 0 to be 99% of the bulk concentration). The effect of particle size (diameter of 0.005 cm–0.05 cm) and adsorption capacity (10–200 mg/g substrate) on glucose yields from enzymatic hydrolysis of dilute acid pretreated poplar (21 min DAP poplar) was predicted for a range of porosity values in separate cases.

4 Results and discussion

4.1 Particle level diffusion and adsorption results

4.1.1 Non-adsorbing enzymes

For non-adsorbing enzymes, the effect of porosity on diffusion is shown in Figure 3A, assuming that all pores are accessible to the enzyme for diffusion into the biomass particle and holding particle radius and adsorption parameters constant. The characteristic time for diffusion decreases by a factor of nearly two with increasing porosity from 0.1–1.0. This trend in the characteristic time for diffusion is due to the effect of porosity on the effective diffusion coefficient according to Equation 3. According to Equation 3, the characteristic time for diffusion of non-adsorbing enzymes is inversely proportional to D _*p*/ε and therefore varies as (2-ε), which over the range 0.1<ε < 1, varies by nearly a factor of 2, in agreement with Figure 4A.

The effect of particle radius alone on characteristic time for diffusion of a non-adsorbing enzyme is shown in Figure 4B. The characteristic time is lowest for particle radius R = 0.0025 cm (16.65 s) and highest for particle radius of R = 0.025 cm (1,665 s). Thus, with a 10x increase in particle radius, there is a 100x increase in characteristic time for diffusion, which is the expected result for purely diffusive transport. Based on the results of this study, if biomass particle radius were to increase much above 0.25 mm, the characteristic time for diffusion of non-adsorbing enzymes into the particle would increase by many thousands of seconds.

4.1.2 Adsorbing enzymes

For adsorbing enzymes, the effect of porosity on diffusion is shown in Figure 5A (while keeping all other parameters constant). The characteristic time for diffusion decreases very strongly with increasing porosity in the range of 0.1< ε < 1.0. For a large adsorption term in Equation 16, as porosity increases, the effective diffusion coefficient in this equation changes by a factor of 20 according to Equation 3. The decrease in diffusion characteristic time with increasing porosity is greater for adsorbing enzymes compared to non-adsorbing enzymes.

For adsorbing enzymes, the adsorption term is much greater than ε, and over the range of 0.1<ε < 1, the term on the left-hand side of Equation 16 involving the Langmuir isotherm is nearly constant. For this case, the characteristic time will vary with changing ε as (2- ε)/ε, which is a factor of 19 over the range of 0.1<ε < 1, in agreement with Figure 5A. However, there are limitations to this analysis of characteristic time for adsorbing enzymes that are worth pointing out. The limitation of this study is that the apparent density ρ a and q m were kept constant throughout the study of the effect of changing ε, but these two parameters should vary with porosity. But due to the lack of data on how these two parameters vary with changing porosity, they were kept constant in our study. With sufficient information, the model could capture the interplay between these parameters and their effect on the characteristic time of diffusion. When the porosity increases, apparent density should decrease based on Equations 1, 11, and maximum adsorption capacity (q _*m*) should increase up to a finite maximum value determined by the number of adsorption sites per unit area on the cellulose surface. However, these counteracting trends in ρ a and q _*m* with increasing ε are not predictable at this time but must rely on experimental measurements to fully understand them.

As shown above for non-adsorbing biomass particles, the other major factor influencing the characteristic time of diffusion is particle size. The effect of particle radius on the characteristic time of diffusion of adsorbing enzymes is shown in Figure 5B. The characteristic time is lowest for particle radius R = 0.0025 cm (526.12 s) and highest for particle radius of R = 0.025 cm (52,613 s); that is, with a 10x increase in particle radius, there is a 100x increase in characteristic time for diffusion, which is also expected for purely diffusive transport. Based on the above model predictions, to reduce the processing time for enzymatic hydrolysis, the particle diameter must be kept to a smaller size range, and the biomass must be pretreated to increase the porosity. We provide more specific interpretations in the multi-scale model results section below.

To further understand the effect of substrate characteristics on the diffusion of enzymes, the maximum adsorption capacity was varied. The effect of maximum adsorption capacity on the diffusion of adsorbing enzymes is shown in Figure 5C. A linear trend between characteristic time and maximum adsorption capacity was predicted. The lowest characteristic time (16.67 s) was at (q m = 0 mg/g substrate) and highest characteristic time (2,925 s) was at (q m = 100 mg/g substrate). This shows that maximum adsorption capacity significantly affects the diffusive penetration of enzymes into a biomass particle. The Langmuir (or affinity) constant K _*a* indicates the extent of interaction between adsorbate and the surface. To understand the effect of the affinity constant on the diffusion of enzymes, the K a value varied from 0 to 10 mg/mL for Avicel PH-101 biomass particle, and results are shown in Figure 5D. When the value of K a was zero, the characteristic time of diffusion was 16.67 s, which means there was no effect of adsorption on the diffusion of enzymes. With increasing affinity constant, the characteristic time also increased to a maximum of 547.4 for K a = 3 mL/mg. After reaching the highest characteristic time, the influence of affinity constant on the characteristic time of diffusion of enzymes is rather modest.

In summary, the diffusion-adsorption model predicted an 18x decrease in characteristic time for diffusion of enzymes when the porosity was varied from 0.1 to 1, a 100x increase in characteristic time for diffusion when the particle size was increased by ten times the typical particle size of Avicel PH-101, a 180x increase in characteristic time of diffusion when the maximum adsorption capacity was increased over an expected range of values from the literature for Avicel PH-101, and a 50x increase in characteristic time when the affinity constant was increased over an expected range of values from the literature. For the range of parameter values obtained from the literature review, we can infer based on these model predictions that the maximum adsorption capacity has the strongest effect on the characteristic time of diffusion, followed by particle size, followed by affinity constant, and followed by porosity (q _*m* > R > K _*a*> ϵ).

4.2 Multi-scale diffusion and adsorption model results

4.2.1 Effect of biomass loading on characteristic time for Diffusion

To understand the effect of biomass loading to enzyme loading ratio on the characteristic time of enzyme diffusion into biomass, the ratio of biomass to enzyme loading was varied from the original ratio (Table 2: W _*biomass*, V, C _*CBo*, C _*CEo*) to 5x the original ratio by maintaining the initial enzyme concentration constant as the biomass loading was increased. The change in enzyme bulk concentration and enzyme concentrations at different radial positions in the biomass particle were predicted using Equations 7, 16 subject to Equations 20–27.

The graphs in Figures 6A,B show the concentration of CE and BG enzymes, respectively, at different radial positions inside the 21 min DAP pretreated biomass particle and in the bulk solution vs. time (secs) for the original biomass to enzyme ratio (Table 2). The concentrations of CE and BG enzymes in the bulk solution (C_CE and C_BG) decrease with time due to the diffusion of enzymes into the biomass particle. The concentration of adsorbing enzyme (C_CE) in the bulk solution in Figure 6A drops from 0.307 mg/mL to 0.26 mg/mL, which is a more significant drop on a percentage basis compared to the drop in C_BG bulk concentration from 0.1107 mg/mL to 0.106 mg/mL shown in Figure 6B due to the effect of adsorption. The concentration of enzymes is initially zero at all radial locations in the particle. As time progresses, the enzymes diffuse into the biomass particle, leading to the increase in enzyme concentrations at different radial positions with a faster increase at the outer radial positions (19, 18) compared to inner radial positions (0, 5, 10). The characteristic time of diffusion of adsorbing enzymes (CE) is around 5,000 s, which is greater than the characteristic time of diffusion of non-adsorbing enzymes (BG), which is around 1,500 . This difference in characteristic times is due to the effect of enzyme adsorption.

Figures 6C,D shows the concentration profiles of adsorbing and non-adsorbing enzymes at different radial positions and in the bulk for biomass loading to enzyme loading ratio equal to 5x the original, respectively. There are a few notable observations compared to the original biomass to enzyme loading case. The bulk concentration of adsorbing enzyme (CE) reduces by nearly half compared to the initial concentration for the 5x case due to accumulation of enzymes inside the biomass particle. The reduction is very rapid at early times and results in an “overshoot” in pore fluid concentration of adsorbing enzymes just inside the biomass particle (positions 19, 18). This same effect is observed for non-adsorbing enzyme (BG), but in a less pronounced manner. The enzyme concentration profiles at different radial positions and in bulk for different biomass loading to enzyme loading ratios (2x to 4x) are available in the SI document in Supplementary Figures S1–S3.

The multi-scale model also predicts that enzyme diffusion characteristic time is affected by the ratio of biomass to enzyme loading, as shown in Figure 7. The model predicted that characteristic time for diffusion of enzymes decreased with an increase in biomass loading for both adsorbing and non-adsorbing enzymes, with a more significant decrease observed for adsorbing enzymes. The characteristic time dropped from 5,107 s to 4,475 with a 5x increase in biomass loading for adsorbing enzymes, a 14.5% decrease. In the case of non-adsorbing enzymes, the characteristic time of diffusion dropped from 1,548 to 1,402 s, which is a 9.4% decrease. The reason for the decrease in characteristic time with increasing biomass to enzyme ratio is related to the decrease in bulk enzyme concentration, which requires less time to achieve 99% of the bulk concentration at r = 0 (definition of characteristic time) when bulk concentration drops more dramatically. Also, it may be interpreted that the “overshoot” in enzyme pore concentration contributes an enhanced diffusion driving force to shorten the characteristic time.

A sensitivity analysis of the characteristic time (CT) threshold value and of the external mass transfer coefficient (k (cm/s)) shown in Table 3 was completed using the model that yielded results in Figure 7 for biomass loading to enzyme ratio of 1. The sensitivity analysis shows that the change in CT for external mass transfer coefficient from 0.01–0.05 cm/s is very minor, indicating that external mass transfer resistance is not a limiting step in the enzyme diffusion process. This conclusion was also arrived at by Rohrbach and Luterbacher (Rohrbach and Luterbacher, 2021) in their in silico analysis of enzyme diffusion and reaction in biomass particles. However, the selection of CT threshold has a big influence on CT for diffusion of enzymes into biomass particles, confirming that the diffusion process reduces in rate significantly as the process approaches equilibrium (Table 3).

TABLE 3. Sensitivity of the predicted characteristic time (CT) to values of k (cm/s) and CT Threshold (%).

4.2.2 Effect of particle size, adsorption capacity, and porosity

The multi-scale diffusion and adsorption model was applied to a range of particle sizes and adsorption capacities for three different porosities. The particle radius (R) range starts from 0.005 cm (typical particle size of Avicel PH-101) to 0.05 cm (typical particle size of DAP poplar biomass), and the range of adsorption capacities (q m) start from 10 to 200 mg/g substrate based on the typical maximum adsorption capacities of different biomass found in different studies (Angarita et al., 2015; Zheng et al., 2016; Mota et al., 2019; Khodaverdi et al., 2012). The porosities used for the model are 0.5, 0.7, and 0.9, based on the porosities reported by Ankathi et al. for dilute-acid pretreated (DAP) hybrid poplar (Ankathi et al., 2019). In these multi-scale model runs, the value of the adsorption affinity parameter, K _*a*, remained constant at 1.0. This is justified because all of the studies in the literature that we found measured K _*a* values greater than 1.0, and the results in Figure 5D show that characteristic time is very insensitive to changes in K _*a* values greater than 1. The model predictions are shown in Figures 8A–C.

The characteristic time (CT) for the diffusion of adsorbing enzyme into 21 min DAP poplar with ε = 0.5 varied between 0.166 h and 233.36 h over the range of R and q _*m* values (Figure 8A). In the case of biomass with ε = 0.7 (Figure 8B), the CT for diffusion varied between 0.069 h and 87.09 h over the range of R and q _*m*. For biomass with ε = 0.9, the CT for diffusion ranged between 0.023 h and 19.52 h over the range of R and q _*m*. Based on the discussion of adsorption kinetic data in Table 1, at a diffusion CT > 4–6 h would result in diffusion limitations on enzyme hydrolysis of biomass. Also, in the in silico study by Rohrbach and Luterbacher, they predicted that for a biomass particle of radius 25 μm there are no external or internal mass transfer resistances to affect overall cellulose enzymatic hydrolysis (Rohrbach and Luterbacher, 2021). Furthermore, for this particle size and for a biomass loading of 2%, which is very close to our study’s 1%, the maximum glucose yield occurs between 5–15 h depending on enzyme loadings from 0.1 to 5.0 of the stoichiometric enzyme adsorption capacity. Based on these prior works, we selected a value of 6 h as a reference diffusion characteristic time to delineate parameter space in Figure 8 between diffusion limited and non-diffusion limited regimes.

In Figures 8A–C, the parameter space of R and q _*m* colored by orange, yellow, green, purple, and deep red are therefore diffusion-limited. Our modeling results shown in Figures 8A–C agree with Rohrbach and Luterbacher (Rohrbach and Luterbacher, 2021), who concluded that with expected inherent enzyme kinetics biomass particles greater than 500 μm in diameter (R > 0.025 cm) would be diffusion-limited for enzymatic hydrolysis. Our results provide additional insights into this conclusion and indicate that the transition to diffusion limitation is sensitive to adsorption strength and particle size. For example, even particles of R < 0.025 cm are expected to be diffusion-limited for the largest values of q _*m*. Porosity also affects this transition, with the observation that diffusion limitations move to larger biomass particle R values with increasing ε.

In reality, porosity does change during enzymatic hydrolysis of cellulosic biomass. The boundary for diffusion limitations would shift from the lower left toward the upper right in the 2-D plots of Figure 8 as porosity increases during enzymatic hydrolysis. An increasing porosity during hydrolysis would accommodate increased maximum adsorption capacity and increasing particle size before diffusion limitations are encountered compared to a particle where porosity remains constant during enzymatic hydrolysis. Therefore, a model with an assumed constant porosity set at the initial value would provide a conservative estimate for characteristic time before diffusion limitations are observed. The results in Figure 8 provide insights into the onset of diffusion limitations in biomass enzymatic hydrolysis with knowledge of particle dimension, porosity, and adsorption strength. With measured porosity and adsorption parameters, these modeling results would help optimize particle size to avoid diffusion limitations). It is recommended to use the results from Figure 8 to identify the biomass particle size reduction that is required to avoid diffusional limitations on enzyme hydrolysis performance. Furthermore, if the initial particle void fraction is selected for this size reduction determination, the choice will be conservative given that the void fraction increases during biomass hydrolysis, which is shown to reduce diffusional limitations as void fraction increases.

5 Conclusion

Diffusion of enzymes into biomass particles is a potential hydrolysis limiting step, and determining the characteristic time for enzyme diffusion as affected by different system parameters helps in understanding the importance of different parameters for the diffusion/adsorption step. The goal of this modeling study is to better understand the effects of system parameters, such as biomass particle size, porosity, adsorption capacity and affinity, and biomass to enzyme loading ratio on the characteristic time of hydrolytic enzyme diffusion into biomass particles. In addressing this goal, we developed both single-particle and multi-scale diffusion models. Key results from the single-particle model are that, over the range of parameters studied, the maximum adsorption capacity has the most substantial effect on the characteristic time for diffusion, then particle size, followed by affinity constant, and then porosity (q _*m* > R > K _*a*> ϵ). The multi-scale diffusion and adsorption model coupled particle diffusion with bulk changes in enzyme concentration. Key findings include the biomass loading effect on diffusion characteristic time, in which characteristic time decreases with increasing biomass to enzyme loading level. In addition, the transition to diffusion limitation is sensitive to both adsorption strength and particle size and is also strongly affected by particle porosity. This study focused only on the CT for diffusion/adsorption of hydrolytic enzymes within biomass particles. As a result, it is now better understood which system parameters affect diffusional CT. Diffusional CT may be short or long compared to enzyme hydrolysis characteristic time. Thus, if enzyme kinetics are very slow, even large biomass particles will not experience diffusional limitations. Ultimately, we must combine both diffusional and enzyme hydrolysis kinetics into system models, as others have done before, yielding important insights. However, those prior studies did not investigate the importance of system parameters in the important diffusional step, which our study does for the first time. Building on this research, future studies to improve predictions of enzymatic hydrolysis kinetics could incorporate more parameters into the modeling space, including changes in accessibility factor during hydrolysis, feedback between changes in ε and q _*m* and ρ _*a*, a thorough review of enzyme adsorption literature to define realizable combinations of R, q _*m*, and ε, as well as advances in modeling enzyme catalytic action at the pore surface. Another potentially promising area of future research is from Figures 4, 5, one can fit mathematical functions to get analytical expressions for the dependence of the characteristic times on parameters like porosity, Radius, etc. With more model complexity, the benefits of more accurate model predictions in reducing costs and environmental impacts of advanced biofuels production can be realized.

Generative AI statement

The author(s) declared that generative AI was not used in the creation of this manuscript.

Any alternative text (alt text) provided alongside figures in this article has been generated by Frontiers with the support of artificial intelligence and reasonable efforts have been made to ensure accuracy, including review by the authors wherever possible. If you identify any issues, please contact us.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fctls.2026.1761691/full#supplementary-material

Supplementary material Supplementary material at the source

Footnotes

Edited by: Putrakumar Balla, Chungnam National University, Republic of Korea

Reviewed by: Partho Sakha De, Heinrich Heine University of Düsseldorf, Germany

Praveen Thomas, Vellore Institute of Technology (VIT), India

Footnotes

AFEX, ammonia fiber expansion; BGβ-glucosid ase enzyme; CE, Cellobiohydrolase and Endoglucanase enzymes; DAP, Dilute acid pretreatment; EH, Enzymatic hydrolysis; WS, Pretreated wheat straw.