Rigid-flexible coupling dynamics for bionic body-caudal-fin propulsion

Abstract

With growing adoption of body-caudal-fin (BCF)-based bionic propulsion robots, optimizing thrust generation has become a critical focus of research. Despite extensive experimental efforts, theoretical studies in this area remain insufficient. This paper attempts to simplify the swing of the caudal fin as a rotation problem of the hub-beam system. Employing the Euler–Bernoulli beam theory to describe the elastic deformation of the beam and the Morison equation to describe the hydrodynamic load, the governing equation of this system is derived based on the principle of virtual work. We then investigate the effects of material properties, motion parameters, and structural parameters on thrust generation. Our model suggests that the flexibility of the beam can greatly increase instantaneous thrust. In addition, the thrust calculated by our model agrees well with the experiment results. Our work can provide guidance on choosing the proper material and rotation mode to improve propulsion performance in the design of BCF bionic robots.

Keywords

Rigid-flexible coupling Euler–Bernoulli beam theory Morison equation propulsion performance

1. Introduction

Underwater flexible bionic robots have become increasingly popular in recent years. Compared with the traditional propeller-driven autonomous underwater vehicle, body-caudal-fin (BCF) bionic propulsion robots are more efficient, mobile, bio-friendly, and concealable [1 –4]. Therefore, BCF bionic propulsion robots have numerous potential applications in maritime research, underwater creature and environment monitoring, and deep-sea resource exploration [5 –7]. However, there are still challenges in optimizing and controlling bionic flexible structures.

The caudal fin is a major source of thrust generation in fish utilizing BCF mode. Using the volume imaging system which is capable of taking three-dimensional instantaneous pictures of wake flow patterns, Lauder [8] investigated the motion patterns of fish and fins, as well as their influence on water flow patterns. Flammang et al. [9] conducted experiments and simulations on the fins of ray-finned fish, investigating how flexible the fin rays are and how flexibility influences their response to environmental disturbances. There are numerous caudal fin constructions in the current bionic structure design, including trapezoidal, rectangular, bio-inspired shapes, or some other geometries. The stiffness of the caudal fin is determined by the involved geometric shapes and the elastic modulus of the material. Experiment reveals that caudal fins with varying stiffness have distinct propulsion properties [10,11] and, within a certain range, there is a positive correlation between thrust and flexibility [12 –15]. Through numerical simulation, several studies were carried out to investigate the influence of caudal fin shape on the propulsion mechanism [16,17]. Zhong et al. [18] provided a simple tuning strategy for the caudal stiffness of the fish-like robot so that the bionic fish can double the swimming efficiency at tuna-like frequencies and speeds.

Nevertheless, an analytical model for the caudal fin swing problem is still missing. This absence significantly hinders the analysis of displacement and stress distribution within the caudal fin during motion. This motivates us to develop a rigid-flexible coupling dynamic model. This new model mimics a rotating beam and is essential for accurately describing the swing and deformation of the caudal fin. Notice that modeling the large deformations and material nonlinearities of such flexible biomimetic structures is grounded in the theoretical framework of nonlinear continuum mechanics, as formulated in finite elasticity theory [19,20].

Several models of in-plane motion of rotating flexible beams were developed by deriving the governing equation based on the Hamilton principle and the Euler–Bernoulli beam theory, with further considering the influence of bending and stretching deformation [21 –24]. The Timoshenko beam theory was also used to study the deformation of flexible beams [25]. Based on this and further considering the influence of inertial effects such as centrifugal and Coriolis forces on bending stiffness, the bending and tensile deformation of the beam were studied [26,27]. Furthermore, for beams of arbitrary shape, the nonlinear rigid-flexible coupling model was developed based on the continuum mechanics containing the finite bending, shearing, stretching, and large-scale motion [28 –30]. In order to study the rotating problem of a beam with tip mass, Zhu and Yuan established models based on the Euler–Bernoulli beam [31] and Timoshenko beam theories [32], and examined the axial deformation and deflection of the beam. As the main deformation in the rotating beam, controlling the transverse deformation is the key factor. It was found that controlling the transverse displacement of the beam was consistent with the dynamic stiffness control [33,34]. By applying the Frobenius approach, Wright et al. [35] obtained the exact frequencies and mode shapes for a rotating beam. Liu and Lu [36] developed a hub-beam system with consideration of three-dimensional motion and the interaction among large-scale overall motion, torsional deformation, and longitudinal deformation. These studies, along with advanced modeling techniques for complex structures [37,38], highlight the sophistication achievable in structural dynamics. However, the above-mentioned investigations have not taken the hydrodynamic load into consideration.

In this study, we endeavor to integrate the Morison equation [39 –42] with the rotating Euler–Bernoulli beam theory [43,44] to formulate a comprehensive rigid-flexible coupling dynamic model for the rotating beam. This model will significantly enhance our capacity for analyzing the instantaneous thrust generated by a fixed-axis rotating beam. Furthermore, it holds substantial potential for guiding the design of BCF biomimetic propulsion robots, offering valuable insights and practical applications in the field of robotic design and engineering.

The organization of the paper is as follows. In section 2, the mathematical model for the rigid beam is developed, considering the hydrodynamic loads during rotation. Section 3 presents the detailed flexible beam model where beam flexibility and its influence on the thrust generation are discussed. Section 4 details the numerical simulations and experimental validation to verify the proposed theoretical model. Section 5 investigates the thrust generation mechanism in depth, and section 6 concludes this study.

2. Rigid beam model

In this section, we propose a rigid model for the hub-beam system, as illustrated in Figure 1. A uniform beam of length L, width b, height h, density $ρ_{solid}$ , and Young’s modulus E is fixed to a rigid hub with radius R and mass moment of inertia I_h. The inertial frame, O-X₁X₂X₃, has its base vectors I , J , and K , aligned with the OX₁-, OX₂-, and OX₃-axes, respectively. The hub is connected to an actuator that applies a torque τ, causing it to rotate about the X₃-axis. Consequently, the beam motion is confined to the horizontal O-X₁X₂ plane, with no gravitational force. The floating frame, o-x₁x₂x₃, moves rigidly with the beam, and its base vectors, i , j , and k , are aligned with the ox₁-, ox₂-, and ox₃-axis, respectively. The x₁-axis is tangential to the neutral axis of the beam at the point of attachment to the hub. The relative rotation angle between the floating frame o-x₁x₂x₃ and the inertial frame O-X₁X₂X₃ is denoted by α.

Figure 1.

A BCF propulsion model based on the hub-beam system.

The position vector r of a point on the neutral axis of the beam, as expressed in the inertial frame O-X₁X₂X₃ can be written as

r = (R + x_{1}) i + x_{2} j,

(1)

where the relationship between the base vectors of the inertial frame and the floating frame satisfies

\begin{matrix} I = \cos α i - \sin α j, \\ J = \sin α i + \cos α j, \end{matrix}

