Research Article  Open Access
Raudel Avila, Yixin Wu, Rinaldo Garziera, John A. Rogers, Yonggang Huang, "Analytical Modeling of Flowrate and Its Maxima in Electrochemical Bioelectronics with Drug Delivery Capabilities", Research, vol. 2022, Article ID 9805932, 13 pages, 2022. https://doi.org/10.34133/2022/9805932
Analytical Modeling of Flowrate and Its Maxima in Electrochemical Bioelectronics with Drug Delivery Capabilities
Abstract
Flowrate control in flexible bioelectronics with targeted drug delivery capabilities is essential to ensure timely and safe delivery. For neuroscience and pharmacogenetics studies in small animals, these flexible bioelectronic systems can be tailored to deliver small drug volumes on a controlled fashion without damaging surrounding tissues from stresses induced by excessively high flowrates. The drug delivery process is realized by an electrochemical reaction that pressurizes the internal bioelectronic chambers to deform a flexible polymer membrane that pumps the drug through a network of microchannels implanted in the small animal. The flowrate temporal profile and global maximum are governed and can be modeled by the ideal gas law. Here, we obtain an analytical solution that groups the relevant mechanical, fluidic, environmental, and electrochemical terms involved in the drug delivery process into a set of three nondimensional parameters. The unique combinations of these three nondimensional parameters (related to the initial pressure, initial gas volume, and microfluidic resistance) can be used to model the flowrate and scale up the flexible bioelectronic design for experiments in medium and large animal models. The analytical solution is divided into (1) a fast variable that controls the maximum flowrate and (2) a slow variable that models the temporal profile. Together, the two variables detail the complete drug delivery process and control using the three nondimensional parameters. Comparison of the analytical model with alternative numerical models shows excellent agreement and validates the analytic modeling approach. These findings serve as a theoretical framework to design and optimize future flexible bioelectronic systems used in biomedical research, or related medical fields, and analytically control the flowrate and its global maximum for successful drug delivery.
1. Introduction
Controlled and targeted drug delivery of pharmacological agents in organs/tissues has helped researchers study local biological responses to drug treatments and determine the drug efficacy while mitigating unwanted side effects often present in systemic drug delivery strategies [1, 2]. The approaches and technologies in drug delivery vary significantly depending on factors related to the type of drug and anatomical target location. A review of drug delivery technologies presented in [3, 4] discusses the evolution of drug delivery approaches that range from tailored nanoparticles that drive the active drug agents through biological membranes to implantable microsystem with injectable probes that enable targeted drug delivery in organs/tissues of interest.
Examples of injectable drug delivery systems include intracerebral [5, 6] and intraarterial [7] injections in the brain, kidneys, eyes, ears, and lymphatic system of small (mice) and mediumsized (cats) animals that are aimed at enhancing the drug effectiveness when treating affected areas. For a timely and successful drug delivery, the flowrate must be accurately controlled as excessively high flowrates can exert dangerous stresses on fragile surrounding tissue [8, 9] and excessively low flowrates can obstruct the fluidic probes and result in an incomplete/ineffective delivery [10]. Table 1 shows a list of reported flowrate ranges used in drug delivery applications that vary depending on the size of the animal and target location. In reported injectable drug delivery systems for small animals, the drug flowrate can vary in the range of tens of nanoliters per minute to hundreds of microliters per minute depending on the drug dose and intended application such as behavioral neuroscience [10, 11] or cancer therapeutics [12–14]. For this reason, accurately controlling the flowrate and its maximum value during the drug delivery process is important to ensure safe delivery of the drug without imposing stresses that can damage the surrounding soft tissues in small animals.

