Coupling CFD and Machine Learning to Assess Flow Properties in Porous Scaffolds for Tissue Engineering

TL;DR

This paper uses computational fluid dynamics and machine learning to study how scaffold structure affects fluid flow properties important for tissue engineering.

Contribution

The novel use of ML to predict flow properties from scaffold topology offers a new design tool for tissue engineering scaffolds.

Findings

ML models achieved R2 above 0.9 for predicting permeability and average wall shear stress from scaffold topology.

Percentile-based wall shear stress metrics provided insights into environments relevant for bone and cartilage differentiation.

Pore shape had no significant effect on permeability or average wall shear stress.

Abstract

Computational fluid dynamics and machine learning (ML) models are employed to investigate the relationships between scaffold topology and key flow parameters, including permeability (k), average wall shear stress (WSSa), and the 25th and 75th percentiles of WSS. Statistical analysis showed that WSSa values are consistent with those found in common scaffold architectures, while percentile-based WSS properties provided insight into shear environments relevant for bone and cartilage differentiation. No significant effect of pore shape was observed on k and WSSa. Correlation analysis revealed that k was positively associated with topological features of the scaffold, whereas WSS metrics were negatively correlated with these properties. ML models trained on six topological and flow inputs achieved a performance of R2 above 0.9 for predicting k and WSSa, demonstrating strong predictive capability based on the topology. Their performance decreased for WSS25% and WSS75%, reflecting the difficulty in capturing more specific shear events. These findings highlight the potential of ML to guide scaffold design by linking topology to flow conditions critical for osteogenesis.

Figures

Full text

1. Introduction

Tissue engineering is an emerging field in biomedical engineering that leverages both life sciences and engineering to repair, maintain, or restore tissue functions, involving biomaterials, cells, and growth factors, among others [1,2,3]. This approach has gained increasing interest in repairing bone defects, as scaffolds can be used to mimic the bone matrix structure, providing an environment for cell attachment, growth, and tissue development. While scaffolds provide mechanical stability to the tissues, they also should exhibit a porous structure to promote osteoconductivity and osteoinductivity, resulting in additional challenges [4,5,6].

Previous studies have shown that the permeability of the scaffold and the wall shear stress (WSS) exerted by fluid flows play a relevant role in osteogenesis, that is, cell adhesion and bone regeneration. However, both properties are influenced by the pore topology, hindering the evaluation of the performance of the scaffold [7]. Computational fluid dynamics (CFD) has emerged as a power computational technique to investigate in detail the fluid flow in scaffolds, such as nutrient diffusion [8], permeability and porosity relantionships [9], and WSS distribution on scaffold surfaces [10], among others. Additionally, Rabiatul et al. explored the relationship between permeability and porosity for different trabecular bone orientations. The authors observed a good correlation, although it was dependent on structural orientation [11]. Fu et al. evaluated different scaffold models for trabecular bones, indicating that bone regeneration was affected by the topology [10]. Ali et al. conducted additional studies focusing on the characterization of permeability and WSS in scaffolds. They found that structures with minimal variation in channel size exhibited increased permeability. Nevertheless, no clear correlation was observed between regular topologies and WSS [9]. Suffo & López-Marín conducted a combined fluid–structure interaction study to evaluate the mechanical performance of scaffolds subjected to fluid flows at high Reynolds numbers, comparing various turbulence models using inlet velocities of 2 m/s [12]. More recently, Manescu et al. performed simulations on Mg-based alloys, known for their biodegradable capability, showing that both the mechanical performance and permeability of these structures fell within the biological range of trabecular bone [13]. Kong et al. explored various scaffold shapes, evaluating flow rates and WSS across different regions of the scaffolds to assess their capacity for efficient nutrient exchange, among other factors [14]. Additional CFD studies have investigated heterogeneous porous scaffolds with Voronoi-based architectures, showing that variations in porosity influence the internal velocity and pressure distributions, as well as permeability, affecting their potential for bone tissue engineering [15].

While numerous studies have explored the flow effects on scaffolds, there is a lack of systematic investigations addressing the relationships between scaffold topology and both permeability and WSS. As reported by Kong et al., biological phenomena such as bone tissue, cartilage tissue, and fibrous tissue differentiation are associated with distinct WSS ranges [14]. This highlights the importance of systematically examining the aforementioned relationships. Here, machine learning (ML) emerges as a formidable approach to establish predictive models for permeability and WSS based on topological properties. Several previous studies have implemented supervised learning models to predict outcomes such as the rupture risk of cerebral aneurysms [16,17], cardiovascular diseases [18,19], and the transmission of airborne infections [20,21]. These studies underscore the predictive power of ML, suggesting its strong potential for advancing scaffold design and optimization.

