Generalizing Murray ’ s law : An optimization principle for fluidic networks of arbitrary shape and scale

Murray's law states that the volumetric flow rate is proportional to the cube of the radius in a cylindrical channel optimized to require the minimum work to drive and maintain the fluid. However, application of this principle to the biomimetic design of micro/nano fabricated networks requires optimization of channels with arbitrary cross-sectional shape (not just circular) and smaller than is valid for Murray's original assumptions. We present a generalized law for symmetric branching that (a) is valid for any cross-sectional shape, providing that the shape is constant through the network; (b) is valid for slip flow and plug flow occurring at very small scales; and (c) is valid for networks with a constant depth, which is often a requirement for lab-on-a-chip fabrication procedures. By considering limits of the generalized law, we show that the optimum daughter-parent area ratio Γ, for symmetric branching into N daughter channels of any constant cross-sectional shape, is Γ=N−2/3 for large-scale channels, and Γ=N−4/5 for channels with a characteristic length scale much smaller than the slip length. Our analytical results are verified by comparison with a numerical optimization of a two-level network model based on flow rate data obtained from a variety of sources, including Navier-Stokes slip calculations, kinetic theory data, and stochastic particle simulations.

Murray's law states that the volumetric flow rate is proportional to the cube of the radius in a cylindrical channel optimized to require the minimum work to drive and maintain the fluid.However, application of this principle to the biomimetic design of micro/nano fabricated networks requires optimization of channels with arbitrary cross-sectional shape (not just circular) and smaller than is valid for Murray's original assumptions.We present a generalized law for symmetric branching that (a) is valid for any cross-sectional shape, providing that the shape is constant through the network; (b) is valid for slip flow and plug flow occurring at very small scales; and (c) is valid for networks with a constant depth, which is often a requirement for lab-on-a-chip fabrication procedures.By considering limits of the generalized law, we show that the optimum daughter-parent area ratio C, for symmetric branching into N daughter channels of any constant cross-sectional shape, is C ¼ N À2=3 for large-scale channels, and C ¼ N À4=5 for channels with a characteristic length scale much smaller than the slip length.Our analytical results are verified by comparison with a numerical optimization of a two-level network model based on flow rate data obtained from a variety of sources, including Navier-Stokes slip calculations, kinetic theory data, and stochastic particle simulations.V C 2015 AIP Publishing LLC.[http://dx.doi.org/10.1063/1.4935288]

I. INTRODUCTION
In 1926, Murray 1 posited that there were two competing factors contributing to the energy cost of blood flow through the arterial system: (1) the energy required to drive the flow, which increases as the vessel radius decreases; and (2) the energy required to metabolically maintain the fluid, which increases with increasing vessel radius.Thus, to minimize the total power requirement, the vessel could be neither too large nor too small.Using the Hagen-Poiseuille law to describe the flow through a cylindrical vessel (i.e., assuming the flow is laminar, Newtonian, steady, and fully developed), the power W f required for the flow to overcome the viscous drag is where DP is the pressure drop over the vessel, l is the dynamic viscosity of the fluid, L is the vessel length, f is the volumetric flow rate, and r is the vessel radius.This is offset by the maintenance "cost of blood" W m which increases linearly with the blood volume: where m is an all-encompassing metabolic coefficient that includes the chemical cost of keeping the blood constituents fresh and functional, and the general cost owing to the weight of the blood and the vessel.The total power requirement where a ¼ 8lL/p and b ¼ mpL.For a constant volumetric flow rate, and given values of a and b, the total power will be a function of a single variable-the vessel radius r-and the minimum power is found by differentiating with respect to r and equating to zero: Rearranging Eq. ( 4) gives where . Thus for any vessel considered independently, Eq. ( 5) describes the optimal relation between volumetric flow rate and vessel radius, such that the power requirement is minimized.If the fluid viscosity and metabolic coefficient are constant throughout a network, then K is constant.It was surmised that in this instance, Eq. ( 5) should hold for all vessels in a network operating at maximum efficiency, which tacitly assumes that the local losses through the junction (due to bends and channel contractions) are negligible compared to the pressure losses over the channel lengths; this limits the applicability of Murray's law to networks of high-aspect-ratio channels.By applying the conservation of mass at a branching point, we retrieve the ubiquitous principle known as Murray's law: where the subscripts p and d i denote the parent and the ith daughter (of N), respectively.In the common case of symmetrical branching, Eq. ( 6) reduces to ][15] Murray's law can also be applied to inorganic systems 5 and has been adapted for networks with rectangular, trapezoidal, 16 and elliptical cross sections. 17The applicability of Murray's law to moderately rarefied gas flows has also been investigated 18 and a departure from Eq. ( 7) has been noted.These developments indicate that some version of Murray's law could be applied as a biomimetic design principle for microfluidic and nanofluidic networks, such as lab-on-a-chip devices for microreactors 19 or tissue engineering, 20 or micro/ nanoscale heat exchangers for high performance fuel cells 21 or the cooling of electronic devices. 22,23ince the original derivation of Murray's law, it has been noted that the application of other optimization principles (not just that of minimum work) result in Eq. ( 7): minimizing the total mass of the network, 24 minimizing volume for a constant pressure drop and flow rate, 25 keeping the shear stress constant in all channels, 26 or minimizing flow resistance for a constant volume. 5,27s noted in Ref. 16, in some circumstances, it is desirable or necessary for the cross-sectional shape to vary between the parent and daughter channels of a branching network.In many lab-on-a-chip fabrication procedures (e.g., photolithograpy, wet or dry etching, or surface micromachining), the depth remains constant throughout the device, and thus the shape (i.e., aspect ratio) of the cross section must vary.A multi-depth approach to fabrication does exist, 28 but it is relatively complex.
Despite many developments to Murray's law, there are three major barriers that prevent it from being relevant to the design of many artificial fluidic networks: (1) it is not demonstrably applicable to cross sections of any arbitrary shape; (2) it is not applicable at the micro/nanoscale, where a fluid can no longer be accurately described as a continuous material; 29 and (3) it is not applicable to networks which maintain a constant depth through branching, wherein the cross-sectional shape changes between the parent and daughter channels.

