The Mesicopter: A Miniature Rotorcraft Concept Phase II ... - CiteSeerX

2 downloads 0 Views 822KB Size Report
best 4-digit section tested at this Re, the NACA 4702. ..... The commercially available rotors specially designed for a 5-gram B-2 motor used in an electric airplane ..... wall features, the structure was judged to build up side down, and two ...
The Mesicopter: A Miniature Rotorcraft Concept Phase II Interim Report Ilan Kroo, Fritz Prinz (P.I.’s) Michael Shantz (Lecturer) Peter Kunz, Gary Fay, Shelley Cheng, Tibor Fabian, Chad Partridge (Ph.D. students) Stanford University July 2000

1. Introduction Work on the Mesicopter Phase II program began in November 1999. This report summarizes work during the first nine months of this project. The research has focused on four main areas including aerodynamic design, fabrication, power systems, and stability and control. During this period we have developed several new analysis, design, and manufacturing techniques, and have completed a number of prototypes and testbeds. In addition, we have participated in numerous discussion with NASA and DARPA groups to identify both near-term and longer-term applications of interest. Results of the research have been published in several technical papers [1-3] and popular press in several countries [4,5]. This interim report serves to document the work to date, providing view of some of the research details for NASA and NIAC reviewers and providing a mechanism for Stanford researchers to document the work to date. The report includes descriptions of rotor development (design and fabrication), power issues (batteries and motors), initial research on stability and control , and summaries of the prototypes developed to date. Future work will include a stronger emphasis on stability and control as the initial results suggest that this is more challenging than expected.

2. Rotor Development 2.1 Aerodynamic Analysis and Design Rotor aerodynamics has been a principal emphasis of the first part of this research. This is because rotor design at very small scales is poorly understood and because an efficient rotor system is necessary before issues of stability and control become critical. This section describes work on section and rotor design at very low Reynolds numbers.

2.1.1 Computational 2-D Aerodynamics The initial computational work in the area of ultra-low Reynolds aerodynamics, carried out in phase one of this program, consisted of a trade-study of the effects of airfoil geometry and Reynolds number on twodimensional airfoil performance. Continuing work in this area has focused on additional geometry considerations caused by manufacturing and structural constraints, and on the development of an optimization based design tool for ultra-low Reynolds number airfoils. Assessment of Constant Thickness Airfoils Manufacturing limitations at small scales place minimum gauge constraints on the design of highly detailed cross-sections. For very small chord-lengths, the choice must be made between a traditional airfoil profile, at the expense of greater maximum thickness, or a constant thickness distribution at the minimum gage. Figure 1 shows a photomicrograph of typical manufactured section. The as-built sections exhibit blunt leading and trailing edges and a nearly constant thickness distribution. In order to ascertain the effects on performance and the possible impact on design, alternative airfoil sections are being analyzed which incorporate the features of the as-built sections.

L.E.

T.E.

1 mm

Figure 1 Photomicrograph of mesicopter rotor section. Several NACA 4-digit airfoils have been compared with 2% constant thickness plates having matching camber distributions. Two leading edge shapes are also considered: a leading edge radius and a nearly blunt leading edge. These two geometries and the trailing edge geometry are pictured in figure 2. All computations utilized the INS-2d incompressible NavierStokes solver [6]. The primary effect of substituting a constant thickness profile is a uniform increase in drag across the operating range of the section. The drag increase ranges from 9% at Re 6000, down to 5% at Re 1000.

Blunt Leading Edge

Radiused Leading Edge

Trailing Edge

Figure 2. Edge shapes studied. Drag polars for the NACA 0002 and 2 plate sections are shown in figure 3. This drag increase is comparable to increasing from a 2% thick airfoil to a 5% thick airfoil, but the plates suffer penalties in L/D that are minimal compared to an actual increase in maximum thickness. An increase in section thickness results in a significant reduction in the lift curve slope and lower maximum steady-state lift coefficients, but the constant 2% thickness plates exhibit no reduction in lift curve slope relative to the 2% thick airfoil. Earlier separation at the leading edge does lower the maximum steady-state lift coefficient, but overall, the use of a constant thickness distribution has little impact on performance for this application.

Drag Comparison of NACA 0002 and 2% Thick Flat Plates Fully Laminar INS2d Results 0.60 NACA 0002 Blunt Plate Plate w/ Radius

0.50

Re 1000 Re 2000

0.40 Re 6000

Cl

0.30 0.20 0.10 0.00 -0.10 -0.20 0.04

0.05

0.06

0.07

0.08

0.09

0.10

0.11

0.12

0.13

Cd

Figure 3. Effect of section edge shape on performance. The principal effect of leading edge shape is on the formation of the leading edge separation bubble. At Re 1000, 2% thick flat plates with either leading edge are fully attached past 4 degrees angle of attack, and once formed, the bubbles grow slowly. The radiused plate gains less than a half degree of angle of attack for equivalent bubble lengths. Leading edge separation bubbles appear on the plates earlier than on the airfoil, but the leading edge stall on the airfoil occurs very quickly. The net effect is a minimal penalty in lift for the plates. At Re 6000, the leading edge separation bubble forms almost immediately on the blunt section, but does not form on the radiused plate until 1.5 degrees later. The blunt plate stalls slightly earlier than the radiused plate, which in turn stalls a full degree earlier than the airfoil. Avoiding a blunt profile is advisable, but as the Reynolds number and maximum section thickness are reduced, the details of the thickness distribution become less relevant and the camber-line becomes the dominant factor in performance. Airfoil Optimization The optimization study makes use of previous results to simplify the problem to its most essential elements. The section thickness is reduced to 2%, representative of the lower manufacturing limit. Although manufacturing at the lower bound results in a constant thickness profile, a NACA 4-digit thickness distribution is used. This should not affect the utility of the results and greatly facilitates automated grid generation. A specified thickness distribution reduces the number of variables considered, but also simplifies the problem by removing minimum thickness constraints. Many optimizers have a tendency to reduce the trailing edge thickness as much as possible. Implementing minimum thickness constraints so as to maintain a smooth thickness distribution is challenging, but without them optimizers often propose negative areas. These infeasible points require a restart or some other method for moving the optimization forward. The added complexity is not warranted here, since variations in the thickness distribution appear to have little effect on thin sections. The free design element is the camber-line. A wide range of variation is achieved by modeling the camberline with an Akima spline [7] anchored at the leading and trailing edges. Akima splines provide a curve that has the benefits of a tensioned spline, avoiding the undulations that can occur with simple cubic splines, but

are generally smoother across the control points and require no tuning to achieve a satisfactory interpolation. Currently, four interior control points, or knots, are used to define the camber line. These are evenly distributed and their chordwise locations are fixed. The knots are free to move perpendicular to the chord line, but the movement of the knots is constrained by upper and lower maximum camber limits. These have been implemented to prevent the optimizer from generating excessively cambered sections that might cause the flow solver to fail during a function evaluation. The optimization study utilizes a constrained simplex optimizer, a modified Nelder-Mead simplex, coupled with the INS2d code and a grid generator. Having only four design variables permits the use of this simple optimization method. In addition, each 2-D steady-state solution of the flow solver is relatively inexpensive. The simplex method is simple to implement and does not require (possibly noisy) gradient calculations. It is likely not the most efficient option, but the small problem size and inexpensive flow calculations make it a good solution. Function evaluations are also designed to increase the robustness of the method. The airfoil metric of interest, the maximum L/D, falls very close to the maximum steady-state angle of attack. The maximum steady-state angle of attack is a strong function of the geometry. If angle of attack was formulated as an independent design variable, it is likely that the optimizer would repeatedly try to evaluate an infeasible point. The flow solver would fail to converge and a restart of some sort would be necessary. The simplest way to avoid this is to remove the angle of attack as a design variable. A coarse performance polar is generated for each function evaluation. Calculations begin at a low angle of attack where convergence is likely, regardless of geometry. The angle of attack is then increased until the L/D passes its maximum or the flow solver fails to converge. In either case, an objective function value is generated for that geometry and the process continues. This is computationally expensive and not particularly efficient, but it is very robust and effective.

Re 6000 Optimized Airfoil Leading Edge

Re 2000 Optimized Airfoil Figure 4. Optimized sections at low Reynolds numbers. Two of the airfoils that have been developed using this approach are designed for Re 6000 and Re 2000. Both sections are shown in figure 4. The optimization runs were initialized with a flat plate airfoil, but the converged solutions have been checked by restarting the optimization with a seed geometry near the upper camber limits. Both airfoils exhibit similar features with a prominent droop near the nose, well-defined aft camber, and distinct hump in the camber distribution that begins near 65% chord and reaches its maximum height at 80% chord. Optimization at Re 6000 resulted in a maximum camber close to 4%, but the Re 2000 solution increases to 6% camber. The increase is a compensation for the larger reductions in effective camber at lower Reynolds number. The Re 2000 airfoil achieves a maximum L/D of 8.2, 5% higher than best 4-digit section tested at this Re, the NACA 4702. The Re 6000 airfoil achieves an L/D of 12.9, 4% better than the 4702 section and a 16% improvement over the 4402. Figure 5 shows the L/D versus geometric angle of attack for these three airfoils. The optimized section begins to show gains past three degrees, increasing until the maximum L/D is reached at 5 degrees. Past this point, the performance of this section falls off sharply.

L/D of Optimized Airfoil and NACA 4-Digit Sections Re 6000, INS2d Fully Laminar 14

12

10

L/D

8

Optimized Airfoil, Re 6000 NACA 4402 NACA 4702

6

4

2

0 -2

-1

0

1

2

3

4

5

6

7

Geo. Angle of Attack (deg.)

Figure 5. Lift-to-drag performance of low Reynolds number sections. A small portion of the drag benefit of this section is due to near-zero skin friction over the middle third of the airfoil without any substantial separation. The droop in the nose and the increase in the ideal angle of attack assist this. The majority of the gains in lift and drag are connected to a reduction in trailing edge separation. The maximum camber is moved to the aft control point and this local feature serves to contain trailing edge separation up to a certain angle of attack. This has been described in the literature as a separation ramp [8]. Up past 4 degrees angle of attack, there is a distinct inflection in the pressure distribution near 80% chord, as shown in figure 6. The reduced adverse gradients ahead of this point do not allow separation to move forward, but the steeper gradients following the inflection point promote separation up to that point. This is the primary source of the drag penalty seen with aft camber at low angles of attack. It appears that the hump in the camber distribution aids this by reducing the pressure gradients locally and improving the health of the boundary layer. Past 5 degrees, the inflection is no longer present. Beyond this point separation moves far forward and performance rapidly degrades. The principal drawback of this section is its abrupt stall behavior.

Pressure Distribution on the Aft 60% of the Re 6000 Optimized Airfoil α=4.0 deg. and α=5.0 deg., INS2d Fully Laminar -0.4 -0.3 -0.2

Cp

-0.1 Beyond α=5.0 deg., separation rapidly moves forward from 80% chord 0.0 0.1 0.2 0.3 0.4 0.4

4.0 Degrees 5.0 Degrees

0.5

0.6

0.7

0.8

0.9

1.0

x/c

Figure 6. Pressure distribution on optimized section showing separation ramp. The optimized design of these two airfoils highlights the ability of small modifications in geometry to be very effective in altering section performance. One goal of future work will be incorporating additional degrees of freedom. Two means of providing additional control over the camber distribution are to add additional knots or to provide an additional lateral degree of freedom to any or all of the existing knots. The ability to modify the thickness of the section, particularly near the leading edge where delaying separation is of interest, also holds promise for additional performance gains.

