test/src/test_mc.cpp

   1 /* Copyright_License {
   2
   3   XCSoar Glide Computer - http://www.xcsoar.org/
   4   Copyright (C) 2000-2013 The XCSoar Project
   5   A detailed list of copyright holders can be found in the file "AUTHORS".
   6
   7   This program is free software; you can redistribute it and/or
   8   modify it under the terms of the GNU General Public License
   9   as published by the Free Software Foundation; either version 2
  10   of the License, or (at your option) any later version.
  11
  12   This program is distributed in the hope that it will be useful,
  13   but WITHOUT ANY WARRANTY; without even the implied warranty of
  14   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  15   GNU General Public License for more details.
  16
  17   You should have received a copy of the GNU General Public License
  18   along with this program; if not, write to the Free Software
  19   Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
  20 }
  21 */
  22
  23 #include "test_debug.hpp"
  24 #include "GlideSolvers/GlideSettings.hpp"
  25 #include "GlideSolvers/GlidePolar.hpp"
  26 #include "GlideSolvers/GlideState.hpp"
  27 #include "GlideSolvers/GlideResult.hpp"
  28 #include "GlideSolvers/MacCready.hpp"
  29 #include "Navigation/Aircraft.hpp"
  30 #include "OS/FileUtil.hpp"
  31
  32 #include <stdio.h>
  33 #include <fstream>
  34 #include <string>
  35 #include <math.h>
  36
  37 const fixed Vmin(5.0);
  38
  39 static void
  40 polar_mc(std::ofstream &ofile, const fixed mc)
  41 {
  42   GlidePolar polar(mc);
  43   ofile << (double)mc << " "
  44         << (double)polar.GetVBestLD() << " "
  45         << (double)polar.GetBestLD() << " "
  46         << (double)polar.GetVMin() << " "
  47         << (double)polar.GetSMin() << " "
  48         << (double)polar.GetVMax() << " "
  49         << (double)polar.GetSMax() << "\n";
  50 }
  51
  52 static void
  53 basic_polar(const fixed mc)
  54 {
  55   char bname[100];
  56   sprintf(bname,"output/results/res-polar-%02d-best.txt",(int)(mc*10));
  57   std::ofstream pfile("output/results/res-polar.txt");
  58   std::ofstream mfile(bname);
  59
  60   GlidePolar polar(mc);
  61   for (fixed V= Vmin; V<= polar.GetVMax(); V+= fixed(0.25)) {
  62     pfile << (double)mc << " "
  63           << (double)V << " "
  64           << -(double)polar.SinkRate(V) << " "
  65           << (double)(V/polar.SinkRate(V))
  66           << "\n";
  67   }
  68
  69   mfile << (double)mc
  70         << " " << 0
  71         << " " << (double)mc
  72         << " " << (double)polar.GetBestLD()
  73         << "\n";
  74   mfile << (double)mc
  75         << " " << (double)polar.GetVBestLD()
  76         << " " << -(double)polar.GetSBestLD()
  77         << " " << (double)polar.GetBestLD()
  78         << "\n";
  79 }
  80
  81
  82 static void
  83 test_glide_alt(const fixed h, const fixed W, const fixed Wangle,
  84                std::ostream &hfile)
  85 {
  86   GlideSettings settings;
  87   settings.SetDefaults();
  88
  89   GlidePolar polar(fixed(0));
  90   polar.SetMC(fixed(1));
  91
  92   AircraftState ac;
  93   ac.wind.norm = fabs(W);
  94   if (negative(W)) {
  95     ac.wind.bearing = Angle::Degrees(fixed(180)+Wangle);
  96   } else {
  97     ac.wind.bearing = Angle::Degrees(Wangle);
  98   }
  99   ac.altitude = h;
 100
 101   GeoVector vect(fixed(400.0), Angle::Zero());
 102   GlideState gs(vect, fixed(0), ac.altitude, ac.wind);
 103   GlideResult gr = MacCready::Solve(settings, polar, gs);
 104   hfile << (double)h << " "
 105         << (double)gr.altitude_difference << " "
 106         << (double)gr.time_elapsed << " "
 107         << (double)gr.v_opt << " "
 108         << (double)W << " "
 109         << (double)Wangle << " "
 110         << "\n";
 111 }
 112
 113 static void
 114 test_glide_stf(const fixed h, const fixed W, const fixed Wangle, const fixed S,
 115                std::ostream &hfile)
 116 {
 117   GlideSettings settings;
 118   settings.SetDefaults();
 119
 120   GlidePolar polar(fixed(0));
 121   polar.SetMC(fixed(1));
 122
 123   AircraftState ac;
 124   ac.wind.norm = fabs(W);
 125   if (negative(W)) {
 126     ac.wind.bearing = Angle::Degrees(fixed(180)+Wangle);
 127   } else {
 128     ac.wind.bearing = Angle::Degrees(Wangle);
 129   }
 130   ac.altitude = h;
 131   ac.netto_vario = S;
 132
 133   GeoVector vect(fixed(400.0), Angle::Zero());
 134   GlideState gs(vect, fixed(0), ac.altitude, ac.wind);
 135   GlideResult gr = MacCready::Solve(settings, polar, gs);
 136
 137   fixed Vstf = polar.SpeedToFly(ac, gr, false);
 138
 139   hfile << (double)h << " "
 140         << (double)gr.altitude_difference << " "
 141         << (double)gr.v_opt << " "
 142         << (double)Vstf << " "
 143         << (double)W << " "
 144         << (double)Wangle << " "
 145         << (double)ac.netto_vario << " "
 146         << "\n";
 147 }
 148
 149 static bool
 150 test_stf()
 151 {
 152   { // variation with height
 153     std::ofstream hfile("output/results/res-polar-s0.txt");
 154     for (fixed h=fixed(0); h<fixed(40.0); h+= fixed(0.1)) {
 155       test_glide_stf(h,fixed(0),fixed(0),fixed(0),hfile);
 156     }
 157   }
 158   { // variation with S, below FG
 159     std::ofstream hfile("output/results/res-polar-s1.txt");
 160     for (fixed S=fixed(-4.0); S<fixed(4.0); S+= fixed(0.1)) {
 161       test_glide_stf(fixed(0), fixed(0),fixed(0),S, hfile);
 162     }
 163   }
 164   { // variation with S, above FG
 165     std::ofstream hfile("output/results/res-polar-s2.txt");
 166     for (fixed S=fixed(-4.0); S<fixed(4.0); S+= fixed(0.1)) {
 167       test_glide_stf(fixed(40), fixed(0),fixed(0),S, hfile);
 168     }
 169   }
 170   { // variation with S, below FG, wind
 171     std::ofstream hfile("output/results/res-polar-s3.txt");
 172     for (fixed S=fixed(-4.0); S<fixed(4.0); S+= fixed(0.1)) {
 173       test_glide_stf(fixed(0), fixed(10.0), fixed(0),S, hfile);
 174     }
 175   }
 176   { // variation with S, above FG, wind
 177     std::ofstream hfile("output/results/res-polar-s4.txt");
 178     for (fixed S=fixed(-4.0); S<fixed(4.0); S+= fixed(0.1)) {
 179       test_glide_stf(fixed(40), fixed(10.0), fixed(0), S, hfile);
 180     }
 181   }
 182   return true;
 183 }
 184
 185 static bool
 186 test_mc()
 187 {
 188   {
 189     std::ofstream ofile("output/results/res-polar-m.txt");
 190     for (fixed mc=fixed(0); mc<fixed(5.0); mc+= fixed(0.1)) {
 191       basic_polar(mc);
 192       polar_mc(ofile, mc);
 193     }
 194   }
 195
 196   {
 197     std::ofstream hfile("output/results/res-polar-h-00.txt");
 198     for (fixed h=fixed(0); h<fixed(40.0); h+= fixed(0.1)) {
 199       test_glide_alt(h, fixed(0), fixed(0), hfile);
 200     }
 201   }
 202
 203   {
 204     std::ofstream hfile("output/results/res-polar-h-50.txt");
 205     for (fixed h=fixed(0); h<fixed(40.0); h+= fixed(0.1)) {
 206       test_glide_alt(h, fixed(5.0), fixed(0), hfile);
 207     }
 208   }
 209
 210   {
 211     std::ofstream hfile("output/results/res-polar-w.txt");
 212     for (fixed w=fixed(-10.0); w<=fixed(10.0); w+= fixed(0.1)) {
 213       test_glide_alt(fixed(50.0), w, fixed(0), hfile);
 214     }
 215   }
 216
 217   {
 218     std::ofstream hfile("output/results/res-polar-a.txt");
 219     for (fixed a=fixed(0); a<=fixed(360.0); a+= fixed(10)) {
 220       test_glide_alt(fixed(50.0), fixed(10.0), a, hfile);
 221     }
 222   }
 223   return true;
 224 }
 225
 226 static void
 227 test_glide_cb(const fixed h, const fixed W, const fixed Wangle,
 228               std::ostream &hfile)
 229 {
 230   GlideSettings settings;
 231   settings.SetDefaults();
 232
 233   GlidePolar polar(fixed(1));
 234
 235   AircraftState ac;
 236   ac.wind.norm = fabs(W);
 237   if (negative(W)) {
 238     ac.wind.bearing = Angle::Degrees(fixed(180)+Wangle);
 239   } else {
 240     ac.wind.bearing = Angle::Degrees(Wangle);
 241   }
 242   ac.altitude = h;
 243
 244   GeoVector vect(fixed(400.0), Angle::Zero());
 245   GlideState gs (vect, fixed(0), ac.altitude, ac.wind);
 246   GlideResult gr = MacCready::Solve(settings, polar, gs);
 247
 248   gr.CalcDeferred();
 249
 250   hfile << (double)W << " "
 251         << (double)Wangle << " "
 252         << (double)gr.vector.bearing.Degrees() << " "
 253         << (double)gr.cruise_track_bearing.Degrees() << " "
 254         << "\n";
 255 }
 256
 257 static bool
 258 test_cb()
 259 {
 260   {
 261     std::ofstream hfile("output/results/res-polar-cb.txt");
 262     for (fixed a = fixed(0); a <= fixed(360.0); a+= fixed(10)) {
 263       test_glide_cb(fixed(0), fixed(10.0), a, hfile);
 264     }
 265   }
 266   return true;
 267 }
 268
 269 int main() {
 270
 271   plan_tests(3);
 272
 273   Directory::Create(_T("output/results"));
 274
 275   ok(test_mc(),"mc output",0);
 276   ok(test_stf(),"mc stf",0);
 277   ok(test_cb(),"cruise bearing",0);
 278
 279   return exit_status();
 280
 281 }