TPMS-based membrane lung with locally-modified permeabilities for optimal flow distribution

Simulation

All simulations were executed using Ansys CFX 19.0 (Ansys, Inc, Canonsburg, PA, USA). Simulation domains were isolated from CAD models using Ansys SpaceClaim, and all computational meshes were created using Ansys Meshing.

TPMS permeability data

Deforming individual TPMS elements facilitated the control of local permeability in the membrane module. In this study, the Schwarz-P (SWP) shape was used. The surface F_SWP of a single SWP element is described by the following implicit equation:

Under normal conditions, the coefficients of each of the cosine terms of the SWP implicit function are equal to 1. Increasing the coefficient preceding a cosine term for a particular axis decreases the smallest cross-sectional area for fluid flowing around the element along that same axis. However, the overall volume within the unit cell enclosed by the surface remains the same, therefore the volume porosity of the element remains the same.

The SWP element was simulated to determine the Darcy permeability for several geometries. The simulation was set up to replicate the experimental procedure described in Schlanstein et al.³⁶. Low flows were spread over theoretical modules created by translationally periodic boundary conditions, and pressure loss was calculated across the thickness of one layer of the module. In subsequent post-processing, the streamwise and transverse permeabilities were calculated for the standard as well as the “most permeable” (c = 0.1) and “most occluded” (c = 1.9) SWP geometries, illustrated in Fig. 6. The degree of occlusion (or inversely, the available cross section) to the flow determined the permeability of the element. For example, for the most permeable streamwise SWP geometry, a higher cross-sectional area for flow is available associated with lower pressure losses than for the most occluded, streamwise configuration.

To obtain a wide range of permissible permeability, a combination of 1 mm and 3 mm SWP elements was conceived for the bundle. Figure 6a shows the range of streamwise permeabilities that were achievable with both element sizes for the different states of deformation. The relationship between the streamwise and transverse permeabilities of 1 mm and 3 mm elements can be seen in Fig. 6b. The linear regression was used in the optimization scheme. As the design intent was to facilitate as many 1 mm elements in the module as possible, the threshold value for the decision to use 3 mm elements instead of 1 mm elements was set to the upper limit of the 1 mm elements’ streamwise permeability range (8.68e−9 m²).

Initial flow field

A geometry corresponding to that of the Novalung Interventional Lung Assist (iLA) was simulated to establish a comparative baseline. The simulation geometry was created in Creo Parametric 4.0 (PTC, Boston, MA, USA). All dimensions for the geometry were measured on a disassembled iLA device (Fig. 1b). The fiber bundle was modeled as a porous domain with a constant porosity of 0.493 and an anisotropic permeability. Streamwise permeability was set to 10.88e−10 m², and transverse permeability was set to 7.71e−10 m² in accordance with previously reported data for permeabilities of stacked fiber mats³⁶.

The computational domain of the TPMS-based prototype device was identical to that of the predicate device except for the removal of the distributor plates and the widening of the membrane module to occupy the resulting empty space. Removing the distributor plates also meant that a symmetrical boundary condition could be assumed along the diagonal plane from bottom to top of the prototype device, reducing computational effort. Blood was modelled as a Newtonian fluid with a density of 1059 kg/m³ and a dynamic viscosity of 3.6 mPas. Flows of 0.5, 1, 2, 3, 4, and 4.5 L/min were tested, corresponding to the range of flow rates for which the iLA is clinically approved.

Iterative flow field optimization

With the range of achievable permeabilities in hand, the goal of the iterative optimization scheme was then to determine how those permeabilities should be distributed in the porous domain in order to satisfactorily distribute blood flow throughout the module, creating an optimized flow field. An optimized flow field in this context meant a flow field with homogeneous flow velocities. Practically speaking, between two iterations, the optimization algorithm would lower the permeability at every point in the fiber bundle where the velocity was too high and raise the permeability at every point where the velocity was too low. This iterative process was stopped once the difference between the variable velocity v_y and the ideal velocity v_ideal was acceptably small or the change in difference between two consecutive iterations stagnated. This meant that the simulation workflow needed to (a) prescribe the value of streamwise and transverse permeability over the porous domain according to data provided from the previous iteration, (b) calculate the desirable value of the streamwise permeability based on the deviation of the real flow to the desired flow, and (c) repeat, using the updated permeability values as the provided data in step (a).

After a simulation solution for the given iteration was found, a correction function was defined that determined the new pointwise permeability value to use for the following iteration. The correction function took the form of a hyperbolic tangent function:

where K_new is the streamwise permeability for the next iteration, K_old that of the current iteration, v_ideal is the ideal flow velocity in streamwise direction, v_y is the streamwise flow velocity in the current iteration, and p is a constant used to adjust the rate of change between iterations. A 1 L/min flow rate was used for all iterations of the optimization process. The membrane module of the prototype device, like that of the iLA, had a footprint of 100 mm × 100 mm. Accounting for the volume porosity of the SWP elements (0.5), this gave an area average velocity v_ideal of 3.33 mm/s. Initially set at 150, this p constant was slowly decreased to 1 as the iterative process proceeded in order to balance rapid progress towards the ideal result with higher precision once the given solution approached the ideal.

In general, it is helpful to establish meaningful initial conditions for the distribution of permeability a priori in order to enhance convergence. Here, a 3D permeability matrix was used to define the streamwise permeability K_init at each point in the computational mesh based on the following equation:

where r_corner (m) is the distance from the lower corner. The dependent transverse permeability was then defined from the streamwise permeability (see Fig. 6b).

Design and manufacturing