2.1.2 Rotor Design Methodology and Computational Tools A low order 3D rotor design code and a closely related rotor analysis code have been created for the development of very low Reynolds number rotors. The tools consist of a rotor performance package coupled with a nonlinear optimizer and 2-D section data from 2-D Navier Stokes analyses. The theoretical development and a preliminaryr comparison of analysis results and test data are provided. This code is currently applicable for the static thrust (hover) case only, but the development of relations for steady climbing rotors is also described. Theoretical Development: The theory applied in these codes is based on the combination of momentum flux concepts and the blade element approach. Thrust and torque expressions are developed for a differential blade element. These expressions are then integrated along the length of the blade. Two expressions are developed for thrust and two for torque. One is based on the conservation of momentum across a differential annulus; the other is based on the forces developed by the bound circulation on a differential blade element and known 2-D section properties.

The momentum theory relations for thrust and torque are based on the conservation of momentum across a differential annulus of the rotor disk. Initially, uniform velocities are assumed across the entire annulus and the flow is assumed inviscid. This is essentially actuator disk theory applied to a differential element. Corrections for these assumptions will be applied later in the formulation. Primes denote quantities per unit length. u = Induced vertical velocity v = Induced horizontal velocity ( + assumed in the direction of blade) U∞ = Freestream Vertical Velocity (U∞= 0 for hover) Thrust’ = T’ = [ 2 ρ u ( u + U∞) (2πr)] Torque’ = Q’ = [ 2 ρ v ( u + U∞) (2πr) r] Viscous effects are incorporated by adding the two-dimensional section drag to each blade element: Thrust’ = T’ = [ 2 ρ u ( u + U∞) (2πr)] - [B q c Cd Sin(φ)] Torque’ = Q’ = [ 2 ρ v ( u + U∞) (2πr) r]+ [B q c Cd r Cos(φ)] Where:

B = Number of Blades c = Local Chord φ = αi + Atan(U∞ / ωr) or

tan(φ) = (u + U∞) / ( ωr - v )

The blade element formulation is a strip theory method applied to rotor blades. At any given differential blade element, the generated forces are assumed determined from two-dimensional section properties. Given the geometry of the blade and the local inflow vector, thrust and torque can be determined. These forces are then integrated over the length of the blade and the total number of blades. In most formulations, the lift force is expressed in terms of the lift curve slope and angle of attack. In this formulation the lift is represented by the Kutta-Joukowski relation: Force per unit length = ρ (V × Γ). Thrust’ = T’ = [ B ρ ( ωr - v ) Γ ] - [B q c Cd Sin(φ)] Torque’ = Q’ = [ B ρ r ( u + U∞ ) Γ ] + [B q c Cd r Cos(φ)] The trigonometric terms are eliminated in all four relations by utilizing the following: L’ Sin(φ) = ρ( u + U∞ ) Γ L’ Cos(φ) = ρ( ωr - v ) Γ Substitution yields: Momentum: Thrust’ = T’ = [ 2 ρ u (u+ U∞) (2πr)] - [B (Cd / Cl ) ( u + U∞ ) ρ Γ ] Torque’ = Q’ = [ 2 ρ v (u+ U∞) (2πr) r ]+ [B (Cd / Cl ) ( ωr - v ) ρ Γ r ] Blade Element: Thrust’ = T’ = [ B ρ ( ωr - v ) Γ ] - [B (Cd / Cl ) ( u + U∞ ) ρ Γ ] Torque’ = Q’ = [ B ρ r ( u + U∞ ) Γ ] + [B (Cd / Cl ) ( ωr - v ) ρ Γ r ] The actuator surface approach used in the momentum equations assumes constant induced velocities across any particular annulus. The presence of a finite number of blades results in a circulation distribution that is markedly different from the infinite blade limit. A simple correction to this assumption is obtained by applying a form of the Prandtl tip loss factor [9]. This is a correction for finite blade numbers. It also has

the desired effect of driving the circulation to zero at the tip as required for a finite span wing. The φtip value represents the tip helix angle. κ = (2/π) Acos ( e -f ) f = (B/2) ( 1 - ( r / R )) ( 1 / Sin (φtip )) This correction is applied the local bound circulation: B Γ = κ Γ∞ blades Applying this to the inviscid portions of the two momentum equations yields: Momentum: Thrust’ = T’ = [ 2 ρ u (u+ U∞) (2πr) κ] - [B (Cd / Cl ) ( u + U∞ ) ρ Γ ] Torque’ = Q’ = [ 2 ρ v (u+ U∞) (2πr) r κ]+ [B (Cd / Cl ) ( ωr - v ) ρ Γ r ] Blade Element: Thrust’ = T’ = [ B ρ ( ωr - v ) Γ κ] - [B (Cd / Cl ) ( u + U∞ ) ρ Γ ] Torque’ = Q’ = [ B ρ r ( u + U∞ ) Γ ] + [B (Cd / Cl ) ( ωr - v ) ρ Γ r ] The last issue is the treatment of the tangential induced velocity term (v). The formulation to this point incorporated only the inviscid induced tangential velocity, also referred to as inviscid swirl. In most conventional large scale, high Reynolds number applications, this is sufficient. For the small scale, very low Reynolds number applications of interest here, the viscous flow entrainment becomes an important consideration. The thick wake regions generated by each blade produce a significant ‘viscous swirl’ effect. The current model used for this phenomena is described later, but for know this term is incorporated by separating the tangential velocity into vinviscid and vviscous . No viscous correction is applied to the vertical induced velocity. This viscous swirl is incorporated into all terms of the formulation except in the inviscid portion of the momentum equation for torque. In this relation the viscous losses are already accounted for in the viscous drag portion of the torque equation. These substitutions results in the final version of the four basic relations: Momentum: Thrust’ = T’ = [ 2 ρ u (u+ U∞) (2πr) κ] - [B (Cd / Cl ) ( u + U∞ ) ρ Γ ] Torque’ = Q’ = [ 2ρ vinviscid (u+ U∞) (2πr) r κ]+ [B(Cd / Cl ) ( ωr - vinviscid - vviscous ) ρ Γ r ] Blade Element: Thrust’ = T’ = [ B ρ ( ωr - vinviscid - vviscous ) Γ κ] - [B (Cd / Cl ) ( u + U∞ ) ρ Γ ] Torque’ = Q’ = [ B ρ r ( u + U∞ ) Γ ] + [B (Cd / Cl ) ( ωr - vinviscid - vviscous ) ρ Γ r ] As currently implemented, these four relations yield two equations for two unknowns (u , vinviscid ) for each differential blade element. The other three ‘unknowns’ ( Γ, κ, and vviscous ) are all dependent variables. The circulation can be expressed as a function of the known geometry, 2-D airfoil characteristics, and the local flow velocity. The tip loss factor is based on the tip helix angle. Since there is no coupling between blade elements, the tip section may be solved first, yielding a closed set. This tip velocity may then be used in solving the rest of the blade elements.

Uncoupled equations for Rotor Induced Velocities and Solution of the Design Problem: After some manipulation of the thrust and torque equations, two uncoupled equations are obtained for the u and vinviscid induced velocities. This formulation allows for a direct solution of the rotor performance problem. The required inputs are the rotor speed, chord distribution, section characteristics, and lift distribution. The incidence distribution is an output of the method along with thrust, torque, and power required. This provides a simple method for design of a rotor from a blank sheet, but analysis of an existing rotor, where the incidence is known, but the lift distribution is not, forces an iterative solution. The equations are still closed for the analysis case, but highly coupled. The determination of the Γ distribution becomes an assortment of inverse trigonometric functions of the induced velocities and the incidence angles. In the design case, Γ may be expressed as function of the input lift distribution and the local velocity vector. After a bit of work the following uncoupled equations for u and vinviscid are obtained: u = (-U∞ + Sqrt [U∞2 + 4 vinviscid (b1 - vinviscid )]) / 2 b1 = (ωr - vviscous) The equation for the tangential induced velocity is a quartic equation for arbitrary U∞ , but for the limiting case of hover (U∞ = 0) the equation reduces to a quadratic with the following solution: vinviscid = ((-b1 (b22)) +( b1 b2 Sqrt [b22 + 4])) / 2 b2 = (B / (8 π r κ)) c Cl With the induced velocities determined, the incidence angle may be calculated using the lift coefficient and lift curve slope of the airfoil. Either set of differential thrust and torque equations may be used to determine the rotor thrust, torque, and power required.

Inputs: C l Distribution Chord Distribution RPM

Calculate v for all viscous blade sections

Solve blade tip induced velocitiesPrandtl and tip loss factor

Solve induced velocities for rest of blade

Calculate blade incidence angles

Re , lC Calculate blade section profile drag

Interpolate across INS2D CFD results for the current airfoil

Cd Calculate differential thrust and torque and sum across blades

Power required = Torque * Angular Rate

Output: Thrust Torque Power Required

Figure 7. Basic structure of the rotor design program.

Viscous Swirl Model: The viscous swirl model is based on computed airfoil wake profiles at very-low Reynolds numbers. The current data is from INS2d calculations at Re 1000, 2000, and 6000. The section is a 2% constant thickness airfoil with a NACA 4402 camberline, a blunt leading edge, and a radiused trailing edge. The airfoil is operating at 4.0 degrees geometric angle of attack. Analysis of a single airfoil provides a reasonable basis for this model. The 2% constant thickness plate represents the current manufacturing minimum gage with SDM. In addition, moderate variations to the camberline would have a small effect on the wake profiles compared to the much larger variations caused by variations in Reynolds number and distance from the trailing edge. Two input parameters determine the viscous swirl correction for any individual blade element, the chord Reynolds number the blade element and the local arc length between the trailing edge of one blade and the leading edge of the next. The model applies the computed two-dimensional wake along this arc. For this first attempt at a viscous swirl correction, the average wake deficit velocity is computed instead of implementing the detailed wake profile. A model using the detailed profile, with the strongest deficit very close to the wake axis, would be very sensitive to the wake trajectory. Currently, our limited understanding of the wake behavior for these very small rotors does not justify this added complexity. One of the near term goals is to further explore the flow field generated by these rotors using both three-dimensional CFD and experimental methods. This work should greatly enhance our ability to develop more detailed and accurate models for rotor design. Figures 8 and 9 illustrate the effects of varying the distance from the trailing edge and the Reynolds number. The pronounced translation of the profiles with increasing distance occurs because the grid is aligned with the chordline and the airfoil is at a positive angle of attack. The mean deficit velocity is calculated over the region of the profile where (u / U∞) drops below 1.0. As expected, the wake diffuses and dissipates as the wakes moves down stream, decreasing the mean deficit velocity for a given Reynolds number. The effect of reducing the Reynolds number is seen to increase the width and intensity of the wake, increasing the mean deficit velocity at a given distance. The results of this analysis, combining the effects of these two parameters, are shown in figure 10. The INS2d computations are for a single airfoil in a constant free-stream flow. Modeling the viscous swirl effect as a two-element system, the forward blade and the trailing blade, neglects the fact the forward blade itself is operating in a non-uniform flow field. This is a large simplification, particularly at low advance ratios, but attempting to include the coupled effects of all the blades would be difficult with this simple model. An iterative approach to solving the system would steadily drive the total velocity to zero as the viscous terms accumulated. For this reason, the local rotational velocity is taken as the normalization velocity for vviscous and as the velocity for the local Reynolds number. This is done only for the calculation of vviscous. While not truly representative of the flow, this model at least captures the first order viscous wake effects, discouraging designs with high local solidity. Power law fitting across the distance parameter and quadratic fitting of the coefficients across Reynolds number results in the following model for vviscous: vviscous = C1 ((Arc Length)^C2) C1 = [-(3.0e-10) (Reωr 2)] + [-(2.0e-6) (Reωr )] + 0.241 C2 = [(3.0e-9) (Reωr 2)] + [-(7.0e-5) (Reωr )] - 0.372

