Automatic Design of Microfluidic Gradient Generators

Concentration gradient generators are common components of microfluidic devices when specified mixing ratios have to be acquired. However, depending on their complexity, the creation of a suitable design is an elaborate task. Although approaches exist that derive such a design, they suffer from low flexibility and usability, which is why most concentration gradient generators are still mainly designed manually thus far—a time-consuming and error-prone task. To tackle this problem, we propose a method for the automatic design of tree-shaped concentration gradient generators and introduce a publicly available online tool that automatically generates a suitable design by adjusting the hydrodynamic resistances within the network. Moreover, the proposed tool is not only able to generate designs complying with user-defined requirements, such as mixing ratios, flow rates, and channel widths but also works for any number of outlets. Simulations using a Computational Fluid Dynamics (CFD) tool as well as evaluations on fabricated designs confirm the quality of the results obtained by the proposed tool.


I. INTRODUCTION
M ICROFLUIDIC concentration gradient generators are important components for biological and chemical research [1,2,3,4,5,6] due to their ability to create spatially resolved concentration values allowing to screen a wide range of concentrations at the same time. Especially tree-shaped gradient generators, comprised of branching, mixing, and recombining channels [7,8,9,10,11,12,13,14,15,16,17,18], are commonly used due to their flexibility in concentration values and their ability to maintain the gradient profile indefinitely [6].
However, their design methods are not yet fully developed. Compared to the semiconductor industry, which heavily relies on automated methods and tools to guarantee a high quality standard of the designed devices, design automation is not commonly used for microfluidic devices [19]. Even though design rules exist for certain cases [13,20,21], there is a lack of available design tools. This is especially critical for gradient generators discussed in this work, since their design requires the consideration of a high number of parameters (such as channel geometries, concentration values, fluid properties, flow rates, etc.), which all affect each other and, hence, make determining a design realizing the desired concentration ratios a highly non-trivial task. Although there are first tools that allow to derive such a corresponding design (discussed later in Sec. II-B), they often lack in flexibility and usability. As a consequence, designing gradient generators is mainly done by hand thus far-resulting in a very time-consuming task where most likely errors are introduced due to the manual approach.
In this work, we address this problem by introducing a publicly available online tool which is capable of automatically designing the layout of concentration gradient generators that satisfy user-defined requirements. To this end, we utilize the tree-shape structure as well as the physical model proposed by Oh et al. [13] to properly balance hydrodynamic resistances of the corresponding channels in order to match the user-defined concentration values as well as flow rates at the outlets. The resulting tool eventually is capable of designing the desired devices in a push-button fashion. Possible applications of this automated design approach lie especially in medical, biological, and chemical research. Changing requirements on the concentration gradient can be met through simple adaption of the input parameters, which effectively allows to shorten turnaround times. The tool is accessible at https://iic.jku.at/eda/research/gradient_generator.
In order to validate the quality of the results generated by the proposed tool, several designs were generated by choosing random values for the concentration values at the outlets. Then, these designs were simulated with a Computational Fluid Dynamics (CFD) tool, which showed that the simulated concentration values at the outlets differ (on average) by no more than 0.4% from the actually desired ones. Moreover, selected designs were fabricated and experimentally evaluated through optical images. Although physical realizations usually introduce additional sources of errors (due to fabrication tolerances, inconsistencies in the measurement setup, etc.), the obtained results also confirm the validity of the proposed tool-showing deviations of the concentration values by no more than 4.24%.
The remainder of this work is structured as follows: Section II first reviews the considered problem as well as the related work. Motivated by that, we describe the proposed solution in Section III. The quality of the results generated by the proposed tool has been validated by simulations and experiments on fabricated designs described in Section IV. Finally, Section V concludes this work.

II. BACKGROUND & MOTIVATION
In this work, we consider the (automatic) generation of microfluidic gradient generators. To this end, we first briefly review the corresponding design task followed by a discussion of the related work. This shows that the state-of-the-art in the design of microfluidic devices generating different mixtures still leaves a substantial amount for improvement-motivating our work.

