## Abstract

New necessary conditions for the nonlinear stability of branched tree-like structures composed of elastic rods are presented. The conditions, which are found by examining the second variation of an energy functional, are established by extending Legendre's classical work on this topic. For the branched tree-like structures of interest in this paper, one of the new necessary conditions is the existence of a bounded solution to a set of Riccati equations, which are coupled at branching points. In addition, a set of new necessary conditions are established in the event that the abscissa of the branching point is not fixed. The results are illustrated using a range of examples featuring planar configurations of branched elastic rods where each of the individual rods is modelled using Euler's theory of the elastica.

## 1. Introduction

Tree-like structures formed by long flexible branches are ubiquitous in nature and the mechanics of these, and other plant structures, have been studied for centuries. Of particular interest in this paper is the use of rod theories to provide the mechanical models. These theories have been used by, among many others, Greenhill (1881) to predict the maximum height of a tree, by McMillen & Goriely (2002) to understand the mechanics of tendril perversion and by Silk *et al.* (1982) to model the evolution of a rice panicle. Many recent works which use rod theories to model tree-like structures have been motivated by the need to develop realistic computer graphics-based images of trees and forests. Incorporating biomechanical aspects of plant development, such as growth and the formation of residual stresses, makes this work particularly challenging. Despite these difficulties, many of the resulting images of trees exhibit extraordinary detail and realism (see the works of Prusinkiewicz and his co-workers (Jirasek *et al*. (2000) and Costes *et al*. (2008)), Fourcaud & Lac (2003) and Fourcaud *et al*. (2003)).

One of the difficulties inherent in computing equilibrium configurations for tree-like structures modelled using elastic rods is the algebraic complexity of solving the coupled boundary-value problems. It is natural to wonder whether multiple configurations are possible, and recently O'Reilly & Tresierras (2011*a*) have shown that this is indeed the case. For example, consider the three equilibrium configurations shown in figure 1. Here, two child branches are joined to a parent branch at a single point. All three branches are subject to terminal loadings, deform under their self-weight, and, for the parameter values selected to construct this example, three possible solutions to the boundary-value problem coexist. Classifying the stability of these configurations is an open question which this paper seeks to address.

By way of additional background, we note that it is easy to imagine buckling a single terminally loaded branch deforming under its own weight by elongating the branch, by increasing the mass density or by reducing its flexural rigidity. Computing estimates for the critical combination of length and flexural rigidity for a straight rod is a classical problem and several equivalent approaches can be used (e.g. Love 1927; Maddocks 1984). Dating to Max Born's seminal dissertation on the stability of a heavy elastic rod (Born 1906), the most popular method for establishing stability criteria is to use Jacobi's analysis of the second variation of the total energy and examine the possible existence of conjugate points. The non-existence of conjugate points can be considered as a necessary condition for the second variation to be non-negative, and consequently a necessary condition for nonlinear stability. We note that subtleties in using Jacobi's analysis have led to a continued effort to broaden its applicability to models featuring elastic rods (e.g. Jin & Bao 2008; Manning 2009 and references therein; Levyakov & Kuznetsov 2010). In the course of attempting to find stability criteria for branched rods, we found a treatment of the second variation due to Legendre in 1786 to be particularly fruitful. Legendre's treatment features prominently in the optimal control literature (e.g. Bryson & Ho 1975) and, as discussed in Bolza (1904) and Gelfand & Fomin (1963), is closely related to Jacobi's treatment from 1837.