Effect of Downstream Distance on Wake Velocity Profiles 2% t/c NACA 4402 Camberline, INS-2d, Re 1000,

α= 4.0 deg.

1.1 1.0 0.9 0.8

u / U∞

0.7 Moving downstream

0.6

from trailing edge 0.5 0.4

0.25c 0.5c 1.0c 1.6c 2.2c 10.3c

0.3 0.2 0.1 0.0 -1.0

-0.8

-0.6

-0.4

-0.2

-0.0

0.2

0.4

0.6

0.8

1.0

y/c

Figure 8. Effect of viscous wake on local velocities.

Effect of Reynold Number on Wake Velocity Profiles 2% t/c NACA 4402 Camberline, INS-2d,

α= 4.0 deg.

1.1 1.0 0.9 0.8

u / U∞

0.7 0.6 0.5

Re 1000 Re 2000 Re 6000

0.4 0.3 0.2 0.1 0.0 -0.5

-0.4

-0.3

-0.2

-0.1

-0.0

0.1

y/c

Figure 9. Effect of Reynolds number on viscous wake.

0.2

0.3

0.4

0.5

Average Wake Deficit Velocity (v

viscous )

α= 4.0 deg.

2% t/c NACA 4402 Camberline, INS-2d, 0.6

0.5

v viscous / U ∞

0.4

0.3 Re 1000 Re 2000 Re 6000

0.2

0.1

0.0 0

1

2

3

4

5

6

7

8

9

10

11

x/c aft of Trailing edge

Figure 10. Model of viscous wake for using in design code.

Rotor Design and Analysis Utilizing Nonlinear Optimization The rotor analysis method just described is only capable of estimating the performance of a given geometry and operating condition. Alone, it is incapable of developing and improving a rotor design for a particular application. The analysis method also does not incorporate the chosen power plant into the analysis. A particular power plant cannot arbitrarily provide any amount of power at any RPM. The rotor operating condition must be matched to the capabilities of the power plant to have a physically realizable system. Complete rotor analysis and design tools have developed by coupling the rotor performance program with a non-linear optimization package. The optimization code being used is the SNOPT package developed by Gill, Murray, and Saunders [10]. In the design mode, the goal is to maximize thrust for a given radius and for a particular motor. The chord distribution, lift distribution, and rotational speed are treated as the free design variables with the primary constraint that the power required matches the motor's power available. As mentioned earlier, the analysis mode is a different problem and requires the ability to input only the geometry of the rotor, chord distribution, incidence, and the rotor speed. Rather than create a second iterative rotor analysis code for this problem, this case is also solved using same rotor analysis code and non-linear optimizer, but the problem posed to the optimizer is modified. The specified incidence angle at each station is treated as an equality constraint. This drives the geometry to the specified incidence angles. The optimal, and only, solution that also satisfies the incidence constraints provides the correct lift distribution for the analysis of rotor performance. There is no constraint on the power required since in the analysis mode it is assumed that the case represents a physically realizable system. The design code produces an optimal solution at a single operating point. The analysis code allows these solutions to be evaluated over a range of operating conditions. It also provides a means of validating the method by comparing predicted performance with experimental results for existing rotors and propellers.

Results Initial results indicate that the both codes do a good job of predicting lifting performance, but the power required is significantly underestimated. Also apparent from initial results is the importance incorporating

viscous swirl effects for very small, low Reynolds number rotors. Figure 11 compares optimal chord distributions with and without the viscous swirl model. Both designs are 2.5 cm diameter, 5-bladed rotors for the 5mm Smoovy motor. Without the swirl correction the stations inboard of 35% span hit the maximum chord constraint. The remainder of the blade quickly tapers down, ostensibly reducing drag at the tip. This results in a rotor with very high solidity inboard. The original predicted thrust at 50,000 RPM was 3.9g, with the swirl correction added the predicted thrust drops 38% to 2.4g. The measured thrust for this rotor was 1.5g. Although the swirl correction does help, the very high solidity of this rotor still poses a problem for the analysis. With the swirl correction present during design the inboard solidity drops significantly and area at the tip increases. The predicted thrust for this design is 3.8 grams with the viscous swirl correction. Figure 12 compares the predicted lift versus RPM, with and without viscous swirl corrections, to experimental data for a 4-bladed, 2.5 cm diameter rotor. The prediction with the viscous swirl correction agrees well with experiment, but the prediction without the swirl correction significantly overestimates the lift for a given rotor speed. Although the lift versus RPM prediction is reasonable, the estimation of the power required appears to be underestimated by as much as a factor of two. It is unlikely that errors in the 2-D section performance predictions could account for this discrepancy. The presence of significant 3-D effects in the flow-field is most likely responsible. A key area of future work is to determine the nature of the large discrepancies between the power required predicted by the current method, and the results of experiment. In conjunction with the Rotorcraft Computational Fluid Dynamics Group at NASA Ames, 3-D Navier-Stokes calculations on several candidate rotors will be undertaken. In parallel with this computational work, experimental flow visualization studies are planned to support and complement the CFD. Between these two studies, an improved understanding of the flow field will be obtained. With this additional insight, it should be possible to improve the lower order models used for design.

Optimal Chord Distributions With and Without Viscous Swirl Correction 5-Bladed Rotor, Diam.= 2.5 cm, 5mm Smoovy 10V 0.35

0.25

0.15

c/R

0.05 0.0 -0.05

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

r/R

-0.15 with viscous swirl without viscous swirl -0.25

-0.35

Figure 11. Optimized blade chord distributions

0.9

1.0

Comparison of Predicted and Experimental Rotor Thrust 4-Bladed, 2.5cm Diameter Rotor 6

5

Thrust (g)

4

3

2 Experiment Predicted w/ Viscous Swirl Predicted w/o Viscous Swirl 1

0 25000

30000

35000

40000

45000

50000

RPM

Figure 12. Computed and measured lift vs. RPM for optimized rotors.

55000

2.2 Fabrication 2.2.1 Small Rotor Fabrication Shape Deposition Manufacturing (SDM) process was selected to fabricate rotors for mesicopter. The basic steps were described in Phase I final report and shown in the following figure (Figure 13). The efforts in the past 9 months were especially made to improve the results by experimenting on material combinations, modifying processes, and brainstorming on other possible solutions to increase rotor strength and stiffness. Moreover, a set of 5-blade rotors with various incidences were built and tested to investigate the effects of incidences of the same type of design.

Figure 13 -- steps for rotor fabrication Two types of materials are involved in SDM process—part material and support material. Epoxy (Adtech EE-501/ 503) and polyurethane (Ciba TDT 205-3) are mainly used as part materials, while wax is the standard support material. Previously, blue wax was mostly used and provided nice machinability. However, only epoxy had nice bonding to the blue wax during the machining. In order to build polyurethane rotors as well, red, yellow, and purple waxes were also tested. Red wax has the worst machinability among three. Both yellow and purple waxes can provide sufficient bonding strength during machining, and polyurethane parts can be easily pulled out of the wax substrate without melting the wax. This gives us a great potential to avoid rotor warping or distortion that may happen in melting wax. Another possible alternative is to use water-soluble support materials. Soldermask and water-soluble wax were experimented with epoxy part material. Soldermask is a UV-cured material with fair machinability, and water-soluble was is a lower temperature wax. The results showed that residual stress from the strong bonding between the water-soluble substrate and epoxy caused severe distortion when the substrate was dissolved. These experiments suggested yellow and purple waxes to be the better support materials for epoxy and polyurethane parts. Three major issues were observed from the finish parts—flash at the edges, wedge shape (see Figure 14) at the leading and trailing edges, and warping/distortion caused by thermal effect in wax removal. Different strategies were applied to resolve these problems. During the machining, over-cutting the edges on the top surface could help to clean flash. However, because of the over-cut, the cambered plate cross-section would have wedge shapes at leading and trailing edges. These wedge shapes removed sharp corners, but in the less favorable side. Therefore, another strategy, changing the building direction (from +Z to –Z), was applied. As a result, the wedge shapes were reversed to the favorable side (Figure 15), and the plate shape became closer to airfoil shape at the leading and trailing edges. Besides, precise temperature control in wax removal process (90°C) and careful handling can reduce warping and distortion.

Figure 14 -- wedge shape at L.E and T.E.

Figure 15 -- L.E. and T.E. of rotor built in the opposite building direction In some lift tests, the measured thrust degraded over several test runs. This may indicate that structural deformation was occurring during the operation. Possible strategies to increase strength and stiffness of rotor were brainstormed and summarized as follows:

Material Combination

Strategy • Place fibers in the wax mold before casting • Sputter thin metal layers (e.g. Ti) on both surfaces

Manufacturing Process



Increase the curing temperature

Structural Design



Increase the thickness as blades merge to the hub

Comments • Difficult handling and placing 5-10µm fibers. • Increase strength and stiffness effectively. • Possible discontinuity of metal layers may lose the benefit of coating and introduce extra loading on the rotor. • Higher curing temperature will strengthen the links between molecules. • The effectiveness of improvement is unknown. • Temperature is limited by the melting point of wax. • Reduce the deformation during the fabrication. • Increase the structural stiffness.

In order to verify the effect of incidence and optimize the design, a set of 5-blade rotors with various incidences were built and tested. The incidence was increased by 4, 6, and 8 degrees based on the original design. The results are summarized in the following table. Lift increases as the incidence increases. When the incidence was increased by 8 degrees, the required power was much higher than others to reach the designed RPM (~45,000 RPM), but the rotor generated more lift at this speed. RPM 46380

Lift (g) 1.325

Voltage (V) 9.1

Current (A) 0.16

48540 44820 35220 47220

1.235 1.836 1.852 2.947

9.0 8.8 8.4 11.4

0.15 0.16 0.16 0.23

Original 4 degrees + 6 degrees + 8 degrees +

2.2.2 Rotor Fabrication for Larger Testbed Vehicles The commercially available rotors specially designed for a 5-gram B-2 motor used in an electric airplane were tested and identified as proper rotors to generate sufficient lift for testbed application. These rotors

can generate about 20 grams at 20,000rpm. Reversed rotors of the same type are also needed in the testbed application. Therefore, the efforts focused on reproducing rotors that match the performance and generating reversed ones. Coordinate Measuring Machine (CMM) was used to measure the coordinates of leading and trailing edges of the sample rotor. This task was completed with the help from Precision Tool Distributors in Sunnyvale. From the measured data, distribution of chords and incidences can be determined. The NACA 4402 airfoil was first adopted and fabricated of polyurethane. However, this version was too thin to retain the blade shape during operation. Therefore, the rotor was modified resembling the shape of cambered plates with thickness of 7.5% of chord. The model was further improved by rounding sharp corners at the junctions to avoid stress concentration. With sanding at trailing edges, the reconstructed rotor can generate 20 grams as expected. By mirroring the solid model in CAD, the reversed rotors can be manufactured without difficulty.

