Optimal control of plant root tip dynamics in soil

This paper aims to propose a novel approach to model the dynamics of objects that move within the soil, e.g. plants roots. One can assume that external forces are significant only at the tip of the roots, where the plant’s growth is actuated. We formulate an optimal control problem that minimises the energy spent by a growing root subject to physical constraints imposed by the surrounding soil at the tip. We study the motion strategy adopted by plant roots to facilitate penetration into the soil, which we hypothesize to be a circumnutation movement. By solving the proposed optimal control problem numerically, we validate the hypothesis that plant roots adopt a circumnutation motion pattern to reduce soil penetration resistance during growth. The proposed formalisation could be applied to replicate such a biological behaviour in robotic systems, to adopt the most efficient strategy for autonomous devices in soil exploration.


Introduction
Plant roots grow into the soil driven by attractive targets (e.g. nutrients or water) and avoiding obstacles [1]. Under non-stressful biological and chemical conditions, the growth of the roots depends mainly on the soil strength and the presence of obstacles at the root tip [2].
Roots have developed morphological and mechanical strategies to reduce the soil strength and handily overcome the resistance of the surrounding environment [3]. Most of these strategies have been deeply investigated. For example, it is well known that roots change the tip diameter [4], move by growing at the tip level [5] and produce mucus at the apex [6] to grow into the soil.
It has also been hypothesised that a characteristic motion in the root's tip, called circumnutation, could be adopted to facilitate soil penetration and improve the seedling of the plant [7]. Circumnutation is an active rhythmical motion pattern (e.g., elliptical, circular, pendulum-like or zigzag-shaped) [8] observed for the first time by Darwin [9]. Several studies have investigated the origin and biological processes that guide this pattern of motion [8,10]. It is well known that the circumnutation is induced by different growth rates obtained by active swelling and deswelling of cells at opposite sides of the organ [10] and is believed to be driven by both gravitropism and internal periodic signals [11][12][13].
Recently, the ability of roots to move into the soil has inspired novel robotic technologies for soil exploration, penetration and monitoring [14,15]. Hence, to accurately replicate the growth dynamics of a plant in robotic devices, it is crucial to fully characterise the motion strategy employed by a root's apex and the role of the circumnutation in order to overcome the soil penetration resistance.
Nevertheless, circumnutation in roots is still a poorly understood process, mainly due to the lack of data. Experimental data collection is limited by two main factors. Firstly, the response of soil to the forces actuated by a growing root is not easy to characterise. Soil is a mixture of organic matter, minerals, gases, liquids and organisms, in which the arrangement of these distinct components defines the soil type, its texture and behaviour in response to external forces [16]. The detailed description of the dynamics of an object moving in soil remains an open problem to date [17]. Secondly, plants have evolved differently in response to the environment, resulting in a wide diversity of genotypes across species and phenotypes over the same species. Consequently, similar growth conditions can produce very different root behaviours [18,19]. Hence, given the complexity of soil composition and plant genotype and phenotype, the characterisation of a specific motion pattern such as root circumnutation remains challenging.
To understand root-soil interactions, previous studies have conducted experiments on plants using homogeneous media, such as agar or phytagel [20][21][22]. However, the results from these studies might not be sufficiently generalizable to natural soil conditions and consequently not readily translatable to plant-inspired robotic systems. In a more recent study [23], it has been proposed an experimental framework to analyse the circumnutation of a robotic tip into topsoil. Here, the experimental setting confirmed that circumnutation, when applied to a robotic tip travelling into the soil, reduces the forces to be exercised compared to straight penetration and, consequently, the energy spent by the plant-inspired robotic system.
Another way to study the root growth in real soil is represented by mathematical modelling. Models can help to investigate the root growth and motion by describing the physical root-soil interactions. For example, in [3], the existing mathematical models of root growth in soils are combined with the macroscopic observations of root behaviours to estimate the magnitude of the forces experienced by the root apex. In the approaches reviewed in [3], the forces acting on the root flank are usually omitted. In the same reference, the authors have also shown that the macroscopic root growth can be investigated by observing the water potential flow inside the cells. However, a quasi-static description of plant growth is not suitable for the study of dynamical interactions between roots and soil, such as the ones described by oscillatory motion patterns. Furthermore, models focussing on cellular growth cannot be directly translated to robotic applications because of the current mechanical and macroscopic nature of engineering systems for soil penetration. In view of the previous considerations, the study of specific movements, such as circumnutation, is still an open question and a great challenge for modelling. Investigations on the patterns of motion during root growth not only can provide effective tools to reveal new insights from the natural system but also it can lead to solutions fundamental in agriculture to improve the uptake, in civil engineering and soil sciences to stabilise slopes and in robotics to design efficient autonomous penetration devices.
In this paper, for the first time, we describe the dynamics of a root penetrating the soil by proposing a model based on optimal control methods. From a control theoretical viewpoint, the root tip is the system to be controlled, while the control is represented by the root elongation. In such a framework, the root moves under the action of a control signal, without a priori prescribing any direction and rate of elongation. As main model outcome, we obtain the con- In the elongation zone, there is fast growth obtained by cell elongation. Into the meristem, the growth is slower and obtained by cell division. Specialised cells are formed into the differentiation zone. The root cap protects the tip, by mucus and cell release. (B) A sequence of root growth into a soil whose particles are smaller than the tip size. The tip penetration affects the soil structure rearranging the soil particles that behave like mobile obstacles. The red arrows inside each particle suggest the resistance of soil particles to be shifted by the tip. The red thicker arrow into the tip describes the averaged action of soil particles on the root growth. The soil resistance is broken by the pushing force exercised by the root during growth.
trol function which minimises the energy dissipated by the elongating apex, taking into account the forces acting at the tip in soils with different densities. Here we propose a general and rigorous model that takes inspiration by plant behaviours, without being based on cellular processes. The model successfully agrees with biological observations [24,25], meets the efficiency requirements introduced in [23] and demonstrates that the optimal motion strategy adopted by plant roots follows an oscillatory pattern similar to circumnutation movements.
In the following, we first introduce the framework in which the system is moving and the forces involved in the motion. Then, we formalize the optimal control problem we aim to solve (section 2). We validate the model (section 3) and study the effects of soil compaction and tip shape on the energy consumption (section 4.1). We investigate the optimal trajectory emerging at different soil densities (section 4.2) and we finally conclude with discussions and open questions (section 5).