In this study, the permeability and WSS of three different types of scaffolds with regular structures were investigated using CFD simulations. Scaffolds with circular, cubical, and hexagonal pores were constructed, with diameters ranging from 0.2 to 0.5 mm, giving a total of 12 different topologies. Supervised machine learning models were implemented to predict the permeability and WSS on the scaffolds based on the inlet velocity, pore surface, porosity, and structure type, resulting in models that captured the underlying relationships. This integrated computational approach provides an efficient framework for evaluating scaffold designs, paving the way for data-driven optimization in tissue engineering.

2. Materials and Methods

2.1. Scaffold Construction

Cubic scaffolds with 6 mm edges and regular pore distributions were constructed in FreeCAD 1.0 software. Three pore shapes were considered, namely circular, square, and hexagonal. For each shape, four pore diameters (0.2–0.5 mm) and three pore numbers on the front face (4, 9, and 16) were modeled, resulting in a total of 36 scaffold configurations. As an example, Figure 1 illustrates scaffolds with 0.5 mm diameter pores in circular, square, and hexagonal topologies. The scaffolds were placed inside a box with dimensions of $[eqn]$ mm^3^ for the CFD simulations, as described below.

2.2. Boundary Conditions and Numerical Schemes

The incompressible and laminar Navier–Stokes equations were used to model the fluid flow through the scaffolds. No-slip boundary conditions were applied to the box and scaffold surfaces, with a constant inlet velocity and zero outlet pressure. A scheme is shown in Figure 2. Four inlet velocities were considered, namely 0.1, 0.2, 0.3, and 0.4 mm/s, as they correspond to typical velocities found in trabecular bones [22,23]. The fluid was modeled as blood, with a constant density of 1065 kg/m^3^, using the non-Newtonian Carreau blood model with parameters $[eqn]$ s, $[eqn]$ , $[eqn]$ Pa·s, and $[eqn]$ Pa·s, as reported in a previous study [24]. The value of $[eqn]$ has been typically used as the constant dynamic viscosity for blood in Newtonian models [25,26,27]. A constant inlet velocity was employed as blood flow, and blood volume in trabecular bone does not vary significantly with the cardiac cycle [28,29]. The equations were solved using the SIMPLE algorithm with second-order schemes for velocity and time. Simulations were performed for a total of 200 steps with a time step size of 0.005 s. The residual criterion of convergence for the continuity and velocity components were set to $[eqn]$ . All CFD simulations were carried out using Ansys Fluent 2025 R1.

2.3. Meshing

All samples consisting of the scaffold and surrounding box were meshed with tetrahedral elements using Ansys meshing. To inspect the effect of mesh refinement, different mesh densities ( $[eqn]$ ) were analyzed, ranging from 862 to 4036 $[eqn]$ . The average WSS ( $[eqn]$ ) over the scaffold surface and the permeability (k) were calculated using the scaffold with a 0.2 mm circular pore diameter. The permeability was calculated by transforming Darcy’s law using the flow rate $[eqn]$ , with A is the scaffold area and v the inlet velocity, leading to

$$[eqn]$$

where v is the inlet velocity, $[eqn]$ is the fluid viscosity within the scaffold, L is the length of the scaffold, and $[eqn]$ is the pressure drop across the scaffold [10]. As observed in Figure 3, both k and $[eqn]$ increased with the mesh density and reached relatively stable values for densities above 3000 $[eqn]$ . Nevertheless, the total variations in both quantities between the lowest and highest mesh densities were below 4.5%, indicating limited sensitivity in this range. Based on these observations, a mesh density of ∼ 3000 m^−3^ was adopted for all simulations in this study. We note here that mesh sensitivity analysis was also conducted on the 0.3 mm, 0.4 mm, and 0.5 mm circular pore scaffolds, resulting in no statistically significant differences for the larger densities discussed here.

2.4. Machine Learning Methods

