CHEMICAL ENGINEERING TRANSACTIONS VOL. 57, 2017 A publication of The Italian Association of Chemical Engineering Online at www.aidic.it/cet Guest Editors: Sauro Pierucci, Jiří Jaromír Klemeš, Laura Piazza, Serafim Bakalis Copyright © 2017, AIDIC Servizi S.r.l. ISBN 978-88-95608- 48-8; ISSN 2283-9216 A Tunable Microfluidic Device to Investigate the Influence of Fluid-Dynamics on Polymer Nanoprecipitation Stefano Cerbelli*a, Alessia Borgognaa, Maria Anna Murmuraa, Maria Cristina Annesinia, Cleofe Paloccib, Marco Bramosantib, Laura Chronopouloub aDip. di Ingegneria Chimica Materiali Ambiente, Sapienza Università di Roma, Via Eudossiana 18, 00184 Roma (Italy) bDip. di Chimica, Sapienza Università di Roma, Piazzale A. Moro 4, 00185 Roma (Italy) stefano.cerbelli@uniroma1.it Polymer drug-embedding nanocapsules are attracting increasing attention as effective tools for the targeted delivery of pharmaceutical molecules on specific biological tissues. Besides, it is well established that an effective selectivity of the delivery dictates that the size of the carrier particles be accurately controlled, thus maintaining the size dispersion of the particle population as low as possible. To this end, microfluidics-assisted precipitation provides a promising alternative to the traditional processes in that the structure of the flow - ultimately controlling the particle size distribution - can be reliably predicted from the solution of Navier-Stokes equations in the laminar regime. Notwithstanding the great potential provided by microfluidics techniques, much about the interaction between fluid-dynamics and polymer transport and precipitation is yet to be understood. In this work, we investigate polymer precipitation in a simple cross-junction inflow-outflow microchannel, which has proven a viable benchmark to gain insight into the physics of nanoprecipitation in that the particle size distribution is sensitively dependent on the flow operating conditions. Specifically, previous experimental work by some of these authors proved that average particle size can vary by an order of magnitude for operating conditions where the solvent flow rate varies by a factor of three, while keeping the non-solvent flow rate constant. The scope of this work is to show that such sensitive dependence on operating conditions finds direct correspondence in the kinematic structure of the flow, which undergoes abrupt qualitative changes in the same range of operating conditions, provided a fully three-dimensional solution of the incompressible Navier-Stokes equation (thus retaining the inertial term in momentum balance) is afforded. 1. Introduction Drug-loaded polymer nanocapsules are attracting increasing research interest thanks to their ability to release pharmaceutical molecules selectively and effectively on target tissues. To this end, biocompatible polymers such as Polylactic Acid (PLA) pol-(lactic-co-yglycolic) acid (PLGA) as well as polyethylene glycol (PEG) have been successfully used. For instance, PLGA nanoparticles entrapping organic actives have been exploited on cell cultures (Chronopoulou et al., 2013). In the case where PLGA is used, the formation of drug encapsulation can be accomplished by nanoprecipitation, i.e. by mixing a solution of a solvent (e.g. acetonitrile, ACN) containing both the polymer and the bioactive molecules dissolved in molecular form with a non-solvent, which is typically water (see Saad & Prud’homme (2016) and therein cited references for a recent account on flash nanoprecipitation). Both the nominal (i.e. average) particle size and the dispersion of the particle size distribution about the mean value depend sensitively on the mixing process involving the background solution of solvent and non-solvent, which defines how fast and how homogeneously supersaturation conditions are reached. The fast and homogeneous supersaturation implies that the nucleation process of NPs prevails over growth, e.g. yielding capsules of order hundred nanometers for mixing times in the order of milliseconds. At laboratory or even industrial scale, these rapid mixing conditions dictate that the flow regime be turbulent. Besides, homogeneous transport conditions are scarcely enforced in turbulent conditions, which often involve a wide range of length- and time-scales. For this reason, driven by swiftly spreading novel micro-fabrication techniques, a shift of attention has been driven towards scaling down the process equipment to characteristic linear dimensions of the order of few tens of µm, where fast yet controlled mixing conditions can be achieved DOI: 10.3303/CET1757143 Please cite this article as: Cerbelli S., Borgogna A., Murmura M.A., Annesini M.C., Palocci C., Bramosanti M., Chronopoulou L., 2017, A tunable microfluidic device to investigate the influence of fluid-dynamics on polymer nanoprecipitation, Chemical Engineering Transactions, 57, 853- 858 DOI: 10.3303/CET1757143 853 (Bally et al. 2012; Capretto et al. 2012). The present study focuses on previous experimental work performed by some of these authors, which showed how controlled particle size in the range of tens of nanometers can be obtained in a simple tunable microfluidic device (Chronopoulou et al. 2014). The core of the device is constituted by a zero dead-volume cross-junction that realizes fast mixing between solvent and non-solvent streams. We target specific ranges of operating conditions for which a sensitive dependence of particle size has been observed. For such conditions, an attempt to correlate particle size to mixing time showed that no simple relationship between these two quantities can be established. In this attempt, mixing time was estimated based on the assumption that the solvent stream undergoes simple hydrodynamic focusing in the device, where the water streams entering the lateral channels push the solvent into a single flux tube towards the centre of the discharge channel. In this work, we challenge this oversimplified kinematic picture and investigate in some detail the fluid dynamic structure arising from the cross-junction mixing geometry. 2. Device geometry and modeling approach. The geometry of the cross-junction constituting the core of the micromixing apparatus is represented in the left panel of Fig. 1. The ACN liquid stream embedding PLGA and the drug enters the equipment from the top channel, whereas pure water is added by the side channels. In the experiments described in Chronopoulou et al. (2014), the length of the outlet tube at the bottom of the cross-junction (not represented in the figure) is equal to 30 cm. Such length ensures that the precipitation process is completed before collecting the nanoparticles for the experimental determination of the particle size distribution. For this configuration, Chronopoulou et al. (2014) performed a series of experiments where the average particle size was determined at fixed flow rate of the non-solvent (water) and variable flow ratio, defined as ACN to water flowrate (see left panel of Fig. 1). As can be observed, the average particle size undergoes a sensitive non-monotonic variation as the ACN flowrate is increased, from a maximum size of order 180 nm to a minimum size of 20 nm. In the case shown, the classical approaches based on simple estimates of the mixing time associated with the background solution water/ACN cannot explain the abrupt change exhibited by the average particle size in such a narrow region of the parameter defining the operating conditions. Motivated by these results, we seek an explanation of this behaviour based on the fluid-dynamic analysis of the flow structure. Figure 1: Left Panel: Three-dimensional structure and two-dimensional projection (red-shaded surface) of the cross-junction micromixer. Axes scales preserve the geometric ratios of the experimental apparatus. Axes units are expressed in meters. The overall length of the cross-junction modelled in the computational domain is equal to 2mm. Right panel: average particle size for fixed water flowrate (800 µl/min) and ACN-to-water flowrate ratio ranging from 0.0065 to 0.150. One notes that, even leaving aside the complex phenomenology associated with polymer precipitation (see, e.g. Di Pasquale et al. (2012); Di Pasquale et al. (2013) for a description of numerically-solved models accounting for full coupling between polymer transport and precipitation and solvent transport momentum balance), the mere prediction of the mixing structures resulting from the transport by advection and diffusion of the solvent and non-solvent species represent a computational challenge in that -(i) the solution of the Navier- 854 Stokes equations defining the convective field is in principle coupled with advective-diffusive transport of the solvent and non-solvent species because both the density and the viscosity of the mixture are composition- dependent; -(ii) the characteristic Peclet number in the targeted experimental conditions is order 104-105, making it necessary to use a highly refined discretization to resolve the sharp boundary layers that are to be expected in the steady-state concentration profiles; -(iii) the lack of azimuthal symmetry peculiar of the geometry suggests that a two-dimensional reduction of the flow domain might be inappropriate for capturing the relevant flow features. In the next Section, we use a two-dimensional reduction of Navier-Stokes and species transport equations to assess the impact of two-way coupling between momentum and mass transport. The conclusions stemming from this analysis are thus exploited to build up the three-dimensional flow model, which we develop in Section 5. Concluding remarks and possible directions of future work are presented in Section 6. All simulation results presented were obtained with a commercial finite-element solver (COMSOL 5.2). 3. Assessment of two-way coupling effects between mass and momentum transport The first issue to be addressed is to determine the impact of two-way coupling between mass and momentum transport, which could be triggered by the viscosity and density dependence on the variable composition of solvent and non-solvent. Figure 2 depicts the comparison between the flow structure obtained by assuming constant density and viscosity and the flow structure as computed by solving the fully coupled momentum and species transport equation. In this regard, density and viscosity dependency on constituents molar fraction have been determined by interpolating available experimental data. Because of the large Pe value at the prevailing operating conditions, a two-dimensional setting for the simulation has been chosen. Figure 2: Comparison between the flow structure obtained by assuming constant density and viscosity and the flow structure as computed by solving the fully coupled momentum and species transport equation. The two- dimensional computational domain is the red-shaded surface represented in the left panel of Fig. 1. As can be gathered from the comparison represented in Fig. 2, the solutions are qualitatively and quantitatively similar, thus suggesting that the constant viscosity and density approximation can be enforced without loss of essential phenomenological fluid-dynamic information. This finding makes it possible to approach the investigation of the flow structure in the device by enforcing the fully three-dimensional incompressible NavierStokes equation. Because the ratio of ACN to water flowrate never exceeds 15%, constant properties of pure water are next assumed for density and viscosity regardless of the mixture composition. 4. Three-dimensional solution of inertial incompressible flow The observations discussed above motivate to approach the description of the flow structure through incompressible Navier-Stokes equations, where density and viscosity dependences on mixture composition are neglected. Because of the lack of azimuthal symmetry, a fully three-dimensional representation of the 855 mixing domain is enforced. The steady-state field equation to be solved can therefore be expressed in dimensionless form as 𝒗 ∙ ∇𝒗 = 1 𝑅𝑒 ∇2𝒗 − ∇𝑃 (1) where variables are made dimensionless by choosing the cross-junction internal radius, R, as reference length, water inlet velocity, U, as reference velocity, and where the dimensionless pressure includes the hydrostatic contribution of the gravitational field. As regards the boundary conditions, parabolic flow profiles with prescribed flowrate are enforced at the inlet cross-sections of the device, whereas uniform pressure is assumed at the outlet cross-section. The presence of a non-vanishing left hand side in Eq.(1) makes the problem nonlinear. Because of the moderate values of the Reynolds number and of the time-independence of the feeding conditions, steady (i.e. nonturbulent) flow is expected. For this reason, a typical approximation enforced in microfluidics is to neglect the nonlinear term, thus enforcing a strictly Stokes (creeping) flow, whose solution is numerically much easier to approach. Besides, the concepts of laminar and creeping flows are not strictly interchangeable, and the case examined here is precisely one example where the solution of the Stokes and of the Navier-Stokes equations differ significantly from each other. Figure 3 provides a clear illustration of this issue by showing the structure of the streamlines (i.e. the trajectories of passively transported fluid elements) starting at the same points of the inlet section at the top of the cross-junction, where ACN is fed to the micromixer. As can be observed, the flow structure is qualitatively different in the two cases, thus suggesting that the action of the nonlinear term cannot be overlooked in the geometry and operating conditions considered. Figure 3: Comparison between the Stokes (creeping) flow and the solution to Navier-Stokes equations for the same boundary conditions in the cross-section geometry. Water flow rate is 800 µl/min and ACN-to-water flow ratio is 0.0065. The streamlines represented (red lines) are started at the same points of the cross-section at the top of the device, where ACN is fed. Based on this observation, we next tackle the flow structure within the micromixer through the numerical solution of the full incompressible Navier-Stokes equation. 5. Qualitative structure of the flow at different values of the flowrate ratio Accurate finite-element numerical simulations were performed to determine the flow structure at fixed water flowrate of 800 µl/min and variable ACN-to-water flowrate ratio. The mesh used - spatially nonuniform to account for the steep gradients near the intersection of the cylinders - corresponds to order 106 degrees of freedom for the values of the primitive variables. A variety of different flow structures was unveiled at increasing flow ratio. Figure 4 depicts two such structures, represented in terms of streamlines. As can be noted, the flow corresponding to the lower ACN flowrate (left Panel in the figure) exhibits vortex structures that 856 are formed just above the water injection. These vortices are connected with two separate flux tubes at the exit cross-section of the device. The right panel of the figure represents the flow structure obtained at a value of ACN flowrate nearly ten times bigger than the previous case. Here, no separate flux tubes are observed. It is clear that these kinematic structures will determine different concentration profiles for the mixture of solvent and non-solvent components at the exit of the cross-junction, where the premixed stream enters the discharge tube. Because the concentration profile defines the features of the mixing process downstream the domain, and therefore ultimately the mixing time, the results found in this work shed some light as to why a simple estimate of the mixing time derived on the basis of a postulated existence of a focused stream of solvent might prove inappropriate to describe mixing in this device. Preliminary results targeting the steady-state concentration profiles (not shown for brevity) suggest that, depending on the value of the flow ratio, the stream of solvent entering the top cross-section of the device can either split into four or two flux tubes that run alongside the channel walls, or, for larger values of the ACN flowrate, attain the expected structure of a hydrodynamical focused stream flowing in the central core of the pipe. These different behaviours are consistent, at least qualitatively, with the different performance of the flash nanoprecipitator as regards the average particle size that has been observed in experiments. Figure 4: Comparison between the numerical solutions to the incompressible Navier-Stokes equations for a water flowrate equal to 800 µl/min and a value of ACN-to-water flow ratio of 0.0065 (left) and 0.0625 (right). The streamlines represented (red lines) originate at the same points of the cross-section at the top of the device, which constitutes the inlet of the solvent stream. 6. Conclusions Motivated by previous experimental studies showing sensitive dependence of average particle size on operating conditions in a cross-junction flash nanoprecipitator, we investigate the flow structure associated with a full three-dimensional steady-state solution of the incompressible Navier-Stokes equation in the assigned geometry. Possible significant two-way coupling effects between momentum transport and advective-diffusive transport of solvent (ACN) and non-solvent (water) species is ruled out on the basis of results obtained in a two-dimensional representation of the simultaneous coupled transport problem. So much established, we approach the flow structure in the three-dimensional setting. Preliminary investigations on the influence of the nonlinear term in Navier-Stokes equation suggest that the full nonlinear structure of momentum balance should be retained. Simulation results conducted at fixed water flowrate and variable ACN-to-water flowrate ratio show that in the same ranges of operating conditions that yield sensitive dependence of particle size, qualitatively different flow features characterize the fluid-dynamics response of the device. This finding provides a clear indication of the fact that simplifying hypotheses typically enforced to estimate the dependence of mixing time on system parameters might be inappropriate in the present context. The results discussed in this work suggest that the flow structure associated with the cross-junction micromixer might control the nanoprecipitation process in the mixing channel downstream the device, thus motivating further quantitative analysis of steady-state mixing in this segment of the flash nanoprecipitator. Because of the large values of the Peclet number that characterize typical flow conditions, this analysis could be based on the assumption that axial diffusion in the discharge channel might be neglected compared to 857 streetwise convection, thus making it possible to use the generalized two-dimensional spectral approach developed in Garofalo et al. (2010), which has been successfully applied up to Peclet values of order 104 -105. References Bally F., Garg D. K., Serra C.A., Hoarau Y., Anton N., Brocon C., Parida D., Vandamme T., Hadziioannou G., 2012, Improved size-tunable preparation of nanoparticles by microfluidic nanoprecipitation, Polymer 53, 5045-5051. Capretto L., Cheng W., Carugo, D., Katsamenis O.L., Hill M., Zhang X., 2012, Mechanisms of co- nanoprecipitation of organic actives and block copolymer in a microfluidic environment, Nanotechnology 23, 375602. Chronopoulou L., Cutonilli A., Cametti C., Dentini M., Palocci C., 2012, Chitosan-Coated PLGA-based nanoparticles: a sustained drug-release strategy for cell cultures, Colloids Surf. B, 103, 117-123. Chronopoulou L., Sparago C., Palocci C., 2014, A modular microfluidic platform for the synthesis of biopolymeric nanoparticles entrapping organic actives, J. Nanopart. Res., 16, 2073-9. Di Pasquale N., Marchisio D.L., Barresi A.A., 2012, Model validation for precipitation in solvent-displacement processes, 2012, Chem. Eng. Sci., 84, 671-683. Di Pasquale N., Marchisio D.L., Carbone P., Barresi A.A., 2013, Identification of nucleation rate parameters with MD and validation of the CFD model for polymer particle nano precipitation, Chem. Eng. Res. Des., 91, 2275-2290. Garofalo F., Adrover A., Cerbelli S., Giona M., 2010, Spectral characterization of static mixers. The S-shaped micromixer as a case study, AIChE J., 56, 318-335. Saad W.S., Prud’homme R.K., 2016, Principles of nanoparticle formation by flash nanoprecipitation, Nano Today 11, 212-227. 858