A. CONSIDERED PROBLEM
Gradient generators are devices that generate particular concentrations of fluids, which are usually predefined by the user. In this work, we focus on concentration generators that rely on diffusive mixing, i.e., which work in regimes of low Peclét numbers [20]. More precisely, the basic concept behind a microfluidic gradient generator is shown in Fig. 1. Here, two fluids are pumped into the microfluidic structure through the inlets on the top-left and top-right and, then, flow towards a single outlet channel where they get finally mixed. In this case, the resulting concentration value at the outlet C o depends on the concentration values of the two fluids (C i1 and C i2 ) as well as on their flow rates (Q i1 and Q i2 ), i.e., Since the concentration values of the two fluids are usually predefined, C o can be specified by properly adjusting the flow rates of the two streams. For example, setting Q i1 = Q i2 would result in a concentration value that is the average between the two inlet concentration values. This simple concept can be expanded to support more sophisticated designs, leading to complex gradient generators that are able to output multiple desired concentration values. Especially tree-shaped gradient generators, are an ideal platform for that [22]. They are characterized by a closed microfluidic network consisting of two inlets at their base that are conjoined through a tree-shaped channel system with multiple outlets at the bottom-level 1 . Its branches are composed of meander channels that ensure the diffusive mixing of the liquids. Additionally, altering the lengths of these meanders affect the flow rates inside the device and, by this, desired concentration values at the outlets can be realized [23]. As an example, Fig. 2 shows such a tree-shaped design of a gradient generator which yields very dedicated concentrations of 100%, 80%, 40%, and 0% in the corresponding outlets.
However, designing such gradient generators is highly non-trivial and constitutes a rather complex task which depends on numerous parameters. This obviously includes the used liquids and the desired concentration values. But, beyond that, also the correct hydrodynamic resistances and flow rates inside the device, the time span for a successful diffusive mixing, as well as various geometrical constraints have to be considered when designing such a gradient generator that is supposed to realize the desired output.

B. RELATED WORK
In the following, we discuss the related work and their corresponding designs and/or corresponding approaches which have been presented in the past in order to realize gradient generators.
The concentration gradient generator design proposed by Lee et al. [23] is highly flexible in the sense that it can generate various concentration values at the outlets, ranging from simple linear gradients to arbitrary monotonically rising functions. This is done not only through variation of the fluidic resistances but also by allowing different input flow rates at the inlets. However, the used design is rather inefficient with respect to the occupied space, since the required meander channels are typically very long caused by the staggered design. Furthermore, while this work shows how to derive the desired channel dimensions, an automatic process to generate the desired design is not presented. Hence, a user is still forced to create the design manually, which can become very tedious considering the long meander channels.
One recent work introduces a concentration gradient generator design tool that is based on the same circuit analysis model [13] that we use in our work [24]. There, a tool was presented that is able to predict the generated concentrations for tree-shaped gradient generators with two different inlet flow rates yet uniform meander lengths. Further, the proposed approach allows to calculate the ideal length of the meanders for ensuring a full mixing based on specific user inputs. Nevertheless, the tool cannot be used for the creation of custom gradient generator designs. Due to the limitation of uniform meander lengths, only specific concentration gradients can be reached which, of course, limits the use-cases of this approach. Thus, the tool does not allow to specify the desired concentration values rather than just outputting them.
The approach proposed by Wang et al. [25] is unconventional and works with randomly generated designs. More precisely, a quadratic grid network of size n×n with two inlets at the top and multiple outlets at the bottom is created. Between each grid point, random connections are drawn which results in a completely random microfluidic network. Then, many of these random networks were simulated (using Computational Fluid Dynamics, i.e., CFD) and the resulting concentration values at the outlets were stored inside a database. Based on that, users can then input their desired concentration values and the database is scanned for the closest match. In order words, this solution is basically resting on randomly generated designs. The obvious limitation of this approach is that it is strongly relying on having a good database that has a closely fitting design. If such a design is not available in the database, only a rough approximation (if at all) is obtained. Additionally, any change in the number of outputs or simply changing to a liquid with different diffusion coefficient requires an extensive extension of the database.
All these related works show that there are approaches how gradient generators can be realized. But they also unveil significant disadvantages and limitations such as low flexibility and usability. As a consequence, most gradient generators are still mainly designed manually thus far-a time-consuming and error-prone task. In this work, we are aiming to address this issue by proposing a tool which is able to automatically derive a design of a gradient generator realizing the desired concentration values based on specifically given user inputs. We decided to focus on tree-shaped concentration gradient generators rather than the design introduced by Lee et. al. [23], since the former produce a more balanced layout with equally long channels for linear gradients, compared to the latter which are more suited for logarithmic gradients. By additionally making this tool easily accessible through a web-based frontend, we believe this allows to overcome the current unsatisfactory situation in which gradient generators are still designed manually.

