Pazy wing time-domain LCO¶

Perform a time-domain limit-cycle oscillation computation for the straight Pazy wing, plotting the deflection of the beam tip over time.

Imports

In [ ]:

from jax import numpy as jnp
from matplotlib import pyplot as plt

from flapjax.aero.flowfields import ConstantFlowField
from flapjax.models.pazy.straight.pazy_wing import generate_pazy_wing

from jax import numpy as jnp from matplotlib import pyplot as plt from flapjax.aero.flowfields import ConstantFlowField from flapjax.models.pazy.straight.pazy_wing import generate_pazy_wing

Case parameters: freestream velocity, angle of attack, tip force/moment perturbation, and total physical simulation time.

In [ ]:

u_inf_mag = 74.0  # freestream velocity
alpha = jnp.deg2rad(0.5)  # angle of attack

# torsional moment to apply to the tip for the initial static solve. This leads to an overprediction in deformation compared to the actual aeroelastic equilbrium, giving the system an initial pertubation
m_tip = 0.5

physical_time = 1.5  # time to run the simulation for, in seconds

u_inf_mag = 74.0 # freestream velocity alpha = jnp.deg2rad(0.5) # angle of attack # torsional moment to apply to the tip for the initial static solve. This leads to an overprediction in deformation compared to the actual aeroelastic equilbrium, giving the system an initial pertubation m_tip = 0.5 physical_time = 1.5 # time to run the simulation for, in seconds

Build the wing at the chosen freestream and AoA.

In [ ]:

case = generate_pazy_wing(
    flowfield=ConstantFlowField(
        u_inf=jnp.array((u_inf_mag, 0.0, 0.0)),
        rho=1.225,
        relative_motion=True,
    ),
    aoa=alpha,
    m=8,
    m_star=80,
    node_multiplier=1,
)

case = generate_pazy_wing( flowfield=ConstantFlowField( u_inf=jnp.array((u_inf_mag, 0.0, 0.0)), rho=1.225, relative_motion=True, ), aoa=alpha, m=8, m_star=80, node_multiplier=1, )

Time-step size is set by the aerodynamic panelling; derive the number of steps required for the requested physical time.

In [ ]:

dt = case.aero.dt
n_tstep = int(physical_time / dt)
print(f"dt={float(dt):.4f} s, n_tstep={n_tstep}")

dt = case.aero.dt n_tstep = int(physical_time / dt) print(f"dt={float(dt):.4f} s, n_tstep={n_tstep}")

Assemble the tip follower force / moment vector used to perturb the static solution.

In [ ]:

f_ext = jnp.zeros((case.structure.n_nodes, 6)).at[-1, 3].set(m_tip)

f_ext = jnp.zeros((case.structure.n_nodes, 6)).at[-1, 3].set(m_tip)

Static equilibrium with the tip perturbation applied. This provides the initial condition for the time-domain solve.

In [ ]:

# wing clamped at root
static_sol = case.static_solve(
    prescribed_dofs=tuple(range(6)), horseshoe=False, f_ext_follower=f_ext
)

# wing clamped at root static_sol = case.static_solve( prescribed_dofs=tuple(range(6)), horseshoe=False, f_ext_follower=f_ext )

Time-domain dynamic solve starting from the perturbed static equilibrium.

In [ ]:

# initialise the dynamic solve with the static solution, prescribing the root DOFs to be fixed
dynamic_sol = case.dynamic_solve(
    init_case=static_sol, prescribed_dofs=tuple(range(6)), n_tstep=n_tstep
)

# initialise the dynamic solve with the static solution, prescribing the root DOFs to be fixed dynamic_sol = case.dynamic_solve( init_case=static_sol, prescribed_dofs=tuple(range(6)), n_tstep=n_tstep )

Uncomment to write files for visualisation in ParaView.

In [ ]:

# dynamic_sol.plot("./out/")

# dynamic_sol.plot("./out/")

Plot beam tip vertical deflection over time.

In [ ]:

tip_z = dynamic_sol.structure.x[:, -1, 2]
t = jnp.arange(n_tstep) * dt

fig, ax = plt.subplots()
ax.plot(t, tip_z)
ax.set_xlabel("Time, s")
ax.set_ylabel("Tip deflection, m")
plt.show()

tip_z = dynamic_sol.structure.x[:, -1, 2] t = jnp.arange(n_tstep) * dt fig, ax = plt.subplots() ax.plot(t, tip_z) ax.set_xlabel("Time, s") ax.set_ylabel("Tip deflection, m") plt.show()