Closed Zcaic closed 11 months ago
I found ""np.argmin()" doesn't support casadi.
Hey @Zcaic !
This is an area of active work-in-progress that I am extremely excited for! :)
Previously, with XFoil-driven aerodynamics, optimizing on airfoils through the asb.Opti()
interface was not possible (as XFoil was in the loop, which is a black box). Instead, one could use the separate airfoil optimizer code here to "pre-optimize" the airfoils before overall aircraft optimization - but this airfoil optimization was done with non-differentiable (Nelder-Mead simplex via scipy.optimize
) optimization that was separate from asb.Opti()
, and quite slow.
Now, with NeuralFoil, optimizing CST parameters of airfoils (i.e., airfoil shape optimization as part of a larger aircraft design optimization problem) should in theory be possible!
However, the nuts and bolts of the code haven't yet been fully implemented - I imagine this will happen in the next few weeks.
In the time between, you may have luck with the following steps:
asb.Opti()
optimization environment.neuralfoil.get_aero_from_kulfan_parameters()
This function is yet-to-be-documented, but the kulfan_parameters
keyword argument is a dictionary of the form described here.Critically, Step 3 bypasses the conversion to/from a coordinate array, eliminating the argmin
needed to identify the LE index.
This approach may need some massaging to work correctly, but in theory it should work. If you give this a shot, please let me know how it goes! I'm extremely excited to see how NeuralFoil does when stress-tested against airfoil aerodynamic shape optimization.
When I have time to take a stab at this and get it working, I will post a Jupyter notebook in AeroSandbox/tutorials/
and update this issue here.
import aerosandbox as asb
import aerosandbox.numpy as np
from aerosandbox.geometry.airfoil.airfoil_families import get_kulfan_coordinates,get_kulfan_parameters
from neuralfoil import get_aero_from_kulfan_parameters
import logging
logging.basicConfig(level=logging.INFO,format='%(levelname)s - %(asctime)s - %(message)s...')
def track_opt(iter):
logging.info(f"optimizer state at iteration {iter}")
track_lp.append(opti.debug.value(lp))
track_up.append(opti.debug.value(up))
oriairfoil=asb.Airfoil(name='rae2822')
initcst=get_kulfan_parameters(coordinates=oriairfoil.coordinates)
opti = asb.Opti()
lp = opti.variable(init_guess=initcst['lower_weights'], upper_bound=-0.03, lower_bound=-0.5)
up = opti.variable(init_guess=initcst['upper_weights'], upper_bound=0.5, lower_bound=0.05)
mycst = dict(lower_weights=lp, upper_weights=up, TE_thickness=0, leading_edge_weight=0)
rsaero = get_aero_from_kulfan_parameters(kulfan_parameters=mycst, alpha=2, Re=6.98886e+06)
ck = -rsaero["CL"] / rsaero["CD"]
opti.minimize(ck)
track_lp=[]
track_up=[]
sol = opti.solve(callback=track_opt,verbose=False)
print(track_lp[-1])
print(track_up[-1])
print(sol.stats()['iterations']['obj'])
the history of obj is:
[-71.97528142783199, -72.13468982455441, -67.62159392513874, -79.97649494975234, -122.52266687532232, -133.51252039590796, -180.17954523397108, -183.93608137653388, -184.81119028819847, -210.3904160677403, -210.82175351449737,
-232.7252212728867, -243.22206547147832, -251.63075279570074, -254.38170153061793, -255.48376530501085, -255.55609826991125, -255.59643616947838, -255.5969481213508, -255.59695565203143]
it seems to work sucessfully\ then I use optimal solution (lp and up) to verify:
lp = np.array([-0.05640185, -0.03, -0.03, -0.03, -0.03, -0.03, -0.03, -0.03]) # the optimal lp
up = np.array([0.13740123, 0.29148949, 0.29450211, 0.27986522, 0.5, 0.5, 0.17367956, 0.3624671]) # the optimal up
myairfoil=asb.Airfoil(coordinates=get_kulfan_coordinates(lower_weights=lp,upper_weights=up))
rsaero=myairfoil.get_aero_from_neuralfoil(alpha=2,Re=6.98886e+06)
ck=rsaero['CL']/rsaero['CD']
print(f"neuralfoil: {ck}")
xf=asb.XFoil(airfoil=myairfoil,Re=6.98886e+06)
rsaero=xf.alpha(2)
ck=rsaero['CL']/rsaero['CD']
print(f"xfoil: {ck}")
output:
neuralfoil: [255.59443472]
xfoil: [248.25301205]
it works well although the CL/CD has a little difference between neuralfoil and xfoil
This is so awesome! 😄 <3% difference in L/D between NeuralFoil and XFoil after optimization is excellent to hear, and better than I was expecting. You may try supplying model_size="xxxlarge"
to the analysis method call, which may tighten this gap even more.
Great work here!
@peterdsharpe hello, Is it possible to use an external optimizer? for example, I use aerosandbox to get objective and gradient values, and pass to the external gradient optimizer.
I use code below to get objective and gradient:
import casadi as ca
import aerosandbox as asb
opti=asb.Opti()
twi=opti.variable()
...
ap=asb.AeroBuildup(
airplane=myairplane,
op_point=asb.OperatingPoint(
velocity=100,
alpha=twi
)
).run()
ck=ap['CL']/ap['CD']
f=ca.Function('ck',[twi],[ck])
dfdt=ca.Function('dfdt',[tw0],[ca.gradient(ck,twi)])
print(f(0))
print(dfdt(0))
it can work as I expect.
Nice work! This is exactly the right way to access the gradient.
This actually reminds me that I should put together an AeroSandbox tutorial on how to connect AeroSandbox to NLOPT, for users who want to try other optimizers (e.g., gradient-free or global). Perhaps even making this a method of asb.Opti
would be useful!
By the way @Zcaic ,
Check out these two brand-new tutorials which show AeroSandbox airfoil optimization through a few examples:
Oh, it's a so amazing work! @peterdsharpe I have learned these two tutorials Airfoil Analysis and Airfoil Optimization. For optimization with many parameters, It saves a lot of time than using evolutionary algorithms.
But I am curious about what happens about optimization at high Mach. I also tried to optimize at Ma=0.8, Re=1.8e7, the optimizer ran 1000 steps and warn that optimization failed. The optimized airfoil also looks a little strange.
I'm glad to hear that the new tutorials are helpful!
I have also noticed some weird behavior in terms of transonic airfoil optimization - it seems to behave worse than subsonic airfoil optimization.
I think that this is because of a few reasons:
I won't have time for the next few weeks to work on these improvements (busy with other projects), but I hope to get around to it within the next few months. This is definitely an area where contributions would be welcome!
Hello, is Aerosandbox.Opti() not supported for optimization of CST parameters of airfoils? I user opti.varibale to define get_kulfan_coordinates()'s upper_weights and lower_weights, then use get_aero_from_neuralfoil() to get CD, then opti.minimize(CD), but it doesn't seem to work.