A total of 144 cases were obtained, corresponding to the three pore shapes, four pore diameters (D), three different numbers of pores in the front face (N), and four inlet velocities (v). Additional topological descriptors were also calculated, including the porosity degree ( $[eqn]$ ), total scaffold surface ( $[eqn]$ ), and pore surface ( $[eqn]$ ). These quantities were considered as input properties. Flow parameters were then calculated to evaluate the scaffold capability for osteogenesis, including permeability (k), average WSS ( $[eqn]$ ), first-quartile WSS ( $[eqn]$ ), and third-quartile WSS ( $[eqn]$ ). These four properties served as response variables. Regression models were implemented to predict k, $[eqn]$ , $[eqn]$ , and $[eqn]$ based on the input properties. Permeability and average WSS were chosen, as they are commonly employed to evaluate the scaffold capability for osteogenesis, while $[eqn]$ and $[eqn]$ were chosen because they provide insights into the localized shear stress environment, which is critical for guiding tissue differentiation processes such as bone and cartilage formation. All variables were standardized for comparison purposes using

$$[eqn]$$

where $[eqn]$ is the original variable, $[eqn]$ is its mean, and $[eqn]$ is its standard deviation. Data handling and transformation were conducted with the Pandas library in Python 3.13.5 [30].

Linear regression, support vector regression (SVR), K-nearest neighbor (KNN) regression, and random forest (RF) regression were explored to predict the flow properties through the scaffolds [31,32,33]. Table 1 shows a summary of the hyperparameters used for each case. The models were trained using 10-fold cross-validation [34], and their performance was evaluated by averaging the coefficient of determination ( $[eqn]$ ) across all cases. $[eqn]$ is defined as

$$[eqn]$$

where n is the total number of observations (144 in this work), $[eqn]$ is the actual value obtained from the CFD simulations, $[eqn]$ is the predicted value by the regression model, and $[eqn]$ is the mean of $[eqn]$ . All calculations were carried out using the Scikit-learn library in Python [35].

3. Results

3.1. Statistical Analysis

Statistical exploration of the k, $[eqn]$ , $[eqn]$ , $[eqn]$ , and $[eqn]$ variables was conducted to gain insights into their behavior. Table 2 shows the mean ( $[eqn]$ ), standard deviation (s), minimum, and maximum values for each variable. The values for k were slightly higher compared to previous studies with similar pore diameters [10,13] but still within the range observed in bulk cancerous bone [36,37] and solid gyroid scaffolds [38]. To further decrease the permeability, the models should be scaled down while maintaining the porosity degree [39]. For $[eqn]$ , the values range from 0.005 Pa to 0.416 Pa, with a mean of 0.148 Pa, in agreement with those reported for other topologies such as lattice diamond, octet, and gyroid [9,14]. The $[eqn]$ and $[eqn]$ were also considered, as they provide insights into the differentiation processes occurring within the scaffold [14,40]. The mean value for $[eqn]$ is 0.004 Pa, which is beneficial to bone tissue differentiation, while the mean value for $[eqn]$ is 0.276, which is beneficial to cartilage tissue differentiation.

The relationship between k and $[eqn]$ and the pore shape was also examined. Table 3 presents the average value for each quantity. Although the values are relatively similar across the groups, a single statistical measure is not enough to determine whether pore shape has a significant effect on the flow properties. To address this matter, the non-parametric Kruskal–Wallis test was performed, yielding p-values of 0.052 for k and 0.11 for $[eqn]$ . These results indicate that there is no statistically significant difference between the groups [41,42] at the conventional 0.05 level, suggesting that pore shape does not have a relevant effect on the permeability and average WSS.

3.2. Correlation Analysis

The relationships between k and the input properties D, N, v, $[eqn]$ , $[eqn]$ , and $[eqn]$ are shown in Figure 4, where each data point corresponds to an individual simulation. Increasing trends are observed in all cases except for v, where a slight decreasing variation occurs. As reported in previous studies, larger pore diameters, surfaces, and porosity degrees promote enhanced flow through the scaffold. In the case of v, the slight decreasing trend observed for a few cases may be attributed to minor deviations from Darcy’s law.

When inspecting the behavior of $[eqn]$ , decreasing trends are observed for all cases except for v, as shown in Figure 5. Previous studies have attributed the inverse relationship between porosity and WSS to flow redistribution and energy dissipation within the porous structure. A higher number of pores results in more interconnected channels that facilitate fluid flow, causing the flow to be more uniformly distributed throughout the scaffold surface. This redistribution reduces the concentration of high velocities near the walls, decreasing velocity gradients and, consequently, the WSS [9,43,44]. Additionally, the increasing trend of $[eqn]$ with v corresponds to the enhanced shear forces exerted by the fluid as velocity rises, leading to greater WSS on the scaffold surfaces.

The Spearman’s correlation coefficient ( $[eqn]$ ) was calculated to evaluate the degree of correlation between the different properties considered in this study. $[eqn]$ is given by

$$[eqn]$$