III. PROPOSED SOLUTION
In this section, we present the proposed method which allows to fully automate the construction of a gradient generator and, thus, aids designers in this, usually, time-consuming and error-prone task. To this end, we first discuss the general idea behind the approach and, afterwards, take a closer look into the implementation details. The method eventually got implemented as an online-tool which is described at the end of this section-including a discussion of the various parameters that can be used to automatically generate the desired design for a specific use-case.

A. GENERAL IDEA
Following the established procedure, also the proposed approach utilizes a tree-shaped structure as shown in Fig. 3, which has two inlets and multiple outlets (in this case, five). At the two inlets, user-defined fluids with certain concentration values are injected into the structure. The two fluids are then mixed multiple times inside the structure while they flow towards the outlets of the device and, by this, realize the desired concentration values. Of course, the number of outlets as well as the corresponding concentration values can be specified by the user.  In order to realize the desired concentration values, the design consists of several cascaded so-called "layers" (cf. Fig. 3). The purpose of each layer is to guide the incoming fluids to corresponding meander channels, where two streams with different concentration values get mixed and form a fluid with a new concentration value. For example, the third outlet of the second layer in Fig. 3 is a mixture between the layer's second and third inlet (or in other words, the second and third outlet of the first layer). As a result, each layer (except the last one) creates intermediate results of concentration values inside the corresponding meanders. The last layer finally mixes these intermediate results one more time and realizes the desired concentration values at the outlets. Additionally, the number of outlets N O defines the number of layers by N L = N O − 2. This is due to the tree-shaped structure of the gradient generator and since the number of outlets is defined by the user, the number of layers is also determined by this. Please note, the higher the number of layers, the larger the structure becomes and, by this, the longer the paths between the inlets and the outlets get. As a result, the total "execution time" increases, i.e., the time interval between the injection of the two fluids and the fully establishment of the desired gradient at the outlets. Similar to Fig. 1, the exact fluid combination (and, by this, the resulting concentration value) inside a meander channel then depends on the flow rates and the corresponding concentration values of the two streams that meet inside the meander. Hence, by properly adjusting these flow rates, the desired concentration value can be realized.
That is, in contrast to tree-shaped structures with uniform meander lengths, we propose to allow for complex concentration profiles even with equal flow rates at the inputs. This can be accomplished by altering the lengths of the meander channels, since the length of a channel is directly proportional to its hydrodynamic resistance and, thus, also affects the flow rates inside the device. Therefore, the meander channels are not only used to ensure a proper diffusive mixing but are also utilized to tune the flow rates inside the structure. Using this as basis, the remaining question to address is how to dimension the channels in order to establish the required flow rates inside the gradient generator and, by this, how to realize the desired concentration values at the outlets.