II. ANALYTICAL SOLUTIONS
Although Murray's original optimization concerns the minimization of work, the principle can be generalized as a maximization of flow conductance per unit volume for a variety of constraint combinations.For a two-level network (consisting of a single parent channel branching into multiple daughter channels), this can be expressed as arg max where Q is the mass flow rate through the parent channel, DP is the total pressure drop (from inlet of the parent to the outlet of the daughters), V is the network volume, and is the daughter-parent cross-sectional area ratio.Note, the three constraint options (pressure-drop minimization, volume minimization, and flow-rate maximization) all lead to the same optimal daughter-parent area ratio.For our optimization, we assume that the channel lengths L are large compared to the size of the parent-daughters junction so that (1) the localized pressure losses through the junction are negligible compared to the pressure drops over individual channels (as in Murray's law) and ( 2) the volume of the network can be considered to be the sum of the channel volumes: This means that the channel lengths are treated as being independent of the branching angle. 30The channel lengths are also treated as being independent of the optimal daughterparent area ratio (as found in Murray's optimization 1 ), i.e., dL/dC ¼ 0. This will be verified later.Inserting Eq. ( 10) into the fitness function of Eq. ( 8), differentiating with respect to C, and equating to zero, gives d dC noting that, for all constraint combinations, d(DP)/dC ¼ 0 and dQ/dC ¼ 0. The pressure drop over the parent and daughter channels can be expressed in terms of the mass flow rate: where L is the channel length and k is flow resistance per unit length, e.g., k p ¼ DP p /(QL p ).The pressure drop over the entire network As d(DP)/dC ¼ 0 and dQ/dC ¼ 0, differentiating Eq. ( 14) with respect to C gives Substituting Eq. ( 15) into Eq.( 11), via the chain rule, gives the generalized optimal area relation for symmetric branching: For brevity, this will be referred to as the generalized law to distinguish it from Murray's law.This expression is not a function of L which suggests that the optimal daughterparent area ratio is independent of the channel lengths (and thus branching angle).Note this generalized law is valid for any Reynolds number (e.g., for turbulent flow) and for any fluid (e.g., non-Newtonian) 31 as long as local losses in the junction remain negligible.In this paper, we restrict our attention to laminar and Newtonian flows, and we now consider some important cases where A can be expressed easily as an analytical function of k.
A. The continuum-flow limit We begin with the steady incompressible Navier-Stokes momentum equation for laminar flow through a long channel with an arbitrary cross-sectional shape, i.e., where u is the streamwise channel velocity.We can nondimensionalize this using DP/L, l and cross-sectional area A, such that where and tilde denotes a dimensionless quantity or operator.We also define y, z as the axes of the cross-sectional plane, and Provided the boundary conditions are fixed (which is the case for the continuum-flow limit, where the no-slip boundary condition applies), the solution of Eq. ( 18), ũðỹ; zÞ, is independent of A, DP, L, and l, and is thus a property of the cross-sectional shape alone.Similarly, so is An expression for the mass flow rate is obtained by integrating the fluid momentum over the cross-sectional area where q is the mass density.Substitution of Eqs. ( 19)-( 21) into (22) gives the mass flow rate for an arbitrary crosssectional shape and the flow resistance per unit length It is assumed that the pressure drop over the network is small such that the viscosity and density are constants.For cylindrical channels, S ¼ 1/(8p) and Eq. ( 23) becomes the Hagen-Poiseuille flow rate.

Constant-shape networks
If the cross-sectional shape is constant throughout the network (i.e., S ¼ const), substituting Eq. ( 24) into the generalized law (Eq.( 16)) gives Equation ( 25) relates the area of the parent channel to the area of the daughter channel in an optimized symmetric network.
Other studies have found that Murray's law holds for some specific cross-sectional shapes (i.e., rectangles and trapezoids, 16 and ellipses 17 ), but this demonstrates that Murray's law for symmetric branching is in fact applicable to any cross-sectional shape, provided the shape is constant through the network.Note this means that if the cross-sectional areas are expressed in terms of a hydraulic radius, Murray's law can be expressed in its original form (i.e., Eq. ( 7)) for any shapehowever, this is only applicable for symmetric branching at the continuum-flow limit.

Constant-depth networks
For a constant-depth network, the optimal branching behavior depends on the change in shape between the parent and daughter channels, which is a function of the aspect ratio a ¼ A/h 2 for rectangular channels.This leads to two aspectratio limits for C (at the continuum-flow limit), as illustrated in Fig. 1.An accurate approximation of the channel shape property S for a rectangle 32 is FIG. 1. Schematic showing the symmetric bifurcation of constant-depth rectangular channels of (a) high aspect ratio (a ) 1) and (b) low aspect ratio (a ( 1).

174302-3
Stephenson et al.J. Appl.Phys.118, 174302 (2015 where h is the constant network depth.Equation ( 26) only includes the first term of an infinite series, but as the denominators of higher terms increase exponentially, the expression is sufficient.It is clear from Eq. ( 26) that at the aspect ratio limits, the shape property reduces to S ¼ A 12h 2 ; for a ( 1: (28)   As S is no longer a constant, the generalized law must consider the change in shape between parent and daughter channels At the aspect ratio limits, Eq. ( 29) reduces to Equations ( 30) and ( 31) can also be derived directly from Eqs. ( 27) and ( 28), respectively.Substituting Eqs. ( 24), (27), and (30) into the generalized law (Eq.( 16)), via the chain rule, gives the high-aspect-ratio continuum-flow limit: Incidentally, this limit is equal to da Vinci's rule of tree branching. 33Similarly, substituting Eqs. ( 24), (28), and (31) into the generalized law gives the low-aspect-ratio continuum-flow limit: B. The plug-flow limit For the other limiting case of scale, the flow becomes dominated by non-continuum and non-equilibrium effects, for liquids and dilute gases, respectively.For more information, see the reviews by Gad-el-Hak, 29 Karniadakis et al., 34 and Conlisk. 35The complexity of channel flow increases dramatically in such conditions, e.g., density layering can occur near to the bounding surface in nanoscale liquid channel flows 36 and Knudsen layers (kinetic boundary layers) can occur in rarefied-gas channel flows. 37Furthermore, the stateof-the art simulation techniques for such flows-molecular dynamics for liquid and dense gas flows, 38 and the direct simulation Monte Carlo method for dilute gas flows 39 although accurate, are extremely computationally expensive.Therefore, for the purpose of generating a simple analytical relation, we make use of the general approximate observation that at very small scales (i.e., small channel cross sections or high degrees of rarefaction), the channel velocity profile assumes a "plug-like" form, wherein the average velocity approximates the velocity of the fluid at the walls (the slip velocity).Navier's slip condition (which is only approximate for these conditions, as discussed in Section II C) allows us to express the wall shear stress s w as a function of the slip velocity u s , which must be nearly constant around the perimeter of the cross section, given the plugflow profile: where g is the slip length determined by the fluid-solid interaction.Note the value of g does not affect the optimal area ratio solution at the plug-flow limit, and thus accurate knowledge of it is not required.As the flow tends to the plug-flow limit, the mass flow rate is simply and a force balance relates the pressure drop to the wall shear, i.e., where P is the perimeter of the section.Combining Eqs. ( 34)- (36) gives the mass flow rate and the flow resistance per unit length where =P is a property of the cross-sectional shape (like S).

Constant-shape networks
If we consider the case where the shape remains the same throughout the network, Eq. ( 38) can be substituted into the generalized law (16) to get the plug-flow limit for an arbitrary cross-sectional shape:

Constant-depth networks
For networks of rectangular cross section with a constant depth, Eq. ( 38) can be rewritten as where the perimeter is Noting again that a ¼ A/h 2 , it can be seen from Eq. ( 40) that at the aspect-ratio limits, k reduces to Substituting Eq. ( 42) into the generalized law gives for the high-aspect-ratio plug-flow limit.Interestingly, the high-aspect-ratio limits are the same for continuum-and plug-flow.Similarly, substituting Eq. ( 43) into the generalized law gives for the low-aspect-ratio plug-flow limit.

C. A slip flow approximation
In certain regimes, between the continuum-and plugflow limits, it is sufficiently accurate to model the channel flow using Navier-Stokes equations with modified boundary conditions accounting for the velocity slip at wall boundaries-such solutions are referred to as slip solutions.There are number of types of slip boundary condition (see, e.g., those in Refs.34 and 40), but the most common is Navier's slip condition (for which Maxwell's slip boundary condition 41 is a special case for gas flows).Slip solutions using Navier's condition (s w ¼ u s l/g) have shown to yield accurate results, as compared to molecular dynamics predictions, for liquid flows in channels as small as $1-2 nm for water 42 and $2-3 nm for Lennard Jones fluids. 43,44For gases, provided the Knudsen number (the ratio of the mean free path to the channel height) is small, the slip solution is also a wellknown and good approximation.
A further approximation to the exact slip solution, for any cross section, can be obtained by assuming that the shear stress, and thus the velocity slip, is constant around the perimeter, i.e., This gives the exact mass flow rate of the slip solution for a circular cross section (which has a uniform shear stress) and accurate approximations for rectangles of any aspect ratio (within 3% of results from a finite-difference slip solver).Note that the slip solution is itself an approximation, and particularly for gas flows must be treated with caution, as discussed later.Equation ( 46) can be rearranged to give the flow resistance per unit length Inserting Eq. ( 47) into the generalized law (Eq.( 16)) for parent and daughter channels produces the general slip solution for the optimal daughter-parent area ratio for all crosssectional shapes across all length-scales (48)

Constant-shape networks
For all networks with a constant shape, S p ¼ S d and dS/ dA ¼ 0. For a circular cross section, A ¼ pR 2 ; S ¼ 1=ð8pÞ; , and dP=dA ¼ ffiffiffiffiffiffiffiffi ffi p=A p . Substituting these values into Eq.( 48) and simplifying produces where rp ¼ r p =g.For rectangular cross sections, P ¼ 2ða þ 1Þh and dP=dA ¼ ða þ 1Þh=A.A constant shape means a variable depth h and a constant aspect ratio a¼A/h 2 .Substituting these values into Eq.( 48) and simplifying gives an optimal area relation of where S is calculated using Eq. ( 26) hp ¼ h p =g and C ¼ 2S(a þ 1) is a constant.

Constant-depth networks
For rectangular cross sections of a constant depth and variable aspect ratio, substituting the expression for dS/dA (from Eq. ( 29)) into Eq.( 48), along the rectangle parameters previously outlined, gives

III. NUMERICAL VERIFICATION AND DISCUSSION
An accurate numerical optimization procedure is now described to demonstrate the following: (a) the generalized law (Eq.( 16)) is applicable to networks with arbitrary crosssectional shape and for all scales; (b) the limits of the generalized law, identified as the continuum-flow limit (Eqs.( 25), (32), and ( 33)) and the plug-flow limit (Eqs.( 39), (44), and ( 45)), are valid and precise for all shapes considered; and (c) the approximate slip solutions to the generalized law (Eqs.( 49)-( 51)) provide reasonable accuracy, even for rarefied gas flows.
The numerical optimization procedure is based on calculations that use non-dimensional mass flow rates, i.e., These mass flow rate are obtained either from published sources, from high-resolution finite-difference slip solutions, from exact analytical expressions, or from stochastic particle calculations (for dilute gases).The numerical model of the network also assumes that channels are sufficiently long such that pressure losses at the parent-daughters junction are negligible compared to the pressure drops over the channels themselves.This allows a model for a two-level network to be constructed from predictions of mass flow rate through individual channels, constrained by a common branching pressure and the requirement for mass continuity.For each shape and physical model, an interpolant is constructed that provides the non-dimensional mass flow rate for any given area, QðAÞ.From this, Eq. ( 52) can be evaluated for parent and daughter channels Combining Eqs. ( 53) and (54), and substituting for the total pressure drop DP ¼ DP p þ DP d , produces a model of the mass flow rate through a two-level network for a particular cross-sectional shape and physical model: For symmetric branching, the volume constraint requires that L p A p ¼ V À NA d L d .The numerical procedure uses an interior-point constrained optimization algorithm 45 in MATLAB V R to find the parent and daughter areas that, in combination, maximize the mass flow rate through the network (as predicted by Eq. ( 55)) for arbitrarily fixed total volume, total pressure drop, channel lengths, and fluid properties.Alternatively, the same result can be obtained via a brute force approach to ensure that the numerical optimization finds a global, not local, maximum.

A. Different shapes and sizes
Murray's original derivation was for circular channel sections at the continuum-flow limit.The first set of optimization results are presented to verify that the generalized law is valid for a variety of cross-sectional shapes across all length scales, with the shape remaining constant through branching.The numerical optimization uses non-dimensional mass flow rates calculated from a standard central-difference solution of the laminar Navier-Stokes equations with a Navier slip boundary condition. 46The results from the analytical and numerical optimizations, for symmetric bifurcations, are presented in Fig. 2.
For large parent areas, relative to the square of the slip length g, the optimum daughter-parent area ratio converges to the continuum-flow limit of the generalized law (Eq.( 25)) for all shapes considered.The same is true for the other extreme of scale: the optimum area ratio for all shapes converges to the plug-flow limit of the generalized law (Eq.( 39)) for small parent areas.In the transition between these limits, the approximate slip solutions of the generalized law are also highly accurate.The difference between the slip approximations of the generalized law and the numerical optimization is less than 0.1% across the entire range of scales for all shapes tested.Clearly, predictions for the optimum dimensions are not particularly sensitive to errors introduced by the approximation in Eq. (46).These results are consistent regardless of the daughter and parent channel length magnitudes, demonstrating that the optimal daughterparent area ratio is independent of the lengths, as was asserted in the analytical solution.
In Fig. 2, it is observed that for networks of channels with rectangular cross sections, the aspect ratio affects the range of areas for which C is in the transition period between the continuum-and plug-flow limits.This occurs because the influence of slip is governed by the size of the smallest cross-sectional length scale (hereon referred to as the characteristic length L) relative to the slip length.So, for the same cross-sectional area, the characteristic length of a rectangle with a high-aspect-ratio is less than the characteristic length of a low-aspect-ratio rectangle.This means that optimal branching of high-aspect-ratio channels will depart from the continuum-flow limit and approach the plug-flow limit at larger parent areas than optimal branching in low-aspectratio channels.FIG. 2. Optimal daughter-parent area ratio C against non-dimensional parent area for networks with a constant shape.Comparison of the approximate slip solution to the generalized law (Eqs.(49) and (50)) and the numerical optimization using data from a Navier-Stokes slip solver.Plotted for circles, squares, and rectangles of aspect ratio a ¼ 5, a ¼ 10, and a ¼ 100 at N ¼ 2.

B. Rarefied gas flow
The numerical optimization is now performed for symmetric branching of rarefied gas flows.The mass-flow-rate data are obtained from a variety of sources, including our own simulations 47 computed using a version of lowvariance deviational simulation Monte Carlo (LVDSMC). 48dditional mass-flow-rate data are procured from solutions of the linearized Boltzmann equation 49,50 to verify the accuracy of the LVDSMC results.For dilute gases, the slip length can be related to the mean free path k via g ¼ C s k, where C s ¼ 1.11 is the first-order slip coefficient for the hard-sphere model of gases with purely diffuse molecular reflection at walls. 51The results for the analytical and numerical optimizations, for circular cross sections, are presented in Fig. 3. To demonstrate that the generalized law is valid for any number of daughter branches (N ! 2), C is plotted for N ¼ 2, N ¼ 3 and N ¼ 5 and, for clarity, is normalized with respect to the continuum-flow limit (N À2=3 ).
Again, the agreement between the numerical optimization and the plug-flow limit of the generalized law is excellent for each case considered.It is perhaps unexpected that a slip solution to the generalized law should converge to the same result as that of kinetic theory and LVDSMC at the free-molecular limit, given the approximate nature of slip boundary condition at such scales.However, as rp !0, C in Eq. ( 49) becomes independent of the slip length g, and is thus unaffected by any inaccuracy in the slip model.Due to computational cost, molecular simulations are not performed for sufficiently large areas to see the solution meet the continuum-flow limit; but, since the results from kinetic theory converge to the solution of the slip model, agreement at the continuum-flow limit is also expected.It is well known that when the Knudsen number Kn ¼ k=L is greater than $ 0.1, slip solutions become inaccurate, explaining the departure in C between the limits.The kinetic theory and LVDSMC results all show a minimum in C beneath the plug-flow/free-molecular limit.This is possibly a manifestation of the Knudsen minimum, 52 a rarefied gas phenomenon that occurs when the diffusive flux starts to dominate the convective flux as the length scale decreases. 53ote although Eq. ( 49) is approximate between the scale limits, the precise result of the numerical optimization can be reclaimed by expressing the interpolated dimensionless mass flow rate, QðAÞ, as a flow resistance per unit length, k(A).This allows the evaluation of the generalized law in (16), but, of course, does not afford an analytical relation.

C. Constant-depth rectangles
Finally, we perform a numerical optimization for a symmetrically bifurcating network of rectangular channels with a constant depth and variable aspect ratio, mimicking the conditions of a micro-fabricated fluidic network.Here, we have used non-dimensional mass-flow-rate data obtained from the two-dimensional Navier-Stokes slip solver.Figure 4 shows the optimum daughter-parent area ratio against dimensionless area (A/h 2 ; incidentally, this is equal to the aspect ratio, a) for rectangles of a constant dimensionless depth h ¼ h=g.
Agreement between the numerical optimization and the approximate slip solution to the generalized law (Eq.( 51)) is excellent for all cases, with the difference being less than 0.5% across the entire range of scales.The results also show convergence to the limits derived in Eqs. ( 32), ( 33), (44), and (45).It can be seen that when the non-dimensional depth h is small, such that the slip length is relatively large, C varies only between the plug-flow limits for high and low aspect ratios (Eqs.( 44) and ( 45), respectively).The highaspect-ratio limits are the same for continuum-and plugflow.When h is large and the aspect ratio decreases, C first tends to the low-aspect-ratio limit continuum-flow (Eq.( 33)), until the area gets sufficiently small that the variable width is comparable to the slip length, at which point FIG. 3. Normalized optimal daughter-parent area ratio against nondimensional parent area.Comparison of the analytical slip solution to the generalized law (Eq.( 49)) and the numerical optimization using data from kinetic theory 49,50 and LVDSMC.Plotted for circles at N ¼ 2, N ¼ 3, and N ¼ 5. FIG. 4. Optimal daughter-parent area ratio C against non-dimensional parent area A/h 2 (equal to aspect ratio a) for rectangles of a constant dimensionless depth h.Comparison of the approximate slip solution to the generalized law (Eq.( 51)) and the numerical optimization using data from a Navier-Stokes slip solver.Plotted for rectangles of depth h ¼ 1; 10; 10 2 ; 10 3 ; 10 4 , and 10 5 , at N ¼ 2.
For constant-depth networks, C behaves very differently compared to constant-aspect-ratio networks.For the most part, when the depth is constant, the optimal daughter-parent area ratio increases with decreasing cross-sectional areathe opposite to the trend found in Figs. 2 and 3.In addition, the range of C values is wider, the minimum C is smaller, and the maximum C is larger for constant-depth networks: 0.5 C 0.71 compared to 0.57 C 0.63.However, there are some similarities.At the low-aspect-ratio limits, the trend for constant-depth networks is the same as for constant-aspect-ratio networks, with the optimal daughter-parent area ratio decreasing as cross-sectional area decreases.

IV. CONCLUSION
We have derived a generalized optimization principle that leads to analytical expressions for the optimum daughter-parent area ratio C for any shape, at any scale, and for any number of daughter branches.We have shown that at the continuum-flow limit, Murray's law is reclaimed (C ¼ N À2=3 ) for networks of any cross-sectional shape, as long as the shape is constant through branching.The generalized law presented is also valid for slip flows and plug flows (C ¼ N À4=5 ) that occur at smaller length scales, as well as for networks with a constant depth, which is often a requirement for lab-on-a-chip devices due to limitations in fabrication procedures.The new optimal design relation we propose will allow the original biomimetic design principle of Murray to be applied to a variety of micro and nanofluidic networks that require non-circular geometry, due to manufacturing constraints, and are designed for increasingly smaller scales in order to achieve a greater degree of control, functionality, and analytical and economic efficiency.Future work will include determining optimal branching angles and optimum channel lengths for arbitrary cross-sectional shapes, across all physical scales.
Generalizing Murray's law: An optimization principle for fluidic networks of arbitrary shape and scale David Stephenson, a) Alexander Patronis, b) David M. Holland, c) and Duncan A. Lockerby d) School of Engineering, University of Warwick, Coventry, CV4 7AL, United Kingdom (Received 17 August 2015; accepted 26 October 2015; published online 6 November 2015)