Model setting
In this section, we will describe the root movement into the soil by using a system of ordinary differential equations (ODEs) and set the analogy with a rootlike artificial device that is supposed to imitate the natural root. To this aim, let us first observe that the soil resistance is relevant just at the tip since the root grows at the apical level by cell elongation and division [26] ( figure 1(A)). Also, we considered that penetration tasks, even in artificial systems, are better accomplished by conical or parabolic probe shapes [17,27], similar to root's apexes. A device (natural or artificial) having a conical or parabolic tip can be described by the relevant parameters, namely the diameter d, the height h and the axis of symmetry a ∈ R 3 . The shape of the tip is assumed constant during all the time evolution, so that only the axis a = a(t) will depend on time, according to the position of the tip into the space.
Let x(t), v(t) ∈ R 3 and u(t) ∈ U ⊂ R 3 be the position, the velocity of the vertex and the acceleration, or control which drives the root growth, respectively (figure 2(A)).
Here, U : is the set of control values and u m , u M > 0 are the maximum deceleration and acceleration values that the device can produce, depending on its physical and mechanical constraints. We use g to denote the vertical downward gravitational acceleration.
Real soil is a complex mixture of particles (clay, silt and sand) and pores. The plant's root might move in a soil whose particles, or clusters of particles, have a volume that is similar or bigger than the size of the tip. In this case, each particle behaves like a possibly fixed obstacle that the tip cannot move. The motion pattern of the tip is mainly driven by the stochastic displacement of the soil pores, in the case in which they are sufficiently large for the tip to pass through. In this case, buckling instabilities can also occur. Buckling is an elastic instability, namely a passive mechanism induced by a compression force (that the root cannot oppose to) acting along the longitudinal axis of the tip [3,24].
A different framework takes place when the root grows in a possibly heterogeneous soil, but with the assumption that soil particles have a volume that is smaller than the size of the tip and are separated by pores that are not sufficiently large for the root to pass through ( figure 1(B)). In this case, each particle behaves like a mobile obstacle that offers a mechanical resistance which can be displaced. The penetration in this environment works only by displacing the particles (or analogously relocating the voids) and by compacting the soil around the tip [3]. We assume the latter to be the framework where our device is moving. Since the root is not able to perceive the fluctuations of the forces from individual particles [28], the action of the forces on the tip can be averaged. This assumption holds for many common soils in agriculture, such as flooded, potting, sand and sandy loam soils. Since the soil density is defined as the mass of soil particles divided by the volume occupied, we can assume the averaged resistance of particles is proportional to the soil density, which we will later denote with k [g cm −3 ]. The averaged resistance of soil particles will affect the root tip motion by reducing the pushing acceleration u(t). In particular, namedĀ(t) [-] the percentage reduction of u(t) experienced by the vertex at x(t) due to the mechanical resistance of the soil particles, the tip can move and penetrate the soil only if: The soil density k and the functionĀ(t) are not intrinsic properties of the soil. They depend on how the soil is handled, namely on the displacement of voids and soil particles. Let k max [g cm −3 ] be the maximum soil density that can be obtained by compacting the soil. From now on, by compacted soil, we mean a soil handled to offer the maximum density k max . For a compacted soil, it corresponds to a maximum percentage reduction A max [-] of u(t). Therefore, one can write: Since we are studying the case of mobile soil particles, we will assume A max < 1 and no buckling can occur. Being the soil a complex mixture of solid and fluid structures [29], in addition to the mechanical resistance, our device experiences a dissipative and nonconservative force, called drag force [30,31]. When an object immersed in a fluid has a speed |v| relatively high with respect to the fluid velocity, the drag force is proportional to |v| 2 [17]. On the other hand, at relatively slow speed, the drag force is linear with respect to v [32]. We will assume a linear dependence of the drag force on the velocity since the root moves relatively slow with respect to the water and air flows into the soil. Furthermore, the experiment proposed in [23], and here used to estimate model parameters, studies a robotic tip that moves with a circumnutating helical pitch comparable to real roots, thus the linearity assumption of the drag force with respect to v is still valid. Therefore, the drag accelerationā d (t) opposing to the tip motion depends on the velocity of the tip v and on the tapering and the lateral surface S l of the device exposed to the interaction with the medium [17]. These properties are expressed in the equation:ā where λ 2 is a parameter estimated in appendix A. The ratio d/(d + h) measures the tapering of the tip that results in a smaller experienced drag force [17]. In addition, for a conical tip of height h and radius r (figure 2(A)), S l = πr √ (h 2 + r 2 ), whereas for a parabolic tip, S l = (πr)((r 2 + 4h 2 ) (3/2) − r 3 )/(6h 2 ). S c l = 2πrh + πr 2 is the surface of a cylinder of height h and radius r that interacts with the soil. The ratio S l /S c l estimates the surface of the device exposed to the interaction with the medium [23]. Also, if the device can direct its tip towards any direction with respect to the vertical axis, the lateral surface interacting with the soil is at its minimum when the tip axis is vertically directed, while it reaches the maximum when the tip axis is horizontal (at 90 • from the vertical axis, ρ = h in figure 2(B)). In the work [23], the authors analysed the effects on the energy dissipated by the system at different inclinations of an artificial device's tip, which penetrates a sandy loam soil at different densities. From the authors' considerations, it emerges that the rotations of the tip might rearrange particles and act as a cavity expansion strategy to lighten the soil resistance around the tip [33]. Due to the random rearrangement of soil particles, it is not possible to accurately provide a complete analytical description of this phenomenon, especially if the tip can perform any possible motion pattern into the soil. To approximate the effects of the particle rearrangement on the soil resistance, we introduce the function R(v(t)) (neglecting the dependences on the parameters to not burden the notation): where is the velocity and λ 1 , μ are parameters whose values are reported in appendix A. The parameter μ 1 is used only to ensure the wellposedness of the evolution equation (4) for velocities such that |v| 1. The angle:  (1) reduces the action of R when either the tip lateral surface which interacts with the soil is maximum (namely if α → π/2 or |cos(α)| → 1), or when the soil is at its highest density (k → k max and no interstitial pores can be resized) or when the soil is extremely light (very low density, namely k → 0) such that the pore resizing does not provide evident benefits. In all these cases, as well as when v 1 = v 2 = 0 (namely straight penetration), one has R ≈ 1. On the other hand, when the motion of the tip is on an horizontal plane (v 3 → 0), one has R ≈ 0. This feature reflects the fact that a tip continuously moving on an horizontal plane encounters only a negligible resistance. The action of R will affect bothĀ(t) andā d (t) as follows: Therefore, the resulting system for the dynamics of the tip in the time interval [t 0 , T f ] is: (4) In the previous equations, the choice of a control u : [t 0 , T f ] → U will determine a specific trajectory x(t) with velocity v(t). Notice that the second equation in the system is a consequence of the Newton's second law: where the force per unit mass F/M * [N kg −1 ] depends on the control u(t).
Following the usual definition of work force × displacement, here we can compute the work performed by the control strategy u(t) as: is required to make the definition of the work W computed by the control u(t) dimensionally consistent (see appendix A for the units of measurement). Furthermore, the absolute value is necessary to take into account the non-conservative nature of the forces acting on the root tip, thus avoiding negative work values.

