Fixing Computation of Optimization Parameters

[manumerous]

Fixing the Symmetric Sigmoid Function

While validating the model I realized, that the scaling factor of the symmetric sigmoid function, used to interpolate between the linear and plate region of the aerodynamics model used to compute lift and drag force in the x-z body plane, was not yet adapted and therefore resulted in both the linear and plate region coefficients contributing to the complete range of present angle of attack.

This can be seen in the plots below. A function output of 1 means purely the plate model is used and 0 means purely the linear/quadratic model for lift/drag respectively is used.

The scale factor is used as follows:

Old function: scale_factor = 1

New function: scale_factor = 30

Worse Force Prediciton min Y-Direction

After implementing the changes (that can be seen here) it was observed that the acceleration prediction in the y-direction actually got quite a bit worse when testing a pure hover flight log. This was expected to be the case in the x and z directions, where the new changes actually reduce the model fidelity. However, since the changes should not affect the features in the y direction i would have not expected a change there.

All following logs are generated using the standard gazebo vtol vehicle from the simulator.

With scale_fac = 1

c_d_wing_y_offset: -0.00188112081996518

With scale_fac = 30

c_d_wing_y_offset: -0.011940344476886905

As can be seen the The second drag coefficient for the wing in y direction is much larger than in the first case and actually performs worse than just setting the coefficient to zero:

This last version contains only the gravity contribution to acceleration in y direction.

[manumerous]

Investigating the Feature Matrix

First i suspected that maybe something in constructing the feature matrix went wrong. Therefore i compared the feature matrix corresponding to the forces in y direction before and after changing the scale_scale factor.

In the current model there is only one force influencing the acceleration in y direction next to gravity:

F_drag_fuselage/m = c_d_wing_y_offset* v_air_y^2

This is the drag force in y direction.

I therefore computed the difference between the features relating to y in the feature matrix X before and after changing the sigmoid function. Hereby the rmse difference was both zero for the column corresponding to c_d_wing_y_offset and all columns for the y direction.

Therefore it was validated, that the feature matrix in that direction are equal. Since the problem is separable (e.g. c_d_wing_y_offset is inflencing only the accel in y direction and vise versa) the two regressions should result in the same optimal coefficient for the same feature matrix.

Therefore, something with the regression must go wrong in one of the cases.

[manumerous]

Indeed when comparing the resulting coefficients we can see some oddity for the new model:

Scale_fac = 1

coefficients:
  I_xx- I_yy: 0.0
  I_yy-I_zz: -11.170911822704054
  I_zz-I_xx: 0.0
  c_d_ail_lin: 1.1546319456101628e-14
  c_d_ail_quad: 1.7763568394002505e-15
  c_d_ele_lin: 0.0
  c_d_ele_quad: 0.0
  c_d_wing_xy_lin: -0.4431122520849753
  c_d_wing_xy_offset: -0.005524389663418377
  c_d_wing_xy_quad: -0.06814651813986927
  c_d_wing_xy_stall: 0.11674522371199962
  c_d_wing_y_offset: -0.00188112081996518
  c_l_ele_lin: 0.0
  c_l_wing_xy_lin: 0.15795022065762762
  c_l_wing_xy_offset: 0.09801699110743775
  c_l_wing_xy_stall: 0.009855032606916492
  horizontal_c_m_drag_z_lin: 0.0
  horizontal_c_m_drag_z_quad: 0.0
  horizontal_c_m_leaver_lin: 0.0
  horizontal_c_m_leaver_quad: 0.0
  horizontal_c_m_rolling: 0.0
  horizontal_rot_drag_lin: 1.9463597400459776e-15
  horizontal_rot_thrust_lin: -2.123301534595612e-15
  horizontal_rot_thrust_quad: -2.78388423424758e-14
  intercept: -0.02638338747066049
  vert_c_m_drag_z_lin: 0.023637085261420863
  vert_c_m_drag_z_quad: 1.0650043171104302
  vert_c_m_leaver_lin: 6.059270976681351
  vert_c_m_leaver_quad: -27.47046425093598
  vert_c_m_rolling: 2.05819064228688
  vert_rot_drag_lin: -0.03862280597347806
  vert_rot_thrust_lin: -0.13515531456125088
  vert_rot_thrust_quad: 6.429530898923768
log_file: logs/2021-05-03/10_32_27.ulg
metrics:
  R2: 0.9483115908333891

Scale_fac = 30