The rapid development of microtechnology in biomedical research has enabled wireless flexible bioelectronics with drug delivery capabilities in miniaturized form factors that integrate chemical, mechanical, and fluidic interfaces to manipulate small fluid volumes in a programmable fashion and realize targeted drug delivery in regions like the brain [7, 9, 10, 15, 16], peripheral nerves [11], and eyes [17]. To circumvent previous tethered delivery strategies that restrict animal motion, flexible bioelectronics utilize a wireless link that harvests energy from nearby electromagnetic sources and powers the corresponding subsystems to deliver the drug in freely moving animals while maintaining negligible device and thermal loads that can affect the animal behavior or alter the chemistry of the drug. For neuroscience research specifically, this aspect is especially important because minimizing the device load on freely moving animal allows to study the drug effects without any interference from the drug delivery actuation mechanisms [18–20].
Figure 1 shows a crosssectional schematic of a typical flexible bioelectronic device used for drug delivery where all the geometric, environmental, fluidic, and mechanic parameters are labeled. The main components of the flexible bioelectronic include a set of interdigitated electrodes, an electrolyte chamber, a flexible membrane, a drug chamber, and microfluidic channels partially implanted in organs/tissues. To deliver the drug, electrical current flows through the electrodes—in direct contact with the electrolyte—to initiate an electrochemical reaction known as water hydrolysis. The byproducts of this electrochemical reaction are hydrogen and oxygen gas that accumulates in the electrolyte chamber and pressurizes the bottom side of a flexible membrane. In response to the accumulated gas pressure, the flexible membrane deforms into the form of a spherical cap and gradually pumps the drug, sitting on the top side of the flexible membrane, through a network of microfluidic channels which outlet is the target organ/tissue. By using electrochemistry, the drug can be pumped from inside the devices in a controlled fashion without generating excessive heat [11, 21, 22] or requiring external moving parts that can increase the complexity of injectable devices operating partially inside the tissue/organs of animals.
Several analytic models that combine more than 12 parameters involved in the drug delivery process related to geometric, environmental, fluidic, electrical, and flexible membrane mechanics into unique combinations of 3 nondimensional parameters related to the initial environmental pressure, initial gas volume, and the microfluidic resistance have been proposed [23]. These 3 nondimensional parameters group all the parameters involved in the drug delivery process and must be carefully selected to control the total delivery time and volume for accelerated drug delivery [24]. However, prior analytical models [23–25] focus on the total volume delivered over time for applications requiring accelerated delivery (e.g., lifesaving medication) and do not satisfy the zeroflowrate initial condition which affects the flowrate over time, particularly the maximum flowrate. For injectable drug delivery systems targeting fragile tissues/organs in small animals, the maximum flowrate always occurs near the beginning of the drug delivery process and is the relevant quantity to control and adhere with reported experimental flowrate guidelines for administering certain types of drugs in experiments like the ones listed in Table 1 and to ensure that the drug is delivered safely without inducing excessive stresses that can damage the surrounding tissue/organs.
Here, we propose an analytical model that separates the drug delivery process into a “slow” variable that satisfies the zerovolume, but not zeroflowrate initial conditions, as used previously in [24, 25] to control the drug delivery time and volume, and a new “fast” variable that corrects the slow variable solution to satisfy both the zerovolume and zeroflowrate initial conditions and is used to control the flowrate and its maximum value without affecting the drug delivery time and volume. The new analytical model featuring a “fast” variable for flowrate control and an explicit formula for the maximum flowrate is validated by showing excellent agreement with the numerical results without using the “slow” and “fast” variables and serves as the theoretical framework to design the flexible bioelectronic devices with flowrate control capabilities that can be adjusted by unique combination of the nondimensional parameters depending on the target organ and drug delivery timeframe.
2. Results
2.1. Flexible Membrane Mechanics
Pumping the drug from inside the flexible device into the target location is achieved by inducing a pressure differential between the bottom and top surfaces of a flexible membrane which causes it to deform into the shape of the drug reservoir (e.g., spherical cap) with a maximum vertical displacement . This deformation is governed by the flexible membrane geometrical (i.e., thickness and radius ) and mechanical properties (i.e., Young’s modulus and Poisson ratio or alternatively a generalized hyperelastic strain energy density function), and the mechanics of deformation can be modeled according to the function that gives the relationship between pressure differential applied to the flexible membrane and the volume it expands.
During the delivery process, the flexible membrane transitions from bendingdominated (i.e., small displacement ) to stretchingdominated deformation (i.e., large displacement ) as it goes from flat into the shape of a spherical cap. Based on experiments for a bioelectronic device that targets the peripheral nerves [11] and numerical models of the drug delivery process [10, 23], the maximum flowrate always occurs near the beginning of the delivery process when the deformation is small, i.e., bendingdominated.
For these bioelectronics with drug delivery capabilities, block copolymers like styreneisoprenestyrene (SIS) have been used as flexible membrane for their soft mechanics and gas permeability properties which allow larger deformations and prevent the gas migration from the electrolyte chamber to the drug chamber [10]. Ideally, the stressstrain relationship of the block copolymers should be obtained from experimental testing (e.g., uniaxial and biaxial), when possible, to accurately model the flexible membrane expansion [32, 33]. Figure 2(a) shows the constitutive stressstrain relationship, obtained from uniaxial testing, of a SIS block copolymer sample. Young’s modulus is calculated from the small strain (linearelastic regime) as . For a representative bioelectronic device (i.e., thickness and radius ) previously used for drug delivery in the mouse brain [10], the function is obtained numerically using finite element analysis (FEA) as shown in Figure 2(b), which considers both bending and stretching effects and uses the Marlow hyperelastic model [33] to build the strain energy density function of the SIS block copolymer based on the stressstrain data shown in Figure 2(a). When the deformation of the flexible membrane is small (i.e., ) such as near the beginning of the drug delivery process when the maximum flowrate occurs, a bendingdominated function is derived from plate theory [34] as which establishes a linear relationship between pressure differential and the deformed volume as shown in our previous work [24, 25] and has excellent agreement with the function obtained from FEA as shown in Figure 2(b) when the deformation of the flexible membrane is small. Generally, the function can be written nondimensionally as , where is the nondimensional volume, which gives for bendingdominated deformation only. It gives a constant derivative with respect to , and for an incompressible membrane material (i.e., ), which will be used in subsequent sections to analyze the maximum flowrate.
Figure 2(b) shows that once the deformation of the flexible membrane becomes large or , the bendingdominated solution does not agree well with FEA as the stretching (and nonlinear deformation) effects in the membrane become relevant. However, for the analysis in this paper, the focus is the maximum flowrate which occurs when the deformation is small as explained in the Supplementary Information (Note 1), and thus, bendingdominated, not the volume or flowrate temporal response, requires the numerical and analytical function to agree over the entire pressurevolume range as shown in our previous work on stretchingdominated deformation [24, 25]. This validates the analytical bendingdominated function in Equation (1) instead of the FEA solution to control only the maximum flowrate in the drug delivery process.
2.2. Governing Equations of Drug Delivery
A typical wireless bioelectronic device features three main components: (1) electrochemical reservoir, (2) flexible membrane, and (3) microfluidic channels where each contributes to the drug delivery dynamics and affects the delivery time, flowrate, and its maxima and must be accounted for explicitly. The physics of wireless bioelectronics relying on electrochemistry as actuation method can be modeled approximately by the ideal gas law: where is the pressure, is volume change inside the electrolyte reservoir, is the initial volume of gas in the electrolyte reservoir, is the ideal gas constant (8.3144 J mol^{1} K^{1}), is the temperature of the electrolyte, and is the number of moles, which is linear to the electrical current in the electrodes and is given by following Nernst equation [35]. Here, is time, is Faraday’s constant (96485 C mol^{1}), the ratio corresponds to water and must be changed if using a different electrolyte according to the redox reaction [36], is the initial amount of gas moles in the electrolyte reservoir and is related to the initial volume of gas in the electrolyte reservoir [25] by , and is the initial value of , i.e., the initial environmental pressure at the target region (e.g., ~117 kPa when considering blood pressure).
The pressures applied at both sides of the flexible membrane can be combined with the pressure differential in the microfluidic channel given by to derive the pressure from equilibrium in the device as [23] where the terms on the righthand side include the resistance to membrane deformation, the initial environmental pressure at the target region, and the microfluidic resistance that includes the following: is the dynamic viscosity of the drug, is the length of the microchannels, and are the width and height of the rectangular microchannel, and is the drug flowrate. The microfluidic resistance term is given for a rectangular cross section with laminar flow and must be modified if using different cross sections; for instance, in the following analysis, it will be simplified to for a square (i.e., ) cross section. At the beginning of the drug delivery process (i.e., ) before any electrical current flows through the electrodes, the initial value of in Equation (5) is which satisfied the presence of an initial amount of gas moles in the electrolyte reservoir . Substituting the pressure and number of moles into the ideal gas law in Equation (4) yields a 1^{st}order ordinary differential equation (ODE) for the drug volume as
The firstorder ODE in Equation (6) is the governing equation for the drug delivery process, and it can be rewritten nondimensionally by normalizing the volume as , introducing a nondimensional time , and writing the general expression for the function nondimensionally as to yield that involves three nondimensional parameters, namely, the initial environmental pressure , the initial volume , and the microfluidic resistance that combine all the dimensional parameters involved in the delivery process into only three nondimensional parameters that can be studied independently to understand how they influence the drug delivery process. In this case, it is important to note that the nondimensional function in Equation (7) is the nonlinear function obtained from FEA as shown in Figure 2(b) because the governing equation is used to obtain the volume temporal profile. The governing equation in Equation (7) can be solved numerically with the initial condition .
2.3. Analytical Model for Flowrate in Drug Delivery: Slow Variable
In general, the microfluidic resistance is small, as compared to the other two nondimensional parameters and , but not zero (), such that the perturbation method [37] can be used to solve the governing equation in Equation (7) analytically [25] or the drug delivery time as which is written explicitly for the total nondimensional delivery time that it takes to deliver a nondimensional volume of drug . The first two terms on the right side of Equation (8) can be regarded as the time required for the flexible membrane to deform and overcome the external environmental pressure and the last term is as the time for the drug to travel through the microfluidic channels where both of these times occur simultaneously; i.e., as the flexible membrane deforms, it pumps the drug through the microchannels. Let denote the solution of the above equation, where the subscript “slow” is used to denote a function of the regular (slow) time (as opposed to the fast time introduced in the next section). The flowrate of the “slow” variable solution is obtained by taking its derivative with respect to time in Equation (8) as
It is important to note that Equation (9) gives a nonzero initial flowrate that is equal to which does not satisfy the zeroflowrate initial condition .
2.4. Analytical Model for Flowrate in Drug Delivery: Fast Variable
The “slow” variable solution presented in the previous section works well when is small except at the initial delivery time because it does not satisfy the zeroflowrate initial condition. Since appears as the highest order derivative in the “slow” variable solution for the flowrate in Equation (9) and is small, the singular perturbation method can be used to introduce a “fast” variable solution of the form , where the is a “fast” changing variable and the presence of ensures a small effect on the volume temporal profile but a large effect initially in the flowrate to satisfy the zeroflowrate initial condition.
The total drug delivery time in Equation (8) then becomes , where is the “slow” variable and is a “fast” changing variable that is relevant near the initial time of the delivery process to ensure a zero initial flowrate. Therefore, for a finite and a very small , the value of is approximately zero such that , , and is a constant value given from Equation (10) and these assumptions can be used to derive as and the derivation details for the “fast” variable are shown in the Supplementary Information (Note 2). The flowrate term of the “fast” variable solution is obtained from Equation (11) by taking a derivative as
The complete expression for the flowrate is . For time , Equation (12) becomes to satisfy the initial condition of zero initial flowrate.
The following terminology is adopted in this section to distinguish the solutions: (i)The term “numerical” is used for solutions of the ODE in Equation (7) described in Section 2.2 when the nondimensional function is obtained from FEA using the Marlow hyperelastic model based on the stressstrain behavior of the SIS polymer(ii)The term “semianalytical slow” is used for the “slow” variable solutions in Equation (9) described in Section 2.3 for the flowrate defined as when the nondimensional function is obtained from FEA using the Marlow hyperelastic model based on the stressstrain behavior of the SIS polymer(iii)The term “semianalytical slow + fast” is used for the solutions in Equations (9) and (12) described in Section 2.4 for the flowrate defined as when the nondimensional function is obtained from FEA using the Marlow hyperelastic model based on the stressstrain behavior of the SIS polymer(iv)The term “analytical slow + fast” is used for the solutions in Equations (9) and (12) described in Section 2.4 for the flowrate defined as when the nondimensional function is derived from plate theory for bendingdominated deformation in Equation (2)
The results in Figure 3 show the flowrate temporal profile for a representative bioelectronic device with a SIS flexible membrane previously used for combined drug and light delivery in the mouse brain [10] with the parameters listed in Table 2 and show the numerical, semianalytical, and analytical solutions of the flowrate and its maximum value. Figures 3(a) and 3(b) both show the numerical solution and “semianalytical slow + fast” solution, and they begin with an initial zero flowrate and increase until reaching a peak value labeled as the maximum flowrate and then gradually decrease as the drug delivery process continues. The main difference between the “semianalytical slow” and “semianalytical slow + fast” solutions in Figure 3(a) is at the beginning of the delivery process (i.e., ) showing that the “semianalytical slow” solution does not satisfy the zero initial flowrate, but the “semianalytical slow + fast” solution does due to the introduction of the “fast” variable which dominates at the beginning of the delivery process. Both semianalytical solutions in Figure 3(a) closely match the numerical solution after the initial time because the function is obtained from FEA based on the Marlow hyperelastic model and it considers both the bending and stretching effects of the deformation. Thus, just like in the previous volume temporal profile models [24, 25], modeling the flowrate temporal profile requires excellent agreement between the FEA and analytical function . However, since the maximum flowrate always occurs at the beginning of the drug delivery process when the bending effects in the flexible membrane are prevalent, the FEA function can be replaced by the linear function given in Equation (1) for bendingdominated deformation to model the drug delivery process up to the point where the maximum flowrate is reached. For this bioelectronic device geometry specifically (i.e., thickness and radius ), the bending effects in the membrane cannot be neglected when controlling the maximum flowrate. Figure 3(b) shows that the “analytical slow + fast” using the in given in Equation (1) satisfies the zero initial flowrate condition due to the presence of the “fast” variable and has excellent agreement with the numerical solution up to the point of the maximum flowrate which is the key quantity of focus in this analysis. It is important to note that the bendingdominated deformation can only be used up to the time point when the maximum flowrate is reached while the deformation remains small; otherwise, a stretchingdominated or the FEA solution is necessary. The reason why the bendingdominated solution is relevant in this particular case is twofold: (1) the maximum flowrate always occurs near the beginning of the drug delivery process when the deformation is small and therefore bendingdominated and (2) the flexible membrane bending effects depend on the nondimensional ratio which is almost five times higher than bioelectronic devices handling larger drug volumes (e.g., 100–1000 μl) that focus on achieving faster drug delivery , where a larger membrane radius ensures stretchingdominated deformation [25] and the control is on the drug delivery time and volume instead of the flowrate and its maximum value. However, for compact bioelectronics with drug delivery capabilities for use in small animals, the key quantity to control is the magnitude of the maximum flowrate to avoid damaging surrounding fragile tissues resulting from excessively high flowrates, not the volume or temporal profile as shown in our previous work [23–25] where the focus was to obtain the total delivery time and volume.

