Thread: [Jsbsim-cvslogs] CVS: JSBSim/utils/aeromatic++/Systems Controls.cpp, 1.48, 1.49

SourceForge Headquarters 225 Broadway Suite 1600 San Diego, CA 92101 +1 (858) 422-6466

Update of /cvsroot/jsbsim/JSBSim/utils/aeromatic++/Systems
In directory sfp-cvs-1.v30.ch3.sourceforge.com:/tmp/cvs-serv10827/Systems

Modified Files:
	Controls.cpp 
Log Message:
Fix a number of issues that arise when the minimum number of input parameters was provided. Especially when no stall speed was specified

Index: Controls.cpp
===================================================================
RCS file: /cvsroot/jsbsim/JSBSim/utils/aeromatic++/Systems/Controls.cpp,v
retrieving revision 1.48
retrieving revision 1.49
diff -C2 -r1.48 -r1.49
*** Controls.cpp	20 Aug 2016 09:00:58 -0000	1.48
--- Controls.cpp	3 Sep 2016 09:45:56 -0000	1.49
***************
*** 113,117 ****
      float CL = W/Q/Sw;
  
!     if (Vs)
      {
          // *** CLmax based on wing geometry and stall speed ***
--- 113,117 ----
      float CL = W/Q/Sw;
  
!     if (Vs > 0.5f)
      {
          // *** CLmax based on wing geometry and stall speed ***
***************
*** 150,154 ****
      float nh = _aircraft->_htail.efficiency;
      float Ee = _aircraft->_htail.flap_ratio;	// elevator
-     float dcgx = -(cg_loc[X] - _aircraft->_aero_rp[X])*INCH_TO_FEET;
      float ch = cbarw*sqrtf(Sh/Sw);
  
--- 150,153 ----
***************
*** 244,253 ****
      M2 = M*M;
  
!     float dsdB = 0.0f;          // ds/dB
!     float CYbeta = -nv*(Sv/Sw)*CLav[0]*(1.0f+dsdB);
!     float CYp_const = -(AR+cosf(sweep))/(AR+4*cosf(sweep))*tanf(sweep);
  
      float CL0 = _aircraft->_CL0;
      float CLalpha = _aircraft->_CLalpha[0];
  
      float zw = -0.0f; // z-pos. wing: positive down
--- 243,258 ----
      M2 = M*M;
  
!     float dsdB = -_aircraft->_vtail.de_da;          // ds/dB
!     float CYbeta = -nv*Sv/Sw*CLav[0]*(1.0f-dsdB);
  
      float CL0 = _aircraft->_CL0;
+     if (Vs > 0.5f)
+     {
+         CL0 = CLaw[0]*(iw - a0w)+(Sh/Sw)*nh*CLah[0]*(ih - E0);
+         _aircraft->_CL0 = CL0;
+     }
+ 
      float CLalpha = _aircraft->_CLalpha[0];
+     _aircraft->_CDalpha.at(0) = CLalpha*(2.0f*CL0)/(PI*AR*Ew);
  
      float zw = -0.0f; // z-pos. wing: positive down
***************
*** 255,263 ****
      float Clbwf = 1.2f*sqrt(AR)*((zw+2.0f*D)/(bw*bw));
      float Clbvt = -(zv/bw)*CLah[0];
- 
-     CL0 = CLaw[0]*(iw - a0w)+(Sh/Sw)*nh*CLah[0]*(ih - E0);
-     _aircraft->_CDalpha.at(0) = CLalpha*(2.0f*CL0)/(PI*AR*Ew);
- 
-     _aircraft->_CL0 = CL0;
      for (int i=0; i<4; ++i)
      {
--- 260,263 ----
***************
*** 265,270 ****
          {
          case 0:
!             CL = 10.0f * _aircraft->_CLmax[0];
!             Vt = Vs;
              break;
          case 1:
--- 265,278 ----
          {
          case 0:
!             if (Vs > 0.5f)
!             {
!                 CL = 10.0f * _aircraft->_CLmax[0];
!                 Vt = Vs;
!             }
!             else	// No stall speed was specified
!             {
!                 CL = CL0;
!                 Vt = sqrtf(W/(0.5f*rho*CL*Sw));
!             }
              break;
          case 1:
***************
*** 288,297 ****
          _aircraft->_alpha.at(i) = alpha;
  
!         float CLmin = CL0 + MIN_ALPHA*CLalpha;
!         float CLmax = CL0 + std::max<float>(alpha,MAX_ALPHA)*CLalpha;
! 
!         _aircraft->_CYp.at(i) = -CL*CYp_const;
          _aircraft->_Cnp.at(i) = -CL/8.0f;
  
          float Cmin, Cmax;
  
--- 296,312 ----
          _aircraft->_alpha.at(i) = alpha;
  
! #if 0
!         float CYp_const = (AR+cosf(sweep))/(AR+4*cosf(sweep))*tanf(sweep);
!         _aircraft->_CYp.at(i) = CL*CYp_const;
          _aircraft->_Cnp.at(i) = -CL/8.0f;
+ #else
+         float l0a = 4*Sw*CL/PI/bw;
+         float av = l0a*(PI/4.0f)*bv/Sv;
+         _aircraft->_CYp.at(i) = 8.0f/(3.0f*PI)*nv*(bv*Sv/(bw*Sw))*av;
+         _aircraft->_Cnp.at(i) = -lv*_aircraft->_CYp[i]/bw;
+ #endif
  
+         float CLmin = CL0 + MIN_ALPHA*CLalpha;
+         float CLmax = CL0 + std::max<float>(alpha,MAX_ALPHA)*CLalpha;
          float Cmin, Cmax;
  
***************
*** 304,312 ****
        _aircraft->_Clbeta.at(i*2+1) = Cmax - Clbwf - alpha*Clbvt;
  
!         float Clr_const = 2.0f*lv*lv/bw/bw*(CYbeta*alpha);
          Cmin = (CLmin/4.0f)-Clr_const;
          Cmax = (CLmax/4.0f)-Clr_const;
          _aircraft->_Clr.at(i*2) = Cmin;
          _aircraft->_Clr.at(i*2+1) = Cmax;
      }
  
--- 319,329 ----
        _aircraft->_Clbeta.at(i*2+1) = Cmax - Clbwf - alpha*Clbvt;
  
!         float Clr_const = 2.0f*lv*zv/bw/bw*CYbeta;
          Cmin = (CLmin/4.0f)-Clr_const;
          Cmax = (CLmax/4.0f)-Clr_const;
          _aircraft->_Clr.at(i*2) = Cmin;
          _aircraft->_Clr.at(i*2+1) = Cmax;
+ 
+         if (Vs <= 0.5f) break;
      }
  
***************
*** 314,336 ****
      _aircraft->_CLadot = _aircraft->_CLq*deda;
  
!     float CLhde = ((CLah[0]/PI)*(acosf(1.0f-2.0f*Ee)+2.0f*sqrtf(Ee-Ee*Ee)));
!     _aircraft->_CLde = (CLhde*Sh/Sw)*2.0f/PI;
  
      // pitch
      if (_aircraft->_user_wing_data > 0)
      {
! //      float dwf = L/_aircraft->_aero_rp[X];
! //      float Kf = 0.033f + 0.538f*dwf + 1.5f*dwf*dwf;
! //      float Cmfus = Kf*D*D*L/Sw/cbarw/CLaw[0];
!         _aircraft->_Cmalpha = CLaw[0]*(dcgx/cbarw) - Vh*CLah[0]*(1.0f-deda);
          _aircraft->_Cmq = -_aircraft->_CLq*(lh/cbarw);
          _aircraft->_Cmadot = -_aircraft->_CLadot*(lh/cbarw);
  
!         float Cmtde = CLah[0]/PI*(1.0f-Ee)*sqrtf(Ee-Ee*Ee);
          _aircraft->_Cmde = (Sh*ch/Sw/cbarw*Cmtde - lh*Sh*CLhde/cbarw/Sw);
      }
  
      // side
!     float Cltdr = (CLav[0]/PI)*(acosf(1.0f-2.0f*Er)+2.0f*sqrtf(Er-Er*Er));
      _aircraft->_CYbeta = CYbeta;
      _aircraft->_CYr = -2.0f*(lv/bw)*(CYbeta);
--- 331,356 ----
      _aircraft->_CLadot = _aircraft->_CLq*deda;
  
!     float CLhde = ((CLah[0]/PI)*(acosf(1.0f-2.0f*Ee)+2.0f*sqrtf(Ee*(1.0f-Ee))));
!     _aircraft->_CLde = Sh*CLhde/Sw; // *2.0f/PI;
  
      // pitch
      if (_aircraft->_user_wing_data > 0)
      {
!         float dcgx = -(cg_loc[X] - _aircraft->_aero_rp[X])*INCH_TO_FEET;
! 
!         float dwf = L/_aircraft->_aero_rp[X];
!         float Kf = 0.033f + 0.538f*dwf + 1.5f*dwf*dwf;
!         float Cmaf = -Kf*D*D*L/Sw/cbarw/CLaw[0];
! 
!         _aircraft->_Cmalpha = CLaw[0]*(dcgx/cbarw) - Vh*CLah[0]*(1.0f-deda) + Cmaf;
          _aircraft->_Cmq = -_aircraft->_CLq*(lh/cbarw);
          _aircraft->_Cmadot = -_aircraft->_CLadot*(lh/cbarw);
  
!         float Cmtde = CLah[0]/PI*(1.0f-Ee)*sqrtf(Ee*(1.0f-Ee));
          _aircraft->_Cmde = (Sh*ch/Sw/cbarw*Cmtde - lh*Sh*CLhde/cbarw/Sw);
      }
  
      // side
!     float Cltdr = (CLav[0]/PI)*(acosf(1.0f-2.0f*Er)+2.0f*sqrtf(Er*(1.0f-Er)));
      _aircraft->_CYbeta = CYbeta;
      _aircraft->_CYr = -2.0f*(lv/bw)*(CYbeta);
***************
*** 339,360 ****
      // roll
      float TRh = _aircraft->_htail.taper;
! //  float TRv = _aircraft->_vtail.taper;
      _aircraft->_Clp = -(CLaw[0]/12.0f)*(1.0f+3.0f*TR)/(1.0f+TR)
                        + (CLah[0]/12.0f)*(Sh/Sw)*(1.0f+3.0f*TRh)/(1.0f+TRh)
!                       + (CLav[0]/12.0f)*(Sv/Sw)*(1.0f+3.0f*TRh)/(1.0f+TRh);
  
      // yaw
      float k0 = 0.075f;
!     float k1 = 1.0f+TC; // 1.256f; // correction factor for wing thickness
!     _aircraft->_Cnbeta = nv*Vv*CLav[0]*(1.0-dsdB);
!     _aircraft->_Cnr = -(k0*CL*CL + k1*CD0w) - 2.0f*nv*Vv*CLav[0]*(lv/bw);
      _aircraft->_Cndr = -Vv*Cltdr;
- 
- #if 0
-  printf("Cma: %f, Cmadot: %f, Cmq: %f, Cmde: %f\n",  _aircraft->_Cmalpha, _aircraft->_Cmadot, _aircraft->_Cmq, _aircraft->_Cmde);
-  printf("CYbeta: %f, CYr: %f, CYp: %f, CYdr: %f\n", _aircraft->_CYbeta, _aircraft->_CYr, _aircraft->_CYp, _aircraft->_CYdr);
-  printf("Cnbeta: %f, Cnr: %f, Cnp: %f, Cndr: %f\n", _aircraft->_Cnbeta, _aircraft->_Cnr, _aircraft->_Cnp, _aircraft->_Cndr);
-  printf("Clbeta: %f, Clr: %f, Clp: %f\n", _aircraft->_Clbeta, _aircraft->_Clr, _aircraft->_Clp);
- #endif
  }
  
--- 359,373 ----
      // roll
      float TRh = _aircraft->_htail.taper;
!     float TRv = _aircraft->_vtail.taper;
      _aircraft->_Clp = -(CLaw[0]/12.0f)*(1.0f+3.0f*TR)/(1.0f+TR)
                        + (CLah[0]/12.0f)*(Sh/Sw)*(1.0f+3.0f*TRh)/(1.0f+TRh)
!                       + (CLav[0]/12.0f)*(Sv/Sw)*(1.0f+3.0f*TRv)/(1.0f+TRv);
  
      // yaw
      float k0 = 0.075f;
!     float k1 = 1.0f+TC;	// correction factor for wing thickness
!     _aircraft->_Cnbeta = nv*Vv*CLav[0]; // *(1.0f-dsdB);
!     _aircraft->_Cnr = -2.0f*nv*Vv*CLav[0]*(lv/bw) - (k0*CL*CL + k1*CD0w);
      _aircraft->_Cndr = -Vv*Cltdr;
  }
  
***************
*** 464,468 ****
  
      float AR = _aircraft->_wing.aspect;
- //  float sweep = _aircraft->_wing.sweep * DEG_TO_RAD;
      float Ew = _aircraft->_wing.efficiency;
  
--- 477,480 ----
***************
*** 482,487 ****
      file << "    <!-- Antennas, struts and wires are not taken into account        -->" << std::endl;
      file << "    <!-- CD for gear (fixed and retractable) is defined below         -->" << std::endl;
!     file << "    <function name=\"aero/force/Drag_zero_lift\">" << std::endl;
!     file << "       <description>Drag at zero lift</description>" << std::endl;
      file << "       <product>" << std::endl;
      file << "          <property>aero/qbar-psf</property>" << std::endl;
--- 494,499 ----
      file << "    <!-- Antennas, struts and wires are not taken into account        -->" << std::endl;
      file << "    <!-- CD for gear (fixed and retractable) is defined below         -->" << std::endl;
!     file << "    <function name=\"aero/force/Drag_minimum\">" << std::endl;
!     file << "       <description>Minimum drag</description>" << std::endl;
      file << "       <product>" << std::endl;
      file << "          <property>aero/qbar-psf</property>" << std::endl;
***************
*** 785,800 ****
      Cnr = _aircraft->_Cnr;
  
!     file << std::setprecision(4) << std::fixed << std::showpoint;
!     file << "    <!-- Stall initiator -->" << std::endl;
!     file << "    <function name=\"aero/moment/Yaw_alpha\">" << std::endl;
!     file << "       <description>Yaw moment due to alpha</description>" << std::endl;
!     file << "       <product>" << std::endl;
!     file << "           <property>aero/qbar-psf</property>" << std::endl;
!     file << "           <property>metrics/Sw-sqft</property>" << std::endl;
!     file << "           <property>metrics/bw-ft</property>" << std::endl;
!     file << "           <property>aero/alpha-rad</property>" << std::endl;
!     file << "           <property>fcs/gear-no-wow</property>" << std::endl;
  
!     file << "           <table>" << std::endl;
          file << "             <independentVar lookup=\"row\">aero/beta-rad</independentVar>" << std::endl;
          file << "             <independentVar lookup=\"column\">aero/Re</independentVar>" << std::endl;
--- 797,814 ----
      Cnr = _aircraft->_Cnr;
  
!     if (_aircraft->_Re[1] != 0.0f)
!     {
!         file << std::setprecision(4) << std::fixed << std::showpoint;
!         file << "    <!-- Stall initiator -->" << std::endl;
!         file << "    <function name=\"aero/moment/Yaw_alpha\">" << std::endl;
!         file << "       <description>Yaw moment due to alpha</description>" << std::endl;
!         file << "       <product>" << std::endl;
!         file << "           <property>aero/qbar-psf</property>" << std::endl;
!         file << "           <property>metrics/Sw-sqft</property>" << std::endl;
!         file << "           <property>metrics/bw-ft</property>" << std::endl;
!         file << "           <property>aero/alpha-rad</property>" << std::endl;
!         file << "           <property>fcs/gear-no-wow</property>" << std::endl;
  
!         file << "           <table>" << std::endl;
          file << "             <independentVar lookup=\"row\">aero/beta-rad</independentVar>" << std::endl;
          file << "             <independentVar lookup=\"column\">aero/Re</independentVar>" << std::endl;
***************
*** 817,824 ****
          file << "             </tableData>" << std::endl;
          file << "           </table>" << std::endl;
! 
!     file << "       </product>" << std::endl;
!     file << "    </function>" << std::endl;
!     file << std::endl;
      file << "    <function name=\"aero/moment/Yaw_beta\">" << std::endl;
      file << "       <description>Yaw moment due to beta</description>" << std::endl;
--- 831,838 ----
          file << "             </tableData>" << std::endl;
          file << "           </table>" << std::endl;
!         file << "       </product>" << std::endl;
!         file << "    </function>" << std::endl;
!         file << std::endl;
!     }
      file << "    <function name=\"aero/moment/Yaw_beta\">" << std::endl;
      file << "       <description>Yaw moment due to beta</description>" << std::endl;
***************
*** 1184,1188 ****
  
      file << std::fixed << std::showpoint;
!     if (C.size() == 1) {
          file << "           <value> " << (C[0]) << " </value>" << std::endl;
      }
--- 1198,1202 ----
  
      file << std::fixed << std::showpoint;
!     if (C.size() == 1 || _aircraft->_Re[1] == 0.0f) {
          file << "           <value> " << (C[0]) << " </value>" << std::endl;
      }

Thread: [Jsbsim-cvslogs] CVS: JSBSim/utils/aeromatic++/Systems Controls.cpp, 1.48, 1.49

An open source flight dynamics & control software library in C++

jsbsim-cvslogs