B. IMPLEMENTATION
In order to properly dimension the design of the gradient generator so that the desired concentration values are realized, a proper physical model is needed that allows to describe the behavior of such a microfluidic device. To this end, we utilize the model proposed by Oh et al. [13], which is applicable at low Reynolds numbers where a laminar flow occurs. As this is the case for the considered gradient generators, this perfectly suits our needs. In fact, the model allows to describe the flow inside a microfluidic channel by Hagen-Poiseuille's [13] law ∆P = Q R, where ∆P is the pressure drop along the channel, Q the volumetric flow rate, and R the hydrodynamic resistance of the channel. The resistance depends on the length l, width w, and height h of the channel as well as the dynamic viscosity of the fluid µ and can be compute as follows (for rectangular channels with a ratio of h w ≤ 1): Hence, the model basically allows to transform hydrodynamic components into their electrical counterparts, i.e., channels become resistances, pressure pumps become voltage sources, and flow rate pumps become current sources. This allows to convert a microfluidic design to an equivalent electrical network with all the benefits that come along with it (e.g., the possibility to apply well-known methods and means from the electrical-domain) 2 . Having such a model, now allows to correctly dimension all the channels inside the gradient generator so that the desired functionality is obtained. More precisely, we break this down into four simpler steps: 1) Defining basic parameters 2) Computing the flow rates in each channel 3) Computing the hydrodynamic resistances of each channel 4) Solving geometrical constraints to generate the design First, we have to define basic parameters that must be provided by the designer and are needed for our implementation. One of the most important parameters here is the number of outlets, because it defines the number of layers and, by this, the size of the underlying "tree-shaped" structure of the equivalent electrical network (i.e., how all the channels are connected). Other important parameters are the general width and height of the channels, the flow rates, as well as the concentration values at the inlets/outlets, and the viscosity of the used fluid. Besides this, there are also parameters which modify the "appearance" of the final design. While all these parameters are discussed in more detail in the next section, we assume in the following that they are already specified.
The second step is to compute the flow rates inside each channel. While certain flow rates are already known (the ones at the inlets and outlets, since they are defined by the user), we only have to obtain the flow rates inside each layer/meander. This is achieved by combining Kirchhoff's current law with the equation for the mixed concentration value (i.e., Eq. (1)) at all relevant points. Note that this is possible without knowing any resistance value of a channel inside the network.
Once all flow rates inside the channels are obtained, the next step is to compute the required resistances of each channel so that these flow rates are established inside the device. Again, only certain resistances have to be computed (namely the ones from the meander channels), since all other resistances are indirectly defined by the user through geometrical and fluid parameters. Here Kirchhoff's current and voltage laws can be applied in order to establish and solve an equation system that allows to obtain these values. Moreover, at this step also a minimal mixing time is considered (specified by the user), to ensure a proper mixing between two streams inside a meander channel. This time basically affects the minimal lengths (or resistances) for each meander.
Having all resistances inside the network, the length of each channel is now computed by rearranging Eq. (2) to

FIGURE 4: Input mask
Basically this only has to be done for the meander channels, since the lengths of all other channels are already specified by user defined parameters. The last step is now to solve the geometrical constraints of all meander channels. More precisely, the number of arcs as well as the width and height of the area a meander channel occupies must be chosen in such a way, that each channel matches the desired length. Once this is accomplished, the design of the gradient generator is fully specified and is ready for further steps such as fabrication.

C. RESULTING TOOL
The method described above has been implemented as an online-tool that is easily accessible at https://iic.jku.at/eda/ research/gradient_generator and can be used by everyone to automatically create the design of a gradient generator realizing the desired concentration ratios. To this end, the user only has to provide the basic parameters into an input mask as shown in Fig. 4.
More precisely, the user has to provide the following parameters:

Geometrical Parameters
• w: General width of all channels. This parameter is one of the most important ones and the user is encouraged to choose a suitable value here by considering different aspects such as the used platform, fabrication limits, required volume of the mixed liquids, etc. Moreover, this value is also the basis for certain constraints (described in the following parameters) that ensure a proper design VOLUME 4, 2016 of the gradient generator. Of course, this value always has to satisfy w > 0. • h: General height of all channels. Similar to the width the user should also consider aspects such as fabrication limits or required volume of the mixed liquids when specifying this value. The inequality w ≥ h > 0 must hold, since Eq. (2) is only valid for these values. • r: Radius of the meander arcs. By increasing this value, the arcs become larger and the meander obviously needs more space to get the same channel length as a meander with a smaller radius. Hence, this parameter primarily affects the appearance of the final gradient generator. In order to prevent too small gaps between the single arcs, this value is also constrained to r ≥ w. • w Meander,max : Maximal width of the area a meander can occupy. This value is also responsible for the distance between two neighbouring inlets/outlets and, thus, allows to control the overall width of the whole gradient generator (i.e., it mainly affects the appearance of the final design). In order to ensure that a meander channel can realize basically every desired hydrodynamic resistance, this value must satisfy w Meander,max ≥ 5w + 8r.

