Three-dimensional evolution of body and fluid motion near a wall

Evolution of three-dimensional body motion within surrounding three-dimensional fluid motion is addressed, each motion affecting the other significantly in a dynamic fluid–body interaction. This unsteady problem is set near a wall. The spatial three-dimensionality present is a new feature. For inviscid incompressible fluid, a basic nonlinear formulation is described, followed by a linearised form as a first exploration of parameter space and solution responses. The problem reduces to solving Poisson’s equation within the underbody planform, subject to mixed boundary conditions and to coupling with integral equations. Numerical and analytical properties show dependence mainly on the normal and pitch motions, as well as instability or bounded oscillations depending on the position of the centre of mass of the body, and a variety of three-dimensional shapes is examined.


Introduction
A dynamic fluid-body interaction involves the unsteady motion of a solid body that is freely moving in a surrounding fluid and affecting the fluid flow substantially, which thereby affects the body motion substantially and produces two-way interplay. In the scenarios of current concern, we have an incident unidirectional fluid flow over a solid wall but with a free finite three-dimensional body located initially above the wall such that part of the oncoming fluid travels through the gap between the underneath of the body and the wall. This yields an unsteady three-dimensional interaction.
Background motivations for focussing on these interactions include the applications to industrial, biomedical and environmental modelling. Dynamic fluid-body interactions have been much studied through direct numerical simulations and experiments. Examples include contributions by Portela et al. [18], Wang and Levy [23], Dehgahan and Basirit Tabrizi [3] concerned with turbulent flow and contributions by Hall [7], Einav and Lee [5], Petrie et al. (1993) [17], Schmidt and Young [19], Loisel et al. [13] on flow transition. Further related studies are by Owen [15], Miller et al. [14], Foucaut and Stanislas [6], Virmavirta et al. [22], Kudrolli et al. [11], Diplas et al. [4], Akoz and Kirkgoz [1], Hewitt et al. [8], Kudrolli et al. [12]. As far as we are aware, comparatively little analysis has ever been performed or presented in three spatial dimensions in the area despite the three-dimensionality of interactions in many practical configurations such as those above.
The works mentioned above are, to repeat, predominantly computational or experimental studies and much of the interest in those works as well as in the background motivations is centred on properties at fairly high Communicated by Vassilios Theofilis.  Reynolds numbers. In that parameter range relatively little exploration of fluid-body interaction has been undertaken even for two spatial dimensions in terms of applied mathematical studies apart perhaps from quite recent papers by Smith and Ellis [20], Smith and Wilson [21], Palmer and Smith [16], Jolley et al. [9] for inviscid fluids concerning bodies in boundary layers and channel flows. This sets the scene for the present study which is focussed on an applied mathematical investigation of a dynamic fluid-body interaction in three spatial dimensions.
The three-dimensional body is taken to be moving in the vicinity of the wall or ground and to be almost aligned with the wall. The body is thin as regards its representative dimension normal to the wall, relative to the typical streamwise and spanwise dimensions, and the spatial variation of the underbody shape is supposed to be of the same order as the thickness of the fluid-filled gap between the wall and the underbody. The fluid is assumed to be inviscid and incompressible and its motion to be laminar. There are connections here with the study in Jones and Smith [10] of three-dimensional flow due to a car underbody with ground effects present, but in the current cases the underbody is freely moving due to interaction with the fluid flow. The current scenario has the coordinate frame of reference to be used being such that the origin moves horizontally with the centre of mass of the body; the wall then appears as non-stationary in general. The oncoming fluid stream is assumed to be uniform. Section 2 describes the formulation of the model, including the scales and the derivation of the unsteady three-dimensional system. See also Figs. 1, 2, 3. Here, the interaction in the gap between the underbody and the wall is found to exert most influence on the combined fluid and body motions. A model simplification with a linearised approach, in Sect. 3, is then adopted as a means of gaining a first appreciation of the parameter space involved. Section 4 describes analytical features and Sect. 5 discusses the numerical method which is applied to the system, while numerical case studies are presented in Sect. 6. Section 7 provides further discussion and conclusions.