Creation of SWP elements

After the optimization process had been conducted, a 3D matrix of point-wise permeability values was obtained with the spatial resolution of the computational mesh. These simulated permeabilities were interpolated in Matlab (version 2019a, MathWorks, Natick, MA, USA) so that a permeability value could be obtained for a value at any point within the spatial domain, not simply the points in the computational mesh. Iteratively, a 3 mm × 3 mm × 3 mm region was queried. If the mean value of permeabilities within that region was above the threshold permeability between 1- and 3-mm elements (8.68e−9 m²), the quadratic regression equation for 3 mm elements’ permeabilities was used to determine what coefficient for the cosine term in the SWP implicit equation was used to produce that permeability. Alternatively, if the mean permeability fell below the threshold value, a cubic grid of nine 1 mm elements was created, using instead the regression equation for 1 mm elements. An STL mesh for each element was created (Fig. 7a,e) using its calculated cosine coefficient (Fig. 7b,f). As they were, SWP elements with disparate coefficient values or sizes would not fit together at their interface. Therefore, within the first and last eighth of the element’s length in each direction, a transition was defined which adjusted the cosine coefficient and size of the current element to match those of the neighboring element (Fig. 7c,g). If an element was not connected to another element in a certain direction, the SWP-surface was edited to include a face sealing the element in that direction (Fig. 7d,h). The SWP-surface meshes were then exported as STL files. This was repeated for every element that would make up the membrane module.

Assembly creation and post-processing

Alongside the individual STL files, the x-, y-, and z-coordinates as well as the size of each element was recorded. This information was then used to import the STL files into a 3-matic (Materialise NV, Leuven, Belgium) project. After translating and scaling each STL file, the single elements were merged into the full membrane module. Due to restrictions regarding the file size of the resulting data, 3-matic’s adaptive remeshing capabilities were used to reduce the mesh complexity before exporting the completed STL file.

Manufacturing of the prototype

After preparation of the CAD files had been completed, manufacturing of the individual components and assembly of the prototype device for feasibility and validation purposes could begin. In order to reach the high spatial resolutions required to model the SWP elements, a stereolithography-based 3D printing technique was used to create the membrane modules (Materialise GmbH, Leuven, Belgium). After printing was completed, the remaining resin was washed away from the structure by a combination of pressurized air, a brief soaking in a 50% isopropanol/water solution, and centrifugation. The flow antechambers on either side of the membrane module were manufactured using a material jetting 3D printing technique (Objet 350 Connex3, Stratasys Ltd., Eden Prairie, USA) and VeroClear material.

Other than the removal of the diffuser plates and the consequent thickening of the membrane module to 24 mm to accommodate the empty space, the prototype device remained identical to the predicate.

Benchtop washout tests

In order to validate the simulation results as well as to provide a real-world comparison between the predicate and prototype devices, a series of washout tests was conducted. For the tests, the predicate and prototype devices were measured at the same flow rates as had been simulated. A flow circuit was designed to determine fluid residence times in each device at different flow rates (Fig. 8). The main circuit consists of a glycerol reservoir containing a transparent water glycerol mixture, a pump (deltastream DP 2, Xenios AG, Heilbronn, Germany), and the test module. An additional side arm connected to this circuit in front of the test module including a reservoir with dye solution (water glycerol solution with the same viscous properties as the circulating fluid) and including the same type of pump running at the same operational point. Behind the test module, another side arm leads to a waste reservoir set to the same hydrostatic level as the reservoir with the glycerol solution. A series of magnetic tube clamps (Fluid Concept GmbH, Stutensee, Germany) were positioned on each tubing branch as illustrated in Fig. 8. Initially, the transparent glycerol solution circulated in a closed loop (state of the magnetic tube clamps during ‘circulation phase’: 1 = closed, 2 = open, 3 = open and 4 = closed). The magnetic tube clamps were controlled remotely such that ink could be injected for a prescribed amount of time (state of the magnetic tube clamps during ‘injection phase’: 1 = open, 2 = closed, 3 = closed and 4 = open). During the measurement phase, the polluted fluid was sent to the waste reservoir to keep the remaining fluid transparent (state of the magnetic tube clamps 1 = closed, 2 = open, 3 = closed and 4 = open). After ink was injected, normal flow was resumed. Photometric color sensors (TCS34725, Taos, Inc., Plano, TX, USA) were used to monitor the flow of the injected ink bolus into and out of the module. The photometric sensor readings were continuously recorded for later analysis, and flow rate was monitored with an ultrasonic flow probe (Transonic Systems Inc., Ithaca, NY, USA). Once the dyed fluid without any remains left the circuit, the magnetic tube clamps were reset to the initial circulation phase.

Prior to testing, the photometric color sensors were calibrated with a series of ink-water/glycerol dilutions to obtain a concentration curve. Additionally, flow rate was shown to remain constant during the time when ink was flowing through the module, therefore the mass flow rate of ink could be calculated. From this, threshold values were defined for the minimum, mean, and maximum residence times. Minimum residence time was defined as the time between when 1% of the total dye mass had passed each sensor, mean residence time as the time between the passage of the largest single amount of ink, and maximum residence time as the time between when 95% had passed each sensor.

The predicate and prototype devices were measured at the same flow rates as had been simulated. All tests were conducted using a η = 3.78 ± 0.13 mPas water/glycerol mixture (n = 5 measurements at 0.01–10 Pa, tested via MCR502 rheometer, Anton Paar GmbH, Graz, Austria).

Ethical approval

This article does not contain any studies with human or animal subjects performed by any of the authors.