Fluid Parameters
• µ: Viscosity of the used fluids needed for Eq. (2). Of course, the value has to satisfy µ > 0. • t min : Minimum time two streams need to flow inside a meander, in order to ensure a successful mixing. Basically, this value affects the length of the meanders and must be t min > 0. Inlets In , c (2) In : Concentration values at the two inlets that have to satisfy 100 ≥ c (1) In > c (2) In ≥ 0. Per default, we assume that c (1) In = 100% and c (2) In : Flow rates at the two inlets, which have to satisfy Q (1) In > 0 and Q (2) In > 0 • l (1) In , l (2) In : Length of the inlets. For design reasons we constrain them to l (1) In ≥ w and l (2) In ≥ w. Outlets • N : Number of outlets, which always has to satisfy N ≥ 3, due to the tree-shaped structure of the gradient generator. • c (2) Out , c   Out are already determined due to physical dependencies. Per default, we assume that all these parameters have the same value and equal a fraction of In )/N . Additionally the constraint Out : Additional hydrodynamic resistances to be considered at the outlets. These parameters are completely optional and allow to consider devices/channels that are attached to the outlets of the gradient generator, but are not included in the final design. Also have to satisfy R (i) After all the parameters are properly defined, a design of a gradient generator can be created by a simple click on a button (using the general idea and the implementation as described above). The resulting design can then be exported as a Scalable Vector Graphics (SVG) file as shown in Fig. 5. Using this file, further steps (such as fabricating the design) can be conducted. Overall, this completely turns a tedious and error-prone manual design step into a fully automatic process in which only basic parameters have to be provided.

IV. QUALITY OF THE RESULTS
Obviously, the quality of the tool proposed in the previous section is only as good as the quality of the designs it generates. Because of this, we additionally conducted several evaluations-both, through simulation as well as physical validation-which allowed for an assessment of this quality. This section summarizes the corresponding endeavors as well as the obtained results.

A. VALIDATION THROUGH SIMULATION
In order to validate the designs generated by the proposed tool, we utilized OpenFOAM [28]-an established Computational Fluid Dynamics (CFD) simulation tool. This tool allows to precisely simulate the flow distributions inside the channels as well as the mixing process between the two fluids and, thus, allows to validate the concentration values at the outlets.
To this end, we first generated multiple designs with 3, 4, 5, 6, and 7 outlets, where each design had to generate (randomly defined) concentration values at the outlets. More precisely, for each distinct number of outlets, we generate 20 designs (100 designs in total), where each design is automatically created by the proposed tool. Afterwards, the designs got converted and meshed so they can be used and simulated with OpenFOAM.
After the simulations have been completed, we compared the desired concentration values at the outlets with the results obtained by the CFD simulations. This comparison is shown in Table 1. Here, each row summarizes the results of 20 simulations (N S ) with the corresponding number of outlets N O . The third column |e| shows the mean absolute error (MAE) of the simulated concentration values compared with the desired ones, while the fourth and fifth column show the minimal and maximal absolute error, respectively. More precisely, the MAE can be calculated as follows, where i S and i O are the indices of the simulation run and outlet, respectively. Please note that we did not consider the concentration values of the first and last outlet since they are identically with the first and second inlet, respectively, and, thus, are not relevant in practical use cases anyway. Due to that, the index i O only runs from 2 to N O − 1 and not from 1 to N O , since this would distort the results. As it can be observed, the simulated concentration values differ on average by not more than 0.4% from the desired ones (regardless of the number of outlets). Even the worst cases (i.e., the maximal errors) remain rather small and never exceed 0.9 %-a negligible error. This clearly confirms the quality of the obtained results and, hence, the validity of the proposed design tool.

B. VALIDATION THROUGH FABRICATION & MEASUREMENT
In addition to the evaluations based on CFD simulation, the quality of the results obtained by the proposed tool has also been evaluated using fabricated devices. To this end, designs with three and five outlets were generated by the proposed tool and, afterwards, fabricated. By measuring the outcomes of the resulting devices, the accuracy of the results has been determined. In the following, we briefly summarize the fabrication and measurement process, before we provide and discuss the obtained results.