Model formulation
The typical scales ofx andỹ are comparable and can be denoted asl say, while the typical scale ofz ish where the ratioh/l 1. The body slopes are typically O(h/l). We therefore write, for the thin gap between the underbody and the wall, with X , Y , Z of order unity and anticipate that in the gap Here, the expansions in (2a) follow from a balancing in the continuity equation as well as in the streamwise and spanwise momentum equations, the pressure expansion in (2b) is also inferred from those momentum balances; the positive constant β in the time scale (2c) is taken to be large, β 1, which stems from the property that for a sufficiently dense body the representative time scale of the body motion is significantly greater than that of the fluid motion. The fluid motion in the thin gap is then governed by a nonlinear thin-layer system in 3D for an inviscid fluid, namely the continuity equation and two momentum equations in X, Y of the form The order unity scaled pressure P = P(X, Y, T ) here is independent of the normal coordinate Z by virtue of the normal momentum balance. The boundary conditions below are written now for a leading edge planform of body with shape which is given in general by X = X 1 (Y ), and F denotes the body surface in contact with the fluid: The condition (4a) corresponds to the requirement of tangential flow at the solid wall and the condition (4b) is the kinematic condition at the moving underbody surface F; we repeat that the scaled time T is slow in terms of the fluid flow. Further, (4c) represents the match with the oncoming flow in the far field and (4d) refers to the expected continuous quantities, but there are also discontinuities in velocity and pressure across the leading edge of the planform as X → X 1 − and X → X 1 +. Here, σ is the angle between the tangent to the curved leading edge and the outer stream direction (X -axis), such that cot σ = X 1 (Y ) [10]. The condition (4d) stems from a thin Euler-like region lying along the leading edge in which the horizontal velocity component normal to the edge adjusts in a quasi-planar manner along with the vertical velocity; the distance coordinate, s, along the edge being of order unity but the tangential component is convected pressure-free along the inward directed streamlines locally. Consequently, the horizontal velocity tangential to this edge is continuous across the above region. In (4d) moreover, U , V are expected to be independent of Z , although usually they depend on s. Constraint (4e) is associated with the Kutta-like requirement that the pressure be continuous across the trailing edge of the planform [10,20]. This last condition leads to the constraint that P must be zero for any trailing edge location at the boundary of the planform, the reason being that the thinness of the body induces only small variations in P except in the gap where P is of order unity. We comment also that as in the 2D case [20,21] there is upstream influence present in spite of the locally parabolic nature of the governing equations (3). The leading edge here corresponds to fluid particles which travel from the surrounding flow onto or beneath the body, whereas the trailing edge refers to particles travelling into the surrounding flow. In the present nonlinear setting, the junctions between the leading and trailing edges are unknown in advance.
The body motion is now considered. Dynamic fluid-body interaction occurs when the body is free to move, its motion being controlled by the forces and moments due to the flow pressure. In the present configuration of a thin body, only three factors are directly interactive with the motion, over the current time scales, namely the scaled lift force L Z acting in the normal direction, the scaled moment or torque M Y acting about the Y -axis and the scaled moment M X about the X -axis. The corresponding drag force L X in the X direction, for instance, is too small to alter the streamwise movements, which consist of the body continuing to translate at its initial velocity (thus, relative to the body, the oncoming stream is uniform as in (4c)), and similarly the spanwise drag force L Y in the Y direction and the moment M Z about the Z -axis are negligible. The translational velocities U = 1, V = 0 in the streamwise and spanwise directions are therefore uniform and likewise for any angular velocity dθ/dT about the Z -axis. On the other hand, the flow pressure is dependent on the body shaping as is evident from the original problem (3)(4). So the fluid and body motions are coupled.
The scaled angles of rotation about the Y, Z , X axes are written Θ 1 (T ), Θ 2 (T ), Θ 3 (T ), denoting the pitch, yaw and roll angles, respectively; the vertical distance moved by the body is H (T ) and X c , Y c denote the X -wise and Y -wise centres of mass in turn. The rates of change of X c , Y c and the yaw angle Θ 3 are taken to be zero as a prime example. We have where the double integrals extend over the XY -planform of the body. Here, m, j 1 , j 2 are, respectively, the scaled mass and moments of inertia with respect to Θ 1 , Θ 2 . Also the leading edge X 1 and trailing edge X 2 are quasi-steady, as is the fluid flow, in contrast with the motion of the body which is fully unsteady with scaled time variable T of order unity. The underbody shape is given by The largeness of the factor β mentioned after (2c) is consistent with the body-motion balances in (5), provided that the density ratio ρ B /ρ F is sufficiently large.