where $[eqn]$ represents the difference between the two ranks of the pair of properties and n is the number of observations, which is 144 in this work. A positive $[eqn]$ indicates a positive correlation between the two pair of variables, while a negative $[eqn]$ corresponds to a negative correlation. A good correlation is inferred when $[eqn]$ and a strong correlation when $[eqn]$ .

The resulting values are shown in Figure 6. Positive correlations are found between topological properties such as the pore number, porosity degree, pore surface, and total scaffold surface, while no correlation is observed for v with these properties. This is expected since the inlet velocity is set as a boundary condition in the CFD simulation, whereas the topological properties are intrinsic to the scaffold design. An inspection of the permeability reveals that relevant positive correlations exist with the pore diameter, porosity degree, total scaffold surface, and pore surface, reflecting the enhanced flow through the scaffold. On the other hand, negative correlations are observed for $[eqn]$ , $[eqn]$ , and $[eqn]$ with the topological properties, consistent with the reduced concentration of high velocities at the scaffold surface discussed previously.

3.3. Predictive Models

Seven quantities were initially considered as input properties to predict flow parameters on the scaffolds:, namely the pore shape, pore diameter (D), pore number in the front face (N), inlet velocity (v), porosity degree ( $[eqn]$ ), pore surface ( $[eqn]$ ), and total scaffold surface ( $[eqn]$ ). However, based on the exploratory statistical analysis, pore shape was discarded, as it showed no correlation with permeability or average WSS. Consequently, the predictive models were constructed as follows:

$$[eqn]$$

where y is the predicted flow parameter, representing either k, $[eqn]$ , $[eqn]$ , or $[eqn]$ , and f is the regression model

In the 10-fold cross-validation process, the coefficient of determination ( $[eqn]$ ) was calculated for each model, yielding ten values per case. The mean $[eqn]$ was then computed for both the training and testing sets, as shown in Figure 7, where black bars represent the standard deviations. The models exhibited excellent performance for permeability, with predictive capability exceeding 90%. Similar results were obtained for $[eqn]$ , indicating that the flow properties commonly used to evaluate scaffold suitability for osteogenesis can be predicted from topological parameters. In contrast, model performance decreased significantly for $[eqn]$ and $[eqn]$ , with larger standard deviations observed in the testing sets. This effect was particularly pronounced for $[eqn]$ , suggesting that predicting localized low-shear and high-shear regions is inherently more challenging due to the higher sensitivity required to capture subtle flow disturbances. Interestingly, the performance across the different ML models is very similar, indicating that the flow parameters are strongly governed by the underlying scaffold topology rather than the choice of the regression algorithm. This robustness suggests that relatively simple models, such as linear regression, can achieve predictive accuracy comparable to more complex approaches, in contrast to other studies where more sophisticated methods, such as RF, achieved superior performance [45,46,47,48]. The dataset generated in this study can be found in the Supplementary Materials.

4. Conclusions

This study combined computational fluid dynamics and machine learning models to investigate the influence of scaffold topology on flow-related properties, including permeability (k), average wall shear stress ( $[eqn]$ ), and the 25th and 75th percentiles of WSS. The analysis demonstrated that topological parameters such as pore diameter, porosity, and surface areas strongly influence k and $[eqn]$ , while pore shape had no statistically significant effect on k or $[eqn]$ . Thus, the effect of shape appeared to be much weaker compared to diameter. The negative correlations between WSS measures and topological parameters reflect the effect of porosity in redistributing flow and reducing shear concentrations. Additionally, this study considered non-Newtonian effects through a Carreau fluid model, revealing differences in magnitudes compared to Newtonian assumptions. This highlights the importance of fluid rheology when modeling the scaffold performance under physiological conditions, a matter that should be further explored in future studies.

Machine learning models trained on topological properties achieved high predictive accuracy for k and $[eqn]$ , demonstrating their capability to predict relevant flow properties of the scaffold for osteogenesis. However, the models delivered reduced performance for predicting $[eqn]$ and $[eqn]$ , indicating that localized shear extremes are more sensitive to subtle geometric and flow variations.

Future studies should extend the analysis to a wider range of scaffold architectures, investigate transient and pulsatile flows and mechanical performance through fluid–structure interaction and non-regular topological structures, and incorporate experimentally measured flow data to further validate the predictive capability of ML approaches. Additionally, dimensionality reduction techniques should be explored to address input variable correlations, reduce model complexity, and further mitigate overfitting. The integration of advanced topology descriptors, realistic boundary conditions, and refined rheological models will enhance the ability to design scaffolds optimized for specific tissue engineering applications.