On parametric representation of the Newton's aerodynamic problem

Newton's problem of finding the surface shape of a rotation body based on the condition of minimal resistance of the body when it moves in a rarefied medium is discussed. The problem is formulated in the form of a classical isoperimetric problem in calculus of variations. The exact solution is given in the class of piecewise differentiable functions. The numerical results of specific calculations of the functional for cone and hemisphere are presented. We prove that the optimization effect is significant by comparison of the results for cone and hemisphere with the value of the optimized functional for the optimal contour.


Introduction
The problem of finding the resistance of an axisymmetric body moving in a rarefied medium was first analyzed by Sir Isaac Newton in his "Mathematical Principles of Natural Philosophy". Newton published the solutions for specific calculations of the resistance of the axisymmetric bodies of a given shape. Also, Newton's minimal resistance problem, i.e., the problem of determining the shape of the surface of a solid of revolution, which experiences a minimum resistance when it moves through a homogeneous fluid with constant velocity in the direction of the axis of revolution, was studied in 1685 [1]. This is the first example of a problem solved in what is now called the calculus of variations, appearing a decade before the brachistochrone problem (see, for example, Ref. [2]).
Let us proceed to the mathematical formulation of the Newton's problem of a body of minimal aerodynamic resistance, or Newton's aerodynamic problem, in terms of the resistance to the movement of a body in an ideally elastic medium of particles of low density, that is , y(0) = 0, y(a) = b (1) where the integrand in Eq. (1) is always positive; then the values of functional J[y] ≥ 0 for arbitrary dependences y(x) for the specified boundary conditions at x = 0 and x = a. The profile is a surface of revolution of a curve y(x) which generates a solid body shape when it is rotated about the x-axis as it is shown by arrow in Fig. 1 in the case of a convex surface, i.e., under the additional conditions that the body is convex and rotationally symmetric.
For an explicit expression of the function y(x) satisfying the boundary conditions in Eq. (1), one can consider where n is an integer. One can note that lim n→∞ J[y n ] = 0 and, thus, the functional J[y n ] reaches its smallest value on contours y n (x) that are non-differentiable in the interval [0, a] (see, for example, Ref. [3]). Thereby, in order to regularize the solutions of variational problem in the formulation of Eq. (1), it was proposed to look for a solution in the class of monotone increasing dependences y(x), i.e., satisfying the condition y ′ (x) ≥ 0. From physics perspective, if the space is filled with a uniformly distributed rarefied fluid, it is assumed that every molecule of the fluid hits the body only once; therefore, it is necessary to introduce a condition that the function y(x), whose graph generates the rotating body shape like in Fig. 1, is monotone increasing. As a result, we arrive at the formulation of the Newton's problem in the form Eq. (2) includes unilateral constraints in the form of coupled differential equations, which make impossible to solve it by using classical algorithms for solving variational problems. A detailed analysis of the solution of problem in the formulation of Eq. (2) by transforming it to an optimal control problem and solving it by using methods for solving optimal control problems is given in Ref. [4], [page 43]. As a result of such analysis, it is shown that the solution of this variational problem is function y(x), which has one corner point at x = x * , wherein One can note various methods of investigating the Newton's aerodynamic problem considering nonsymmetric cases within the framework of the Newtonian sine-squared pressure law [5], and an extension to the Newton solution on the shape of the profile of an axisymmetric rigid body which, under prescribed length and caliber, has the least aerodynamic resistance in the presence of friction [6], as well as the existence of a body of minimal resistance with prescribed volume was proved in the class of profiles which satisfy a two sided bound on the mean curvature [7]. Also, the Newton's problem was generalized to the bodies which are not rotating, to determine a convex body of the given length and of the given maximal cross section, which has the minimal resistance when moving in the rarefied medium [8]. The Newton's aerodynamic problem of finding the local properties of convex surfaces near ridge point was solved in the case of convex bodies [9], while giving insights into the structure of singular points of a solution to Newton's least resistance problem [10]. We also note Ref. [11] in which, from a unified point of view, the well-known classical solutions of the Newton's aerodynamic problem are discussed using the Legendre transform and Hessian measures.
We analyze in this paper the variational problem in the formulation of Eq. (2), which is transformed into an equivalent isoperimetric variational problem, and, in such a way, we arrive at the parametric representation for the Newton's aerodynamic problem for the solution of which classical methods for solving variational problems can already be applied. For convenience and completeness of the presentation in Section 2, a detailed transformation of the Newton's problem to an isoperimetric variational problem is given and the necessary conditions for an extremum are formulated. Parametric representation and solutions for the Newton's aerodynamic problem, including their asymptotic expressions are analyzed in Section 3. Based on these analytical results and numerical data obtained in Section 4, we formulate the main conclusions in Section 5.