Model simplifications and linearised interactions
Returning to the fluid dynamical part of the model, we seek a solution in which U, V are uniform in Z inside the gap. The conditions (4a, 4b) then establish that the velocity components must satisfy while (3b, 3c) with (4c) imply that the vertical vorticity is zero, leading to the forms U = ∂Φ/∂ X, V = ∂Φ/∂Y . Here, Φ(X, Y, T ) is the scaled vertical potential function, it is governed by the elliptic equation on account of (7), with the pressure then given by The planform-edge boundary conditions (4d, 4e) written in terms of Φ thus become The fluid flow is thus controlled by (8-10) but with the unknown moving shape F(X, Y, T ) in (8) coupling the fluid flow with the body motion through (5-6). Some ideas on solving the fluid-flow part (8)(9)(10) alone are presented in Sect. 5.2 of Jones and Smith [10], especially equations (5.4)-(5.6) there. The presence of the body-motion part (5)(6) in the current overall system makes it more difficult however. So we turn to an analytical investigation next.
A linearised version of the problem described in the previous section is found to yield valuable insight. We write: For this almost-constant case F = 1 + f with the parameter being small, we can perform a linearisation such that the potential function Φ and pressure P expand as Substitution into (8) shows that at leading order we are then left with solving for the perturbation potential φ(x, y, t). The boundary conditions at the given planform edge of the body are now from (10). For prescribed underbody f (x, y, t) and planform shape the formulation (13)-(14) constitutes a closed problem for φ(x, y, t), once φ is found then the pressure p follows from Further, the rates of change of x c , y c and the yaw angle θ 3 referred to in the previous section are taken to be zero as a prime example. It is notable also that the junctions between inflow and outflow edges (leading and trailing edges) are known in advance for the current linearised setting. The underbody shape perturbation is given, from (6), by and coupled with the above we have, from (5): where, to repeat, the double integrals are over the x y-planform of the body. The governing equations and boundary conditions of the interaction are now thus (13)(14)(15)(16)(17); initial conditions on h, θ 1 , θ 2 and their temporal derivatives h , θ 1 , θ 2 are assumed known at time zero. Concerning solution properties, first we see that the roll angle θ 2 has no effect on the total interaction. This is a perhaps surprising feature. It happens because, in (16), θ 2 contributes to f only a term independent of x and has no influence on φ due to the right-hand side of (13), and hence no influence on pressure p and no influence on the double-integral forces in (17). Therefore, the body's roll angle is able to evolve simply as θ 2 (0) + θ 2 (0)t (where the constants are determined by the initial conditions), giving a uniform rolling movement of the body. This makes physical sense since the rolling takes place in a plane perpendicular to the oncoming stream. Second, it is similarly clear that the body's vertical displacement h(t) in (16) has no effect on the evolution of the pitch angle θ 1 (t) by virtue of the x-derivative on the right-hand side of (13). Hence, the linearised interaction reduces to (13-17b) in effect. Third, a numerical treatment of this linearised interaction is called for in general, following the analysis immediately below.

Analysis of special model cases
This section provides analytical features of several special model shapes, which can serve as helpful comparisons for the numerical work to be discussed in the next section, as well as shedding some extra light on the overall solution properties. The special cases considered here are (a) the circular planform, (b) the almost 2D case associated with a large span in the y direction and (c) the thin-planform case of a relatively small span in the y direction.

On the circular planform when g is negligible
This corresponds to the special case where f = θ 1 x (which is also concerned with exponential growth, see below), whose visual illustration is given in Fig. 4a. The basic problem in the real plane is then to solve the Poisson equation inside the semi-circle bounded by r = 1 and y = 0, together with the boundary conditions In this subsection, we are working for convenience in terms of x = r cos θ, y = r sin θ as the planar Cartesian coordinates and ω = x + iy is the associated complex coordinate. The condition (19c) stems from symmetry about the x-axis in the original circle case. The analytical solution via conformal mapping can be found in "Appendix 7", where the pressure along the boundary of the flow domain is obtained and shown in Fig. 4. The result is found to agree closely with the direct numerical solution of this special case.

