1*fc37ad8cSJames Wrightimport pandas 2*fc37ad8cSJames Wrightimport matplotlib.pyplot as plt 3*fc37ad8cSJames Wrightimport seaborn as sns 4*fc37ad8cSJames Wrightimport numpy as np 5*fc37ad8cSJames Wrightimport scipy.signal as sig 6*fc37ad8cSJames Wright 7*fc37ad8cSJames Wright 8*fc37ad8cSJames Wrightdef coeff(force, rho=1, u=1, D=1, zspan=0.2): 9*fc37ad8cSJames Wright S = np.pi * D * zspan # surface area 10*fc37ad8cSJames Wright return 2 * force / (rho * u**2 * S) 11*fc37ad8cSJames Wright 12*fc37ad8cSJames Wright 13*fc37ad8cSJames Wrightdef shedding_period(df): 14*fc37ad8cSJames Wright sample = df[df["Time"] > 70] # once the initial transient has passed 15*fc37ad8cSJames Wright peaks, _ = sig.find_peaks(sample["ForceY"]) 16*fc37ad8cSJames Wright period = np.diff(sample["Time"].iloc[peaks]) 17*fc37ad8cSJames Wright return period.mean() 18*fc37ad8cSJames Wright 19*fc37ad8cSJames Wright 20*fc37ad8cSJames Wrightdf = pandas.read_csv("force.csv") 21*fc37ad8cSJames Wrightdf["Drag Coefficient"] = coeff(df["ForceX"]) 22*fc37ad8cSJames Wrightdf["Lift Coefficient"] = coeff(df["ForceY"]) 23*fc37ad8cSJames Wrightperiod = shedding_period(df) 24*fc37ad8cSJames Wright 25*fc37ad8cSJames Wrightsns.set_theme(style="ticks") 26*fc37ad8cSJames Wrightpalette = sns.color_palette() 27*fc37ad8cSJames Wrightfig, ax_drag = plt.subplots() 28*fc37ad8cSJames Wrightax_lift = ax_drag.twinx() 29*fc37ad8cSJames Wright 30*fc37ad8cSJames Wrightsns.lineplot( 31*fc37ad8cSJames Wright data=df, 32*fc37ad8cSJames Wright x="Time", 33*fc37ad8cSJames Wright y="Drag Coefficient", 34*fc37ad8cSJames Wright ax=ax_drag, 35*fc37ad8cSJames Wright color=palette[0]) 36*fc37ad8cSJames Wrightsns.lineplot( 37*fc37ad8cSJames Wright data=df, 38*fc37ad8cSJames Wright x="Time", 39*fc37ad8cSJames Wright y="Lift Coefficient", 40*fc37ad8cSJames Wright ax=ax_lift, 41*fc37ad8cSJames Wright color=palette[1]) 42*fc37ad8cSJames Wrightax_drag.set_title(f"Shedding period {period}") 43*fc37ad8cSJames Wrightax_drag.set_ylim(0.41, 0.49) 44*fc37ad8cSJames Wrightax_drag.tick_params(axis="y", colors=palette[0]) 45*fc37ad8cSJames Wrightax_drag.yaxis.label.set_color(palette[0]) 46*fc37ad8cSJames Wrightax_lift.tick_params(axis="y", colors=palette[1]) 47*fc37ad8cSJames Wrightax_lift.yaxis.label.set_color(palette[1]) 48*fc37ad8cSJames Wright 49*fc37ad8cSJames Wright# plt.savefig("vortexshedding.svg") 50*fc37ad8cSJames Wrightplt.show() 51