2.2.3 Geometry Verification Two types of analysis were conducted—geometry characterization, and structural analysis. Rotor was laser scanned to provide overall picture of fabrication deviations, and blade cross-sections were cut and observed under optical microscope. These Two techniques helped to characterize features after manufacturing. Structural analysis investigates the effects of lift and pitching moment in terms of bending and torsion of a single blade. • Laser Scanning Compared to mechanical measuring and coordinate measuring (CMM), laser scanning is a better solution in our targeting dimension range (blade thickness of 50-80 µm). A company in Michigan, GKS Inspection Services, conducted measurements for us. The manufactured rotor was laser scanned with 40,000 points from the top surface and compared with the CAD model of original design. The rotor sent for laser scanning was 5-blade rotor with 4 more degrees of incidence from the original design. The compared result is shown in Figure 16. Red values represent positive side of material, meaning excess material or bent areas of the parts. Conversely, blue values represent the negative side of material. As observed in the figure, there was excess material in the leading edges and deficient in trailing edges. These substantial variations indicated that blades were twisted clockwise, resulting in more incidences than those expected and perhaps explained some of the excessive torque values measured. We further analyzed the deviations of incidence at 50% and 75% of radius for each blade based on the figure. By measuring the difference between leading and trailing edges and taking chord values from the model, the twist angles (in degrees) were calculated and shown in the below table. Blade 1 50%R 75%R Twist Angle (degrees)

1.41

8.62

Blade 2 50%R 75%R 10.27

2.28

Blade 3 50%R 75%R 3.05

2.07

Blade 4 50%R 75%R 4.03

11.38

Blade 5 50%R 75%R 2.72

0.98

In order to reduce the twist introduced during manufacturing, the blade thickness was increased at 25% of radius and inboard. The modified rotor was laser scanned (Figure 17) for comparison. The anomalous twist has been substantially reduced.

Figure 16 – laser scanning result of 4 degrees+ rotor

Figure 17 – laser scanning result of modified rotor



Cross-section Examination To characterize cross-section geometry after fabrication, the rotor was embedded in another material and cut at desired positions (Figure 18). Pictures taken under microscope provided cross-section geometry information of the fabricated rotors. In particular, shape of leading and trailing edges are verified. This also gave us great feedback on process improvement regarding building directions (Figure 19).

Figure 18 – cutting cross-section at desired positions

Figure 19 – cross-section pictures taken under optical microscope for rotors with different building directions



Bending Analysis The blade was analyzed as a cantilever beam subject to various distributed loads (lift) along the radius direction (Figure 20). Both finite element and analytical solutions were calculated and compared. For 5blade rotor, the beam was approximated by 19 elements, 20 stations. The moment of inertia, I, at each station was calculated in the CAD software, and was assumed to be liner in each element. The values of distributed loads were taken directly from the design sheet, and the force density in each element was assumed to be constant. The deflections of epoxy rotor obtained from finite element method (FE) and analytical calculation are compared in Figure 21. Moreover, the deflections (from finite element analysis) of epoxy and polyurethane rotors are compared in Figure 22.

f

Figure 20 -- rotor blade modeled as a cantilever beam

analytical-ep

Epoxy Rotor

fem-ep-result

0.009

deflection (mm)

0.008 0.007 0.006 0.005 0.004 0.003 0.002 0.001 0.000

0.00

2.00

4.00

6.00

8.00

10.00

12.00

r (mm)

Figure 21 -- 5-blade FE and analytical solution

epoxy-analytical epoxy-fem pu-analytical pu-fem

Deflection (mm)

Bending 0.05 0.04 0.03 0.02 0.01 0 0.00

2.00

4.00 6.00 r (mm)

8.00

10.00

Figure 22 -- 5-blade bending for epoxy and PU rotors The deflections associated with bending are very small for epoxy rotors. However, with the same loading condition and rotor geometry, deflections of PU rotors are about 4-5 times of those of epoxy ones. This may explain why the PU rotors performed a significant lift reduction during spinning. For big-hub 4-blade epoxy rotors, similar analyses were conducted with 18 elements and 19 stations. The deflections calculated by FE were compared with 5-blade rotor (Figure 23). The 4-blade rotor is less stiff because of smaller chords and larger diameter. The deflection at tip of 4-blade rotor is about 1.6 times of that of 5-blade rotor.

4-blade

Deflection due to bending

5-blade

deflection (mm)

0.015 0.0125 0.01 0.0075 0.005 0.0025 0 0.0

2.0

4.0

6.0

r (mm)

8.0

10.0

12.0

Figure 23 -- bending for epoxy 4- and 5- blade rotors •

Torsion Analysis The torsion analysis here only considered the effect from pitching moment and lift. The pitching moment coefficients, Cm, at different cross-sections were given. The moments calculated from Cm and lift act at the quarter chord of the cross-section. Therefore, in order to find the total moment at the center of each cross-section, the moment caused by lift should be also considered. Figure 24 shows the torsion results of 5-blade epoxy and PU rotors, while Figure 25 is the comparison between 4-blade and 5-blade epoxy rotors. At the tip of 5-blade epoxy rotor, the torsion angle (.27°) is about 3% of the incidence (8.9°). On the other hand, at the tip of 4-blade rotor, the torsion angle (0.46°) is about 5.8% of the incidence (7.9°).

epoxy

Torsion

polyurethane

angle (degrees)

1.5 1.25 1 0.75 0.5 0.25 0 -0.25 0.0

2.0

4.0

6.0

r (mm)

8.0

10.0

Figure 24 -- 5-blade torsion for epoxy and PU rotors

4-blade

angle (degrees)

Torsion

5-blade

0.50 0.40 0.30 0.20 0.10 0.00 -0.10 0.00

2.00

4.00

6.00

r (mm) Figure 25 -- torsion for epoxy 4- and 5- blade rotors

8.00

10.00

12.00

2.2.4 Airframe Fabrication The airframe is a structural component of the mesicopter integrating all subsystems. It includes the rotor shroud as well as structural connections. In the process of mesicopter assembly batteries and motors with controllers are connected to the airframe. Additionally the airframe can accommodate shrouds around propellers to first protect the fragile propellers and second to improve the lift generated by the propellers. Besides the above-mentioned functionality the airframe has to satisfy size, weight and manufacturability constraints Figure 26 shows the evolution of the airframe design. The manufacturability of every particular design was the main driving force. The first three models, see 26a-c, are based on a printed circuit board (PCB) as a structural material for the airframe. This choice allows elegant integration of electronics circuitry on the airframe. On the other side the PCB is limited to planar structures and therefore inappropriate for the shroud. Figure 26a shows the initial design - the first iteration. The second design, see 26b, included modifications based on stability analysis. To move the center of gravity and to reduce the moment of inertia the batteries were moved inboard also the rotor spacing was reduced. This design showed following drawbacks: still a lot of unused frame area, overall size too large and components in the airflow of the propellers. The third design, see 26c, addressed these issues, and included the computed rotor tilt angles to maximum damping of the system dynamics. The last three models shown in Figure 26d-f implement the shrouds. These models cannot be built using PCB and therefore require different manufacturing approach. The model shown in Figure 26d consists of three parts (shrouds, connecting pins, and motor holders) that are later assembled. The model shown in Figure e consists of a single-part airframe. The structure is therefore stiffer and the alignment between shrouds and motor holders is more precise. The final design is shown in Figure 26f. The overall size of this design is 65.6 mm. The calculated weight of the airframe is 3.24 g with polyurethane as the airframe material. •

Material Investigation

a

b

c

d

e

f

Figure 26: Evolution of the airframe design.

Printed circuit board was considered as an airframe material in the early design. The advantage is to integrate motor controllers and sensor components with the frame structure. However, as we included 3D shrouds in design to increase rotor performance and provide protection, it’s getting harder to shape PCboards into precise 3D features in this small scale and then assemble. Alternatively, other possible materials that are compatible with SDM process were considered. Three polymers, polyurethane, white polyurethane, and epoxy, were proved in the earlier research [11] as suitable materials for SDM process. The material selection for air-frame was based on these three materials. The properties are compared in the following table. Ciba TDT 205-3 Adtech LUC 4180 Adtech EE-501/ Polyurethane (PU) Polyurethane (white PU) 503 Epoxy Material Density (g/cm^3) Tensile Strength (MPa) Curing time

1.12

1.09

1.6

70

78-80

86-89

2-3 hours

12 hours

24 hours

In airframe application, the structural strength were calculated in advance to ensure these three materials are all qualified to provide sufficient strength. The second important issue is the weight due to the very limited weight budget for air-frame. As shown in the table above, epoxy has larger density and is considered as less suitable material, even though it has better strength. Therefore, we focused on the other two (PU and white PU), and investigated the possibilities to reduce the density further. In addition, a rigid polyurethane foam material available in hobby shops was also evaluated. A feasible solution to reduce density is to add glass bubbles, Micro-Balloons (MB), into PU or white PU. Three materials were explored and cast into a block to determine density. The results are compared as follows: Material Density (g/cm^3) Comment PU w/ MB 0.64 Short mixing and handling time before cured PU foam < 0.15 Expansion is not uniform and not controllable White PU w/ MB 0.43~0.51 Long handling time before cured White PU became more attractive than PU because of longer handling time and lower density before and after mixing with MB. However, the longer curing hours of white PU will increase the total fabrication cycle time. The machinability of both materials remains fair after mixing with MB. • Manufacturing The airframe was built by SDM process. Multiple layers and different materials can be manufactured sequentially onto the previous layers without additional assembly. For the final air-frame design, the building direction and required layers were identified. Due to the accessibility of milling tools and thin wall features, the structure was judged to build up side down, and two principal layers are required (see Figure 27, 28).

+Z

Building direction Figure 27 -- design and building direction

Layer 1

Layer 2

Final Part

Figure 28 -- two principal layers A simplified air-frame design, which is similar to final design without 15-degree tilt, was built first by using 3-axis CNC machine in order to evaluate feasibility. The first test implemented the idea of using white PU w/ MB for the main structure and PU foam for shrouds. Since the shroud structure contributes about 2/3 of the total weight, we expected to reduce weight significantly by using PU foam. The result is shown below (Figure 29). The total weight of this version is 1.1 grams. Two materials were bonded nicely without any compatible problem. However, the PU foam is not rigid enough to retain the required shape when the thickness is small. Another experiment was performed to build the whole structure out of white PU w/ MB (Figure 30). The whole air-frame weights 1.7 grams. The shrouds are stiffer to retain the shape, but are still flexible when they are very thin at the bottom. Besides, some cavities were observed in the air-frame due to some air bubbles trapped inside while casting. These cavities may help to reduce weight as long as they don’t appear at the critical places that endanger the structure.

Figure 29 -- 3-axis CNC machined air-frame with PU foam shrouds

Figure 30 -- 3-axis CNC machined airframe

This simplified airframe helped us to identify material and manufacturing issues. White PU w/ MB with effective weight reduction and sufficient stiffness was determined to be the suitable air-frame material. General fabrication sequences and machining strategies are mostly applicable to the final version that needs 5-axis CNC machine for manufacturing. Moreover, in order to improve the stiffness at the bottom of shrouds, a thin layer of PU was added during the fabrication. The fabricated airframe by 5-axis CNC machine is shown in Figure 31. The final part is in one piece without further assembly. The total machining time for the whole air-frame is about 6-8 hours, while the material curing time for each layer is 12 hours.

