• <tr id="yyy80"></tr>
  • <sup id="yyy80"></sup>
  • <tfoot id="yyy80"><noscript id="yyy80"></noscript></tfoot>
  • 99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

    Biomimetic flexible plate actuators are faster and more efficient with a passive attachment

    2016-12-05 07:39:48PeterYehAlexanderAlexeev
    Acta Mechanica Sinica 2016年6期

    Peter D.Yeh·Alexander Alexeev

    Biomimetic flexible plate actuators are faster and more efficient with a passive attachment

    Peter D.Yeh1·Alexander Alexeev1

    ?The Chinese Society of Theoretical and Applied Mechanics;Institute of Mechanics,Chinese Academy of Sciences and Springer-Verlag Berlin Heidelberg 2016

    Using three-dimensional computer simulations, we probe biomimetic free swimming of an internally actuated flexible plate in the regime near the first natural frequency. The plate is driven by an oscillating internal moment approximating the actuation mechanism of a piezoelectric macro fiber composite(MFC)bimorph.We show in our simulations that the addition of a passive attachment increases both swimming velocity and efficiency.Specifically,if the active and passive sections are of similar size,the overall performance is the best.We determine that this optimum is a result of two competing factors.If the passive section is too large,then the actuated portion is unable to generate substantial deflection to create sufficient thrust.On the other hand,a large actuated section leads to a bending pattern that is inefficient at generating thrust especially at higher frequencies.

    Flexible plate·Biomimetic·Lattice Boltzmann· Fluid–structure interaction·Robotic fish·Piezoelectric actuator

    DOI 10.1007/s10409-016-0592-0

    1 Introduction

    The development of efficient aquatic machines that mimic the swimming motion of biological fish is an ongoing challenge.Fish leverage the distributed flexibility of their body and fins in order to swim quickly and with great agility [1].Experimental studies have revealed that individual finstrokes are highly complex and mostly three-dimensional [2–5].Conventional aquatic machines employ servomotors or hydraulic actuators in conjunction with a series of linked rigid sections in order to mimic body flexibility in a single bending plane[6–11].However,generating a large range of three-dimensional fin strokes requires a complicated and cumbersome system of links and pulleys.Although well-studied,these motor-based actuators still cannot match the swimming performance of biological fish.

    On the other hand,internally actuated(active)smart materials can be adapted to perform complex fish-like motions and controlled to follow different strokes without the physical complexity of motor-based designs.Typical smart materials used for biomimetic swimming devices include ionic polymer–metal composites[12–18],shape memory alloys [19–24],magnetostrictive thin films[25–27],and others [28–30].Among other active smart materials,piezoelectric actuators are attractive because of their geometric scalability,high efficiency,and noiseless performance[31].Despite these advantages,the current understanding and usage of piezoelectric actuators in swimming devices is still limited [32–35].

    Macro fiber composites(MFCs)[36–38]are a new class of piezoelectric actuators that have been recently used and tested in small-scale aquatic propulsors[39–41].MFCs offer large dynamic stresses in bending and high performance at low and high frequencies.The thrust performance of MFCs as piezoelectric biomimetic actuators was recently investigated by Erturk and Delporte[41].In this study,the authors fabricated MFC bimorph actuators,which are two MFC laminates bonded together by a thin layer of epoxy.When oppositely signed electric fields are applied to the different layers,one stretches and the other contracts,creating bending.Thus, these flexible MFC bimorph actuators are driven by oscil-lating internal moments.The thrust generation of the MFC bimorph cantilevers was subsequently measured experimentally.Furthermore,these MFC bimorphs were also recently used to build an untethered swimming device that resembled a robotic fish[39].The MFC bimorph acted as the fin,providing oscillating flexible propulsion.A waterproof“body”housed the electronics that powered the bimorphs,allowing for unrestricted movement within a fluid.

    Fig.1 Schematic of internally actuated flexible plate with passive attachment.a Geometry of plate showing the dimensions and size of the active portion(darker shade)and passive portion(lighter shade).b Oscillating internal moment in the active section leads to asymmetric bending patterns that generate propulsion

    In addition,Erturk and Delporte[41]probed the effect of a passive caudal fin attachment by comparing thrust with and without the passive fin.It was found that the passive fin attachment introduced an additional bending mode in the same frequency range,leading to wideband thrust generation.Furthermore,the thrust amplitude was higher for the cantilever beam with a passive attachment.These results suggest that active,internally actuated fins with and without passive attachments behave differently and that a particular combination of active and passive sections may lead to optimal swimming performance.Thus,the purpose of our study is to use computational simulations to systematically study the performance of biomimetic swimmers with active(via internal moment)actuation with a passive attachment and to identify the physical mechanisms that enhance the swimming performance of the composite active–passive fins.

    Simulations facilitate the fundamental understanding of the physics behind biomimetic oscillating propulsion.Computational studies of flexible fins have initially started with simple two-dimensional models with torsional springs[42] or elastic filaments[43–47].Recent advances in computing power have allowed for more complex three-dimensional simulations.Current models are capable of simulating the three-dimensional interaction between the oscillating plate and the surrounding viscous fluid[48–51].Recent computational[52–55]models and also experiments[56–64]have focused on the thrust production and free swimming of flexible plates in order to explore fast and efficient swimming regimes.However,these modeling efforts have exclusively focused on passive fins,whose behavior is qualitatively different from that of internally actuated plates.To the best of our knowledge,the computational study of free swimming of oscillating plates with internal actuation has not yet been undertaken.

    Here,we use computational simulations to study the swimming performance of internally actuated flexible plates with passive attachments as a suitable fin actuator for aquatic swimming devices.Experimentally,the active fins,such as those made of MFCs,could be attached to a power housing enclosing electronic components necessary to control fin beating and to enable unrestrained locomotion of the device. To understand the role of the active fin in the swimmer propulsion,we focus on the hydrodynamics of an internally actuated plate(or fin)with a clamped boundary condition imposed at the plate leading edge to mimic its attachment to a power housing.We do not include the housing in our computational model to isolate the effects governing the plate propulsion. We expect that the power housing would impose a drag force resisting the plate propulsion and decreasing the swimming speed.The magnitude of this resistive drag will strongly depend on the particular geometry of the power housing and, therefore,can greatly vary between different experimental realizations.We expect,however,that the basic trends found in this study will hold,thereby providing practical guidelines for designing efficient biomimetic devices that employ active flexible propulsors.

    Our model actuator is shown in Fig.1a.The active fin is modeled as a flexible plate of uniform bending stiffness with total length L,width w,and thickness b,which is small compared to the length(b? L).The aspect ratio L/w is 2.5.The plate is segmented into two sections.First,an internally actuated section,called the active portion(shown in a darker shade),extends a distance d from the leading edge.Beyond that,a passively responding tail section(the passive portion,lighter shade)is attached.Both sections have the same isotropic material properties.The active portion actuates the swimmer using a sinusoidally oscillating internal moment,M(t)=M0sin(ωt),uniformly distributed within the active portion.Here,ω is the driving frequency, and M0is the amplitude of the applied moment.This actu-ation pattern mimics that of an oscillating bimorph MFC composite[39].

    The leading edge vertical displacement and slope are both kept at zero(clamped in the z-direction),but the plate is allowed to swim forward horizontally.The clamped boundary condition mimics an oscillating fin attached to a power housing with a large enough mass to prevent leading edge deflection.The goal of our study is to investigate the effects of adding the passive elastic attachment on the steady state swimming velocity.To isolate these effects,we hold constant the total length of the swimmer and vary the size of the active portion.The plate is driven at constant frequency,but its stiffness is varied.Depending on the proximity to resonance,different bending patterns emerge that lead to faster or slower swimming velocity.

    2 Computational model

    Our computational model is a fully coupled solver,capturing both the elastic deformation response and the hydrodynamics.The fluid mechanics is modeled using a three dimensional lattice Boltzmann model(LBM),while the solid mechanics is captured by a lattice spring model(LSM).These two models are fully coupled by boundary conditions.This fluid–structure interaction model has been used previously to study the free swimming of elastic plates undergoing plunging actuation[52,53].

    The LBM is a particle based mesoscale method that simulates an incompressible Newtonian fluid.The computational fluid domain is discretized into a cubic lattice of equally spaced“fluid nodes”.At each node,continuous velocity distribution functions,fi(r,ci,t),characterize the flow.The distribution functions represent the mass density of“fluid particles”at position r that propagate with velocity c in specifically defined direction i.The distribution function evolves in time according to the discrete Boltzmann equation [65].We employ a D3Q19 model,meaning that a three-dimensional cubic lattice with 19 distribution functions is used.Hydrodynamic fields such as density ρ,momentum j=ρu,and stresses Π are computed as moments of fias follows To capture the solid mechanics of the plate,we employ an LSM[66,67].In this simple model,the continuous elastic thin plate is discretized into a network of masses connected together by harmonic springs arranged on a regular lattice. As shown in Fig.2,we use a triangular lattice,leading to isotropic behavior and Poisson ratio ν=1/3.The stretching springs that connect nodes have stiffness ks.Note that ksis relatively large so the plate is nearly inextensible.In order to model bending,every collinear triplet of nodes defines one bending unit.At the center node of each bending unit is a torsional spring with stiffness kb,which resists out-of-plane bending.On a triangular lattice,this bending spring configuration is related to the one-dimensional beam bending rigidity EI as the following[67,68]

    The different types of nodes are also depicted in Fig.2a. On the leading edge(left),the open squares indicate nodes that are prescribed to be fixed in the vertical direction.Having two fixed rows leads to a zero slope boundary condition. The darker circles indicate nodes within the active section, while the lighter circles represent nodes forming the passive section.The division between active and passive sections isindicated by the dotted line.All the nodes to the left are within the active section,while those to the right are within the passive section.Both the active and passive section nodes have the same equilibrium spacing Δs=2.325 and identical masses,yielding mass per area ρsb.

    Fig.2 a LSM triangular lattice of flexible plate.Open squares indicate nodes at the boundary that are restricted in the vertical direction.Darker and lighter circles indicate mass nodes in the active and passive portions,respectively.b Distribution of external forces on single bending unit to model a local couple at the center node

    In order to model an internal bending moment within the LSM framework,we apply forces judiciously to bending units parallel to the x-direction within the active section.Figure 2b illustrates the appropriate distribution of the forces on a single bending unit.As shown,a bending unit is composed of three nodes,i?1,i0,and i+1.The distance between nodes is given by r±1,0=|r±1,0|,where r±1,0=r±1?r0is the difference between the nodal position vectors.We apply the bending forces to i+1and i?1in the direction perpendicular to both r±1,0and?y,which is the direction of the imposed internal moment(into the page).We emphasize that the bending units concerned are only the ones parallel to the x-direction. The force directions can be computed by the following

    Using this procedure,a positive bending moment will lead to positive curvature in the x z-plane,presupposing that the entire bending unit does not rotate entirely around so that node i?1and i+1have switched positions.The bending nodal forces are thus given by

    Here,Nr=11 is the number of horizontal rows in the lattice. The opposite force on node i0ensures that the total bending force on a bending unit is zero,representing a local couple at i0.These bending forces are applied to all bending units within the active region.

    The mass nodes experience forces due to the springs,the fluid at the interface,and forces due to the applied moment for active nodes.The displacement of the mass nodes is computed by integrating Newton’s equations of motion using the velocity Verlet algorithm[69,70].Forevery single LBM time step,we perform the LSM integration eight times to ensure numerical stability.

    The computational fluid domain size is 16 L×6 L×8 L with grid spacing Δc=2.In order to capture more accurately the fluid flow,the grid is refined near the plate surface. This refined grid has size 4 L×2 L×3 L and spacing Δf=1 is appropriately coupled with the coarser grid at their boundaries[71].On the outer boundaries of the coarse grid,periodic boundary conditions are applied,modeling an unbounded fluid.

    The solid nodes are independent from the fluid lattice nodes and form a boundary surface that moves within the fluid domain.Our fluid–solid coupling procedure has been used and extensively tested in previous studies[52,53,72–80].Briefly,a local set of three solid nodes define a triangular surface adjacent to fluid nodes.Distribution functions from these adjacent LBM nodes propagate onto the surface and reflect,transferring momentum.This momentum transfer is associated with fluid forces that act on the three surface nodes.At the same time,an interpolated bounce-back condition at the LSM surface is applied in the LBM distribution functions,leading to no-slip and no-penetration boundary conditions[81].Details of the fluid–solid coupling procedure are described in a prior study[76].

    To simulate free locomotion,we impose a uniform flow throughout the fluid domain and let the plate translate freely in the x-direction.The initial flow velocity is chosen so that free swimming velocity relative to the fluid domain is nearly zero so that the plate stays within the refined grid.The large computational domain also ensures that the far field velocity does not change appreciably as the plate moves forward.The simulations are performed for50 oscillation periods to ensure that the plate locomotion reaches a periodic steady state.

    3 Model validation

    In order to validate our LSM model,we calculate the static deflection of the internally actuated plate under the application of a constant moment.We consider a cantilevered plate with active and passive sections.In this case the plate bending deformation is in one dimension and can be modeled by the moment–curvature relationship for beam bending.Mathematically,the curvature κ(s)along the plate length in terms of the arc length coordinate s is expressed as

    In these expressions,M is a constant moment that is applied on the active section for a length d and EI is the one-dimensional bending rigidity.The equation states that the curvature is constant within the active section,while it is zero in the passive section—the deformation extends in a straight line.The plate deformation curve z(x)can be expressed parametrically as follows

    Here,R=EI/M is the radius of curvature.

    Fig.3 Static deflection of an internally actuated plate with active and passive sections.Symbols represent LSM simulations,while lines represent the analytical solution good agreement is found between simulations and theory

    We use our LSM with internal actuation model as described in the previous section to compute the deformation for R/L=9.6 and d/L values of 0.2,0.6,and 1.These bending curves are shown in Fig.3 with symbols.As shown, the plate maintains constant curvature with its slope increasing until the active–passive boundary at x=x(d)≈d.In the passive section,the plate deformation continues in a straight line until its free end.Consequently,the plate with the largest active section exhibits the largest deflection.

    Our results are compared to the analytical solution, depicted using the straight and dotted lines as shown in Fig.3. We find good agreement between the LSM results and the analytical solution,differing by no more than 2%,which is related to the grid resolution of our LSM model.Therefore, the LSM moment actuation model accurately captures the beam bending behavior.

    The LBM model and the coupling procedure between LBM and LSM has been previously validated and used for studies in flexible flapping aerodynamic and free swimming of plunging elastic plates[52,53,72].To examine the effect of grid density on the solution accuracy,we have performed a grid convergence study.We found that doubling the grid size did not change our results by more than 1%.

    4 Results and discussion

    As previously discussed in the introduction,we investigate the free locomotion of sectioned elastic plates that are internally actuated in the section closer to the leading edge.The goal of the study is to understand how the swimming performance is affected by the size ratio between the active and passive sections while keeping total length of the plate constant.In addition,we investigate the bending patterns that arise from varying its bending rigidity.We focus on the regime near the first natural frequency where resonance amplification results in the larger deflections and faster swimming.

    We first identify the parameters that characterize the system and list the corresponding values used in our simulations.All dimensional quantities are given in LBM units. The plate dimensions are w=20 and L=50,and the plate is actuated at a constant driving frequency with period τ=2π/ω=2000.This defines the characteristic velocity Uc=L/τ,which represents the number of body lengths per period.We set the Reynolds number to be a constant in our simulations R e=ρUcL/μ=400,where ρ=1 is the fluid density.This determines the dynamic viscosity μ in the simulations.In addition,the mass ratio χ=ρL/(ρsb)is set to 2.5.The mass ratio characterizes the ratio between added and apparent masses,and this value is appropriate for artificial aquatic robots that propel using oscillating flexible plates. The mass per area ρsb thus can be determined.Moreover, the reduced moment

    The dimensionless parameter characterizing bending rigidity is treated specially.An oscillating flexible plate in fluid has a natural frequency ω1,fthat depends on bending rigidity EI.Consequently,when the bending rigidity is varied,so does the natural frequency.We define the frequency ratio φ=ω/ω1,fto be the ratio between the driving frequency and first natural frequency in a fluid[82],where ω1,fis computed using linear theory for thin,large aspect ratio plates using the following expression

    Here,λ1is the smallest root satisfying 1+cosλcoshλ=0 andΓ is the hydrodynamic function,whose values depend on aspect ratio and recursively on ω1,f.The values of Γ can be found elsewhere[82].Our previous study has shown that for similar plate parameters ω1,festimates the natural frequency in fluid to within 2%[53].The use of the frequency ratio allows us to characterize the current oscillation state by its proximity to resonance.

    For plates with different values of d/L,we sample the range of φ between 0.8 and 2.0 and compute the period-averaged velocity U after the swimmer reaches a steady state. We plot our results of U/Ucversus φ in Fig.4,where different lines indicate different values of d/L.We find that the swimming velocity is maximized in the vicinity of the first natural frequency for each case of d/L.Indeed,the resonance oscillations near the first natural frequency are linked to faster swimming[53].Furthermore,the maximum swimming velocity is found to be roughly the same for all cases, except for that of d/L=0.2.This is surprising,as the fully actuated plate is not the fastest and even slightly slower than the plate with d/L=0.6.At a higher frequency ratio ofφ=2.0,the discrepancy is even more pronounced with the velocity at d/L=0.6 several times faster than that of d/L=1.At higher frequency ratios,having a passive section appears to dramatically enhance the swimming performance. We conclude that at all post-resonance frequencies,the plates with passive attachments(of size ratios d/L=0.4–0.6) achieve either equal or better performance than a fully actuated plate with the same total length.

    Fig.4 Dimensionless swimming velocity as a function of frequency ratio velocity is maximized near the first natural frequency for all active section sizes.At higher frequencies,the plate with either d/L=0.4 or 0.6 achieves the fastest swimming speeds

    To identify the physical mechanism that allows the passive flap to enhance the swimming velocity,we first examine the trailing edge deflection amplitude,δt,0.A larger trailing edge deflection is usually correlated to faster swimming [52,53].In Fig.5a,we plot δt,0/L as a function of frequency ratio φ.The data shows that,similar to swimming velocity, δt,0/L is also maximized near the first natural frequency,and the deflection amplitudes are nearly the same other than a lower amplitude for d/L=0.2.The similar trends between trailing edge deflection and velocity in the resonance regime suggest that a larger deflection leads to faster swimming. This does not hold,however,at higher frequencies when the trailing edge displacement is still nearly the same,but the plates with passive attachments are significantly faster than the fully actuated one.

    To explore the correlation between trailing edge deflection and swimming velocity,we graph the velocity as a function of deflection in Fig.5b.We find that in general,the larger the trailing edge deflection,the larger the swimming velocity.For plates with spatially similar bending patterns,one with larger trailing edge amplitude would create a larger momentum flux behind the swimmer and,thus,generate more thrust.

    Fig.5 a Trailing edge deflection as a function of frequency.Despite higher velocity at higher frequencies for d/L≥ 0.4,trailing edge deflection remains nearly the same for these swimmers and fully actuated swimmer with d/L=1.b Velocity versus trailing edge deflection showing that in general a larger trailing edge deflection implies faster swimming

    This explains the poorer performance of the d/L=0.2 plate.The moment is applied near the leading edge,far away from the trailing edge,and as such is unable to oscillate the larger passive portion with sufficiently large amplitude compared to plates with larger active sections.Despite the general correlation of larger displacement with faster swimming,Fig.5b shows that for the same trailing edge amplitude, the velocity can also vary significantly and can change up to three times.The slower velocities correspond to fully active plates at higher frequency ratios.These observations suggest that the addition of a passive section affects the bending patterns of the fin which significantly alters the swimming performance.

    In order to characterize the differences in the bending patterns,we examine two cases:d/L=0.6 and 1,both at the frequency ratioφ=2.The maximum trailing edge displacement differs by no more than 20%,but surprisingly,the fully actuated plate is nearly five times slower.In Fig.6a,we plot the time evolution during one period of the non-dimensional internal moment M′(t),trailing edge deflection,δt(t)/L, and trailing edge angle,αt(t)/L.We find that the deflection curves(the solid lines with symbols)for both cases are in phase,which is expected because the frequency ratio is the same.The deflection amplitudes are confirmed to be about 20%apart.Also,the deflection curves are nearly out of phase with the applied moment,which again is expected because the frequency ratio is higher than the natural frequency.

    Fig.6 a Time history of trailing edge kinematics for plates with d/L= 1 and 0.6 at φ=2.Note that the trailing edge deflection and angle are out of phase for d/L=1.Snapshots of deflection curve are shown for b d/L=1 and c d/L=0.6.The negative slope at the trailing edge for d/L=1 contributes to loss of swimming performance

    The most pronounced differences between the two cases are in the trailing edge angle curves.For d/L=0.6,the faster plate,the trailing edge angle(dotted line)is found to be in phase with the corresponding deflection.In contrast, the trailing edge angle for d/L=1,the poorer performer, is out of phase with the deflection.The bending patterns are visualized in Fig.6b for d/L=1 and Fig.6c for d/L=0.6 as snapshots of the instantaneous deflection curves during one period.In Fig.6b,the trailing edge angle and deflection are out of phase,so when the deflection is maximized,the angle is negative.The result is a shape that“cups”the fluid, limiting the momentum flux that leads to forward thrust.In contrast,Fig.6c depicts the bending pattern for a plate with a passive attachment.Here,the angle and deflection are in phase,so at the maximum extent of the trailing edge,the angle is positive,leading to a slope that effectively pushes fluid backwards to enhance thrust.

    The physical mechanism leading to these bending patterns can be explained by analyzing the actuation mechanism.An internal moment acts to change the local curvature,which in turn affects the slope of the deflection curve.When the plate is fully actuated,the trailing edge exhibits a significant local curvature because it lies within the actuation region.Thus, we expect that the angle at the trailing edge responds quickly to the applied moment,i.e.,the angle and moment have a small phase lag.Indeed,a smaller phase lag is seen in Fig.6a between the moment and angle time histories for d/L=1 compared to d/L=0.6.Thus,for fully actuated plates,at the time of maximum deflection,the internal moment,which is out of phase with deflection,creates the negative trailing edge slope.In turn,this leads to a loss of swimming performance.

    Fig.7 Phase between trailing edge deflection and angle as a function of frequency ratio.The passive flap causes the phase to stay nearly zero in a wide frequency range for d/L=0.6 and smaller,which contributes to better swimming performance

    The aforementioned effect does not occur with the passive attachment because the internal moment is confined to the active section and does not directly affect the curvature of the passive attachment if the passive section is long enough.Furthermore,this trailing edge curving effect does not occur with lower frequency ratios closer to the natural frequency because the moment,deflection,and angle are all in phase regardless of the size of the active section.Thus,the swimming performance at low frequencies is all rather similar among all swimmers(with the exception of d/L=0.2 whose poor performance is related to small deflection).

    To corroborate further our analysis,we plot in Fig.7 the phase Δθ between the trailing edge deflection and angle for multiple simulations.For active sections sized d/L=0.6 and smaller,the deflection and angle stay nearly in phase for all frequency ratios.This bending pattern is correlated with better performance.However,for d/L=0.8 and 1,phase increases with frequency ratio.This suggests that a too small passive attachment allows the internal moment to affect the trailing edge curvature and change the trailing edge angle detrimentally.

    Thus,we conclude that the addition of the passive attachment leads to a better swimming performance at higher frequencies because it suppresses the detrimental curvature at the trailing edge.At lower frequencies,the passive attachment does not affect the swimming performance unless it is too large,in which case the reduced trailing edge deflection suppresses the generation of adequate thrust.

    The previous discussion was centered on the swimming speed.Another important parameter for characterizing the swimmer performance is the power consumption.The total input power is calculated in our LSM model as the dot product between the external nodal forces FM,iand the nodal velocity vi,summed over all nodes:

    Fig.8 a Power coefficient and b swimming economy as a function of frequency ratio.With decreasing d/L,power decreases while swimming economy increases.Plates with passive flaps exhibit better swimming economy as well as faster velocity

    Here,the brackets indicate period-averaging.The power is normalized by characteristic powerw L to yield the power coefficient CP=P/Pc.In Fig.8a,we plot the characteristic power as a function of frequency ratio for the different values of d/L.We find that the input power is maximized near the first natural frequency,similar to the trends found for the swimming velocity(Fig.4)and the trailing edge deflection(Fig.5a).We also find that the power increases with increasing d/L,which is expected since the deflection in the active section(and by extension vi)is smaller when the active section is smaller.

    To quantify how efficient swimming is,we introduce the swimming economy ε=(U/Uc)/CP,defined as the ratio between dimensionless swimming velocity and power coefficient.This quantity represents the distance traveled by the plate per unit applied work.In Fig.8b,we plot ε as a function of frequency ratio and find that the swimming economy decreases with increasing d/L.Thus,the swimming economy is the best with the smallest active section,but in this case the swimming velocity is slower.To swim fast at a wide range of frequencies,the swimmer should have an active section equal to about a half of its full length,in which case the swimming is more economical than the fully active plate. Thus,an addition of the passive flap is advantageous for both swimming speed and economy.

    Fig.9 Swimming economy as a function of center of mass deflection of internally actuated swimmers.Swimming economy generally decreases with increasing center of mass deflection with the exception of the poorer cases corresponding d/L=1 and 0.8 at high frequencies

    We have previously showed for a plunging passive swimmer that a higher swimming economy is correlated with smaller center of mass displacement[53].Large vertical deflections in the bending pattern lead to strong side vortices and contribute to vortex induced drag,hindering the swimming efficiency[52,63].Here,we plot the swimming economy of internally actuated swimmers as a function of the center of mass deflection in Fig.9 and find a general trend in which a smaller center of mass displacement leads to higher swimming economy.

    Interestingly,we find for the cases showing poor economy at higher frequencies corresponding to d/L=1 and 0.8 that the economy slightly increases with increasing center of mass deflection.This is likely because the detrimental trailing edge angle dominates the swimming performance, whereas the increased viscous losses due to a greater center of mass displacement have a secondary effect.In this case,a reduction of the trailing edge angle can lead to slightly better economy even for an increased center of mass displacement.

    5 Summary

    MFC bimorphs and other active materials are attractive oscillating propulsor alternatives to traditional cumbersome motor-based actuators[6–11]and have been used recently in untethered aquatic devices[39].Our study demonstrates that the addition of passive flaps enhances the swimming speed and efficiency of moment actuated flexible fin propulsors. This conclusion agrees well with recent experimental observations that the addition of a passive flap to active fins leads to increased swimmer thrust and propulsion speed[41].

    Specifically,we used three-dimensional,fully coupled computer simulations to investigate free underwater locomotion of an internally actuated elastic plate with a passive elastic attachment.The internal sinusoidal moment applied to the swimmer modeled piezoelectric actuation of an MFCbimorph.We isolated the hydrodynamics of the plate representing a swimmer fin and tested the swimming performance with different active to passive section size ratios.This allowed us to determine the effect of a passive attachment while holding the total length of the plate constant.

    Our results indicated that the addition of a passive elastic flap improved both swimming speed and swimming efficiency.By probing the trailing edge kinematics and deflection curves,we determined that if the active and passive sections are of similar size,the overall performance is optimal.A large passive attachment was found to hinder overall swimmer deflection leading to slow swimming velocity.When the active section was too large,the internal moment actuation led to swimmer profiles that are not conducive to generating thrust at higher(post-resonance)frequency ratios.Specifically,the internal moment changed the angle of the trailing edge so that it is out of phase with the deflection.

    We characterized efficiency using the swimming economy,and found that the economy increases with the size of the passive attachment.The plate with the largest passive section was the most economical,but also the slowest.For active section sizes in the range d/L=0.4 and 0.6(i.e.,the active and passive sections are similar sized),the swimming velocity was overall the fastest and efficiency was better compared to fully active swimmers.In other words,we showed that by simply adding a passive attachment of similar material properties to the smart plate actuator better swimming performance can be achieved.We correlated this effect with the changes in the swimmer bending pattern.Thus,our results will be useful for the design of faster and more efficient biomimetic underwater propulsion devices.

    Acknowledgments The authors would like to thank L.Cen and A. Erturk for the stimulating and insightful discussions.

    References

    1.Sfakiotakis,M.,Lane,D.M.,Davies,J.B.C.:Review of fish swimming modes for aquatic locomotion.IEEE J.Ocean.Eng.24, 237–252(1999)

    2.Flammang,B.,Lauder,G.:Caudal fin shape modulation and control during acceleration,braking and backing maneuvers in bluegill sunfish,Lepomis macrochirus.J.Exp.Biol.212,277–286(2009)

    3.Esposito,C.J.,Tangorra,J.L.,Flammang,B.E.,et al.:A robotic fish caudal fin:effects of stiffness and motor program on locomotor performance.J.Exp.Biol.215,56–67(2012)

    4.Flammang,B.E.,Lauder,G.V.:Speed-dependent intrinsic caudal fin muscle recruitment during steady swimming in bluegill sunfish, Lepomis macrochirus.J.Exp.Biol.211,587–598(2008)

    5.Lauder,G.V.,Tangorra,J.L.:Fish Locomotion:Biology and Robotics of Body and Fin-based Movements.In:Robot Fish,25–49. Springer,New York(2015)

    6.Barrett,D.S.:Propulsive efficiency of a flexible hull underwater vehicle.[Ph.D.Thesis],Massachusetts Institute of Technology, Boston(1996)

    7.Anderson,J.M.,Chhabra,N.K.:Maneuvering and stability performance of a robotic tuna.Integr.Comp.Biol.42,118–126(2002)

    8.Hirata,K.:Development of experimental fish robot.In:Sixth International Symposium on Marine Engineering 2000,235–240

    9.Ding,R.,Yu,J.,Yang,Q.,et al.:CPG-based dynamics modeling and simulation for a biomimetic amphibious robot.In:2009 IEEE International Conference on Robotics and Biomimetics(ROBIO), 1657–1662.IEEE(2009)

    10.Liu,F.,Lee,K.M.,Yang,C.J.:Hydrodynamics of an undulating fin for a wave-like locomotion system design.IEEE/ASME Trans. Mechatron.17,554–562(2012)

    11.Kopman,V.,Laut,J.,Acquaviva,F.,et al.:Dynamic modeling of a robotic fish propelled by a compliant tail.IEEE J.Ocean.Eng.40, 209–221(2015)

    12.Shahinpoor,M.:Conceptual design,kinematics and dynamics of swimming robotic structures using ionic polymeric gel muscles. Smart Mater.Struct.1,91–95(1992)

    13.Mojarrad,M.,Shahinpoor,M.:Noiseless propulsion for swimming robotic structures using polyelectrolyte ion-exchange membrane. Proc.SPIE 2716,183–192(1996)

    14.Chen,Z.,Um,T.I.,Zhu,J.:Bio-inspired robotic cownose ray propelled by electroactive polymer pectoral fin.In:ASME 2011 International Mechanical Engineering Congress and Exposition, 1–8(2011)

    15.Chen,Z.,Um,T.I.,Bart-Smith,H.:Ionic polymer–metal composite enabled robotic manta ray.Proc.SPIE 7976,797637(2011)

    16.Chen,Z.,Shatara,S.,Tan,X.B.:Modeling of biomimetic robotic fish propelled by an ionic polymer–metal composite caudal fin. IEEE/ASME Trans.Mechatron.15,448–459(2010)

    17.Aureli,M.,Kopman,V.,Porfiri,M.:Free-locomotion of underwater vehicles actuated by ionic polymer metal composites.IEEE/ASME Trans.Mechatron.15,603–614(2010)

    18.Shen,Q.,Wang,T.,Kim,K.J.:A biomimetic underwater vehicle actuated by waves with ionic polymer–metal composite soft sensors.Bioinspir.Biomim.10,055007(2015)

    19.Cho,K.J.,Hawkes,E.,Quinn,C.,et al.:Design,fabrication and analysis of a body-caudal fin propulsion system for a microrobotic fish.In:2008 IEEE International Conference on Robotics and Automation,706–711(2008)

    20.Wang,Z.L.,Hang,G.R.,Li,J.A.,et al.:A micro-robot fish with embedded SMA wire actuated flexible biomimetic fin.Sens.Actuators A 144,354–360(2008)

    21.Shinjo,N.,Swain,G.W.:Use of a shape memory alloy for the design of an oscillatory propulsion system.IEEE J.Ocean.Eng. 29,750–755(2004)

    22.Wang,Z.,Hang,G.,Wang,Y.,et al.:Embedded SMA wire actuated biomimetic fin:a module for biomimetic underwater propulsion. Smart Mater.Struct.17,025039(2008)

    23.Rossi,C.,Colorado,J.,Coral,W.,et al.:Bending continuous structures with SMAs:a novel robotic fish design.Bioinspir.Biomim. 6,045005(2011)

    24.Rossi,C.,Coral,W.,Colorado,J.,et al.:A motor-less and gear-less bio-mimetic robotic fish design.In:2011 IEEE International Conference on Robotics and Automation,3646–3651 (2011)

    25.Zhang,Y.S.,Liu,G.J.:Design,analysis and experiments of a wireless swimming micro robot.In:2005 IEEE International Conference on Mechatronics and Automations,946–951 (2005)

    26.Zhang,Y.S.,Liu,G.J.:Wireless micro biomimetic swimming robot based on giant magnetostrictive films.In:2005 IEEE International Conference on Robotics and Biomimetics,195–200(2006)

    27.Zhang,Y.S.,Liu,G.J.:Wireless swimming microrobot:design, analysis,and experiments.J.Dyn.Syst.Meas.Control131,011004 (2009)

    28.Zhang,Z.,Philen,M.,Neu,W.:A biologically inspired artificial fish using flexible matrix composite actuators:analysis and experiment.Smart Mater.Struct.19,094017(2010)

    29.Zhang,Z.G.:Design and control of a fish-like robot using an electrostatic motor.In:2007 IEEE International Conference on Robotics and Automation,974–979(2007)

    30.Philen,M.,Neu,W.:Hydrodynamic analysis,performance assessment,and actuator design of a flexible tail propulsor in an artificial alligator.Smart Mater.Struct.20,094015(2011)

    31.Erturk,A.,Inman,D.J.:Piezoelectric Energy Harvesting.Wiley, New York(2011)

    32.Fukuda,T.,Kawamoto,A.,Arai,F.,et al.:Steering mechanism and swimming experiment of micro mobile robot in water.In:IEEE Proceedings,Micro Electro Mechanical Systems,300–305(1995)

    33.Wiguna,T.,Heo,S.,Park,H.C.,et al.:Design and experimental parametric study of a fish robot actuated by piezoelectric actuators. J.Intell.Mater.Syst.Struct.20,751–758(2009)

    34.Heo,S.,Wiguna,T.,Park,H.C.,et al.:Effect of an artificial caudal fin on the performance of a biomimetic fish robot propelled by piezoelectric actuators.J.Bionic Eng.4,151–158(2007)

    35.Ming,A.G.,Park,S.,Nagata,Y.,et al.:Development of underwater robots using piezoelectric fiber composite.In:2009 IEEE International Conference on Robotics and Automation,3435–3440(2009)

    36.Wilkie,W.K.,Bryant,R.G.,High,J.H.,et al.:Low-cost piezocomposite actuator for structural control applications.Proc.SPIE 3991, 323–334(2000)

    37.High,J.W.,Wilkie,W.K.:Method of Fabricating NASA-Standard Macro-fiber Composite Piezoelectric Actuators.National Aeronautics and Space Administration,Langley Research Center, Hampton(2003)

    38.Bryant,R.G.:Overview of NASA Langley’s Piezoelectric Ceramic Packaging Technology and Applications.National Aeronautics and Space Administration,Langley Research Center,Hampton(2007)

    39.Cen,L.,Erturk,A.:Bio-inspired aquatic robotics by untethered piezohydroelastic actuation.Bioinspir.Biomim.8,016006(2013)

    40.Erturk,A.,Inman,D.J.:An experimentally validated bimorph cantilever model for piezoelectric energy harvesting from base excitations.Smart Mater.Struct.18,025009(2009)

    41.Erturk,A.,Delporte,G.:Underwater thrust and power generation using flexible piezoelectric composites:an experimental investigation toward self-powered swimmer-sensor platforms.Smart Mater. Struct.20,125013(2011)

    42.Zhang,J.,Liu,N.S.,Lu,X.Y.:Locomotion of a passively flapping flat plate.J.Fluid Mech.659,43–68(2010)

    43.Hua,R.-N.,Zhu,L.,Lu,X.-Y.:Locomotion of a flapping flexible plate.Phys.Fluids 25,121901(2013)

    44.de Sousa,P.J.S.A.F.,Allen,J.J.:Thrust efficiency of harmonically oscillating flexible flat plates.J.Fluid Mech.674,43–66(2011)

    45.Li,G.J.,Zhu,L.D.,Lu,X.Y.:Numerical studies on locomotion performance of fish-like tail fins.J.Hydrodyn.24,488–495(2012)

    46.Huang,W.-X.,Shin,S.J.,Sung,H.J.:Simulation of flexible filaments in a uniform flow by the immersed boundary method.J. Comput.Phys.226,2206–2228(2007)

    47.Lee,J.H.,Huang,W.-X.,Sung,H.J.:Flapping dynamics of a flexible flag in a uniform flow.Fluid Dyn.Res.46,055517(2014)

    48.Facci,A.L.,Porfiri,M.:Analysis of three-dimensional effects in oscillating cantilevers immersed in viscous fluids.J.Fluid.Struct. 38,205–222(2013)

    49.Huang,W.-X.,Chang,C.B.,Sung,H.J.:An improved penalty immersed boundary method for fluid–flexible body interaction.J. Comput.Phys.230,5061–5079(2011)

    50.Engels,T.,Kolomenskiy,D.,Schneider,K.,et al.:Numerical simulation of fluid–structure interaction with the volume penalization method.J.Comput.Phys.281,96–115(2015)

    51.Li,G.-J.,Liu,N.-S.,Lu,X.-Y.:Dynamic performance and wake structure of flapping plates with different shapes.Acta Mech.Sin. 30,800–808(2014)

    52.Yeh,P.D.,Alexeev,A.:Effect of aspect ratio in free-swimming plunging flexible plates.Comput.Fluids 124,220–225(2016)

    53.Yeh,P.D.,Alexeev,A.:Free swimming of an elastic plate plunging at low Reynolds number.Phys.Fluids 26,053604(2014)

    54.Dai,H.,Luo,H.X.,de Sousa,P.J.S.A.F.,et al.:Thrust performance of a flexible low-aspect-ratio pitching plate.Phys.Fluids 24,101903(2012)

    55.Huang,W.-X.,Sung,H.J.:Three-dimensional simulation of a flapping flag in a uniform flow.J.Fluid Mech.653,301–336(2010)

    56.Moored,K.,Dewey,P.,Boschitsch,B.,et al.:Linear instability mechanisms leading to optimally efficient locomotion with flexible propulsors.Phys.Fluids 26,041905(2014)

    57.Quinn,D.B.,Lauder,G.V.,Smits,A.J.:Maximizing the efficiency of a flexible propulsor using experimental optimization.J.Fluid Mech.767,430–448(2015)

    58.Dewey,P.A.,Boschitsch,B.M.,Moored,K.W.,et al.:Scaling laws for the thrust production of flexible pitching panels.J.Fluid Mech. 732,29–46(2013)

    59.Quinn,D.B.,Lauder,G.V.,Smits,A.J.:Scaling the propulsive performance of heaving flexible panels.J.Fluid Mech.738,250–267 (2014)

    60.Alben,S.,Witt,C.,Baker,T.V.,et al.:Dynamics of freely swimming flexible foils.Phys.Fluids 24,051901(2012)

    61.Pi?eirua,M.,Godoy-Diana,R.,Thiria,B.:Resistive thrust production can be as crucial as added mass mechanisms for inertial undulatory swimmers.Phys.Rev.E 92,021001(2015)

    62.Marais,C.,Thiria,B.,Wesfreid,J.E.,et al.:Stabilizing effect of flexibility in the wake of a flapping foil.J.Fluid Mech.710,659–669(2012)

    63.Raspa,V.,Ramananarivo,S.,Thiria,B.,et al.:Vortex-induced drag and the role of aspect ratio in undulatory swimmers.Phys.Fluids 26,041701(2014)

    64.Shelley,M.J.,Zhang,J.:Flapping and bending bodies interacting with fluid flows.Annu.Rev.Fluid Mech.43,449–465(2011)

    65.Ladd,A.J.C.,Verberg,R.:Lattice-Boltzmann simulations of particle-fluid suspensions.J.Stat.Phys.104,1191–1251(2001)

    66.Ostoja-Starzewski,M.:Lattice models in micromechanics.Appl. Mech.Rev.55,35–60(2002)

    67.Buxton,G.A.,Verberg,R.,Jasnow,D.,et al.:Newtonian fluid meets an elastic solid:coupling lattice Boltzmann and lattice-spring models.Phys.Rev.E 71,056707(2005)

    68.Mao,W.:Modeling particle suspensions using lattice Boltzmann method.[Ph.D.Thesis],Georgia Institute of Technology,Atlanta, USA(2013)

    69.Verlet,L.:Computer“experiments”on classical fluids.I.Thermodynamical properties of Lennard-Jones molecules.Phys.Rev.159, 98–103(1967)

    70.Verlet,L.:Computer“experiments”on classical fluids.II.Equilibrium correlation functions.Phys.Rev.165,201–214(1968)

    71.Chen,H.,Filippova,O.,Hoch,J.,et al.:Grid refinement in lattice Boltzmann methods based on volumetric formulation.Physica A 362,158–167(2006)

    72.Mills,Z.G.,Aziz,B.,Alexeev,A.:Beating synthetic cilia enhance heat transport in microfluidic channels.Soft Matter 8,11508–11513(2012)

    73.Masoud,H.,Bingham,B.I.,Alexeev,A.:Designing maneuverable micro-swimmers actuated by responsive gel.Soft Matter 8,8944–8951(2012)

    74.Alexeev,A.,Balazs,A.C.:Designing smart systems to selectively entrap and burst microcapsules.Soft Matter 3,1500–1505(2007)

    75.Ballard,M.,Mills,Z.G.,Beckworth,S.,et al.:Enhancing nanoparticle deposition using actuated synthetic cilia.Microfluid. Nanofluid.17,317–324(2014)

    76.Alexeev,A.,Verberg,R.,Balazs,A.C.:Modeling the motion of microcapsules on compliant polymeric surfaces.Macromolecules 38,10244–10260(2005)

    77.Masoud,H.,Alexeev,A.:Resonance of flexible flapping wings at low Reynolds number.Phys.Rev.E 81,056304(2010)

    78.Mao,W.,Alexeev,A.:Hydrodynamic sorting of microparticles by size in ridged microchannels.Phys.Fluids 23,051704(2011)

    79.Mao,W.B.,Alexeev,A.:Motion of spheroid particles in shear flow with inertia.J.Fluid Mech.749,145–166(2014)

    80.Kilimnik,A.,Mao,W.,Alexeev,A.:Inertial migration of deformable capsules in channel flow.Phys.Fluids 23,123302 (2011)

    81.Bouzidi,M.,Firdaouss,M.,Lallemand,P.:Momentum transfer of a Boltzmann-lattice fluid with boundaries.Phys.Fluids 13,3452–3459(2001)

    82.Van Eysden,C.A.,Sader,J.E.:Frequency response of cantilever beams immersed in viscous fluids with applications to the atomic force microscope:arbitrary mode order.J.Appl.Phys.101,044908 (2007)

    ? Alexander Alexeev alexander.alexeev@me.gatech.edu http://cfms.gatech.edu

    1George W.Woodruff School of Mechanical Engineering, Georgia Institute of Technology,Atlanta,USA

    7 December 2015/Revised:4 March 2016/Accepted:6 May 2016/Published online:7 September 2016

    美女国产高潮福利片在线看| 手机成人av网站| 久久香蕉国产精品| 成人国语在线视频| 99国产精品一区二区蜜桃av| 两个人免费观看高清视频| 久久久久久免费高清国产稀缺| 国产片内射在线| 少妇裸体淫交视频免费看高清 | 后天国语完整版免费观看| 亚洲专区中文字幕在线| 国产熟女xx| 99在线视频只有这里精品首页| 男女做爰动态图高潮gif福利片| 亚洲精品一卡2卡三卡4卡5卡| 精品熟女少妇八av免费久了| videosex国产| 波多野结衣av一区二区av| 成年版毛片免费区| 精品欧美一区二区三区在线| 日韩欧美一区二区三区在线观看| 成人国产综合亚洲| 男女那种视频在线观看| 麻豆成人午夜福利视频| 亚洲国产精品sss在线观看| 中文字幕久久专区| 亚洲九九香蕉| 97碰自拍视频| 校园春色视频在线观看| 熟妇人妻久久中文字幕3abv| 欧美日韩福利视频一区二区| 午夜视频精品福利| 老司机福利观看| 日韩av在线大香蕉| 1024香蕉在线观看| 身体一侧抽搐| 国产成人精品久久二区二区免费| 国产视频内射| 亚洲一码二码三码区别大吗| 欧美黑人精品巨大| 欧美日韩黄片免| 精品福利观看| 成人国产综合亚洲| www国产在线视频色| 少妇裸体淫交视频免费看高清 | 在线观看舔阴道视频| 国产亚洲精品av在线| 国产精品永久免费网站| 欧美丝袜亚洲另类 | 日韩有码中文字幕| 99国产综合亚洲精品| 久久国产乱子伦精品免费另类| 动漫黄色视频在线观看| 亚洲一区二区三区色噜噜| 欧美一区二区精品小视频在线| 久久欧美精品欧美久久欧美| 亚洲真实伦在线观看| 成年版毛片免费区| 99国产精品一区二区三区| 在线播放国产精品三级| 久久热在线av| 黄色毛片三级朝国网站| 老司机深夜福利视频在线观看| 亚洲电影在线观看av| 日本成人三级电影网站| 日本免费一区二区三区高清不卡| 欧美+亚洲+日韩+国产| 亚洲专区国产一区二区| www日本在线高清视频| 亚洲专区字幕在线| 50天的宝宝边吃奶边哭怎么回事| 黄色a级毛片大全视频| 好看av亚洲va欧美ⅴa在| 老熟妇乱子伦视频在线观看| 久久久久国产精品人妻aⅴ院| 亚洲中文日韩欧美视频| 99久久99久久久精品蜜桃| 亚洲国产欧美网| 法律面前人人平等表现在哪些方面| 免费人成视频x8x8入口观看| 免费看美女性在线毛片视频| 免费观看人在逋| 久久精品国产亚洲av高清一级| 国产亚洲精品久久久久5区| 国产不卡一卡二| 中文亚洲av片在线观看爽| 国产激情偷乱视频一区二区| 国产精品 欧美亚洲| 桃红色精品国产亚洲av| 久久久精品欧美日韩精品| 在线免费观看的www视频| 欧美激情 高清一区二区三区| 国产欧美日韩一区二区三| 亚洲国产欧美日韩在线播放| 成人永久免费在线观看视频| 伦理电影免费视频| aaaaa片日本免费| 久久性视频一级片| 免费人成视频x8x8入口观看| 日本a在线网址| 色精品久久人妻99蜜桃| 国产色视频综合| 久久久久久久午夜电影| 婷婷六月久久综合丁香| 男人舔女人下体高潮全视频| 国产精品一区二区三区四区久久 | 免费女性裸体啪啪无遮挡网站| 亚洲一码二码三码区别大吗| 制服丝袜大香蕉在线| 亚洲av熟女| 久久性视频一级片| 成年女人毛片免费观看观看9| 亚洲欧洲精品一区二区精品久久久| 一级片免费观看大全| 久久久久久久精品吃奶| 国产精品爽爽va在线观看网站 | 丝袜美腿诱惑在线| 大型av网站在线播放| 18美女黄网站色大片免费观看| 久久午夜亚洲精品久久| 在线观看日韩欧美| 禁无遮挡网站| 18禁裸乳无遮挡免费网站照片 | 精品一区二区三区四区五区乱码| 人成视频在线观看免费观看| 免费搜索国产男女视频| 亚洲五月婷婷丁香| 好男人在线观看高清免费视频 | 国产伦在线观看视频一区| 熟妇人妻久久中文字幕3abv| 精品国产超薄肉色丝袜足j| 亚洲精品中文字幕在线视频| 色哟哟哟哟哟哟| 午夜福利免费观看在线| 国产精品,欧美在线| 国产欧美日韩一区二区三| 亚洲国产精品999在线| 久久国产精品男人的天堂亚洲| 国产成年人精品一区二区| a级毛片a级免费在线| 日本a在线网址| 99热只有精品国产| 国内揄拍国产精品人妻在线 | 亚洲熟妇熟女久久| 久久香蕉国产精品| 哪里可以看免费的av片| 国产欧美日韩精品亚洲av| 中文字幕人妻丝袜一区二区| 丰满的人妻完整版| 香蕉av资源在线| 亚洲真实伦在线观看| 人人妻,人人澡人人爽秒播| 久久性视频一级片| 99国产综合亚洲精品| 国产一区二区在线av高清观看| 一区二区日韩欧美中文字幕| 伊人久久大香线蕉亚洲五| 少妇粗大呻吟视频| 亚洲成人久久性| 色哟哟哟哟哟哟| 欧美亚洲日本最大视频资源| 欧美成人午夜精品| 少妇 在线观看| 成年版毛片免费区| 欧美不卡视频在线免费观看 | 九色国产91popny在线| 伦理电影免费视频| 中文字幕最新亚洲高清| www.熟女人妻精品国产| 免费看日本二区| 亚洲人成伊人成综合网2020| 夜夜夜夜夜久久久久| 中文字幕最新亚洲高清| 精品久久久久久成人av| 成年女人毛片免费观看观看9| 成人三级做爰电影| 日日夜夜操网爽| 亚洲精品一卡2卡三卡4卡5卡| 2021天堂中文幕一二区在线观 | 欧美成狂野欧美在线观看| 国产精品美女特级片免费视频播放器 | 久久午夜亚洲精品久久| 亚洲男人的天堂狠狠| 日韩欧美国产在线观看| 日本在线视频免费播放| 欧美一级毛片孕妇| 伦理电影免费视频| 免费一级毛片在线播放高清视频| 久久久久国产精品人妻aⅴ院| 12—13女人毛片做爰片一| 波多野结衣高清无吗| 亚洲 国产 在线| 午夜视频精品福利| 老熟妇仑乱视频hdxx| avwww免费| 欧美性猛交╳xxx乱大交人| 亚洲,欧美精品.| 久久人妻福利社区极品人妻图片| 一边摸一边抽搐一进一小说| 一边摸一边抽搐一进一小说| xxxwww97欧美| 亚洲 欧美一区二区三区| 亚洲欧美精品综合一区二区三区| 少妇粗大呻吟视频| 国产又色又爽无遮挡免费看| 听说在线观看完整版免费高清| 后天国语完整版免费观看| 成人18禁在线播放| 天天一区二区日本电影三级| 人人妻,人人澡人人爽秒播| 91在线观看av| 国产伦一二天堂av在线观看| 国产精品98久久久久久宅男小说| 国产男靠女视频免费网站| 美女高潮喷水抽搐中文字幕| 日本免费a在线| 免费搜索国产男女视频| 久久性视频一级片| 免费在线观看黄色视频的| 桃色一区二区三区在线观看| 欧洲精品卡2卡3卡4卡5卡区| 99精品在免费线老司机午夜| 两性午夜刺激爽爽歪歪视频在线观看 | a在线观看视频网站| 精品午夜福利视频在线观看一区| 亚洲av片天天在线观看| 亚洲欧美日韩无卡精品| 欧美日韩乱码在线| 国产成人欧美在线观看| 好看av亚洲va欧美ⅴa在| 日韩欧美国产一区二区入口| 在线天堂中文资源库| 欧美日韩福利视频一区二区| 精品久久久久久久久久免费视频| 在线观看www视频免费| 国产精品二区激情视频| 久久亚洲真实| 色av中文字幕| 午夜激情av网站| 香蕉av资源在线| 欧美乱妇无乱码| 久久青草综合色| 老熟妇仑乱视频hdxx| 女性生殖器流出的白浆| 老司机午夜十八禁免费视频| 麻豆一二三区av精品| 亚洲专区中文字幕在线| 成在线人永久免费视频| 久久这里只有精品19| 精品卡一卡二卡四卡免费| 男人舔奶头视频| 黄色视频不卡| 久久性视频一级片| 国产精品自产拍在线观看55亚洲| 婷婷丁香在线五月| 国产精品一区二区精品视频观看| 深夜精品福利| 中亚洲国语对白在线视频| 黑人欧美特级aaaaaa片| 日本撒尿小便嘘嘘汇集6| 又黄又粗又硬又大视频| av天堂在线播放| 两个人视频免费观看高清| 国产精品 欧美亚洲| 亚洲 欧美 日韩 在线 免费| 日韩三级视频一区二区三区| 欧美人与性动交α欧美精品济南到| 99热这里只有精品一区 | 成年版毛片免费区| 国产成人系列免费观看| 日日干狠狠操夜夜爽| 亚洲第一电影网av| 精品一区二区三区视频在线观看免费| 19禁男女啪啪无遮挡网站| 成熟少妇高潮喷水视频| 妹子高潮喷水视频| 国产精品二区激情视频| 在线观看免费午夜福利视频| 免费一级毛片在线播放高清视频| 亚洲精品国产区一区二| 国产精品久久久人人做人人爽| 日本 欧美在线| 老司机午夜福利在线观看视频| 丁香六月欧美| 久久香蕉精品热| 久久狼人影院| 国产精品爽爽va在线观看网站 | 日本一本二区三区精品| 久久久久国产精品人妻aⅴ院| 男女午夜视频在线观看| 国产精品av久久久久免费| 最新在线观看一区二区三区| 亚洲va日本ⅴa欧美va伊人久久| 国产高清视频在线播放一区| 精品国产美女av久久久久小说| 视频在线观看一区二区三区| 国产精品野战在线观看| 婷婷精品国产亚洲av在线| 嫩草影院精品99| 久久精品国产99精品国产亚洲性色| 成年版毛片免费区| 欧美黑人欧美精品刺激| 美女 人体艺术 gogo| 亚洲一区中文字幕在线| 久久中文看片网| 国产成人一区二区三区免费视频网站| 一本综合久久免费| 日韩精品免费视频一区二区三区| av超薄肉色丝袜交足视频| 日韩大尺度精品在线看网址| 日韩欧美国产一区二区入口| 成在线人永久免费视频| 制服丝袜大香蕉在线| 久久精品91无色码中文字幕| 狂野欧美激情性xxxx| 久久香蕉精品热| 女性生殖器流出的白浆| 久久久久精品国产欧美久久久| 男人舔女人的私密视频| 啪啪无遮挡十八禁网站| 看片在线看免费视频| 在线视频色国产色| 午夜免费激情av| 欧美日本视频| 香蕉av资源在线| 亚洲成a人片在线一区二区| 国产精品1区2区在线观看.| 在线免费观看的www视频| 可以在线观看的亚洲视频| 久久久久国产精品人妻aⅴ院| 国产成年人精品一区二区| 午夜精品久久久久久毛片777| 国产亚洲精品综合一区在线观看 | 午夜视频精品福利| 91av网站免费观看| 老司机午夜福利在线观看视频| 亚洲五月色婷婷综合| 久久中文看片网| 十八禁人妻一区二区| 国产99久久九九免费精品| 国产精品影院久久| 这个男人来自地球电影免费观看| 成人一区二区视频在线观看| 男人舔女人下体高潮全视频| 色播亚洲综合网| 亚洲在线自拍视频| 村上凉子中文字幕在线| 欧美乱妇无乱码| 欧美激情久久久久久爽电影| 久久精品国产亚洲av香蕉五月| 俺也久久电影网| 久久久久久久久中文| 男人舔女人的私密视频| 黄色片一级片一级黄色片| 黄片播放在线免费| 动漫黄色视频在线观看| 精品久久久久久久毛片微露脸| 欧美另类亚洲清纯唯美| 免费在线观看影片大全网站| 哪里可以看免费的av片| 日本一区二区免费在线视频| 一区二区三区国产精品乱码| 美女午夜性视频免费| 大香蕉久久成人网| 精品卡一卡二卡四卡免费| 久久国产精品男人的天堂亚洲| 成人手机av| 亚洲国产欧美网| 亚洲国产看品久久| 国产欧美日韩一区二区精品| 脱女人内裤的视频| 亚洲专区国产一区二区| 色在线成人网| 国产成人影院久久av| 成人精品一区二区免费| 黑人巨大精品欧美一区二区mp4| 男女午夜视频在线观看| 黄片大片在线免费观看| 国内久久婷婷六月综合欲色啪| 欧美绝顶高潮抽搐喷水| 宅男免费午夜| 国产真实乱freesex| 不卡av一区二区三区| 国产精品二区激情视频| 免费在线观看黄色视频的| 天天一区二区日本电影三级| 欧美激情 高清一区二区三区| 国产麻豆成人av免费视频| 搞女人的毛片| 欧美午夜高清在线| 最新在线观看一区二区三区| 午夜激情av网站| 中文字幕人妻丝袜一区二区| 热99re8久久精品国产| 不卡av一区二区三区| 国产成年人精品一区二区| 久9热在线精品视频| 精品久久久久久久人妻蜜臀av| 国产av在哪里看| а√天堂www在线а√下载| www.精华液| 国产视频一区二区在线看| 男女床上黄色一级片免费看| 三级毛片av免费| 国产精品亚洲美女久久久| 久久亚洲精品不卡| 桃红色精品国产亚洲av| 9191精品国产免费久久| 国产高清激情床上av| 香蕉国产在线看| 久久国产亚洲av麻豆专区| 国产男靠女视频免费网站| av视频在线观看入口| 一卡2卡三卡四卡精品乱码亚洲| 白带黄色成豆腐渣| 久久婷婷人人爽人人干人人爱| 中国美女看黄片| 观看免费一级毛片| 999精品在线视频| 亚洲av日韩精品久久久久久密| 亚洲免费av在线视频| 在线十欧美十亚洲十日本专区| 国产精品香港三级国产av潘金莲| 国产精品久久视频播放| 精品久久久久久久久久免费视频| 国产高清视频在线播放一区| or卡值多少钱| 色综合欧美亚洲国产小说| 亚洲专区国产一区二区| 十八禁网站免费在线| 99热6这里只有精品| 国产免费男女视频| 亚洲精品久久成人aⅴ小说| 午夜福利视频1000在线观看| 可以在线观看毛片的网站| 亚洲精品一区av在线观看| 亚洲aⅴ乱码一区二区在线播放 | 69av精品久久久久久| 久久国产精品人妻蜜桃| 18禁观看日本| 亚洲九九香蕉| 亚洲性夜色夜夜综合| 十分钟在线观看高清视频www| 日韩国内少妇激情av| 岛国视频午夜一区免费看| 欧美一级毛片孕妇| 日韩欧美 国产精品| 成人午夜高清在线视频 | 无遮挡黄片免费观看| 人成视频在线观看免费观看| 男女那种视频在线观看| 亚洲真实伦在线观看| 夜夜爽天天搞| 男男h啪啪无遮挡| 非洲黑人性xxxx精品又粗又长| 国产单亲对白刺激| 午夜激情av网站| 日韩视频一区二区在线观看| 日本一本二区三区精品| 日本 av在线| 99久久精品国产亚洲精品| 身体一侧抽搐| 免费高清视频大片| 国产精品永久免费网站| 神马国产精品三级电影在线观看 | 90打野战视频偷拍视频| 国产亚洲欧美精品永久| 国产99白浆流出| 国产精品爽爽va在线观看网站 | 亚洲精品一卡2卡三卡4卡5卡| 成人免费观看视频高清| 成人18禁高潮啪啪吃奶动态图| 国产精品久久久av美女十八| 日本 欧美在线| 国产精品一区二区免费欧美| 一a级毛片在线观看| 97超级碰碰碰精品色视频在线观看| 女同久久另类99精品国产91| 天天添夜夜摸| 国内精品久久久久精免费| 午夜福利一区二区在线看| 亚洲男人的天堂狠狠| 无遮挡黄片免费观看| 国产极品粉嫩免费观看在线| 国产黄片美女视频| 成人亚洲精品一区在线观看| 国内久久婷婷六月综合欲色啪| 精品久久久久久,| 热99re8久久精品国产| 国产一区二区三区视频了| 久久性视频一级片| 国产精品一区二区精品视频观看| 国产亚洲精品久久久久5区| 国产亚洲精品综合一区在线观看 | 国产一级毛片七仙女欲春2 | 脱女人内裤的视频| 久久天堂一区二区三区四区| 丰满人妻熟妇乱又伦精品不卡| 午夜福利在线在线| a在线观看视频网站| 又黄又爽又免费观看的视频| 久久久久久大精品| 亚洲精品av麻豆狂野| 国产成年人精品一区二区| 在线播放国产精品三级| 午夜福利欧美成人| 欧美色欧美亚洲另类二区| 熟女电影av网| 啪啪无遮挡十八禁网站| 搡老熟女国产l中国老女人| 老汉色av国产亚洲站长工具| 欧美国产日韩亚洲一区| 日韩成人在线观看一区二区三区| 亚洲成人久久爱视频| 一区二区三区激情视频| 人成视频在线观看免费观看| 国产亚洲精品一区二区www| 长腿黑丝高跟| 国产高清videossex| 巨乳人妻的诱惑在线观看| 黄色片一级片一级黄色片| 国产视频一区二区在线看| 免费人成视频x8x8入口观看| 欧美国产精品va在线观看不卡| 少妇裸体淫交视频免费看高清 | 大型av网站在线播放| 亚洲午夜理论影院| 亚洲精品美女久久久久99蜜臀| 亚洲男人的天堂狠狠| 精品人妻1区二区| 麻豆国产av国片精品| 欧美激情久久久久久爽电影| 亚洲人成77777在线视频| 国产成人av激情在线播放| 日韩欧美一区视频在线观看| 97超级碰碰碰精品色视频在线观看| 黄色 视频免费看| 又大又爽又粗| 手机成人av网站| 免费观看人在逋| 亚洲av成人一区二区三| 长腿黑丝高跟| 午夜久久久久精精品| 国产一区在线观看成人免费| 日本撒尿小便嘘嘘汇集6| 日韩欧美 国产精品| 国产成人影院久久av| 精品福利观看| 一夜夜www| 女人爽到高潮嗷嗷叫在线视频| 国产一区二区激情短视频| 久久热在线av| 一区二区日韩欧美中文字幕| 夜夜夜夜夜久久久久| ponron亚洲| 91国产中文字幕| 在线观看66精品国产| 亚洲黑人精品在线| 亚洲最大成人中文| 欧美激情高清一区二区三区| 色播在线永久视频| 美女高潮到喷水免费观看| 精品久久久久久久久久久久久 | 琪琪午夜伦伦电影理论片6080| 国产精品一区二区免费欧美| 欧美黑人精品巨大| 亚洲精品国产一区二区精华液| 国产精品香港三级国产av潘金莲| 亚洲av电影不卡..在线观看| 精华霜和精华液先用哪个| 欧美中文日本在线观看视频| 国产一区在线观看成人免费| 日本撒尿小便嘘嘘汇集6| 国产一区二区三区视频了| 俺也久久电影网| 久热爱精品视频在线9| 啪啪无遮挡十八禁网站| 成年版毛片免费区| 母亲3免费完整高清在线观看| 亚洲五月色婷婷综合| 久久九九热精品免费| 欧美亚洲日本最大视频资源| 亚洲人成电影免费在线| 男人舔女人的私密视频| 国产视频内射| 搡老岳熟女国产| 视频区欧美日本亚洲| 国产精品一区二区免费欧美| 国产精品免费一区二区三区在线| 国产aⅴ精品一区二区三区波| 在线免费观看的www视频| 亚洲欧美日韩高清在线视频| 黄色视频不卡| 在线国产一区二区在线| 制服丝袜大香蕉在线| 亚洲午夜精品一区,二区,三区| 午夜免费激情av| 每晚都被弄得嗷嗷叫到高潮| 99热6这里只有精品| 亚洲美女黄片视频| 中出人妻视频一区二区| av欧美777| 精品欧美一区二区三区在线| 免费看a级黄色片| 久久久久久亚洲精品国产蜜桃av| 国产av在哪里看| 国产精品亚洲av一区麻豆| 国产一区在线观看成人免费| 亚洲国产精品久久男人天堂| 欧美激情 高清一区二区三区| 国产99久久九九免费精品| 久久久国产精品麻豆| 国产精品乱码一区二三区的特点| 一级片免费观看大全| 精品一区二区三区视频在线观看免费|