(2a)

\begin{matrix} i = \cos α I + \sin α J, \\ j = - \sin α I + \cos α J . \end{matrix}

(2b)

Consequently, the velocity v and acceleration of the arbitrary point r on the beam in the inertial frame are given by

\begin{matrix} v = (R + x_{1}) \frac{d i}{dt} + x_{2} \frac{d j}{dt} \\ = - x_{2} \overset{\cdot}{α} i + (R + x_{1}) \overset{\cdot}{α} j, \end{matrix}

(3a)

and

\begin{matrix} \overset{\cdot}{v} = (- x_{2} \overset{\cdot}{α}) \frac{d i}{dt} + [(R + x_{1}) \overset{\cdot}{α}] \frac{d j}{dt} \\ = [- x_{2} \overset{\cdot\cdot}{α} - (R + x_{1}) {\overset{\cdot}{α}}^{2}] i + [- x_{2} {\overset{\cdot}{α}}^{2} + (R + x_{1}) \overset{\cdot\cdot}{α}] j, \end{matrix}

(3b)

where $\overset{\cdot}{α}$ and $\overset{\cdot\cdot}{α}$ denote the first- and second-order time derivatives of the rotation angle α, respectively.

Similarly, the virtual displacement of an arbitrary point on the beam can be written as follows:

δ r = - x_{2} δ α i + (R + x_{1}) δ α j .

(4)

The hydrodynamic forces acting on any point of the beam can be expressed using the differential form of the Morison equation [39], including both the drag force d F _D and the added mass force d F _A. These forces are given by

\begin{matrix} d F_{D} = - \frac{1}{2} ρ_{fluid} C_{D} ‖ v_{n} ‖ v_{n} d S, \\ d F_{A} = - ρ_{fluid} C_{A} \overset{\cdot}{v} d V, \end{matrix}

(5)

where $ρ_{fluid}$ is the fluid mass density per unit volume, C_D and C_A are the dimensionless drag and added mass coefficients, respectively, and $v_{n} = v_{2} n$ denotes the normal velocity at the given point on the beam. Specifically, $‖ v_{n} ‖$ represents the magnitude of this velocity, n is the unit vector normal to the mid-surface of the beam, dS and dV represent the differential elements of the mid-surface, and the volume of the beam, respectively. The cross-sectional area of the beam is A = bh.

In addition to the hydrodynamic forces, the inertial forces due to the beam motion and mass distribution must also be accounted for. These inertial forces are expressed as follows:

d F_{I} = - ρ_{solid} \overset{\cdot}{v} d V .

(6)

Therefore, the total virtual work $δ W_{Total}$ acting on the hub-beam system is given by

\begin{matrix} δ W_{Total} = δ W_{D} + δ W_{A} + δ W_{I} + δ W_{H} + δ W_{E} \\ = - \int_{S} \frac{1}{2} ρ_{fluid} C_{D} ‖ v_{n} ‖ v_{n} d S δ r - \int_{V} ρ_{fluid} C_{A} \overset{\cdot}{v} d V δ r - \int_{V} ρ_{solid} \overset{\cdot}{v} d V δ r - I_{h} \overset{\cdot\cdot}{α} δ α + τ δ α \\ = - \int_{S} \frac{1}{2} ρ_{fluid} C_{D} ‖ v_{n} ‖ v_{n} d S δ r - \int_{V} ρ_{eff} \overset{\cdot}{v} d V δ r - I_{h} \overset{\cdot\cdot}{α} δ α + τ δ α \\ = - \int_{0}^{L} \frac{1}{2} ρ_{fluid} C_{D} b | (R + x_{1}) \overset{\cdot}{α} | (R + x_{1})^{2} \overset{\cdot}{α} d x_{1} δ α \\ - ρ_{eff} b \int_{0}^{L} \int_{- \frac{h}{2}}^{\frac{h}{2}} [x_{2}^{2} \overset{\cdot\cdot}{α} + (R + x_{1}) x_{2} {\overset{\cdot}{α}}^{2} - x_{2} (R + x_{1}) {\overset{\cdot}{α}}^{2} + (R + x_{1})^{2} \overset{\cdot\cdot}{α}] d x_{2} d x_{1} δ α - I_{h} \overset{\cdot\cdot}{α} δ α + τ δ α \\ = - {\frac{1}{8} ρ_{fluid} C_{D} b [(R + L)^{4} - R^{4}] | \overset{\cdot}{α} | \overset{\cdot}{α} + ρ_{eff} IL \overset{\cdot\cdot}{α} + \frac{1}{3} ρ_{eff} A [(R + L)^{3} - R^{3}] \overset{\cdot\cdot}{α} + I_{h} \overset{\cdot\cdot}{α} - τ} δ α, \end{matrix}

(7)

where $δ W_{D}$ denotes the virtual work due to the drag force acting on the beam, $δ W_{A}$ the virtual work from the added mass force, $δ W_{I}$ the virtual work arising from the inertial forces on the beam, $δ W_{H}$ the virtual work due to the inertial forces on the hub, and $δ W_{E}$ corresponds to the virtual work induced by external forces. The combined term $ρ_{solid} + ρ_{fluid} C_{A}$ is defined as an effective density, denoted by $ρ_{eff}$ , I_h is the moment of inertia of the hub about the OX₃-axis, and I represents the cross-sectional moment of inertia around the x₃-axis on the beam, defined as follows:

I = \int_{- h / 2}^{h / 2} \int_{- b / 2}^{b / 2} x_{2}^{2} d x_{3} d x_{2} .

(8)

Using the principle of virtual work $δ W_{Total} = 0$ [45 –47], the dynamic equation governing the motion of the rigid beam is derived as follows:

{I_{h} + ρ_{eff} IL + \frac{1}{3} ρ_{eff} A [(R + L)^{3} - R^{3}]} \overset{\cdot\cdot}{α} + \frac{1}{8} ρ_{fluid} C_{D} b [(R + L)^{4} - R^{4}] | \overset{\cdot}{α} | \overset{\cdot}{α} = τ .

(9)

Furthermore, based on equations (2a), (3), (5), and (6), the instantaneous thrust force along OX₁-axis generated by the rotational motion of the beam around OX₃-axis can be calculated as follows:

\begin{matrix} F_{X} = \frac{1}{2} ρ_{fluid} C_{D} \sin α \int_{S} | v_{2} | v_{2} d S - ρ_{eff} \cos α \int_{V} {\overset{\cdot}{v}}_{1} d V + ρ_{eff} \sin α \int_{V} {\overset{\cdot}{v}}_{2} d V, \\ = \frac{1}{6} C_{D} ρ_{fluid} b [(R + L)^{3} - R^{3}] | \overset{\cdot}{α} | \overset{\cdot}{α} \sin α + \frac{1}{2} ρ_{eff} A [(R + L)^{2} - R^{2}] {\overset{\cdot}{α}}^{2} \cos α \\ + \frac{1}{2} ρ_{eff} A [(R + L)^{2} - R^{2}] \overset{\cdot\cdot}{α} \sin α . \end{matrix}