Figure 31 -- Airframe fabricated by 5-axis CNC machine

2.3 Experimental Tests of Rotor Performance Experimental work has focused on three areas, the evaluation of candidate electric motors, thrust testing of rotor/motor packages, and the low-order validation of the two-dimensional CFD airfoil analyses. Proper characterization of candidate electric motors is critical to optimal rotor design. It is also important for rotor testing where thrust may be easily measured directly, but torque is more simply derived from RPM and motor curves. The experimental work in the area of airfoil performance is driven by the dearth of experimental data at ultra-low Reynolds numbers and the difficulties associated with traditional testing methods under these conditions. The goal of this work is not to provide a high-fidelity validation of the CFD analyses, but to ascertain if the results are accurate enough to discern the effects of significant geometry variations and if the accuracy is sufficient to be applied to the current rotor design methods. 2.3.1 Motor and Rotor Testing The detailed aerodynamic design of the rotors is highly dependent on the performance characteristics of a given motor. Manufacturers commonly provide some form of performance data, but it is often either insufficient, or based on theoretical models rather than detailed experimental data. Test apparatus to accomplish these tasks has been designed and fabricated. The test stand has the capability of evaluating rotor lift for all currently considered motors with the 5mm Smoovy motor as the smallest motor of current interest. Rotor torque is also measurable for motors starting from the 5mm Smoovy. The evaluation of motor performance is accomplished by coupling the shafts of two motors. The one being tested is mounted to the test stand. The second motor acts as a torque load and is mounted in a fixture such that the drive shafts are co-linear. The motor being tested is powered by a constant voltage power supply. The torque generated by the motor is translated to a force on the balance by a load bar that is free to rotate about an axis parallel to the motor shaft. Torque, RPM, and current data are recorded as the load torque is gradually increased until the rated stall current is achieved. From this data, power and efficiency curves are calculated. The testing fixture and electronic balance are pictured below.

Figure 32. Experimental rotor test stand.

Several different test rigs are used for calculating thrust. Two are based on the torque measuring setup. The torque load bar can be replaced by a balance beam that can accept a variety of motor mounts. Different beam lengths allow the testing of rotors for both the mesicopter prototypes and larger test vehicles. This arrangement allows the rotor to be operated out of ground effect and eliminates rotor wash on the balance. The beam can also be counterbalanced to accommodate a wide range of rotor/motor packages. Test data includes thrust, RPM, motor voltage and current. Rotor torque and power required are derived from motor data. Thrust testing in a hostile or environment requires another approach. Preliminary tests of rotors for applications on Mars were carried out inside a Martian atmospheric test chamber at JPL in Pasadena. The low pressures involved and the inability to directly read a balance eliminated the testing method just described. Time and budgetary constraints for this test also did not allow for a more complex solution for direct thrust measurement. The test stand fabricated for this situation is a pendulum setup. This is similar to the approach used for initial thrust testing with the 3mm Smoovy motors during the phase1 research effort, but scaled up considerably to accept the 20 cm diameter Mars rotors. This thrust stand is pictured below, inside the test chamber at JPL. Thrust is estimated by observing the deflection angle of the pendulum during testing.

Figure 33. Mars rotor and test stand in JPL’s Mars environment chamber.

2.3.2 Airfoil Validation via Micro-Glider Flight Testing Classical theory and experimental work at higher Reynolds numbers indicates that the INS-2d results for ultra-low Reynolds number airfoils are reasonable, but detailed experimental verification of airfoil performance estimates is needed. The very low Reynolds numbers of interest make wind tunnel testing essentially impossible. Tunnels or tanks using water, oil, or glycerin, offer the best solution, but such facilities are rare. Even in these fluids the model chords become quite small and direct force determination is difficult. The approach taken here may not provide the resolution and accuracy attainable in a static test facility, but it is relatively inexpensive, allows for the testing of many different airfoils in a simple manner, and should provide a good indicator of the accuracy of the computational performance estimates. Very small gliders have been developed for flight testing at very low Reynolds numbers. The aircraft is a conventional aft-tail configuration with an empty mass of 0.16 grams. The constant chord wing has a span of 5-cm and a 4.5mm chord. The measurement of only three quantities, the mass of the glider, the steady-state glide velocity, and the glide slope, will provide all the information necessary for verification of a particular airfoil’s performance estimates. The gliders are manufactured from polyurethane resin using shape-deposition manufacturing [12]. Fuselage and tail assemblies are standardized, but each wing assembly is designed and manufactured with a desired airfoil geometry. A completed glider is shown in figure 34. Once assembled and trimmed, a range of flight

Reynolds numbers are achievable by simply ballasting the glider to higher wing loadings. Flying unballasted, the prototype design achieves Re 1400 on the wing at 4.25 m/s with a predicted glide ratio of 3.3, ballasting will allow operation up to around Re 3000. Flight at higher Re will require modifications to the wing design to bring the glide speed down.

Figure 34. Microglider for section testing. Predicted performance of the glider is closely coupled to the two-dimensional airfoil analyses. Section data, calculated across a range of Reynolds numbers, is compiled into a database. This database is used by an aircraft performance program in a strip-theory viscous drag buildup. A modified lifting-line method is used to calculate the lift distribution and the induced drag. The performance code has been developed and validated for use with high-performance sailplanes [13]. The database can easily be modified, by simple scaling or some other manner, to match the predicted aircraft performance with the measured steady-state glide ratio. The differences between the modified database and the original will indicate of the accuracy of the two-dimensional computational results. Several prototype gliders have been assembled and flown. After careful trimming, steady-state glides are possible with minimal deviations from a straight flight path. Current work is on developing a system for data acquisition and a method for launching the gliders at a precise velocity. This should greatly reduce the initial transients present during a hand-launch and provide longer periods of steady-state flight. The preliminary flights have already proved useful. The current rotor design tools significantly under-predict the power required for a given RPM. Simple observation of glider test flights indicates that the prototype glider L/D is within 30% of the estimate. This amount of error in the airfoil performance prediction cannot account for the errors in rotor power prediction, so emphasis is being placed on possible power losses from three-dimensional viscous effects and unsteady flow phenomena.

3. Motor and Power Systems 3.1 Motor system The motor system is electro mechanical part of the mesicopter transforming the electrical energy (controller voltage * current) into mechanical energy (shaft torque * angular speed). The weight budget estimations and comparison with currently available micro-motors suggested that the best motor candidate for the self-powered centimeter sized flying vehicle (mesicopter) will be commercially available brushless DC motor from RMB miniature bearings Inc (RMB). RMB offers micro motors in different sizes. Motor diameter 3 mm 5 mm Weight

326 mg

1400 mg

4V

6V

Maximum output power Pmax at Un

94 mW

340 mW

Current at Pmax

19 mA

100 mA

84 %

56 %

55 000 rpm

40 000 rpm

Nominal voltage Un

Efficiency at Pmax rpm at Pmax

Table 1: Weight and performance comparison of RMB micro-motors Table 1 shows weight and the maximum output power operating point parameters at nominal voltage for RMB motors with 3 and 5 mm outer diameter [14,15]. These motors are characterized by high efficiency1 and low weight. Unfortunately a DC brushless motor controller is needed for every motor. Therefore the weight of the controller has to be added to the weight of the motor when compared with DC motors. There are two types of controllers for the DC brushless motors: open- and closed-loop controllers. The maximum efficiency of the DC brushless motor is achieved when the angle between the rotating magnetic field of the permanent magnet rotor and the magnetic field built by current flowing through stator coils is app. 90 degrees. The open loop controller supplies the stator coils with periodic voltage pattern independently from the position of the rotor. This results in a simple circuitry but inefficient and not optimal control behavior of the motor. The closed loop controller uses the angular position of the rotor to generate the stator field that is shifted by app. 90 degrees. There are two approaches how to sense the position of the rotor. The first one uses external sensors mostly Hall sensors, the second one drives the motor at every instant only with two coils. The induced voltage (back EMF), in the third coil, caused by the magnetic field of the rotor, is used to determine the rotor position. The closed loop controller requires more complicated circuitry but guarantees most efficient motor operation. In general the commercially available controller circuits consists of at least two integrated circuits (IC). The first one is the logic IC responsible for the periodic signal generation. The second (IC) is usually a driver stage capable to supply the high currents into the stator coils. RMB offers both open and close loop controllers for their micro motors. Figure 35a shows the closed loop controller from RMB. It is a single IC back EMF DC brushless controller. The IC is a commercially available Phillips TDA 5145. Figure 35b shows Philips TDA 5144 IC in a surface mount device (SMD) package in a minimal configuration with three external SMD capacitors and interconnections. The weight of this controller is about 700 mg. This controller solution is acceptable for the larger (5 mm) micro-motor but too heavy for the smallest (3 mm) micro motor.

1

Only copper losses were taken into account to calculate the motor efficiency.

Figure 35: Back EMF closed loop controller from RMB on right. Back EMF controller using Philips TDA 5144 IC in a SMD package on left.

3.2 Power System Power system is component of the mesicopter supplying electrical energy for the motor system. The general specification for the power system was to supply enough energy for a self-powered mesicopter take-off and short flight. After measuring the performance of motor system (RMB micro motor with available propellers) the mesicopter power requirements were first estimated to 6 V and 600 mA with a weight limit of app. 6g! This corresponds to power density of 600 mW/g.

3.2.1 Batteries The most efficient approach is to supply required voltage and currents directly with batteries. Commercially available batteries with different chemistries were tested and compared see Table 2. The parameters of a specific battery depend on the battery chemistry and on the battery design. Only the NiCd Sony N50AAA can deliver current levels above 600 mA. New lithium based batteries are besides high energy density properties also capable of high drain currents but currently there are no ultra light (< 3 g) lithium batteries available. The cell voltages of tested batteries vary according to the used chemistry between 1.2 V for NiCd batteries up to 3.6 V for lithium based batteries. Therefore several batteries have to be stacked together to get the required voltage level. Cell stacking by putting single cells together is inefficient because the battery package contributes between 40 to 60 % to the total battery weight. Stacking two Sanyo NiCd N50AAA batteries together we already reach the weight 7g. This stack provides only 2.4 V. Series combination of currently available batteries cannot satisfy our power requirement.

Mass Chemistry (g)

Nom Rat Rated Capacity Voltage (V) (Ah)

Estimated Max Power Energy (W) (Wh)

Power Density (W/g)

Energy Density (Wh/g)

N50AAA Sanyo

NiCd

3.52

0.050

1.20

1.961

0.06

0.557

0.02

V20HR

Varta

NiMH

0.89

0.020

1.20

0.004

0.02

0.004

0.03

2L76

Energizer

Li/MnO2

3.28

0.190

3.00

0.787

0.57

0.240

0.17

675 ZA

Rayovac

ZnO2

2.11

0.190

1.55

0.005

0.29

0.002

0.14

G13

Panasonic

AgO2

2.18

0.190

1.55

0.283

0.29

0.130

0.14

Model

Make

Table 2: Performance comparison of commercially available batteries with weight less than 4g.

