peterdsharpe / AeroSandbox

Aircraft design optimization made fast through modern automatic differentiation. Composable analysis tools for aerodynamics, propulsion, structures, trajectory design, and much more.
https://peterdsharpe.github.io/AeroSandbox/
MIT License
687 stars 111 forks source link

Inconsistent rotation convention in operating_point #89

Closed cdhainaut closed 1 year ago

cdhainaut commented 1 year ago

Hi,

I spotted an inconsistency in rotations conventions of operating_point. Or maybe there is something I do wrong, but if I do:

def test_operating_points():
    op_point = OperatingPoint(velocity=1.0, alpha=1.0, beta=-1.0)

    value_1 = op_point.compute_freestream_direction_geometry_axes()
    value_2 = np.array(
        [
            *op_point.convert_axes(
                x_from=-1,
                y_from=0,
                z_from=0,
                from_axes="wind",
                to_axes="geometry",
            )
        ]
    )
    print(value_1)
    print(value_2)

I get the following output:

[ 0.99969541 -0.01745241  0.01744975]
[0.99969541 0.01745241 0.01744975]

I digged a bit into the code, and it seems that the problem is in the convert_axes function. The current rotation from wind to body is done the following way:

x_b = (cb * ca) * x_from + (-sb * ca) * y_from + (-sa) * z_from
y_b = (sb) * x_from + (cb) * y_from  # Note: z term is 0; not forgotten.
z_b = (cb * sa) * x_from + (-sb * sa) * y_from + (ca) * z_from

Whereas if we assume wind to body rotation should be a rotation around y by -alpha and by z of -beta, the expression should be:

x_b = (cb * ca) * x_from + (sb * ca) * y_from + (-sa) * z_from
y_b = (-sb) * x_from + (cb) * y_from  # Note: z term is 0; not forgotten.
z_b = (cb * sa) * x_from + (sb * sa) * y_from + (ca) * z_from

It gives consistent results with the modifications

Am I wrong?

Thank you

peterdsharpe commented 1 year ago

Hey @carlitador,

Thank you for the report!

I agree that this mismatch is unexpected behavior, and that I would expect value_1 and value_2 to be equal.

However, I think that value_2 is correct and value_1 is incorrect, rather than the other way around.

An airplane with alpha = 1 deg and beta = -1 deg results in an airplane with its nose vector pointed slightly up and to the right of its velocity vector. (Or, in other words, oncoming wind is flowing from the front-left side of the airplane to the back-right.)

Thus, the wind-axes vector [-1, 0, 0]_wind should be roughly equivalent to the geometry-axes vector [0.99, 0.01, 0.01]_geom. This would mean that .convert_axes is correct and instead that the issue lies in .compute_freestream_direction_geometry_axes().

Please let me know if this logic agrees with your thinking or if you're thinking about it a different way - always good to corroborate with a second opinion before making a big update like this!

In the meantime, I'll start digging through the .compute_freestream_direction_geometry_axes() method and see if I can find where the error there is.

Best, Peter

peterdsharpe commented 1 year ago

Update: fixed by commit https://github.com/peterdsharpe/AeroSandbox/commit/05ed3d7770a0f6b77ab90c68e374919b2617ef4d .

This is on develop; will go live on master on next patch release (likely today or tomorrow). Will close issue once released.

Verification of fixed sign convention: This is an asb.VortexLatticeMethod analysis with asb.OperatingPoint(..., beta=10). The resulting side force coefficient CY is now (correctly) negative; previously it was positive.

image

peterdsharpe commented 1 year ago

ASB v4.0.7 published, which includes this fix. Closing.