Formulation of the isoperimetric variational problem
To remove the unilateral constraints in Eq. (2), we introduce a real function u(x) according to the condition and considering the boundary condition y(0) = 0, the function y(x) can be represented as where 0 ≤ x ≤ a, and the boundary condition at x = a can be expressed accordingly by the equation Therefore, constant b is the upper limit of the function y(x) as defined by Eq. (1) and depicted in Fig. 1. As a result, the variational problem (3) can be formulated with respect to the function u(x) as an isoperimetric variational problem, that is , Note that the functionals J[u], I [u] do not contain derivatives u ′ (x), which leads to significant simplification in analysis of the isoperimetric variational problem (6).

Necessary conditions for an extremum
Let us further use the method of indefinite Lagrange multipliers to obtain the necessary conditions for an extremum in the problem under discussion, according to which we introduce the extended Lagrange functional where λ is Lagrange multiplier, and one can proceed to the mathematical formulation of the variational problem For convenience of analysis, let us represent the Lagrange functional in the form Using standard calculations in the expression for the first variation of functional (8), that is we arrive at the necessary condition for an extremum, i.e., the variation term vanishes ( Thus, the extremals of the variational problem (6) are given by function u(x) defined by the relations First of all, we note that the implementation of each of the solutions (9) on the entire interval [0, a] is not possible, which means that there are no solutions to the problem (6)  Regarding the possibility of realizing the second solution in Eq. (9) on the entire interval [0,a], we note that it is necessary that the Lagrange multiplier λ be nonzero (λ ∕ = 0) and negative. Therefore, we introduce a positive constant c > 0 according to the equality λ = − c, and rewrite the second of the identities in Eq. (9) in the form Integrating both parts of Eq. (10) in the interval [0, a] and, taking into account the integral constraint in Eq. (6), we arrive at the equality the fulfillment of which is not possible due to the divergence of the last integral in Eq. (11). As a result, we conclude that the solution of variational problem (6) is given in the class of piecewise continuous functions, while on each of the continuity intervals for u(x) the solutions are described by relations represented in a convenient form for the subsequent analysis and which are equivalent to Eq. (9), that is The structure of the Lagrange function in Eq. (7) and the necessary extremum conditions (12) allow us to draw an important conclusion about the properties of solutions u(x) for the variational problem described by Eqs. (6) and (7).

Properties of solutions and relations between parameters
Since u(x) has an even power in all relations (7)-(12), when u(x) is replaced on this interval by − u(x), then all Eqs. (7)-(12), including the value of functional in Eq. (7), do not change. This property allows us to make significant simplifications in the subsequent analysis and representation of the variational problem solutions. Taking this remark into account, we solve the variational problem for the class of non-negative functions, i.e., u(x) ≥ 0, without of any loss of generality in the analysis, and thus one can further consider Let us denote by x * the coordinate of discontinuity points of the function u(x) or one of the points in the case of the existence of several or even an infinite number of discontinuity points. Let us also introduce the notations u − = lim ε→0>0 u(x * − ε) for the left-hand and u + = lim ε→0>0 u(x * +ε) for the right-hand values of the solution u(x) at the discontinuity point (in this case u + ∕ = u − ), and let us formulate the conjugation condition for solutions at the discontinuity point. For this purpose, we use the expression for the first variation of the extended functional, due to the variation of the discontinuity point x * . We have δΛ[u] = (F(x * , u − ) − F(x * , u + )) δx * and, based on the condition that this variation vanishes, we obtain the conjugation condition for solutions at the discontinuity point F(x * ,u − ) = F(x * ,u + ), which takes the form where the definition of F(x, u) in Eq. (8) is taking into account. Using standard calculations in Eq. (14) and considering the expressions u + ∕ = u − and λ = − c, one can write Using Eq. (12) and taking into account Eq. (15) for u + ∕ = 0 and u − ∕ = 0, we arrive at two additional relations for the parameters x * , u + , u − , c, that is Considering Eq. (16) and the fulfillment of condition u ± > 0 (see Eq. (13)), after excluding parameters x * and c, we arrive at the relation In particular, if one of the parameters u + , u − vanishes, conjugation condition (15) takes the form Also, considering expressions (16) and (18), and the properties of solutions (13), we arrive at the values A.A. Barsuk and F. Paladi Thus, in the case u + ∕ = 0 and u − ∕ = 0 the solution parameters x * , u + , u − , and c satisfy Eqs. (15) and (16) at each discontinuity point of the solutions for the variational problem under discussion. When one of the parameters u + , u − vanishes, the corresponding values of parameters are given by Eq. (19) and, in accordance with Eq. (18), we arrive at an expression for the coordinate of the discontinuity point of the solution One can note that, in accordance with Eq. (20) for the cases under consideration, there is a single discontinuity point of the solution with the coordinate determined by this expression. The only admissible solutions of the variational problem are characterized by the values of solutions at discontinuity points for which u − = 0, u + ∕ = 0 and u − ∕ = 0, u + = 0 (see Eq. (19)). However, for this class of solutions, the coordinates of the discontinuity points are determined by a single expression (20), and thus the solution of the variational problem (6) has only one discontinuity point.
Let us now analyze the possibility of obtaining the solution of variational problem with a discontinuity point at which u + ∕ = u − ∕ = 0. In this case, four parameters x * , u + , u − , and c satisfy conditions defined by Eqs. (15) and (16). For convenience of the subsequent analysis, we introduce notation γ = u − /u + and, taking into account Eq. (17), one can conclude that the following parametric representations for the values of u + , u − can be further used in solving the problem Let us analyze in more detail Eqs. (15) and (16). We exclude the parameters x * and c in these equations, and arrive at the relation In addition, taking into account Eq. (21), the equation with respect to the parameter γ can be represented as 1 + 2γ 2 = γ + γ 2 + γ 3 or in its final form A similar analysis of Eq. (15) and the second of relations in Eq. (16) also leads to Eq. (22) with respect to the parameter γ. Thus, the fulfillment of Eqs. (15) and (16) is possible only if the value of parameter γ = 1, i.e., for u − = u + , which corresponds to the case of continuous solutions. Therefore, solutions of this variational problem with discontinuity points at which u − ∕ = 0, u + ∕ = 0, and u + ∕ = u − cannot be obtained.
One can also note that the interval section for which the solution u(x) ≡ 0 is adjacent to the left end of the interval [0,a]. As a result, the solution of variational problem (15) is described by the relations while u + = u(x * ) = 1 and the parameter value c is determined from the isoperimetric condition I[u] = b.

Second variation of the optimized functional
At the end of this section, let us proceed to the analysis of the extremum reached by the functional of the variational problem in the formulation of Eq. (8) on the given solutions. In this regard, we calculate the second variation of the optimized functional and represent it in the form Considering the definition of F(x, u) in Eq. (8), the expression for the Lagrange multiplier λ in Eq. (9), and the representation of the optimal solution by Eq. (23), we arrive at a representation for the second variation in the form In accordance with the above analysis u(x) ≥ 1 for x ≥ x * . For these values the expression − 1 + 3u 4 > 0, and, as a result, we conclude that δ 2 Λ[u] > 0. Thus, the optimized functional in the class of functions under consideration reaches its minimum on the obtained solution of the variational problem.

Parametric and asymptotic representations and analysis of solutions
First of all, we note that y(x) ≡ 0 in the interval 0 ≤ x ≤ x * due to the solutions u(x) ≡ 0 for the variational problem (3) and the representation (4) for y(x). In order to obtain the solution y(x) on the interval x * ≤ x ≤ a, one can use Eq. (3) in the form dy(x) = u 2 dx, and the corresponding boundary condition defined by Eq. (5). Let us further consider the representation for solution u(x) in Eq. (12) as a parametric assignment of the coordinate x = c 2 (1+τ 2 ) 2 τ (τ = u 2 > 0). Taking this into account, we obtain the expression for the differential dy = u 2 dx = τdx(τ), after the integration of which we arrive at the parametric definition of the desired contour y(x), that is where C is the integration constant. Note that the value x = x * = 2c corresponds to the value of parameter τ = τ + = 1 and, as a result, we have x ′ (τ) > 0 at τ > τ + . Indeed, x ′ (τ) = 0 at τ = τ * = 1/3 and for the interval τ > τ + > τ * we have x ′ (τ) > 0. Therefore, as the parameter τ increases, the function x(τ) monotone increases to its limiting value x = a at τ = τ 1 . On the other hand, we have y(x * ) = y(τ + ) = 0 at x = x * , from which we arrive at the value of the integration constant C = − 7 8 c and the representations which are used to determine the values of c and τ 1 .
In particular, excluding the parameter c from Eq. (25), we arrive at the equation Eq. (26) determines the dependence τ 1 (β), which can be also represented in the form For the asymptotic values of the parameter, β ≪ 1 and β ≫ 1, the solutions of Eq. (26) have the following asymptotic representations The values of optimized functional for the variational problem in the formulation of Eq. (2) are obtained in accordance with the above parametric representations, i.e., Eq. (24), and Eq. (20), that is where β = b a and the dependence τ 1 (β) is solution of Eq. (27).

Numerical results
This section presents the numerical results for dependences τ 1 (β) and y(x). First of all, in order to estimate the efficiency of the asymptotic representations (28) and for the asymptotic values of the parameter, taking into account Eq. (28), it is described by the relations Finally, Fig. 3 shows the shape of the optimal contour y(x) corresponding to the values a = 1, b = 1 (β = 1). For these values we get τ 1 = 1.916801245641185, x * = 0.3509425720484109, and c = x * /2.
The asymptotic expressions for the values of optimized functional, that is Eq. (29), in the case of optimal contours corresponding to the asymptotic values of parameter, i.e., β ≪ 1 and β ≫ 1, and the value a = 1 in Eqs. (28) and (30) take the form It follows from the second expression in Eq. (31) that as the ratio β = b a increases, the value of optimized functional decreases and tends to zero. From a practical point of view, this means that an axisymmetric body with a blunt head-form moving in the rarefied medium experiences less resistance.
We also give the value of functional in the formulation of Eq. (2) corresponding to the optimal contour y(x) for a = b = 1(β = 1), that is J[y(x)] = 0.18741. To evaluate the efficiency of the obtained optimal solution, we note that the functional value is J[y 1 (x)] = 0.25 for the dependence y 1 (x) = x, where 0 ≤ x ≤ 1, which corresponds to the generatrix of the lateral surface of a cone with the vertex placed at the origin of coordinates and the x-axis directed along the cone axis. We can also calculate the value of this functional for the contour y 2 (x) determined by the relations (x − 1) 2 + y 2 = 1, where 0 ≤ x ≤ 1 and y ≥ 0, so that in the considered interval of variation of the variable x we have y 2 (x) = , which is describing a hemisphere during rotation around the x-axis, that is J[y 2 (x)] = 0.41667. One can note that y 1 (1) = 1 and y 2 (1) = 1 for both contours y 1 (x) and y 2 (x), which corresponds to the value of β = 1. In passing, one should further note that the models for cone and convex hemisphere are important as these types of surfaces form the structural elements of the real flying objects (see, for example, Ref. [12]).

Conclusions
In this paper, we present the solution of the Newton's aerodynamic variational problem characterized by the presence of unilateral constraints in its mathematical formulation. This property does not allow the use of classical methods for solving variational problems, and, in connection with this, the methods of the optimal control theory have been recently applied to solve the Newton's aerodynamic problem. Transformation of the Newton's aerodynamic variational problem with unilateral constraints to the classical isoperimetric variational problem is given in this paper, followed by the exact solution of this problem. The obtained solution coincides with the solution given in Ref. [3] by transforming the variational problem in the formulation of Eq. (2) to the optimal control problem. Coincidence of these solutions confirms the correctness and efficiency of the method adopted in the article for solving this variational problem. The value of the optimized functional for the optimal contour is given. Using standard procedures of numerical computation, e.g., MATLAB (MathWorks) or Wolfram Mathematica, is extremely efficient in verifying numerical results, while analytical methods are desirable when it is possible to obtain exact formulas. Finally, to evaluate the efficiency of the optimal solution, the values of the functional for cone and hemisphere are also given. From the comparison of these numerical results with the optimal contour value it follows an increase of their resistance as the ratio tends to 1.33 for cone and 2.22 for hemisphere, i.e., the optimization effect is significant.

Author contribution statement
Alexandr A. Barsuk: Conceived and designed the experiments; Performed the experiments; Analyzed and interpreted the data; Wrote the paper.
Florentin Paladi: Conceived and designed the experiments; Performed the experiments; Analyzed and interpreted the data; Contributed reagents, materials, analysis tools or data; Wrote the paper.

Data availability statement
Data will be made available on request.

Declaration of competing interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.