coefficients:
  I_xx- I_yy: 0.0
  I_yy-I_zz: -9.81127245558316
  I_zz-I_xx: 0.0
  c_d_ail_lin: 6.407499313354492e-07
  c_d_ail_quad: 0.0
  c_d_ele_lin: 0.0
  c_d_ele_quad: 0.0
  c_d_wing_xy_lin: 571398767.529881
  c_d_wing_xy_offset: 251200055.67065364
  c_d_wing_xy_quad: 324885678.99236184
  c_d_wing_xy_stall: 0.007205955345232133
  c_d_wing_y_offset: -3.9264559745788574e-05
  c_l_ele_lin: 0.0
  c_l_wing_xy_lin: -5028208.640236674
  c_l_wing_xy_offset: -5910828.535789433
  c_l_wing_xy_stall: -0.004750092635262541
  horizontal_c_m_drag_z_lin: 0.0
  horizontal_c_m_drag_z_quad: 0.0
  horizontal_c_m_leaver_lin: 0.0
  horizontal_c_m_leaver_quad: 0.0
  horizontal_c_m_rolling: 0.0
  horizontal_rot_drag_lin: -1894.5095026329705
  horizontal_rot_thrust_lin: 5.229480857909259
  horizontal_rot_thrust_quad: 455.25045994998084
  intercept: -0.08124219738458427
  vert_c_m_drag_z_lin: -0.732857391244094
  vert_c_m_drag_z_quad: 9.352136732260952
  vert_c_m_leaver_lin: 6.370684764162895
  vert_c_m_leaver_quad: -29.838398466629076
  vert_c_m_rolling: 2.1771723811935977
  vert_rot_drag_lin: -0.035044292948547016
  vert_rot_thrust_lin: -0.1565760156516847
  vert_rot_thrust_quad: 5.954966841284679
log_file: logs/2021-05-03/10_32_27.ulg
metrics:
  R2: 0.9409315346913837

As can be seen the following coefficients get extremely large in the second case:

c_d_wing_xy_lin: 571398767.529881
c_d_wing_xy_offset: 251200055.67065364
c_d_wing_xy_quad: 324885678.99236184

c_l_wing_xy_lin: -5028208.640236674
c_l_wing_xy_offset: -5910828.535789433

Hereby the regression tries to fit the data by making the coefficients extremly large because the features corresponding to the linear aero region get really small for high angles of attack ( e.g. 1- sym_sigmoid(AoA) ~~ 0.0001).

I suspect that moving towards this really large coefficients in the gradient based methods takes so much time that a stopping criteria is triggered before we are close to the optimal value for c_d_wing_y_offset.

I will try to verify this and add a cropped version of the symmetric sigmoid function that only uses the sigmoid when the values differ by more than e.g. 0.025 from 0 or 1.

[manumerous]

I forgot the second force influencing the acceleration in y direction which is the drag of the rotors:

F_drag_rotor = c_d_rotor u_actuator v_air_perpendicular

Therefore, it makes sense, that the rotor drag constant, which also influences the x direction changes due to other forces in the x direction changing.

@tstastny it seems that you were right about the sign in the end. However it makes little sense to me why the regression would choose that parameter to be negative when the result looks visually much worse than before across the board.

Also I am still not sure why there is such a large offset from 0 right from the beginning of the predicted forces. I think that could have to do with the intercept and the linear regression not knowing, that we regress over 3 different linear and rotational directions. After a short exercise it seems to choose the same intercept for all of them.

[manumerous]

So eliminating the intercept certainly helped to get rid of that static offset mentioned above.

Furthermore, i manually set the fuselage drag coefficient to 0 instead of the physically impossible negative value after fitting all the parameters. Now It can easily bee seen why the regression puts it at a negative value to improve the fit:

It seems that the rotor drag is off to compensate some model error in x direction and is then compensated by a negative fuselage drag. Even though this looks visually much worse, this actually results in a slightly better R2 value than if the regression would choose the fuselage drag coefficient c_d_wing_y_offset to be zero.

c_d_wing_y_offset = 0:

R2: 0.9387121475255327

c_d_wing_y_offset < 0:

R2: 0.9397419425512858

I think this hole issue is a good motivation to switch the optimization framework from a linear regression to a QP with constraints. As could be seen, the cross influence of all those features already present can lead to very counterintuitive behaviors. In that sense it would be beneficial to being able to set sanity boundaries on the coefficients.

On the other hand it looks like there is still some other flaws in the model from looking at x and z predictions. Those also still need to be fixed.

[manumerous]

This graphs were build up using the log logs/2021-05-03/10_32_27.ulg

[manumerous]

I realized that the simulator showed some very quirky behavior here with logging positive acceleration in z direction while the collective thrust should be far than high enough to compensate gravity.

The issue i found was documented under PX4/PX4-SITL_gazebo.

With the reference log the adaptations documented in this issues look fine: