Computational hemodynamic optimization predicts dominant aortic arch selection is driven by embryonic outflow tract orientation in the chick embryo

In the early embryo, a series of symmetric, paired vessels, the aortic arches, surround the foregut and distribute cardiac output to the growing embryo and fetus. During embryonic development, the arch vessels undergo large-scale asymmetric morphogenesis to form species-specific adult great vessel patterns. These transformations occur within a dynamic biomechanical environment, which can play an important role in the development of normal arch configurations or the aberrant arch morphologies associated with congenital cardiac defects. Arrested migration and rotation of the embryonic outflow tract during late stages of cardiac looping has been shown to produce both outflow tract and several arch abnormalities. Here, we investigate how changes in flow distribution due to a perturbation in the angular orientation of the embryonic outflow tract impact the morphogenesis and growth of the aortic arches. Using a combination of in vivo arch morphometry with fluorescent dye injection and hemodynamics-driven bioengineering optimization-based vascular growth modeling, we demonstrate that outflow tract orientation significantly changes during development and that the associated changes in hemodynamic load can dramatically influence downstream aortic arch patterning. Optimization reveals that balancing energy expenditure with diffusive capacity leads to multiple arch vessel patterns as seen in the embryo, while minimizing energy alone led to the single arch configuration seen in the mature arch of aorta. Our model further shows the critical importance of the orientation of the outflow tract in dictating morphogenesis to the adult single arch and accurately predicts arch IV as the dominant mature arch of aorta. These results support the hypothesis that abnormal positioning of the outflow tract during early cardiac morphogenesis may lead to congenital defects of the great vessels due to altered hemodynamic loading.


Introduction
Defects of the great vessels represent more than 20% of congenital heart disease in the United States (Roger et al. 2011) and can require multiple surgeries and substantial risk to partially restore normal arterial structure and function. While the etiologies of such defects are not yet fully understood, most can be traced back to morphogenetic insults that occur during the growth and early patterning of the embryonic aortic arches, the precursors to the great vessels. These errors in morphogenesis are likely due to the combination of genetic events and environmental or epigenetic insults and embryo survival depends upon adaptive morphogenesis and remodeling. During embryonic development, up to six pairs of symmetric arch arteries (numbered I-VI) emerge consecutively in a cranio-caudal fashion to surround the anterior foregut (Kardong 2009). Initially, these vessels connect the heart to the embryonic circulation, serving as the conduits for the entire cardiac output. As the embryo develops, this immature symmetric system undergoes large-scale asymmetric morphogenesis to form species-specific adult great vessel patterns.
Given the similarities between avian and human cardiovascular development, as well as ease of access and optical transparency at early stages, the chick embryo has become a standard model organism for the study of embryonic cardiovascular morphogenesis, function, and biomechanics. The embryonic chick heart originates as a pair of non-pulsatile endocardial tubes formed on opposite sides of the embryo which then fuse ventrally to form the straight heart tube located at the midline of the embryo and align parallel to the cranio-caudal axis (Patten 1920). The heart tube is formed at Hamburger Hamilton stage 10 (33 h) and begins to beat almost immediately (Hamburger and Hamilton 1951;Martinsen 2005). From stages 10-12 (2 days), a process termed "cardiac looping" begins whereby the heart tube elongates, rotates along its longitudinal axis, and bends to the right (Manner 2000). Peristaltic contractions drive blood flow (caudal to cranial) through the sinus venosus, atrium, ventricle, outflow tract, and primitive aortic arch by stage 12. Cardiac looping continues to stage 24 (4 days) as the heart twists to form a looped tube, shifting the atrium dorsally and cranially and the ventricle and outflow tract ventrally and caudally (Manner 2000). During the initial phase of cardiac looping, the aortic arches are emerging and connect the outflow tract to the dorsal aorta (Hiruma and Hirakow 1995). Cranial-most arches I and II appear at stage 14, but neither persist, disappearing at stages 19 and 21, respectively. Arch pair III forms by stage 15 and both laterals persist as portions of the brachiocephalic and common carotid arteries. Arch pair IV arises at stage 18, and the left lateral involutes by stage 34 while the right forms a segment of the vital systemic adult aortic arch. Arch pair V either appears transiently or never completely forms. The caudal-most arch pair VI emerges at stage 21 and contributes to the segments of the central pulmonary arteries and ductus arteriosus. This sequence of aortic arch patterning is summarized in Fig. 1.
The emergence of arch pairs I-VI coincides with the formation of the cervical flexure, during which the cranio-caudal axis of the chick embryo transforms from a straight line to a U-shape, bending at the level of the pharyngeal arches (Flynn et al. 1991). Previous studies have shown that this process is closely related to cardiac looping (Flynn et al. 1991;Manner et al. 1993Manner et al. , 1995. Additionally, in the final looping events (stage 18-24), the outflow tract undergoes a rotation and migration as it moves from a more rightward to centrally located position (Manner 2000). Together, these events may change the orientation of the outflow tract relative to the aortic arches over the course of cardiac looping. The orientation and position of the outflow tract continue to shift as the spiral aortopulmonary septum septates the common trunk into the mature pulmonary artery and aorta (de Fig. 1 Schematic representation and timeline of the derivatives of the embryonic chick aortic arch vessels. Schematics are viewed from a ventral perspective. The left sketch depicts a virtual six arch pattern, and the right sketch is the mature great vessel configuration. Gray lines represent vessels that involute. The timeline at the center provides a simplified duration that individual arches maintain a connection to the descending aorta. Note that at any individual time point, there may be 2-4 arches present. Right lateral bars are darker. Schematics are modeled after Backer and Mavroudis (1997), OT outflow tract la Cruz et al. 1977;Dor and Corone 1985;Thompson et al. 1987). Arrested migration and rotation of the outflow tract by mechanical (Bremer 1928;Dor and Corone 1985;Gessner 1966) and genetic (Bajolle et al. 2006;Liu et al. 2002;Yashiro et al. 2007) intervention produces both outflow tract and a aortic arch abnormalities and has been implicated as a possible cause for several congenital defects in humans Goor and Edwards 1973;Lomonico et al. 1988). These studies suggest that the positioning of the outflow tract may be a critical factor in the growth and selection of the aortic arches.
Extensive studies using the chick embryo model have shown that epigenetic perturbations in blood flow can also lead to altered cardiovascular morphogenesis and congenital defects, implicating hemodynamics as a major environmental, epigenetic factor (de Almeida et al. 2007;Gessner 1966;Hogers et al. 1999;Lucitti et al. 2005;Reckova et al. 2003;Sedmera et al. 1999;Tobita et al. 2005). This relationship has been validated for the global growth of the aortic arches by several studies, where changes in aortic arch perfusion by occlusion of a vitelline vein (Hogers et al. 1999;Rychter and Lemez 1965), left atrial ligation (Hu et al. 2009), manipulation of the outflow tract (Gessner 1966), or ligation of individual aortic arches (Rychter 1962) generated abnormal arch growth and defects in arch derivatives. Altered flow patterns resulting from a venous clip model were also shown to affect the temporal and spatial expression of genes such as endothelin-1, Krüppel-like factor-2, and endothelial nitric oxide synthase within the aortic arches (Groenendijk et al. 2004;Poelmann et al. 2008). These studies have established that hemodynamic loading such as wall shear stress (τ w ), which is sensed by endothelial cells, acts to initiate signal transduction pathways leading to the regulation of vessel growth (Culver and Dickinson 2010;Groenendijk et al. 2004;Poelmann et al. 2008). Though the exact molecular mechanisms involved in the relationship between aortic arch hemodynamics, gene expression, and tissue remodeling are not yet fully understood, τ w is believed to be the primary mechanical stimulus for vascular remodeling (Culver and Dickinson 2010;Rodbard 1975), clinically referred to as the "flow-dependency principle" (Kamiya and Togawa 1980;Lu et al. 2001;Rodbard 1975). The temporal and spatial events associated with flow sensitive altered aortic morphogenesis (at the first aortic arch) in the alk1 mutant zebrafish model can be more complex as illustrated by recent studies Corti et al. 2011;Patrick et al. 2011).
Prolonged exposure to altered blood flow results in the normalization of τ w via an increase or decrease in vessel caliber (Bayer et al. 1999;Girerd et al. 1996;Langille and O'Donnell 1986). This phenomenon has been explained by the mechanical hyper-restoration theory, which postulates that biological systems tend to grow and remodel to return to a homeostatic state (Beloussov 2008;Beloussov and Grabovsky 2006;Taber 2009). A variety of numerical studies have employed this theory to formulate restorative stressbased growth laws to model vascular growth and remodeling (Alford et al. 2008;Fung 1990;Rodriguez et al. 1994;Taber and Humphrey 2001;Valentin et al. 2011). These studies have been primarily performed in adult, diseased states (Figueroa et al. 2009;Watton et al. 2009) and postnatal growth of the aorta (Wagenseil 2010). Limited embryonic applications (Taber and Eggers 1996) were decoupled from the hemodynamic environment and only examined a straight, axisymmetric vessel segments.
Lack of sufficient data on the micro-structural and mechanical properties of the embryonic aortic arches has limited the direct application of available growth rate laws (Rodriguez et al. 1994;Taber 1998a) and micro-structurally motivated phenomenological models (Gleason and Humphrey 2004;Humphrey 1999) to arch morphogenesis. An alternative approach is to employ an optimization-based growth strategy to model vessel growth in complex cardiovascular morphologies in response to secondary anatomical alterations and sustained flow and pressure perturbations. Theoretical (Sherman 1981;Taber 1998b;Zamir 1977) and experimental (Gruionu et al. 2005;Kassab and Fung 1995;LaBarbera 1990) studies have provided strong evidence that vascular morphology is governed by a set of physiological principles that tend to optimize cardiovascular performance. The most widely regarded is Murray's law, which states that minimization of the total energy to drive blood flow and maintain the metabolic cost of blood volume determines vascular morphology (Murray 1926). Previous research has shown that minimization of energy is achieved through restoration of an optimal τ w (Kassab and Fung 1995;Zamir 1977), which is dependent on local pressure (Pries et al. 1995;Taber 1998b). This relationship suggests that perturbed hemodynamic loading drives vascular remodeling, potentially resulting in diseased conditions. While Murray's law generally holds for the patterning of vascular trees involved in bulk transport, other models have been proposed for the optimal design of microvasculature, which seek to balance energy utilization with diffusive capacity (Kamiya et al. 1987(Kamiya et al. , 1993Khanin and Bukharov 1994). Although these models have hinted the feasibility of a predictive optimization-based growth numerical framework, their zero-or one-dimensional formulation inherently lacked the capability for accurate evaluation of hemodynamic loading, effect of dynamic anatomical variations, and functional requirements of neighboring vasculature, which may also dictate the local vessel growth (Zamir 1977).
In our previous work, the flow patterns and hemodynamic loading of developing aortic arches were generated through micro-computer tomography (micro-CT) 3D reconstructions of stage 18 and 24 chick embryos and computational fluid dynamics (CFD) analysis ( Fig. 2)  ). The results demonstrated that the distribution of cardiac output to the individual arch vessels shifts significantly between stages 18 and 24. Qualitative assessments based on the 3D reconstructions reveal a change in the angle of the outflow tract relative to the aortic arches across developmental time points, consistent with the studies described above (de la Cruz et al. 1977;Flynn et al. 1991;Manner 2000;Manner et al. 1993Manner et al. , 1995. These results suggested that the change in position of the outflow tract may redistribute cardiac output, leading to differential growth of the aortic arches due to hemodynamic loading. In the present study, we tested the hypothesis that the angle of the outflow tract relative to the aortic sac dictates the local hemodynamic environment of the individual arch vessels and drives differential growth and selection by an optimization principle. We performed in vivo measurements of the intraluminal angle between the outflow tract and aortic sac in stage 18, 21, and 24 chick embryos. We then investigated the influence of this angle on the growth of the stage 18 aortic arches using an optimization-based vascular growth model based on idealized anatomy of the three right lateral arch vessels. To accomplish this goal, we developed an inverse CFD-coupled shape optimization framework to predict the great vessel vascular morphology in silico by evaluating multiple objective functions: hydraulic energy budget and particulate diffusive capacity. Integrating the in vivo arch morphometry data and computational modeling, we forecasted the progression of the embryonic arch development subsequent to stage 18 and compared with the normal progression of arch development. We then evaluated the predicted arch geometries with reference to observed in vivo data and explored the implications of our model for investigating normal and aberrant cardiac morphogenesis. We finally proposed a theoretical model in which the change in outflow tract orientation and a shift in objective function during development can guide the selection and growth of the aortic arches. The model geometry used in the CFD-coupled optimization. a 3D polymeric cast of the stage 18 aortic arches, viewed from a dextro-ventral perspective. Right lateral aortic arches are labeled II, III, and IV. Color indicates the wall shear stress magnitude at peak flow velocity, ranging from 0 (deep blue) to 36 Pa (red). Adapted from Wang et al. (2009). b Parameterization of the stage 18 right lateral aortic arch model geometry. A 2D representation of the right lateral aortic arches was generated from fluorescent dye microinjection and micro-CT reconstructions (see text for details). Arch diameters (∅ II , ∅ III , ∅ IV ) served as the optimization shape variables, with θ OT as an additional anatomical parameter. The inset is a right lateral perspective of the 3D reconstruction for comparison with the 2D model. c A flow chart of the CFD-coupled optimization framework used in this study. The framework returns the optimal diameter set for the θ OT specified at the start. Abbreviations are defined in the text

Fluorescent dye micro-injection
Fertilized white Leghorn chick eggs (Eicher's Family Farm, Wexford, PA) were incubated blunt end up at 37 • C and 50-60% relative humidity to stage 18, 21, and 24 (3, 3.5, and 4 days of incubation, respectively) (Hamburger and Hamilton 1951). Between stages 18 and 24, the cranial arch pair II gradually disappears and the more caudal arch pair VI emerges (Hiruma and Hirakow 1995;Wang et al. 2009). Additionally, between these stages, cardiac looping progresses and the outflow tract moves from a more rightward to a central position (Manner 2000). As the outflow tract is oriented cranially and toward the right of the embryo during the investigated stages, Fig. 3 The angle of the outflow tract relative to the aortic sac was measured using fluorescent dye injection. a Sketch of a stage 24 chick embryo. Arterial vasculature is red and venous vasculature is blue. Green represents neuronal structures, yellow digestive. b-d Representative fluorescent dye injection and angle measurement of stage 18, 21, and 24 chick embryos, respectively. Embryos are viewed from the right lateral; cranial is oriented toward the top. The proximal end of outflow tract shifts cranially during these stages of embryonic development, increasing the angle. a Adapted from Patten (1920). Scale bars correspond to 1 mm only right side up embryos were utilized and all measurements were performed for the right lateral aortic arches. Eggs were positioned under a stereomicroscope (Leica MZ16F and M165FC, Leica Microsystems GmbH, Germany), and the embryo was exposed by windowing the shell and removing the overlying membranes. Embryos that were dysmorphic or exhibited overt bleeding were rejected. Fluid microinjection was performed using a pulled and beveled glass capillary pipette (35 µm inner diameter, 30 • bevel) connected via Silastic tubing to a disposable 22-gauge needle and a 10 µl glass Hamilton syringe, mounted on a three-axis micromanipulator. Fluorescent dye was injected at the ventricular apex and visualized using an external fluorescent light source (Leica EL6000, Leica Microsystems GmbH, Germany). Water soluble dyes Cy3.18 (excitation wavelength 550 nm, emission wavelength 565 nm) and Cy5-diethyl (excitation wavelength 646 nm, emission wavelength 666 nm) diluted in phosphate buffered saline were chosen as most suitable to visualize flow streams (Mujumdar et al. 1993). Typical microinjection volumes ranged from 0.4 to 2 µl. Time-lapse movies of the flow streams through the outflow tract and aortic arches were recorded using a low light monochrome digital camera (Fig. 3, Leica DFC 350FX, Leica Microsystems GmbH, Germany).

Measurement of the outflow tract-aortic arch angle
Still images were extracted from the fluorescent injection movies to measure the angle between the centerline of the distal region of the outflow tract and the plane of the aortic sac.
Using an ad hoc image processing tool developed in MAT-LAB (Mathworks, Natick, MA), the centerline of the outflow tract was established by tracing the course of dye and the plane of the aortic sac was determined by fitting a line to the points of brightest intensity at the branching of the aortic arches (Fig. 3). The outflow tract angle was defined as the caudal-most angle formed by the two intersecting lines. A total of 63 injections and measurements were performed (n = 23 stage 18, n = 23 stage 21, and n = 17 stage 24). A student t test was performed to determine the statistical significance of the variation in the outflow tract position between the stages.
We performed a detailed analysis of micro-CT 3D reconstructions to justify the use of fluorescent injection and validate our methodology used to measure the outflow tract angle. All fluorescent dye injections were performed in vivo, and although embryos were viewed from the right lateral, the orientation of the dorso-ventral and cranio-caudal axes may not be parallel to the microscope stage. The position of the outflow tract, aortic sac, and aortic arches was therefore not uniform across all samples. To determine an acceptable method for measuring the outflow tract angle relative to the aortic arches, we performed an analysis using 3D reconstructions of the stage 18 arch manifold. A 3D polymeric cast of the stage 18 aortic arches was created and scanned using micro-CT as previously described ). The 3D cast was imported into Geomagics (Geomagic Inc., Durham, NC), and snapshots were obtained at various dorso-ventral and cranio-caudal orientations. A reproducible measurement of the outflow tract to aortic arch angle could be obtained by representing the aortic sac as a plane and calculating the relative angle of the outflow tract. The aortic sac plane could be successfully established by connecting the branch points of the aortic arches. Using this measurement method, we found that the range of measurement was within one standard deviation of the in vivo angle, supporting the independence of our measurement technique from the in ovo orientation of the embryo (Online Resource 1).

Model construction and parametrization
Our CFD-coupled shape optimization framework was applied to an idealized 2D representation of the right lateral stage 18 aortic arches (II, III, and IV). The bulk shape and dimensions were extracted from the previous fluorescent dye injections and 3D micro-CT reconstructions and represented as a parametric geometry, extending from the distal region of the outflow tract through to the dorsal aorta ( Fig. 2)  ). The outflow tract was modeled as a straight tube with length (L OT ) and diameter (∅ OT ) of 250 and 50 µm, respectively. The dorsal aorta was modeled a straight, tapered tube with length (L DA ), cranial diameter (∅ CRA ), and caudal diameter (∅ CAU ) of 2,195, 112.5, and 225 µm, respectively, with the ventral wall at an angle of 60 • to the horizontal. From the distal most point of the centerline of the outflow tract, the branch point of each aortic arch was defined by a vector with length (L II = 150 µm, L III = 135 µm, L IV = 134 µm) and angle (α II = 18 • , α III = 80 • , α IV = 143 • ), measured with respect to the vertical. The branch point is the most proximal point on the centerline of the aortic arch, and together the three arch branch points define the plane of the aortic sac, as in the in vivo measurement. The end point of each arch, defined as the distal most point of the arch centerline, was chosen such that the distance from branch point to end point of arch III was 610 µm, and the line connecting these two points formed an angle of 100 • with respect to the plane of the aortic sac. These values were determined based on measurements from the fluorescent dye injections and micro-CT reconstruction. The end points of arches II and IV were defined by their distance from that of arch III, measured along the line of the ventral side of the dorsal aorta.
Third-order Bezier curves were used to connect arch branch points to end points, creating the centerline of each arch. The proximal anastomosis angle (ϕ p ) and length (ψ p ) defined the tangent vector at the aortic sac (ϕ p, , and the distal anastomosis angle (ϕ d ) and length (ψ d ) defined the tangent vector at the dorsal aorta (ϕ d,II = 130 • , ϕ d,III = 95 • , ϕ d,IV = 95 • , ψ d,II = 450 µm, ψ d,III = 400 µm, ψ d,IV = 400 µm). Angle ϕ p was measured with respect to the line connecting the distal most point of the centerline of the outflow tract to the arch branch point, and angle ϕ d was measured with respect to the line of the ventral side of the dorsal aorta. The arch lumen was formed by extruding the vessels walls in the normal directions along the centerline. Third-order Bezier curves were used to make smooth connections to the outflow tract, dorsal aorta, and between the arches. Values for the fixed shape parameters described here were manually determined by trial and error to obtain a qualitative agreement with right lateral views from the fluorescent dye injections and micro-CT reconstructions. ∅ OT was chosen to maintain an inlet Reynolds number (Re) of 4, as observed in our previous stage 18 CFD simulations ). The 2D model in Fig. 2 can be compared to the inset right lateral view of the micro-CT reconstruction in the same figure and to the images from the fluorescent dye injections presented in Fig. 3 to observe its similarity to the right lateral arches.
During the shape optimization, the individual arch diameters, ∅ II , ∅ III , and ∅ IV , served as the inverse shape design variables. The angle of the outflow tract relative to the aortic sac (θ OT ) was an additional input parameter. The generation of the geometric models was automated using an inhouse subroutine that accepts the design variables as inputs.
The geometry was created using FEMLAB (COMSOL Inc., Burlington, MA) internal shape commands, which allowed for seamless interaction between the model construction, meshing, and flow simulation. To translate the 3D tubular geometry to 2D, the input diameters were treated as hydraulic diameters, and the corresponding width of a rectangular channel was assigned as the model diameter. Model parameterization and translation of a 3D geometry to 2D has been implemented previously by our group . Limitations of our model related to the simplification to 2D are discussed later (see Sect. 4.6).

Numerical simulation of blood flow and species transport
The 2D incompressible Navier-Stokes and passive-scalar transport equations were modeled over the right lateral arch domain using the commercial finite element (FE) solver FEMLAB (COMSOL Inc., Burlington, MA). This solver uses the Petrov-Galerkin formulation and has been validated for 2D flow (Hanke 2004). A steady-state solution was chosen as our previous pulsatile CFD simulations revealed that the Womersley number was less than 1 within the aortic arches, suggesting that a steady solution is sufficient ). Flow was modeled using a rigid walls assumption, and blood was treated as a non-Newtonian power-law fluid, μ = mγ n−1 , where μ is the apparent viscosity,γ is shear rate, m = 0.017 Pa n , and n = 0.708 (Shibeshi and Collins 2005). Although the arches may distend during flow, our past experience indicates that wall motion is minimal and would be a negligible factor in this study. We use viscosity parameters for adult blood as there is limited data for the embryo. A steady parabolic flow profile was specified as the inlet boundary condition, with the average velocity corresponding to one-half the average flow rate over the cardiac cycle for the stage 18 chick (0.44 mm 3 /s) (Yoshigi et al. 2000), since only one side of the symmetric aortic arches was modeled. The outflow tract diameter and velocity profile maintain an inlet Re of 4, which was measured from our previous work ). A pressure boundary condition was specified at both the cranial and caudal dorsal aorta outlets. The reference pressure at the caudal end was 120 Pa (0.9 mm Hg), equal to measured dorsal aorta pressures (Hu and Clark 1989;Yoshigi et al. 2000). The cranial pressure was then adjusted during the each inner iteration automatically to provide the desired physiological caudal/cranial flow split of 90/10. Diffusion of a passive species (generic paracrine factors or oxygen) was subsequently modeled, using the velocity field from the flow simulation to solve the convective diffusion equation. A unit normal concentration of 1 mm −3 was supplied as a boundary condition along the walls of the arch vessels, with zero concentration at the outflow tract and a zero diffusive flux boundary condition on all other boundaries. A diffusion coefficient (D) of 1.0 ×10 −4 cm 2 /s was set as constant. This large diffusivity was chosen so that both energy and diffusion were of the same order of magnitude to simplify the optimization routine. The geometry was discretized using an unstructured tetrahedral mesh, consisting of approximately 4,000 elements for grid independency with a minimum quality > 0.5. Six levels of grid refinement were used to verify grid independency of the objective functions.
The model construction, meshing, CFD simulation, and post-processing were all automated in a single modularized function, which accepts inputs of the design variables and outputs the objective function. The utility of coupling CFD to inverse shape optimization was demonstrated previously by our group as a pre-surgical planning tool to improve patient-specific coronary artery bypass graft designs ).

Multi-objective optimization problem
The automated CFD framework described above was coupled to a multi-objective shape optimization algorithm to achieve the optimum set of arch diameters over a range of θ OT . The two-element objective vector consisted of the total power to sustain blood flow (P T ) and the total diffusive exchange capacity (J ) of the arch domain. Minimization of energy is an established principle of vascular form (Gruionu et al. 2005;Kassab and Fung 1995;LaBarbera 1990;Murray 1926;Pries et al. 1995;Sherman 1981;Taber 1998b;Zamir 1977), while maximizing diffusive capacity is derived from the homology of the aortic arches in lower order vertebrates where terminal arch pairs are the precursor of gill vasculature (see Sect. 4.4).
The total power to sustain blood flow in the system (Eq. 1) is comprised of the linear combination of three components: the power to (1) drive the flow of blood (P f ), (2) maintain the metabolic cost of blood volume (P b ), and (3) maintain the metabolic cost of the vessel wall (P w ) (Taber 1998b). Formulas for these three components are given in Eqs. (2-4), respectively, where u is the velocity field, p is the pressure, ρ is the density, ∅ is the vessel diameter, L is the vessel length, and T is the wall thickness. The subscripts ot, cra, and cau stand for the outflow tract inlet and cranial and caudal dorsal aorta outlets, respectively. (1) In Eq. (3), α b is a metabolic constant for blood, while α w and β w in Eq. (4) are the passive and active metabolic constants for the wall, respectively. Values for these constants vary significantly among species and vessel types; in Zamir (1977), α b has a range of 5-5,000 W/m 3 . For all calculations in this study, values of 1,700, 4,000 W/m 3 , and 0.433 s −1 were used for α b , α w , and β w , respectively. In Eq. (4), f , the active fraction of the total fiber stress, was assumed to be 1 and T /∅, the ratio of wall thickness to diameter, was fixed at 0.05. Diffusive capacity was quantified as the total diffusive uptake between the inlet at the outflow tract and dorsal aorta outlets (Eq. 5), where c is the concentration and other notations are as above. As a unit normal concentration was applied along the arch walls, J has a maximum value of 1.
Optimization was performed using the sequential quadratic programming (SQP) method (MATLAB, Mathworks, MA). SQP achieves efficient evaluation of constrained optimization problems by attempting to solve for the Lagrange multiplier directly (Schittowski 1985). At each major iteration, the search direction is obtained by solving a quadratic sub-problem based on quadratic approximation of the Lagrangian. The following optimization problem was posed: where: subject to: In the above problem, s denotes the shape design variables (the arch diameters) and G (s) and F(G(s), s) refer to the CFD solution (where u is the velocity field, p is pressure and c is concentration) and multi-objective functions, respectively. Weighted global criterion scalarization (Collette and Siarry 2003;Marler and Arora 2004) was used to solve the multi-objective problem (Eq. 7), where U (F) represents the scalar objective function and the reference point P o T , J o was chosen such that P T is minimized and J maximized. Optimization is subject to upper and lower bound constraints on the design variables, denoted as s u and s l , respectively. A lower bound of 10 µm arch diameter was enforced rather than implementing arch disappearance in order to maintain a continuous objective function and an upper bound of 180 µm was given to avoid errors in overlapping arch vessels. We interpret an arch converging to the minimum allowed diameter as disappearing.

Framework algorithm
The framework algorithm is summarized in Fig. 2. First, θ OT and the initial state s are used to generate the starting geometry. Model construction, meshing, FE simulation, and post-processing are run in series to evaluate the objective functions at the current diameter set s. This function evaluation loop is executed for sequential changes in each diameter to determine the gradients of the scalar objective function U . Gradients are used to compute the step size, and the updated diameters are fed to the function evaluation loop, repeating the above process. The optimization loop terminates when a minimum is located, returning the optimal diameters for the θ OT defined at the start. For the present study, a ∇U of less than 10 −7 was sufficient for optimization. The entire process was repeated for each θ OT .
To assess the effects of the relative contributions of energy and diffusion on arch morphology, optimizations were performed for three P T /J weighting ratios: 50/50, 65/35, and 100/0. The final ratio, 100/0, refers to a mono-objective optimization, where only total energy utilization is minimized. In the remaining text, we refer to the multi-objective cases as "power + diffusion optimization" and the mono-objective as "power-only optimization". We also tested the effect of the caudal/cranial outlet flow split by conducting the 65/35 case with a 50/50 prescribed split. The simulations with different flow splits did not produce significantly different results (data not shown).
We performed preliminary simulations to examine the behavior of each objective function and determine our initial condition. The diameter of one arch was systematically incremented while keeping the other two at a fixed, equal value. Five diameters were tested, ranging from 10 to 150 µm, and the process was repeated for the same five diameters applied to the two fixed arches (Online Resource 2). The P T gradient was steep for small diameters but quickly decayed to approach a shallow minimum. This minimum is the result of the inverse relation between the gradients for P f and P b and P w ; the larger P f gradient at small diameters drives an increase in arch caliber, which is eventually suppressed due to the larger blood and vessel wall volume. Over the [s l , s u ] interval in this study, the effects of P f dominate and arch vessels approach or reach the maximum allowed diameter to achieve minimum energy. The gradient for J behaved similarly, though not as steep as P T . The maximum value of J is limited to 1, given the unit normal concentration at the arch walls, which is approached as all three arch diameters near zero. Thus, no maximum is obtained on the given [s l , s u ] interval, and J acts to reduce vessel caliber. We ran a single-variant optimization, with all arch diameters equal, for θ OT of 97 • and P T /J of 50/50 and used the resulting diameter (95 µm) as the initial condition for all arches. We tried several other initial conditions, including two small caliber arches and one large and two large arches and one small, and found that starting with three arches of equal diameter consistently converged to a global optimum and that the 95 µm uniform arch diameter tended to converge with the fewest function evaluations. In addition, several tolerances for the ∇U required for termination of the optimization were tried, and we found that stricter values (i.e., <10 −7 ) did not significantly influence the final diameters.

In vivo angle measurements
Embryonic heart looping can significantly alter the distribution of blood flow through the multiple pairs of aortic arch vessels. Quantification of the orientation of the outflow tract relative to the aortic arches using fluorescent dye injection revealed an increase in the angle as the proximal end shifted cranially (Fig. 3). Mean angles ± standard deviation were calculated as 96 • ± 3 • for stage 18, 116 • ± 6 • for stage 21, and 117 • ± 6 • for stage 24. A student t test indicated a statistically significant difference ( p < 0.001) between stages 18 and 21, but not between stages 21 and 24, consistent with the fact that the cardiac looping process is largely completed by stage 21.
Morphogenetic events from stage 18 to 24 include movement of the outflow tract from a right lateral to a more centrally located position (Manner 2000;Patten 1920). Our measurements indicate that this motion causes a further shift in the outflow tract-aortic arch angle as the proximal end of the outflow tract moves toward the cranial pole of the embryo (Fig. 3), as was qualitatively observed earlier in polymeric casts of the embryonic arch system ).

Aortic arch patterning for power and diffusive capacity objective functions
Adaptive growth of the aortic arches was explored using the power + diffusion optimization based on the ratio between P T and J ; analyses were conducted for two P T /J ratios: 50/50 and 65/35. Both cases indicated that three patent arch vessels provide a balance between reducing energy consumption and increasing diffusion. For the P T /J = 65/35, θ OT was a determining factor in selecting which of the three arches had the largest diameter, i.e., the dominant arch. For smaller θ OT , the dominant arch is predicted as arch II, whereas for larger θ OT arch IV was dominant ( Fig. 4; Online Resource 5). When the contribution of P T and J were equivalent, arch configuration showed minor dependence on θ OT , with similar arch diameters predicted for each angle (Online Resource 3). The dependence of arch patterning on θ OT suggests that the orientation of the outflow tract has a significant role in directing the growth of that specific arch that is in line with the centerline of the outflow tract.

Aortic arch patterning for power minimization only
We also examined the case where minimization of power acted as the sole objective function (P T /J = 100/0). Unlike the multi-arch configurations resulting from a balance between P T and J , minimization of P T alone transformed the multi-arch manifold to a single arch geometry, where the remaining two arches converged to the minimum allowed diameter (10 µm). Further, the selection of the single arch displayed strong dependence on θ OT (Fig. 5; Online Resource 6), similar to the effect of θ OT in the power + diffusion case (Sect. 3.2). Total energy utilization was reduced by 10% compared to the 65/35 P T /J ratio. To verify that global minima of the total energy expenditure was achieved at the single aortic arch state, independent of the initial arch configuration, two single-variant optimization problems were conducted with (1) three arches of the same diameter, and (2) a single arch. For the single arch configuration, the design parameters and variables remained unchanged, with the exception that the diameters of arches II and IV were held to 0. Thus, the geometry of arch III was used to model the single arch configuration. A single variable, s = ∅ arch in (Eq. 13), was used to model these two cases. Fig. 4 Power + diffusion optimization with P T /J of 65/35 converged to a multiple arch configuration; however, the dominant arch was highly dependent on θ OT . a The optimal diameter set for each θ OT for P T /J of 65/35. Arch dominance changed from arch II for θ OT < 58 • , to arch III for 58 • < θ OT < 95 • , to arch IV for θ OT > 95 • . All arches remained patent. b Behavior of P T and J over all optimizations. Labels indicate the mean (SD). c and d The velocity and concentration fields, respectively, for θ OT of 102 • and P T /J of 65/35 Minimizing total energy (P T /J = 100/0) demonstrated that, while both cases led to the same total cross sectional area, the optimal single arch case reduced energy utilization by 12% compared to the optimal multi-arch case (Online Resource 4). This relation is consistent with fluid mechanics, which indicates that for equal cross sectional area, the total resistance of three vessels in parallel is greater than that of a single vessel, leading to larger pressure drops.
The optimization histories for sample 65/35 and 100/0 P T /J ratios are shown in Fig. 6. In both cases, there was an initial period of large diameter changes, lasting for 30-40 iterations. After this period, the diameters tended to stabilize, remaining relatively constant. The initial period can be characterized as a search for the dominant arch, with the later iterations determining its optimal diameter. Although the final minimum is shallow, the initial search is associated with steeper gradients, which only level off once the dominant arch is found. This trend suggests that selection of the dominant arch is a strong impulse, while subsequent changes to its diameter are less well defined. Given the SQP method used, it is possible that the dominant arch selected by our model was a saddle point and not the true global minimum. Based on our results, we expect three such saddle points, with each arch as being dominant. To verify independence from the initial condition we also performed optimizations for a subset of θ OT , where we began with an initial condition of one arch at the maximum diameter and the remaining two at the minimum. The largest arch was varied so that all three possible initial conditions were tested. In most cases, the geometry converged to the same as our original results. There were a few cases, however, where the initialized dominant arch did not change. Comparing to our original results revealed these to be saddle points, and in all cases, our original initial condition converged to the global optimum.
The results of the power + diffusion and power-only optimizations demonstrate that maximizing diffusive capacity drives smaller and smaller diameters, while minimizing power results in one large vessel. Imposing both objectives produces a compromise configuration as the influence of diffusive capacity prevents any arch from constant enlargement and energy expenditure inhibits diameter reduction. As the relative importance of one objective increases, its effects are amplified, as seen in the 50/50 and 65/35 P T /J cases. We also performed power + diffusion optimizations for 70/30, 75/25,  b and d). Once the optimal dominant arch is found, arch diameters stabilize as the minimum is approached and 80/20 ratios for a subset of θ OT . We found that, for the vast majority, the single arch configuration was optimal. These additional ratios demonstrate that a single arch is required to reduce the total power below a certain value. We performed optimizations only for these specific P T /J ratios to maintain a manageable number of FE simulations. The 65/35 ratio was reached by manual trial and error to match our previous stage 18 diameter data ).
Both the power + diffusion and power-only optimizations converged to within a small range of P T and J for each θ OT (Figs. 4,5). This trend suggests that a global optimum may exist and that the arch morphology attempts to reach this value by re-patterning the aortic arch diameters as θ OT changes. The adaptive growth as a response to hemodynamic perturbations in order to restore a preferred (optimal) state is previously described as the hyper-restoration theory (Beloussov 2008;Taber 2009). These findings therefore support the idea that a change in the angle of the outflow tract can provide a stimulus for arch development.