3.2.2 DC/DC converter One possibility to increase the voltage from batteries is to use voltage multiplier also called DC/DC converter. For currents up to 250 mA capacitor charge pump converters can be used. They typically consist of a single integrated circuit (IC) and two capacitors. Inductor based converters are used for higher currents. They typically consist of a single IC, diode, coil and some filter capacitors. For higher current levels extra power transistor is needed. Table 3 compares two DC/DC converters based on IC’s from company called MAXIM. The first one uses MAX 1703 IC. It is capable to deliver up to 1A but only at 5.5 V. The MAX 608 can supply about 1 A at 10 V from input voltage as low as 2 V. Unfortunately the required coil weights about 2.7 g. Combining the Sanyo N50AAA NiCd 50mAh battery with these converters we can get following power supply configurations: one Sanyo N50AAA cell with MAX 1703 delivering 5.5 V x 1 A at ≈ 5.5 g or two Sanyo N50AAA cells in series with MAX 608 delivering 10 V x 1 A at 11 g. DC/DC converter IC Input Voltage Output Voltage Output current Weight

Maxim 1703 0.9 V 5.5 V ≈1A ≈2g

Maxim 608 1.8 V 10 V ≈1A ≈4g

Table 3: Comparison of DC/DC converters [16,17]

3.2.3 Supercapacitors In the recent years, due to the advance in the capacitor development, the double layer capacitors (DLC) or “supercapacitors” have become a viable alternative solution to the rechargeable batteries. The rechargeable batteries usually weight much more than the DLC and also it takes much longer to recharge them than to recharge the capacitor. Also the capacity of a rechargeable battery decreases gradually with the increasing number of recharge cycles. After 500 to 1000 life cycles the rechargeable batteries lose most of their storage capacity. The DLC could be charged and discharged for more than 1 million times without any change in its energy storage capacity. An import advantage of the DLC is its high current charge and discharge capacity. After being fully charged the terminals of the capacitor could be even shorted without

damaging the capacitor. On the other hand the DLC shows a linear discharge characteristics compared to a flat one of a rechargeable battery. Additionally for any given size, the energy storage capacity of the rechargeable batteries is still many times greater than the double layer capacitors of the same size, [18] Table 4 compares supercapacitors from several manufacturers [19-23] maximum voltage ratings vary between 2.5 and 5 V. Therefore a series combination of several DLC is needed to satisfy our power requirement. Supercapacitors with the highest energy density come from Maxwell Technologies and then from Elna America Inc. The lightest 10 V rated series combination from Maxwell Technologies weights 13.6 g with a capacitance of 1 F. The lightest 10 V rated series combination from Elna America Inc. has a capacitance of 0.25 F and weights 6.32 g. Company

Rated voltage [V]

Capacitance [F]

Weight [g]

Energy [J]

Energydensity [J/g]

5.5 5.5 5.5

0.047 0.100 0.220

3.9 5 9.5

0.71 1.51 3.33

0.18 0.30 0.35

5.5 11

0.047 0.022

6.2 7.5

0.71 1.33

0.11 0.18

2.5 2.5

1.000 3.300

1.58 3

3.13 10.31

1.98 3.44

2.5 2.5

4.000 10.000

3.4 6

12.50 31.25

3.68 5.21

CapXX

5 2.5

0.500 0.330

12 12

6.25 1.03

0.52 0.09

Powerstore

2.75 2.75

0.470 1.000

2.1 3.6

1.78 3.78

0.85 1.05

BestCap -AVX

5.5 5.5 5.5

0.200 0.120 0.060

15 7 5

3.03 1.82 0.91

0.20 0.26 0.18

Tokin

Elna

Maxwell

Table 4: Performance overview of commercially available supercapacitors.

4. Initial Stability and Control Studies The long term goal for the mesicopter program is to autonomously fly a swarm of vehicles in order to carry out cooperative tasks. This requires, first, that a single mesicopter can be flown autonomously. The current mesicopter prototypes are inherently unstable, so the immediate goal is to stabilize a single mesicopter either passively or actively. Dynamic models of the mesicopter with varying design parameters have been created and tested in simulation. Larger scale prototypes have also been built as a controls test bed to evaluate control laws as they are developed.

4.1 Kinematics The mesicopter is envisioned to be a small unmanned rotor craft vehicle utilizing four equally spaced rotors to provide thrust. Alternative designs incorporate ducts around the rotors to increase survivability. The additional weight is assumed to be offset by an increase in thrust gained by the ducted fan arrangement. The ducts do not adversely effect the dynamic stability of the mesicopter, and it is believed they may actually help. However, for the purposes of the initial dynamic analysis, the ducts have been ignored. A top view of the assumed layout for the mesicopter is drawn in the figure below. The larger scale prototypes follow this layout.

Motor Rotor

Frame Speed Controller Battery Pack Micro-Receiver

Figure 36. Top View of Prototype

The mesicopter has two planes of symmetry as seen in the top view. The only control inputs available will be the four voltages across the motors. This in turn varies the speed and therefore the thrust produced by each motor. This control method was chosen over the more conventional rotor craft control method which utilizes articulated rotor hubs in order to save weight and complexity. These four input degrees of freedom translate into four output degrees of freedom. This can be seen by varying groups of motors differentially or collectively as shown in the diagram below. Assuming the "front" of the vehicle lies along the positive X-axis, then differential control of the two motors along the X-axis cause pitch moments about the Y-axis, while differential control of the motors along the Y-axis cause roll moments about the X-axis. Collective

control of all four motors changes the vertical speed and position of the mesicopter. Finally, differential control of pairs of motors produces yawing moments about the Z-axis. Control of the horizontal velocity and position cannot be controlled directly, and are achieved indirectly through roll or pitch to tilt the thrust vector of the mesicopter producing a horizontal thrust component.

X

X Y

Y Z

Z Pitch Control

Roll Control

X

X Y

Y

Z

Z

Yaw Control

Collective

Figure 37. Control degrees of freedom.

4.2 2-D Dynamic Model In order to effectively design a control system for the mesicopter, a dynamic model must be created. The first step is to stabilize the mesicopter in hover, so only the hover dynamics were considered. As mentioned above, the mesicopter has two planes of symmetry, so the initial model was created based on the assumption that the dynamics decouple along these planes for the hover case. Aerodynamic drag forces on the body are also negligible in hover. Therefore, the initial model used for dynamic analysis was a 2-D model aligned with one of the symmetry planes. Figure 38 describes the geometry used for this initial dynamic model.

T1

T2

r

θΒ

X

H2 H1

h l

β

Z

Figiure 38. 2-D Dynamic Model of the Mesicopter In 2-D, the moments produced by the motors are out of plane, therefore they can be ignored. This does not pose a problem in three dimensions since all four rotors are assumed to be operating at the same speed in hover, therefore the torques of the counter-rotating motors cancel. The forces produced by the rotors lie tangential and horizontal to the rotor tip path plane. These forces will be described in detail later. Some design variables were added to determine their affect on the stability of the mesicopter. A tilt angle, β, was

also included in the model to determine its effect on stability. Also, the position of the rotor hubs above the center of gravity can be adjusted in the model. The only other force acting on the mesicopter is gravity. The rigid body dynamics of the mesicopter in the inertial frame hover dynamics are:

M ⋅& x&= −T1 ⋅ sin(θ B + β ) − T2 ⋅ sin(θ B − β ) − H 1 ⋅ cos(θ B + β ) − H 2 ⋅ cos(θ B − β ) M ⋅& z&= −T1 ⋅ cos(θ B + β ) − T2 ⋅ cos(θ B − β ) + H 1 ⋅ sin(θ B + β ) + H 2 ⋅ sin(θ B − β ) + M ⋅ g I y ⋅ θ B = (T1 − T2 ) ⋅ l ⋅ cos( β ) + ( H 1 + H 2 ) ⋅ (h − l ⋅ sin( β )) The rotor forces, H and T, are derived theoretically using classical momentum theory and actuator disc analysis. Assuming non-flexible hingeless blades with no cyclic actuation, the expressions for the tangential and vertical rotor forces in the rotor reference frame can be written as:

T = ρ ⋅ A ⋅ (Ω ⋅ R ) 2 ⋅ CT H = ρ ⋅ A ⋅ (Ω ⋅ R ) 2 ⋅ CH where the tangential and horizontal thrust coefficients are given by:

CT 1 = Tθ ⋅ θ 0 − ⋅ Tλ ⋅ ( z&r + wG ) σ ⋅a 2 2 ⋅ CH = −( H µ + Rµ ) ⋅ ( − x&r − uG ) σ ⋅a The subscript, r, represents velocities in the rotor's frame of reference. The subscript, G, indicates disturbance terms due to gusts that modify the local velocity seen by each rotor. The remaining coefficients are defined as follows:

−1 1 , and Tθ = , Tλ = 6 4 θ 3 ⋅ CT 3 3 ⋅ cd + ⋅ λHP − TW ÿ + H µ + R µ = λHP ⋅ σ ⋅a 4 8 4⋅a CT whereλHP = C 2⋅ µ2 + T 2 These aerodynamic equations were derived using numerous assumptions that hold for rotor craft designed to carry passengers. Future aerodynamic modeling and testing of the rotors will determine if these equations due indeed hold at for the smaller scales we are using where the Reynolds numbers are much smaller. For now, the assumption is that they should hold to at least first order.

4.3 Stability Analysis and 2-D Simulation The model defined above can now be used to analyze the stability of the mesicopter as well as simulate the time response history. The stability of the mesicopter is studied using root locus analysis. Unfortunately, the equations above are highly non-linear which poses problems for classic root locus analysis. In order to use root locus, the dynamic equations must be linearized. Non-linear and linear time response simulations are then used to verify that the linearized equations sufficiently represent the true dynamics of the mesicopter. The linearization is performed by first linearizing the angular deflections in θ about zero for the hover case. Then x&and z&are assumed to be small. Finally the exponents on various state variables are eliminated using binomial expansion. The linearized equations are then used to create root loci versus the cant angle, β, and versus the length of the rotor shaft, h. The latter parameter affects the height of the rotor hub above the center of gravity of the mesicopter. The root loci are plotted below:

Figure 39. Root Locus plot for the 2-D Mesicopter

The root locus versus the rotor shaft length, h, provides the most insight. The black * symbols indicate the operating conditions of the prototype that has been built. The first observation is that varying h has a large affect on the damping, but little affect on natural frequency while varying β has a large affect on natural frequency, but little affect on the damping. In particular, increasing β increases the natural frequency, and increasing h decreases the damping. Also shown on this plot is the fact that for small h and β, the system poles lie in the right half plane along the real axis. The system poles also head into the right half plane if h is increased to a very large value. However, even though the root locus indicates a statically stable system with some rigid body poles at the origin, the natural frequency is 3.5 rad/sec or roughly 1 Hz. It is nearly impossible for a human to operate the vehicle at this bandwidth. The conclusion is that while the mesicopter is stable and formally controllable, with a pilot in the loop, it becomes unstable.

Figure 40. Root Locus plot for the 2-D Mesicopter

The root locus versus h is less useful, but it confirms the relationships of the parameters h and β with the natural frequency and damping. The prototype was designed using reasonable parameters of h and β from the stability analysis above. Simulations of the time response of the prototype are then run in two dimensions using reasonable parameters for h and β chosen above. The goal was to confirm that the system is indeed stable and that the non-linear simulation agrees with the linearized simulation. The results are shown in figure 41. These plots correspond to the time response of the mesicopter to a small initial horizontal velocity. The linear solution reasonably captures the horizontal behavior as well as the angular rotation corresponding to θ. The vertical dynamics which happen to be mostly second order are lost through the linearization. The plots indicate that these motions are on the order of a couple of millimeters, so this does not pose a great problem. These plots show that the system is indeed stable about θ=0 as shown in the root loci plots.