(10)

We point out again that equation (10) is the expression for the horizontal thrust force in the rigid hub-beam system when it is under the rotation. In the next section, we derive the solution for the same horizontal thrust force by assuming that the beam is flexible.

3. Flexible beam model

Based on the rigid beam model from the previous section, flexible deformations are introduced using Euler–Bernoulli beam theory, leading to the governing equations for the flexible hub-beam system. The equations are then nondimensionalized to facilitate computation and highlight the key dimensionless parameters that govern the system dynamics.

3.1. Governing equations of the flexible beam

The flexible deformations of the beam are described in the floating frame o-x₁x₂x₃. A point on the neutral axis ox₁ of the beam in its undeformed state moves to a new position under deformation. According to Euler–Bernoulli theory, for slender beams, deformations are primarily assumed to be bending, with the shear effect being neglected. In terms of this theory, the displacement field of the beam is defined as a function of the longitudinal position x₁ as follows [23]

\begin{matrix} u_{1} (x_{1}, x_{2}, t) = u (x_{1}, t) - x_{2} \frac{\partial v (x_{1}, t)}{\partial x_{1}}, \\ u_{2} (x_{1}, x_{2}, t) = v (x_{1}, t), \end{matrix}

(11)

where u and v (distinguished from the velocity vector v ), respectively, represent the axial (along x₁) and transverse (along x₂) displacements at (x₁, x₂) on the neutral axis.

Following the same approach as in the rigid beam model, we modify the expression in equation (1) using equation (11) to account for the deformation. Therefore, the position vector of a point on the deformed beam, as expressed in the inertial frame O-X₁X₂, can be written as

\begin{matrix} r = (R + x_{1} + u_{1}) i + (x_{2} + u_{2}) j \\ = (R + x_{1} + u - x_{2} v^{'}) i + (x_{2} + v) j, \end{matrix}

(12)

where the prime notation $()^{'}$ denotes partial differentiation with respect to x₁.

The velocity and acceleration at a point on the deformed beam can be expressed as the time derivatives of the position vector r [48], namely

\begin{matrix} v & = \frac{\partial}{\partial t} [(R + x_{1} + u_{1}) i + (x_{2} + u_{2}) j] \\ = (R + x_{1} + u_{1}) \frac{d i}{dt} + {\overset{\cdot}{u}}_{1} i + (x_{2} + u_{2}) \frac{d j}{dt} + {\overset{\cdot}{u}}_{2} j \\ = [{\overset{\cdot}{u}}_{1} - (x_{2} + u_{2}) \overset{\cdot}{α}] i + [(R + x_{1} + u_{1}) \overset{\cdot}{α} + {\overset{\cdot}{u}}_{2}] j \\ = [\overset{\cdot}{u} - x_{2} {\overset{\cdot}{v}}^{'} - (x_{2} + v) \overset{\cdot}{α}] i + [(R + x_{1} + u - x_{2} v^{'}) \overset{\cdot}{α} + \overset{\cdot}{v}] j, \end{matrix}

(13a)

and

\begin{matrix} \overset{\cdot}{v} & = \frac{\partial}{\partial t} {[{\overset{\cdot}{u}}_{1} - (x_{2} + u_{2}) \overset{\cdot}{α}] i + [(R + x_{1} + u_{1}) \overset{\cdot}{α} + {\overset{\cdot}{u}}_{2}] j} \\ = [{\overset{\cdot}{u}}_{1} - (x_{2} + u_{2}) \overset{\cdot}{α}] \frac{d i}{dt} + [{\overset{\cdot\cdot}{u}}_{1} - {\overset{\cdot}{u}}_{2} \overset{\cdot}{α} - (x_{2} + u_{2}) \overset{\cdot\cdot}{α}] i \\ + [(R + x_{1} + u_{1}) \overset{\cdot}{α} + {\overset{\cdot}{u}}_{2}] \frac{d j}{dt} + [{\overset{\cdot}{u}}_{1} \overset{\cdot}{α} + (R + x_{1} + u_{1}) \overset{\cdot\cdot}{α} + {\overset{\cdot\cdot}{u}}_{2}] j \\ = [{\overset{\cdot\cdot}{u}}_{1} - 2 {\overset{\cdot}{u}}_{2} \overset{\cdot}{α} - (x_{2} + u_{2}) \overset{\cdot\cdot}{α} - (R + x_{1} + u_{1}) {\overset{\cdot}{α}}^{2}] i \\ + [{\overset{\cdot\cdot}{u}}_{2} + 2 {\overset{\cdot}{u}}_{1} \overset{\cdot}{α} - (x_{2} + u_{2}) {\overset{\cdot}{α}}^{2} + (R + x_{1} + u_{1}) \overset{\cdot\cdot}{α}] j \\ = [\overset{\cdot\cdot}{u} - x_{2} {\overset{\cdot\cdot}{v}}^{'} - 2 \overset{\cdot}{v} \overset{\cdot}{α} - (x_{2} + v) \overset{\cdot\cdot}{α} - (R + x_{1} + u - x_{2} v^{'}) {\overset{\cdot}{α}}^{2}] i \\ + [\overset{\cdot\cdot}{v} + 2 (\overset{\cdot}{u} - x_{2} {\overset{\cdot}{v}}^{'}) \overset{\cdot}{α} - (x_{2} + v) {\overset{\cdot}{α}}^{2} + (R + x_{1} + u - x_{2} v^{'}) \overset{\cdot\cdot}{α}] j . \end{matrix}

(13b)

Similarly, the virtual displacement at a point on the beam can be expressed as

\begin{matrix} δ r = [δ u_{1} - (x_{2} + u_{2}) δ α] i + [(R + x_{1} + u_{1}) δ α + δ u_{2}] j \\ = [δ u - x_{2} δ v^{'} - (x_{2} + v) δ α] i + [(R + x_{1} + u - x_{2} v^{'}) δ α + δ v] j . \end{matrix}

(14)

Based on equations (5), (6), (13), and (14), the virtual works associated with the drag force, the added mass force, and the inertial forces are given by

\begin{matrix} δ W_{D} = - \int_{S} \frac{1}{2} ρ_{fluid} C_{D} ‖ v_{n} ‖ v_{n} d S δ r \\ = - \int_{0}^{L} [\frac{1}{2} ρ_{fluid} C_{D} b | (R + x_{1} + u) \overset{\cdot}{α} + \overset{\cdot}{v} | ((R + x_{1} + u) \overset{\cdot}{α} + \overset{\cdot}{v})] δ v d x_{1} \\ - \int_{0}^{L} [\frac{1}{2} ρ_{fluid} C_{D} b | (R + x_{1} + u) \overset{\cdot}{α} + \overset{\cdot}{v} | ((R + x_{1} + u) \overset{\cdot}{α} + \overset{\cdot}{v})] (R + x_{1} + u) δ α d x_{1}, \end{matrix}