Orientation of the outflow tract determines the dominant aortic arch in the developing embryo
We found that, for all optimal arch configurations, one arch vessel must be larger in diameter than the other two, which we designate as the dominant aortic arch (Figs. 4, 5). Further, the selection of the dominant arch was dependent on θ OT , suggesting a significant role of vessel orientation in arch patterning. According to our previous ) and current in vivo measurements and optimization-based growth predictions, a P T /J ratio of 65/35 allows good correlation between the model and observed arch morphogenetic events, with predicted arch diameters within one standard deviation of the in vivo data reported in our previous work ). Arch dominance changed from arch II to arch III and finally to arch IV as θ OT increased (Fig. 4). The dominant arch III predicted for θ OT near the measured angle of stage 18 (96 ± 3 • ) is consistent with our previous in vivo measurements ). Minor deviations from quantitative agreement between measured diameters and angles may be associated with modeling the 3D symmetric anatomy in one lateral (see Sect. 4.6). As θ OT increased during development, the dominant arch shifted to arch IV, which is supported by previous studies on aortic arch growth (Hiruma and Hirakow 1995) ( Fig. 4; Online Resource 5). Transitions in arch dominance occurring across the stage 18, 21, and 24 angles is consistent with the dynamic arch patterns observed at stage 21, the intermediate stage between the II, III, IV configuration of stage 18 and III, IV, VI state of stage 24. The excellent correlation between the model and experimental data for the P T /J ratio of 65/35 suggests that at stage 18, the aortic arches balance energy utilization and diffusive capacity, similar to the optimized capillary tree network structure of microcirculation (Kamiya et al. 1987(Kamiya et al. , 1993. When minimizing energy alone (P T /J = 100/0), the change in arch dominance as a function of angle change was sharper, creating a single arch geometry as the remaining two vessels regressed completely (Fig. 5). In this trend, two configurations were observed: arch II dominance and arch IV dominance, with the transition occurring at around θ OT = 84 • (Fig. 5; Online Resource 6). The dependence of the dominant arch on θ OT seems to be influenced mainly by the minimization of energy, consistent with the hypothesis that energy dissipation at the aortic sac drives the selection of the dominant arch. From the fluid mechanics perspective, the rapid expansion of the narrow outflow tract into the aortic sac generated a jet flow regime that propagates several vessel diameters downstream toward the aortic arches (Figs. 4, 5) ). The angle of the outflow tract may dictate the pathway of the jet flow within the aortic sac and direct cardiac output preferentially toward one of the aortic arches. This jet direction can result in selection of the dominant aortic arch that lies on the flow pathway and thus dissipates the minimum energy. Although flow should be directed toward arch III for θ OT ∼ 90 • , a dominant arch III is never observed when minimizing energy alone. This may be due to the low Reynolds flow (inlet Re = 4) that creates a tendency for blood to adhere to the walls of the outflow tract and aortic sac and thus flowing toward the outermost arches II and IV. Studies at higher Re numbers are needed to verify this idea.

Change in objective function and competition among parallel arches may drive convergence to a single aortic arch
The results of the power + diffusion and power-only optimization pathways reveal that embryonic arch selection can be characterized similar to competing collateral arterial vessels. Examining a simple network of two parallel vessels, Keenan and Rodbard demonstrated that vascular adaption to τ w alone causes one vessel to enlarge and the other to disappear (Keenan and Rodbard 1973). In their model, vessel caliber increased or reduced to restore "optimal" τ w conditions, effectively minimizing total energy. Small perturbations in the distribution of blood flow leads to the dilation of one vessel due to increased τ w and constriction of the other due to reduced τ w . The change in resistance caused more flow to the larger vessel and less flow to the smaller, creating positive feedback loop, which eventually led to a single vessel. This dynamic vessel competition has also been demonstrated in the adaptation of larger microvascular networks, where τ w maintenance led to rarefaction to a single vessel (Hacking et al. 1996;Hudetz and Kiani 1992). Our methodology for the power-only optimization follows these previous models, as the number of aortic arches reduces to one. In our model, the change in θ OT is the source of the flow perturbation, preferentially selecting the dominant arch.
Other dynamic models of microvascular adaptation have shown that applying metabolic controls in addition to τ w can result in stable, multi-vessel networks that closely match the in vivo topology. Using a 1D mathematical model, Pries et al. introduced metabolic demand for O 2 as an additional parameter for the control of vessel diameter (Pries et al. 1998). When applying both τ w and metabolic controls, a stable network was achieved, whereas τ w alone again resulted in instability and rarefaction. The target τ w level in this model was not constant as in Keenan and Rodbard (1973), but varied according to transmural pressure, as suggested by Pries et al. (1995) and incorporated into our model through vessel wall maintenance (Taber 1998b). Applying their model to more complex networks demonstrated that convergence to the in vivo topology could only be obtained by enforcing both τ w and metabolic controls (Pries et al. 2001). Their results indicated that the functional requirements of vascular systems are of equal, and possibly greater, importance as τ w in regulating vessel diameter. This multi-stimuli approach is similar to our power + diffusion model, which generated a stable, multi-arch configuration. Our CFD-coupled optimization also revealed that the higher dimensional geometry and flow field can influence vessel morphology. Minimization of energy accurately predicted that the multi-arch geometry will converge to a single vessel, derived from arch IV (Barry 1951;Hiruma and Hirakow 1995). While the reduction to a single vessel to reduce energy is a wellestablished phenomenon (Keenan and Rodbard 1973), our novel CFD-coupled optimization-based model demonstrates that θ OT is an important factor in selecting the dominant arch due to its role in distributing blood flow. Our power+ diffusion simulations suggest that the multi-arch configuration is due to functional objective functions such as metabolic demand. Based on these results, we propose a theoretical model in which the change in outflow tract orientation and a shift in objective function during development guide the selection and growth of the aortic arches (Fig. 7). However, it is also possible that the emergence of multiple arches is the result of phylogenetic selection, driven by the genetic code. Later convergence to the single vessel may then simply be the adaptation of the aortic arch network to minimize energy. As new arches arise due to inherent genetic expression, those that are placed in a more advantageous position due to θ OT (i.e., arch IV) enlarge at the expense of the already existing arches. This type of dynamic instability is characteristic of a competition among the arches, with θ OT dictating the "winner". Additional research into the formation of new arches will be required to address these issues.

Abnormal orientation of the outflow tract may lead to defects in aortic arch growth and selection
Classical descriptive morphologic studies have long suggested that the migration and rotation of the outflow tract is a critical factor in great vessel development (Bremer 1928;Congdon and Wang 1926) and errors in outflow tract location and orientation are frequently associated with congenital heart disease. In his experiments on chick embryos, Bremer speculated that the "descent" of the heart into the thorax results in significant torsion of the outflow tract, which is translated to the aortic arches, causing obliteration of some arches due to mechanical force (Bremer 1928). This torsion was explicitly implicated in the growth of arch pair IV, where the counter-clockwise axial rotation of the outflow tract causes twisting and bending of the left lateral arch IV, leading to its obliteration. Disruption of the normal leftward rotation of the embryo and dextral looping of the heart by various surgical manipulations resulted in abnormal patterning of the arches, including a persistent left arch IV and aberrant derivatives of arch III. Similar results were obtained by Gessner (1966), where the outflow tract in stage 19 chick embryos was displaced and made immobile by inserting a surgical wire. Arch defects involved deviations from the normal arch derivatives, such as disappearance of arch IV and the definitive arch of aorta formed by arch III (both laterals were observed), absence of the left lateral of arch VI with the left pulmonary artery formed from the right lateral, and a double aortic arch. Our proposed model is consistent with these findings in the literature and further establishes the role of outflow tract orientation in selecting the dominant arch of aorta. We anticipate such surgical manipulations of outflow tract position to alter the hemodynamic environment and selection of the optimal dominant arch, generating these abnormal patterns.
While substantial progress has been made in defining critical genes and proteins associated with normal outflow tract and aortic arch morphogenesis, the molecular mechanisms that determine changes in cell proliferation, differentiation, and matrix remodeling involved in asymmetric arch selection and remodeling have yet to be identified. Targeted gene studies in mice have shown that splotch and Pitx2 are necessary for normal outflow tract rotation (Bajolle et al. 2006;Liu et al. 2002). Defects in these genes result in arrested rotation and produce arch phenotypes similar to those observed by Bremer, including double aortic arch and randomization in laterality of the dominant arch of aorta. However, Pitx2 positive cells are not found in the aortic arches themselves, indicating that it does not have a direct role in asymmetric arch morphogenesis (Liu et al. 2002;Yashiro et al. 2007). Rather, its is believed that the arrested outflow tract rotation in Pitx2 mutants disrupts the distribution of cardiac output to the aortic arches, leading to abnormal hemodynamics that alter shear stress sensitive signaling pathways such as PDGFR and VEGFR2 (Yashiro et al. 2007). While these studies have primarily examined the lateral asymmetry of arch pair IV, our model suggests that a similar mechanism may exist in the selection of the dominant arch from the set of all six arch pairs.
Though not directly related to outflow tract orientation, several other studies have validated the effects of altered flow patterns on aortic arch malformations and growth. Unilateral vitelline vein ligation in the chick embryo results in abnormal intracardiac flow patterns, incomplete looping resulting in double outlet right ventricle, and redistributes flow to the individual arch vessels (Hogers et al. 1999;Rychter and Lemez 1965). These changes in blood flow are associated with a variety of aortic arch defects, including hypoplastic right brachiocephalic artery, interrupted aortic arch, double aortic arch, and hypoplastic pulmonary artery (Hogers et al. 1999). Reducing left ventricle volume by left atrial ligation also causes abnormal perfusion of the aortic arch vessels (Hu et al. 2009). Arch anomalies such as absent arches III and IV and arch hypoplasia have been described. These studies provide ample evidence that irregularities in the outflow tract-aortic arch angle generate altered flow patterns in the embryonic aortic arches and cause defects in the formation of the great vessels.
Many of the aberrations observed in these studies are relevant to known congenital defects. Double aortic arch, caused by persistence of both arch IV laterals and right aortic arch, occurring when the dominant arch of aorta forms from the right-IV arch (normally formed from left-IV in humans), result in a vascular ring compressing the trachea and esophagus, causing upper airway obstruction and dysphagia (Backer and Mavroudis 1997). Interrupted aortic arch is induced by unilateral obliteration of the left arch IV in left arch dominant aortas, is often associated with loss of the right arch IV resulting in an aberrant origin of the right subclavian artery, and is always associated with a ventricular septal defect, in which blood from both ventricles flows through the pulmonary artery and is shunted to the systemic circulation through a patent ductus (Sadler and Langman 2006). It is worth noting that tetralogy of Fallot, transposition of the great arteries, and double outlet right ventricle associated with heterotaxy syndrome can be associated with an abnormal aortic-to-pulmonary valve axis, which occurs as the result of morphogenetic errors in outflow tract orientation Lomonico et al. 1988). Aortic arch hypoplasia is normally associated with hypoplastic left heart syndrome, a severe defect requiring multiple palliative surgeries over the patient's lifetime, and has the highest mortality rate of any birth defect in patients less than 1 year of age (Roger et al. 2011). Fetal intervention in humans has shown that restoration toward normal biomechanical loading can improve the growth and remodeling of left heart structures (mitral valve, aortic valve, and aorta), though recovery of fetal myocardial growth after fetal intervention remains suboptimal (McElhinney et al. 2010;Pekkan et al. 2008). Therefore, the optimal timing for intervention to recover normal growth and remodeling of the developing fetal heart and vasculature remains unknown. 4.4 Diffusive capacity as an objective in early development is supported by comparative anatomy of the aortic arches Imposing diffusive capacity as an early primary objective function for aortic arch growth may not be apparent as all diffusion activities occur through the skin during early embryo development. Additionally, the size of the arch arteries, each on the order of 100 µm at stages 18, 21, and 24 ), is too large to allow substantial diffusion and implicates these vessels as conduits for bulk transport. However, our analysis of diffusive capacity as an optimization objective relates to the morphogenesis of the gill arches in lower order vertebrates rather than to suggest a diffusion function for the avian aortic arches at stage 18. The ontogeny of the aortic arches has long been a hallmark of comparative anatomy, demonstrating clear phylogenetic patterns and illustrating the evolution of the arterial system as the gas exchange medium transitioned from water to air (Kardong 2009;Liem et al. 2001). Unlike the chick and human, in fish, the embryonic symmetric pattern of the aortic arches persists into adulthood as these vessels form the afferent and efferent branchial arteries that feed and drain the gills (Evans et al. 2005;Hughes and Morgan 1973). The resulting pattern of multiple branchial arteries, each feeding several hundred filaments, creates a gas exchange organ highly efficient in extracting oxygen from the surrounding water (Hughes and Morgan 1973;Piiper 1982). It is speculated that hypoxic environments drove the transition from water to air as the external gas exchange medium, a transition, which produced changes in aortic arch growth (Evans et al. 2005;Kardong 2009;Liem et al. 2001). It is worth noting that the number of persistent arches as well as gill surface area is reduced in fishes that are obligate air breathers (Brauner et al. 2004;Hughes and Morgan 1973). In addition to gas exchange, the gills provide other vital functions including ionoregulation, acid-base regulation, and the removal of metabolic waste (Evans et al. 2005).
These diffusion-dependent functions of piscine aortic arch derivatives suggest that the pattern of these vessels is meant to maximize diffusive capacity by increasing the number of capillary beds between the heart and dorsal aorta. Although the arch vessels themselves do not actively participate in diffusive exchange, the stimulus to grow and remain patent may be driven by an optimization principle. This same stimulus would then be present in the early embryonic chick arch vessels, decreasing during the course of development as the foundations for aerial respiration form. Thus, the emergence of multiple arches early in development is driven by a balance between energy budget and diffusive capacity, while the later development of the single systemic arch is the result of minimization of energy alone (Fig. 7). This evidence from comparative anatomy initiated our choice of diffusive capac-ity as an objective function and seems to support our model of optimization-based arch growth.

