src/FDM/UIUCModel/uiuc_coefficients.cpp

   1 /**********************************************************************
   2
   3  FILENAME:     uiuc_coefficients.cpp
   4
   5 ----------------------------------------------------------------------
   6
   7  DESCRIPTION:  computes aggregated aerodynamic coefficients
   8
   9 ----------------------------------------------------------------------
  10
  11  STATUS:       alpha version
  12
  13 ----------------------------------------------------------------------
  14
  15  REFERENCES:
  16
  17 ----------------------------------------------------------------------
  18
  19  HISTORY:      01/29/2000   initial release
  20                02/01/2000   (JS) changed map name from aeroData to
  21                             aeroPart
  22                02/18/2000   (JS) added calls to 1Dinterpolation
  23                             function from CLfa and CDfa switches
  24                02/24/2000   added icing model functions
  25                02/29/2000   (JS) added calls to 2Dinterpolation
  26                             function from CLfade, CDfade, Cmfade,
  27                             CYfada, CYfbetadr, Clfada, Clfbetadr,
  28                             Cnfada, and Cnfbetadr switches
  29                04/15/2000   (JS) broke up into multiple
  30                             uiuc_coef_xxx functions
  31                06/18/2001   (RD) Added initialization of Alpha and
  32                             Beta.  Added aileron_input and rudder_input.
  33                             Added pilot_elev_no, pilot_ail_no, and
  34                             pilot_rud_no.
  35                07/05/2001   (RD) Added pilot_(elev,ail,rud)_no=false
  36                01/11/2002   (AP) Added call to uiuc_bootTime()
  37                12/02/2002   (RD) Moved icing demo interpolations to its
  38                             own function
  39
  40 ----------------------------------------------------------------------
  41
  42  AUTHOR(S):    Bipin Sehgal       <bsehgal@uiuc.edu>
  43                Jeff Scott         <jscott@mail.com>
  44                Robert Deters      <rdeters@uiuc.edu>
  45                Ann Peedikayil     <peedikay@uiuc.edu>
  46 ----------------------------------------------------------------------
  47
  48  VARIABLES:
  49
  50 ----------------------------------------------------------------------
  51
  52  INPUTS:       -V_rel_wind (or U_body)
  53                -dyn_on_speed
  54                -ice on/off
  55                -phugoid on/off
  56
  57 ----------------------------------------------------------------------
  58
  59  OUTPUTS:      -CL
  60                -CD
  61                -Cm
  62                -CY
  63                -Cl
  64                -Cn
  65
  66 ----------------------------------------------------------------------
  67
  68  CALLED BY:    uiuc_wrapper
  69
  70 ----------------------------------------------------------------------
  71
  72  CALLS TO:     uiuc_coef_lift
  73                uiuc_coef_drag
  74                uiuc_coef_pitch
  75                uiuc_coef_sideforce
  76                uiuc_coef_roll
  77                uiuc_coef_yaw
  78                uiuc_controlInput
  79
  80 ----------------------------------------------------------------------
  81
  82  COPYRIGHT:    (C) 2000 by Michael Selig
  83
  84  This program is free software; you can redistribute it and/or
  85  modify it under the terms of the GNU General Public License
  86  as published by the Free Software Foundation.
  87
  88  This program is distributed in the hope that it will be useful,
  89  but WITHOUT ANY WARRANTY; without even the implied warranty of
  90  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  91  GNU General Public License for more details.
  92
  93  You should have received a copy of the GNU General Public License
  94  along with this program; if not, write to the Free Software
  95  Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307,
  96  USA or view http://www.gnu.org/copyleft/gpl.html.
  97
  98 **********************************************************************/
  99
 100 #include "uiuc_coefficients.h"
 101
 102 void uiuc_coefficients(double dt)
 103 {
 104   static string uiuc_coefficients_error = " (from uiuc_coefficients.cpp) ";
 105   double l_trim, l_defl;
 106   double V_rel_wind_dum, U_body_dum;
 107   static bool ap_pah_on_prev = false;
 108   int ap_pah_init = 1;
 109   static bool ap_alh_on_prev = false;
 110   int ap_alh_init = 1;
 111
 112   if (Alpha_init_true && Simtime==0)
 113     Alpha = Alpha_init;
 114
 115   if (Beta_init_true && Simtime==0)
 116     Beta = Beta_init;
 117
 118   // calculate rate derivative nondimensionalization factors
 119   // check if speed is sufficient to compute dynamic pressure terms
 120   if (dyn_on_speed==0)
 121     {
 122       uiuc_warnings_errors(5, uiuc_coefficients_error);
 123     }
 124   if (nondim_rate_V_rel_wind || use_V_rel_wind_2U)         // c172_aero uses V_rel_wind
 125     {
 126       if (V_rel_wind > dyn_on_speed)
 127         {
 128           cbar_2U = cbar / (2.0 * V_rel_wind);
 129           b_2U    = bw /   (2.0 * V_rel_wind);
 130           // ch is the horizontal tail chord
 131           ch_2U   = ch /   (2.0 * V_rel_wind);
 132         }
 133       else if (use_dyn_on_speed_curve1)
 134         {
 135           V_rel_wind_dum = dyn_on_speed_zero + V_rel_wind * (dyn_on_speed - dyn_on_speed_zero)/dyn_on_speed;
 136           cbar_2U        = cbar / (2.0 * V_rel_wind_dum);
 137           b_2U           = bw /   (2.0 * V_rel_wind_dum);
 138           ch_2U          = ch /   (2.0 * V_rel_wind_dum);
 139           Alpha_dot      = 0.0;
 140         }
 141       else
 142         {
 143           cbar_2U   = 0.0;
 144           b_2U      = 0.0;
 145           ch_2U     = 0.0;
 146           Alpha_dot = 0.0;
 147         }
 148     }
 149   else if(use_abs_U_body_2U)      // use absolute(U_body)
 150     {
 151       if (fabs(U_body) > dyn_on_speed)
 152         {
 153           cbar_2U = cbar / (2.0 * fabs(U_body));
 154           b_2U    = bw /   (2.0 * fabs(U_body));
 155           ch_2U   = ch /   (2.0 * fabs(U_body));
 156         }
 157       else if (use_dyn_on_speed_curve1)
 158         {
 159           U_body_dum = dyn_on_speed_zero + fabs(U_body) * (dyn_on_speed - dyn_on_speed_zero)/dyn_on_speed;
 160           cbar_2U    = cbar / (2.0 * U_body_dum);
 161           b_2U       = bw /   (2.0 * U_body_dum);
 162           ch_2U      = ch /   (2.0 * U_body_dum);
 163           Alpha_dot  = 0.0;
 164         }
 165       else
 166         {
 167           cbar_2U   = 0.0;
 168           b_2U      = 0.0;
 169           ch_2U     = 0.0;
 170           Alpha_dot = 0.0;
 171         }
 172     }
 173   else      // use U_body
 174     {
 175       if (U_body > dyn_on_speed)
 176         {
 177           cbar_2U = cbar / (2.0 * U_body);
 178           b_2U    = bw /   (2.0 * U_body);
 179           ch_2U   = ch /   (2.0 * U_body);
 180         }
 181       else if (use_dyn_on_speed_curve1)
 182         {
 183           U_body_dum = dyn_on_speed_zero + U_body * (dyn_on_speed - dyn_on_speed_zero)/dyn_on_speed;
 184           cbar_2U    = cbar / (2.0 * U_body_dum);
 185           b_2U       = bw /   (2.0 * U_body_dum);
 186           ch_2U      = ch /   (2.0 * U_body_dum);
 187           Alpha_dot  = 0.0;
 188           beta_flow_clean_tail = cbar_2U;
 189         }
 190       else
 191         {
 192           cbar_2U   = 0.0;
 193           b_2U      = 0.0;
 194           ch_2U     = 0.0;
 195           Alpha_dot = 0.0;
 196         }
 197     }
 198
 199   // check if speed is sufficient to "trust" Alpha_dot value
 200   if (use_Alpha_dot_on_speed)
 201     {
 202       if (V_rel_wind     < Alpha_dot_on_speed)
 203           Alpha_dot = 0.0;
 204     }
 205
 206
 207   // check to see if icing model engaged
 208   if (ice_model)
 209     {
 210       uiuc_iceboot(dt);
 211       uiuc_ice_eta();
 212     }
 213
 214   // check to see if data files are used for control deflections
 215   if (outside_control == false)
 216     {
 217       pilot_elev_no = false;
 218       pilot_ail_no = false;
 219       pilot_rud_no = false;
 220     }
 221   if (elevator_step || elevator_singlet || elevator_doublet || elevator_input || aileron_input || rudder_input || flap_pos_input)
 222     {
 223       uiuc_controlInput();
 224     }
 225
 226   if (ap_pah_on && icing_demo==false)
 227     {
 228       double V_rel_wind_ms;
 229       V_rel_wind_ms = V_rel_wind * 0.3048;
 230       ap_Theta_ref_rad = ap_Theta_ref_deg * DEG_TO_RAD;
 231       if (ap_pah_on_prev == false)
 232         ap_pah_init = 0;
 233       elevator = pah_ap(Theta, Theta_dot, ap_Theta_ref_rad, V_rel_wind_ms, dt,
 234                         ap_pah_init);
 235       if (elevator*RAD_TO_DEG <= -1*demax)
 236         elevator = -1*demax * DEG_TO_RAD;
 237       if (elevator*RAD_TO_DEG >= demin)
 238         elevator = demin * DEG_TO_RAD;
 239       pilot_elev_no=true;
 240     }
 241
 242   if (ap_alh_on && icing_demo==false)
 243     {
 244       double V_rel_wind_ms;
 245       double Altitude_m;
 246       V_rel_wind_ms = V_rel_wind * 0.3048;
 247       ap_alt_ref_m = ap_alt_ref_ft * 0.3048;
 248       Altitude_m = Altitude * 0.3048;
 249       if (ap_alh_on_prev == false)
 250         ap_alh_init = 0;
 251       elevator = alh_ap(Theta, Theta_dot, ap_alt_ref_m, Altitude_m,
 252                         V_rel_wind_ms, dt, ap_alh_init);
 253       if (elevator*RAD_TO_DEG <= -1*demax)
 254         elevator = -1*demax * DEG_TO_RAD;
 255       if (elevator*RAD_TO_DEG >= demin)
 256         elevator = demin * DEG_TO_RAD;
 257       pilot_elev_no=true;
 258     }
 259
 260   CD = CX = CL = CZ = Cm = CY = Cl = Cn = 0.0;
 261   CLclean_wing = CLiced_wing = CLclean_tail = CLiced_tail = 0.0;
 262   CZclean_wing = CZiced_wing = CZclean_tail = CZiced_tail = 0.0;
 263   CXclean_wing = CXiced_wing = CXclean_tail = CXiced_tail = 0.0;
 264   CL_clean = CL_iced = 0.0;
 265   CY_clean = CY_iced = 0.0;
 266   CD_clean = CD_iced = 0.0;
 267   Cm_iced = Cm_clean = 0.0;
 268   Cl_iced = Cl_clean = 0.0;
 269   Cn_iced = Cn_clean = 0.0;
 270
 271   uiuc_coef_lift();
 272   uiuc_coef_drag();
 273   uiuc_coef_pitch();
 274   uiuc_coef_sideforce();
 275   uiuc_coef_roll();
 276   uiuc_coef_yaw();
 277
 278   //uncomment next line to always run icing functions
 279   //nonlin_ice_case = 1;
 280
 281   if (nonlin_ice_case)
 282     {
 283       if (eta_from_file)
 284         {
 285           if (eta_tail_input) {
 286             double time = Simtime - eta_tail_input_startTime;
 287             eta_tail = uiuc_1Dinterpolation(eta_tail_input_timeArray,
 288                                             eta_tail_input_daArray,
 289                                             eta_tail_input_ntime,
 290                                             time);
 291           }
 292           if (eta_wing_left_input) {
 293             double time = Simtime - eta_wing_left_input_startTime;
 294             eta_wing_left = uiuc_1Dinterpolation(eta_wing_left_input_timeArray,
 295                                                  eta_wing_left_input_daArray,
 296                                                  eta_wing_left_input_ntime,
 297                                                  time);
 298           }
 299           if (eta_wing_right_input) {
 300             double time = Simtime - eta_wing_right_input_startTime;
 301             eta_wing_right = uiuc_1Dinterpolation(eta_wing_right_input_timeArray,
 302                                                   eta_wing_right_input_daArray,
 303                                                   eta_wing_right_input_ntime,
 304                                                   time);
 305           }
 306         }
 307
 308       delta_CL = delta_CD = delta_Cm = delta_Cl = delta_Cn = 0.0;
 309       Calc_Iced_Forces();
 310       add_ice_effects();
 311       tactilefadefI = 0.0;
 312       if (eta_tail == 0.2 && tactile_pitch && tactilefadef)
 313         {
 314           if (tactilefadef_nice == 1)
 315             tactilefadefI = uiuc_3Dinterp_quick(tactilefadef_fArray,
 316                                            tactilefadef_aArray_nice,
 317                                            tactilefadef_deArray_nice,
 318                                            tactilefadef_tactileArray,
 319                                            tactilefadef_na_nice,
 320                                            tactilefadef_nde_nice,
 321                                            tactilefadef_nf,
 322                                            flap_pos,
 323                                            Alpha,
 324                                            elevator);
 325           else
 326             tactilefadefI = uiuc_3Dinterpolation(tactilefadef_fArray,
 327                                             tactilefadef_aArray,
 328                                             tactilefadef_deArray,
 329                                             tactilefadef_tactileArray,
 330                                             tactilefadef_nAlphaArray,
 331                                             tactilefadef_nde,
 332                                             tactilefadef_nf,
 333                                             flap_pos,
 334                                             Alpha,
 335                                             elevator);
 336         }
 337       else if (demo_tactile)
 338         {
 339           double time = Simtime - demo_tactile_startTime;
 340           tactilefadefI = uiuc_1Dinterpolation(demo_tactile_timeArray,
 341                                                demo_tactile_daArray,
 342                                                demo_tactile_ntime,
 343                                                time);
 344         }
 345       if (icing_demo)
 346         uiuc_icing_demo();
 347     }
 348
 349   if (pilot_ail_no)
 350     {
 351       if (aileron <=0)
 352         Lat_control = - aileron / damax * RAD_TO_DEG;
 353       else
 354         Lat_control = - aileron / damin * RAD_TO_DEG;
 355     }
 356
 357   // can go past real limits
 358   // add flag to behave like trim_case2 later
 359   if (pilot_elev_no)
 360     {
 361       if (outside_control == false)
 362         {
 363           l_trim = elevator_tab;
 364           l_defl = elevator - elevator_tab;
 365           if (l_trim <=0 )
 366             Long_trim = l_trim / demax * RAD_TO_DEG;
 367           else
 368             Long_trim = l_trim / demin * RAD_TO_DEG;
 369           if (l_defl <= 0)
 370             Long_control = l_defl / demax * RAD_TO_DEG;
 371           else
 372             Long_control = l_defl / demin * RAD_TO_DEG;
 373         }
 374     }
 375
 376   if (pilot_rud_no)
 377     {
 378       if (rudder <=0)
 379         Rudder_pedal = - rudder / drmax * RAD_TO_DEG;
 380       else
 381         Rudder_pedal = - rudder / drmin * RAD_TO_DEG;
 382     }
 383
 384   return;
 385 }
 386
 387 // end uiuc_coefficients.cpp