Figure 41. Time Response of Mesicopter with Small Initial Horizontal Velocity

The above simulations are conducted assuming identical motors and rotors with the center of gravity exactly centered along the X-axis. In reality, these assumptions cannot be made due to irregularities in the manufacturing of the rotors and motors, and due to uncertainty in the center of gravity position. Additional simulations were conducted removing these assumptions and allowing the motor/rotor characteristics and the center of gravity position to vary. In this case, the mesicopter will still be dynamically stable, however it will have a constant horizontal and vertical velocity. Under these circumstances, hover is still possible using computer control since the system is still dynamically stable. Finally, analysis using noisy components and disturbance inputs due to gusts are planned for the future. The next goal is to use a 3-D dynamic model to design simple speed controllers that can stabilize the mesicopter in hover. It is unclear right now which state variables must be sensed in order to accomplish this task. Therefore, one use of the 3-D simulator is to determine if the mesicopter can be using a minimal amount of information, for example, by using rate gyros alone. Because the envisioned missions of the mesicopter include meteorological data collection in the upper atmosphere as well as planetary exploration, it is assumed that all the sensors needed for control should be carried on-board the mesicopter. This research all leads to the final goal of autonomous flight of the mesicopter without continual human intervention.

4.4 Closed loop Control Development is underway to automate the flight of the mesicopter. We present the work to date and work in progress. The current goal of our application work is to construct the smallest, autonomous vehicle possible, which possesses a full range of motion and additional environment sensing capability. To meet this goal, the mesicopter must be designed with the necessary sensors and actuators to ensure stability and must still be as small, lightweight, and low power as possible. This mesicopter would be an ideal 3dimensional testbed for studies in distributed control and would be suited for various environmental sensing applications. The minimal space requirements of the vehicle make the mesicopter a desirable candidate for in-lab flight and unrestricted experimentation.

4.4.1 Initial Demonstration In June, we completed our first autonomous demonstration that shows stabilization of the vertical position of a mesicopter translating along a wire. This application was a proof of concept showing that the rotor/motor pairs could easily provide enough lift for vertical stabilization even with the weight of motors, velocity controllers, a receiver, a battery, and an infared LED. The mesicopter’s position was tracked using an external vision system. The control effort was calculated off-board and the velocity command desired was provided through a radio transmitter. We describe the components making up the demonstration in more detail. The feedback process is also described. The demonstration is confined within a 6 foot tall steel frame. Within the frame a wire runs from top to bottom. Figure 42 shows this frame and the wire within. The mesicopter slides along the wire.

Figure 42. Test setup for mesicopter vision-based control.

This mesicopter is shown in Figure 43. The copter weighs 60 grams when all the demonstration specific components are connected. The 4 motors supply 80 grams of lift, sufficient for control purposes. The distance from each motor to the center of the mesicopter is 3 inches (longest dimension is 6 inches). Each motor/rotor pair is connected to an individual velocity controller using lightweight magneto copper wire. The remainder of the wiring is also done with the lightweight magneto copper wire. On top of the mesicopter is a receiver module, which receives a multi-channel pulse width modulated signal and provides the velocity controllers with a single pulse corresponding to it’s respective channel. For external sensing purposes, an infared LED is attached to the front of the copter. A 6 volt battery is attached to the bottom of the copter and powers all on-board circuitry. A metal shaft, in which the wire runs through, is embedded in middle of the copter which restricts the motion of the copter.

Figure 43. Stability and control testbed used for 1-DOF control tests. Off-board, a feedback system monitors and updates the velocity command sent to the mesicopter receiver. This feedback system consists of a vision system for position tracking, a VME based processing system for calculating velocity desired, and a microprocessor, which transforms the velocity command into a suitable pulse width signal. This signal is transmitted via a RC radio. The camera uses CCD technology for image acquisition and outputs an analog signal conforming to the RS-170 standard (monochrome). The infared LED, when in the vision plane of the camera, charges the camera’s CCD array, and provides a small region of high intensity which is necessary for the vision tracking scheme. The camera signal is processed by a VME point tracking vision processing board developed by members of the the Stanford Aerospace Robotics Lab (ARL). This vision board is seated within a VME cage along with two Motorola MV147 processing boards (68020 processors). Each of these processing boards has ControlShell, software for real-time control, running on them. These applications are executed over the operating system VxWorks. One board is used as the vision server for hosting the point tracking board, and the other board runs the mesicopter control loop. The mesicopter control application is a client to the vision server. The server provides the client with the location of a tracked point in the camera’s vision plane. The vision server finds these points using an intensity thresholding scheme. The vision server returns the location and estimated velocity of tracked

points, in this case points from the LED. The client receives this information and uses scale factors to calculate the actual position and velocity of the mesicopter. Position and velocity are updated at 30 Hz. The position and velocity information is fed into a simple PD control law. A bias term is added to this control effort. The bias added is the force required to compensate for gravity induced by the 60 gram vehicle. The control effort is translated into a message that will be sent across a serial port to the pulse width producing microcontroller. A PIC-17C56a processor was designed to provide an 8 channel pulse width modulation command at 40Hz. The output pulse width signal has negligible jitter. The processor updates the command as it receives update messages from a serial input. The output of the processor is sent into the trainer port of a Futaba RC radio. The pulse width signal is sent transparently from the radio to the mesicopter. This system completes a feedback loop which allows the mesicopter to maintain a desired position on the wire of the steel frame. Using the simple control law, we were able to see the mesicopter nearly maintain the desired position on the wire with small oscilliations near the desired posotion. The oscillations were in the form of a stable limit cycle. This limit cycle was likely induced because we did not account for the nonlinearities in the actual plant model.

4.4.2 Planned Activities We will enhance this demonstration by to removing the current VME system and replacing it with a Linux based PC and multiple PCI framegrabbers. There is no longer any need for the slow, restrictive, and specialized VME hardware. The processing power of the system will be greatly enhanced allowing a wider range of applications. Using multiple framegrabbers, we will be able to track multiple LEDs in 3dimensions. Using vision data and gyro data from the mesicopter (to be added), we should be able to calculate the 6 dimensions of position/orientation information necessary to produce a controller for stable flight. The controller used will be a product of work done by the stability control group of the mesicopter team. A team of students from the Flight Research group and Aerospace Robotics Lab are currently developing a Control Shell component based vision solution which will be modular and reusable. This vision solution will consist of selectable acquisition, display, image processing, state, and network components simplifying the process of creating a vision processing system. The new vision system will be based on the specification designed from this group. We estimate two months for the completion of the components necessary for use in with the mesicopter.

5. Prototype and Testbed Development Development of the mesicopter has relied on a variety of prototypes, of different scales, each intended to address particular problems that need to be solved before a true meso-scale vehicle is successful. The following sections summarize the characteristics and objectives of some of these vehicles.

5.1 Initial 3g Vehicle The initial mesicopter work was undertaken to determine whether a very small scale vehicle would be at all feasible. The device was constructed using 3mm Smoovy motors with 1.5cm diameter rotors and is shown in figure 44. With power supplied from an external source the device was used to measure thrust, power, and efficiency of the motors and rotors.

Figure 44. Initial 3g mesicopter and rotor tests.

5.2 Nearer-Term Mesicopter (15g) Since the motor speed controllers and batteries that are currently available commercially were too heavy to provide energy and control for the 3g mesicopter, a larger device was constructed. Figure 45a shows the assembled mesicopter prototype. The total weight of the prototype is currently 17 grams, but we expect that this can be reduced to approximately 15g. Figure 45b shows the comparison of single propeller lift versus mesicopter lift as a function of controller voltage. The shrouds produce a negligible effect on the lift performance due to the current large gap between rotor tips and the shroud. Using this lift information from it can be estimated that at least 15 V at the controller input are needed to take off. Tests with external power supply confirmed this estimation. The mesicopter lifts off at 13 V and is able to hover out of ground effect with 16 V draining almost 1 A. The prototype at 17g currently cannot fly with the power available from the onboard supercapacitors. This will addressed in the next few months of development. 18 16 14 4*cos(15)*(singlepropeller lift)

lift [g]

12 10 8 lift of the mesicpter

6 4 2 0 5

6

7

8

9 10 11 controller voltage [V]

12

13

14

15

Figure 45: (from left) a) Assembled mesicopter, b) Comparison of single propeller lift versus mesicopter lift as a function of controller voltage

5.3 Stability and Control Testbed The small mesicopters have been developed primarily to study aerodynamics and power issues. For development of stability and control strategies, it is more effective to start with larger devices with larger lifting capacity and slower dynamics. A series of stability and control testbeds have been designed and constructed Photos of these prototypes are included below.

Figure 46. Six-inch (rotor tip to rotor tip) stability and control testbed.

Figure 47. Side View Illustrating Rotor Cant Angle

The device weighs approximately 60g with component weights as follows: Frame:

3.35 gm

Motor:

5.65 gm x 4

Rotor:

0.6 gm x 4

Speed Controller:

1.08 gm x 4

Micro-Receiver:

2.34 gm

Batteries:

26.18 gm (Seven cell NiCd pack)

TOTAL:

61.19 gm

The testbeds were built using mostly off-the-shelf components available from hobby retailers. Control is sent to the mesicopter via a standard programmable transmitter with channel mixing described previously. A five channel micro-receiver then sends the command signals to the four micro-scale speed controllers. The speed controllers convert the PWM signal generated by the receiver into a variable voltage to the motors. Power to the motors is supplied using either lithium ion or NiCd cells. The motors are brushed DC motors that are commercially available. Finally, the rotors are manufactured at Stanford since clockwise and counter-clockwise versions are needed. The current frame is balsa. The rotors on this testbed are capable of producing more than 20 gm of thrust each, meaning that the vehicle can hover using only 75% of its maximum thrust. This provides a sufficient margin for maneuvers and limited payload or additional sensors. The maximum thrust may be increased further using optimized rotors as described in the aerodynamics section. It has been used successfully for the vision-based closedloop control tests.

5.4 Payload-Carrying Testbed A larger testbed vehicle was also constructed using larger, coreless DC motors of much higher efficiency. This testbed weighs approximately 100g but the rotors provide approximately 160g of lift making it possible to carry significant payloads. The total dimension is about twice that of the smaller testbed. Rotors were constructed from carbon-fiber.

Figure 48. 100g testbed

5.5 Systems Integration Testbed A systems integration testbed has been constructed for exploring various approaches to communication, sensing, stabilization, and control. Exploration of algorithms and techniques for flight stability and semiautonomous navigation requires a flexible and programmable prototyping platform. Such a platform is needed to validate physical simulations of both the flight characteristics and sensor functions. In order to move rapidly on prototypes, the testbed should be built from readily available chips and components to avoid the time and expense of custom ASICs. Completely custom silicon could, of course, produce a much smaller system but would not be a good investment of time and money at this stage of the project. We chose to design with off the shelf surface mount chips giving us a good tradeoff between effort and small size. The components were chosen for their support of digital communications, sensor support, and reconfigurability. The system consists of two main parts, a custom transmitter on the ground for radio communications and control, and a flyer with onboard micro-controller, radio, and other peripherals. Without motors or battery the flyer weighs 9 grams.