Utility of optimization-based growth models
While principles such as Murray's law of energy minimization during vasculogenesis seem to correlate well with observed anatomy, the fundamental relations that govern embryonic morphogenesis, or more generally biological systems, have not been completely determined (Taber 2009). With regards to the vascular system, it is well established that the local mechanical environment is sensed by multiple cell lineages including endothelial cells, but the translation of this local information with respect to complex, global 3D vascular lumen shape adaptation and cellular responses remains under investigation. Both biomolecular pathways and mechanical factors play an important role in this process, likely governed through complex mechano-genetic feedback. The underlying cardiovascular biology that drives the topological and structural alterations is significantly different in disease, normal embryonic growth/development, and acute/chronic adaptation. For each case, once the governing biology is well understood, the computational framework to predict fluid-tissue interaction can quantify the relative degree of mechanical and genetic contributions to the observed modes of vascular shape alterations. Recently, a number of such growth models with varying degrees of complexity have been applied to predict diseased cerebral aneurysm development (Li and Robertson 2009;Watton et al. 2009) and angiogenesis (McDougall et al. 2002;Sun et al. 2005), but limited studies are devoted to normal embryonic growth response. Local growth rules proportional to the various fluid-induced loading indices, including the τ w sensitivity, need to be tested and compared to in vivo measured 3D surface reconstructions. Linear vessel diameters can be statistically correlated with τ w indices ) and would be a starting point to quantify the complex 3D association.
Through the present optimization-based growth modeling approach presented here and substituting CFD-coupled applied shape optimization for integrative growth laws (Alford et al. 2008;Taber 1998a;Taber and Eggers 1996), we are able to determine gross morphologic variations. Unlike dynamic models of vascular adaption (Pries et al. 1998) and solid mechanics-based growth and remodeling simulations (Gleason and Humphrey 2004;Taber 1998a), our method cannot determine the intermediate steps of the growth pathway. In vivo information of microstructural architecture or employing constitutive biomechanical relationships (Holzapfel et al. 2000;Holzapfel and Ogden 2010) would augment the present approach. The present modeling strategy can also be adapted to incorporate other structurally motivated growth and remodeling principles by posing new objective functions and can incorporate variations in cardiac function through the fluid boundary conditions. Based on established physiologic optimization principles (Baba et al. 1995;Kamiya et al. 1987Kamiya et al. , 1993Murray 1926;Taber 1998b), our model was able to predict aortic arch lumen diameters within one standard deviation of reported values ) and forecast future morphogenetic events including dominant arch selection remarkably well. The performance of our model represents an encouraging first step toward an alternative paradigm in fluid-coupled arterial growth and remodeling.