2.5. Analytical Model for Maximum Flowrate
The flowrate temporal profiles in Figure 4(a) show that the maximum flowrate occurs near the beginning of the drug delivery process. Currently, the “analytical slow + fast” solution for the flowrate is divided into the “slow” and “fast” terms given by , and the exact time when the maximum flowrate occurs can be obtained from , which, however, is difficult to yield an explicit formula for the exact time when the maximum flowrate occurs. When time is small, such as in the beginning of the delivery process (i.e., ), the term can be approximated by the constant because in the “slow” variable solution, the maximum flowrate always occurs near , but in the “fast” variable solution is not zero and must be determined. Then, the approximate time when the maximum flowrate occurs can be rewritten as where the first term in the righthand side of Equation (13) is a constant and the second term is only a function of . This gives an explicit formula for at which the maximum flowrate occurs, and its substitution into the flowrate expression gives an explicit formula for the maximum flowrate as where all the derivation details are shown in the Supplementary Information (Note 3). Since the bending effects are relevant at the beginning of the drug delivery process when the maximum flowrate occurs, the expression for can be derived from Equation (1) as . The explicit nondimensional expression for the maximum flowrate becomes that depends on the ratio (due to the bending effects) and the three nondimensional parameters , , and described in Section 2.3. Typically, in these bioelectronic devices with the representative parameters listed in Table 2, the electrical current can be modulated to control (increase or decrease) the maximum flowrate as shown in Figure 4(a) where the electrical current is changed from 0.10 mA to 1.00 mA which in turn increases the maximum flowrate from 0.60 μl/min to 5.20 μl/min. The corresponding nondimensional parameters for this example are calculated and shown in Figure 4(a) for comparison, where increases from 0.0002 to 0.0018 when the electrical current changes and the other two nondimensional parameters and are fixed since they do not depend on the electrical current.
The “analytical slow + fast” solution in Figure 4(a) shows excellent agreement up to the time when the maximum flowrate is reached; after, the agreement deteriorates for larger current values (e.g., 0.75 and 1.00 mA) due to (1) the magnitude of which increases and (2) the differences between the function obtained from FEA and bendingdominated deformation in Equation (1). Most notably, the magnitude of the maximum flowrate shows excellent agreement between the numerical and “analytical slow + fast,” which is the key quantity to control during the drug delivery process to avoid damaging fragile surrounding tissues. Figure 4(b) shows the excellent agreement of the value for the maximum flowrate obtained from the explicit analytical formula in Equation (15) with the numerical values computed from the peaks in the flowrate temporal profile. When the electrical current is less than 0.50 mA, the agreement between the numerical and explicit formula is excellent. As the electrical current increases to 1.00 mA, the explicit formula overpredicts the maximum flowrate by ~10.8%, which is still a reasonable agreement, which validates the explicit analytical expression in Equation (15) for the maximum flowrate.
2.6. Parametric Study of the Maximum Flowrate
The maximum flowrate in Equation (15), like the scaling law in Equation (7), depends on the three nondimensional parameters , , and . Therefore, understanding how the maximum flowrate scales with each of the three nondimensional parameter is important to design and optimize the bioelectronic device geometry and ensure safe and successful drug delivery. To explore the influence of , the crosssectional area of the microfluidic channel is reduced from 50 μm to 18 μm which in turn increases to 0.027, while the other two nondimensional parameters are fixed to and . This crosssectional reduction is relevant when targeting smaller areas (or cells) within the tissues to ensure that the drug is being delivered only in a specific region. Figure 5(a) shows that the maximum nondimensional flowrate decreases approximately nonlinearly with ; this nonlinearity is clearly captured from the flowrate temporal profile peaks in the “analytical slow + fast” model as shown by the excellent agreement with the numerical solution. The explicit formula in Equation (15) has a linear dependence on as shown in Figure 5(a) and provides an excellent agreement when is small (key assumption when using the perturbation method) and a reasonably well analytical approximation to the nondimensional maximum flowrate when increases over the relevant range of microchannel cross sections.
To study the influence of , which is relevant in the refill ability aspect of the bioelectronic device, the other two nondimensional parameters were fixed to and . Figure 5(b) shows that both the “analytical slow + fast” and the explicit formula agree extremely well with the numerical solution in the range , which corresponds to the electrolyte reservoir being completely full and partially full (i.e., 50%) that introduces the presence of an initial gas volume . Here, the explicit formula slightly overpredicts within ~6% the magnitude of the nondimensional maximum flowrate.
To understand the influence of , the other two nondimensional parameters are fixed to and . Significantly changing the value of is difficult as it depends on the initial environmental pressure of the target organ which might vary only by a few kilopascals in humans. Figure 5(c) shows that the “analytical slow + fast” and the explicit formula agree very well with the numerical solution.
3. Discussion
Wireless drug delivery technologies have attracted significant attention for their ability to precisely deliver small drug volumes in a programmable fashion to different target tissues/organs without restricting the animal’s ability to move. Designing these drug delivery technologies requires careful consideration of the geometrical, electrical, flexural, fluidic, and environmental parameters to optimize the device layout to ensure safe and complete delivery within the required timeframe based on the flowrate and delivery time. Specifically, controlling the maximum flowrate during the drug delivery process can help to reduce any damage to the surrounding fragile tissues caused by stresses generated from excessively high flowrates. The introduction of a “fast” variable analytical model in Section 2.4 and the explicit formula for the maximum flowrate in Section 2.5 provides a theoretical framework to control the maximum flowrate by changing any of the three nondimensional parameters , , and . Future flexible bioelectronic systems with drug delivery capabilities can be designed by explicitly considering the influence of geometric dimensions, electronics, flexible membrane mechanics, and microfluidic geometries. For example, wireless drug delivery devices in neuroscience experiments are aimed at being small and lightweight such as to not influence the behavior of the small animal during experiments. This design goal of achieving a small and functional device can be achieved by carefully studying the unique combinations of the nondimensional parameters to yield the optimal configuration for the bioelectronic device. After the bioelectronic device is fabricated, only the initial gas volume and electric current can be modified. The initial gas volume can be changed during the refilling process between experiments, and the electric current can be modified by operating at different duty cycles. Restrictions in size and geometric dimensions affect the device radius and available electric current, and these effects can be modeled using the nondimensional parameters to identify if the bioelectronic device is able to reach the desired maximum flowrate range during delivery before fabricating the device. Further, the influence of the flexible membrane mechanics can be explicitly considered to use a block copolymer with Young’s modulus that deforms quickly with the applied pressure (i.e., soft) and helps achieve the desired maximum flowrate. Similarly, the influence of the microfluidic channel geometry and dimensions can be considered explicitly such as to not impose excessive fluid resistance that can delay the delivery or cause a blockage inside the device while still targeting drug delivery in specific locations. For example, Figure 5(a) shows that changing the cross section of the microchannel from 50 μm (delivery area 2500 μm^{2}) to 18 μm (delivery area 324 μm^{2}) will decrease the maximum flowrate during delivery which can be important to consider especially when interested in drug delivery to cells with dimensions in the tens of micrometers or less. The separation of the flowrate solution into a “slow” variable (which is relevant to determine the total delivery time) and a “fast” variable (which is relevant to satisfy the zeroflowrate initial condition at the beginning of the drug delivery process) allows to prioritize which dimensional parameters (and consequently nondimensional parameters) to change depending on the quantity to control, e.g., delivery time or flowrate and its maxima, and most importantly how to change them to increase or decrease the maximum flowrate based on the linear dependence on and inverse linear dependence on and shown in Equation (15) which can be used to control the maximum flowrate when changing the microchannel cross section to target a smaller region in the organs, initial volume of gas during refill and reuse process, and the physiological pressure of the target region organ as shown in Figure 5. For example, the proposed analytical model can be used to design the bioelectronic device to comply with the maximum flowrate applications listed in Table 1 and tailor specific experiments, target locations, drug delivery timeframe, and animal size. The physics of the maximum flowrate are presented as follows: the first term in Equation (15) is the maximum flowrate achieved while the flexible membrane is overcoming the external pressure to deform, and the second term (which is a negative term of the nondimensional parameter ) is the delayed effect caused by the drug traveling through the microchannels. The relevance of the “slow + fast” variable analytical model presented in Section 2.4 in the context of timesensitive experiments in freely moving animals (1) provides control over the rates of drug delivery which are important in many behavioral neuroscience studies and the total time to deliver the drug and (2) the “fast” variable model allows to determine the time required to reach the maximum flowrate to enable safe wireless pharmacology experiments. The benefits of employing the “analytical slow + fast” variable over the numerical model and FEA to design this emerging class of bioelectronic devices are as follows: (1) the iterative design and optimization process can be done in minutes to properly tune the maximum flowrate by studying only the unique combinations of the nondimensional parameters, subject to practical limits in the fabrication process, rather than each individual dimensional parameter (e.g., electric current, device, and microchannel geometries) and help to study only optimized geometries that need to be characterized experimentally using microparticle tracking velocimetry, a confocal microscope technique used to examine flowrate characteristics in microfluidic devices with drug delivery capabilities. Noting that the “slow + fast” variable analysis presented here is for a bioelectronic device where the bending effects in the polymer flexible membrane are not negligible since the maximum flowrate occurs at the beginning of the delivery process, a similar analysis can be performed for bioelectronic devices mainly experiencing large deformation (i.e., stretchingdominated deformation) to determine the total drug delivery time, flowrate, and its maximum value. Although further experimental testing is required to scale these bioelectronics from small animals to medium and large animals for targeted drug delivery studies, the proposed “slow + fast” variable analytical model for the flowrate and its maxima provides a scalable understanding to control the flowrate via , , and during the drug delivery process.
4. Materials and Methods
4.1. FEA of the Flexible Membrane
The commercial software ABAQUS was used to calculate the function for a flexible membrane with a thickness and radius . A pressure load of was applied uniformly to the bottom surface of the flexible membrane to deform it into the shape of a spherical cap. The flexible membrane was meshed with 3D stress elements (C3D8H), and the total number of 3D elements used in the model is ~50,000. The displacement and rotational degrees of freedom of the element nodes located at the circumference of the flexible membrane were fixed to zero; the remaining elements were free to deform because of the applied pressure resulting in a spherical cap shape. To calculate the deformed volume of the flexible membrane, the changes in the element volume (EVOL) and vertical displacement (U) were output for 250 evenly spaced time intervals. The hyperelastic material properties of the flexible membrane were defined from the uniaxial stressstrain relationship data shown in Figure 2(a) using the Marlow hyperelastic model to build the strain energy function of the SIS block copolymers.
4.2. Numerical Model for Drug Delivery
The numerical solver ode45 was used in MATLAB to solve the 1^{st}order ODE in Equation (6) with the initial condition The step size was set to 1/10,000^{th} of the total time, and the function was obtained from FEA. All the parameters used in the model are listed in Table 2.
Data Availability
The data sets generated and supporting the findings of this article are obtainable from the corresponding author upon reasonable request. The authors attest that all data for this study are included in the paper.
Conflicts of Interest
The authors declare no conflict of interest.
Authors’ Contributions
Raudel Avila is responsible for the conceptualization, methodology, validation, formal analysis, investigation, and visualization; wrote the original draft; and wrote, reviewed, and edited the manuscript. Yixin Wu is responsible for the conceptualization, validation, investigation, and data curation and wrote, reviewed, and edited the manuscript. Rinaldo Garziera is responsible for the conceptualization, methodology, and formal analysis. John A. Rogers is responsible for the conceptualization and resources; wrote, reviewed, and edited the manuscript; and supervised the study. Yonggang Huang is responsible for the conceptualization, methodology, validation, formal analysis, investigation, resources, visualization, supervision, and project administration; wrote the original draft; and wrote, reviewed, and edited the manuscript.
Acknowledgments
R.A. acknowledges support from the National Science Foundation Graduate Research Fellowship (NSF grant number DGE1842165) and from the Ford Foundation Predoctoral Fellowship.
Supplementary Materials
Supplementary Note 1: derivation of the upper bound time for bendingdominated deformation. Supplementary Note 2: derivation of the analytical model for drug delivery. Supplementary Note 3: derivation of the analytical maximum flowrate. (Supplementary Materials)
References
 A. B. Bußmann, L. M. Grünerbel, C. P. Durasiewicz, T. A. Thalhofer, A. Wille, and M. Richter, “Microdosing for drug delivery application—a review,” Sensors and Actuators, A: Physical, vol. 330, p. 112820, 2021. View at: Publisher Site  Google Scholar
 A. Kumar and J. Pillai, “Implantable drug delivery systems: an overview,” Nanostructures for the Engineering of Cells, Tissues and Organs, pp. 473–511, 2018. View at: Publisher Site  Google Scholar
 A. M. Vargason, A. C. Anselmo, and S. Mitragotri, “The evolution of commercial drug delivery technologies,” Nature Biomedical Engineering, vol. 5, no. 9, pp. 951–967, 2021. View at: Publisher Site  Google Scholar
 M. L. Laracuente, M. H. Yu, and K. J. McHugh, “Zeroorder drug delivery: state of the art and future prospects,” Journal of Controlled Release, vol. 327, pp. 834–856, 2020. View at: Publisher Site  Google Scholar
 S. Spieth, A. Schumacher, T. Holtzman et al., “An intracerebral drug delivery system for freely moving animals,” Biomedical Microdevices, vol. 14, no. 5, pp. 799–809, 2012. View at: Publisher Site  Google Scholar
 B. Mathon, M. Nassar, J. Simonnet et al., “Increasing the effectiveness of intracerebral injections in adult and neonatal mice: a neurosurgical point of view,” Neuroscience Bulletin, vol. 31, no. 6, pp. 685–696, 2015. View at: Publisher Site  Google Scholar
 N. Takeda and M. Diksic, “Relationship between drug delivery and the intraarterial infusion rate of SarCNU in C6 rat brain tumor model,” Journal of NeuroOncology, vol. 41, no. 3, pp. 235–246, 1999. View at: Publisher Site  Google Scholar
 M. Y. Chen, R. R. Lonser, P. F. Morrison, L. S. Governale, and E. H. Oldfield, “Variables affecting convectionenhanced delivery to the striatum: a systematic examination of rate of infusion, cannula size, infusate concentration, and tissuecannula sealing time,” Journal of Neurosurgery, vol. 90, no. 2, pp. 315–320, 1999. View at: Publisher Site  Google Scholar
 P. F. Morrison, M. Y. Chen, R. S. Chadwick, R. R. Lonser, and E. H. Oldfield, “Focal delivery during direct infusion to brain: role of flow rate, catheter diameter, and tissue mechanics,” American Journal of PhysiologyRegulatory Integrative and Comparative Physiology, vol. 277, no. 4, pp. R1218–R1229, 1999. View at: Publisher Site  Google Scholar
 Y. Zhang, D. C. Castro, Y. Han et al., “Batteryfree, lightweight, injectable microsystem for in vivo wireless pharmacology and optogenetics,” Proceedings of the National Academy of Sciences of the United States of America, vol. 116, no. 43, pp. 21427–21437, 2019. View at: Publisher Site  Google Scholar
 Y. Zhang, A. D. Mickle, P. Gutruf et al., “Batteryfree, fully implantable optofluidic cuff system for wireless optogenetic and pharmacological neuromodulation of peripheral nerves,” Science Advances, vol. 5, no. 7, pp. 1–12, 2019. View at: Publisher Site  Google Scholar
 D. Rosenblum, N. Joshi, W. Tao, J. M. Karp, and D. Peer, “Progress and challenges towards targeted delivery of cancer therapeutics,” Nature Communications, vol. 9, no. 1, p. 1410, 2018. View at: Publisher Site  Google Scholar
 S. Senapati, A. K. Mahanta, S. Kumar, and P. Maiti, “Controlled drug delivery vehicles for cancer treatment and their performance,” Signal Transduction and Targeted Therapy, vol. 3, no. 1, pp. 1–19, 2018. View at: Publisher Site  Google Scholar
 A. Cobo, R. Sheybani, H. Tu, and E. Meng, “A wireless implantable micropump for chronic drug infusion against cancer,” Sensors and Actuators, A: Physical, vol. 239, pp. 18–25, 2016. View at: Publisher Site  Google Scholar
 R. Hunt Bobo, D. W. Laske, A. Akbasak, P. F. Morrison, R. L. Dedrick, and E. H. Oldfield, “Convectionenhanced delivery of macromolecules in the brain,” Proceedings of the National Academy of Sciences of the United States of America, vol. 91, no. 6, pp. 2076–2080, 1994. View at: Publisher Site  Google Scholar
 R. Qazi, A. M. Gomez, D. C. Castro et al., “Wireless optofluidic brain probes for chronic neuropharmacology and photostimulation,” Nature Biomedical Engineering, vol. 3, no. 8, pp. 655–669, 2019. View at: Publisher Site  Google Scholar
 P. Y. Li, J. Shih, R. Lo et al., “An electrochemical intraocular drug delivery device,” Sensors and Actuators, A: Physical, vol. 143, no. 1, pp. 41–48, 2008. View at: Publisher Site  Google Scholar
 D. Fan, D. Rich, T. Holtzman et al., “A wireless multichannel recording system for freely behaving mice and rats,” PLoS One, vol. 6, no. 7, p. e22033, 2011. View at: Publisher Site  Google Scholar
 R. Qazi, C. Yeon Kim, I. Kang, D. Binazarov, J. G. McCall, and J. W. Jeong, “Implantable optofluidic systems for wireless in vivo photopharmacology,” ChemPhotoChem, vol. 5, no. 2, pp. 96–105, 2021. View at: Publisher Site  Google Scholar
 Y. Chen, N. J. Rommelfanger, A. I. Mahdi et al., “How is flexible electronics advancing neuroscience research?” Biomaterials, vol. 268, article 120559, 2021. View at: Publisher Site  Google Scholar
 K. N. Noh, S. I. Park, R. Qazi et al., “Miniaturized, batteryfree optofluidic systems with potential for wireless pharmacology and optogenetics,” Small, vol. 14, no. 4, article 1702479, 2018. View at: Publisher Site  Google Scholar
 L. Vermaak, H. W. J. P. Neomagus, and D. G. Bessarabov, “Recent advances in membranebased electrochemical hydrogen separation: a review,” Membranes, vol. 11, no. 2, pp. 127–132, 2021. View at: Publisher Site  Google Scholar
 R. Avila, C. Li, Y. Xue, J. A. Rogers, and Y. Huang, “Modeling programmable drug delivery in bioelectronics with electrochemical actuation,” Proceedings of the National Academy of Sciences of the United States of America, vol. 118, no. 11, p. e2026405118, 2021. View at: Publisher Site  Google Scholar
 R. Avila, Y. Wu, J. A. Rogers, and Y. Huang, “A mechanics model for injectable microsystems in drug delivery,” Journal of the Mechanics and Physics of Solids, vol. 156, article 104622, 2021. View at: Publisher Site  Google Scholar
 R. Avila, J. Ciatti, A. VázquezGuardado et al., “Electrochemical bioelectronics in drug delivery  effect of the initial gas volume,” Journal of Applied Mechanics, vol. 89, no. 1, 2022. View at: Publisher Site  Google Scholar
 F. Forouzandeh, N. N. Ahamed, X. Zhu et al., “A wirelessly controlled scalable 3Dprinted microsystem for drug delivery,” Pharmaceuticals, vol. 14, no. 6, 2021. View at: Publisher Site  Google Scholar
 J. S. Speed and K. A. Hyndman, “_In vivo_ organ specific drug delivery with implantable peristaltic pumps,” Scientific Reports, vol. 6, no. 1, pp. 1–7, 2016. View at: Publisher Site  Google Scholar
 Z. V. Taylor, B. Khand, A. Porgador, A. Monsonego, and E. Eremenko, “An optimized intracerebroventricular injection of CD4^{+} T cells into mice,” STAR Protocols, vol. 2, no. 3, article 100725, 2021. View at: Publisher Site  Google Scholar
 J. W. Jeong, J. G. McCall, G. Shin et al., “Wireless optofluidic systems for programmable in vivo pharmacology and optogenetics,” Cell, vol. 162, no. 3, pp. 662–674, 2015. View at: Publisher Site  Google Scholar
 H. Fujii, S. Horie, K. Takeda, S. Mori, and T. Kodama, “Optimal range of injection rates for a lymphatic drug delivery system,” Journal of Biophotonics, vol. 11, no. 8, pp. 1–7, 2018. View at: Publisher Site  Google Scholar
 K. Moussi, A. Bukhamsin, and J. Kosel, “Implantable 3D printed drug delivery system,” in 2019 20th International Conference on SolidState Sensors, Actuators and Microsystems & Eurosensors XXXIII (TRANSDUCERS & EUROSENSORS XXXIII), pp. 2243–2246, Berlin, Germany, 2019. View at: Publisher Site  Google Scholar
 Dassault Systèmes, Abaqus 6.11 Theory Manual, Dassault Systèmes Simulia, Inc., 2013.
 R. S. Marlow, “A General FirstInvariant Hyperelastic Constitutive Model,”,” in Constitutive Models For Rubber – Proceedings Constitutive Models for Rubber, pp. 157–160, London, UK, 2003. View at: Google Scholar
 G. C. Pardoen and R. L. Hagen, “Symmetrical bending of circular plates using finite elements,” Computers and Structures, vol. 2, no. 4, pp. 547–553, 1972. View at: Publisher Site  Google Scholar
 A. S. Feiner and A. J. McEvoy, “The Nernst equation,” Journal of Chemical Education, vol. 71, no. 6, pp. 493494, 1994. View at: Publisher Site  Google Scholar
 D. E. Lee, S. A. Soper, and W. Wang, “Fabrication and mathematical analysis of an electrochemical microactuator (ECM) using electrodes coated with platinum nanoparticles,” Microsystem Technologies, vol. 16, no. 3, pp. 381–390, 2010. View at: Publisher Site  Google Scholar
 M. Kumar and Parul, “Methods for solving singular perturbation problems arising in science and engineering,” Mathematical and Computer Modelling, vol. 54, no. 12, pp. 556–575, 2011. View at: Publisher Site  Google Scholar
Copyright
Copyright © 2022 Raudel Avila et al. Exclusive Licensee Science and Technology Review Publishing House. Distributed under a Creative Commons Attribution License (CC BY 4.0).