OOne-dimensional reduction of viscous jets. I. Theory
Cyril Pitrou
1, 2, ∗ Institut d’Astrophysique de Paris, CNRS-UMR 7095, Universit´e Pierre & Marie Curie - Paris VI,Sorbonne Universit´es, 98 bis Bd Arago, 75014 Paris, France Saint-Gobain Recherche, 39 Quai Lucien Lefranc, 93300 Aubervilliers, France (Dated: 31 March 2017)We build a general formalism to describe thin viscous jets as one-dimensional objects with aninternal structure. We present in full generality the steps needed to describe the viscous jets aroundtheir central line, and we argue that the Taylor expansion of all fields around that line is convenientlyexpressed in terms of symmetric trace-free tensors living in the two dimensions of the fiber sections.We recover the standard results of axisymmetric jets and we report the first and second corrections tothe lowest order description, also allowing for a rotational component around the axis of symmetry.When applied to generally curved fibers, the lowest order description corresponds to a viscous stringmodel whose sections are circular. However, when including the first corrections we find that curvedjets generically develop elliptic sections. Several subtle effects imply that the first corrections cannotbe described by a rod model, since it amounts to selectively discard some corrections. However, in afast rotating frame we find that the dominant effects induced by inertial and Coriolis forces shouldbe correctly described by rod models. For completeness, we also recover the constitutive relationsfor forces and torques in rod models and exhibit a missing term in the lowest order expression ofviscous torque. Given that our method is based on tensors, the complexity of all computations hasbeen beaten down by using an appropriate tensor algebra package such as xAct , allowing us to obtaina one-dimensional description of curved viscous jets with all the first order corrections consistentlyincluded. Finally, we find a description for straight fibers with elliptic sections as a special case ofthese results, and recover that ellipticity is dynamically damped by surface tension. An applicationto toroidal viscous fibers is presented in the companion paper [Pitrou, Phys. Rev. E 97, 043116(2018)].
Contents
I. Introduction II. Geometry
III. Fields expansion on sections
IV. Kinematics ∗ Electronic address: [email protected]
G. Boundary kinematics 11H. Slenderness perturbative expansion 11I. Gauge fixing 12J. Shape restriction 12K. Velocity shear rate 12
V. Dynamics
VI. Application to axisymmetric jets
VII. Application to curved fibers
VIII. Conclusion a r X i v : . [ phy s i c s . f l u - dyn ] A p r Acknowledgments References A. STF formalism
B. Alternate shape representation C. Velocity of the coincident point D. Cartan structure relation E. Second set of corrections foraxisymmetric viscous fibers F. Higher order constraints for curved fibers G. Higher order corrections for curved fibers I. INTRODUCTION
Solving exactly the non-linear fluid equations for longviscous jets is extremely complicated and one needs toresort to an approximation scheme to study the dynam-ics of these systems. Due to the elongated shape, thereis an obvious simplification which consists in consideringa one-dimensional description. A body is considered asbeing slender if its radius R is typically much smallerthan the inverse size of velocity gradients L , that is if thevelocity field changes on length scales which are largerthan the fiber radius. Hence the one-dimensional reduc-tion induces naturally an expansion in the slendernessparameter (cid:15) R ≡ R/L . Given that at lowest order, asolid is approximated by a point particle, then we expectthat a slender jet is approximated at lowest order bysome type of string. Furthermore, as extended objectsare described as point particles with an internal struc-ture encoded in various moments (e.g. in the momentof inertia), the internal structure of the one-dimensionalobject which approximates a viscous jet is encoded insome moments which vary continuously along the one-dimensional fiber. From the perturbative expansion inthe small parameter (cid:15) R we show how this series of mo-ments must be truncated at a given order of correctionsaround the string description. The various moments de-scribing the viscous jet happen to separate naturally intomoments which evolve dynamically and moments whichare related by constraints to the former ones.For simplicity, we restrict our analysis to incompress-ible Newtonian fluids whose internal forces are capturedentirely by a constant viscosity parameter, and we allowfor surface tension effects. These ingredients are suffi-cient to describe the dynamics of drop formation fromthe Rayleigh-Plateau instability (Plateau 1873, Rayleigh1878, Eggers 1997). However, concerning the global shape of the viscous jet, our aim is to remain as general aspossible, allowing for curved fibers (that is curved centrallines) with possibly non-circular cross sections. Indeedthere are a series of geometrical simplifications which areusually performed given the symmetries of specific prob-lems. From the most restrictive to the most general, wefind the axisymmetric case, the straight fiber case withnon-circular sections, the curved fiber case with circu-lar sections, and the curved fiber case with non-circularsections. • Axisymmetric fibers: the fiber central line is astraight line and the cross sections around that cen-tral line are disks only. The one-dimensional re-duction of viscous jets for this geometry has beenextensively studied in previous literature with sev-eral non-equivalent methods. A first method con-sists in using the Cosserat theory (Bogy 1979),and it has been shown that this method is in factequivalent to expanding the velocity fields alonga suitable basis of functions (Eggers 1997, Eg-gers and Villermaux 2008). The second methodis based on a radial expansion (mathematically aTaylor expansion) of velocity fields and it has beendeveloped in, e.g., Garc´ıa and Castellanos (1994),Eggers and Dupont (1994), or Bechtel et al. (1995)when allowing for a possible angular rotationaround the axis of symmetry. The validity of thesemethods has been studied in details in the subse-quent literature, e.g. in Perales and Vega (2010),Ga˜n´an-Calvo et al. (2011), Montanero et al. (2011)or Vincent et al. (2014). In § VI, we recover thestandard lowest order results plus first correctionsusing the radial expansion method. We also reporta general method to obtain recursively its correc-tions up to any order and report the second set ofcorrections. In this geometry, once the constraintsfrom the stress tensor on the fiber side have beenused, the fundamental dynamical variables appearto be the velocity along the axis v , the local rota-tion rate around that axis ˙ φ , and the radius R . • Straight fibers: the fiber central line is still astraight line, but the cross sections can have moregeneral shapes. We find that the section shape ismost conveniently expanded into shape multipoleswhich are symmetric trace free tensors. The lowestmultipole describes for instance the elliptical modu-lation of the cross sections (Bechtel et al. 1988a,b).Under this description, the shape multipoles areadditional fundamental variables. • Curved fibers with circular sections: the fiber cen-tral line can have any general shape as long as thecurvature radius remains larger than the typicalextension of the cross sections. A formalism wasinitially developed in Entov and Yarin (1984) andfurther summarized in Yarin (1993, 2011). Curvedfibers were considered with surface tension effectsin Dewynne et al. (1992), Cummings and Howell(1999) and also in Arne et al. (2009, 2015) to studyrotational spinning processes such as those used inthe production of glass wool or candy floss. A sim-ilar viscous rod model, based on curvilinear coordi-nates adapted to the problem, has been developedby Ribe (2004), Ribe et al. (2006, 2012) to studythe coiling of viscous jets, and numerical methodswere developed by e.g. Audoly et al. (2013), Bergouet al. (2010) to obtain general solutions. In thisarticle, we develop a formalism based on a 2 + 1splitting (Miyamoto 2010) of equations, that is, aseparation between the two-dimensional fiber sec-tions and its one-dimensional central line, to reducecurved viscous jets as one-dimensional object. Wefirst describe in full generality the central line alongwhich the jet is described by following essentiallythe method developed by Ribe. The tangential di-rection of this central line naturally determines afiber direction and a fiber section which is orthogo-nal to it, along which our 2+1 splitting of equationsis performed. Then, using the irreducible repre-sentation of SO(2), we build an expansion of thevelocity field. It is based on symmetric trace-freetensors which are lying in the fiber sections andwe show that these tensors are the moments whichnaturally take into account the internal structureof the fiber. Eventually, the fundamental variablesare the same as for straight fibers with non-circularsections (velocity along the axis v , local rotationrate around that axis ˙ φ , fiber radius R , and shapemultipoles), since all other velocity moments can beobtained as constraints from these variables. Thesefundamental variables must also be supplementedby the fiber central line position and velocity. Thedifference with the straight case lies mainly in thefact that circular sections are only compatible withthe lowest order description, that is with the vis-cous string model. Indeed, as soon as correctionsare included, the shape multipoles are necessarilysourced. For instance, terms which are quadraticin the central line curvature generically source theellipticity.Our consistent description of elongated but possiblycurved viscous fibers allows to find a number of qual-itative results on the structure of the one-dimensionalmodels which contrast with past literature. Most im-portantly, we find that the first corrections for curvedfiber geometries cannot be encompassed by the rod modelof Arne et al. (2009), Ribe (2004), Ribe et al. (2006).These methods are based on the observation that, whenconsidering extended solid objects instead of point parti-cles, we must supplement the momentum balance equa-tion by an angular momentum equation. It is thus ex-pected that to go beyond the string approximation whichcan also be obtained from a momentum balance equation,we should use some form of angular momentum balanceequation. However, this method inspired from solids fails for viscous fluids for a number of reasons which are ab-sent in non-deformable solids.At lowest order in (cid:15) R , the rotation of the fiber sectionfollows the rotation of the fluid on the fiber central line.Since the dynamics of the central line is determined fromthe momentum balance equation, its local rotation rateis also derived from it, implying that the angular momen-tum method cannot bring any new dynamical informa-tion about fiber section rotation. In fact, it is preciselybecause the rotation of sections is determined from therotation of the central line at lowest order that the an-gular momentum method is instead a constraint on thesectional component of the viscous forces which appearas corrections to the lowest order description. Eventu-ally, the coupling between the momentum balance andangular momentum balance equation amounts to select-ing only some corrections and discarding the other onesas several order (cid:15) R effects are missing. For instance, thesectional component of the velocity on the central line differs from the sectional component of the velocity ofthe central line by corrections of order (cid:15) R . Additionally,the longitudinal velocity develops a Hagen-Poiseuille pro-file (that is a parabolic profile in terms of the radial dis-tance), which blurs the notion of solid displacement offiber sections.For these various reasons we find that we should notbuild a one-dimensional reduction of viscous fibers fromthe usual methods which have been developed to describethe continuous deformation of solids, but we should in-stead start from a Taylor expansion of the velocity fieldand find a consistent truncation at any given order.When deriving corrections to a viscous string model, thisrequires to abandon the hypothesis of circular sectionsand to derive the dynamical equations for the shape mo-ments as well.
1. Outline In § II we review the general formalism to describe thefiber central line and fiber sections. We introduce a coor-dinates system and a vector basis which are adapted tothe description of curved fibers. In § III, we review in de-tails how scalar fields and vector field such as velocity canbe expanded in multipoles using an adapted 2+1 decom-position. This leads us to introduce irreducible represen-tations of SO(2) according to which these multipoles areclassified. With this formalism clearly established, it isthen possible to explore in full generality all the kinemati-cal relations of fluid fibers in § IV and then the dynamicallaws in § V. The formalism is then applied for the twomain geometries of interest. First, in § VI we rederivethe axisymmetric results up to first corrections, includ-ing the interplay between the velocity along the axis ofsymmetry, and rotation around that same axis. Secondcorrections of the axisymmetric case are also reportedin Appendix E. In § VII we then tackle the problem ofcurved fibers, first deriving the lowest order string modeland then discussing the general method to obtain correc-tions. We report the first corrections in the curved caseand show that elliptical shapes are necessarily sourced atthat order. We also discuss the special case of straightbut non-circular sections in § VII D. Finally, we com-pare our method with the rod models where dynamicalequations are usually obtained from momentum/angularmomentum balance and we discuss the range of valid-ity of these methods. Several technical developments aregathered in the Appendices, among which the symmetrictrace-free tensors in Appendix A which we use through-out the article. The formalism is applied to study toroidalviscous fibers in Pitrou (2018).
2. Notation
Assuming incompressibility, the mass density ρ is con-stant. Hence, in the remainder of this article we will usethe notation µ/ρ → µ , ν/ρ → ν , P/ρ → P, (1.1)where µ is the viscosity of the fluid, ν the surface tensionparameter, and P the pressure.We use Einstein summation convention whenever indicesare placed in pairs with one index up and one index downas, e.g., in x µ x µ or X i Y i . II. GEOMETRY
For axisymmetric fibers, it is natural to use cylindricalcoordinates. The third coordinate ( z ) is naturally asso-ciated with the axis of symmetry, and the other coordi-nates ( r, θ ) parameterze the two-dimensional space whichis orthogonal to this axis. However, for generally curvedfibers, we must use fiber adapted coordinates. They areclosely related to cylindrical coordinates in the sense thatwe choose a central line inside the fiber and we use thenatural coordinate of this line as the third coordinate.The two other coordinates are then used to describe theplanes which are orthogonal to this central line. In thissection, we construct fiber adapted coordinates and theorthonormal basis which is naturally associated with it. A. Description of the fiber central line
Throughout this article, we use Cartesian coordinates x µ , with the corresponding canonical basis of vector andco-vectors (forms) e µ and e µ = d x µ . Greek indices referto components in this canonical basis. The scalar andwedge products of two vectors X = X µ e µ and Y = Y µ e µ are simply X · Y ≡ X µ Y µ , [ X × Y ] α ≡ ε αµν X µ Y ν , (2.1)where ε αµν is the totally antisymmetric tensor with ε = 1. We assume that it is possible to define a fiber centralline (FCL) as an approximation to the fiber shape. Wepostpone to § IV I the discussion about the geometricalconstruction of this line. The position of this FCL andits tangent vector T = T µ e µ are given by R µ ( s, t ) , T µ ≡ ∂ s R µ , T µ T µ = 1 . (2.2)The last condition ensures that the coordinate s can alsobe used to measure lengths along the FCL, and t is theabsolute time. A parallel projector, also named longitu-dinal projector, and an orthogonal projector, also namedsectional projector, can be defined as P (cid:107) µν = T µ T ν , P ⊥ µν = δ µν − T µ T ν = δ µν − P (cid:107) µν , (2.3)and can be used to project any vectorial quantity alongand orthogonally to the fiber tangential direction T µ .The velocity of the FCL is given by U ≡ ∂ t R (2.4)from which we deduce that the evolution of the tangentvector obeys ∂ s ∂ t R = ∂ t ∂ s R ⇒ ∂ t T = ∂ s U . (2.5)From the normalization condition (2.2) of T µ , we deducethat T · ∂ t T = 0 ⇒ T · ∂ s U = 0 . (2.6)This last relation means that the FCL velocity can onlyrotate the FCL direction T , but without stretching it.Indeed, since it is only a geometrical line there wouldbe no physical information contained in such stretchingand we thus chose the normalization (2.2). The FCLcan be understood as a non-extensible but flexible wirethat would be inside the bulk of the viscous jet, and thecoordinate s can be thought of as the distance along thatwire from a reference point s = 0.The curvature of the FCL is defined as the rate ofchange of the tangential direction along the line for afixed time, that is κ ≡ T × ∂ s T , ∂ s T = κ × T , (2.7)and it depends on ( s, t ). Note that we have chosen de-liberately κ · T = 0 and our convention thus differs fromArne et al. (2009), Ribe (2004), Ribe et al. (2006) wherethe longitudinal component of curvature does not neces-sarily vanish.Finally, to alleviate the notation, for any vector X wewill use the notation (cid:102) X ≡ T × X . (2.8)Note that the vector (cid:101) κ points toward the exterior of theFCL curvature. B. Local orthonormal basis
On each fiber section we can consider an orthonormalbasis ( d , d , d ≡ T ) where the vectors of the basis de-pend only on ( s, t ). We use coordinates i, j, k, . . . to referto components along this basis. It is oriented such that d i × d j = ε ijk d k , (2.9)where ε ijk = ε ijk = ε ijk with ε = 1 is the alternat-ing symbol which is fully antisymmetric. This orthonor-mal basis is related to the canonical Cartesian basis by achange of basis d i = d iµ e µ . (2.10)Similarly, the co-basis which depends also only on ( s, t ) isrelated to the canonical co-basis by a change of co-basis d i = d iµ e µ . (2.11)The change of basis and co-basis satisfy the basic prop-erties d iµ = δ ij δ µν d jν = d iµ (2.12a) d µ = T µ d µ = T µ = δ µν T ν . (2.12b)Any vector on the central line, is decomposed as X ( s, t ) = X i ( s, t ) d i ( s, t ) . (2.13)The fiber central line and its associated orthonormal ba-sis are illustrated in Fig. 1. Note that for components inthe orthonormal basis, the position of indices does notmatter. Indeed from (2.12), we find that for any vector X i = X i . In the orthonormal basis, the scalar productsand wedge products of two vectors X and Y are X · Y = X i Y i , [ X × Y ] i = ε ijk X j Y k . (2.14) FIG. 1: Fiber central line position R µ , and the associatedorthonormal basis d µi . A point lying in a fiber section is thenlabeled by x µ defined in (2.29). Following the evolution of the tangent vector along thefiber (2.7), we choose to transport this basis along theFCL according to ∂ s d i = κ × d i ⇔ κ i = (cid:15) ijk ( ∂ s d j ) · d k , (2.15)which is consistent with (2.7) since T = d . C. Rotation of the orthonormal basis
The rotation rate of the orthonormal basis is definedby ∂ t d i = ω × d i ⇔ ω i = ε ijk ( ∂ t d j ) · d k . (2.16)In particular ∂ t T = ω × T and the sectional projectionof rotation satisfies P ⊥ ( ω ) = T × ∂ t T = T × ∂ s U . (2.17)From the definitions (2.7) and (2.16) of curvature androtation, we get the structure relation ∂ t κ − ∂ s ω = ω × κ . (2.18)This relation is a consequence of the flatness of the clas-sical structure of spacetime as explained in Appendix D.When projected on T , it also implies( ∂ s ω ) · T = 0 . (2.19)We use the indices a, b, · · · = 1 , κ is projectedsince P ⊥ ( κ ) = κ , but the rotation is not a projectedvector, so their decompositions in components read as κ = κ a d a , ω = ω i d i , P ⊥ ( ω ) = ω a d a . (2.20)The notation (2.8) in components is simply (cid:102) X a = ( T × X ) a = X b ε ba , (2.21)where ε ab = ε ab ≡ ε ab is the two-dimensional alternatingsymbol ( ε = 1). Since ∂ s ε ab = ∂ t ε ab = 0 the tildeoperation and the ∂ s or ∂ t derivatives commute.From (2.18) the relations between the derivatives ofcomponents of curvature and rotation read as ∂ t κ i − ∂ s ω i = ε ijk κ j ω k . (2.22)Separating the sectional and longitudinal componentsthese relations are just ∂ t κ a − ∂ s ω a = ε ab κ b ω = − (cid:101) κ a ω , (2.23a) ∂ s ω = ε ab ω a κ b = − ω a (cid:101) κ a = (cid:101) ω a κ a . (2.23b)This last relation (2.23b), which is the component versionof (2.19), states that once the curvature and the sectionalprojection of the rotation are fixed at a given time alongthe FCL, then the longitudinal component of rotation ω is also determined at that time along the FCL, as aconsequence of the structure relation (2.18). We call thistype of relation a constraint equation . D. Essential relations for components
For a vector X ( s, t ), the component of the derivatives[e.g., ( ∂ t X ) i ≡ d i · ∂ t X ] are not the derivatives of thecomponents [e.g., ∂ t X i ≡ ∂ t ( d i · X )]. Indeed, they arerelated by ( ∂ t X ) i = ∂ t X i + ( ω × X ) i = ∂ t X i + ε ijk ω j X k , (2.24a)( ∂ s X ) i = ∂ s X i + ( κ × X ) i = ∂ s X i + ε ijk κ j X k . (2.24b)In particular, for the projected (or sectional) indices, thisreads as ( ∂ t X ) a = ∂ t X a − (cid:101) ω a X + (cid:102) X a ω , (2.25a)( ∂ s X ) a = ∂ s X a − (cid:101) κ a X . (2.25b)If a vector X is projected, that is orthogonal to the tan-gential direction ( X = 0), then in that special case,( ∂ s X ) a = ∂ s X a but note that we still have ( ∂ t X ) a (cid:54) = ∂ t X a . From the derivatives of the components, that is,from ∂ t X i , ∂ s X i , we can conversely obtain the deriva-tives of the vector using (2.24) with the simple relations ∂ s X = ( ∂ s X ) i d i = ( ∂ s X i ) d i + κ × X , (2.26a) ∂ t X = ( ∂ t X ) i d i = ( ∂ t X i ) d i + ω × X . (2.26b)Let us report in coordinates some relations previouslyobtained in a covariant form. Eq. (2.6) reads as in com-ponents ∂ s U + ε ab κ a U b = ∂ s U + (cid:101) κ b U b = 0 . (2.27)As for the property (2.17) for the projection of the rota-tion, it reads as simply ω a = ∂ s ‹ U a + κ a U . (2.28)This relation which is also a constraint equation statesthat once the sectional projection of the central line ve-locity ( U a ) is known at a given time along the centralline, then the projection of rotation ( ω a ) is determined.Since the longitudinal component of rotation ω is de-termined from (2.23b), we can then determine the timeevolution of the components of curvature from (2.23a). E. Fiber adapted coordinates
We use Cartesian coordinates y a = y , y inside thefiber section labeled by ( s, t ) so as to parametrize pointswhich do not lie exactly on the FCL. With y ≡ s , thefiber adapted (FA) coordinates are the set of ( y i ) =( y a , s ). The FA coordinates of a point in the fiber arerelated to the Cartesian coordinates x µ by x µ ( y i ) = R µ ( s, t ) + y a d aµ ( s, t ) . (2.29) These coordinates are illustrated in Fig. 1. The canonicalbasis and co-basis associated with the FA coordinates are e i = e iµ e µ , e i = e iµ e µ , e iµ = ∂x µ ∂y i , e iµ = ∂y i ∂x µ (2.30)and they are related to the orthonormal basis by e = h T = h d , e = 1 h d , (2.31a) e a = d a , e a = d a , (2.31b)where we defined h ≡ (cid:101) κ a y a . (2.32)Given the relations (2.31b), we will ignore in the rest ofthis article the difference between indices referring to e a (resp. e a ) and d a (resp. d a ) since they are equal andrefer to unit vectors in both cases. The distinction isonly meaningful for the third components, that is, thecomponents along the tangential direction, since d isnormalized whereas e is not.The relation between the FA coordinates canonical ba-sis and the orthonormal basis can also be written in theform d i = d ij e j , d i = d ij e j (2.33)with d = h − , d a = 0 , d a = 0 , d ab = δ ba (2.34a) d = h, d a = 0 , d a = 0 , d ab = δ ab . (2.34b)The metric and its inverse in the canonical basis of theFA coordinates are simply g ij ≡ e i · e j = Ü h ê , (2.35a) g ij ≡ e i · e j = Ü h − ê , (2.35b)We notice that it looks very much like a metric of cylin-drical coordinates. The main difference lies in the factthat h depends on all coordinates ( y , y , s ) whereas forcylindrical coordinates h is replaced by the radial coor-dinate [ h → r in (2.35)].Finally, we define the unit directional vector n a , andthe unit orthoradial vector (cid:101) n a in the section plane by y a ≡ rn a , (cid:101) n a ≡ − ε ab n b , r ≡ y a y a . (2.36)Since n a and (cid:101) n a are unit vectors which are mutually or-thogonal and projected, then we find the identities P ⊥ ab = δ ab = n a n b + (cid:101) n a (cid:101) n b (2.37) (cid:15) ab = n a (cid:101) n b − (cid:101) n a n b . (2.38)Using ( r, n a ) rather than ( y , y ) in the section planeamounts simply to using polar coordinates instead ofCartesian coordinates. With the variables ( r, n a , s ), theFA coordinates can be understood as a cylindrical coor-dinates system that one would have deformed so that the z axis is curved and onto the FCL. The factor h definedin (2.32) which is larger that unity in the exterior of cur-vature and smaller than unity in the interior ( (cid:101) κ a pointstoward the exterior of curvature) is the local amount ofstretching which had to be applied to perform this defor-mation. F. decomposition As mentioned in § II A, any vector can be decomposedinto its part along T (its longitudinal part) and a sec-tional part according to X µ = XT µ + X ⊥ µ X ⊥ µ ≡ P ⊥ µν X ν , X ≡ X · T (2.39)or in components X = X , X ⊥ a = X a . (2.40)Given this last property, we will omit the ⊥ subscriptwhen dealing with the components of the projection of atensor. This decomposition is easily extended to tensorsas each index needs to be decomposed into a longitudinaland a sectional part.From now on, we thus use the notation U = U and ω = ω for the longitudinal components of the fiber ve-locity and of the rotation rate. G. Spatial derivatives in fiber adapted coordinates
The components of the derivative of a vector in thecanonical basis of the FA coordinates are given by e iµ e jν ( ∂ µ v ν ) = ∂ i v j − Γ kij v k , (2.41a) e iµ e jν ( ∂ µ v ν ) = ∂ i v j + Γ jik v k , (2.41b) For simplicity we have chosen to use the Cartesian coordinates x µ in the ambient space. If we were to choose general curvilinearcoordinates, e.g., spherical coordinates, then in all our equationswe should promote the partial derivative to a covariant deriva-tive and perform the replacement ∂ µ → ∇ µ , ∂ i → ∇ i and δ µν → g µν , δ µν → g µν everywhere. We would also have to usethe rule ∂ s → e µ ∇ µ . The indices µ, ν . . . could even be givenan abstract meaning and not refer to a particular system of coor-dinates. This is a standard notation for general relativity (Wald1984) in which the background space is unavoidably curved, butfor classical physics this extra layer of abstraction is not neces-sary. We have thus chosen the simplest and most transparentnotation based on an ambient set of Cartesian coordinates withits associated partial derivative, but we must bear in mind thatthe results are more general. where the Christoffel symbols are defined asΓ jik ≡ e iµ e jν ( ∂ µ e kν ) . (2.42)They are related to the components of the metric in thecanonical basis of the FA coordinates byΓ ijk = 12 g il ( ∂ j g lk + ∂ k g lj − ∂ l g jk ) . (2.43)The only non-vanishing Christoffel symbols areΓ = ∂ s hh , Γ a = (cid:101) κ a h , Γ a = − h (cid:101) κ a . (2.44)The derivative of the orthonormal basis is then de-duced from (2.41) and the relations (2.31). We find d jν (cid:0) ∂ µ d iν (cid:1) = d jν e µ ε αβν κ α ( d i ) β = e µ ε aij κ a , (2.45)where from (2.31) we must use e µ = T µ /h . Finally,the divergence of a vector has a simple expression in FAcoordinates. For a general vector X µ it reads as ∂ µ X µ = 1 h ∂ i ( hX i ) = 1 h (cid:2) ∂ a ( hX a ) + ∂ s X (cid:3) (2.46)= 1 h (cid:2)(cid:101) κ a X a + ∂ s X (cid:3) + ∂ a X a . In particular, this allows to obtain the Laplacian of ascalar function S by using X µ = ∂ µ S , and we get∆ S = ∂ a ∂ a S + ∂ s Sh − ( ∂ s S ) y a ∂ s (cid:101) κ a h + (cid:101) κ a ∂ a Sh . (2.47)
III. FIELDS EXPANSION ON SECTIONS
With the FA coordinates we have an appropriatemethod to describe the position of a fluid particle in-side the viscous jet. However in order to describe thedynamics of the viscous fluid inside the fiber, we alsoneed to find an appropriate description for the fluid ve-locity itself, and this is the goal of this section. Clearly,in the slender approximation the velocity field is neces-sarily very close to the velocity U µ on the FCL. It isthus natural to Taylor expand the velocity field aroundthe velocity on the FCL. However, we cannot keep allorders of this expansion and we need general principlesto guide us in truncating such expansion. Since the fibersections are two-dimensional, it is natural to classify thevarious orders of the Taylor expansion according to theirtransformation property under SO(2), that is the localgroup of rotations around the tangential direction. Inthis section we give the essential steps to find the irre-ducible representations of SO(2) and their tensorial ex-pressions. Further details for this method are reportedin Appendix A. A. Taylor expansion
Any 2-field, that is a field on the fiber section withonly sectional indices, can be Taylor expanded in thevariables ( y , y ) for each ( s, t ). For a scalar field, thisTaylor expansion is of the form S ( y a , s, t ) = ∞ (cid:88) (cid:96) =0 S L ( s, t ) y L (3.1)where we use the multi-index notation L = a . . . a (cid:96) onwhich the Einstein summation convention applies. The S L are necessarily symmetric rank- (cid:96) tensors since the y L are symmetric. They are obtained by successive applica-tions of ∂/∂y a on S . If there is no index, that is, for (cid:96) = 0in (3.1), then we use the notation S ∅ . For a 2-vector theexpansion is instead of the form V a ( y i , t ) = ∞ (cid:88) (cid:96) =0 V aL ( s, t ) y L , (3.2)and this is obviously extended to higher order tensors.The tensors V aL are symmetric in the (cid:96) indices L butthere is no particular symmetry involving the index a .If we decompose the total fluid velocity V µ into itslongitudinal part V and its projected part V a as in § II F,then the former is a scalar function whereas the latter isa 2-vector and they are Taylor expanded, respectively, asin (3.1) and (3.2).
B. Irreducible representations of
SO(2)
It proves useful to decompose the y a dependence in theTaylor expansion of (3.1) by separating the dependencein the radial coordinate r and the direction n a . Thedependence in the direction n a can be further decom-posed onto the irreducible representations (irreps here-after) of SO(2). This method follows a method longknown in three dimensions where the direction vector liesin the two-sphere, and the irreps considered are those ofSO(3) (Courant and Hilbert 1953, Thorne 1980, Blanchetand Damour 1986, Blanchet 1998).The irreps of SO(2) are given by the functions e i nθ .More precisely the irreps are two-dimensional as theyare represented by e ± i nθ , except for D which is one di-mensional and which is just the set of constants. Wenote these irreps D n . Any function depending on n a =(cos θ, sin θ ) is indeed expanded in Fourier series as f ( θ ) = ∞ (cid:88) n = −∞ f n e i nθ = f + ∞ (cid:88) n =1 f n e i nθ + f − n e − i nθ . (3.3)The symmetric trace-free (STF) 2-tensors are also irrepsof the rotation group SO(2) just like STF 3-tensors areirreps of the group SO(3). The corresponding expansion of f reads as simply f ( n i ) = ∞ (cid:88) (cid:96) =0 f L n L . (3.4)where this time the f L are symmetric but also tracelesstensors of rank (cid:96) (we recall the multi-index notation L = a . . . a (cid:96) ). The two expansions are related thanks to (for m > Y ma ...a m ≡ ( d + i d ) a . . . ( d + i d ) a m (3.5a) Y − ma ...a m ≡ ( d − i d ) a . . . ( d − i d ) a m . (3.5b)Indeed it is easily found that f L = (cid:88) m = ±| (cid:96) | Y mL f m , f m = 12 (cid:96) f L Y mL . (3.6)The STF tensors of rank (cid:96) have indeed only two degrees offreedom. For instance for a symmetric rank-4 tensor, theonly independent degrees of freedom are f and f ,since from the traceless condition the other componentsare related through f = − f , f = − f = f . We are led to decompose all tensors appearingin Taylor expansions in STF tensors so as to obtain adecomposition in terms of irreps. C. Irreps of scalar functions
For scalar functions in the fiber section plane, the ex-pansion (3.1) is already made in terms of symmetric ten-sors and we only need to remove the traces. For anysymmetric tensor S L , the traceless part can be extractedas S (cid:104) a ...a (cid:96) (cid:105) = [ (cid:96)/ (cid:88) n =0 a n(cid:96) δ ( a a . . . δ a n − a n S a n +1 ...a (cid:96) ) b ...b n b ...b n (3.7) a n(cid:96) ≡ ( − n (2 (cid:96) − n − n − (cid:96) − Ç (cid:96) n å . (3.8)In the expression above, the symmetric part of a ten-sor T a ...a (cid:96) which is not initially symmetric is denoted T ( a ...a (cid:96) ) , and the STF part is denoted with angle brack-ets T (cid:104) a ...a (cid:96) (cid:105) . If the symmetrization ranges on the indicesof a product of tensors, this product has to be consideredas a single tensor for which the indices are symmetrized.We also use the notation n !! = n ( n − n − . . . .From the Taylor expansion (3.1), we see that by re-moving the traces in the tensors S L , we get factors of r ≡ y a y a . The expansion in terms of STF tensors onlyis thus of the form S ( y i , t ) = ∞ (cid:88) (cid:96) =0 ∞ (cid:88) n =0 S ( n ) L ( s, t ) y L r n , (3.9)where the S ( n ) L are STF. The dependence in the directionand the radial coordinates have now been clearly sepa-rated. The directional dependence is decomposed ontoSTF tensors, and the dependence in the radial distanceis an even polynomial as it is a polynomial in r . If noambiguity can arise, we can omit the sums over (cid:96), n so asto alleviate the notation. The STF multipoles of the de-composition (3.9) can be obtained from angular integralsas explained in Appendix A.Finally, in the case (cid:96) = 0, that is for the monopoleof the directional dependence, the coefficients are noted S ( n ) ∅ . For instance, if we consider the longitudinal part ofthe velocity field V , then V (0) ∅ corresponds to a uniformflow, and V (1) ∅ corresponds to a parabolic velocity profile,also known as a Hagen-Poiseuille (HP) flow. D. Irreps of -vector fields It is shown in Appendix A 2 that the expansion in ir-reps of a 2-vector is necessarily of the form V a = “ V ( n ) aL y L r n + y a Ù V ( n ) L y L r n − ε ab y b ◦ V ( n ) r n , (3.10)where the “ V ( n ) L and Ù V ( n ) L are STF tensors and where werecall that there is an implicit sum on (cid:96) and n . At firstsight, this form is very cumbersome if we are familiarwith the irreps of 3-vectors [that is, irreps of SO(3)]. In-deed for 3-vectors, the directional dependence can alwaysbe expanded in terms of electric type and magnetic typemultipoles (Thorne 1980, Pitrou 2009). The result for2-vectors is necessarily very different because the anti-symmetric tensor has only two indices in two dimensions,hence the expressions from the three-dimensional results,which also involve the antisymmetric tensor in three di-mensions, cannot be exported directly to two dimensions.Note also that the tensors “ V ( n ) L have no monopole sincethey must have at least one index. Conversely, the func-tions ◦ V ( n ) are purely monopolar since they do not haveany tensorial index. So, we can interpret (3.10) by sayingthat we have also two sets of STF tensors. The first setis made of the Ù V ( n ) L , (cid:96) ≥
0, and the second set consists inthe ◦ V ( n ) and the “ V ( n ) L , (cid:96) ≥ ◦ V (0) corresponds to a solid rotation aroundthe fiber tangential axis (axial rotation hereafter) and Ù V (0) ∅ corresponds to a radial infall of the fluid.Finally, note that we could also decompose a vectorfield by projecting its free index so as to obtain a scalarfunction. For instance, we can select the radial and ortho-radial components by projecting along n a and (cid:101) n a . Theneach scalar function can be expanded as in (3.9). The re-lation between the two methods is (with an implied sum on (cid:96) and n ) rV r ≡ y a V a = “ V ( n ) L y L r n + Ù V ( n ) L y L r n +1) , (3.11a) rV θ ≡ (cid:101) y b V b = ε a b “ V ( n ) a L − b y a L r n + ◦ V ( n ) r n +1) . (3.11b)This points to a simple method to extract the STF com-ponents of a vector. We first extract the STF componentsof its radial and orthoradial projections as for scalar func-tions (see Appendix A) and then deduce the STF multi-poles of the decomposition (3.10) from the relations (cid:2) rV r (cid:3) ( n ) L = “ V ( n ) L + Ù V ( n − L (3.12a) (cid:2) rV θ (cid:3) ( n ) L = ε (cid:104) a b “ V ( n ) a L − (cid:105) b + δ (cid:96) ◦ V ( n − , (3.12b)which are inverted as ◦ V ( n ) = (cid:2) rV θ (cid:3) ( n +1) ∅ (3.13a) “ V ( n ) L = − ε (cid:104) a b (cid:2) rV θ (cid:3) ( n ) a L − (cid:105) b (3.13b) Ù V ( n ) L = − “ V ( n +1) L + 2 (cid:2) rV r (cid:3) ( n +1) L . (3.13c)As a final comment on STF tensors, we must stress thatthese are much more adapted to abstract tensor manipu-lation than the usual representation (3.3), and we choseto perform all the tensor manipulations of this articlewith xAct (Mart´ın-Garc´ıa 2004). IV. KINEMATICS
We are now equipped with all the necessary formalismto study in details the kinematics and dynamics of vis-cous fibers. In this section we present all the relationsneeded for the kinematics, and in the next section weshall focus on the description of dynamics, so as to ob-tain the one-dimensional reduction of viscous fibers fromits intrinsic physical laws.
A. Velocity parameterization
We separate the total fluid velocity V µ into the velocityof the FCL ( U µ ) and the small difference V µ as V µ = U µ + V µ . (4.1) V is decomposed into a longitudinal part V and a sec-tional part V a which are decomposed as in (3.9) and(3.10). In order to avoid cluttering of indices when thedecomposition (3.10) is used for the velocity field, weshall use the short notation v ( n ) L ≡ V ( n ) L v ( n ) ≡ v ( n ) ∅ v ≡ v (0) ∅ , (4.2a)˙ φ ( n ) ≡ ◦ V ( n ) ˙ φ ≡ ◦ V (0) (4.2b) u ( n ) L ≡ Ù V ( n ) L u ( n ) ≡ u ( n ) ∅ u ≡ u (0) ∅ . (4.2c)0 B. Incompressibility
The incompressibility translates into a condition on thevelocity field as it implies that its divergence vanishes.This incompressibility condition is a scalar equation C E ≡ [ ∂ µ V µ = 0] , (4.3)which is a constraint for the velocity field. Since (2.6) or(2.27) imply that the fiber central line velocity U µ is alsodivergenceless ( ∂ µ U µ = 0), we deduce that ∂ µ V µ = 0 andfrom (2.46) it reads as in terms of the FA coordinates h∂ a V a + (cid:101) κ a V a + ∂ s V = 0 . (4.4)As this constraint is scalar, it can be expanded into irrepsjust like (3.9) and each STF tensor of this expansion mustvanish identically. Using the property (A13), we get0 = (cid:2) C E (cid:3) ( n ) L = ∂ s v ( n ) L + 2( n + 1) “ V ( n +1) L + δ (cid:96) κ a ˙ φ ( n ) +2( n + 1) “ V ( n +1) (cid:104) L − κ a (cid:96) (cid:105) + ( n + 1) (cid:101) κ a “ V ( n ) aL (4.5)+ Å (cid:96) n + 1 ã Å (cid:101) κ a u ( n − aL + u ( n ) L + u ( n ) (cid:104) L − (cid:101) κ a (cid:96) (cid:105) ã . In general, the tensors u ( n ) L can always be expressed interms of the other tensors using the incompressibilityconditions (4.5). In particular, from [ C E ] (0) ∅ , we get u (0) ∅ = − “ V (0) a (cid:101) κ a − ∂ s v . (4.6)This relation has a simple physical interpretation. In-deed, the right-hand side is (minus) the stretching rate ofthe velocity field on the central line T · ∂ s V given that onthe central line the velocity is just V = “ V (0) a d a + v (0) T .It thus states that a radial infall necessarily appears whenthe fiber is stretching, so as to ensure volume conserva-tion. C. Coordinate velocity
Initially, the total velocity V µ is defined as the rate ofcoordinate change, that is as d x µ / d t for the fluid elemen-tary particles. However, if we take the components of thevelocity in the basis associated with FA coordinates, thatis, the V i , we do not get the rate of change of the FA co-ordinates d y i / d t . We thus need to infer the non-trivialrelation between V i and d y i / d t . For this we define thespeed of the coincident point as V µC ≡ ∂x µ ∂t (cid:12)(cid:12)(cid:12)(cid:12) y i V iC ≡ − ∂y i ∂t (cid:12)(cid:12)(cid:12)(cid:12) x µ = e iµ V µC , (4.7)that we gave in both the Cartesian and the FA system ofcoordinates (see Appendix C for details). The speed ofthe coincident point is the speed of a point that would have constant FA coordinates y , y , s . The total velocityis related by V µ = d x µ d t = ∂x µ ∂y i d y i d t + V µC (4.8)that is, in FA coordinates by V i = V iR + V iC , V iR ≡ d y i d t . (4.9) V iR is thus the coordinate velocity for the FA coordi-nates since it equals the rate of change of FA coordinates.From the above definition of V µC and using the parame-terization (2.29) together with the property of the rota-tion (2.16), we get the expressions V µC = U µ + [ ω × ( y a d a )] µ = U µ + ( (cid:101) ω a y a ) T µ + ω (cid:101) y a d aµ (4.10a) V µR = V µ − [ ω × ( y a d a )] µ = V µ − ( (cid:101) ω a y a ) T µ − ω (cid:101) y a d aµ . (4.10b)On the expression (4.10a) it appears that the coincidentpoint velocity is the velocity of the FCL, on which isadded the solid rotation of the system of FA coordinateswhich does not vanish when the point considered is notlying exactly on the FCL. Part of the rotation is dueto the projected part of the rotation rate ω a and corre-sponds to the rotation of section planes, and the rest ofthe rotation is due to the longitudinal part of rotation ω and corresponds to a rotation of the basis vectors d and d around the fiber tangential direction, that is, toa rotation inside the section plane itself. D. Section shape description
The section shape can be characterized by its radiusas a function of the direction in the section, that is, byits curve in polar coordinates R ( s, t, n a ) . (4.11)As any function depending on the direction inside thefiber section, it can be decomposed in STF multipoles R L as R ( s, t, n a ) = R ( s, t ) (cid:32) ∞ (cid:88) (cid:96) =1 R L ( s, t ) R (cid:96) n L (cid:33) . (4.12)We call R ( s, t ) the radius of the fiber and the lowest mul-tipole R ab can be interpreted as an elliptic elongation ofthe fiber section shape. Instead of working with the mul-tipoles R L which have a dimension 1 /L (cid:96) , we can definedimensionless STF moments as “ R L ≡ R L R (cid:96) . (4.13)There exist alternate ways for describing the shape of thefiber section, and we give an example of another methodin Appendix B and relate it to the description (4.12).1 E. Normal vector
For a function depending on the direction f ( n a ), wecan define an orthoradial derivative by DDn a f ( n c ) ≡ ⊥ ba ∂f ( y c ) ∂y b (cid:12)(cid:12)(cid:12)(cid:12) y c = n c (4.14)where the orthoradial projector is defined by ⊥ ab ≡ δ ab − n a n b . (4.15)This projector satisfies ⊥ ab n b = 0, ⊥ ab (cid:101) n b = (cid:101) n a , and ⊥ aa = 1. Note that it can also be written as ⊥ ba = (cid:101) n a (cid:101) n b = r∂ b n a . (4.16)When applied on the fiber radius, this yields D R Dn a = r ∂ R ∂y a = R ∞ (cid:88) (cid:96) =1 ( (cid:96) + 1) ⊥ ba “ R bL n L (4.17)= R ∞ (cid:88) (cid:96) =1 ( (cid:96) + 1) (cid:101) n a Ä“ R bL (cid:101) n b n L ä . (4.18)A normal co-vector (not necessarily unity) to theboundary surface is N = N µ d x µ with N µ = ∂ µ Φ | Φ=0 (4.19)where the boundary function Φ is defined asΦ( s, t, y a ) = r − R ( s, t, n a ) = √ y a y a − R ( s, t, n a ) . (4.20)The components of the normal vector are in our case N a = n a − R D R Dn a (4.21a)= n a − (cid:101) n a ( (cid:96) + 1) “ R bL (cid:101) n b n L “ R L n L (4.21b) N = − ∂ s R , N = 1 h N , (4.21c)where sums over (cid:96) are implied. This vector can be nor-malized and the unit normal vector is just (cid:99) N µ ≡ N µ » N a N a + ( N ) . (4.22) F. Scalar extrinsic curvature
In order to use Young-Laplace law for surface tension,we need the general expression for the scalar part of theextrinsic curvature of the fiber boundary, which by def-inition is the divergence of the unit normal vector. Us-ing (2.46), it is thus obtained as
K ≡ ∂ µ (cid:99) N µ = (cid:101) κ a (cid:99) N a h + ∂ a (cid:99) N a + ∂ s (cid:99) N h . (4.23) G. Boundary kinematics
Since the boundary must follow the velocity field, theconstraint Φ = 0 must propagate with the velocity, thatis, ï dΦd t ò Φ=0 = ï ∂ Φ ∂t + d y i d t ∂ Φ ∂y i ò Φ=0 = 0 . (4.24)Since ∂ s r = ∂ t r = 0, then ∂ t Φ = − ∂ t R and given thedefinition (4.19) for the non-unit normal vector and thedefinition (4.9) of the coordinate velocity, this constraintis simply rewritten as R E ≡ (cid:2) ∂ t R = V iR N i = V R N + V aR N a (cid:3) Φ=0 . (4.25)We must be careful with the fact that it is V R = d s/ d t which appears and not V R , and we must thus use V R = d V R = h − V R to relate them.Despite its apparent simplicity, this equation is actu-ally rather complicated. First we stress that it is thecoordinate velocity which appears since it is the velocitywhich gives the rate of change for FA coordinates. Butmore importantly, all quantities must be evaluated on thefiber side, meaning that every occurrence of y a must bereplaced by R n a where R itself has a directional depen-dence given by (4.12). This equation is thus in generalextremely non-linear in the STF multipoles R L . We areforced to realize that it is hopeless to solve the generalproblem of curved viscous jets without major simplifica-tions, which leads to consider a perturbative scheme onwhich we should perform a consistent truncation. H. Slenderness perturbative expansion
If the typical length for variations in the velocity field is L and is much larger than the radius R , then the approx-imation that the body considered is elongated holds andwe can hope to find a coherent one-dimensional reduc-tion. The small parameter of our perturbative expansionis thus (cid:15) R ≡ R/L . (4.26)For instance, the moment v (1) comes originally from aTaylor expansion of the velocity, as can be checked fromthe dimension [ v (1) ] = [ v (0) ] /L . A term like v (1) R isthus of order (cid:15) R compared to the lowest order velocity v (0) , meaning that the former is really a correction tothe latter. In general, higher order multipoles correspondto higher orders in (cid:15) R because they primarily come fromgradients of the velocity fields which bring inverse powersof L . As an example, v (0) ab y a y b is of order (cid:15) R comparedto v = v (0) ∅ .When including the effect of curvature, we also assumethat the scale 1 / | κ | (the curvature radius) is also of theorder of L at most. Indeed if there is curvature, the ve-locity flow must adapt on scales which are commensurate2with the curvature radius. As a consequence, terms of thetype κ a y a must also be of order (cid:15) R . In any case it wouldbe impossible to consider sections which have a sectionradius larger than the fiber curvature radius. Indeed, inthat case fiber sections would intersect, hence | κ | R ∝ (cid:15) R must be small to obtain a satisfactory one-dimensionalapproximation. I. Gauge fixing
Since there are three degrees of freedom in the positionof the FCL inside the viscous jet, we can fix two of theseby asking that there is no dipole in R . This correspondsto the intuitive requirement that the fiber central lineshould be in the middle of the fiber section. Formally,this means that we fix the gauge by setting R a = 0 , (4.27)and it leads to a constraint equation when consideringthe dipole of (4.25).The gauge restriction (4.27) is an order (cid:15) R expressionas it is automatically satisfied at lowest order. It readsas indeed “ V (0) a = R Å − “ V (1) a + 6 H v (0) a − H v (cid:101) κ a + 2 κ a ˙ φ − H (cid:101) ω a − (cid:101) κ a ∂ s v + 2 ∂ s v (0) a ä + O ( (cid:15) R )(4.28)where we have defined the stretching factor H ≡ ∂ s ln R . (4.29)The gauge restriction fixes the global velocity shift withrespect to the fiber central line velocity U a which is en-coded by the velocity moment “ V (0) a . And, since thisshift is of order (cid:15) R , we can state that the projectedcomponent of the velocity on the central line is nearlythe projected velocity of the central line itself, that is,[ V a (cid:39) U a ] y = y =0 .With (4.27), we have fixed only two of the three gaugedegrees of freedom which arise from the fact that we arefree to choose any curve as the FCL. The third degree offreedom corresponds to the possible reparameterizationof the fiber inside the same curve, that is to the replace-ment s → s + f ( s, t ). Given the choice of normaliza-tion for the tangent vector in (2.2), this freedom is onlya global but time-dependent reparameterization freedom s → s + f ( t ). For every problem considered, there is anatural way to fix unambiguously the affine parameter s ,for instance setting it to s = 0 at the boundary. J. Shape restriction
With the gauge choice of the previous section, we havemanaged to cancel the shape dipole. However, we can-not assume that in general higher order shape multipoles vanish. Indeed, if we have curvature, then terms of thetype ∝ R κ (cid:104) a κ b (cid:105) which are of order (cid:15) R would source thesection shape quadrupole R ab . Typically, we expect tofind that shape terms like R L y L or R (cid:96) R L are of order (cid:15) (cid:96)R . In § VII, we discuss this scaling and show that cir-cular sections are compatible with the string descriptionof curved fibers, which is the lowest order description inwhich order (cid:15) R effects are ignored.Of course, if we restrict to axisymmetric jets as weshall do in § VI, then κ a = 0 and circular sections areconsistent throughout even though one might still wantto consider straight jets with non-circular sections. K. Velocity shear rate
In order to describe the dynamics of the viscous fluidinside the fiber, we will need to consider the gradient ofthe velocity field. Let us define the non-symmetric tensor S µν ≡ ∂ µ V ν . (4.30)Its components are given by S ab = ∂ a V b (4.31a) S a = ∂ a V (4.31b) hS a = ∂ s V a − (cid:101) κ a V (4.31c) hS = ∂ s V + (cid:101) κ a V a . (4.31d)These components are obtained either from the compo-nents in the canonical basis of the FA coordinates, that isthe S ab , S a , S a and S which we compute from (2.41),that we then project on the orthonormal basis, or usingdirectly (2.45) in the velocity decomposition V = V i d i .This tensor can be decomposed as S µν = 12 σ µν + ω µν . (4.32)where we used that the velocity shear rate is (twice) thesymmetric part σ µν ≡ S ( µν ) = S µν + S νµ , (4.33)and the vorticity is the antisymmetric part (cid:36) µν ≡ S [ µν ] = 12 ( S µν − S νµ ) . (4.34)We also define the vorticity (Hodge) dual vector by (cid:36) α ≡ ε αµν (cid:36) µν ⇒ (cid:36) µν = ε µνα (cid:36) α . (4.35)From (4.31), we find that the components of the shearin the orthonormal basis are then given by σ ab = ∂ a V b + ∂ b V a , (4.36a) hσ a = h∂ a V + ∂ s V a − (cid:101) κ a V , (4.36b) hσ = 2( ∂ s V + (cid:101) κ a V a ) . (4.36c)3Similarly, the components of the vorticity vector are (cid:36) a = ε ab (cid:0) S b − S b (cid:1) (4.37a)= ε ab ï ∂ b V − h (cid:0) ∂ s V b − (cid:101) κ b V (cid:1) ò (cid:36) = ε ab ∂ [ a V b ] = ε ab ( ∂ a V b − ∂ b V a ) . (4.37b) V. DYNAMICS
The dynamics of viscous fluids is well known and arisesfrom the Navier-Stokes equation. However, in order toachieve a one-dimensional reduction for viscous jets, wemust find a way to get rid of the physics on the fiberboundary. Enforcing the boundary conditions on thestress tensor on the fiber boundary leads to a set of threeconstraints which can be conveniently used to reduce thenumber of free fields in our one-dimensional reduction.This section is dedicated to the general construction ofthis method and we then apply it in the subsequent sec-tions for axisymmetric and curved fibers.
A. Total stress tensor and viscous forces
The total stress tensor is decomposed as τ µν = τ ( P ) µν + τ ( µ ) µν , (5.1)where the two components arise from pressure forces andshear viscosity. For a Newtonian fluid, they are simplygiven by τ ( P ) µν = − P g µν , τ ( µ ) µν = µσ µν . (5.2)The pressure is a scalar and it is thus decomposed as in(3.9) P = ∞ (cid:88) (cid:96) =0 ∞ (cid:88) n =0 P ( n ) L y L r n , (5.3)where the P ( n ) L are STF tensors.From the total stress tensor (5.1), we can define a vis-cous force per unit area on the fiber sections F µ ≡ τ µ . (5.4)As any vector it has a longitudinal part F which corre-sponds to viscous traction or compression on fiber sec-tions, and a sectional part F a . B. Boundary conditions
The boundary conditions for the stress tensor is thevector constraint C µ ≡ î τ µν (cid:99) N ν + ν K (cid:99) N µ = 0 ó (5.5) ≡ î τ ( µ ) µν (cid:99) N ν + ν K (cid:99) N µ − P (cid:99) N µ = 0 ó . As any vector, it can be decomposed into its longitudinalcomponent C and its sectional component C a . The lattercan be further decomposed into a radial contribution andan orthoradial contributions as r C ≡ C a n a , θ C ≡ C a (cid:101) n a . (5.6)However, these are not fields on the fiber sections sincethey are defined only at the boundary. They dependon the position on the FCL s , on time t , but their sec-tional dependence is only a dependence in the directionvector n a . They also depend on the fiber radius R andon the shape multipoles R L . So, in general they can beexpanded in STF components as C ( s, t, n a ) = ∞ (cid:88) (cid:96) =0 C L ( s, t, R ) R (cid:96) n L (5.7a)= ∞ (cid:88) (cid:96) =0 (cid:88) n C ( n ) L ( s, t ) R n + (cid:96) n L , (5.7b)where in the second line we have also expanded the de-pendence of the STF multipoles in powers of R , and withsimilar expansions for the radial and orthoradial bound-ary constraints θ C and r C . To be precise, for the latter,the expansion takes the form (5.7) for r C /R . Note againthe difference with (3.9) as the powers of r are replacedby powers of R because constraints are only defined onthe boundary. In practice, it is simpler to consider theconstraint (5.5) with the non-normalized normal vector N µ instead of (cid:99) N µ since they are equivalent.We realize that the total sum involves all orders of theform (cid:15) mR . Since we are eventually interested in resultswhich are valid up to a given order in (cid:15) R , we define themoments of the constraints up to a given order by C ( ≤ n ) L = m = n (cid:88) m =0 C ( m ) L R m , (5.8)with similar definitions for the radial and orthoradial con-straints.As we shall detail in two examples in § VI and § VII,we will deduce general relations from the vector con-straint (5.5), or, more precisely, from its three scalarcomponents (longitudinal, radial, and orthoradial).
C. Volumic forces
Once the stress tensor is computed, it is straightfor-ward to get the volumic forces f µ since they are expressedas f µ = ∂ ν τ νµ + g µ = µ ∆ V µ − ∂ µ P + g µ , (5.9)where g µ are long range volumic forces such as grav-ity. Its components are related to the components of4the stress tensor thanks to f = ∂ a τ a + 1 h (cid:0) ∂ s τ + 2 (cid:101) κ a τ a (cid:1) + g (5.10) f a = ∂ b τ ba + 1 h (cid:0) ∂ s τ a + (cid:101) κ b τ ba − (cid:101) κ a τ (cid:1) + g a . D. Navier-Stokes equation
The dynamics of a viscous fluid is governed by theconservation equation and the Navier-Stokes equation.The conservation equation has already been used sincefrom incompressibility it implied the divergenceless con-dition (4.3). Having developed all the tools to expressthe volumic forces, we are now in position to write theNavier-Stokes equation. It is of the form D µ ≡ [ A µ = f µ ] (5.11)where the acceleration vector is A µ ≡ ∂ t V µ | y i + V Rν S νµ + G µ . (5.12)If the dynamics is considered in a (constantly) rotatingframe, the fictitious or geometrized forces (inertial andCoriolis forces) gathered in G µ are expressed as G µ ≡ ε µαβ Ω α V β + (Ω ν x ν )Ω µ − Ω x µ (5.13)where Ω µ is the rotation of the frame with respect to aGalilean (inertial) frame.Just as any vector, the Navier-Stokes equation (5.11)is decomposed into a longitudinal part D and a sectionalpart D a which we can further decompose in irreps asin § III D. The components of the first term of (5.12) aresimply obtained from( ∂ t V ) ≡ d µ ∂ t V µ = ∂ t V + (cid:101) ω a V a (5.14a)( ∂ t V ) a ≡ d aµ ∂ t V µ = ∂ t V a − V (cid:101) ω a + (cid:101) V a ω , (5.14b)where we recall the notation V ≡ V for the longitudinalcomponent of the velocity. E. Secondary incompressibility constraint
The incompressibility implies the divergenceless con-dition (4.3), and as such it can be considered as a pri-mary constraint. Then by ensuring that this constraintremains satisfied when time evolves, we obtain a sec-ondary constraint. Indeed taking the divergence of theNavier-Stokes equation (5.11) we obtain an incompress-ibility constraint ∂ µ A µ = ∂ µ f µ = − ∆ P , (5.15)where we have assumed that the long range forces areconstant (e.g., for gravity) or have at least no divergence ( ∂ µ g µ = 0). If the long range forces have a divergence,then this should be included in (5.15).From the expression (5.12) of the acceleration, thisconstraint is D C ≡ (cid:2) S µν S νµ + 2 ε αβµ Ω β S αµ − = − ∆ P (cid:3) . (5.16)Using the decomposition (4.32) of the velocity gradienttensor into velocity shear rate and vorticity, it can berecast nicely as14 σ µν σ νµ − µ + (cid:36) µ )(Ω µ + (cid:36) µ ) = − ∆ P. (5.17)The right-hand-side of this equation can be computedby using (2.47). We will use this secondary incompress-ibility constraint to get constraints on the various mo-ment of the pressure field. F. Dimensionless reduction
So far, all physical quantities have a physical dimen-sion. It is, however, possible to build dimensionless quan-tities. Usually for viscous fluids, this is done by notingthat the dimensions of viscosity and surface tension [re-calling that we have divided them by the mass densityin (1.1)] are [ µ ] = L T − [ ν ] = L T − . (5.18)It is thus possible to define a length scale µ /ν and atime scale µ /ν from which we can define dimensionlessquantities for all physical quantities in the problem athand. However, we want to allow for the possibility ofhaving no surface tension, and we do not wish to use ν todefine dimensionless quantities. Instead, we decide thatthere is a natural length scale L in our problem whichcorresponds to the typical length of velocity variations.The time scale is then obtained from L /µ . The mainquantities in the problem are simply adimensionalizedusing these length scales and timescales. For instance,we define dimensionless quantities with t = L µ t , s = L s , v = µL v , g = µ L g , ˙ φ = µL ˙ φ , R = LR , ν = µ L ν ,κ a = 1 L κ a , ω a = µL ω a , (5.19)where we recall the notation v = v (0) ∅ and ˙ φ = ˙ φ (0) ∅ .For higher order multipoles this construction of dimen-sionless variables is straightforwardly performed as thedimension of the multipoles is read from the expansionfrom which it is defined. For instance from (3.9) we get v ( n ) L = µ/L (cid:96) +2 n v ( n ) L . Note also that by construction (cid:15) R = R so that in practice it is by expanding in pow-ers of R that we identify for any expression the variousorders in powers of (cid:15) R .5The dimensionless ratios of fluid mechanics which arerelevant for viscous jets in a rotating frame are theReynolds number, the Froude number, the Rossby num-ber, and the Weber number. If we consider a typicalreference reduced velocity v r , then they are simply ex-pressed asRe = v r , Fr = ( v r ) g , Rb = v r Ω , We = ( v r ) ν . (5.20)In the remainder of this article, we use the dimensionlessphysical quantities rather than the dimensionless num-bers, but using the dictionary (5.20), all expressions canbe recast with these dimensionless numbers.In the next section, we present two main applicationsof our formalism and we shall assume that we have per-formed such a dimensionless reduction for all quantities.In order to avoid cluttering of notation, we will omit inthe remainder of this article to specify that the quan-tities are dimensionless. In practice, the dimensionlessreduction amounts simply to replacing µ → ν , and it is thus used as a consis-tency check.Note that several other schemes would have been possi-ble to build dimensionless quantities, just by using otherphysical quantities that might be present in the problem.For instance if we use gravity, then we can build a lengthscale and a time scale without resorting to a choice of L simply by L = µ / g / , T = µ / g − / . (5.21)After building dimensionless quantities with these scales,the gravity vector g µ is replaced by a unit vector while µ → T = Ω − , L = » µ/ Ω . (5.22)After building dimensionless quantities, the rotation vec-tor Ω µ is replaced by a unit vector while µ → VI. APPLICATION TO AXISYMMETRIC JETS
Using the formalism developed so far in the particularcase of a straight viscous jet ( κ a = 0) would be equiv-alent to kill a fly with a sledgehammer, especially if wealso require axisymmetry of the fiber around the FCL. In-deed, in that case the FA coordinates are just Cartesiancoordinates, with the third coordinate s correspondingto the axis of symmetry of the problem. This problemhas been studied already in the literature (Garc´ıa andCastellanos 1994, Eggers and Dupont 1994) and we shallrederive, recover and extend these standard results to in-clude higher order corrections, and a possible rotation of the viscous fluid around the axis of symmetry (calledtorsion by Bechtel et al. (1995)).From the assumed rotational symmetry of the prob-lem, only the monopolar moments are non-vanishing andwe need only to consider the v ( n ) = v ( n ) ∅ , ˙ φ ( n ) = ˙ φ ( n ) ∅ and u ( n ) = u ( n ) ∅ . We recall that we use the short nota-tion introduced in the definitions (4.2), and in particularwe emphasize the notation v = v (0) for the fundamentallowest order component of the longitudinal velocity. Itis also possible to further restrict the problem to non-rotational flows, and require that all the ˙ φ ( n ) vanish, asdone in e.g. Garc´ıa and Castellanos (1994), Eggers andDupont (1994) but we will not perform such simplifica-tion and allow for a rotation of the fluid around the axisof symmetry as in Bechtel et al. (1995).From the symmetry of the problem we also deduce that U a = 0 and the rotation rate of the central line satisfiesnecessarily ω a = 0. Then from (2.23b), the longitudinalcomponent of rotation ω is constant along the fiber so itis reasonable to choose ω = 0. Eventually, the only pos-sible velocity components for the central line is U = U .However it is natural to also set U = 0 which amounts totaking a non-moving FCL. As a consequence from (4.10)the total velocity V µ is also the coordinate velocity V µR .The sledgehammer comes from the fact that the FA co-ordinates ( y , y , y = s ) are just Cartesian coordinatesso they can be chosen to be the Cartesian coordinates( x , x , x ) and all the machinery of FA coordinates isnot used in this case.Finally, the long range volumic forces need also to re-spect the symmetry and g a = 0, so if we consider theeffect of gravity the fiber needs to be vertical and wewill write simply g = g = g . Let us also mention thatwe discard the possibility of considering a rotating frame(Ω a = Ω = Ω µ = 0). Indeed, even though it is possiblein principle to also consider the problem in a frame ro-tating around the axis of symmetry (that is Ω (cid:54) = 0), thiswould be of very limited interest as it can be obtainedfrom the replacement ˙ φ (0) → ˙ φ (0) + Ω (see the discussionin § VII E).
A. The lowest order viscous string model
The divergenceless condition (4.5) leads simply to theset of relations u ( n ) = − n + 1 ∂ s v ( n ) , (6.1)so that we need only to consider the v ( n ) and the ˙ φ ( n ) .At lowest order ( n = 0), (6.1) has a very simple inter-pretation. A gradient in the longitudinal velocity ∂ s v (0) leads to a radial inflow u (0) because of incompressibility.Indeed, if the flow stretches, the radius shrinks to ensurevolume conservation and thus incompressibility.The lowest contributions of the Navier-Stokes equation6are D (0) for the longitudinal part and ◦ D (0) for the rota-tional part. At lowest order they lead to ∂ t v + v∂ s v = − ∂ s P (0) + 4 v (1) + ∂ s v + g (6.2a) ∂ t ˙ φ + v∂ s ˙ φ − ˙ φ∂ s v = ∂ s ˙ φ + ˙ φ (1) . (6.2b)The evolution of the radius is easily obtained from (4.25)and it reads as, at lowest order, ∂ t ln R = −H v − ∂ s v + O ( (cid:15) R ) . (6.3)In order to find a closed system of equations from (6.2),we need to find P (0) , v (1) , and ˙ φ (1) from the boundaryconstraint equation (5.5). First, it turns out that in theaxisymmetric case, the contribution of surface tensioncan only be in P (0) , and we can separate the pressurefield as P ≡ P ν + p , P ν ≡ ν K , (6.4)where p is the contribution coming from viscous forcesand where the extrinsic scalar curvature K does not de-pend on r . It is also more convenient to combine the lon-gitudinal part and the radial part of the boundary con-straint (5.5) to obtain a constraint which gives directly p and remove the pressure dependence in the longitudinalconstraint. Indeed, the components of the normal vectortake exactly the form N a = n a N = −H R , (6.5)and we find that the three scalar constraints can be ex-pressed in the form p C ≡ (cid:2)(cid:0) − α (cid:1) p = µ (cid:0) σ ab n a n b − α σ (cid:1)(cid:3) (6.6a) C ≡ (cid:2)(cid:0) − α (cid:1) σ a n a = α ( σ − σ ab n a n b ) (cid:3) (6.6b) θ C ≡ î σ ab n a (cid:101) n b = ασ a (cid:101) n a ó , (6.6c)where α ≡ H R .Using that constraints are expanded according to (5.7),then from the lowest order of the pressure constraint( p C (0) ∅ ) we obtain the pressure monopole at lowest order p (0) = u (0) + O ( (cid:15) R ) = − ∂ s v + O ( (cid:15) R ) . (6.7)The longitudinal constraint C at lowest order (that is C (0) ∅ ) gives the HP profile encoded by v (1) v (1) = H ∂ s v + ∂ s v + O ( (cid:15) R ) . (6.8)Finally the orthoradial constraint θ C does not vanishidentically if we have rotation. Instead, if we consider θ C (1) ∅ , we get ˙ φ (1) = H ∂ s ˙ φ + O ( (cid:15) R ) . (6.9)The lowest order of the incompressibility constraint D C is not needed for the lowest order dynamics and is only required when considering higher order corrections as weshall see in § VI B.Using (6.7), (6.8), and (6.9) replaced in (6.2) we arenow able to obtain ∂ t v = g − ν∂ s K − v∂ s v +6 H ∂ s v + 3 ∂ s v + O ( (cid:15) R ) (6.10a) ∂ t ˙ φ = ˙ φ∂ s v + 4 H ∂ s ˙ φ − v∂ s ˙ φ + ∂ s ˙ φ + O ( (cid:15) R ) . (6.10b)Together with (6.3), it forms a closed set of equations.Since rotation does not couple to the longitudinal velocityin (6.10a) it is reasonable to consider that the dynamicalequation (6.10b) should be considered only when includ-ing the first corrections. In fact in order to obtain it wehad to consider θ C (1) ∅ and not θ C (0) ∅ which vanishes iden-tically, so we realize that we have been using a higherorder constraint to be able to close (6.2b).We note finally that if surface tension is ignored, theevolution of velocity and rotation does not depend on theradius. However, as soon as we consider surface tension,the dependence in R appears of course through K , andwe thus need (6.3) to complement the dynamical equa-tions for v . We intentionally did not replace in (6.10) theexpression of the extrinsic scalar curvature K since eventhough we might use our perturbative expansion in pow-ers of R , it proves often useful to keep its most generalexpression when considering the axisymmetric geometry.Instead of using K = 1 /R , which is the lowest order ex-pression obtained from (4.23), we obtain a much betterdescription if we use instead the exact expression (Eggersand Villermaux 2008) K = 1 R (cid:112) ∂ s R ) − ∂ s R [1 + ( ∂ s R ) ] / , (6.11)which is easily obtained from the general expression ofthe extrinsic scalar curvature (4.23). In (6.10a), oneshould thus use (6.11). It amounts to resumming allhigher order contributions from surface tension effects,and this is made possible thanks to the decoupling prop-erty (6.4).Equation (6.10a) together with (6.11) for surface ten-sion induced pressure, and the boundary kinematic re-lation (6.3) constitute the lowest order model for an ax-isymmetric jet. Note that the expression for ∂ t v involvessecond order derivative with respect to the affine param-eter s . The factor 3 in front of the last term of (6.10a) isthe famous Trouton factor (Trouton 1906) and we reviewits origin in § VII F 1. As for the factor 6 H in the previousterm, it is easily understood from the longitudinal com-ponent of viscous forces per unit area (5.4) F (0) ∅ = 3 µ∂ s v .Indeed, this implies that the longitudinal viscous forcesintegrated on a circular section are F (cid:39) µπR ∂ s v , (6.12)which implies that the lineic density of longitudinal forcesis µ∂ s (3 πR ∂ s v ).7In the next section we detail the general method toobtain corrections up to any order in the small parameter (cid:15) R and report the detailed expressions of the first setof corrections (corrections up to order (cid:15) R ). The secondset of corrections (that is up to order (cid:15) R ) is reported inAppendix E. B. General method for higher order corrections
The equations (6.2) are formally unchanged when con-sidering higher orders because they are exact. Howeverthey involve quantities that we have determined from theside constraints up to order (cid:15) R contributions. We thusneed to determine these quantities ( v (1) , ˙ φ (1) , p (0) ) withgreater precision, that is also taking into account contri-butions of order (cid:15) R . To this end, we need to consider p C ( ≤ , C ( ≤ , and θ C ( ≤ . After a straightforward com-putation, we can show that the longitudinal constraint C reads as ∞ (cid:88) n =0 R n î ∂ s u ( n ) + 4( n + 1) v ( n +1) − H ∂ s u ( n − − n H v ( n ) + (8 n + 6) H u ( n ) ó = 0 , (6.13)and the truncated constraint C ( ≤ is found by keepingonly n = 0 and n = 1 in this sum. However we noticethat if we want to deduce v (1) from C ( ≤ , that is, keepingcorrections of order (cid:15) R , then we need an expression for v (2) at lowest order. Similarly, if we want corrections upto order (cid:15) R , we need v (3) at lowest order and v (2) up tocorrections of order (cid:15) R and v (1) up to corrections of order (cid:15) R .It is straightforward to show that the pressure con-straint p C has the general form (cid:88) n p ( n ) R n = (cid:88) n R n ï (2 n + 1)(1 + H R )1 − H R ò u ( n ) (6.14)and if we want to deduce p (0) up to order corrections oforder (cid:15) R we need the truncation p C ( ≤ in which we needto replace an expression for p (1) at lowest order. Sim-ilarly, if we want p (0) up to corrections (cid:15) R , then from p C ( ≤ we need p (2) at lowest order and p (1) up to correc-tions (cid:15) R , and so on.Finally, the orthoradial constraint θ C reads as in fullgenerality ∞ (cid:88) n =1 R n (cid:104) n ˙ φ ( n ) − H ∂ s ˙ φ ( n − (cid:105) = 0 , (6.15)and in particular, from θ C ( ≤ we deduce that if we need˙ φ (1) up to corrections of order (cid:15) R , then we need an ex-pression for ˙ φ (2) at lowest order. The structure of these Equation Variable Dependence p C (6.13) p (0) v (0) and p (1) R , p (2) R . . . C (6.14) v (1) v (0) and v (2) R , v (3) R . . . θ C (6.15) ˙ φ (1) ˙ φ (0) and ˙ φ (2) R , ˙ φ (3) R . . . TABLE I: Structure of dependencies from the pressure, radial,and orthoradial boundary constraints. recursive dependencies in the three boundary constraintsis summarized in Table I.This problem is solved if we now also consider highermoments of the Navier-Stokes equation (5.11) togetherwith the incompressibility constraint (5.15) so as to findexpressions for the missing moments. The key is to noticethat these equations contain a Laplacian, e.g., ∆ P for theincompressibility constraint or ∆ V µ for the Navier-Stokesequation. Since for any scalar S expanded as (3.9), thecoefficients of the expansion of ∆ S are[∆ S ] ( n ) L = 4( n + 1)( n + 1 + (cid:96) ) S ( n +1) L + ∂ s S ( n ) L , (6.16)then in the axisymmetric case we can use this prop-erty for (cid:96) = 0, and from the incompressibility con-straint (5.15) we can express P ( n +1) in terms of ∂ s P ( n ) .Similarly from the longitudinal part of the Navier-Stokesequation (5.11) we can express v ( n +1) in terms of ∂ s v ( n ) since it contains ∆ V .Indeed, the general expansion of the longitudinal partof the Navier-Stokes equation in powers of r n (the D ( n ) )is ∂ t v ( n ) + (cid:88) m u ( n − m ) ∂ s v ( m ) + (cid:88) m mv ( m ) u ( n − m ) (6.17)= 4( n + 1) v ( n +1) + ∂ s v ( n ) − ∂ s p ( n ) + δ n ( g − ∂ s P ν ) . For instance, using D (1) , we can obtain v (2) at lowestorder, in function of ∂ t v (1) and also ∂ s p (1) . Let us ignorethis latter dependence for the sake of simplicity. Giventhat we already know the lowest order expression of v (1) in terms of v = v (0) from (6.8), then using it we obtain v (2) as a function of v (0) . The time derivatives on v (0) can be further replaced with the lowest order dynamicalequation of v (0) (6.10a). In the end, we have obtained v (2) in terms of v (0) and its derivatives with respect to s ,thus having the structure of a constraint equation.As for the rotational part, its expansion in powers of r n , that is ◦ D ( n ) , leads to ∂ t ˙ φ ( n ) + (cid:88) m (cid:16) v ( n − m ) ∂ s ˙ φ ( m ) + ( m + 1) ˙ φ ( m ) u ( n − m ) (cid:17) = ∂ s ˙ φ ( n ) + 4( n + 1)( n + 2) ˙ φ ( n +1) . (6.18)Similarly, from ◦ D ( n ) we see that we can obtain ˙ φ ( n +1) asa function of ∂ s ˙ φ ( n ) , and the time derivatives are eventu-ally eliminated in the same manner once all replacements8 Equation Variable Dependence D ( n ) (6.17) v ( n +1) p ( n ) , v ( m ) m ≤ n ◦ D ( n ) (6.18) ˙ φ ( n +1) ˙ φ ( m ) , v ( m ) , m ≤ n D C ( n ) (6.19) p ( n +1) p ( n ) , v ( m ) , ˙ φ ( m ) , m ≤ n TABLE II: Structure of dependencies for the constraints de-duced from the Navier-Stokes equation and the incompress-ibility constraint. with lowest order relations are performed. And finally,using the expansion in powers of r n of the incompress-ibility constraint, that is using the D C ( n ) , we get (cid:88) m ï n + 52 + 2 m ( n − m ) ò u ( m ) u ( n − m ) + (cid:88) m (cid:104) − n + 1) ˙ φ ( m ) ˙ φ ( n − m ) + 2 mv ( m ) ∂ s u ( n − m ) (cid:105) = − ∂ s p ( n ) − n + 1) p ( n +1) . (6.19)We see that from D C ( n ) we can obtain p ( n +1) in termsof ∂ s p ( n ) and this time it is directly in form of a con-straint since it does not involve any time derivatives.The structure of these recursive dependencies deducedfrom the Navier-Stokes equation and the incompressibil-ity constraint is summarized in Table II.We understand that the general procedure is very re-cursive but simple, and this motivates the use of abstractcalculus packages (such as Mathematica ) to circumventthe complexity of these tedious abstract computations.1. Initially from the lowest order of the constraints p C , C and θ C (6.13,6.14,6.15) we get p (0) , v (1) and ˙ φ (1) at lowest order in function of v (0) and ˙ φ (0) , namely(6.7), (6.8) and (6.9)2. Then, if we know the p ( q ) , v ( q ) and ˙ φ ( q ) with 0 ≤ q ≤ n up to order (cid:15) mR , then from (6.17,6.18,6.19)we can deduce p ( n +1) , v ( n +1) and ˙ φ ( n +1) up to order (cid:15) mR as summarized in Table II.3. By using the constraints p C , C and θ C (6.13,6.14,6.15) we see that we can find p (0) , v (1) and ˙ φ (1) up to order (cid:15) mR if we know the p ( q ) , v ( q +1) and ˙ φ ( q +1) up to order (cid:15) m − q ) R , and this issummarized in Table I.4. From the second point, we see that if we know the p (0) , v (1) and ˙ φ (1) up to order (cid:15) mR then we alsoknow all p ( q ) , v ( q +1) and ˙ φ ( q +1) to the same order (cid:15) mR and from the third point this is what we needto know the p (0) , v (1) and ˙ φ (1) up to order (cid:15) m +1) R ,which validates this recursive method. 5. Eventually, time derivatives on the fundamentalvariables v and ˙ φ which appear in corrective terms,can be replaced by using their dynamical equationsat a lower order and it is thus possible to obtainthe corrections to the lowest order model up to anyorder only in terms of derivatives with respect to s .Finally, the evolution of the radius is easily found upto any given order in (cid:15) nR . Indeed given that ∂ t ln R = (cid:88) m Ä u ( m ) − H v ( m ) ä R m , (6.20)then in order to obtain the evolution up to (cid:15) nR correc-tions, we need the u ( m ) and v ( m ) (with m ≤ n ) up to (cid:15) n − m ) R corrections. C. First corrections
Implementing the procedure described in the previoussection, we first get the constraints p (1) = ˙ φ − ( ∂ s v ) − ν ∂ s K + ∂ s v + O ( (cid:15) R ) (6.21a) v (2) = H ∂ s H ∂ s v − H ( ∂ s v ) + ˙ φ∂ s ˙ φ − ν H ∂ s K + H ∂ s v − ∂ s v∂ s v − ν ∂ s K + H ∂ s v + ∂ s v + O ( (cid:15) R )(6.21b)˙ φ (2) = − ˙ φ∂ s H ∂ s v + H ∂ s H ∂ s ˙ φ − ∂ s ˙ φ∂ s H− H ˙ φ∂ s v + H ∂ s ˙ φ − ∂ s H ∂ s ˙ φ − ˙ φ∂ s v + O ( (cid:15) R ) . (6.21c)We also find the corrections to the constrained quanti-ties which were already computed at lowest order whenderiving the lowest order string model. Indeed, the ex-pressions (6.7), (6.8), and (6.9) need to be corrected withthe terms p (0) ⊃ R (cid:104) − ˙ φ − H ∂ s v − ∂ s H ∂ s v + ( ∂ s v ) + ν ∂ s K − H ∂ s v − ∂ s v (cid:3) (6.22a) v (1) ⊃ R (cid:2) H ∂ s v + H ∂ s H ∂ s v + H ( ∂ s v ) − ˙ φ∂ s ˙ φ + ∂ s v∂ s H + ν H ∂ s K + H ∂ s v + ∂ s H ∂ s v + ∂ s v∂ s v + ν ∂ s K + H ∂ s v − ∂ s v (cid:3) (6.22b)˙ φ (1) ⊃ R î ˙ φ∂ s H ∂ s v + H ∂ s H ∂ s ˙ φ + ∂ s ˙ φ∂ s H + H ˙ φ∂ s v + H ∂ s ˙ φ + ∂ s H ∂ s ˙ φ + ˙ φ∂ s v ó . (6.22c)9Once replaced in (6.2), we finally get the corrections to(6.10) which read as ∂ t v ⊃ R (cid:104) H ˙ φ + 12 H ∂ s v + H ∂ s H ∂ s v + H ( ∂ s v ) + ˙ φ∂ s ˙ φ + 3 ∂ s v∂ s H + ν H ∂ s K + H ∂ s v + 6 ∂ s H ∂ s v + ∂ s v∂ s v + H ∂ s v + ∂ s v (cid:3) (6.23a) ∂ t ˙ φ ⊃ R î ˙ φ∂ s H ∂ s v + H ∂ s H ∂ s ˙ φ + ∂ s ˙ φ∂ s H + H ˙ φ∂ s v + H ∂ s ˙ φ + ∂ s H ∂ s ˙ φ + ˙ φ∂ s v ó (6.23b) ∂ t ln R ⊃ R (cid:0) − H ∂ s v − ∂ s H ∂ s v − H ∂ s v − ∂ s v (cid:1) . (6.23c)We report in Appendix E the next order correctionswhich are of order (cid:15) R . Note also that surface tensioneffects do not enter explicitly the dynamical equation foraxial rotation (6.23b) nor the dynamical equation for theradius (6.23c), but they still matter due to the couplingsbetween v , ˙ φ and R .In general, if we keep terms of order (cid:15) nR in the expres-sion giving ∂ t v , that is, if we consider the n -th correction,then it involves terms which have 2( n + 1) order deriva-tives in the affine parameter s , typically from terms ofthe form R n ∂ n +1) s v , and the differential complexity isincreased. If we consider for instance a steady regimein which all time derivatives vanish, this can lead to arather stiff differential system as the coefficients in frontof the highest derivatives are typically the smallest. D. Comparison with the Cosserat model
As shown by Garc´ıa and Castellanos (1994) theCosserat model can be viewed as an averaged model.Indeed, given that the longitudinal velocity is given by v (0) + r v (1) + r v (2) + . . . it could be natural to consideran averaged velocity and derive the dynamical equationfor this variable and not for v (0) which should be con-sidered as a derived variable. The average longitudinalvelocity is simply defined as v ≡ R (cid:90) R V r d r = 2 R (cid:90) R (cid:88) n v ( n ) r n r d r = (cid:88) n R n n + 1 v ( n ) . (6.24)We obtain that from (6.20) and (6.1) the kinematic equa-tion for the time evolution of the fiber radius reads ex-actly (Garc´ıa and Castellanos 1994) ∂ t ln R = −H v − ∂ s v . (6.25)Obviously the expression (6.24) needs to be truncatedat a given power of R which is the order at which the equations are considered. Once truncated, we can invertit because from our algorithm, the v ( n ) are expressed interms of (the derivatives of) v (0) and ˙ φ (0) . After invert-ing, we obtain v = v (0) as a power series in v and itsderivatives. For instance up to the first corrections weget v = v + R (cid:0) H ∂ s v + ∂ s v (cid:1) + O ( (cid:15) R ) (6.26a) v = v − R (cid:0) H ∂ s v + ∂ s v (cid:1) + O ( (cid:15) R ) . (6.26b)At lowest order, the form of (6.10a) is the same if it is ex-pressed with v or with v because both velocities are equalat lowest order. Differences appear only when includingthe first set of corrections. Eventually, we find ∂ t v + R Å − H ∂ t ∂ s v − ∂ t ∂ s v ã (6.27)= g − ∂ s P (0) ν + 6 H ∂ s v − v∂ s v + 3 ∂ s v + R (cid:16) H ˙ φ − H ∂ s v − H ( ∂ s v ) + ˙ φ∂ s ˙ φ + ∂ s v∂ s H − H ∂ s v + H v∂ s v + ∂ s H ∂ s v + v∂ s v (cid:1) which we have written in a form which matches Eq. (55)of Garc´ıa and Castellanos (1994). However it does notmatch the Cosserat equation Eq. (53) of Garc´ıa andCastellanos (1994) [which is also Eq. (74) of Eggers(1997) derived with a Galerkin approximation method].As explained in details by Garc´ıa and Castellanos (1994),this is because some terms involving v (1) are removed.In the case of the Cosserat model of Eggers (1997), theseterms would probably be recovered when including higherorders of the Galerkin method, given that at lowest orderin the Galerkin approximation the information that thenext corrections are parabolic in nature has not been putin. In a sense, the lowest order of the Galerkin approx-imation also amounts to ignoring some parabolic termsof the type v (1) r . In our method, the contributions v (1) and v (2) are taken into account from the constraints (6.8)[corrected by (6.22b) and (6.21b)].Furthermore, our results are extended to include a pos-sible axial rotation from the inclusion of the ˙ φ ( n ) . We findthat when including the first corrections to the viscousstring approximations, the corrected dynamics of v cou-ples with the lowest order axial rotation ˙ φ = ˙ φ (0) whoseevolution needs then to be computed from the lowest or-der dynamical equation (6.10b). This interplay betweenlongitudinal velocity and axial rotation which takes placewhen the first corrections are included has already beendescribed in Bechtel et al. (1995), but our formalism al-lows already for more compact and geometrically moremeaningful expressions. However, it is only when con-sidering curved fibers that our formalism based on STFmultipoles appears to be powerful as we shall now see indetails in the following section.0 VII. APPLICATION TO CURVED FIBERSA. Overview of curved fiber specificities
Whenever we consider curved fibers, we must take intoaccount the property that the FCL curvature κ a does notvanish anymore, nor does the rotation rate ω i of the or-thonormal basis. As stressed in § IV J, even if we restrictto a stationary regime for which all time derivatives van-ish and thus ω a = 0, we can still generically form STFproducts of the type κ (cid:104) a . . . κ a (cid:96) (cid:105) (7.1)which might source the STF moments of order (cid:96) . Forinstance, terms proportional to κ (cid:104) a κ b (cid:105) would source theshear part of the velocity field “ V ab from the Navier-Stokesequation, and this would in turn induce a deformation ofthe fiber shape of the type R ab from the kinetic conditionon the fiber side (4.25). Since the lowest order descrip-tion corresponds to a string approximation for which thesection size and shape are irrelevant, these terms are ex-pected to arise only when including the first corrections.For instance, the combinations κ (cid:104) a κ b (cid:105) y a y b , or κ (cid:104) a κ b (cid:105) R (7.2)are of order (cid:15) R . For the lowest order description, onemight consider circular sections, but as soon as we con-sider refinements to this description, we must abandonthe circular shape assumption. As we shall find in theremainder of this section, the fundamental dynamicalvariables are the same as for the axisymmetric case( v = v (0) ∅ , ˙ φ = ˙ φ (0) ∅ , and R ) on which we add the variousshape multipoles R L and also the position and velocityof the FCL. If we consider a description up to order n inpowers of (cid:15) R , then we shall find that we must include atleast the multipoles R L with (cid:96) ≤ n . As soon as we leavethe realm of straight fibers, we open Pandora’s box andwe cannot obtain all STF devils from constraints, as theshape multipoles become dynamical.This departure from circular fiber sections is even moreobvious when considering the motion in a steady rotatingframe. Indeed, the inertial forces will bring typical con-tributions of the form Ω (cid:104) a Ω b (cid:105) , which are similar to tidalforces and induce an elliptic elongation. Such contribu-tion arises naturally in the Navier-Stokes equation (5.11)as can be seen from the general expression of the iner-tial forces (5.13), and they typically source the velocityshear moments “ V ( n ) ab , which in turns source the ellipticdeformation R ab from the boundary kinematics (4.25).We first derive the lowest order viscous string modelin the curved case in § VII B and discuss correctionsin § VII C. In § VII D we also consider the special caseof straight but non axisymmetric fibers with mild ellipticshapes.
B. Viscous string model
1. Incompressibility conditions
The incompressibility conditions on moments (4.5) areused to replace the moments u ( n ) L . In particular, from[ C E ] (0) ∅ , [ C E ] (1) ∅ , [ C E ] (0) a , [ C E ] (1) a , and [ C E ] (0) ab we obtainthe general conditions u (1) ∅ = − “ V (1) a (cid:101) κ a − (cid:101) κ a Ù V (0) a − ∂ s v (1) ∅ , (7.3a) u (0) a = − “ V (1) a − “ V (0) ab (cid:101) κ b − (cid:101) κ a u (0) ∅ − κ a ˙ φ − ∂ s v (0) a , (7.3b) u (1) a = − “ V (2) a − “ V (1) ab (cid:101) κ b − (cid:101) κ a u (1) ∅ − (cid:101) κ b Ù V (0) ab − κ a ˙ φ (1) − ∂ s v (1) a , (7.3c) u (0) ab = − “ V (1) ab − “ V (1) (cid:104) a (cid:101) κ b (cid:105) − (cid:101) κ (cid:104) a Ù V (0) b (cid:105) − ∂ s v (0) ab − “ V (0) abc (cid:101) κ c , (7.3d)that are used extensively together with (4.6) throughout § VII.
2. Normal vector and curvature
At lowest order, the unit normal vector componentsare simply (cid:99) N = −H R + O ( (cid:15) R ) (7.4a) (cid:99) N a = n a + O ( (cid:15) R ) (7.4b)and from (4.23) the scalar extrinsic curvature reads as K = 1 R + (cid:101) κ a n a + O ( (cid:15) R ) . (7.5)If we expand R H as in (5.7b), then [ R K ] (0) ∅ = 1 and[ R K ] (0) a = (cid:101) κ a . Note that we cannot separate the pressurecontribution from viscous and surface tension effects asin (6.4). The scalar extrinsic curvature must be useddirectly inside the boundary constraint (5.5).
3. Navier-Stokes equation
Using the components expression for the velocity gra-dient (4.31), the velocity shear (4.36), and the velocityderivatives (5.14), the lowest multipoles of the longitudi-nal and sectional components of acceleration and of thevolumic forces are [reminding that we are using the no-1tation (4.2) throughout] A (0) ∅ = ( ∂ t + v∂ s ) v + ∂ t ¯ U + U a (cid:101) ω a + I + 2 U a (cid:101) Ω a , (7.6a) f (0) ∅ = g + 4 v (1) + v (0) a (cid:101) κ a − κ a ω a − ∂ s P (0) ∅ (7.6b) − vκ a κ a + 2 (cid:101) κ a ∂ s “ V (0) a + “ V (0) a ∂ s (cid:101) κ a + ∂ s v, (cid:99) A (0) a = − v (cid:101) κ a + 2 ‹ U a Ω + I a − U (cid:101) Ω a − v (cid:101) Ω a − ¯ U (cid:101) ω a − v (cid:101) ω a + ∂ t U a + ‹ U a ω, (7.7a) “ f (0) a = g a − P (0) a + “ V (1) a − κ a ˙ φ − (cid:101) κ a ∂ s v − ∂ s v (0) a − v∂ s (cid:101) κ a − ∂ s (cid:101) ω a , (7.7b) ◦ A (0) = − v (0) a (Ω a + ω a + κ a v ) − Ä ˙ φ + Ω ä ∂ s v + ( ∂ t + v∂ s ) ˙ φ, (7.8a) ◦ f (0) = − “ V (1) a κ a − κ a κ a ˙ φ + 8 ˙ φ (1) − κ a ∂ s v (0) a − v (0) a ∂ s κ a + v (cid:101) κ a ∂ s κ a + ω a ∂ s (cid:101) κ a + (cid:101) κ a ∂ s ω a + ∂ s ˙ φ . (7.8b)In these expressions, we have defined the inertial force onthe FCL by I ≡ Ω × ( Ω × R ) . (7.9) I and I a are as usual its longitudinal and sectional com-ponents. Similarly, we recall that Ω a and g a are the sec-tional components of the steady frame rotation (if any)and of long range forces, whereas Ω and g are their longi-tudinal projections. Since Ω and g are constant vectors,then from (2.15) their components vary along the FCLaccording to ∂ s g = − (cid:101) κ a g a , ∂ s g a = κ a g , (7.10)with similar expressions for the components of Ω . How-ever I is not constant and we must use ∂ s I = − Ω − (cid:101) κ a I a , ∂ s I a = ΩΩ a + (cid:101) κ a I . (7.11)In reality, the expressions (7.6), (7.7) and (7.8) areformally more complex, as they also involve terms whichcontain “ V (0) a and “ V (0) ab . The former vanishes at lowestorder from the gauge condition (4.28), and when study-ing the structure of corrections in § VII C, we will showthat from the orthoradial boundary condition the lattervanishes as well at lowest order. The missing terms aregathered in § VII C 3.From (7.6) we can infer the lowest moment of the lon-gitudinal part of Navier-Stokes equation, that is D (0) ∅ ,which could also be obtained from the longitudinal pro-jection of a momentum balance equation. From (7.7)we can infer the lowest order multipoles of the sectional part of the Navier-Stokes equations D (0) a . It could also beobtained from the sectional projection of the momentumbalance equation. Finally, from (7.8), we infer ◦ D (0) whichgoverns the dynamics of the axial rotation rate ˙ φ = ˙ φ (0) .In the axisymmetric case, it was not strictly speakingpart of the lowest order string description since axial ro-tation decoupled at lowest order. For curved fibers, itseems at first sight when examining (7.7) that axial rota-tion retroacts on the dynamical equations for U a . How-ever, when including the boundary constraints it will ap-pear that this is not the case. As in the axisymmetriccase, the dynamical equation for the axial rotation is infact already a first correction and not part of the viscousstring model.In order to obtain closed dynamical equations from(7.6), (7.7), and (7.8), we need to find the lowest orderexpressions for v (0) a , P (0) = P (0) ∅ , v (1) = v (1) ∅ , “ V (1) a , and˙ φ (1) , which as in the axisymmetric case are going to beobtained from the lowest order of the boundary condi-tion (5.5).
4. Constraint equations
We recall the notation introduced in § V B for the vari-ous constraints obtained on the boundary and their mul-tipoles expansions. • First, from the lowest order monopole and dipoleof the longitudinal constraint, that is C (0) ∅ and C (0) a we get v (1) = H ∂ s v + ∂ s v + O ( (cid:15) R ) (7.12a) v (0) a = v (cid:101) κ a + (cid:101) ω a + O ( (cid:15) R ) . (7.12b)The constraint (7.12a) is exactly the same as theone obtained in the axisymmetric case (6.8) and wewill find that the effect of curvature appears onlyat higher orders. However, the latter constraint(7.12b) deserves a thorough comment. From theexpression of the total velocity (4.1) and the de-composition (3.10) for the relative velocity V , thevelocity field on the FCL (that is when y = y = 0)is V Cen = U + v T + “ V (0) a d a = U + v T + O ( (cid:15) R ) (7.13)since from the gauge condition (4.28) “ V (0) a is onlyan order (cid:15) R quantity. The rotation rate of the fluidon the central line is thus approximately ω a Cen ≡ [ T × ∂ s ( U + v T )] a = ω a + κ a v . (7.14)Furthermore, the dipolar component v (0) a corre-sponds to the sectional part of the solid rotationof the fluid contained in a fiber section. Indeed,2if we focus on terms which are linear in y a in thedecomposition (3.10), and if we ignore radial infall( u (0) ∅ ) and axial rotation ( ˙ φ (0) ), we realize that therelative velocity V contains V ⊃ v (0) a y a T = [ ω a Sec d a × y b d b ] , ω Sec a ≡ − (cid:101) v (0) a . The constraint (7.12b) thus states that the solidrotation of the fiber section ω a Sec is equal to therotation of the fluid on the FCL ω a Cen . It meansthat once the central line velocity U a is determined(together with the curvature and the longitudinallowest order velocity v ), then the rotation of thefluid on the FCL is determined, and the fiber sec-tion rotation must follow exactly the same rotationrate. This fact can be rephrased more rigorouslyby replacing the constraint (7.12b) in the expres-sion (4.37) for the vorticity, as we obtain (cid:36) a = ω a Cen + O ( (cid:15) R ) = ω a + κ a v + O ( (cid:15) R ) . (7.15)It is yet another way to see that the rotation of thefluid in fiber sections (vorticity) is guided by therotation of the fluid on the FCL.The consequence is that the fluid contained in agiven section, which by construction is orthogonalto the tangential direction T , remains always in ageometrically defined fiber section. Or, said differ-ently, the fluid particles belonging to different fibersections are not mixed by the fluid velocity. Wemust stress again that this result is only valid atlowest order in (cid:15) R .Hence, just by contemplating of (7.12b) we canunderstand why the lowest order approximation iscalled a string approximation. It is because thefiber sections are slaves of the FCL, as they are de-termined from it without retroaction at lowest or-der. Furthermore, and this is even more important,this type of behavior also corresponds to a form offlexible rod model since the fiber sections are notmixed and remain orthogonal to the FCL. If therewas no longitudinal velocity ( v = 0), a good anal-ogy would be the spinal column with the vertebrabeing the sections. If we want to consider a longi-tudinal velocity, a good analogy would be a collarmade of beads. The beads have a cylindrical holethrough which the collar string is passed, and if thebeads can slide along the string of the collar, theirorientation with respect to the string tangential di-rection is necessarily fixed thanks to the cylindricalhole. We understand already at that point that it ishopeless to try to find consistently the correctionsof the string model if we start from a flexible rodmodel, since the lowest order description is alreadya form of flexible rod model. • Second, from the lowest order monopole and dipoleof the orthoradial constraint, that is θ C (1) ∅ and θ C (0) a , we obtain˙ φ (1) = ( H (cid:101) κ a ω a + H ∂ s ˙ φ ) + O ( (cid:15) R ) , (7.16a) “ V (1) a = (2 κ a ˙ φ − (cid:101) κ a ∂ s v + 2 v∂ s (cid:101) κ a + 2 ∂ s (cid:101) ω a )+ O ( (cid:15) R ) . (7.16b) • And third, from the lowest order monopole anddipole of the radial boundary constraint, that isfrom r C (0) ∅ and r C (0) a , and using the gauge condi-tion (4.28) and the previous constraints we get P (0) ∅ = νR − ∂ s v + O ( (cid:15) R ) (7.17a) P (0) a = − κ a ˙ φ + (cid:101) κ a ∂ s v − v∂ s (cid:101) κ a − ∂ s (cid:101) ω a + ν (cid:101) κ a R + O ( (cid:15) R ) . (7.17b)
5. Dynamics of the string model
Inserting the constraints of the previous section in theexpressions (7.6), (7.7), and (7.8), we finally obtain thesystem of equations ∂ t v = − ∂ t ¯ U − U a (cid:101) ω a + g + ν H R − v∂ s v + 6 H ∂ s v + 3 ∂ s v + O ( (cid:15) R ) , (7.18a)( ∂ t U ) a = g a + v (cid:101) κ a + 2 v (cid:101) ω a − (cid:101) κ a ∂ s v − ν (cid:101) κ a R + O ( (cid:15) R ) , (7.18b) ∂ t ˙ φ = 4 H (cid:101) κ a ω a + ˙ φ∂ s v + ω a ∂ s (cid:101) κ a + 4 H ∂ s ˙ φ − v∂ s ˙ φ + ∂ s ˙ φ + O ( (cid:15) R ) , (7.18c)where using (2.25a) we used the compact expression( ∂ t U ) a = ∂ t U a + ω ‹ U a − ¯ U (cid:101) ω a . (7.19)Note also that from (2.24a) we also get ( ∂ t U ) = ∂ t U + (cid:101) ω a U a allowing to rewrite (7.18a) in a slightly more com-pact form if desired. Finally, when surface tension effectsare included, we also need to determine the dynamics ofthe fiber radius as it couples to (7.18b) and (7.18a), andat lowest order it reads exactly as in the axisymmetriccase, that is ∂ t ln R = −H v − ∂ s v . (7.20)Several comments are in order for this viscous stringmodel. • We can check that there is no retroaction of ˙ φ onthe lowest order dynamical equations for v and U a .The axial rotation dynamical equation (7.18c) is infact part of the first corrections and not part of thelowest order string model. • The equations are at most linear in the curvature κ a even though we have not linearized in this variable.3 • To compare (7.18a) and (7.18b) with the result ob-tained from a momentum balance equation by Arneet al. (2009, 2015), Ribe (2004), Ribe et al. (2006),we must first use that the average velocity insidea fiber is approximately the velocity on the FCLgiven at lowest order by (7.13). Then, from theproperties( ∂ t V Cen ) = ∂ t v + ( ∂ t U ) (7.21a)( ∂ t V Cen ) a = ( ∂ t U ) a − v (cid:101) ω a (7.21b) v ( ∂ s V Cen ) = v∂ s v (7.21c) v ( ∂ s V Cen ) a = − v (cid:101) κ a − v (cid:101) ω a , (7.21d)and defining a convective derivative by D t ≡ ∂ t + v∂ s , the string equations are simply recast as( D t V Cen ) (cid:39) g + ν H R + 6 H ∂ s v + 3 ∂ s v (7.22a)( D t V Cen ) a (cid:39) g a − (cid:101) κ a ∂ s v − ν (cid:101) κ a R . (7.22b) • We can in particular check that the contributionsfrom the surface tension are exactly those foundin Eqs. (37) (49) and (57) of Ribe et al. (2006).Since the vector (cid:101) κ a points toward the exterior ofthe FCL curvature, the contribution − ν (cid:101) κ a /R tendsto unfold the viscous jet as strongly curved regionwill be accelerated toward the center of curvature. • If we ignore the surface tension contributions,and using the covariantization relations (2.26), theequation (7.22) can be recast in the compact form D t V Cen (cid:39) g + 1 πR ∂ s F , (7.23) F ≡ πR F (0) ∅ T , F (0) ∅ = 3( ∂ s v )(7.24)where the link with the momentum balance has nowbeen made obvious. The expression of the viscousforce F which is purely longitudinal has the sameorigin as in the axisymmetric case [see Eq. (6.12)].
6. Rotating frame
If we now consider that the problem is studied in asteady rotating frame, then Eqs. (7.18) need to be sup-plemented with the contributions ∂ t v ⊃ − I − U a (cid:101) Ω a (7.25a)( ∂ t U ) a ⊃ − ‹ U a Ω − I a + 2( U + v ) (cid:101) Ω a (7.25b) ∂ t ˙ φ ⊃ − (cid:101) Ω a ( ω a + vκ a ) + Ω ∂ s v . (7.25c)The first two equations arise naturally from the longitudi-nal and sectional components of the inertial and Coriolisforces.
7. Longitudinal central line velocity
We notice that the dynamical equation for v (7.18a)contains ∂ t U . However, we must remember that the lon-gitudinal velocity of the central line contains a remain-ing gauge degree of freedom. Indeed, we have fixed theposition of the central line inside the fiber by askingthat there should be no shape dipole. However, as ar-gued in § IV I, we can still displace the fiber inside thatcurved central line, equivalent to a reparameterization s → s + f ( t ) which changes the velocities as U → U − ∂ t f and v → v + ∂ t f .Eq. (7.18a) is in fact an equation for ∂ t ( v + U ) and wemust find a unique way to determine U independently.Let us fix the value of the central line longitudinal veloc-ity U for a given affine parameter (say s = 0 at a fiberboundary) at all times. Typically the fiber is attached atthe boundary so we choose simply U ( s = 0) = 0. Thenfrom the condition (2.6) we obtain ∂ s U = (cid:101) κ a U a (7.26)and thus U must be determined everywhere on thefiber at all times, effectively breaking the degeneracy in(7.18a).
8. Full-set of equations for the string model
We are now in position to gather all the equationswhich are required for the curved viscous string model.First there is a set of equations which are constraint equa-tions and which must be solved at a given initial time,since they are ordinary differential equations in s .1. Once R ( s, t ) is known at a given time, e.g. at initialtime, then the tangential direction of the FCL isobtained from (2.2).2. The curvature κ of the FCL is then determined atthat same given time from its definition (2.7).3. The orthonormal basis d i , can be determined ev-erywhere on the fiber at that same given timefrom (2.15), provided some choice is made on a fiberboundary.4. With this orthonormal basis we can extract thecomponents κ a of the curvature.5. When the FCL velocity components U a and U areknown at a given time, and the curvature compo-nents κ a as well, then the rotation components ω i can be found from (2.28) for the sectional compo-nent and from (2.23b) for the longitudinal compo-nent at the same given time.6. The longitudinal part of the FCL velocity U is notdynamical but it is instead constrained by (7.26).4Then we have dynamical equations which give the timeevolution of variables from initial conditions, but theydepend on partial derivatives ∂ s so they are partial dif-ferential equations.1. The position in space of the FCL R ( s, t ) is modifieddue to the FCL velocity U , as seen on (2.4).2. The curvature components evolve in time thanksto (2.23a).3. The orthonormal basis evolves in time with (2.16).4. The longitudinal velocity v of the fluid inside thefiber evolves in time according to Eq. (7.18a).5. The sectional part of the FCL velocity U a evolvesaccording to Eq. (7.18b).6. The fiber radius R evolves according to (7.20).When considering corrections to the viscous string model,this structure between constraints and dynamical equa-tions is preserved.
9. Covariant expressions
If we prefer to work fully in the Cartesian canonicalbasis, that is, if for a vector we prefer using the covariantform X than the X i , then from (2.13) this is immedi-ate. However, we need to use the covariantization rela-tions (2.26) to recast the derivatives in the desired form ∂ t X or ∂ s X . The resulting equations take a more trans-parent form if we separate all vectors in sectional and lon-gitudinal parts according to (2.39) and write equationsfor these components. Since the longitudinal and sec-tional projections involve only the tangential direction T ,then the dependence in the section basis d a disappears.The FA coordinates and the orthonormal basis were in-troduced to perform all intermediary computations, butthe final results need not be expressed in this language.This justifies a posteriori why we have chosen the specialchoice (2.7) for the curvature and discarded the possibil-ity of having a non-vanishing longitudinal component forthe FCL curvature. Indeed, this would lead to the samefinal equations when expressed in a covariant form. How-ever, all intermediate computations would be more in-volved because the expression of the metric (2.35) wouldbe much more complicated and non-diagonal, and theChristoffels (2.44) would also be more complex. In par-ticular, as a consequence of this choice for curvature, thepartial derivatives with respect to s are easily written ina covariant form since for any vector X we deduce from(2.25b) the property P ⊥ ( ∂ s X ⊥ ) = ( ∂ s X a ) d a . (7.27)This allows to read the covariant form of a given equationwritten in terms of sectional components nearly instantly. We gather in covariant form the complete set of equa-tions described in § VII B 8. First, the vectors ω and U are decomposed in sectional parts and longitudinal partsas in (2.39). The constraint equations are ∂ s R = T , (7.28a) ∂ s T = κ × T = − (cid:101) κ ⇔ κ = T × ∂ s T , (7.28b) ω ⊥ = T × ∂ s U = ∂ s ‹ U + κ U , (7.28c)0 = T · ∂ s ω ⇔ ∂ s ω = − ω · (cid:101) κ , (7.28d)0 = T · ∂ s U ⇔ ∂ s U = − U · (cid:101) κ . (7.28e)As for the dynamical equations, they are recast as ∂ t R = U , (7.29a) ∂ t κ = ∂ s ω + ω × κ ⇔ P ⊥ ( ∂ t κ ) = P ⊥ ( ∂ s ω ⊥ ) ,∂ t T = ω × T , (7.29b) ∂ t ln R = −H v − ∂ s v, (7.29c) ∂ t v = − ∂ t U − U · (cid:101) ω + g + ν H R − v∂ s v (7.29d)+6 H ∂ s v + 3 ∂ s v − I − Ω × U ) · T ,P ⊥ [ ∂ t U ⊥ ] = T × (cid:104) ( U + 2 v ) ω + (cid:16) v − ∂ s v − νR (cid:17) κ (cid:105) + [ g − I ] ⊥ + 2 T × (cid:2) ( v + U ) Ω − Ω U (cid:3) , (7.29e)where we should use that for any vector ∂ t X ⊥ = P ⊥ [ ∂ t X ⊥ ] + ( X ⊥ · (cid:101) ω ) T (7.30a) ∂ s X ⊥ = P ⊥ [ ∂ s X ⊥ ] + ( X ⊥ · (cid:101) κ ) T . (7.30b)
10. Stationary regime
In the stationary regime, all partial time derivativesvanish and the viscous fiber model takes a simpler form.The velocity of the fiber center also necessarily vanishesand U a = U = 0 as well as the rotation ( ω i = 0). In thatcase, it becomes an ordinary differential equation in theFCL parameter s . For completeness, we report here theset of stationary equations, and they read as ∂ s R = T , (7.31a) ∂ s ln R (cid:39) − ∂ s ln v, (7.31b)3 ∂ s v + 6 H ∂ s v − v∂ s v (cid:39) − ν H R + I − g, (7.31c) (cid:101) κ a (cid:16) v − νR − ∂ s v (cid:17) (cid:39) I a − g a − v (cid:101) Ω a , (7.31d) ∂ s d i = κ × d i . (7.31e)Eq. (7.31d) is used to determine the curvature κ a (but itcan become singular) and (7.31c) is used to integrate v along the FCL. Eq. (7.31b) is the statement that vR isconstant in a stationary regime, due to incompressibility.If written in covariant form, the last three equations5of (7.31) read as3 ∂ s v + 6 H ∂ s v − v∂ s v (cid:39) − ν H R + T · ( I − g ) , (7.32a) κ (cid:16) v − νR − ∂ s v (cid:17) (cid:39) T × [ g − I ] − v Ω , (7.32b) ∂ s T = κ × T = − (cid:101) κ , (7.32c)which is the standard form for the stationary curvedstring model in the literature. C. Beyond the string model
1. Limitations of the string model
Apart for surface tension effects, the spatial exten-sion of sections, that is the radius R , does not play anyrole in the dynamical equations, meaning that the inter-nal structure of the fiber has no impact on the dynam-ics. Furthermore, at lowest order the sectional part ofthe viscous forces which reduces to the components “ F (0) a vanishes. The viscous forces are purely longitudinal, aswould be the case in a string, thus justifying the nameof the approximation. If we want to consider a model forviscous fibers in which the size of spatial sections plays arole, we must necessarily consider higher order terms inthe parameter (cid:15) R .Since fiber sections rotate at the same angular velocityas the fluid located on the FCL [see Eq. (7.15)], it meansalso that a rod model in which sections are not mixedand remain orthogonal to the fiber tangential directioncannot be part of this higher order model. Hence, whenconsidering the lowest order of the angular momentumbalance equation, we do not obtain a dynamical equa-tion which gives the time evolution of the fiber sectionrotation as a function of viscous forces, but rather ob-tain a constraint on the viscous forces (more precisely ontheir sectional component) given that the fiber sectionrotation is already determined by the string model.For higher order models, the corrections of order (cid:15) R inEq. (7.15) will imply that fluid particles of different sec-tions will be mixed as a result of time evolution. Higherorder models must also necessarily take into account thatthe velocity of the fluid on the central line , is not exactlythe velocity of the central line . Indeed, Eq. (4.28) im-plies that there is a tiny shift between the two which isan order (cid:15) R correction.Finally, terms of the type (7.2) typically source theshape quadrupole R ab and a model restricted to circularsections cannot be sufficient when considering correctionsto the string limit. As we shall explain in this section,when including order (cid:15) nR effect, we must include all mul-tipoles R L with (cid:96) ≤ n . Since we are interested in the firstcorrections which are of order (cid:15) R , we consider quadrupo-lar shape moments thereafter.
2. Normal vector and extrinsic curvature
When including a first set of corrections, the compo-nents of the normal vector and the unit normal vectorare approximately given by N = −H R + H R ( n a (cid:101) κ a ) + O ( (cid:15) R ) (7.33a) N a = n a − R ⊥ ba R bc n c − R ⊥ ba R bcd n c n d + O ( (cid:15) R ) (7.33b) (cid:99) N = −H R + H R ( n a (cid:101) κ a ) + O ( (cid:15) R ) (7.33c) (cid:99) N a = n a − R Å H n a − ⊥ ba R bc n c ã (7.33d)+ R Ä − ⊥ ba R bcd n c n d + H n a n b (cid:101) κ b ä + O ( (cid:15) R ) . From (4.23), the extrinsic curvature can then be ob-tained. Expanding R K in moments as in (5.7a), the firstmoments which are used to compute the first set of cor-rections to the string model are[ R K ] ∅ = 1 − R (cid:0) H + κ a κ a + ∂ s H (cid:1) + O ( (cid:15) R ) (7.34a)[ R K ] a = (cid:101) κ a (cid:2) R (cid:0) H + κ b κ b + 2 ∂ s H (cid:1)(cid:3) + R Ä H ∂ s (cid:101) κ a − R ab (cid:101) κ b ä + O ( (cid:15) R ) (7.34b)[ R K ] ab = 3 R ab + κ (cid:104) a κ b (cid:105) + O ( (cid:15) R ) . (7.34c)
3. Fundamental dynamical equations
The general method to build higher order correctionsis similar to the axisymmetric case. The only differenceis that now a given equation will not just give a con-straint on the monopole but also on other moments. Aswe restrict to second order, we shall need to consider themonopole, dipole, and quadrupole of equations only. Theincompressibility conditions (7.3) are also used through-out to express any dependence in multipoles u ( n ) L in termsof other types of multipoles.We start from lowest moments of the Navier Stokesequation components ( D (0) ∅ , “ D (0) a ), which give the funda-mental dynamical equations for the longitudinal velocity v = v (0) ∅ , and the FCL sectional velocity U a . The evo-lution of the axial rotation ˙ φ = ˙ φ (0) has already beenfound in (7.18c) and we have argued that it should beconsidered as part of the first corrections. Eqs (7.6) and(7.7) were not given in full generality as we had removedthe contributions from “ V (0) a and “ V (0) ab which are order (cid:15) R quantities. The Navier-Stokes components (7.6) and(7.7) must be supplemented by the contributions A (0) ∅ ⊃ “ V (0) a Ä (cid:101) Ω a + (cid:101) ω a + v (cid:101) κ a + v (0) a ä (7.35a) f (0) ∅ ⊃ (cid:101) κ a ∂ s “ V (0) a + “ V (0) a ∂ s (cid:101) κ a (7.35b)6 Equation Variable Essential dependence D (0) ∅ (7.6) and (7.35) ∂ t v (0) v (1) ∅ , P (0) ∅ , v (0) a , (cid:98) V (0) a (cid:98) D (0) a (7.7) and (7.36) ∂ t U a (cid:98) V (1) a , P (0) a , v (0) a , (cid:98) V (0) a , (cid:98) V (0) ab R E ( ≤ ∅ (4.25) ∂ t R v (0) a , (cid:98) V (1) a , v (0) ∅ , v (1) ∅ R R E ( ≤ ab (4.25) ∂ t R ab (cid:98) V (0) ab , (cid:98) V (1) ab TABLE III: Fundamental dynamical equations, with the cor-responding variables and the main variables on which theirevolution depends. The dependencies for the last two equa-tions giving the evolution of radius and shape quadrupole aregiven only when order (cid:15) R terms are included. (cid:99) A (0) a ⊃ “ V (0) b “ V (0) ab − “ V (0) a Å ∂ s v + “ V (0) b (cid:101) κ b ã + v∂ s “ V (0) a − ε ab “ V (0) b Ä ˙ φ + ¯ ω − ä + ∂ t “ V (0) a (7.36a) “ f (0) a ⊃ “ V (0) b κ a κ b − “ V (0) a κ b κ b + “ V (0) ab (cid:101) κ b + ∂ s “ V (0) a to be fully general.We also need to consider the dynamics of the shape.As in the string model, we need to consider the monopoleof the boundary kinematics (4.25) to determine the evo-lution of the radius. The dipole of this equation hasalready been considered in (4.28) to fix the gauge and de-termine “ V (0) a . Furthermore, we also need to consider thequadrupole of the boundary kinematics equation (4.25)so as to determine the evolution of R ab . These dynami-cal equations need to be truncated at the required order.In Table III we summarize the essential dependencies ofthe fundamental dynamical equations which need to bedetermined from constraints.
4. General structure of constraint equations
The method follows essentially the same steps as inthe axisymmetric case. The constraints which were ob-tained at lowest order in § (VII B 4) need to be extendedto include order (cid:15) R corrections, so as to be replaced inthe fundamental dynamical equations. We follow thesame procedure except that all constraints are consid-ered up to a higher order. For instance, when derivingthe string model we considered the monopole of the ra-dial constraint at lowest order r C (0) ∅ to constrain P (0) ∅ ,and now we must consider r C (1) ∅ .However, just like when finding the corrections of theaxisymmetric case, the price to pay is that we introducenew dependencies. For instance from r C (1) ∅ we can obtaincorrections for the constraint which determines P (0) ∅ , butit involves P (0) ∅ R . This set of dependencies is summa-rized in Table IV.The solution to this problem follows exactly themethod used in the axisymmetric case. We use the higher Equation Variable Essential dependence r C ∅ P (0) v (0) and P (1) R , P (2) R . . . r C a P (0) a v (0) a and P (1) a R , P (2) a R . . . r C ab P (0) ab v (0) ab and P (1) ab R , P (2) ab R . . . C ∅ v (1) v (0) and v (2) R , v (3) R . . . C a v (0) a v (1) a R , v (2) a R . . . C ab v (0) ab v (1) ab R , v (2) ab R . . . θ C ∅ ˙ φ (1) ˙ φ (0) and ˙ φ (2) R , ˙ φ (3) R . . . θ C a (cid:98) V (1) a (cid:98) V (2) a R , (cid:98) V (3) a R . . . θ C ab (cid:98) V (0) ab (cid:98) V (1) ab R and (cid:98) V (2) ab R . . . Gauge fixing (cid:98) V (0) a v (0) a , (cid:98) V (1) a , v (0) . . . TABLE IV: Structure of boundary constraints. Each multi-pole of each constraint is used to determine a variable as afunction of other variables. In each case, we specify whichvariable is determined and report the essential variables onwhich it depends to emphasize the structure of the recursivemethod. All these constraints need to be truncated at a givenpower of R . order moments of the Navier-Stokes, which we consider asconstraints and not dynamical equations, together withthe incompressibility constraint (5.15) to remove thesenewly introduced variables. This is made possible since,as in the axisymmetric case, these equations have Lapla-cians which allow us to use the property (6.16).We thus need to follow a recursive algorithm which isvery similar to the one used in the axisymmetric case,but which is more involved since it involves the (cid:96) -th or-der STF moments when considering order (cid:15) (cid:96)R corrections.Furthermore, the structure of the recursive algorithm isslightly different for the dipole components since i) thegauge constraint already determines one dipolar moment( “ V (0) a ) and ii) one equation is used to determine the exter-nal variable U a . The corresponding set of dependenciesis summarized in Table V.
5. Quadrupoles of constraints
Let us first examine the quadrupoles of the constraintsfor which we need only the lowest order expressions.From the quadrupole of the orthoradial constraint at low-est order, that is, from θ C (0) ab we get that “ V (0) ab = R ï − “ V (1) ab − κ (cid:104) a (cid:101) κ b (cid:105) ˙ φ (7.37) − κ (cid:104) a κ b (cid:105) ∂ s v + κ (cid:104) a ( v∂ s κ b (cid:105) + ∂ s ω b (cid:105) ) (cid:3) + O ( (cid:15) R )7 Equation Variable Essential dependence D ( n +1) ∅ v ( n +2) ∅ P ( n +1) ∅ , v ( n +1) ∅ D ( n ) a v ( n +1) a v ( n ) a , P ( n ) a D ( n ) ab v ( n +1) ab v ( n ) ab , P ( n ) ab ◦ D ( n +1) ˙ φ ( n +1) P ( n +1) ∅ , ˙ φ ( n ) (cid:98) D ( n +1) a (cid:98) V ( n +2) a (cid:98) V ( n +1) a (cid:98) D ( n ) ab (cid:98) V ( n +1) ab (cid:98) V ( n ) abD C ( n ) ∅ P ( n +1) ∅ P ( n ) ∅ D C ( n ) a P ( n +1) a P ( n ) aD C ( n ) ab P ( n +1) ab P ( n ) ab TABLE V: Structure of dependence for the constraints ob-tained either from the higher moments of the Navier-Stokesequation or from the incompressibility constraint. We haveindicated only the essential dependence, so as to emphasizeclearly the structure of the recursive method, but it dependsin general on the full set of lower order variables. where we recall that the notation (cid:104) a . . . a n (cid:105) means theSTF part and the notation ( a . . . a n ) means the sym-metric part. Actually, the property that surface tensionterms start at order (cid:15) R and viscous terms start at order (cid:15) R arises for all “ V (0) L with (cid:96) ≥ r C (0) ab we find that P (0) ab depends on “ V (0) ab /R . Once thisdependence is replaced, we get P (0) ab = νR (cid:0) R ab + κ (cid:104) a κ b (cid:105) (cid:1) + κ (cid:104) a (cid:101) κ b (cid:105) ˙ φ (7.38)+ κ (cid:104) a κ b (cid:105) ∂ s v − κ (cid:104) a ( v∂ s κ b (cid:105) + ∂ s ω b (cid:105) ) + O ( (cid:15) R ) . Finally from the quadrupole of the longitudinal con-straint C (0) ab we get v (0) ab = 0 + O ( (cid:15) R ) . (7.39)
6. Monopoles and dipoles of constraints
We now examine the monopoles and dipoles of theconstraints. The lowest orders have been found alreadyin § (VII B 4) when deriving the viscous string model,and they need to be used to replace variables appearingin the order (cid:15) R terms. Once this is done, then from theconstraints C ( ≤ ∅ and C ( ≤ a we get that the constraints (7.12) should be supplemented by v (1) ⊃ R î − v (2) − v (1) a (cid:101) κ a + 6 H κ a κ a ∂ s v +84 H ∂ s H ∂ s v − H vκ a ∂ s κ a − H κ a ∂ s ω a +6 ∂ s v∂ s H + 84 H ∂ s v + 12 ∂ s H ∂ s v +48 H ∂ s v + 20 H ∂ s v + ∂ s v (cid:3) (7.40a) v (0) a ⊃ R î − v (1) a + 10 H κ a ˙ φ − H (cid:101) κ a ∂ s v +10 H v∂ s (cid:101) κ a + 10 H ∂ s (cid:101) ω a − (cid:101) κ a ∂ s v (cid:3) . (7.40b)From the constraint θ C ( ≤ ∅ we get an expression for ˙ φ (1) but we do not need any higher order contribution sincethe dynamical equation for ˙ φ (7.18c) is already part of thefirst corrections. However, from θ C ( ≤ a we get the firstcorrections for the constraint on “ V (1) a . Finally, from theconstraints r C ( ≤ ∅ and r C ( ≤ a we obtain the correctionson the constraints for P (0) ∅ and P (0) a . These are ratherlarge expressions, and we gathered them in Appendix F.As in the symmetric case, the constraint for v (1) nowdepends on v (2) when including order (cid:15) R corrections, andwe shall thus need another constraint to replace it. In factwe also need constraints for v (1) a , “ V (2) a , and “ V (1) ab . As inthe axisymmetric case, these additional constraints willcome from higher moments of the Navier-Stokes equa-tions and the incompressibility constraint (see Table V).Since the expressions of these constraints can be ratherlarge, we report them in Appendix F.
7. Corrections for dynamical equations
The evolution of the radius with the first correctionsincluded is given by ∂ t ln R ⊃ R ( − H ∂ s v − ∂ s H ∂ s v − H ∂ s v − ∂ s v ) . (7.41)It is striking that once all the constraints are properlyreplaced, we reach the same expression as in the axisym-metric case (6.23c), even though the original expressiondeduced from the monopole of (4.25) is formally morecomplex. In principle, we could apply the same methodas in § VI D and use the average longitudinal velocityinstead of v (0) , thus changing the fundamental variable.Concerning the corrections of the dynamical equationsfor v and U a , we report them in Appendix G. They areobtained from the replacement of all constraints, andthen the replacement of the lowest order dynamical equa-tion to remove the time derivatives appearing in the cor-rective terms. Several comments are in order here. • As for the lowest order string model, these equa-tions could be recast in a covariant form, that isusing the canonical Cartesian basis, following themethod described in § VII B 9. This is especiallystraightforward when using the property (7.27), sowe do not write it explicitly.8 • Odd powers of R correspond to surface tension ef-fects ( R − for the lowest order and R for the cor-rections), whereas even powers of R correspond toinertial and viscous effects ( R for the lowest orderand R for the first corrections). • We remark that the quadrupoles of shape R ab which appear from surface tension effects do notretroact on the dynamical equations for v and U a . • In practice, it proves easier to multiply the equa-tions considered [e.g., boundary constraint (5.5) orNavier-Stokes (5.11)] by h or h , so as to avoid un-necessary factors h − = 1 / (1 + (cid:101) κ a y a ) or h − whichwould need to be expanded in an infinite series in (cid:101) κ a y a . For instance, K given by (4.23) involves h − ,but it is not the case for h K . All products of tensorsfully contracted with vectors y a are then handledthanks to (A13). Hence, it is in principle possibleto find very general relations between multipoles,as we did for instance with the incompressibilitycondition (4.5).Note that the dynamical evolution of ˙ φ needs to beobtained from Eq. (7.18c). Similarly, the dynamical evo-lution of the quadrupole now needs to be determined in-dependently. From the quadrupole of (4.25) and usingthe constraint (7.37), we find the simple lowest order dy-namical equation( ∂ t + v∂ s ) R ab = − “ V (1) ab + R ab ∂ s v − ε ( ac R b ) c ( ˙ φ − ω ) , (7.42)in which one should replace the constraint (F7).Note that since the evolution of R ab is sourced by “ V (1) ab ,then from (F7) we see that it contains typical terms ofthe type κ (cid:104) a κ b (cid:105) ∂ s v . If we were to consider the dynam-ics of higher multipoles such as R L , it would be sourcedby terms of the type κ (cid:104) a . . . κ a (cid:96) (cid:105) ∂ s v among other terms.For comparison, from (7.20) we see that ln R is typicallysourced by ∂ s v . Hence, the higher the shape multipoleis, the more spatial derivatives are involved in its dy-namics. This justifies why when considering correctionsof order (cid:15) nR we need only to consider the multipoles R L with (cid:96) ≤ n . It also justifies a posteriori why we areworking with the shape multipoles R L (more preciselytheir dimensionally reduced variables R L ) as defined in(4.12) and not the “ R L of (4.13), which are better suitedto describe relative shape perturbations. D. Straight fibers with elliptic sections
It is now easy to consider the case of straight fibers butwith non-circular sections. We need only to consider thespecial case κ a = ω a = ω = U a = U = 0. We recover im-mediately the dynamical equations found in the axisym-metric case with the first corrections included [(6.10), (6.23c) and corrections (6.23a)]. However, we also ob-tain the dynamical evolution of the shape quadrupole,which evolves as an independent equation. Indeed, toretroact on the dynamics of v we would need terms ofthe type R ab R ab which would appear only when correc-tions of order (cid:15) R are included.Eq. (7.42) restricted to straight fibers takes formallythe same form, except that the constraint (F7) now needsto be also considered in that restriction when replaced.The last term of Eq. (7.42) is expected, as it just statesthat axial rotation ˙ φ will rotate the ellipticity, but thesectional rotation of the orthonormal basis ( ω ) must alsobe taken into account and subtracted. However it is onlyreally an effect of the choice of basis to measure com-ponents. Indeed we can rewrite it in a manifestly co-variant form following the method of § VII B 9. Defin-ing R µν ≡ R ab d aµ d bν the shape quadrupole evolution isgiven by P ⊥ αµ P ⊥ βν [( ∂ t + v∂ s ) R αβ ] = − “ V (1) ab d aµ d bν (7.43)+ R µν ∂ s v − φT β ε β ( µα R ν ) α , and we can check that the contribution from the or-thonormal basis rotation ω has disappeared.Let us now compute the evolution of the relative el-lipticity, that is, the evolution of the dimensionless mo-ments “ R ab = R R ab . By using the lowest order of thefiber radius evolution (7.41), we finally find that the di-mensionless quadrupole evolves according to( ∂ t + v∂ s ) “ R ab = − νR “ R ab − H v “ R ab (7.44) − ε ( ac “ R b ) c ( ˙ φ − ω ) . In a stationary regime ( ∂ t “ R ab = 0) we see that thereis a competition between surface tension effects whichtend to decrease ellipticity and stretching ( H ≤
0) whichtends to increase the relative contribution of ellipticity.Indeed if we assess the evolution of Q ≡ “ R ab “ R ab fromthe previous equation, its second line which is a purelyrotational effect does not contribute since ε ac “ R cb “ R ba =0, and we get simply( ∂ t + v∂ s ) Q = − νR Q − H vQ . (7.45)If it is true that stretching increases the ellipticity, thesurface tension effects encompassed in the first term onthe right hand side would eventually dominate and dampellipticity, as R is reduced by stretching. E. Alternative method for rotating frames
In our formalism we have allowed for the possibilityto be working in a rotating frame. This was taken intoaccount in the Navier-Stokes equation, by adding the fic-titious forces term (5.13). There is, however, a simpler9method to recover all expressions in a rotating frame.First we derive the results in a non-rotating frame, al-lowing to cut by approximately half the number of termsin the final results, and then we relate all variables inthe non-rotating frame with their counterparts in the ro-tating frame. We denote by G ω i the components of therotation rate of the orthonormal basis [defined by (2.16)]in the non-rotating frame, and by R ω i its counterpart inthe rotating frame. Similarly, we define G V i and R V i asthe velocity in the non-rotating frame and rotating framerespectively, and adopt a similar notation for the FCLvelocity components U i . These quantities are related by G V i = R V i + [ Ω × x ] i (7.46a) G U i = R U i + [ Ω × R ] i (7.46b) G ω i = R ω i + Ω i , (7.46c)where we recall that the wedge products are performedaccording to (2.14). In particular, when considering therelative velocity with respect to the FCL [ V µ definedin (4.1)], we note that the previous relations imply thatthe moments of its longitudinal and sectional part areunchanged except for G ˙ φ (0) = R ˙ φ (0) + Ω , G v (0) a = R v (0) a + (cid:101) Ω a . (7.47)In order to replace all variables referring to the non-rotating frame in terms of the variables referring to therotating frame, we also need to be able to relate the timederivatives of the FCL velocity components. Specifically,we need ∂ t G U i = ∂ t R U i (7.48a) − (cid:2) R ω × ( Ω × R ) (cid:3) i + [ Ω × R U ] i ∂ t Ω i = [ Ω × R ω ] i , (7.48b)where we emphasize that in these expressions, we areconsidering time derivatives of components. The rela-tion (7.48a) is obtained from the definition (2.4) (whichreads as here ∂ t R µ = G U µ ) and (7.46b), by using thedefinition [ Ω × R ] i ≡ d iµ ( Ω × R ) µ and the property(2.16). The relation (7.48b) is obtained from the defi-nition Ω i ≡ d iµ Ω µ and the property (2.16). Finally, wealso need the derivatives ∂ s (cid:0) [ Ω × R ] i (cid:1) = − (cid:101) Ω i − [ κ × ( Ω × R )] i . (7.49)Using (7.46), (7.47), and (7.48), we were able to checkthat starting from the expressions found in a non-rotatingframe (Ω a = Ω = 0), we recover the expressions valid ina rotating frame, for all constraints and all dynamicalequations. It is thus a healthy consistency check for thevalidity and correctness of the method. F. Physical insights on constraints
In this section we give physical interpretations for thevarious velocity components. For simplicity, we neglect the effect of curvature, so as to emphasize the physicaleffects more clearly.
1. Pressure constraint and the Trouton ratio
Let us first analyze how the constraint for v (1) (7.12)arises at lowest order. If we have a longitudinal velocitygradient ( ∂ s v (cid:54) = 0), then we have a radial infall whichis constrained by incompressibility as we get u = − ∂ s v from (4.6). As the radial infall is proportional to the ra-dial distance V a ⊃ uy a /
2, there is a radial gradient inthis radial velocity. Thanks to viscous forces, this createsa component τ ( µ ) ab ∝ µδ ab u . From the boundary condition(where we ignore surface tension for simplicity), there isnecessarily a pressure appearing to compensate for this, P (cid:39) µu . Hence on the fiber section, the longitudinalforce per unit area is F = τ = 2 µ∂ s v − µu = 3 µ∂ s v ,and we recover the standard Trouton enhancement ra-tio (Trouton 1906) 3 /
2. Indeed, in the end we get afactor 3 instead of the factor 2 that would have beenfound if we had forgotten the pressure contribution. Tosummarize, in order to satisfy the boundary constraint,a pressure must appear when we have a gradient in thelongitudinal velocity, and it transforms a 2 into a 3 inthe longitudinal viscous forces because this pressure alsoacts on the sections.
2. The Hagen-Poiseuille profile
If we now consider a gradient in the gradient of thelongitudinal velocity ( ∂ s v (cid:54) = 0), it will induce a gradientof radial infall ( ∂ s u = − ∂ s v ). This gradient of radial in-fall induces viscous forces per unit area applied on fibersections, of the form F a = τ a (cid:39) µ ( ∂ s u/ y a . We mustrealize that the component τ a gives the sectional com-ponents of the force per unit area applied onto the fibersections, but τ a n a gives also the longitudinal forces ap-plied on the boundary. A way to have a longitudinalcomponent of viscous forces on the fiber side, is to havea HP profile, that is a parabolic profile. Indeed, a ve-locity profile V = v (1) r creates τ a (cid:39) µ v (1) y a and if v (1) = − ∂ s u/ ∂ s v/
4, then the boundary constraint issatisfied. If the fiber radius is not constant along the fiber(
H (cid:54) = 0), this reasoning is slightly altered, thus explainingthe corresponding contribution in (7.12). To summarize,when we have a gradient of stretching, then we have agradient of infall velocity which creates a component τ a on the fiber side, which in turn needs to be canceled by aHP profile so as to satisfy the boundary conditions, andeventually when everything is taken into account, thereis no sectional component for the forces per unit area onsections ( F a (cid:39)
3. Dipolar longitudinal velocity
Let us consider a simple case in which there is no lon-gitudinal velocity v = v (0) ∅ , but we allow for a gradient inthe sectional components of the FCL velocity ( ∂ s U a (cid:54) = 0).This induces a force per unit area applied on sections ofthe form F a = τ a (cid:39) µ∂ s U a . Following the same rea-soning as in the previous section, we realize that the flowmust adapt in a way which creates an additional contri-bution for τ a in order to satisfy the boundary constraint.If we consider a dipolar modulation of the longitudinalvelocity V = v (0) a y a , then it creates a stress τ a = µv (0) a .If v (0) a = − ∂ s U a , then the boundary constraint is satis-fied, and given that we have considered v = v (0) ∅ = 0,this is equivalent to the constraint (7.12b) once we use(2.28). Physically, if neighbor sections slide along eachother (that is if they have a relative velocity which issectional), the viscous forces force them to rotate so thatthey do not slide along each other. To conclude, we mightthink that a gradient in the sectional velocity would cre-ate a sectional force per unit area on sections, but in factthe boundary condition would ensure that this is not thecase, just like in the previous section. Only when con-sidering higher order constraints will there be a sectionalviscous force per unit area on sections, and this is whythe lowest order model is a string model where viscousforces per area are necessarily longitudinal.
4. Dipolar pressure constraint
If we have a dipole modulation of the radial infall( u (0) a (cid:54) = 0), then following the reasoning of § VII F 1, itwill induce a dipole in the pressure so as to satisfy bound-ary constraints. Indeed the sectional velocity contains V a = y a / u (0) b y b ) and, thus, τ ab ⊃ µδ ab ( u (0) c y c ). In or-der to compensate for this component, we need a pressuregradient P (0) a = µu (0) a , which once all other constraintsare used gives (7.17b). To be fully correct, we should alsomention that u (0) a also contributes as τ ab ⊃ µy ( a u (0) b ) , butthere is also a parabolic sectional component V a = “ V (1) a r which contributes as τ ab ⊃ µy ( a “ V (1) b ) , and from the or-thoradial constraint they must cancel. Again, the bot-tom line is that at lowest order the system adapts sothat τ a (cid:39)
0, but this property also implies that thereis no sectional component for the viscous forces per unitarea on sections.
G. Comparison with rod models
The methods based on rod models follow a slightly dif-ferent logical route in order to obtain a one-dimensionalreduction. Indeed, in our method we solve for the con-straints and then use them to obtain the volumic forces. In rod models, we solve instead only some of the bound-ary constraints, those needed to get the pressure mo-ments, and then we use them to compute the forces perunit area applied on fiber sections. These are then in-tegrated to get total forces and total torques applied onsections, allowing to establish a momentum balance equa-tion and an angular momentum balance equation on aslice of fluid contained between two infinitesimally closesections. The information of the boundary constraintswhich was not used explicitly is then used implicitly be-cause we use that no force is applied on the sides of theinfinitesimal slice, so the balance equations involve onlythe forces and torques on the fiber sections. If surfacetension effects are considered, then they are of courseadded to the side of the infinitesimal slice. In this sectionwe collect the expressions for forces and torques that wefind with our formalism, so as to facilitate comparisonswith existing literature using rod models.
1. Viscous forces
The total force on sections is F (cid:39) (cid:90) r = Rr =0 F r d r d θ , (7.50)where we recall the definition (5.4) for the forces per unitarea on sections. This relation is approximate because wehave neglected the effect of non-circular sections. If wewere to take into account correctly the effect of a non-circular shape, then we would need to perform a changeof variables as in Appendix B. At lowest order the totallongitudinal force is given by F (cid:39) πR F (0) ∅ = πR (cid:16) − νR + 3 µ∂ s v (cid:17) . (7.51)Note that in order to write a momentum balance equationon an infinitesimal slice, we should also consider i) longdistance forces g in the bulk of the slice, and ii) theeffect of surface tension on the side of the slice [e.g., Eq.(49) of Ribe et al. (2006)]. This can also be computedby integrating all the volumic forces on the infinitesimalslice, as the total lineic force is obtained fromd F d s tot ≡ (cid:90) r = Rr =0 f hr d r d θ , (7.52)where we recall that the total volumic forces are givenby (5.9). The factor h , which mathematically is √ g ij with the metric (2.35), takes into account the fact thatif the fiber is curved, then for an infinitesimal slice thereis more fluid in the exterior of curvature and less in theinterior of curvature. At lowest order we get simplyd F d s tot = πR g + 2 ∂ s ( πνR T ) + ∂ s ( F T ) . (7.53)The first term is the effect of gravity, the second is theeffect of surface tension on the side of the infinitesimal1slice, and the last is the net effect of viscous forces onsections. Using (7.51), the total lineic force is simplyd F d s tot = πR g + ∂ s ( πνR T ) + 3 µ∂ s ( πR ∂ s v T ) , (7.54)in agreement with the r.h.s. of (7.22). If we had ignoredthe effect due to the induced pressure on the sections[given by the first term of (7.51)], we would have overes-timated the total lineic force by a factor 2 as seen whencomparing (7.54) with (7.53). So we can state that iffor viscous forces the effect of the constrained pressuredis to enhance by a factor 3 /
2. Viscous torques
Still neglecting the effect of non-circular sections, thetotal torque applied on a fiber section is given by Γ (cid:39) (cid:90) r = Rr =0 ( y a d a ) × F d θr d r (7.55)or in components (that is, using Γ = Γ a d a + Γ T )Γ a (cid:39) − (cid:90) r = Rr =0 (cid:101) y a F d θr d r (7.56)Γ (cid:39) (cid:90) r = Rr =0 (cid:101) y b F b d θr d r . (7.57)In order to obtain the lowest order expressions for thetorque, we only need to keep contributions which are lin-ear in y a in the components of F , and we findΓ a (cid:39) πR ε ab F (0) b = − πR ‹ F (0) a , Γ (cid:39) πR ◦ F (0) . (7.58)Note that “ V (0) ab does not contribute to the torque as itcorresponds to a shear flow inside the section.The expression of torques takes a simple form whenexpressed in terms of fluid vorticity. At lowest order,the vorticity components obtained from (4.37) with thelowest order constraints replaced are given by (cid:36) = ˙ φ + O ( (cid:15) R ) (7.59a) (cid:36) a = ω a + κ a v + O ( (cid:15) R ) . (7.59b)The longitudinal component of the torque, is then foundfrom ◦ F (0) (cid:39) µ î ∂ s ˙ φ − κ a Ä (cid:101) ω a + v (0) a äó = µ Ä ∂ s ˙ φ + (cid:101) κ a ω a ä = µ ( ∂ s (cid:36) ) (7.60)where in the second equality we have used the lowest or-der constraint (7.12b) and in the third we used the prop-erty (2.24b). This component of the torque is induced by twisting (longitudinal difference of vorticity) and itsphysical origin is thus obvious.The sectional torque is found from − ‹ F (0) a (cid:39) − νR κ a + 3 µ î ∂ s ( ω a + κ a v ) − κ a ∂ s v − ˙ φ (cid:101) κ a ó = − νR κ a + 3 µ ( ∂ s (cid:36) ) a − µκ a ∂ s v . (7.61)The physical origin of the second term is simple. Thesectional component of vorticity corresponds to a rota-tion around a sectional axis. If we consider two neigh-bor sections which have different sectional vorticities, aswould happen if the fiber is bent, then the fluid locatedinside will be squeezed on one side and stretched on theother side, that is, there will appear a dipole of stretch-ing. Then from the viscous forces induced, this createsa sectional torque. With this naive view we would geta factor 2 and not a factor 3, but as in § VII F 1, thereis also a dipolar pressure which is induced to satisfy theboundary conditions, and it implies again the appearanceof the Trouton factor enhancement 3 / stretching implies the viscous force (7.51), twisting im-plies the longitudinal torque (7.60), and bending inducesthe second term of the sectional torque (7.61), and thesegeometries have been considered separately in these ref-erences even though they can have mixed effects. How-ever, the coupling of stretching with curvature has beenignored, and it happens to have an effect on the torquewhich is contained in the last term in (7.61). We stressthat this effect arises at the same order, and is not anorder (cid:15) R correction.In a simple case, the physical origin of this last termcan also be understood. Let us ignore surface tensionand consider a stationary regime with no axial rotation( ω a = ˙ φ = 0) and constant curvature ( ∂ s κ a = 0). Theexpression of the sectional force (7.61) is simply − ‹ F (0) a (cid:39) − µκ a ∂ s v . (7.62)As the fiber is stretched ( ∂ s v > u = − ∂ s v . Since the fiber is curved, the particles inthe exterior of curvature ( (cid:101) κ a y a >
0) are compressed whilethey move closer to the FCL, and conversely the particlesin the interior of curvature ( (cid:101) κ a y a <
0) are stretched whilethey move toward the FCL. As a result of viscous forces,this creates a contribution to the sectional torque. And,as usual, we get a factor 3 / (cid:101) κ a points in the ex-terior of curvature (or − (cid:101) κ a points toward the center ofcurvature of the FCL). This means that extrinsic cur-vature is increased in the exterior (for points such that (cid:101) κ a y a >
0) and decreased in the interior (for points suchthat (cid:101) κ a y a < i) the torque of long distance forces g in the bulk of theslice and ii) the contribution of surface tension on theside of the slice which also creates a torque [e.g., Eq. (50)of Ribe et al. (2006)]. Just as for the momentum bal-ance equation, the angular momentum balance equationis best computed by integrating the torques of volumicforces on the infinitesimal slice, as the total lineic torqueis obtained fromd Γ d s tot ≡ (cid:90) r = Rr =0 ( y a d a ) × f h d θr d r . (7.63)We find that it can be expressed asd Γ d s tot = T × F + d γ d s , (7.64)where the first term involves the sectional part of thetotal viscous forces defined in (7.50), and the second termis given at lowest order byd γ d s (cid:39) ∂ s Γ + πR (cid:16)(cid:101) κ × g + 4 νR H κ (cid:17) . (7.65)It appears clearly that the last term on the right-handside comes from surface tension effects on the side of theinfinitesimal slice, whereas the first term is the net ef-fect of torques applied on the sections of the infinitesimalslice. The middle term is the torque induced by gravity,which comes from the fact that when the fiber is curved,then the center of mass of an infinitesimal slice is notexactly on the FCL, but is instead offset by R (cid:101) κ /
4. Incomponents (7.65) reads as simply ï d γ d s ò a (cid:39) ∂ s Γ a − (cid:101) κ a Γ + νπR H κ a (7.66a) ï d γ d s ò (cid:39) ∂ s Γ + (cid:101) κ a Γ a , (7.66b)where the expressions of the torques applied on fiber sec-tions are given by (7.58) with (7.60) and (7.61). Forcompleteness, we report the explicit result which is ï d γ d s ò a (cid:39) πR (cid:16) gκ a + ν H κ a R − µ H (cid:101) κ a ˙ φ − µκ b κ b ω a +2 µκ a κ b ω b − µ H κ a ∂ s v − ν∂ s κ a R + 12 µ H v∂ s κ a + µ∂ s v∂ s κ a − µ ˙ φ∂ s (cid:101) κ a − µ (cid:101) κ a ∂ s ˙ φ + 12 µ H ∂ s ω a − µκ a ∂ s v + 3 µv∂ s κ a + 3 µ∂ s ω a (cid:17) (7.67) ï d γ d s ò (cid:39) πR (cid:16) − g a κ a − µκ a κ a ˙ φ + 8 µ H (cid:101) κ a ω a +3 µv (cid:101) κ a ∂ s κ a + 2 µω a ∂ s (cid:101) κ a + 8 µ H ∂ s ˙ φ +5 µ (cid:101) κ a ∂ s ω a + 2 µ∂ s ˙ φ (cid:17) . (7.68)
3. Sectional forces and dipolar HP profile
The component of the longitudinal velocity v (1) b r y b can be considered as a dipolar HP profile, as it is aparabolic profile with a dipolar modulation. This ve-locity component induces a force per area on sections F a = τ a (cid:39) µ y a ( v (1) b y b ) + µv (1) a r and it remains unde-termined by the boundary constraints. When averagedover directions it contributes to the forces per unit areaas F a (cid:39) µv (1) a r , and after integration over the wholesection, it contributes to the sectional part of the totalviscous force. Indeed, the sectional components of thetotal viscous force applied on sections [defined in (7.50)]are F a (cid:39) πR (cid:16) − µv (1) a + 3 µ H κ a ˙ φ + µκ b (cid:101) κ a ω b (7.69) − µ H (cid:101) κ a ∂ s v + µ ˙ φ∂ s κ a + 3 µ H v∂ s (cid:101) κ a + µ∂ s v∂ s (cid:101) κ a + µκ a ∂ s ˙ φ + 3 µ H ∂ s (cid:101) ω a − µ (cid:101) κ a ∂ s v + µv∂ s (cid:101) κ a − µκ b κ b (cid:101) ω a + µ∂ s (cid:101) ω a (cid:17) . As the rotation of sections is constrained by (7.12b),then an angular momentum balance equation would infact determine the value of v (1) a because part of the torquebalance equation (7.64) comes from T × F . In our methodit is determined from D (0) a (see Table. V), that is, fromthe dipole of the longitudinal part of the Navier-Stokesequation. Indeed it determines the rate of change of thelocal dipolar longitudinal velocity v (0) a which is related tothe vorticity and the local rotation rate of the fluid on theFCL [see the discussion which follows (7.12)], and it thuscontains the same information as the angular momentumbalance equation used in rod models. The expressionobtained for v (1) a is reported in Appendix F, and oncereplaced in Eq. (7.69), the sectional components of thetotal viscous force applied on sections read as F a (cid:39) πR (cid:16) − ν R ∂ s (cid:101) κ a + 3 µ H κ a ˙ φ + µκ b (cid:101) κ a ω b (7.70) − µ H (cid:101) κ a ∂ s v + µ ˙ φ∂ s κ a + 3 µ H v∂ s (cid:101) κ a − µ∂ s v∂ s (cid:101) κ a + µκ a ∂ s ˙ φ + 3 µ H ∂ s (cid:101) ω a − v ˙ φκ a − µ (cid:101) κ a ∂ s v + µv∂ s (cid:101) κ a − µκ b κ b (cid:101) ω a + µ∂ s (cid:101) ω a − ˙ φω a − µ (cid:101) κ b ω b κ a + (cid:101) κ a v∂ s v + (cid:101) ω a ∂ s v (cid:17) .
4. Rod models and their validity
As explained in the previous sections, the constitutiverelations of rod models are based on the determination offorces per unit area on fiber sections, and boundary con-straints are only used explicitly to determine the pressureprofile. The boundary constraints are then used implic-itly in balance equations.At lowest order, only the momentum balance equationis used, and we find (7.23) which is also exactly what is3found in the viscous string model. This equation deter-mines the motion of the FCL, and thus the rotation rate ω a can be inferred from it. Given that the fiber vorticityis constrained to match the rotation rate of the fluid lo-cated on the FCL (see § VII B 4), the angular momentumbalance equation is in fact used to determine the sectionalcomponents of the viscous forces per unit area F a , as ex-plained in the previous section, and it corresponds to theaddition of an order (cid:15) R correction.To summarize, the rod model amounts to usingEq. (7.23), formally exactly like the viscous string model,but it differs from it in the expression of the force usedon the right hand side which has sectional components.Hence the rod model corresponds to D t V Cen (cid:39) g + 1 πR ∂ s F , (7.71)with the force given by F ≡ ( πR ∂ s v + νπR ) T + F a d a , (7.72)and where we recall the definition of the convectivederivative D t ≡ ∂ t + v∂ s . The components of the left-hand side of Eq. (7.71) are obtained from Eqs. (7.21)exactly like for the string model. The components ofthe first term on the right-hand side are obtained fromEqs. (7.22) as in the string model. However, the rodmodel differs from the string model thanks to the sec-tional components of the total viscous force which arereported in Eq. (7.70). Note that from the relations of § II D we must use ∂ s ( F a d a ) = ∂ s F a d a + T ˜ κ a F a (7.73)so as to evaluate the last term in Eq. (7.71). Furthermore,the dynamics of ˙ φ is given by the lowest order dynamicalequation (7.18c) as in our model, and it is required since ˙ φ appears in Eq. (7.70). In the rod model of Ref. Ribe et al.(2006), this is equivalently found from the longitudinalpart of the momentum balance equation, even though itdoes not appear so explicitly as the physical case studiedis stationary.Finally, we already mentioned that for straight fibersthe model is improved by considering the full expressionof the boundary curvature Eggers and Dupont (1994)given by Eq. (6.11) instead of the lowest order 1 /R . Wecan use a similar ansatz for curved fibers by noting thatin Eq. (7.72), the term νπR is in fact νπ K ∅ R at lowestorder in (cid:15) R . Hence from Eq. (7.34a) an improved rodmodel is obtained by the replacement νπR → νπR [1 − R ( H + κ a κ a + ∂ s H )] (7.74)in Eq. (7.72). This improved rod model is necessary tocompute in Ref. Pitrou (2018) the Rayleigh-Plateau in-stability of a viscous fibers.The validity of rod models is limited by the fact thatwe are considering one correction while discarding othersources of corrections, whereas in our approach we sys-tematically consider all corrections of order (cid:15) R . Amongthe effects ignored there are the following: • the difference between the velocity of the centralline and the velocity of the fluid on the central line(see § IV I); • the HP profile induced by the constraint (7.12a)which mixes the fluid particles belonging to neigh-boring sections; • the shape moments which are sourced and invali-date the assumption that the sections remain cir-cular.Furthermore, even though it is computationally in-volved, our approach allows to find the corrections up toany order when rod models would fail because of the im-possibility to deal with the mixing sections. In principle,we have a clear recursive algorithm made of constraintsreplacements in fundamental dynamical equations.However, there are cases in which rod models capturethe essential corrections. First, shape moments appearin all intermediary expressions but do not appear in thefinal dynamical equations (G3) and (G4), hence, if we arenot interested in the fiber sections shape but only on thecentral line, they can be ignored. Furthermore, if we areconsidering the steady motion in a rotating frame, and ifthe Rosby number is very small, that is, for high rotationrates, then it does not matter if we have ignored most ofthe corrective effects. Indeed, the boundary constraintsdo not involve the frame rotation, and frame rotationenters essentially only in the determination of v (1) a from D (0) a , or equivalently in the determination of the totalsectional forces since it is related through (7.69). In thisregime, the fast rotation induces a sectional force andits expression should be captured correctly by an angu-lar momentum balance equation thanks to (7.64), pro-vided the last term in (7.61) is correctly included. In theend, rod models take only some corrective terms, but ifthe system is considered in a fast-rotating frame, theseretained corrective terms should also be the most im-portant ones, and rod models should lead to a reliableextension of the viscous string model. VIII. CONCLUSION
We have developed all the theoretical tools which arerequired to obtain a general one-dimensional descriptionof curved fibers. A concrete application for toroidal vis-cous fibers is presented separately in Pitrou (2018). Froma theoretical point of view, our 2 + 1 splitting, our useof fiber adapted coordinates, and more importantly ourparametrization of the velocity field in terms of STF ten-sors allow for a clear discussion about constraints anddynamical equations. It avoids the cluttered componentby component expressions which are usual in such con-text (Dewynne et al. 1992), and it bears a more trans-parent geometrical meaning since all quantities are an-alyzed in terms of their monopole, dipole, quadrupole,4and higher order multipoles. We find that it is the nat-ural language which allows to overcome the complexityof equations for curved viscous fibers. Indeed, the use ofan adapted formalism is the key to understand in depthapparently complex problems. From a practical or com-putational point of view, this STF based approach isvery powerful as it is possible to handle tensors withappropriate abstract tensor packages, and to this endwe used xAct (Mart´ın-Garc´ıa 2004). The correspondingnotebooks are available upon request from the author.The main results of this article are the following: • We have recovered the standard results of axisym-metric fibers at lowest order in (6.3) and (6.10),including also axial rotation. • The first corrections for this model are collected in(6.23). • We extended these results to include a second set ofcorrections and these can be found in Appendix E. • The main purpose of this article was to develop aformalism for curved fibers and we first rederivedthe viscous string model whose central equationsare (7.18), (7.25) and (7.20). Its covariant formu-lation is summarized in (7.28) and (7.29). • We found the first corrections for curved fibers infull generality in (G3), (G4), and (7.41), and theseare relevant when (cid:15) R is not so small since they areof order (cid:15) R . • Elliptic shape perturbations are sourced at that or-der and their dynamics is governed at that orderby (7.42) with (F7) replaced. • In particular, when restricting to straight fibers,the dynamical equation for the evolution of ellipticshape perturbations takes the simple form (7.44). • Finally, when comparing with rod models methods,we have exhibited a missing term in the expression(7.61) for the torque applied on fiber sections.
Acknowledgments
I would like to thank G. Faye for his help on the ir-reducible representations of SO(2), and R. Gy and F.Vianey for their encouragement to write this article. Ialso thank J. Eggers and N. Ribe for comments on earlierversions of this article. This research was initiated whenthe author was working for Saint-Gobain Recherche.
REFERENCES
C. Pitrou, Phys. Rev. E , 043116 (2018), 1801.10127. J. Plateau, Statique exp´erimentale et th´eorique des liquidessoumis aux seules forces mol´eculaires (Paris, Gauthier-Villars, 1873).J. W. S. Rayleigh, Proc. R. Soc. London , 4 (1878).J. Eggers, Reviews of Modern Physics , 865 (1997).D. B. Bogy, Annual Rev. of Fluid Mech. , 207 (1979).J. Eggers and E. Villermaux, Reports on Progress in Physics , 036601 (2008).F. J. Garc´ıa and A. Castellanos, Physics of Fluids , 2676(1994).J. Eggers and T. F. Dupont, J. of Fluid Mech. , 205(1994), physics/0110081.S. E. Bechtel, K. D. Bolinger, J. Z. Cao, and M. G. Forest,SIAM J. Applied Math.
55 (1) , 58 (1995).J. M. Perales and J. M. Vega, Physics of Fluids , 112114(2010).A. M. Ga˜n´an-Calvo, C. Ferrera, and J. M. Montanero, J. ofFluid Mech. , 427 (2011).J. M. Montanero, M. A. Herrada, C. Ferrera, E. J. Vega, andA. M. Ga˜n´an-Calvo, Physics of Fluids , 122103 (2011).L. Vincent, L. Duchemin, and S. Le Diz`es, J. of Fluid Mech. , 220 (2014).S. E. Bechtel, M. G. Forest, D. D. Holm, and K. J. Lin, J. ofFluid Mech. , 241 (1988a).S. E. Bechtel, K. J. Lin, and M. G. Forest, J. of Non-Newtonian Fluid Mech. , 87 (1988b).V. M. Entov and A. L. Yarin, J. Fluid. Mech. , 91 (1984).A. L. Yarin, Free liquid jets and films: hydrodynamics andrheology (Longman, Wiley, 1993).A. L. Yarin,
Springer handbook of atomization and sprays, ch.2, pp55-73 (Springer, Heidelberg, 2011).J. N. Dewynne, J. R. Ockendon, and P. Wilmott, J. of FluidMech. , 323 (1992).L. J. Cummings and P. D. Howell, J. of Fluid Mech. , 361(1999).W. Arne, N. Marheineke, A. Meister, and R. Wegener,Berichte des Fraunhofer ITWM (2009).W. Arne, N. Marheineke, A. Meister, and R. Wegener, J. ofComput. Physics , 20 (2015).N. M. Ribe, Royal Society of London Proceedings Series A , 3223 (2004).N. M. Ribe, M. Habibi, and D. Bonn, Physics of Fluids ,084102 (2006).N. M. Ribe, M. Habibi, and D. Bonn, Annual Review of FluidMechanics , 249 (2012).B. Audoly, N. Clauvelin, P.-T. Brun, M. Bergou, E. Grinspun,and M. Wardetzky, J. of Comput. Physics , 18 (2013),1202.4971.M. Bergou, B. Audoly, E. Vouga, M. Wardetzky, and E. Grin-spun, ACM Transactions on Graphics , 116 (2010).U. Miyamoto, JHEP , 011 (2010), 1007.4302.R. M. Wald, General Relativity (Chicago Univ. Press, 1984).R. Courant and D. Hilbert,
Methods of mathematical physics- Vol.1; Vol.2 (New York: Interscience Publication, 1953,1953).K. S. Thorne, Rev. Mod. Phys. , 299 (1980).L. Blanchet and T. Damour, Phil. Trans. Roy. Soc. Lond. A320 , 379 (1986).L. Blanchet, Class. Quant. Grav. , 1971 (1998), gr-qc/9801101.C. Pitrou, Class. Quant. Grav. , 065006 (2009), 0809.3036.J. M. Mart´ın-Garc´ıa, (2004), URL . F. T. Trouton, Proc. Royal Soc. Lond. A , 426 (1906).J. D. Jackson, Classical Electrodynamics (Wiley, 1998).M. Nakahara,
Geometry, topology and physics (Taylor andFrancis, USA, 2003).
Appendix A: STF formalism1. Extraction of STF tensors
Integrals on directions are simply (cid:90) d θ π n a . . . n a n = (2 n − n )!! δ ( a a . . . δ a n − a n ) (cid:90) d θ π n a . . . n a n +1 = 0 . (A1)Here, δ ( a a . . . δ a n − a n ) means that the indices need tobe fully symmetrized among the (2 n − (cid:90) d θ π n a n b = 12 δ ab (A2) (cid:90) d θ π n a n b n c n d = 18 (cid:0) δ ab δ cd + δ ac δ bd + δ ad δ bc (cid:1) . (A3)We recall that in general for STF tensors, we can usethe multi-index notation K ≡ a . . . a k or L ≡ b . . . b (cid:96) .However, in order to avoid confusion on multi-indices inthis section and the next one, we use the weaker multi-index notation a K ≡ a . . . a k or b L ≡ b . . . b (cid:96) . Let usdefine I a L b L ≡ δ a (cid:104) b . . . δ a (cid:96) b (cid:96) (cid:105) , J a L b L ≡ (cid:15) a (cid:104) b δ a b . . . δ a (cid:96) b (cid:96) (cid:105) , (A4)where we recall that when indices are enclosed in (cid:104) . . . (cid:105) wemust take the symmetric trace-free part. Clearly, thesequantities are just related by (cid:15) a c I ca L − b L = J a L b L , (cid:15) a c J ca L − b L = − I a L b L . (A5)They have the interesting properties I ca L − cb L − = I a L − b L − I aa = 2 I ab (cid:15) ab = 0 (A6a) J ca L − cb L − = J a L − b L − J aa = 0 J ab (cid:15) ab = 2 , (A6b) I a L b L I b L a L = 2 , J a L b L J b L a L = 2 , I a L b L J b L a L = 0 . (A7)The tensors (A4) are used when computing the followingintegrals on direction vectors2 (cid:96) (cid:90) d θ π n a L n (cid:104) b K (cid:105) = δ (cid:96)k I a L b L (A8a)2 (cid:96) (cid:90) d θ π (cid:15) ac n c n a L − n (cid:104) b K (cid:105) = δ (cid:96)k J a L b L . (A8b) It is then immediate to show that for a scalar functionexpanded in STF tensors as S = (cid:88) (cid:96) S L n (cid:104) L (cid:105) (A9)then the STF moments can be extracted through S L = 2 − (cid:96) (cid:90) d θ π n (cid:104) L (cid:105) S . (A10)This type of integral is very well suited for a tensor com-puter algebra system such as xAct (Mart´ın-Garc´ıa 2004)since we need only to implement the rules (A1).
2. Products of STF tensor
As explained in § III B, STF tensors in two-dimensionsare irreducible representations of SO(2). If the tensorhas (cid:96) indices, then it is in the representation D (cid:96) . Whenwe have a product of two STF tensors of rank (cid:96) and (cid:96) (cid:48) ,it means we have the tensor product of the representa-tions D (cid:96) ⊗ D (cid:96) (cid:48) . This is not irreducible, but it can bedecomposed in irreducible representations. Given thatthe dimension of D (cid:96) is 2 (except for D for which thedimension is 1), then this tensor product is of dimension2 × D | (cid:96) − (cid:96) (cid:48) | ⊕ D (cid:96) + (cid:96) (cid:48) if (cid:96) (cid:54) = (cid:96) (cid:48) , or D (cid:96) + (cid:96) (cid:48) ⊕ D ⊕ D if (cid:96) = (cid:96) (cid:48) . When counting the dimensions,it is the statement that 2 × × A K and B L .If we assume first that k < (cid:96) , then A a K B b L = A (cid:104) a K B b L (cid:105) + A c K B c K (cid:104) b L − K I a K b K (cid:105) . (A11)Under this form, we have indeed decomposed the prod-uct into two irreducible parts D | k − (cid:96) | and D | k + (cid:96) | (that is2 × A c K B c K (cid:104) b L − K (cid:105) and A (cid:104) a K B b L (cid:105) .However, if the tensors are of equal rank ( k = (cid:96) ), thisgets slightly different since A a L B b L = A (cid:104) a L B b L (cid:105) + 12 A c L B c L I a L b L + 12 (cid:15) cd A cc L − B dc L − J a L b L . (A12)In that case we have decomposed the product as a sum oftwo scalar functions (both corresponding to the represen-tation D ) and an element of D (cid:96) (that is 2 × A c L B c L , (cid:15) cd A cc L − B dc L − and theSTF tensor A (cid:104) a L B b L (cid:105) .In both cases we get (with k ≤ (cid:96) ) A K B L y K y L = A (cid:104) K B L (cid:105) y K y L + r k k A K B Kc L − K y c L − K . (A13)6If we consider the expansion (3.2), we can first removethe traces in the L indices to recast the expansion as V a ( y i , t ) = ∞ (cid:88) (cid:96) =0 ∞ (cid:88) n =0 V ( n ) a (cid:104) L (cid:105) ( s, t ) y L r n . (A14)Each V ( n ) a (cid:104) L (cid:105) can be handled exactly as a product of tensors A a and B L . If (cid:96) > (cid:96) = 1 we must use (A12). Combining all the termsin the sum (3.2) we finally conclude that the decomposi-tion of a 2-vector field in terms of irreps is necessarily ofthe form (3.10). Appendix B: Alternate shape representation
We can consider the following STF moments M L ≡ (cid:96) (cid:90) y (cid:104) L (cid:105) ¯ ρ ( y b )d y b , (B1)which are integrals on the fiber section which should cap-ture its shape. ¯ ρ is a step function which is unity if thereis a fluid particle and vanishes otherwise. These momentsare built just like the material moments of constant den-sity extended objects, or like the electric moments of uni-formly charged extended objects (Jackson 1998). Thesemoments can be related to the moments of radial dimen-sions R L defined in (4.12), but the relation is non linear.To see this, let us change variables and define rescaledcoordinates z a ≡ y a “ R L n L y a = z a (1 + “ R L n L ) . (B2)The Jacobian of the transformation isd y b = J d z b , J = Ä “ R L n L ä (B3)and the integrals (B1) are recast as M L = 2 (cid:96) πR (cid:96) +2 (cid:96) + 2 (cid:90) Ä “ R M n M ä (cid:96) +2 n (cid:104) L (cid:105) d θ π . (B4)For the monopole we obtain M ∅ = πR (1 + ∞ (cid:88) (cid:96) =1 − (cid:96) “ R L “ R L ) (cid:39) πR , (B5)which is just the area of the section. If (cid:96) > M L (cid:39) πR (cid:96) +2 “ R L = 2 πR (cid:96) +1) R L , (B6)where (A8) was used to compute the integral (B4).The kinematic equation, giving the evolution in time ofthese moments is found from the conservation equationof the density function ¯ ρ∂ t ¯ ρ + V iR ∂ i ¯ ρ = ∂ t ¯ ρ + V R ∂ s ¯ ρ + V aR ∂ a ¯ ρ = 0 . (B7) Indeed integrating over directions as in (B1), and afterintegrations by parts, we get simply ∂ t M L = − (cid:90) y (cid:104) L (cid:105) ( V R ∂ s ¯ ρ − ¯ ρ∂ a V aR ) d y + (cid:96) (cid:90) y (cid:104) a L − V a (cid:96) (cid:105) R ¯ ρ d y . (B8) V aR needs to be expressed in terms of its multipoles. Tothis end, we first express it in terms of V a from the re-lation (4.10b), and then use the expansion (3.10) for V a .As for V R , we should first use that V R = V R /h , expand1 /h = 1 / (1 + (cid:101) κ a y a ) which brings increasing powers of (cid:101) κ a y a , and then we should relate it to V from (4.10b) soas to use the expansion (3.9) for V . The angular inte-grals can then be performed with (A8). The resultingdynamical equations for these shape multipoles M L arerather complicated because of the high powers in the cur-vature vector κ a , but they are linear in both the velocitymultipoles and the shape multipoles M L . Instead, thedynamical equation (4.25) was still linear in velocity mul-tipoles, but extremely non-linear in the shape multipoles R L . Hence, it is not surprising that the relation betweenthe two types of multipoles (B4) is very non-linear. Sincethe normal vector and thus the extrinsic curvature aremore easily expressed with the multipoles R L as seenon (4.21), we chose to work with these multipoles so asto be able to include surface tension effects. Appendix C: Velocity of the coincident point
Let us define the space-time Cartesian coordinates X ˆ µ = ( t, x µ ) with ˆ µ = 0 , , , Y ˆ ı = ( t, y ˆ ı ) with ˆ ı = 0 , , ,
3. Theneach coordinate system is a function of the other one,that is we have the functions X ˆ µ ( Y ˆ ı ) and Y ˆ ı ( X ˆ µ ) whichare related by ∂X ˆ µ ∂Y ˆ ı ∂Y ˆ ı ∂X ˆ ν = δ ˆ µ ˆ ν , ∂Y ˆ ı ∂X ˆ µ ∂X ˆ µ ∂Y ˆ = δ ˆ ı ˆ . (C1)In particular, using the first relation for ˆ µ = µ and ˆ ν = 0,and the second relation with ˆ ı = i and ˆ = 0 we get0 = ∂x µ ∂t (cid:12)(cid:12)(cid:12)(cid:12) y + ∂x µ ∂y ∂y i ∂t (cid:12)(cid:12)(cid:12)(cid:12) x = ∂x µ ∂t (cid:12)(cid:12)(cid:12)(cid:12) y + d iµ ∂y i ∂t (cid:12)(cid:12)(cid:12)(cid:12) x (C2a)0 = ∂y i ∂t (cid:12)(cid:12)(cid:12)(cid:12) x + ∂y i ∂x µ ∂x µ ∂t (cid:12)(cid:12)(cid:12)(cid:12) y = ∂y i ∂t (cid:12)(cid:12)(cid:12)(cid:12) x + d iµ ∂x µ ∂t (cid:12)(cid:12)(cid:12)(cid:12) y (C2b)and we recover (4.7). Appendix D: Cartan structure relation
In this appendix, we build on the four dimensional per-spective of the previous section. Let us define the space-time tetrad d ˆ ı = ( d , d i ) . (D1)7It is made from the spatial orthonormal basis on whichwe have added a time directed vector d ˆ µ = δ ˆ µ , ⇒ ∂ s d = ∂ t d = 0 . (D2)Let us define the infinitesimal rotation matrices[ J i ] jk ≡ η ijk ⇒ [ J i , J j ] = − η ijk J k (D3)where η ijk is the permutation symbol with η = 1 andwhere the sum on the index k is implied. We can definean operator valued (rotation valued) one-form in the fourdimensional classical space-time by Ω ≡ κ i J i d s + ω i J i d t . (D4)The components of this form in the Cartesian canonicalbasis are Ω = Ω ˆ µ d x ˆ µ Ω ˆ µ ≡ κ i J i δ µ + ω i J i δ µ . (D5)This form is clearly the connection form of the Cartanformalism since (still with the sum on k implied) ∂ s d j = [Ω ] jk d k ⇔ ∂ s d j = κ i [ J i ] jk d k ⇔ ∂ s d j = κ × d j ∂ t d j = [Ω ] jk d k ⇔ ∂ t d j = ω i [ J i ] jk d k ⇔ ∂ t d j = ω × d j which in a four dimensional perspective reads exactly asthe first Cartan structure equation ∂ ˆ µ d ˆ = [Ω ˆ µ ] ˆ ˆ k d ˆ k . (D6)Then, since the classical space-time is flat, the secondCartan structure equation reads as (Nakahara 2003)d Ω + Ω ∧ Ω = 0 . (D7)Expressed explicitly in a basis of two-forms d y ˆ ı ∧ d y ˆ =d y ˆ ı ⊗ d y ˆ − d y ˆ ⊗ d y ˆ ı , the terms of this equation read asd Ω = (cid:2) ∂ t κ i − ∂ s ω i (cid:3) J i d t ∧ d s Ω ∧ Ω = ω i κ j [ J i , J j ]d t ∧ d s = [ κ × ω ] i J i d t ∧ d s and therefore we recover the structure relation (2.18). Appendix E: Second set of corrections for axisymmetric viscous fibers
As discussed in § VI A, the dynamical equation (6.10b) is part of the first set of corrections, implying that (6.23b) isin fact part of the second set of corrections, so we only report the second corrections for ∂ t v and ∂ t ln R . The secondcorrections to the dynamical equation for v are ∂ t v ⊃ R (cid:104) − H ˙ φ − H ˙ φ ∂ s H + 48 H ∂ s v + H ∂ s H ∂ s v + H ( ∂ s H ) ∂ s v + H ( ∂ s v ) + H ∂ s H ( ∂ s v ) + H ( ∂ s v ) − H ˙ φ∂ s ˙ φ − ˙ φ∂ s H ∂ s ˙ φ − H ( ∂ s ˙ φ ) + H ∂ s v∂ s H + 12 ∂ s H ∂ s v∂ s H + ( ∂ s v ) ∂ s H + H ∂ s v + ˙ φ ∂ s v + H ∂ s H ∂ s v + ( ∂ s H ) ∂ s v + H ∂ s v∂ s v + ∂ s H ∂ s v∂ s v + ( ∂ s v ) ∂ s v + H ∂ s H ∂ s v + H ( ∂ s v ) − H ˙ φ∂ s ˙ φ − ∂ s ˙ φ∂ s ˙ φ + H ∂ s v∂ s H + ∂ s v∂ s H + H ∂ s v + H ∂ s H ∂ s v + H ∂ s v∂ s v + ∂ s H ∂ s v + ∂ s v∂ s v − ˙ φ∂ s ˙ φ + ∂ s v∂ s H + H ∂ s v + ∂ s H ∂ s v − ∂ s v∂ s v + ∂ s v + ν (cid:16) H ∂ s K + H ∂ s H ∂ s K + H ∂ s v∂ s K + ∂ s H ∂ s K + ∂ s v∂ s K + H ∂ s K − H v∂ s K + ∂ s H ∂ s K + ∂ s v∂ s K − H ∂ t ∂ s K + H ∂ s K − v∂ s K − ∂ t ∂ s K + ∂ s K (cid:17)(cid:105) . (E1)As for the radius evolution, it should be corrected at that order by ∂ t ln R ⊃ R (cid:104) − H ∂ s v − H ∂ s H ∂ s v − ( ∂ s H ) ∂ s v − H ( ∂ s v ) − ∂ s H ( ∂ s v ) + H ˙ φ∂ s ˙ φ + ( ∂ s ˙ φ ) − H ∂ s v∂ s H − H ∂ s v − H ∂ s H ∂ s v − H ∂ s v∂ s v − ∂ s H ∂ s v − ( ∂ s v ) + ˙ φ∂ s ˙ φ − ∂ s v∂ s H− H ∂ s v − ∂ s H ∂ s v − ∂ s v∂ s v − H ∂ s v + ∂ s v − ν (cid:16) H ∂ s K + ∂ s H ∂ s K + H ∂ s K + ∂ s K (cid:17)(cid:105) . (E2) Appendix F: Higher order constraints for curved fibers
From the boundary constraint (5.5) we get the additional contributions “ V (1) a ⊃ R (cid:16) − “ V (2) a − H v (1) a − ε ac “ V (1) bc κ b + H κ a ˙ φ + κ a κ b κ b ˙ φ + H κ b (cid:101) κ a ω b − H κ b κ b (cid:101) ω a + κ a ˙ φ∂ s H − H (cid:101) κ a ∂ s v − κ b κ b (cid:101) κ a ∂ s v − (cid:101) κ a ∂ s H ∂ s v − ∂ s v (1) a + H ˙ φ∂ s κ a + H v∂ s (cid:101) κ a + vκ b κ b ∂ s (cid:101) κ a + v∂ s H ∂ s (cid:101) κ a + H ∂ s v∂ s (cid:101) κ a + H κ a ∂ s ˙ φ + H ∂ s (cid:101) ω a + κ b κ b ∂ s (cid:101) ω a + ∂ s H ∂ s (cid:101) ω a − ∂ s (cid:101) κ a ∂ s v + H v∂ s (cid:101) κ a + H ∂ s (cid:101) ω a − H (cid:101) κ a ∂ s v − (cid:101) κ a ∂ s v (cid:17) (F1) P (0) ∅ ⊃ − νR (cid:0) H + κ a κ a + ∂ s H (cid:1) + R (cid:16) − P (1) ∅ − H ∂ s v − κ a κ a ∂ s v − ∂ s H ∂ s v + vκ a ∂ s κ a + κ a ∂ s ω a − H ∂ s v − ∂ s v (cid:17) (F2) P (0) a ⊃ νR (cid:16) − ε ac R bc κ b + H (cid:101) κ a + κ b κ b (cid:101) κ a + 2 (cid:101) κ a ∂ s H + H ∂ s (cid:101) κ a (cid:17) (F3)+ R (cid:16) − P (1) a − “ V (2) a + 4 H v (1) a + ε ac “ V (1) bc κ b − H κ a ˙ φ − κ a κ b κ b ˙ φ + H κ b (cid:101) κ a ω b − H κ b κ b (cid:101) ω a − κ a ˙ φ∂ s H + H (cid:101) κ a ∂ s v + κ b κ b (cid:101) κ a ∂ s v + (cid:101) κ a ∂ s H ∂ s v + ∂ s v (1) a − H ˙ φ∂ s κ a − vκ b (cid:101) κ a ∂ s κ b − H v∂ s (cid:101) κ a − vκ b κ b ∂ s (cid:101) κ a − v∂ s H ∂ s (cid:101) κ a + H ∂ s v∂ s (cid:101) κ a − H κ a ∂ s ˙ φ − κ b (cid:101) κ a ∂ s ω b − H ∂ s (cid:101) ω a − κ b κ b ∂ s (cid:101) ω a − ∂ s H ∂ s (cid:101) ω a + H (cid:101) κ a ∂ s v + ∂ s (cid:101) κ a ∂ s v − H v∂ s (cid:101) κ a − H ∂ s (cid:101) ω a + (cid:101) κ a ∂ s v (cid:17) . From the higher order of the Navier-Stokes equation that we consider as constraints, we get v (2) = νR (cid:0) − H − H κ a κ a − κ a ∂ s κ a + ∂ s H (cid:1) − H ( ∂ s v ) − gκ a κ a + H ( v ) κ a κ a + H vκ a Ω a + H Ω a Ω a + U a κ b (cid:101) κ a Ω b + κ a κ a I + κ a ˙ φ (cid:101) Ω a − κ a Ω (cid:101) Ω a + U a κ a κ b (cid:101) Ω b + H vκ a ω a − (cid:101) κ a ˙ φω a + (cid:101) κ a Ω ω a + H Ω a ω a + H ω a ω a + (cid:101) κ a ω a ¯ ω + U a κ b κ b (cid:101) ω a − H κ a κ a ∂ s v + vκ a κ a ∂ s v + κ a Ω a ∂ s v + κ a ω a ∂ s v + H ∂ s H ∂ s v + ( v ) κ a ∂ s κ a + (cid:101) κ a ˙ φ∂ s κ a + v Ω a ∂ s κ a + vω a ∂ s κ a − κ a ∂ s v∂ s κ a + v∂ s κ a ∂ s κ a + ˙ φ∂ s ˙ φ + Ω ∂ s ˙ φ + vκ a ∂ s ω a + Ω a ∂ s ω a + ω a ∂ s ω a + ∂ s κ a ∂ s ω a + κ a κ a ∂ t ¯ U − κ a ∂ t ω a + H ∂ s v − κ a κ a ∂ s v − ∂ s v∂ s v + vκ a ∂ s κ a + κ a ∂ s ω a + H ∂ s v + ∂ s v (F4) v (1) a = ν R ∂ s (cid:101) κ a + vκ a ˙ φ + vκ a Ω + ˙ φ Ω a + ΩΩ a + ˙ φω a + Ω ω a − κ b (cid:101) κ a ω b + κ a (cid:101) κ b ω b + κ b κ b (cid:101) ω a + H (cid:101) κ a ∂ s v − v (cid:101) κ a ∂ s v − (cid:101) Ω a ∂ s v − (cid:101) ω a ∂ s v − ˙ φ∂ s κ a + ∂ s v∂ s (cid:101) κ a − κ a ∂ s ˙ φ + (cid:101) κ a ∂ s v − v∂ s (cid:101) κ a − ∂ s (cid:101) ω a (F5) “ V (2) a = νR (cid:0) − ε ac R bc κ b + H (cid:101) κ a + κ b κ b (cid:101) κ a − (cid:101) κ a ∂ s H − H ∂ s (cid:101) κ a + ∂ s (cid:101) κ a (cid:1) − ε ac “ V (1) bc κ b − Ω (cid:101) κ a − ( ˙ φ ) (cid:101) κ a − ( ∂ s v ) (cid:101) κ a − ( v ) κ b κ b (cid:101) κ a κ a κ b κ b ˙ φ − (cid:101) κ a ˙ φ Ω + (cid:101) κ a Ω − vκ b (cid:101) κ a Ω b − (cid:101) κ a Ω b Ω b + vκ b κ b (cid:101) Ω a + κ b Ω b (cid:101) Ω a + H κ b (cid:101) κ a ω b − vκ b (cid:101) κ a ω b − (cid:101) κ a Ω b ω b + κ b (cid:101) Ω a ω b − (cid:101) κ a ω b ω b − H κ b κ b (cid:101) ω a + vκ b κ b (cid:101) ω a + κ b Ω b (cid:101) ω a + κ b ω b (cid:101) ω a − H v (cid:101) κ a ∂ s v + κ b κ b (cid:101) κ a ∂ s v − κ a ˙ φ∂ s v + κ a Ω ∂ s v − H (cid:101) Ω a ∂ s v − H (cid:101) ω a ∂ s v + (cid:101) κ a ∂ s H ∂ s v + v ˙ φ∂ s κ a + v Ω ∂ s κ a − vκ b (cid:101) κ a ∂ s κ b − (cid:101) κ a ω b ∂ s κ b + κ b (cid:101) ω a ∂ s κ b + vκ b κ b ∂ s (cid:101) κ a − κ b ω b ∂ s (cid:101) κ a + H ∂ s v∂ s (cid:101) κ a − v∂ s v∂ s (cid:101) κ a + H κ a ∂ s ˙ φ − vκ a ∂ s ˙ φ − Ω a ∂ s ˙ φ − ω a ∂ s ˙ φ − ∂ s κ a ∂ s ˙ φ + ˙ φ∂ s ω a + Ω ∂ s ω a − κ b (cid:101) κ a ∂ s ω b + κ b κ b ∂ s (cid:101) ω a − ∂ s v∂ s (cid:101) ω a + H (cid:101) κ a ∂ s v − v (cid:101) κ a ∂ s v − (cid:101) Ω a ∂ s v − (cid:101) ω a ∂ s v + ∂ s (cid:101) κ a ∂ s v − ˙ φ∂ s κ a + ∂ s v∂ s (cid:101) κ a − κ a ∂ s ˙ φ + (cid:101) κ a ∂ s v − v∂ s (cid:101) κ a − ∂ s (cid:101) ω a (F6) “ V (1) ab = νR (cid:0) R ab + κ (cid:104) a κ b (cid:105) (cid:1) + κ (cid:104) a κ b (cid:105) (cid:2) ( v ) + ∂ s v (cid:3) + κ (cid:104) a (cid:101) κ b (cid:105) ˙ φ − vκ (cid:104) a ∂ s κ b (cid:105) − κ (cid:104) a ∂ s ω b (cid:105) + (cid:0) ω (cid:104) a + Ω (cid:104) a (cid:1) (cid:0) ω b (cid:105) + Ω b (cid:105) (cid:1) + 13 vκ (cid:104) a (cid:0) ω b (cid:105) + Ω b (cid:105) (cid:1) . (F7)From the incompressibility constraint (5.15) we also obtain P (1) ∅ = ν R (cid:0) −H − κ a κ a + ∂ s H (cid:1) (F8)+ Ω + ( ˙ φ ) − ( ∂ s v ) + ( v ) κ a κ a + ˙ φ Ω + vκ a Ω a + vκ a ω a + Ω a ω a + ω a ω a − κ a κ a ∂ s v + vκ a ∂ s κ a + κ a ∂ s ω a + ∂ s vP (1) a = νR (cid:0) − ε ac R bc κ b + H (cid:101) κ a + κ b κ b (cid:101) κ a − (cid:101) κ a ∂ s H + H ∂ s (cid:101) κ a − ∂ s (cid:101) κ a (cid:1) (F9) − Ω (cid:101) κ a − ( ˙ φ ) (cid:101) κ a + ( ∂ s v ) (cid:101) κ a − ( v ) κ b κ b (cid:101) κ a − κ a κ b κ b ˙ φ − (cid:101) κ a ˙ φ Ω − vκ b (cid:101) κ a Ω b − vκ b (cid:101) κ a ω b − (cid:101) κ a Ω c ω c − (cid:101) κ a ω c ω c + vκ b κ b (cid:101) ω a + κ b Ω b (cid:101) ω a + κ b ω b (cid:101) ω a + H v (cid:101) κ a ∂ s v + κ b κ b (cid:101) κ a ∂ s v − κ a ˙ φ∂ s v − κ a Ω ∂ s v + H (cid:101) Ω a ∂ s v + H (cid:101) ω a ∂ s v + v ˙ φ∂ s κ a + v Ω ∂ s κ a − vκ b (cid:101) κ a ∂ s κ b − vκ b κ b ∂ s (cid:101) κ a − v∂ s v∂ s (cid:101) κ a − vκ a ∂ s ˙ φ − Ω a ∂ s ˙ φ − ω a ∂ s ˙ φ + ∂ s κ a ∂ s ˙ φ + ˙ φ∂ s ω a + Ω ∂ s ω a − κ b (cid:101) κ a ∂ s ω b − κ b κ b ∂ s (cid:101) ω a − ∂ s v∂ s (cid:101) ω a + v (cid:101) κ a ∂ s v + (cid:101) Ω a ∂ s v + (cid:101) ω a ∂ s v − ∂ s (cid:101) κ a ∂ s v + ˙ φ∂ s κ a + ∂ s v∂ s (cid:101) κ a + κ a ∂ s ˙ φ − (cid:101) κ a ∂ s v + v∂ s (cid:101) κ a + ∂ s (cid:101) ω a . Appendix G: Higher order corrections for curved fibers
Once the boundary constraints of § VII C 6 and Appendix F are replaced, the dynamical evolution for the longitu-dinal velocity and the FCL velocity are given by ∂ t v = g + ν H R − I − U a (cid:101) Ω a − U a (cid:101) ω a + 6 H ∂ s v − v∂ s v − ∂ t ¯ U + 3 ∂ s v + νR ( H + H κ a κ a + H ∂ s H + κ a ∂ s κ a + ∂ s H ) + R (cid:16) Ω H − v (2) ∅ + H v κ a κ a − v (1) a (cid:101) κ a + H ˙ φ + 2 H ˙ φ Ω + 2 H vκ a Ω a + κ a ˙ φ (cid:101) Ω a + 2 H vκ a ω a + (cid:101) κ a ˙ φω a + (cid:101) κ a Ω ω a + 2 H Ω a ω a + H ω a ω a + 12 H ∂ s v − H κ a κ a ∂ s v + vκ a κ a ∂ s v + κ a Ω a ∂ s v + κ a ω a ∂ s v + 21 H ∂ s H ∂ s v − H ( ∂ s v ) + H vκ a ∂ s κ a + v κ a ∂ s κ a + (cid:101) κ a ˙ φ∂ s κ a + v Ω a ∂ s κ a + vω a ∂ s κ a + κ a ∂ s v∂ s κ a + ˙ φ∂ s ˙ φ + Ω ∂ s ˙ φ + H κ a ∂ s ω a + vκ a ∂ s ω a + Ω a ∂ s ω a + ω a ∂ s ω a + 3 ∂ s v∂ s H + 18 H ∂ s v − κ a κ a ∂ s v + 6 ∂ s H ∂ s v − ∂ s v∂ s v + vκ a ∂ s κ a + κ a ∂ s ω a
0+ 6 H ∂ s v + ∂ s v (cid:17) + O ( (cid:15) R ) (G1)( ∂ t U ) a = g a − ν (cid:101) κ a R + ( v ) (cid:101) κ a − U a Ω − I a + 2 ¯ U (cid:101) Ω a + 2 v (cid:101) Ω a + 2 v (cid:101) ω a − (cid:101) κ a ∂ s v + νR ( ε ac R bc κ b − H (cid:101) κ a − κ b κ b (cid:101) κ a − (cid:101) κ a ∂ s H − H ∂ s (cid:101) κ a − ∂ s (cid:101) κ a )+ R (cid:16) − “ V (2) a − H v (1) a − ε ac “ V (1) bc κ b − Ω (cid:101) κ a − v κ b κ b (cid:101) κ a + 12 H κ a ˙ φ + κ a κ b κ b ˙ φ − (cid:101) κ a ˙ φ − (cid:101) κ a ˙ φ Ω − vκ b (cid:101) κ a Ω b − vκ b (cid:101) κ a ω b + H κ a (cid:101) κ b ω b − (cid:101) κ a Ω b ω b − (cid:101) κ a ω b ω b + κ b Ω b (cid:101) ω a + κ b ω b (cid:101) ω a + 3 κ a ˙ φ∂ s H − H (cid:101) κ a ∂ s v + H v (cid:101) κ a ∂ s v − κ b κ b (cid:101) κ a ∂ s v − κ a ˙ φ∂ s v − κ a Ω ∂ s v + H (cid:101) Ω a ∂ s v + H (cid:101) ω a ∂ s v − (cid:101) κ a ∂ s H ∂ s v + (cid:101) κ a ( ∂ s v ) − ∂ s v (1) a + 4 H ˙ φ∂ s κ a + v ˙ φ∂ s κ a + v Ω ∂ s κ a + vκ b (cid:101) κ a ∂ s κ b − vκ a (cid:101) κ b ∂ s κ b + 12 H v∂ s (cid:101) κ a + 3 v∂ s H ∂ s (cid:101) κ a + H ∂ s v∂ s (cid:101) κ a − v∂ s v∂ s (cid:101) κ a + H κ a ∂ s ˙ φ − vκ a ∂ s ˙ φ − Ω a ∂ s ˙ φ − ω a ∂ s ˙ φ + ∂ s κ a ∂ s ˙ φ + ˙ φ∂ s ω a + Ω ∂ s ω a + ¯ ω∂ s ω a + κ b (cid:101) κ a ∂ s ω b − κ a (cid:101) κ b ∂ s ω b + 12 H ∂ s (cid:101) ω a + 3 ∂ s H ∂ s (cid:101) ω a + ∂ s v∂ s (cid:101) ω a − ∂ s (cid:101) κ a ∂ t v − κ a ∂ t ˙ φ − H (cid:101) κ a ∂ s v + v (cid:101) κ a ∂ s v + (cid:101) Ω a ∂ s v + (cid:101) ω a ∂ s v − ∂ s (cid:101) κ a ∂ s v + ˙ φ∂ s κ a + 4 H v∂ s (cid:101) κ a − v ∂ s (cid:101) κ a + ∂ s v∂ s (cid:101) κ a + κ a ∂ s ˙ φ + 4 H ∂ s (cid:101) ω a − v∂ s (cid:101) ω a + (cid:101) κ a ∂ s ∂ t v − ∂ s ∂ t (cid:101) ω a − (cid:101) κ a ∂ s v + v∂ s (cid:101) κ a + ∂ s (cid:101) ω a (cid:17) + O ( (cid:15) R ) . (G2)If the higher moments of Navier-Stokes equation and the incompressibility constraint are then used to replace theunknown variables as summarized in Table V, and also if time derivatives in corrective terms are replaced using thelower order dynamical equations, we obtain the closed and final results with order (cid:15) R corrections included ∂ t v = g + ν H R − I − U a (cid:101) Ω a − U a (cid:101) ω a + 6 H ∂ s v − v∂ s v − ∂ t ¯ U + 3 ∂ s v + νR ( H + H κ a κ a + H ∂ s H + κ a ∂ s κ a + ∂ s H ) + R (cid:16) Ω H + H v κ a κ a + H ˙ φ + 2 H ˙ φ Ω + H vκ a Ω a − H Ω a Ω a + κ a ˙ φ (cid:101) Ω a + κ a Ω (cid:101) Ω a + H vκ a ω a − (cid:101) κ a ˙ φω a − (cid:101) κ a Ω ω a + H Ω a ω a + H ω a ω a + 12 H ∂ s v − H κ a κ a ∂ s v + 2 vκ a κ a ∂ s v + 2 κ a Ω a ∂ s v + 2 κ a ω a ∂ s v + H ∂ s H ∂ s v + H ( ∂ s v ) + H vκ a ∂ s κ a − v κ a ∂ s κ a + (cid:101) κ a ˙ φ∂ s κ a − v Ω a ∂ s κ a − vω a ∂ s κ a − κ a ∂ s v∂ s κ a − v∂ s κ a ∂ s κ a + ˙ φ∂ s ˙ φ + Ω ∂ s ˙ φ + H κ a ∂ s ω a − vκ a ∂ s ω a − Ω a ∂ s ω a − ω a ∂ s ω a − ∂ s κ a ∂ s ω a + 3 ∂ s v∂ s H + H ∂ s v − κ a κ a ∂ s v + 6 ∂ s H ∂ s v + ∂ s v∂ s v + vκ a ∂ s κ a + κ a ∂ s ω a + H ∂ s v + ∂ s v (cid:17) + O ( (cid:15) R ) (G3)( ∂ t U ) a = g a − ν (cid:101) κ a R + ( v ) (cid:101) κ a − U a Ω − I a + 2 ¯ U (cid:101) Ω a + 2 v (cid:101) Ω a + 2 v (cid:101) ω a − (cid:101) κ a ∂ s v + νR ( − H (cid:101) κ a − κ b κ b (cid:101) κ a − (cid:101) κ a ∂ s H − H ∂ s (cid:101) κ a − ∂ s (cid:101) κ a )+ R (cid:16) v κ b κ b (cid:101) κ a + 12 H κ a ˙ φ − H vκ a ˙ φ − (cid:101) κ a ˙ φ − H vκ a Ω − (cid:101) κ a ˙ φ Ω − (cid:101) κ a Ω − H ˙ φ Ω a − H ΩΩ a + vκ b (cid:101) κ a Ω b + (cid:101) κ a Ω b Ω b − κ b Ω b (cid:101) Ω a − vκ a κ b (cid:101) Ω b − H ˙ φω a − H Ω ω a + vκ b (cid:101) κ a ω b + 2 H κ a (cid:101) κ b ω b + (cid:101) κ a Ω b ω b − κ b (cid:101) Ω a ω b + (cid:101) κ a ω b ω b + 3 κ a ˙ φ∂ s H − H (cid:101) κ a ∂ s v + H v (cid:101) κ a ∂ s v − κ b κ b (cid:101) κ a ∂ s v − κ a ˙ φ∂ s v − κ a Ω ∂ s v + H (cid:101) Ω a ∂ s v + H (cid:101) ω a ∂ s v − (cid:101) κ a ∂ s H ∂ s v + (cid:101) κ a ( ∂ s v ) + 6 H ˙ φ∂ s κ a − v ˙ φ∂ s κ a − v Ω ∂ s κ a + vκ b (cid:101) κ a ∂ s κ b + (cid:101) κ a ω b ∂ s κ b − κ b (cid:101) ω a ∂ s κ b + 12 H v∂ s (cid:101) κ a + κ b ω b ∂ s (cid:101) κ a + 3 v∂ s H ∂ s (cid:101) κ a − H ∂ s v∂ s (cid:101) κ a + v∂ s v∂ s (cid:101) κ a + 8 H κ a ∂ s ˙ φ − vκ a ∂ s ˙ φ − Ω a ∂ s ˙ φ − ω a ∂ s ˙ φ + 2 ∂ s κ a ∂ s ˙ φ − ˙ φ∂ s ω a − Ω ∂ s ω a + κ b (cid:101) κ a ∂ s ω b + κ a (cid:101) κ b ∂ s ω b + 12 H ∂ s (cid:101) ω a + 3 ∂ s H ∂ s (cid:101) ω a + ∂ s v∂ s (cid:101) ω a − H (cid:101) κ a ∂ s v + v (cid:101) κ a ∂ s v + (cid:101) Ω a ∂ s v + (cid:101) ω a ∂ s v − ∂ s (cid:101) κ a ∂ s v + ˙ φ∂ s κ a + 6 H v∂ s (cid:101) κ a − ∂ s v∂ s (cid:101) κ a + κ a ∂ s ˙ φ + 6 H ∂ s (cid:101) ω a − (cid:101) κ a ∂ s v + v∂ s (cid:101) κ a + ∂ s (cid:101) ω a (cid:17) + O ( (cid:15) R ) ..