In our earlier work (O'Reilly & Peters 2011), we used Legendre's treatment to establish stability criteria for terminally loaded elastic struts, and, in this paper, we show how the treatment can be extended to tree-like structures composed of elastic rods. The resulting necessary condition, which we denote by L1, for the minimization of the second variation features a series of Riccati equations, one for each branch, which are coupled at the branching points. If a solution to each individual Riccati equation can be found, then L1 is satisfied. The condition includes the appealing result that a necessary condition for a branched structure to be stable is that each branch in the structure should be stable. Another set of necessary conditions, which we denote by conditions (B3), is also established. These conditions pertain to the case where the abscissa of the branching point is variable. This variability is present in situations where the branches are held together by an adhesive at a branching point. Sufficient conditions for stability of branched tree-like structures remain to be found.

An outline of this paper is as follows: in §2, several details on the functional *I* of interest are presented and results on its first variation are recalled from O'Reilly & Tresierras (2011*a*). These results are based on the methods of Ivanov & Tuzhilin (2001). The necessary conditions for the extremization of *I* feature *N*+1 Euler–Lagrange differential equations and a pair of conditions, (B1) and (B2), at the branching point. Section 3 contains a discussion of the second variation and establishes a necessary condition, which is denoted by L1, for minimization of the functional *I*. To accommodate the case where the branching point is allowed to vary, additional necessary conditions, which are labelled conditions (B3), are established. In §4, the analyses of §§2 and 3 are applied to a tree-like structure composed of elastic rods. To enable a tractable analysis, we model the rods using Euler's theory of the elastica. With the help of a series of examples in §§5–7 we illustrate how the condition L1 can be used to determine the instability of planar tree-like structures. As regards the configurations shown in figure 1, we show that the configurations labelled A and C satisfy L1, while B is unstable because L1 is violated. This paper concludes in §8 with a discussion of potential extensions of the stability criterion to other rod theories.

## 2. The first variation of a functional defined on the branched structure

Although the calculus of variations has been in existence for well over 200 years, it is only in recent years that the systematic application of this calculus to branched structures has been explored. One example is the work of Ciarlet and his co-workers on multi-structures (see Ciarlet *et al.* 1989; Le Dret 1989). The state-of-the-art on variational principles for branched structures can be found in Ivanov & Tuzhilin (2001) and the vast majority of the developments discussed by them pertain to the first variation of a functional.

Recently, O'Reilly & Tresierras (2011*a*) used some of Ivanov & Tuzhilin's developments to establish necessary conditions for a set of functions to extremize a functional defined on a branched structure. These conditions were then applied to determine the equilibria of tree-like structures composed of elastic rods. In this section, we recall some of the developments in O'Reilly & Tresierras (2011*a*). These preliminary results are needed to establish the new necessary conditions in §3 and are framed using classical (indirect) methods of the calculus of variations.

### (a) Background and notation

We consider tree-like structures formed by the graphs of piecewise smooth functions (see figure 2). It suffices to consider the case where one branching point (or vertex) exists and to restrict attention to situations featuring the following set of piecewise smooth real-valued functions:
2.1
The domains of *x* and *x*_{1},…,*x*_{N} are
2.2
For our purposes, we assume that
2.3
It is convenient to define complementary functions for *f* and *f*_{K},
2.4
Referring to figure 2, we choose the coordinates *x*_{K} such that at a branching point *x*=*β* and *x*_{K}=*β*. Compatibility of *y* and *y*_{K} is assumed
2.5
where we have introduced the abbreviations
2.6
In this paper, we consider functionals of the form
2.7
We seek a set of extrema of *I* which satisfy the fixed boundary condition at *x*=0 and *N* compatibility conditions,
2.8
We have elected to leave the *N* boundary conditions at *x*_{J}=*L*_{J} unspecified.

To compress several lengthy expressions, we define the following notation for a pair of functions *a*(*x*) and *a*_{J}(*x*_{J}) evaluated at a branching point:
2.9
We also employ three distinct abbreviations for jumps in functions,
2.10

### (b) Variations and compatibility

We are now in a position to consider changes to *I* which arise when the functions *y**(*x*) and and the branching point *β* are varied,
2.11
where
2.12
A set of 2*N* compatibility conditions is obtained by taking the first and second derivatives of equation (2.12)_{2} with respect to *ϵ* and then letting *ϵ*→0,
2.13
Conditions of the form (2.13)_{1} feature in adhesion problems where *N*=1 (e.g. Majidi & Adams (2009) or Seifert (1991)).

### (c) Necessary conditions for an extremal

Following O'Reilly & Tresierras (2011*a*), we can establish necessary conditions for an extremal by examining the first variation of *I*. We start by substituting equation (2.11) into *I*. Next we compute d*I*/d*ϵ* with the help of the Leibniz rule, take the limit as *ϵ*→0 and perform an integration by parts. At the conclusion of these manipulations, we find that
2.14
With the help of equation (2.13), we can simplify the right-hand side of equation (2.14) by substituting for *η*^{+}_{J}. Then, using equation (2.10)_{2}, we can write equation (2.14) in another form,
2.15
It is illuminating to note that the first two terms on the right-hand side of this equation are obtained by varying *y*(*β*) and *β*, respectively, of the branching point. For those branches where *y*_{J}(*x*_{J}=*L*_{J}) are prescribed, the third term on the right-hand side of equation (2.15) will vanish because *η*_{J}(*x*_{J}=*L*_{J})=0 for these branches.

In order for to be extremals, the right-hand side of equation (2.15) needs to vanish for all *μ*,*η*^{−},*η*(*x*) and *η*_{J}(*x*_{J}). The arbitrariness of *η*(*x*) and *η*_{J}(*x*_{J}) leads to the respective Euler–Lagrange equations,
2.16
For those branches where *y*_{J}(*x*_{J}=*L*_{J}) is unspecified, the vanishing of the right-hand side of equation (2.15) also yields the natural boundary conditions
2.17
The peculiar feature of the branched system arises when we consider the first two terms on the right-hand side of equation (2.15). We suppose that we can vary *η*^{−} and *μ* in an arbitrary manner. This leads to a pair of conditions, which we refer to as the branching point conditions (B1) and (B2), respectively:
B1
and
B2
These conditions are similar in form to the classical Erdmann–Weierstrass conditions, which hold at a corner point of an extremal. The condition [[*f*_{r}]]_{B}=0 can also be related to the junction conditions discussed in works on multi-structures and continuity of force and moments at a branching point (see Le Dret 1989; O'Reilly & Tresierras 2011*a*). Furthermore, the relationship between the condition [[*g*]]_{B}=0 and material (or configurational) forces is discussed in Faruk Senan *et al.* (2008) and O'Reilly & Tresierras (2011*a*).

The assumptions (2.3) on the functions *f*,*f*_{1},…,*f*_{N} imply that the *N*+1 strengthened Legendre necessary conditions for the set of extremals are identically satisfied. By fixing all but one of the functions *y*,*y*_{1},…,*y*_{N}, we can establish *N*+1 Weierstrass necessary conditions (featuring his excess function) and 2*N*+2 Erdmann–Weierstrass corner conditions. In the applications considered in the sequel, these 3*N*+3 conditions are satisfied and, in the interests of brevity, are not discussed further.

In summary, if constitutes a set of extremals for the variational problem, then (as expected) each element satisfies the appropriate Euler–Lagrange necessary condition (cf. (2.16)). If one end of an extremal is not specified, then the natural boundary condition (2.17) is assumed to hold there. The peculiar nature of branched structures manifests at the branching point. Here, we have the pair of conditions (B1) and (B2). As emphasized in O'Reilly & Tresierras (2011*a*), apart from the condition [[*g*]]_{B}=0, the developments in this section can be considered as special cases of the results established by Ivanov & Tuzhilin (2001).

## 3. A stability criterion from the second variation

To establish a stability criterion for the set , we next examine the second variation *δ*^{2}*I* of *I*. We shall assume that *δ*^{2}*I*≥0 is necessary for stability of the set of extremals , while positive definiteness is sufficient for stability of this set. A typical development of a criterion involving *δ*^{2}*I* features the Jacobi condition and the search for conjugate points. Although some work on generalizing the condition to networked structures has been presented in Pronin (2001), we find it easier in dealing with the conditions at branching points to extend a treatment originally proposed by Legendre in 1786 for a single extremal to branched structures. This Legendre-inspired treatment leads to a series of *N*+1 Riccati equations: one for each branch. The equations are coupled in a simple manner by a boundary condition at the branching point. The existence of a bounded solution to these equations is one of the sought-after necessary conditions for *δ*^{2}*I*≥0. We were unable to establish a condition for *δ*^{2}*I*>0 for a branched structure.

We consider the same variations as those used in §2*b*. Starting with the expression (2.7) for *I* and substituting the variations (2.11), we compute the second derivative of *I* with respect to *ϵ* and then set *ϵ* to zero,
3.1
In writing equation (3.1), we have used the standard abbreviations
3.2
and
3.3
The functions *R*,*R*_{1},…,*R*_{N} are strictly positive (cf. (2.3)).

With the help of the 2*N* compatibility conditions (2.13), we can simplify the right-hand side of equation (3.1) by substituting for *η*^{+}_{J} and (*η*′_{J})^{+}. With some rearranging we find that
3.4
and
3.5
With the help of condition (B1), we substitute equations (3.4) and (3.5) back into equation (3.1) to find that the second variation can be decomposed into two terms,
3.6
Here, *I*_{β} is entirely associated with varying the branching point *x*=*β* and the components of *I*_{2} have a familiar form,
3.7
where
3.8
We now turn to examining necessary conditions for the minimization of *I*.

### (a) Necessary conditions for the minimization of *I* by

To establish necessary conditions it is convenient to first consider the case where the branching point is fixed: *μ*=0 and *I*_{β}=0. It suffices to seek conditions for *I*_{2}≥0. Following Legendre, we introduce the functions *w*(*x*),*w*_{1}(*x*_{1}),…,*w*_{N}(*x*_{N}) and add the following identity to the right-hand side of equation (3.7)_{2}:
3.9
After performing the usual manipulations of *I*_{2}, it becomes evident that we can dramatically simplify *I*_{2} provided solutions to the following *N*+1 Riccati equations exist,
3.10
These solutions are subject to the conditions^{1}
3.11
which were obtained by examining the underbraced terms on the right-hand side of equation (3.9). We also note that, for many applications, equation (3.11)_{3} reduces to the condition [[*w*]]_{B}=0. This condition couples the solutions of the Riccati equations at the branching point and is surprisingly easy to accommodate.

If bounded solutions *w*,*w*_{1},…,*w*_{N} to equations (3.10) and (3.11) can be found, then the resulting expression for is non-negative:
3.12
This expression is an expected generalization of a part of the second variation for problems involving multiple functionals. Conversely, by examining each of the branches separately, one can conclude (using standard arguments from the calculus of variations) that the existence of the bounded solutions *w*,*w*_{1},…,*w*_{N} is also necessary for *I*_{2}≥0. We are now in a position to state the following.

### Condition L1

*If the set of extremals* *minimize I then the solutions* {*w*(*x*),*w*_{J}(*x*_{J})} *to equation (3.10) which satisfy the boundary conditions (3.11) must remain bounded*.

We now examine the case where the branching point is variable. We again proceed by adding equation (3.9) to equation (3.6). After some rearranging and repeated use of the compatibility equation (2.13)_{1}, , we find that
3.13
Noting that the variations *η*^{−} and *μ* are independent, we can conclude from equation (3.13) that a necessary condition for d^{2}*I*/d*ϵ*^{2}|_{ϵ=0}≥0 is the satisfaction of the condition L1 and, in addition, the satisfaction of the following set of conditions, which we label B3:
B3
Clearly, it is necessary to first solve the set of Riccati equations (3.10) in order to compute the boundary values featured in conditions (B3). We emphasize that the conditions (B3) only apply to the case where the branching point is variable,^{2} while the condition L1 applies in all cases. Finally, it is worth mentioning that while L1 is intimately related to the buckling stability of the individual branches, conditions (B3) pertain to the stability of the adhesion mechanism at the branching point.

### (b) A set of Jacobi transformations

The necessary condition L1 can also be expressed using an equivalent set of Jacobi equations. To see this, we recall that each of the Riccati equations in (3.10) can be transformed to a linear second-order ordinary differential equation using a Jacobi transformation. Thus, the Jacobi transformations
3.14
transform equation (3.10) to a set of *N*+1 Jacobi differential equations
3.15
As stated in Reid (1972, theorem 2.1) bounded solutions to the Riccati equation for *w*_{K}(*x*_{K}) on a given interval exist if, and only if, a solution *u*_{K}(*x*_{K}) for the corresponding Jacobi differential equation exists on the same interval where *u*_{K}(*x*_{K})≠0 and *w*_{K}=−(*Q*_{K}*u*_{K}+*R*_{K}*u*′_{K})/*u*_{K}. This theorem can be used to establish the well-known equivalence of conjugate points for the solutions to equation (3.15) to blow-up of the solutions to equation (3.10).

With the help of equation (3.14), we conclude that equation (3.15) could be investigated in place of equation (3.10). The coupling condition (3.11)_{3} now features *u*,*u*_{1},…,*u*_{N} and their derivatives,
3.16
Without loss in generality, one can choose *u*^{−}=1 and equation (3.16) can then be used to prescribe (*u*′)^{−}.

## 4. Application to a planar-branched structure

We now consider the application of the stability criteria to a branched structure composed of elastic rods. For concreteness, we consider the situation shown in figure 3, in which a tree-like structure consists of a base rod which supports two rods and . Each of the rods is modelled using Euler's theory of the elastica (see Love (1927) for a discussion of this theory). In addition, to accommodate rod-like models for plants, the rods can be tapered and possess intrinsic curvature. The latter is used to model the ‘branch-shape memory’ exhibited by plant stems.^{3} Fortuitously, for the tree-like structure of interest, the elastic rod model, boundary conditions and loading allow us to establish an expression for the total energy *V* in terms of three angles, *θ*, *θ*_{1} and *θ*_{2}, and their derivatives (see equation (4.7)).

Referring to figure 4, the position vector **r** of a material point labelled *ξ*∈[0,*L*], where *ξ* is the arclength parameter, on the inextensible centreline of has the representation **r**=*X***E**_{1}+*Y* **E**_{2}. Similarly, material points on the inextensible centreline of where *K*=1,2 are identified using a coordinate *ξ*_{K}∈[*L*,*L*_{K}] and their position vectors have the representations **r**_{K}=*X*_{K}**E**_{1}+*Y* _{K}**E**_{2}. The unit tangent vectors to the centrelines of and have the respective representations
4.1
The curvature *κ* of the centreline of is defined by the relation *κ*=∂*θ*/∂*ξ*. When this centreline is unloaded it relaxes into a centreline with a curvature *κ*_{g}. The curvature *κ*_{g} is known as the intrinsic curvature. Related remarks pertain to the two other rods.

In addition to a gravitational force per unit mass, , we assume that the tips of the rods and are subject to the respective terminal dead loadings −**F**_{1} and −**F**_{2}. These forces have the representations
4.2
We assume that the rod has a flexural stiffness *EI*_{K}=*EI*_{K}(*ξ*_{K}), a length *L*_{K}−*L*, an intrinsic curvature *κ*_{gK}=*κ*_{gK}(*ξ*_{K}), and a mass per unit length *ρ*_{0K}=*ρ*_{0K}(*ξ*_{K}), where *K*=1,2. The rods and are connected at *ξ*_{K}=*L* to the rod , which has a length *L*, flexural rigidity *EI* and mass density per unit length *ρ*_{0}=*ρ*_{0}(*ξ*). The bending moment in each elastica is prescribed by a classic constitutive equation
4.3
More general constitutive equations are possible, but equations (4.3) suffice to illustrate the stability criterion.

The total energy of the rods consists of the sum of the strain energies, gravitational potential energies and the potential energy of the terminal loads,
4.4
This representation of energy is not particularly useful. Thus, with the assistance of equation (4.1), we proceed to express the components of **r**, **r**_{1} and **r**_{2} in an integral form. For example,
4.5
It is also convenient to define the masses
4.6
Note that *m*(0),*m*_{1}(*L*) and *m*_{2}(*L*) are the total masses of the rods and , respectively.

With the help of equations (4.5) and (4.6) and some additional manipulations, we arrive at an equivalent expression for the total energy,
4.7
where **W** is the combined weight of the rods and ,
4.8
We have expressed equation (4.7) in a form which allows it to be identified with the functional *I*: *y*=*θ*, *y*_{1}=*θ*_{1}, *y*_{2}=*θ*_{2}, *x*=*ξ*, *x*_{K}=*ξ*_{K}, etc. The resulting expression also makes it transparent that the rod must support the weight and terminal loading of the branches and . It thus facilitates the extension of our work to situations featuring an arbitrary number of heavy terminally loaded rods.

### (a) First variation

The Euler–Lagrange equations (2.16) associated with the functional *V* are
4.9
The solutions *θ** and (i.e. the extremals) to these equations satisfy the boundary conditions
4.10
Finally, as the branching angles are prescribed, we have the auxiliary conditions
4.11
Because the branching angles are prescribed, the branch conditions (B1) and (B2) are not applicable to this problem.

### (b) Second variation

Prior to examining the second variation of *V* , we recall the definitions (3.2) and (3.3). For the branched rod structure of interest, it is straightforward to show that
4.12
and
4.13
The Riccati equations we need to solve can be inferred from equation (3.10),
4.14
The desired solutions to these equations need to satisfy the boundary conditions (from equation (3.11))^{4}
4.15
To see whether equation (4.14) has solutions, we integrate equation (4.14)_{2,3} backwards from *x*_{K}=*L*_{K} using the conditions (4.15)_{1,2}. At the branching point *ξ*_{1}=*L*^{+} and *ξ*_{2}=*L*^{+}, we use equation (4.15)_{3} to deduce the boundary condition for the third Riccati equation (4.14)_{1}. If solutions to all three Riccati equations can be found in this manner, then we can say that L1 is satisfied and, thus, so too is a necessary condition for stability.

## 5. Buckling of a heavy strut under a terminal load

For illustrative purposes we examine the classic problem of a heavy strut which is terminally loaded. The uniform strut has a mass density per unit length of *ρ*_{0}, a flexural rigidity *EI*, and a length *L*. As can be seen in figure 5*a*, the rod is clamped at *ξ*=0 and is subject to a vertical force −*F***E**_{1} at the free end *ξ*=*L*.

This problem has a celebrated history (see Love 1927) and it is illuminating to examine how the stability analysis we are using applies to this case. The analysis is similar to that for a terminally loaded strut in the absence of gravity that we presented in O'Reilly & Peters (2011), and we use the necessary condition L1 to classify unstable configurations. In addition, for this problem, 5.1 Consequently, we can adapt the arguments presented in Gelfand & Fomin (1963, §26) to establish a sufficient condition, which we label LS1, for stability of a single branch.

### Condition LS1

*A sufficient condition for the second variation to be positive definite in the case where there is a single branch and (5.1) hold is the existence of a solution w*(*x*) *for all x*∈[0,*L*] *to the Riccati equation for w*(*x*): *w*′+*P*−*w*^{2}/*R*=0 *with w*(*L*)=0.

For the heavy strut, the Euler–Lagrange and Riccati equations for this problem can be deduced from equations (4.9)_{1} and (4.14)_{1},
5.2
The solutions to these equations are subject to the boundary conditions
5.3
In the sequel, we shall fix the dimensionless weight parameter and vary the terminal load parameter *γ*=*FL*^{2}/*EI*.

Examining the solution to the Euler–Lagrange equation (5.2)_{1}, we observe that the straight strut (i.e. *θ**=0) is a solution for all *F* and . Assuming that *α* is sufficiently small, we find that the Riccati equation for the straight strut,
5.4
has a bounded solution provided *γ* is smaller than a critical value *γ*_{crit}.^{5} A representative sample of *w* for varying *γ* is shown in figure 5*b*. When *γ*≥*γ*_{crit}, then (5.4) does not have a solution. We conclude with the help of LS1 that the straight strut is stable for *γ*<*γ*_{crit} and, using L1, that the straight strut is unstable for *γ*≥*γ*_{crit}. For *γ*>*γ*_{crit}, the Euler–Lagrange equation (5.2)_{1} admits two non-trivial solutions (or buckled states) which are mirror images of each other. The evolution of these solutions as *F* is increased is shown in figure 5*a*. We note that, as *γ* is increased, the strut shows considerable deflection from the vertical. It suffices to examine a single Riccati equation (5.2)_{2} to determine the stability of both buckled solutions. Referring to figure 5*c*, we find that equation (5.2)_{2} possesses bounded solutions for each one of the pair of buckled states and conclude, with the help of the sufficient condition LS1, that the buckled states are stable.

## 6. Application to a tree-like structure

We next consider a system of three rods which are connected at a branching point (see figure 3). We assume that the base of the parent branch is clamped while the tips of the child branches, and , are free from external loading. To explore the wealth of solutions possible for a branched structure of this type we vary the angle *θ*_{b} at the base of and see how the resulting strain *κ*_{b} at the base changes as a function of *θ*_{b},
6.1
As emphasized in O'Reilly & Tresierras (2011*a*) the resulting function *κ*_{b}(*θ*_{b}) clearly indicates the presence of multiple solutions.^{6}

The parameter values for the example discussed in the present section are
6.2
The branching angles for the example are
6.3
The strain *κ*_{b} as a function of *θ*_{b} is shown in figure 6. Observe that, for *θ*_{b}∈(−30.37^{°},30.37^{°}), three distinct solutions to the Euler–Lagrange equations (4.9) are possible.

To characterize the stability of the solutions to the boundary-value problem, the Riccati equations (4.14) are examined. An example of such a calculation is shown in figure 7 for the case where *θ*_{b}=20^{°}. As expected, two of the solutions (A and C) satisfy the necessary condition L1, while a bounded solution to the Riccati equation (4.14)_{1} for *w* is non-existent for B. Combining the results for all values of *θ*_{b}, we conclude that the necessary condition L1 is satisfied except for the dashed portion of the graph of *κ*_{b}(*θ*_{b}) in figure 6*a*. For the portion of the graph of *κ*_{b}(*θ*_{b}) where L1 is not satisfied, we can conclude that the tree-like structures for these cases are unstable. In these unstable cases, the parent branch has buckled under the weight of the child branches.

## 7. A tree of Riccati equations

It is illuminating to apply the condition L1 to branched structures such as those shown in figure 8. This structure consists of seven branches joined at three branching points. The branches have respective lengths, arclength parameters and flexural rigidities of , and . We assume that the equilibrium configuration of the structure has been computed. Then, in order to determine whether the total energy of the structure satisfies the necessary condition for a minimum, we need to find a bounded set of solutions to a set of Riccati equations. The boundary conditions for and the jump conditions that are used to prescribe , and are easily inferred from our earlier developments (cf. equations (3.11) and (4.15)). As a result, a set of seven Riccati equations can be established,
7.1
subject to the boundary conditions
7.2
We refer to the set of Riccati equations (7.1) as a *tree of Riccati equations*. Alternatively, using the transformations (3.15), a *tree of Jacobi equations* of the form (3.16) could also be established. With the help of L1, either tree of equations will be useful in distinguishing unstable solutions in computer-generated images of trees which are based on rod models.

## 8. Concluding comments

The work presented in this paper provides the first criteria for *δ*^{2}*I*≥0 of a functional *I* defined on a branched tree-like structure. One of these criteria, which we denote by L1, is an extension of both a classical necessary condition for a single functional which was first developed by Legendre in 1786 and our application of this condition to a single elastic rod in O'Reilly & Peters (2011). For a tree-like structure, the condition L1 is satisfied if bounded solutions to a set of Riccati equations can be found. For each branch of this structure, we have shown that the boundary conditions for the functions *w*,*w*_{1},…,*w*_{N} are intimately related to the topology of the structure. If one branch in a tree-like structure is unstable, then we find that we are unable to compute a bounded solution to the Riccati equation associated with this branch. This feature has the appealing interpretation that none of the individual branches of a stable tree-like structure can be unstable. Our analysis has also taken advantage of the fact that the tree-like structures of interest have no kinematically closed loops.

One restriction in our work is that we used the simplest nonlinear rod theory to develop applications of the condition L1. We believe that it is possible to extend L1 to other more elaborate rod theories such as the KirchhoffLove rod theory (which is also known as Kirchhoff's theory of rods) and some preliminary work on this extension can be found in Peters (2011). The importance of this extension lies in the fact that Kirchhoff's theory forms the basis for the works by Prusinkiewicz and his co-workers (Jirasek *et al.* 2000; Costes *et al.* 2008), and this group has developed efficient algorithms to compute the non-planar equilibrium configurations of branched tree-like structures. The development of a criterion, which is equivalent to L1, for the rod theories these authors use would allow one to eliminate some physically unrealistic solutions in a systematic fashion.

In the interests of brevity, we have not explored applications of the set of conditions (B3). Closely related examples featuring elastic rods adhering to rigid surfaces are discussed in a forthcoming paper. The stability criteria for these examples are based on the arguments used to establish the conditions (B3).

## Acknowledgements

This research was partially supported by grant no. CMMI-0726675 from the U. S. National Science Foundation. The authors are also grateful to Prof. Ivanov for his assistance with the literature on variational methods for branched structures, Timothy Tresierras for his advice and assistance in computing the results shown in §6, and the two anonymous referees for their constructive criticisms and suggestions on an earlier version of this paper.

## Footnotes

↵1 In writing equation (3.11) we have employed the simplification that at a branching point when

*μ*=0 (cf. equation (2.13)_{1}).↵2 An example of such a situation can be found in O'Reilly & Tresierras (2011

*a*).↵3 For further details on the morphogenesis of branches, see Faruk Senan

*et al.*(2008), Goriely & Goldstein (2006), Jirasek*et al.*(2000), O'Reilly & Tresierras (2011*b*) and references therein.↵4 The condition

*η*^{2}(0)*w*(0)=0 is identically satisfied because*η*(0)=0 owing to the fixed boundary condition at*ξ*=0.↵5 An analytic expression, featuring Airy functions, for

*w*(*ξ*) can be established, but we do not pursue this here.↵6 In O'Reilly & Tresierras (2011

*a*), a branched structure of three rods (similar to that shown in figure 3) with nine possible configurations is presented.

- Received May 7, 2011.
- Accepted August 8, 2011.

- This journal is © 2011 The Royal Society