5.5.1 Radio Transmitter and Controller The custom transmitter shown below uses a Great Planes Real Flight Futaba controller with a joystick port output plugged into a board containing a Microchip PIC17C756A micro-controller and a Linx Technologies 418MHz RF digital transmitter. This radio transmits serial data at 4800 baud with a 300 ft range, and can easily be upgraded to a higher bandwidth two-way digital communications link. The PIC micro-controller samples the joystick port, converts the joystick values for roll, pitch, yaw, and throttle, to motor controls for front, left, back, and side motors. It then transmits this data plus button settings to the flyer.

Figure 49. Controller and transmitter for systems integration testbed.

The following equation expresses the conversion. roll pitch yaw throttle

=

0 d

d 0 -d 0 -d 0

front left

−k k −k k 1/4 1/4 1 / 4 1 / 4

back right

where d represents the distance from the center of mass to each motor and k expresses a linear (approximate) relationship between lift and drag. The matrix is inverted in order to solve for the power to the front, left, back, and right motors given the joystick settings for roll, pitch, yaw, and throttle.

5.5.2 Rotary Wing Flyer The flyer, shown below, is essentially a flying printed circuit board. This is in line with an emphasis on multipurpose components. The thin, light, 20mil circuit board also serves as the frame to which the motors are attached with CA glue. The connector for the battery is also the battery support post. The PIC17 micro-controller can easily be extended with additional sensory and communication components as well as additional control software to incorporate their functionality. The surface mount PIC17 reads the incoming motor controls from the Linx Technologies 418MHz receiver via a serial digital interface and one of the PIC17’s two UARTs. The PIC17 uses its 3 on-chip PWM outputs plus a fourth PWM, implemented in software using timers, to drive the motors. L293DD power chips input the PWM and deliver 1.2A (max) of current to each motor. The vehicle with motors and battery weighs 65 grams and generates 80 g of total thrust. The PIC17 software is written in C giving the capability for expansion through use of the 12 10bit A/D inputs, I2C interface, dual UARTs, pulse width modulators, and many robust I/O pins.

Figure 50. PCB forms the structure of this prototype.

5.5.3 Vision Stabilization The flyer has been tested and it lifts off and flies. However, it is very difficult to control without some form of sensor based stabilization. An investigation is underway into the use of a CMOS image sensor and fixed lens as input into a vision stabilization system. The Photobit PB0100 CMOS camera chip has been selected for its flexible register control, dynamic on-chip exposure control, and convenient I2C interface to

the PIC17 micro-controller. The PIC17 and 4800 baud radio are not fast enough to communicate video to an image processor on the ground but initial calculations indicate that the PIC17 is fast enough to sample small regions of the camera image and to compute the optical flow vectors for these regions. These vectors will be used to compute an estimation of the six degree-of-freedom spatial vector describing the vehicle motion. This sensed motion will be used to stabilize the vehicle. The PB0100 chip is shown below on the PCB for the new vision stabilized flyer. The lens in its black plastic lens holder and the PCB for the camera are also shown.

Figure 51. PC board for next version of vision-stabilized prototype. Note lens and CMOS camera chip.

5.5.4 Simulator for Vision Stabilization A simulator is under construction that will facilitate the testing of various algorithms and parameters in vision stabilization. The simulator will include a rigid body dynamic simulation of the flyer, a geometric model of the camera, an anti-aliased ray tracer for sampling the image that the camera sees, a visual environment consisting of a texture mapped ground plane, and a motion estimation and feedback control system. The rigid body dynamic model of the flyer is given by the following Newton Euler equation of motion using spatial vector notation. i

I f i = I i a i = ~i - hi i

~ hi

α

Mi

a

~ where h = m i ~si

mi

0

0

and M i = 0 0

mi 0

0 mi

Modeling the motors as point masses on the x and z axis of the body fixed coordinate frame having its origin at the center of mass, and modeling the battery as a point mass on the –y axis, gives the following inertia tensor. Ii =

2d 2 mk + l2bk

0

0 0

2

0

4d mk 0 0 2d 2 mk + l2bk

where mi is total flyer mass, si is the vector to the center of mass and is equal to zero, d is the distance to each motor, l is the distance to the battery, mk is the mass of a motor, bk is the mass of the battery, and a is the spatial acceleration. The spatial force fi is the sum of the forces acting on the flyer including rotor lift forces, drag forces (torques), gravity, and noise.

0 0 fi =

df1cos − f1sin

− df 4cos 0

0 0 +

f1cos 0

df3cos f3sin

0 0

+

f3cos 0

f 4cos − f 4sin

df 2cos 0 +

0 0

0 f1 + f3 − f 2 − f 4 +k

f 2cos − f 2sin

0 0 0 0

+

0 0 0 Gx Gy

+

ÿn fn

Gz

Here θ is the inward tilt angle of the motor thrust line of force, and k relates drag torques to thrust as a linear approximation. The spatial acceleration is obtained by multiplying the inverse of the constant 6x6 diagonal inertia matrix II times fi. Spatial velocity is computed by integrating the acceleration. Spatial position, using a quaternion to represent angular orientation, is computed by integrating the linear velocity to get the position, and by computing a difference quaternion for the orientation. The difference quaternion being þq e = 0.5ÿþ t and the associated incremental rotation matrix is ∆A = E + 2(~qe ~qe - q e3~qe ) where ω is angular velocity, qe = , E is the identity matrix, and | qe0 qe1 qe2 qe3 | = 1. The camera model comprises a center of projection COP expressed in the body fixed frame, a direction vector, an up vector and a region in the view plane. This geometry determines the rays to be traced to obtain the pixel data. Anti-aliasing is done by casting rays through the corners of pixels, finding their intersections with environment textures, then averaging the texture regions using a quad-tree. The pixel data is used in the following optical flow equation due to Berthold Horn 1986 and made rigorous by Oliver Faugeras 1996. ∇I • V +

∂I =0 ∂t

This equation expresses the assumption that in any small region of the image and within a small step of time, only a translation of the image will occur. Here I is the image as a function of x, y, and t, and V= is the velocity of the translation in image space. Since we can compute the spatial gradient of the image and the time derivative from the pixels, this gives one equation with two unknowns. Additional information needed to solve this comes from a smoothness assumption, giving the following optical flow algorithm for computing V. u = uav - Ix P/D P = Ixuav + Iyvav + It v = vav - Iy P/D D = λ2 + Ix2 + Iy2 where Ix is ∂I/∂x, and uav is a spatial or temporal running average of u. The simulator will explore several control and feedback mechanisms for stabilization using these optical flow vectors. The resulting algorithm will be tested in the prototype flyer.

5.5.5 Communication We are also tracking developments in high bandwidth, short range radio communication. The Bluetooth radios operating in the 2.4GHz ISM band with up to 10 Mbps data rates are especially interesting and are now becoming available. The following table shows bandwidth requirements for various types of data. MPEG2 1900x1200 MPEG2 720x480 MPEG1 360x240 Video Conference MP3 Audio Intercom misc sensors

~20 Mbps 5 Mbps 1-1.5 Mbps 128-350 Kbps 128 Kbps 33 Kbps 0-50 Kbps

5GHz radios ? Bluetooth ? Bluetooth Bluetooth Bluetooth Linx Technologies SC Linx Technologies LC

6. Future Work Work to date has illustrated both the feasibility of meso-scale rotorcraft and the challenges still ahead before practical implementations are possible. As described in our original proposal, we plan to continue work in each of the areas reported here, with special emphasis on stability and control. Many of researchers with whom we have talked, have suggested that a device which is capable of carrying 10-50g of payload would be useful today. While we expect that advances in sensor technologies will provide interesting missions for vehicles carrying less than 1g of payload in 10-40 years, such systems are not presently available. Thus, in addition to our efforts directed at the smallest mesicopters, we will continue to work with the larger vehicles which serve as useful testbeds while providing greater near-term payload capabilities. Work in the second half of the Phase II program will also include increased interaction with the science community to identify current and future opportunities for such vehicles.

7. References 1.

2. 3.

4. 5. 6. 7.

8. 9. 10.

11. 12.

13.

14. 15. 16. 17. 18. 19. 20. 21. 22. 23.

Kroo, I., Kunz, P., “Development of the Mesicopter: A Miniature Autonomous Rotorcraft,” presented at the American Helicopter Society International Vertical Lift Aircraft Design Specialists’ Meeting, San Francisco, CA, Jan. 2000. Kroo, I., Kunz, P., “Meso-scale Flight and Miniature Rotorcraft Development,” Proceedings of a Workshop on Fixed and Flapping Flight at Low Reynolds Numbers, Notre Dame, June 2000. Kunz, P., Kroo, I., “Analysis, Design, and Testing of Airfoils for Use at Ultra-Low Reynolds Numbers,” Proceedings of a Workshop on Fixed and Flapping Flight at Low Reynolds Numbers, Notre Dame, June 2000. Henry Bortmann, “Whirlybugs,” New Scientist, June 5, 1999. Neil Gross, ed., “Developments to Watch: Tiny Helicopters that Go Boldly Where No Man…,” Business Week, June 21, 1999. Rogers, S. E. and Kwak, D., “An Upwind Differencing Scheme for the Steady-state Incompressible Navier-Stokes Equations,” NASA TM 101051, November 1988. Akima, Hiroshi, “A Method of Univariate Interpolation that Has the Accuracy of a Third-Degree Polynomial,” ACM Transactions on Mathematical Software, Vol. 17, No. 3, September 1991, Pages 341-366. Maughmer, Mark D., Somers, Dan M., “Design and Experimental Results for a High-Altitude, Long Endurance Airfoil,” Journal of Aircraft, Vol. 26, No. 2, February 1989, Pages 148-153. McCormick, Barnes W. Jr., Aerodynamics of V/STOL Flight, Dover Publications Inc., Mineola, New York, 1999, pp. 82-84. Gill, Murray, Saunders, “User’s Guide for SNOPT 5.3: A Fortran Package for Large-scale Nonlinear Programming,” Technical Report SOL 98-1, Department of Engineering Economic Systems and Operations Research, Stanford University, Stanford, California, May 1998. John Kietzman, "Rapid Prototyping Polymer Parts via Shape Deposition Manufacturing", Ph.D. Thesis, Stanford University, February 1999. A.G. Cooper, S. Kang, J. W. Kietzman, F. B. Prinz, J. L. Lombardi and L. Weiss, “Automated Fabrication of Complex Molded Parts Using Mold SDM,” Proceedings of the Solid Freeform Fabrication Symposium, University of Texas at Austin, Austin, Texas, August 1998. Kunz, Peter J., “Development of a Software Package for the Assessment of High-Performance Sailplanes,” Master of Science Thesis, Department of Aerospace Engineering, Pennsylvania State University, University Park, Pennsylvania, August 1997. RMB Miniature Bearings, Inc. SYH59002 Datasheet; www.smoovy.com RMB Miniature Bearings, Inc. SYH39001 Datasheet; www.smoovy.com Maxim Integrated Products, Inc; MAX608 Datasheet; www.maximic.com Maxim Integrated Products, Inc; MAX1703 Datasheet; www.maximic.com Elna America Inc.; www.elna-america.com Tokin America Inc., www.tokin.com Maxwell Technologies Inc., www.maxwell.com Cap-XX Pty Ltd.; www.cap-xx.com PowerStor Corporation; www.powerstor.com AVX Corporation, www.avxcorp.com