1) Chip Fabrication & Concentration Measurement
The microfluidic devices were fabricated from poly(methyl methacrylate) (PMMA) and patterned using a commercial computer numerical control (CNC) milling machine with a 300 µm square endmill [29]. The channel widths and heights were predefined as 300 µm and 200 µm, respectively, and used as inputs to the proposed design tool. The concentrations at the outlets of the fabricated devices were determined using a colorimetric method, which is based on the evaluation of mixing of different colored liquids in the HSV (hue, saturation, value) color space. This method was introduced by Rezk et al. [30] and later also successfully utilized by other groups [31,24]. To perform the measurements, two differently colored liquids were prepared by diluting yellow and blue food coloring (Lianyungang Xinai Food Technology Co., Guannan County, China) with deionized (DI) water. The test liquids were prepared with a 8:2 and 2:8 ratio of yellow:blue dye, and represent 100% and 0% concentration at the left and right inlet of the gradient generators, respectively. From that, a monotonically decreasing hue value can be expected with increasing concentration at the outlets. Optical images were obtained using a digital SLR camera (EOS 760D, Canon, Tokyo, Japan) with a macro lens (EFS 60 mm, Canon, Tokyo, Japan).
In order to be able to relate the outlet hue values to the concentration, first a calibration curve was established. To this end, handmade mixtures of the two liquids were prepared and filled into dedicated channels on the chips. Fig. 6 shows the hue values of the handmade mixtures and a corresponding quadratic fit. Note that the colors for the 0% and 100% concentrations were specifically chosen to be able to fit a simple quadratic function.
Afterwards, the two colored test liquids were injected at the inlets using syringe pumps (ExiGo Pump, Cellix, Dublin, Irland) at a flow rate of 5 µL min −1 (also used as inputs to the proposed design tool). Again, optical images were taken and the hue value of each outlet was determined in the same fashion as during calibration. Using the calibration   curve from before, the measured hue was mapped to the corresponding concentration. 8 show a snapshot of the fabricated chips and the measurement results compared to the designed outlet concentrations for the 1-layer and 3-layer design, respectively. The measured concentrations are in good agreement with the desired values that the proposed design tool was supposed to generate. In fact, the maximum deviation of the fabricated devices from the desired value is 4.24% in the worst case. Again, this is a very good accuracy-particularly considering that many deviations may not be caused by the design, but can be attributed to tolerances in the fabrication process (since the width and height of the channels play a significant role in the overall accuracy of the gradient generator) and the photo-setup (which may introduce small errors through inconsistencies between individual pictures). Overall, also the evaluation using fabricated devices and measurement confirm the quality of the results obtained by the proposed approach.

V. CONCLUSIONS
In this work, we considered the automatic design of microfluidic gradient generators. We introduced an online tool that balances the hydrodynamic resistances such that complex concentration profiles at the outlets can be achieved by means of varying the meander lengths. Using this, designers only have to provide basic input parameters and can then obtain a functional design in a push-button fashion. This greatly enhances the flexibility and usability compared to previous approaches, which are rather limited in this regard.
The quality of the results generated by the proposed tool has been validated through CFD simulations confirming that the absolute differences between the desired concentration values and simulated results is negligible (less than 0.9% for all tested designs with an average error of 0.4%). Even experiments with fabricated designs (which are additionally subject to further deviations caused by fabrication processes) confirm the quality of the obtained results.
Even though this work presents an easy tool for the automatic design of concentration gradient generators, it is still the user's responsibility to check physical limitations related to the fabrication itself or the applicability of the underlying simplifications. Whereas a suitable user interface can aid the designer in the latter aspect, the limits of the fabrication cannot be easily considered in a design tool, since they are highly dissimilar and dependent on the method of choice.
Future work can focus not only on further increasing usability by highlighting restrictions in the user interface, but also on implementing different layout options for the user to choose from. Whereas the presented tool focusses on tree-shaped designs, alternative design options that provide a more optimal result in terms of channel length, fluid volume, or footprint can be implemented.
While concentration gradient generators served as proof-of-principle, one can imagine to use the same approach for the design of other microfluidic components such as droplet generators, mixers, or filters. A successive combination of individual tools developed for various components might then even allow the automated design of more complex composite microfluidic chips. Finally, with increasing complexity, the value of automatic design tools for microfluidics is becoming more and more evident.
The resulting tool is available online and can be accessed on https://iic.jku.at/eda/research/gradient_generator. In these areas, he published more than 350 papers in journals and conferences and served in editorial boards and program committees of numerous journals/conferences such as TCAD, ASP-DAC, DAC, DATE, and ICCAD. For his research, he was awarded, e.g., with an ERC Consolidator Grant, Best Paper Awards, e.g., at TCAD and ICCAD, a DAC Under-40 Innovator Award, a Google Research Award, and more. VOLUME 4, 2016