(15a)

and

\begin{matrix} δ W_{A} + δ W_{I} & = - \int_{V} ρ_{eff} \overset{\cdot}{v} d V δ r \\ = - ρ_{eff} \int_{0}^{L} {[A \overset{\cdot\cdot}{u} - 2 A \overset{\cdot}{v} \overset{\cdot}{α} - Av \overset{\cdot\cdot}{α} - A (R + x_{1} + u) {\overset{\cdot}{α}}^{2}] δ u \\ + {[- I {\overset{\cdot\cdot}{v}}^{″} + {Iv}^{″} {\overset{\cdot}{α}}^{2} + A \overset{\cdot\cdot}{v} + 2 A \overset{\cdot}{u} \overset{\cdot}{α} - Av {\overset{\cdot}{α}}^{2} + A (R + x_{1} + u) \overset{\cdot\cdot}{α}]} δ v \\ + [- A \overset{\cdot\cdot}{u} v + I {\overset{\cdot\cdot}{v}}^{'} + 2 Av \overset{\cdot}{v} \overset{\cdot}{α} + I \overset{\cdot\cdot}{α} + A v^{2} \overset{\cdot\cdot}{α} + A (R + x_{1} + u) \overset{\cdot\cdot}{v} \\ + 2 A (R + x_{1} + u) \overset{\cdot}{u} \overset{\cdot}{α} + 2 {Iv}^{'} {\overset{\cdot}{v}}^{'} \overset{\cdot}{α} + A (R + x_{1} + u)^{2} \overset{\cdot\cdot}{α} + I {v^{'}}^{2} \overset{\cdot\cdot}{α}] δ α} d x_{1} . \end{matrix}

(15b)

In addition, due to the deformation of the flexible beam, the virtual work done by the internal forces must also be considered. The relationship between the infinitesimal strain $ε_{11}$ and the deformation displacement u₁ is given by