Energetic cost of metabolic maintenance
When computing the total energy expenditure, calculating both the metabolic cost of blood volume and the metabolic cost of the vessel wall involves the use of metabolic constants (Eqs. 3,4). Reported values for these constants vary greatly for different species (Taber 1998b;Zamir 1977). Although some trends do exist, such as an inverse relation with animal size (Taber 1998b), the choice for the values of these constants was hindered by available literature data. Limited information on the composition and metabolic requirements of embryonic blood and energetic cost of vascular growth prohibit sufficient estimation of these parameters. The values used in our model are within the reported physiologic range; however, future studies would be required to determine precise values.

2D geometry versus 3D aortic arch manifold
Although the model geometry used in this study is based on reconstructions of the stage 18 aortic arches, the simplification to 2D neglects the complex 3D flow and pressure fields of the true anatomical shape ). For example, in the embryo, the aortic arches emerge as paired structures surrounding the foregut, with large arch curvatures, which are not captured in our planar model. Additionally, while we have examined the effects of the cranio-caudal shift of the outflow tract, the lateral movement and axial rotation cannot be investigated in 2D and may play a similar role in shaping the aortic arches (Bremer 1928). Extending our optimization framework to 3D would improve the predictive capabilities and allow further investigation into arch development. However, this will also increase the computational cost quite significantly.

Conclusions
We have shown that hemodynamic changes associated with outflow tract orientation due to cardiac looping can influence the growth and selection of the embryonic aortic arches. Measurement of the angle of the outflow tract relative to the aortic sac revealed a statistically significant shift between stages 18 and 21, a critical point during cardiovascular development. The optimization-based growth model demonstrated that this change in angle re-patterns the individual arch diameters and accurately predicts the mature arch IV dominance over development. A multi-objective optimization that balances energy utilization with diffusive capacity modeled the physiologic stage 18 arch diameters well, while a minimization of energy alone converged to the adult single systemic aortic arch. Additional research will require extension to 3D modeling and investigation of the cellular and molecular mechanisms responsible for maintaining/achieving optimal function.