Optimal control
The main contributions of this paper are achieved by estimating the optimal trajectory for the tip that minimises the dissipated energy W. To be more precise, consider the control u : [t 0 , T f ] → U, the angle α : [t 0 , T f ] → [0, π/2] and the solution y(t) = (x, v)(t) of (4) for a given initial position and velocity y( t ∈ [t 0 , T f ] to emphasise the dependence of y(t) on the choices of the functions u(t), α(t). Furthermore, one can regard W as a function whose value depends on the choice of the control u(t), the angle α(t) and the related trajectory y(t) = y(u, α)(t). Accordingly, we write W = W(y, u, α, T f ). Let us also introduce the closed set T ⊆ R 3 , which represents the target that the device should reach during the soil penetration. For example, T could be a desired depth or, for roots, an underground pool of nutrients. The aim is to find both the optimal controlū(t) and the optimal slopē α(t) such that the trajectory solutionȳ =ȳ(t,ū,ᾱ) of (4) reaches the target T minimising the cost function W in the time interval [t 0 , T f ].
It is equivalent to solve the optimal control problem (P), defined as: (4), In the next sections, we will focus on the numerical solution of problem (P) by looking for the optimal couple (ū,ᾱ)(t) that minimises W (in reaching a given depth in the time [t 0 , T f ]) only in the case of a parabolic tip.
To numerically solve the optimal control problem (P), we will use the direct method, being more stable with respect to other methods [34]. More precisely, the time interval [t 0 , T f ] is divided in n subintervals. In any subinterval, the control is assumed to be a constant function. Starting from an initial guess of the control for any subinterval, the constrained optimisation method fmincon of Matlab is used to fix values of the control in each subinterval [35]. Both the dynamics and the cost function are converted in algebraic equations so that the resulting nonlinear programing problem can be solved by well established methods [36]. The optimisation variables are both the control u and the angle α. To speed up the convergence of the method, the initial guess can be chosen similar to a circumnutating control with a constant slope α (since by evidence in section 4, the straight penetration is not expected to be the optimal solution). Namely, u 0 ∈ R 3,n , α 0 ∈ R n are such that:

Parameter estimation
To validate the model, we verify if the dynamics described by equation (4) is able to capture motion behaviours of root-like intruders already observed and described in the literature. Motivated by the results obtained in [20], we aim to investigate the role of the root circumnutation in reducing soil penetration resistance. In [23], the dissipated work of a parabolic robotic tip (with height h = 3.3 cm and diameter d = 2 cm) is measured. The tip is reaching a depth of x M = 30 cm in a real soil (topsoil of sandy loam type) under three different densities: k = 0.38, 0.4, 0.42 g cm −3 , obtained, respectively, without pressing the soil, with a gentle press and with a stronger press on top of the soil volume (see appendix A for the values of all variables). In the experiments produced in [23], the tip is forced to perform either a straight penetration or a circular motion (circumnutation) with a fixed slope (i.e. α introduced in section 2) with respect to the vertical direction. We refer to [23] for details about the experimental setup, the protocols and the acquisition of data.
Here, we will use the experiment and the known density coefficients in [23], to calibrate the free parameters of the dynamic equation (4), while in section 4.1 we will investigate the effects of both soil compaction and shape of the tip (section 4.1). Finally, we will derive the optimal control arising at different soil density values (section 4.2).
Simulations in this section and section 4.1 are obtained by integrating the ODE system (4) with the Euler method (step of Δx = 10 −3 ) and computing the dissipated work W with the trapezoidal rule. It remains to evaluate the control u(t) to perform the wished trajectory. At each step of the Euler method, the control u(t) is computed by solving a nonlinear equation, as follows. For the straight penetration, the control u = u(t) solves:   Indeed, the distinct combinations of circumnutating periods and tip slopes that both are in the previous range and provide the same high lead angle are close to each other. As a consequence, there are no evident benefits in the energy consumptions of the resulting tip trajectories. On the other hand, many different combinations of (T, α) in the previous range provide lower lead angles with different energy consumptions for the resulting tip trajectories. e 3 . Therefore: .
The authors in [23] conclude that the circumnutation could save up to the 33% of energy with respect to the straight penetration. Figure 3 shows the similar behaviour obtained by the two approaches (left from [23], right from our simulations) for the case of straight penetration (with α = 0 • ) and for the case of different kinds of circumnutations (α ∈ {10 • , 20 • } and T ∈ {30 s, 60 s, 120 s, 240 s}).
Moreover, it is possible to characterise the helix resulting from the circular circumnutation by the lead angle ε which is defined as: where P 0 = v 3 T is the pitch of the helix and v 3 is the downward velocity of the tip. In [23], it is shown that the circumnutation is more efficient than the straight penetration if, regardless of the soil density, the lead angle is in the range ε ∈ [45 • , 63 • ] ( figure 4(A)). To verify this condition, we simulated W at different lead angles ε in the range [20 • , 90 • ] within the three different soil densities as in [23]. The simulations (figures 4(B)-(D)) are in agreement with the experimental results in [23] and confirmed the presence of an optimal value for the lead angle such that the energy is minimised. In section 2, we have estimated the strength of the rotation by the term v 2 1 + v 2 2 /|v| (using the relation μ 1). In [23], the authors have shown that the reduction of forces due to the rotation  is directly affected by the term cos(ε). The formulation here proposed is more general, since we are comparing any possible motion in the soil and, thus, the lead angle ε is not always easy to characterise. Nevertheless, for the circular circumnutation (as in [23]), the velocity is v = (−ρω sin(ωt), ρω cos(ωt), −v 3 ) and the two terms coincide: In particular, in [23], the authors couple the cos(ε) with two varying parameters to take into account the effects of different soils and different inclinations of the tip with respect to the vertical direction. Here, we estimate this effect for any inclination of the tip by introducing the exponent in the formulation of R.

Effects of soil compaction
It has been conjectured that the circumnutation movement could rearrange the soil particles and reduce the soil penetration resistance [23]. Indeed, as it has already been noted in [24], the inclination of the root apex plays a fundamental role in facilitating soil penetration.
When the density increases up to a maximal threshold value, the energy to rearrange soil particles increases and the circumnutation may be no longer efficient. To investigate this condition, we compared the instantaneous energy between the straight penetration and the penetration performed by the root with circumnutation, in the case of a high soil density (k = 0.42) and a compacted soil (with density arbitrarily assumed in the model as k max = 0.6) (see figure 5(A)). Results confirmed the expected behaviour, showing the presence of an advantage of implementing circumnutations with a soil density of k = 0.42, which does no longer exist in the compacted soil. We additionally wanted to verify the behaviour of the energy in the presence of a medium having progressively increasing strength. Figure 5(B) shows the instantaneous energy of both the straight penetration and the circumnutation with the medium composed of different soil layers: a topsoil (up to −1.5 cm) with density progressively increasing from k = 0.38 to k = 0.4, a second layer (up to −3.5 cm) with density increasing from k = 0.4 to k = 0.42 along with depth, a third layer (up to −4.5 cm) with density going from k = 0.56 to k = 0.6 and the latter layer having a uniform k = 0.6. It is interesting to notice that, when the soil density increases, the resizing of pores becomes more difficult (if not impossible) to achieve, making the circumnutation alone a less efficient strategy in highly dense soils. These results are in agreement with biological investigations demonstrating a decrease in penetration capability of roots with the increasing of medium density [21] which requires the activitation of complementary strategies, like the promotion of radial growth [37].

Effects of tip shape
The shape and the dimension of a robotic tip can affect the efficiency of a penetrometer [17,27,38]. According to [17], a conical body that is moving into a granular medium feels a smaller drag force in the cases in which it is more tapered. In figure 6, the total work in a low dense soil (k = 0.38) is plotted when the tip becomes more and more tapered (i.e., the height-diameter ratio h/d increases). For each h/d ratio, both the straight penetration (lead angle ε = 90 • ) and the circumnutation at different periods and amplitudes are performed. The simulations show a behaviour in agreement with results in [17], where it has been observed that the drag acceleration and, consequently, the energy dissipated to perform the motion, decrease with an increment of the h/d ratio, but with less evident benefits at high ratios.

Optimal control
In this section, we will investigate the best strategy for an intruder to penetrate in a medium such that W is minimised. Following the work in [23], in the previous section 4.1 we have imposed and compared two different strategies: a straight penetration and a penetration performed with circumnutations. In this section, we extend previous results by evaluating which strategy emerges from the simulation without fixing a priori any path for the trajectory x(t) or any value for the inclination angle α of the tip with respect to the vertical axis −e 3 . We use here the direct method (see section 2.2) to compute the best strategy that minimises the dissipated work W in three soils having densities k = 0.38, 0.4, 0.42, respectively. Figure 7(A) shows that the optimal penetration strategy obtained is an oscillatory trajectory for each soil. In particular, in section 3, in the case of soil density k = 0.42, the most efficient motion has been observed to be a circumnutation with α = 10 • and T = 60 s. The optimal control performed by the direct method allows to save the 21.3% of W with respect to the previous circumnutation when the depth to reach is fixed to 5 cm, even if better results could be obtained by improving the accuracy of the direct method.
Since with the direct method also the slope α is a variable to optimise, we evaluated the behaviour for this parameter in each soil density ( figure 7(B)). It is worth noting that: • In each soil, the optimal angle α tends to continuously oscillate around an average value, without reaching a constant. • This average value decreases in more dense soils (with the increase of soil penetration resistance). This behaviour is in agreement with the experimental results in [20]. • Oscillations are wider (smaller) in the low (high) dense soil and this observation is in line with results of other investigations [24] demonstrating a strong relationship between penetration success, soil strength and the angle of inclination of roots. Indeed, the study in [24] assessed the tendency of roots to more likely penetrate in denser soils when oriented almost perpendicularly with respect to the interface from lower to higher dense soils, thus suggesting the inefficacy of performing wide oscillations in denser media. Furthermore, in [3], it has been reviewed that circumnutation is sometimes viewed as a way to reorient the tip and explore the environment. Therefore, since a common agricultural soil can be viewed as a succession of horizontal planar interfaces whose densities increase with the depth, the root tip has to manage with two conflicting demands: the need to increase its slope (to explore the soil and reorient itself) and the need to grow vertically (to improve the probability of penetrating the soil while growing downward). Figure 7(B) shows that the optimal management of the two conflicting demands can be satisfied just by looking at the minimisation of the energy dissipation.
To summarise, two main results are here obtained: (I) an oscillatory trajectory (circumnutation-like motion) as the optimal, without imposing a prescribed behaviour to the control, and (II) a dynamics of α in agreement with the experimental observations available in the literature. These results strongly validate the model and its accuracy in describing the root tip-soil interaction; while at the same time, they strengthen the hypothesis that circumnutation might be adopted by plant roots to reduce soil penetration resistance [20,23]. From these results, one can deduce that circumnutation can emerge from the mechanical interaction between the root tip and the impeding medium, as the outcome of an adapting strategy of roots to reduce the dissipated work.

Discussions
In this paper, we investigated the role of circumnutation in plant roots. Following the idea in [23], we hypothesised that circumnutation can reduce the resistance experienced by a system during penetration of a medium. We proposed a mathematical model to describe the movement of a root tip-like intruder, taking into account the soil resistance experienced during the penetration.
Previous experiments [23] showed that circumnutation could reduce the energy in vertical soil penetration. However, the complexity of the soil, the genetic diversity in plants and the lack of biological data make difficult to generalise these results to plants in real soils.
Here, we described the dynamics of a root tip-like object moving into a real soil whose particles have a diameter that is smaller than the tip size at different densities, modelling the physical phenomenon and considering only the resistance experienced at the tip. This assumption was made possible by the ability of plant roots to move in soil by growing at their apical level [26].
Firstly, we estimated model parameters using data in [23], and obtained results in agreement with the existing literature (see section 4.1). By fixing different circumnutation amplitudes, we demonstrated that below a maximal soil compaction level, the circumnutation can reduce the mechanical resistance offered by the soil with respect to a straight penetration.
Borrowing tools from the optimal control theory, we then evaluated the optimal penetration strategy for an intruder in soils which reduces the mechanical work. The resulting optimal trajectory confirmed the presence of oscillatory patterns matching circumnutation movements performed by the tip (section 4.2). Moreover, looking at the optimal slope α of the tip (the amplitude of the circumnutation motion), we observed greater amplitude for less dense soils and smaller amplitude in more dense soils, a behaviour already verified by experimental results [24].
The work showed that the circumnutation might be a behaviour modulated by the root-soil interaction, and developed by plants to save energy during growth. This suggestion is also in agreement with several biological investigations showing the influence of mechanical stimulation in the oscillatory growth patterns performed by plant root's tip [39,40]. In particular, the optimal control-based approach has successfully demonstrated that the circumnutation is the most efficient penetration strategy among an infinite combination of possible motions. Since the dimension of the system penetrating the soil can be easily changed, the model can also provide the most efficient strategy for soil exploration in autonomous systems forecasting the optimal control that minimise the mechanical work.
As far as we know, this is the first attempt to model the root dynamics by adopting an optimal control theory approach. We believe this mathematical tool could be successfully applied to investigate many other plant growth processes as well as the interactions among roots. For example, the dynamics (equation (4)) can be generalised to more than one root tip to investigate the complex plant root system development in real soil. Furthermore, if the cost function W is modified to take into account the metabolic needs of a plant, then the optimal controlbased approach could help to unveil growth mechanisms and resource allocation patterns in plants. This approach would also facilitate the analysis of a combination of strategies adopted by plants to adapt to the environment and negotiate with physical constraints. For instance, our analysis could combine circumnutation movements with the radial expansion of roots [37].
Besides the biological implications, such an approach can suggest optimal strategies to design efficient robotic devices for soil exploration. Indeed, the model can help to estimate soil resistance, forecast the requirements for the robots and identify the most convenient trajectory to follow. In particular, the framework presented in [41] can be used to generalise the results here proposed. Indeed, if the soil particles are still mobile but have a diameter that is bigger or similar to the tip size, then the differences among the resistance of each particle cannot be neglected and the soil resistance should be weighted all around the tip surface, obtaining an integrodifferential equation whose well-posedness is investigated in [41]. Furthermore, if the robotic tip moves so fast that the drag accelerationā d may not depend anymore linearly on the speed, the formulation proposed in this paper for the drag acceleration could be modified with the following: where f(|v|) is a nonlinear function to be defined. It is worthy to note that this framework is not anymore Lipschitz and the existence of a minimum solution is studied in [41]. As a future step, we mean to increase the complexity of the proposed model by adding more variables to be optimised. For example, the optimal shape could be investigated as well as the interaction among multiple tropisms (directional growth responses to environmental stimuli [39]). Moreover, as a further step, we aim to insert into the dynamics a stochastic term to simulate the random search of root tips and unveiling how the root's need to explore and exploit the soil affects the motion pattern.