On the nearly-2D case
Here, the parameter q is large; |y| is large compared with |x| and the y derivatives become relatively small. With leading and trailing edge positions written suitably as x = x 1 (y), x = x 2 (y), respectively, equation (13) for the special case f = θ 1 x (see Fig. 5c) gives us to leading order and so since we require φ = 0 at x = x 1 and ∂φ/∂x = 0 at x = x 2 we have Hence, the pressure solution is simply The torque τ is, therefore, the integral from Here, The tipping point is thus 1/3 of the distance from the leading edge to the trailing edge. Because of (17b), a CoM (centre of mass) lying ahead of the tipping point leads to oscillations of θ 1 and a CoM behind the tipping point leads to exponential growth.

On the very thin case
This is for q 1 in the ellipse notation where the y scale the x scale. So the governing equation (13) becomes to leading order. If the body surface function f is independent of y throughout then we obtain, for x < 0, after applying the conditions that φ = 0 at the planform edges which are given by y = Ω(x) and y = −Ω(x) say. The pressure here is therefore (26)     The body surface function is given as with the initial values of θ 1 and h taken to be −1 and 0, respectively For x > 0 application of the requirement ∂φ/∂x = 0 at the edges leads to the pressure being In the special case of f = θ 1 x, these pressures reduce to: P − = −θ 1 ΩΩ x , P + = 0. Thus, when the scaled angle θ 1 is positive the scaled pressure is negative under the front half of the body but zero under the rear half.
Additionally for the example of the ellipse where ( we have in the special case that The torque τ follows from (5) as: This implies that the critical position of the centre of mass, say x k , takes the value x k = −2 A/3, corresponding to criticality at 1/6 of the chord measured from the leading edge of the planform. If x c < x k then τ/θ 1 is negative from (6) and so the angular momentum equation (A7) yields (d 2 θ 1 /dt 2 )/θ 1 being negative, implying oscillatory behaviour. If x c > x k on the other hand then exponential growth is implied. Such behaviour can be seen in our direct numerical solutions in Fig. 6; we also demonstrate the agreement in pressure solution implied in (28) in Fig. 6b.

Numerical methods
Here, we present a numerical method of solving the Poisson equation with mixed boundary conditions in (13,14), which involves at each time step encoding their finite difference discretisations as coefficients of a sparse linear equation matrix. So long as our model is well posed and correctly discretised, this encoding matrix should be well conditioned and invertible, and its solution yields our desired flow potential over the discretisation grid. This method can be repeated and marched forward in time to obtain the body's motion. We begin by establishing the shape of the flow channel. Suppose the vertical projection of the body onto the ground is in the shape of an ellipse, which can be configured as a reasonable approximation to a wide range of model objects. The flow gap of interest is the space between this ground area and the lower surface of the body. The lower body's surface need not be smooth, instead it can have various "bumps" and "grooves" prescribed via the underbody surface function (16); various surface configurations will be demonstrated as part of our solution analysis.
Since the 3D flow problem is reduced to characteristically 2D as outlined in Sect. 3's linearisation analysis and has a domain of elliptical shape, we first translate the model equations into elliptical coordinates before applying any numerical treatments. Suppose the flow domain centred at the origin is given by x 2 /A 2 + y 2 /B 2 = 1, where 2 A and 2B are the lengths of the ellipse's two axes, or equivalently: where q = B/A; the parameter q is a constant and prescribed based on the model shape. In elliptical coordinates, we write x = r cos(θ ) and y = qr sin(θ ) for r ∈ (0, A] and θ ∈ (0, 2π]; thus, the Poisson model equation (13) becomes: (cos 2 θ + q −2 sin 2 θ) ∂ 2 φ ∂r 2 + r −2 (sin 2 θ + q −2 cos 2 θ) Now introducing the finite difference scheme, we take A to be unity without loss of generality and impose a M × N discretisation grid over this elliptical domain such that MΔr = 1, N Δθ = 2π, where r = iΔr ,   (31) can be discretised using a 4th-order central finite difference scheme: where We leave the Right-Hand-Side (R.H.S.) of the equation undiscretised for presentation simplicity: a body surface function f can be discretised using an appropriate finite difference scheme or, in many cases, its derivatives are explicitly known. Since there are in total M × N grid points and thus as many unknown fluid potential data points, this "interior discretisation" offers (M − 2) × N linear equations with On the exterior boundary, there are two boundary conditions. Firstly at a = 1 and θ ∈ [ π 2 , 3π 2 ], we have the Dirichlet condition: secondly on the rest of the curve we have the Neumann condition for θ ∈ (− π 2 , π 2 ): ∂φ ∂ x = cos θ ∂φ ∂r − sin θ ∂φ ∂θ = 0, which can be discretised using a 4-th order backward difference scheme: These two boundary conditions together offer another N linear equations for the unknowns of φ M, j for j ∈ [1, N ]. Finally, to discretise the innermost curve we find an approximation to the fluid potential at the origin φ 0, j , which itself is off the discretisation grid. To do so we apply a Taylor expansion of the fluid potential at r = 0: Note that we opted to write the fluid potential at the origin as φ 0, j as opposed to φ 0 ; this should be the same value for all angle indexers j at this point of origin (at least in theory). We do so to keep track of the discretisation angle as will become clear below. The first-order derivative ∂φ ∂r (0, j) can be discretised via a central difference scheme as: whereθ =jΔθ is the 180 • reflection of the angle θ = jΔθ over the point of origin, that isθ = |θ + π, 2π|; likewise for the second-order derivative ∂ 2 φ ∂r 2 (0, j): Substituting the two expressions into Taylor expansions of (35) gives the following discretisation of the fluid potential at the point of origin: We now can discretise the Poisson equation on the innermost curve via a central difference scheme. Substituting in the origin discretisation of (36) where necessary we have all the grid points we need to apply our scheme: This, therefore, delivers another set of N equations for φ 1, j for j ∈ [1, N ], and we thus have altogether M × N linear equations for the same number of unknowns φ i, j .
To solve this finite difference system, we encode the equations of (32, 33, 34, 36) into a (M × N )×(M × N ) This sparse matrix, which consists of M × N linearly independent equations, is of full rank and invertible: it can be solved using a linear algebra software package such as NeRSC's SuperLU. At each time step, we solve this linear system of fluid potentials φ, based on which the pressure p can be deduced via (15), and we subsequently obtain the body's positional parameters h, θ 1 and θ 2 through (17), where the double integrals are evaluated via the composite trapezoidal rule. This scheme can be repeated and marched forward in time to obtain the combined evolution of the body and fluid motions. Figure 7 demonstrates this numerical solver over a circular domain for a test case where the fluid potential is φ = r 2 (sin θ + cos θ) and r ∈ (0, 1]. In this test case, q = 1 in equation (30) and the left-hand-side of equation (31) becomes 3r 2 (sin θ + cos θ); consequently, the Dirichlet and Neumann conditions become φ = sin θ + cos θ and cos θ ∂φ ∂r − sin θ r ∂φ ∂θ = sin θ cos θ + cos 2 θ + 1, respectively. The test case shows the solution of our numerical solver agrees with the real solution closely.

Numerical case studies
In this section, we apply the aforementioned numerical scheme to bodies with various shapes, particularly by adding a "thickness" function g(x, y) to the body shape f , i.e. f (x, y, t) = g(x, y) + xθ 1 (t). We assume the body's centre of mass is at the origin throughout this section and θ 1 (t) is the angle that controls the pitch of the body. In the first case, the body resembles an inverted cone whose shape is given by f = x 2 + (y/q) 2 − 1 + xθ 1 , that is g(x, y) = x 2 + (y/q) 2 − 1. See Fig. 8b for a surface profile plot where q is taken to be unity, that is the body has the same length span in the x and y dimensions. The fluid flows in the direction of the positive x-axis and slows down as it comes in contact with the body's leading surface. This is associated with an increase in flow pressure in a region at the body's leading edge; Fig. 8a shows the pressure of the fluid and Fig. 8b depicts such pressure exerted on the body surface in the form of a heat map. Past this leading contact region the flow increases in speed as the flow channel narrows, the speed reaches a maximum at the tip of the cone where the channel is at its narrowest. Such an increase in flow speed is marked by a drop in fluid pressure which is visualised by the heat map: a negative pressure region can be seen surrounding the tip of the cone and extending to the body's trailing edge, where the atmospheric pressure condition applies. Figure 8c and d contrasts the fluid pressure for a body with narrow frontal profile facing the oncoming flow (q = 0.1) against one with a wide profile (q = 10), as can be seen by the different scales of the y-axis in the two figures. The fluid pressure for a flow past a body with a narrow frontal profile is comparatively small in magnitude and has less variation throughout the flow channel, but, on the other hand, the flow past a body with wide contact profile causes a large pressure variation. It is also interesting to note that a body with narrow profile where q = 0.1 has approximately 13% of its body surface subject to a positive pressure force, whereas for a wide profile case q = 10 less than 1% of its body surface experiences positive pressure lift. Figure 9 shows the lift and torque forces exerted on the body for a range of span parameter q values (from 0.1 to 10). The body experiences a net suction as well as negative torque in all cases. The suction force on the body grows larger as the body increases in span q, while the torque force has two distinct behaviours: for a narrow body where 0 < q < 1 it experiences an increasingly larger torque about the y-axis as q increases; for a wide body where q ≥ 1 on the other hand, the torque gradually decreases as the body's span becomes wider. Figure 10 demonstrates the evolution of the body's position over time where q is taken to be unity. Figure      angle on the other hand becomes more negatively acute over time due to the torque resulting from the pressure distribution acting on the body surface. Figure 10b shows pressure profile under the body at simulation end time (t = 6): over time the positive pressure region grows in size as well as magnitude near the leading edge, which contributes to the net positive lift on the body as well as a negative torque about the y-axis.

Case 2
In our second case, the body shape is given by f = cos(π x/4) cos(π y/(2q)) + xθ 1 and loosely resembles the chassis of an automobile or the under-carriage of a catamaran, as seen in Fig. 11b surface plot where q is taken to be unity. The interactive behaviour for the present relatively complicated shape of case 2 and for the shape in case 1 can be explained in terms similar to those for the earlier more basic shapes. The fluid flows through a (vertically narrow) channel entrance at the leading edge into a wide-open region under the mid-section of the body; in contrast to the inverted cone case, the fluid first increases in speed as it flows past this leading edge and then slows down at the wide mid channel section, this flow behaviour being marked by a low pressure region near the body's leading edge and a correspondingly high pressure region in its mid flow section, where the body's surface height is at its greatest. As the flow approaches the trailing edge it increases in speed once again and the fluid pressure gradually drops to be atmospheric at this edge. See Fig. 11a for a 3D plot of the pressure solution and Fig. 11b for the corresponding pressure heat map on the body. Figure 11c and d contrasts the fluid pressure for a body with narrow frontal profile facing the oncoming flow (q = 0.1) against one with a wide profile (q = 10). Consistent with the previous case study, the fluid pressure for a flow past a body with a narrow frontal profile is comparatively small in magnitude and has less variation throughout the flow channel, whereas the flow past a body with wide span has a comparatively greater pressure variation. We also note that a body with narrow profile where q = 0.1 has approximately 13% of its      body surface near the leading edge subject to a negative pressure force, whereas for a wide profile case q = 10 less than 1% of its leading surface region experiences a negative pressure force. Figure 12 shows the effect of the span parameter q on the lift and torque forces exerted on the body. We first note that the body experiences a net positive pressure force for all configurations of span (q ranges from 0.1 to 10), this contrasting with the case of an inverted cone where a net suction force is produced by the flow, thus reinforcing the importance of the surface shape in controlling pressure force effects. The suction force on the body grows larger as the body increases in span q, while the torque force has two distinct behaviours: for a narrow body where 0 < q < 1 it experiences an increasingly larger torque about the y-axis as q increases; for a wide body where q ≥ 1 on the other hand, the torque gradually decreases as the body's span becomes wider. Figure 13 demonstrates the evolution of the body's position over time where q is taken to be unity. Figure 13a shows changes in the body's pitch angle against the oncoming stream and its vertical centre of mass (CoM) position: its vertical CoM initially climbs higher from its starting position due to positive pressure lift under the body, but as the pitch angle evolves this pressure force drops to below zero as the region of negative pressure grows from the leading edge. As a consequence, the body's vertical CoM drops lower; its pitch angle on the other hand becomes more positive over time due to the torque from the pressure distribution. Figure 13b shows pressure profile under the body at the simulation end time (t = 6): over time the negative pressure region grows in size as well as magnitude near the leading edge, which contributes to the gradual reduction in lift on the body as well as a positive torque about the y-axis.

Further discussion and conclusions
This study which has addressed numerically and analytically the evolution of three-dimensional body motion within, and interacting with, surrounding three-dimensional fluid motion is felt to have revealed some interesting features of fluid-body interaction. The unsteady interaction problem here is set near a solid wall and  supposes the body to be solid and the fluid to be inviscid; flow separation is assumed to be absent. As regards the main findings, a dependence predominantly on the normal and pitch motions of the body has been shown in the formulation of Sects. 2 (nonlinear) and 3 (linearised), while instability or bounded oscillations are found from Sects. 4, 5 to depend on the position of the centre of mass of the body. The issue of instability or oscillations here is associated broadly with exponential growth of the normal and pitch deviations in which the underbody shape becomes negligible. Various cases of three-dimensional shapes have also been considered.
There is a useful comparison to be made between the present three-dimensional findings and those for two-dimensional interactions addressed in previous works. Foremost among these is probably the effect of moving the centre of mass of the body forward towards the leading edge. When the body is two-dimensional, or has a relatively wide span as discussed in Sect. 5, then the positioning x c of the centre of mass has a critical value at the 1/3 chord station. If x c is positioned ahead of this station then the fluid-body interaction tends to be stable, inducing oscillations, whereas an x c value greater than 1/3 gives instability corresponding to exponential growth of small disturbances. By contrast, when the body is "very" three-dimensional, having a narrow span, then the critical positioning is at the 1/6 station, with stability (instability) being present upstream (downstream) of the critical value. The stabilising effect here is physically sensible. Another notable feature is that more variations in the underbody pressure responses arise in the three-dimensional interactive results.
A number of aspects that follow on from the present first study of unsteady three-dimensional interplay between a body and the surrounding fluid are equally interesting. First, the nonlinear version of the interactive system as originally formulated in Sect. 2 would be of considerable value. The roll motion of the body represented by variations in the angle θ 2 is then expected to come into play. Next, allowance could be made for an incident stream which is angled, with constant velocity components U = U w, W = W w instead of the incoming stream in (4c). This yields a situation not dissimilar to that considered by Jones & Smith (2003) and a similar mathematical approach can be expected to hold. Third, the centre of mass positions x c , y c might vary, travelling at constant nonzero speed, and there may be a rotation of the body planform, as for a thrown discus or Frisbee, implying that the underbody shape can rotate. Finally here, other cases of dynamic fluid-body interaction such as in boundary layers and channel flows (Smith and Ellis [20], Smith and Wilson [21], Palmer and Smith [16], Jolley et al. [9]) would seem to call for spatially three-dimensional investigations. Although rather little progress has been seen previously, we believe there are many intriguing facets to explore in the future within this area of three-dimensional freely moving objects.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.

Data availability
The datasets generated during and/or analysed during the current study are available from the corresponding author on reasonable request.

Appendix: Analytical solution of a body with circular planform
This is the special case where f = θ 1 x whose visual illustration is given in Fig. 4a. The basic problem in the real plane is then to solve the Poisson equation inside the semi-circle bounded by r = 1 and y = 0, together with the boundary conditions This corresponds to the special case where f = θ 1 x (which is also concerned with exponential growth, see below), whose visual illustration is given in Fig. 4a. The basic problem in the real plane is then to solve the Poisson equation inside the semi-circle bounded by r = 1 and y = 0, together with the boundary conditions given in (19).
In this subsection, we are working for convenience in terms of x = r cos θ, y = r sin θ as the planar Cartesian coordinates and ω = x + iy is the associated complex coordinate. The condition (19c) stems from symmetry about the x-axis in the original circle case.