ε_{11} = u^{'} - x_{2} v^{″} .

(16)

Using the stress-strain constitutive relation, the stress $σ_{11}$ in the beam can be expressed as follows:

σ_{11} = E (u^{'} - x_{2} v^{″}) .

(17)

The strain energy associated with the beam deformation is then given by

\begin{matrix} U = \frac{1}{2} \int_{V} σ_{11} ε_{11} d V \\ = \frac{1}{2} \int_{0}^{L} \int_{- h / 2}^{h / 2} \int_{- b / 2}^{b / 2} E (u^{'} - x_{2} v^{″})^{2} d x_{3} d x_{2} d x_{1} \\ = \frac{1}{2} EA \int_{0}^{L} {u^{'}}^{2} d x_{1} + \frac{1}{2} EI \int_{0}^{L} {v^{″}}^{2} d x_{1} . \end{matrix}

(18)

Therefore, the virtual work done by the internal forces can be expressed as follows:

\begin{matrix} δ U = \int_{0}^{L} {EAu}^{'} δ u^{'} d x_{1} + \int_{0}^{L} {EIv}^{″} δ v^{″} d x_{1} \\ = - \int_{0}^{L} {EAu}^{″} δ u d x_{1} + \int_{0}^{L} {EIv}^{″ ″} δ v d x_{1} . \end{matrix}

(19)

Substituting equations (15) and (19) into the principle of minimum total potential energy $δ U - δ W_{Total} = 0$ [46,47], and taking the variations with respect to the elastic deformations u and v, and the rotation of angle α, the governing equations for the hub-beam system and the corresponding boundary conditions are derived as follows

ρ_{eff} A [\overset{\cdot\cdot}{u} - 2 \overset{\cdot}{v} \overset{\cdot}{α} - v \overset{\cdot\cdot}{α} - (R + x_{1} + u) {\overset{\cdot}{α}}^{2}] - {EAu}^{″} = 0,

(20a)

\begin{matrix} ρ_{eff} A [\overset{\cdot\cdot}{v} + 2 \overset{\cdot}{u} \overset{\cdot}{α} - v {\overset{\cdot}{α}}^{2} + (R + x_{1} + u) \overset{\cdot\cdot}{α}] + ρ_{eff} I (- {\overset{\cdot\cdot}{v}}^{″} + v^{″} {\overset{\cdot}{α}}^{2}) \\ + {EIv}^{″ ″} + \frac{1}{2} ρ_{fluid} C_{D} b | (R + x_{1} + u) \overset{\cdot}{α} + \overset{\cdot}{v} | [(R + x_{1} + u) \overset{\cdot}{α} + \overset{\cdot}{v}] = 0, \end{matrix}

(20b)

\begin{matrix} I_{h} \overset{\cdot\cdot}{α} + \int_{0}^{L} {ρ_{eff} A [- \overset{\cdot\cdot}{u} v + (R + x_{1} + u) \overset{\cdot\cdot}{v}] + ρ_{eff} A [2 v \overset{\cdot}{v} + 2 (R + x_{1} + u) \overset{\cdot}{u}] \overset{\cdot}{α} \\ + ρ_{eff} A [v^{2} + (R + x_{1} + u)^{2}] \overset{\cdot\cdot}{α} + ρ_{eff} I [\overset{\cdot\cdot}{α} + {v^{'}}^{2} \overset{\cdot\cdot}{α} + 2 v^{'} {\overset{\cdot}{v}}^{'} \overset{\cdot}{α} + {\overset{\cdot\cdot}{v}}^{'}] \\ + \frac{1}{2} ρ_{fluid} C_{D} b (R + x_{1} + u) | (R + x_{1} + u) \overset{\cdot}{α} + \overset{\cdot}{v} | [(R + x_{1} + u) \overset{\cdot}{α} + \overset{\cdot}{v}]} d x_{1} = τ, \end{matrix}

(20c)

with the corresponding boundary conditions

\begin{matrix} {u |}_{x_{1} = 0} = 0, {v |}_{x_{1} = 0} = 0, {v^{'} |}_{x_{1} = 0} = 0, \\ {u^{'} |}_{x_{1} = L} = 0, {v^{″} |}_{x_{1} = L} = 0, {v^{' ″} |}_{x_{1} = L} = 0 . \end{matrix}

(21)

Similar to the rigid beam model, and based on equations (2a), (5), (6), and (13), the instantaneous thrust generated by the rotational motion of the flexible beam along OX₁-direction is given by

\begin{matrix} F_{X} = \int_{0}^{L} {- ρ_{eff} A [\overset{\cdot\cdot}{u} - 2 \overset{\cdot}{v} \overset{\cdot}{α} - v \overset{\cdot\cdot}{α} - (R + x_{1} + u) {\overset{\cdot}{α}}^{2}] \cos α \\ + ρ_{eff} A [\overset{\cdot\cdot}{v} + 2 \overset{\cdot}{u} \overset{\cdot}{α} - v {\overset{\cdot}{α}}^{2} + (R + x_{1} + u) \overset{\cdot\cdot}{α}] \sin α \\ + \frac{1}{2} ρ_{fluid} C_{D} b | (R + x_{1} + u) \overset{\cdot}{α} + \overset{\cdot}{v} | [(R + x_{1} + u) \overset{\cdot}{α} + \overset{\cdot}{v}] \sin α} d x_{1} . \end{matrix}

(22)

Equation (22) represents the horizontal thrust in the hub-beam system when the beam is flexible beam. We point out that when Young’s modulus E approaches infinity, the deformation (u and v) tends to zero. Then, equation (22) reduces to equation (10), which describes the horizontal thrust in the rigid hub-beam system.

3.2. Nondimensional governing equations

To account for the size effects in the solid, the governing equations of the fluid-structure interaction, as presented in equations (20) and (21), need to be nondimensionalized. This process eliminates unnecessary physical parameters and introduces dimensionless groups, making the analysis more efficient and scalable. Nondimensionalization is a widely used numerical technique for solving partial differential equations, and it is applied here to obtain the system response.

The following dimensionless parameters are used for the nondimensionalization of the governing equations

\begin{matrix} \bar{t} = ω t, \\ \partial / \partial t = ω \partial / \partial \bar{t}, \partial / \partial t^{2} = ω^{2} \partial / \partial {\bar{t}}^{2}, \\ [\bar{u}, \bar{v}, {\bar{x}}_{1}, \bar{R}, \bar{L}, \bar{b}, \bar{h}] = [u, v, x_{1}, R, L, b, h] / L, \\ \frac{\partial}{\partial x_{1}} = \frac{1}{L} \frac{\partial}{\partial {\bar{x}}_{1}}, \frac{\partial}{\partial {x_{1}}^{2}} = \frac{1}{L^{2}} \frac{\partial}{\partial {\bar{x}}_{1}^{2}}, \\ \bar{A} = A / L^{2}, \bar{I} = I / L^{4}, \\ {\bar{ρ}}_{fluid} = ρ_{fluid} / ρ_{solid}, {\bar{ρ}}_{e ff} = ρ_{eff} / ρ_{solid}, \\ \bar{E} = \frac{E}{ρ_{solid} ω^{2} L^{2}}, \end{matrix}

(23)

where ω is the angular frequency of the hub actuator motion, which is closely related to the variation in α.

Using these dimensionless parameters, the governing equations in equations (20a-b) and (21) are defined as follows

\begin{matrix} {\bar{ρ}}_{e ff} \bar{A} [\frac{\partial^{2} \bar{u}}{\partial {\bar{t}}^{2}} - 2 \frac{\partial \bar{v}}{\partial \bar{t}} \frac{\partial α}{\partial \bar{t}} - \bar{v} \frac{\partial^{2} α}{\partial {\bar{t}}^{2}} - (\bar{R} + {\bar{x}}_{1} + \bar{u}) {(\frac{\partial α}{\partial \bar{t}})}^{2}] - \bar{E} \bar{A} \frac{\partial^{2} \bar{u}}{\partial {\bar{x}}^{2}} = 0, \\ {\bar{ρ}}_{e ff} \bar{A} [\frac{\partial^{2} \bar{v}}{\partial {\bar{t}}^{2}} + 2 \frac{\partial \bar{u}}{\partial \bar{t}} \frac{\partial α}{\partial \bar{t}} - \bar{v} {(\frac{\partial α}{\partial \bar{t}})}^{2} + (\bar{R} + {\bar{x}}_{1} + \bar{u}) \frac{\partial^{2} α}{\partial {\bar{t}}^{2}}] + {\bar{ρ}}_{e ff} \bar{I} [- \frac{\partial^{4} \bar{v}}{\partial {\bar{x}}_{1}^{2} \partial {\bar{t}}^{2}} + \frac{\partial^{2} \bar{v}}{\partial {\bar{x}}_{1}^{2}} {(\frac{\partial α}{\partial \bar{t}})}^{2}] \\ + \bar{E} \bar{I} \frac{\partial^{4} \bar{v}}{\partial {\bar{x}}_{1}^{4}} + \frac{1}{2} {\bar{ρ}}_{fluid} C_{D} \bar{b} | (\bar{R} + {\bar{x}}_{1} + \bar{u}) \frac{\partial α}{\partial \bar{t}} + \frac{\partial \bar{v}}{\partial \bar{t}} | [(\bar{R} + {\bar{x}}_{1} + \bar{u}) \frac{\partial α}{\partial \bar{t}} + \frac{\partial \bar{v}}{\partial \bar{t}}] = 0, \end{matrix}

(24a)

with

\begin{matrix} {\bar{u} |}_{({\bar{x}}_{1} = 0)} = 0, {\bar{v} |}_{({\bar{x}}_{1} = 0)} = 0, {\partial \bar{v} / \partial {\bar{x}}_{1} |}_{({\bar{x}}_{1} = 0)} = 0, \\ {\partial \bar{u} / \partial {\bar{x}}_{1} |}_{({\bar{x}}_{1} = 1)} = 0, {\partial^{2} \bar{v} / \partial {\bar{x}}_{1}^{2} |}_{({\bar{x}}_{1} = 1)} = 0, {\partial^{3} \bar{v} / \partial {\bar{x}}_{1}^{3} |}_{({\bar{x}}_{1} = 1)} = 0 . \end{matrix}

(24b)

Subsequently, the driving torque required for the rotational motion of the beam in equation (20c) and the resulting thrust in equation (22) are expressed as follows

\begin{array}{l} τ = I_{h} \frac{\partial^{2} α}{\partial {\bar{t}}^{2}} + ρ_{solid} L^{5} ω^{2} \int_{0}^{1} {{\bar{ρ}}_{e ff} \bar{A} [- \frac{\partial^{2} \bar{u}}{\partial {\bar{t}}^{2}} \bar{v} + (\bar{R} + {\bar{x}}_{1} + \bar{u}) \frac{\partial^{2} \bar{v}}{\partial {\bar{t}}^{2}}] + {\bar{ρ}}_{e ff} \bar{A} [2 \bar{v} \frac{\partial \bar{v}}{\partial \bar{t}} + 2 (\bar{R} + {\bar{x}}_{1} + \bar{u}) \frac{\partial \bar{u}}{\partial \bar{t}}] \frac{\partial α}{\partial \bar{t}} \\ + {\bar{ρ}}_{e ff} \bar{A} [{\bar{v}}^{2} + {(\bar{R} + {\bar{x}}_{1} + \bar{u})}^{2}] \frac{\partial^{2} α}{\partial {\bar{t}}^{2}} + {\bar{ρ}}_{eff} \bar{I} [\frac{\partial^{2} α}{\partial {\bar{t}}^{2}} + {(\frac{\partial \bar{v}}{\partial {\bar{x}}_{1}})}^{2} \frac{\partial^{2} α}{\partial {\bar{t}}^{2}} + 2 \frac{\partial \bar{v}}{\partial {\bar{x}}_{1}} \frac{\partial^{2} \bar{v}}{\partial {\bar{x}}_{1} \partial \bar{t}} \frac{\partial α}{\partial \bar{t}} + \frac{\partial^{3} \bar{v}}{\partial {\bar{x}}_{1} \partial {\bar{t}}^{2}}] \\ + \frac{1}{2} {\bar{ρ}}_{fluid} C_{D} \bar{b} (\bar{R} + {\bar{x}}_{1} + \bar{u}) | (\bar{R} + {\bar{x}}_{1} + \bar{u}) \frac{\partial α}{\partial \bar{t}} + \frac{\partial \bar{v}}{\partial \bar{t}} | [(\bar{R} + {\bar{x}}_{1} + \bar{u}) \frac{\partial α}{\partial \bar{t}} + \frac{\partial \bar{v}}{\partial \bar{t}}]} d {\bar{x}}_{1}, \end{array}

(25a)

and

\begin{matrix} F_{X} = ρ_{solid} L^{4} ω^{2} \int_{0}^{1} {- {\bar{ρ}}_{e ff} \bar{A} [\frac{\partial^{2} \bar{u}}{\partial {\bar{t}}^{2}} - 2 \frac{\partial \bar{v}}{\partial \bar{t}} \frac{\partial α}{\partial \bar{t}} - \bar{v} \frac{\partial^{2} α}{\partial {\bar{t}}^{2}} - (\bar{R} + {\bar{x}}_{1} + \bar{u}) {(\frac{\partial α}{\partial \bar{t}})}^{2}] \cos α \\ + {\bar{ρ}}_{e ff} \bar{A} [\frac{\partial^{2} \bar{v}}{\partial {\bar{t}}^{2}} + 2 \frac{\partial \bar{u}}{\partial \bar{t}} \frac{\partial α}{\partial \bar{t}} - \bar{v} {(\frac{\partial α}{\partial \bar{t}})}^{2} + (\bar{R} + {\bar{x}}_{1} + \bar{u}) \frac{\partial^{2} α}{\partial {\bar{t}}^{2}}] \sin α \\ + \frac{1}{2} {\bar{ρ}}_{fluid} C_{D} \bar{b} | (\bar{R} + {\bar{x}}_{1} + \bar{u}) \frac{\partial α}{\partial \bar{t}} + \frac{\partial \bar{v}}{\partial \bar{t}} | [(\bar{R} + {\bar{x}}_{1} + \bar{u}) \frac{\partial α}{\partial \bar{t}} + \frac{\partial \bar{v}}{\partial \bar{t}}] \sin α} d x_{1} . \end{matrix}

(25b)

In the next section, we present both experimental and numerical results of this force and discuss the new features involved.

4. Numerical simulations and experimental validation

In this section, we first obtain the hydrodynamic coefficients based on experimental measurements. These coefficients are then utilized to validate the model introduced in the previous sections by comparing the thrust generated by the rotational beam under various operating conditions, and for both rigid and flexible models.

First, the hydrodynamic loads on the thin plate were measured using a force sensor at different flow velocities. Based on the Morison equation, the hydrodynamic coefficients C_D and C_A were determined. These constant coefficients were found to be 1.5499 and 0.1284, respectively.

As shown in Figure 2, the beam is mounted on a fixed hub, with an actuator used to oscillate the beam around the hub. The angle α between the OX₁- and ox₁-axes is related to the angular amplitude B and frequency f, which can be expressed as follows:

α (t) = B \sin (2 π ft) \equiv B \sin (ω t) .

(26)

Figure 2.

Experimental setup diagram.

In our experimentation, the beam is with a length of L = 200 mm, width b = 20 mm, height h = 2 mm, and the radius of the hub R = 40 mm. The material selected is polypropylene (PP), with a density of 910 kg/m³ and a Young’s modulus E = 890 MPa. By setting the angular amplitude B and frequency f in equation (26), and substituting it into the boundary-value problem described by equation (24), the axial displacement u and deflection v at any point on the beam and at any given time can be calculated. Subsequently, by substituting u, v, and the rotational angle α into equation (25), the torque τ required to maintain the beam’s rotation and the instantaneous thrust F_X generated by the beam swing can be determined.

Figure 3 presents the instantaneous thrust curves calculated using the flexible beam model (equation (25b)) at an angular amplitude of B = π/3, and for frequencies f = 0.4 Hz, f = 0.6 Hz, and f = 0.8 Hz, along with the experimental instantaneous thrust measured by the force sensor.

Figure 3.

Comparison of theoretical and experimental thrust in the rotating flexible hub-beam for fixed rotation angular amplitude B = π/3; and frequencies (a) f = 0.4 Hz, (b) f = 0.6 Hz, (c) f = 0.8 Hz. The flexible beam is made of polypropylene (PP).

It can be observed from Figure 3 that the horizontal thrusts calculated and measured exhibit strong consistency. The calculated force agrees well with the experimentally measured one, although the latter one has slightly varying amplitudes due to the fact that the fluid used was not an ideal one. Figure 3 further shows that with increasing frequency (from Figure 3(a)–(c)), the thrust force gradually shifts its amplitude toward its positive side, providing an opportunity to enhance the directional thrust force for propulsion. As such, this comparison between theory and experiment indicates that the proposed rigid-flexible coupling dynamic hub-beam model, which is based on the Euler–Bernoulli beam theory and the Morison equation, is efficient in the simulation of bionic BCF propulsion. Therefore, in the next section, we study the propulsion mechanisms in BCF based on our theoretical solution.

5. Investigation of propulsion mechanisms

This section investigates the propulsion mechanisms of the hub-beam system by analyzing the thrust generated by both the rigid and flexible beam models introduced in previous sections. The analysis focuses on the effects of driving frequency, amplitude, and sectional moment of inertia on the thrust, providing a deeper understanding of how model characteristics, driving parameters, and structural properties influence the performance of the proposed hub-beam system.

We first study the effect of beam rigidity on the thrust. When the beam material approaches rigidity, the thrust calculated using the rigid model from section 2 (equation (10)) is compared with the thrust obtained using the flexible model from section 3 (equation (25b)), where the material is polyamide 66 (PA66) with a density of 1150 kg/m³ and a Young’s modulus of 8300 MPa. The other, even softer material is PP, the one used in Figure 3. This comparison aims to investigate how variations in material properties influence the propulsion mechanism, while keeping the driving parameters and other material characteristics constant.

As shown in Figure 4, the thrust curve for the rigid model, displayed over a dimensionless motion period, exhibits periodic symmetric fluctuations, indicating that the thrust is solely driven by the beam motion with no net thrust. In contrast, the thrust curve for the flexible beam model exhibits asymmetric fluctuations, with the average thrust significantly greater than zero, indicating that the beam deformation contributes additional net thrust. As the flexibility of the beam increases, the net thrust also increases, indicating a potential channel for creating more propulsion. Furthermore, the thrust curve shifts to the right, revealing a time-delay effect due to deformation. Again, it is this deformation which not only affects the beam dynamic response but also alters the interaction between the beam and the surrounding fluid, resulting in an asymmetric thrust distribution. This demonstrates that the deformation of flexible materials is the key factor in generating net thrust. Overall, within a certain range, flexible materials exhibit superior forward propulsion performance, with PP being the best among our two testing materials (PP and PA66).

Figure 4.

Comparison of thrust for rigid and flexible beam models at B = π/3 and f = 0.8 Hz. Note that PP is more flexible than PA66.

Second, we study the effect of the two dynamic parameters, frequency and amplitude, on the thrust, for the fixed beam material made of the soft PP. For the beam made of PP, the thrust generated at different frequencies f and amplitudes B in the flexible hub-beam system during one nondimensional vibration period is shown in Figure 5. These plots distinctly illustrate the significant dependence of thrust on both frequency and amplitude. Specifically, Figure 5(a) demonstrates that for fixed amplitude B = π/3, an increase in frequency results in a corresponding rise in net thrust. Similarly, Figure 5(b) shows that, for fixed frequency at 0.8 Hz, an increase in amplitude leads to a rise in net thrust. This clearly highlights the positive correlation between both frequency and amplitude with thrust. This behavior suggests that the deformation of the flexible beam, influenced by both frequency and amplitude, plays a key role in enhancing the thrust output.

Figure 5.

Thrust variation in the flexible hub-beam (with the soft PP beam) at different frequencies and amplitudes: (a) frequency effect for fixed B = π/3 and (b) amplitude effect for fixed f = 0.8 Hz.

Finally, we investigate the effect of the flexible beam dimensions on thrust (for the beam made of PP). The beam length is fixed at L = 200 mm, while the width b varies as 20, 30, and 40 mm. The beam height h is also adjusted to study the effect of the sectional moment of inertia on thrust generation. For ease analysis, the average thrust ${\bar{F}}_{X}$ over one motion period was used as the evaluation metric, with the driving parameters being fixed at amplitude B = π/3 and frequency f = 0.8 Hz.

Figure 6 illustrates the variation of the average thrust with respect to the beam height h for different beam width b. It can be observed that as h increases, the beam flexibility decreases, resulting in a gradual reduction in average thrust. On the other hand, thrust increases significantly with the beam width b, due to the increased flow area which can help enhance the hydrodynamic response. Notably, for different values of b, when h reaches a certain threshold, the average thrust approaches zero, indicating that the beam behaves increasingly like a rigid body. This further confirms that increasing stiffness suppresses thrust generation.

Figure 6.

Variation of average thrust (over one motion period) with beam height for different beam widths at fixed amplitude B = π/3 and frequency f = 0.8 Hz.

Results in Figure 6 indicate that thrust variation is influenced not only by the driving frequency and amplitude but also by the geometric parameters of the beam. In general, increasing the flow area and/or decreasing the beam height could effectively enhance the thrust output of the flexible beam, offering valuable insights for the application of flexible materials in propulsion systems.

6. Conclusion

In this paper, we delve into the rotational dynamics of a rigid-flexible coupling beam subjected to hydrodynamic loads. Utilizing Euler–Bernoulli beam theory, we accurately model the elastic deformation of the beam and integrate this with the rigid displacement of the hub-beam system to accurately solve the displacement on the beam at any given time. The Morison equation is applied to find the fluid load, leading to a comprehensive model based on the principle of virtual work.

The comparison between the rigid and flexible models shows that the beam deformation is essential for generating net thrust. A more flexible beam model corresponds to a higher propulsion efficiency, with thrust results further being consistent with our experimental measurements. This validation demonstrates the accuracy and efficiency in potential future applications in bionic BCF propulsion.

As numerical examples, we show that both driving frequency and amplitude can positively influence thrust generation. Namely, an increase in the driving frequency and/or amplitude will increase the thrust value. We further show that reducing the beam height will increase its flexibility, leading to a greater thrust output under consistent swing conditions. These findings offer valuable insights for optimizing beam design and for improving propulsion performance in applications such as bionic robots.

Footnotes

ORCID iDs

Ernian Pan

Yilin Qu

Funding

The authors disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: This work is supported by the National Natural Science Foundation of China (grant no. 52401391 (Y.Q.)), the Young Elite Scientists Sponsorship Program by CAST (grant no. YESS20240453 (Y.Q.)), the Ningbo Yongjiang Talent Program (grant no. 2024A-179-G (Y.Q.)), the “Innovation Yongjiang 2035” Key R&D Program (Grant No. 2024T001 (Y.Q.)), the China Postdoctoral Science Foundation (grant no. 2025T181165, No. 2023M732863 and GZC20233516 (Y.Q.)), the Open Funds of KuiYuan Laboratory (grant no. KY202436 (Y.Q.)), the Fundamental Research Funds for the Central Universities (grant no. ZDGC03 (Y.Q.)), and the Yushan Fellow Program and the National Science and Technology Council of Taiwan (grant no. NSTC 114-2221-E-A49-017 (E.P.)).

Declaration of conflicting interests

The authors declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

References

Meng

Dong

, et al. Toward a novel robotic manta with unique pectoral fins. IEEE Trans Syst Man Cybern Syst 2020; 52: 1663–1673.

Wen

Ren

Di, Santo

, et al. Understanding fish linear acceleration using an undulatory biorobotic model with soft fluidic elastomer actuated morphing median fins. Soft Robot 2018; 5: 375–388.

Tian

Wang

, et al. Robust optimization design for path planning of bionic robotic fish in the presence of ocean currents. J Mar Sci Eng 2022; 10: 1109.

Cao

Pan

, et al. Effect of cross-joints fin on the thrust performance of bionic pectoral fins. J Mar Sci Eng 2022; 10: 869.

Katzschmann

DelPreto

MacCurdy

, et al. Exploration of underwater life with an acoustically controlled soft robotic fish. Sci Robot 2018; 3: aar3449.

Liu

Jiang

Optimum curvature characteristics of body/caudal fin locomotion. J Mar Sci Eng 2021; 9: 537.

Picardi

De Luca

Chimienti

, et al. User-driven design and development of an underwater soft gripper for biological sampling and litter collection. J Mar Sci Eng 2023; 11: 771.

Lauder

GV.

Fish locomotion: recent advances and new directions. Annu Rev Mar Sci 2015; 7: 521–545.

Flammang

Alben

Madden

, et al. Functional morphology of the fin rays of teleost fishes. J Morphol 2013; 274: 1044–1059.

10.

Kopman

Porfiri

Design, modeling, and characterization of a miniature robotic fish for research and education in biomimetics and bioinspiration. IEEE/ASME Trans Mech 2012; 18: 471–483.

11.

Kancharala

Philen

MK.

Study of flexible fin and compliant joint stiffness on propulsive performance: theory and experiments. Bioinspir Biomim 2014; 9: 036011.

12.

Lucas

Thornycroft

Gemmell

, et al. Effects of non-uniform stiffness on the swimming performance of a passively-flexing, fish-like foil model. Bioinspir Biomim 2015; 10: 056019.

13.

Cui

Jiang

Design and implementation of thunniform robotic fish with variable body stiffness. Int J Robot Autom 2017; 32: 109–116.

14.

Heathcote

Wang

Gursul

Effect of spanwise flexibility on flapping wing propulsion. J Fluids Struct 2008; 24: 183–199.

15.

Quinn

Lauder

Smits

AJ.

Scaling the propulsive performance of heaving flexible panels. J Fluid Mech 2014; 738: 250–267.

16.

Song

Zhong

, et al. Tail shapes lead to different propulsive mechanisms in the body/caudal fin undulation of fish. Proc Inst Mech Eng C J Mech Eng Sci 2021; 235: 351–364.

17.

Feilich

Lauder

GV.

Passive mechanical models of fish caudal fins: effects of shape and stiffness on self-propulsion. Bioinspir Biomim 2015; 10: 036002.

18.

Zhong

Zhu

Fish

, et al. Tunable stiffness enables fast and efficient swimming in fish-like robots. Sci Robot 2021; 6: abe4088.

19.

Steigmann

DJ.

Finite elasticity theory. Oxford: Oxford University Press, 2017.

20.

Steigmann

Bîrsan

Shirani

Lecture notes on the theory of plates and shells. Cham: Springer, 2023.

21.

Choura

Jayasuriya

Medick

MA.

On the modeling, and open-loop control of a rotating thin flexible beam. J Dyn Syst Meas Contr 1991; 113: 26–33.

22.

Hanagud

Sarkar

Problem of the dynamics of a cantilevered beam attached to a moving base. J Guid Control Dynam 1989; 12: 438–441.

23.

Yang

Jiang

Chen

DC.

Dynamic modelling and control of a rotating Euler–Bernoulli beam. J Sound Vib 2004; 274: 863–875.

24.

Yigit

Scott

Ulsoy

AG.

Flexural motion of a radially rotating beam attached to a rigid body. J Sound Vib 1988; 121: 201–210.

25.

Naganathan

Soni

AH.

Coupling effects of kinematics and flexibility in manipulators. Int J Rob Res 1987; 6: 75–84.

26.

Kane

Ryan

Banerjeer

AK.

Dynamics of a cantilever beam attached to a moving base. J Guid Control Dynam 1987; 10: 139–151.

27.

Simo

Vu-Quoc

The role of non-linear theories in transient dynamic analysis of flexible structures. J Sound Vib 1987; 119: 487–508.

28.

Tang

C-M.

Nonlinear modeling and partial linearizing control of a slewing Timoshenko-beam. J Dyn Syst Meas Control 1996; 118: 75–83.

29.

Simo

Vu-Quoc

. On the dynamics in space of rods undergoing large motions—a geometrically exact approach. Comput Methods Appl Mech Eng 1988; 66: 125–161.

30.

Simo

Vu-Quoc

. On the dynamics of flexible beams under large overall motions—the plane case: part II. J Appl Mech 1986; 53: 855–863.

31.

Zhu

Mote

Jr . Dynamic modeling and optimal control of rotating Euler–Bernoulli beams. J Dyn Syst Meas Control 1997; 119: 802–808.

32.

Yuan

C-M

. Nonlinear modeling and partial linearizing control of a slewing Timoshenko-beam. J Dyn Syst Meas Control 1996; 118: 75–83.

33.

S-Q

Zou

H-L

Deng

Z-C

, et al. Dynamic analysis on hub–beam system with transient stiffness variation. Int J Mech Sci 2019; 151: 692–702.

34.

Dangarwala

Gopal

KN.

Coupled free vibration analysis of rotating non-uniform cantilever beams by an element-wise Ritz method using local hierarchical functions. Comput Struct 2023; 288: 107133.

35.

Wright

Smith

Thresher

, et al. Vibration modes of centrifugally stiffened beams. J Appl Mech 1982; 49: 197–202.

36.

Liu

J-Y

Rigid-flexible coupling dynamics of three-dimensional hub-beams system. Multibody Syst Dyn 2007; 18: 487–510.

37.

Qian

Batra

Chen

LM.

Static and dynamic deformations of thick functionally graded elastic plates by using higher-order shear and normal deformable plate theory and meshless local Petrov–Galerkin method. Compos Part B Eng 2004; 35: 685–697.

38.

Batra

RC.

Finite plane strain deformations of rubberlike materials. Int J Numer Methods Eng 1980; 15: 145–156.

39.

Morison

Johnson

Schaaf

SA.

The force exerted by surface waves on piles. J Pet Tech 1950; 2: 149–154.

40.

Morison

Johnson

O’Brien

MP.

Experimental studies of forces on piles. Coast Eng Proc 1953; 1: 25.

41.

Kopman

Laut

Acquaviva

, et al. Dynamic modeling of a robotic fish propelled by a compliant tail. IEEE J Ocean Eng 2014; 40: 209–221.

42.

Shahab

Erturk

Electrohydroelastic Euler–Bernoulli–Morison model for underwater resonant actuation of macro-fiber composite piezoelectric cantilevers. Smart Mater Struct 2016; 25: 105007.

43.

Baruh

Tadikonda

SSK

. Issues in the dynamics and control of flexible robot manipulators. J Guid Control Dynam 1989; 12: 659–671.

44.

Liang

Liu

Zhou

Quasi-variational principles of single flexible body dynamics and their applications. Sci China Ser G Phys Mech Astron 2009; 52: 775–787.

45.

Chen

Zhang

, et al. Virtual work principle for piezoelectric semiconductors and its application on extension and bending of ZnO nanowires. Crystals 2023; 13: 1368.

46.

Pan

Zhu

, et al. Modeling thermoelectric effects in piezoelectric semiconductors: new fully coupled mechanisms for mechanically manipulated heat flux and refrigeration. Int J Eng Sci 2023; 182: 103775.

47.

Pan

Zhu

, et al. A thermodynamically consistent theory for flexoelectronics: interaction between strain gradient and electric current in flexoelectric semiconductors. Int J Eng Sci 2025; 208: 104165.

48.

Xie

Zhang

, et al. A rigid-flexible coupling dynamic model for robotic manta with flexible pectoral fins. J Mar Sci Eng 2024; 12: 292.