home *** CD-ROM | disk | FTP | other *** search
/ io Programmo 41 / IOPROG_41.ISO / soft / c++ / NUMCPP11.ZIP / test.cpp < prev    next >
Encoding:
C/C++ Source or Header  |  1999-05-17  |  15.9 KB  |  513 lines
  1. //===================================================================
  2. // test.cpp
  3. //
  4. // Version 1.
  5. //
  6. // Written by:
  7. //   Brent Worden
  8. //   WordenWare
  9. //   email:  Brent@Worden.org
  10. //
  11. // Copyright (c) 1998-1999 WordenWare
  12. //
  13. // Created:  May 17, 1999
  14. // Revised:  
  15. //===================================================================
  16.  
  17. #include <iostream>
  18. #include <string>
  19. #include <list>
  20.  
  21. #include "algorthm.hpp"
  22. #include "contdist.h"
  23. #include "correlat.hpp"
  24. #include "descript.hpp"
  25. #include "deviate.h"
  26. #include "discdist.h"
  27. #include "eigensys.hpp"
  28. #include "linalg.hpp"
  29. #include "matrix.hpp"
  30. #include "rootfind.hpp"
  31. #include "vector.hpp"
  32.  
  33. static const int SIZE = 5;
  34.  
  35. int randomValue(int n)
  36. {
  37.     return (int)((double)rand() / ((double)RAND_MAX + 1) * n);
  38. }
  39.  
  40. template<class T>
  41. void initArray(T first, T last)
  42. {
  43.     while(first != last){
  44.         *first = randomValue(SIZE);
  45.         ++first;
  46.     }
  47. }
  48. template<class T>
  49. void outputArray(T array[], int n)
  50. {
  51.     std::cout << array[0];
  52.     for(int i = 1; i < n; ++i){
  53.         std::cout << ' ' << array[i];
  54.     }
  55.     std::cout << std::endl;
  56. }
  57.  
  58. template<class T>
  59. std::ostream& operator<<(std::ostream& os, num::Vector<T>& vec)
  60. {
  61.     for(int i = 0; i < vec.size(); ++i){
  62.         os << vec[i] << ' ';
  63.     }
  64.     os << std::endl;
  65.     return os;
  66. }
  67.  
  68. template<class T>
  69. std::ostream& operator<<(std::ostream& os, num::Matrix<T>& mat)
  70. {
  71.     for(int i = 0; i < mat.rows(); ++i){
  72.         os << mat[i];
  73.     }
  74.     return os;
  75. }
  76.  
  77. void testAlgorithm()
  78. {
  79.     double x[SIZE];
  80.     num::Vector<double> y(SIZE);
  81.  
  82.     initArray(x, x + SIZE);
  83.     initArray(y.begin(), y.end());
  84.  
  85.     std::cout << "x: "; outputArray(x, SIZE);
  86.     std::cout << "y: " << y;
  87.  
  88.     std::cout << "average(x[0], x[1]): " << num::average(x[0], x[1]) << std::endl << std::endl;
  89.  
  90.     double pair[SIZE*SIZE];
  91.     num::pairwdiff(x, x + SIZE, y.begin(), y.end(), pair);
  92.     std::cout << "pairwdiff(x, y): "; outputArray(pair, SIZE*SIZE); std::cout << std::endl;
  93.     
  94.     double walsh[SIZE*(SIZE+1)/2];
  95.     num::walshavg(x, x + SIZE, walsh);
  96.     std::cout << "walshavg(x): "; outputArray(walsh, SIZE*(SIZE+1)/2); std::cout << std::endl;
  97. }
  98.  
  99. void testContDist()
  100. {
  101.     double p;
  102.     double x;
  103.  
  104.     p = num::betap(.5, 1, 2); x = num::betav(p, 1, 2);
  105.     std::cout << "betap(.5, 1, 2):     " << p << "\t betav(" << p << ", 1, 2):     " << x << std::endl;
  106.  
  107.     p = num::cauchyp(.5, 1, 2); x = num::cauchyv(p, 1, 2);
  108.     std::cout << "cauchyp(.5, 1, 2):   " << p << "\t cauchyv(" << p << ", 1, 2):   " << x << std::endl;
  109.  
  110.     p = num::dblexpp(.5, 1, 2); x = num::dblexpv(p, 1, 2);
  111.     std::cout << "dblexpp(.5, 1, 2):   " << p << "\t dblexpv(" << p << ", 1, 2):   " << x << std::endl;
  112.  
  113.     p = num::expp(.5, 1); x = num::expv(p, 1);
  114.     std::cout << "expp(.5, 1):         " << p << "\t expv(" << p << ", 1):         " << x << std::endl;
  115.  
  116.     p = num::gammap(.5, 1, 2); x = num::gammav(p, 1, 2);
  117.     std::cout << "gammap(.5, 1, 2):    " << p << "\t gammav(" << p << ", 1, 2):    " << x << std::endl;
  118.  
  119.     p = num::logisticp(.5, 1, 2); x = num::logisticv(p, 1, 2);
  120.     std::cout << "logisticp(.5, 1, 2): " << p << "\t logisticv(" << p << ", 1, 2): " << x << std::endl;
  121.  
  122.     p = num::lognormp(.5, 1, 2); x = num::lognormv(p, 1, 2);
  123.     std::cout << "lognormp(.5, 1, 2):  " << p << "\t lognormv(" << p << ", 1, 2):  " << x << std::endl;
  124.  
  125.     p = num::rayleighp(.5, 1); x = num::rayleighv(p, 1);
  126.     std::cout << "rayleighp(.5, 1):    " << p << "\t rayleighv(" << p << ", 1):    " << x << std::endl;
  127.  
  128.     p = num::uniformp(1.5, 1, 2); x = num::uniformv(p, 1, 2);
  129.     std::cout << "uniformp(1.5, 1, 2): " << p << "\t uniformv(" << p << ", 1, 2):  " << x << std::endl;
  130.  
  131.     p = num::weibullp(.5, 1, 2); x = num::weibullv(p, 1, 2);
  132.     std::cout << "weibullp(.5, 1, 2):  " << p << "\t weibullv(" << p << ", 1, 2):  " << x << std::endl;
  133. }
  134.  
  135. void testCorrelation()
  136. {
  137.     double x[SIZE];
  138.     num::Vector<double> y(SIZE);
  139.  
  140.     initArray(x, x + SIZE);
  141.     initArray(y.begin(), y.end());
  142.  
  143.     std::cout << "x: "; outputArray(x, SIZE);
  144.     std::cout << "y: " << y;
  145.  
  146.     std::cout << "covariance(x, y):  " << num::covariance(x, x + SIZE, y.begin()) << std::endl;
  147.     std::cout << "correlation(x, y): " << num::correlation(y.begin(), y.end(), x) << std::endl;
  148.     std::cout << "kendall(x, y):     " << num::kendall(x, x + SIZE, y.begin()) << std::endl;
  149.     std::cout << "spearman(x, y):    " << num::spearman(y.begin(), y.end(), x) << std::endl;
  150. }
  151.  
  152. void testDescript()
  153. {
  154.     double x[SIZE];
  155.     num::Vector<double> y(SIZE);
  156.     std::list<double> z(SIZE);
  157.  
  158.     initArray(x, x + SIZE);
  159.     initArray(y.begin(), y.end());
  160.     initArray(z.begin(), z.end());
  161.  
  162.     std::cout << "centmon(x, 3):        " << num::centmom(x, x + SIZE, 3) << std::endl;
  163.     std::cout << "coeffvar(y):          " << num::coeffvar(y.begin(), y.end()) << std::endl;
  164.     std::cout << "geomean(z):           " << num::geomean(z.begin(), z.end()) << std::endl;
  165.     std::cout << "harmean(y):           " << num::harmean(y.begin(), y.end()) << std::endl;
  166.     std::cout << "iqr(x):               " << num::iqr(x, x + SIZE) << std::endl;
  167.     std::cout << "kurtosis(y):          " << num::kurtosis(y.begin(), y.end()) << std::endl;
  168.     std::cout << "mean(z):              " << num::mean(z.begin(), z.end()) << std::endl;
  169.     std::cout << "meandev(y):           " << num::meandev(y.begin(), y.end()) << std::endl;
  170.     std::cout << "meddev(x):            " << num::meddev(x, x + SIZE) << std::endl;
  171.     std::cout << "medain(y):            " << num::median(y.begin(), y.end()) << std::endl;
  172.     std::cout << "quantile(z):          " << num::quantile(z.begin(), z.end(), .1) << std::endl;
  173.     std::cout << "quartile1(y):         " << num::quartile1(y.begin(), y.end()) << std::endl;
  174.     std::cout << "quartile3(x):         " << num::quartile3(x, x + SIZE) << std::endl;
  175.     std::cout << "range(y):             " << num::range(y.begin(), y.end()) << std::endl;
  176.     std::cout << "rms(z):               " << num::rms(z.begin(), z.end()) << std::endl;
  177.     std::cout << "skewness(y):          " << num::skewness(y.begin(), y.end()) << std::endl;
  178.     std::cout << "stddev(x):            " << num::stddev(x, x + SIZE) << std::endl;
  179.     std::cout << "stderrmean(y):        " << num::stderrmean(y.begin(), y.end()) << std::endl;
  180.     std::cout << "sum(x):               " << num::sum(x, x + SIZE) << std::endl;
  181.     std::cout << "sum2(y):              " << num::sum2(y.begin(), y.end()) << std::endl;
  182.     std::cout << "sump(x):              " << num::sump(x, x + SIZE, y.begin()) << std::endl;
  183.     std::cout << "sumsquare(y):         " << num::sumsquare(y.begin(), y.end()) << std::endl;
  184.     std::cout << "trimmean(x, .05):     " << num::trimmean(x, x + SIZE, .05) << std::endl;
  185.     std::cout << "trimmean(y, .05, .1): " << num::trimmean(y.begin(), y.end(), .05, .1) << std::endl;
  186.     std::cout << "variance(z):          " << num::variance(z.begin(), z.end()) << std::endl;
  187. }
  188.  
  189. void testDeviate()
  190. {
  191.     int n = SIZE;
  192.     int i;
  193.     long seed = 10203040L;
  194.  
  195.     std::cout << "ran0: ";
  196.     for(i = 0; i < SIZE; ++i){ std::cout << num::ran0(&seed) << ' '; } std::cout << std::endl;
  197.     
  198.     std::cout << "ran1: ";
  199.     for(i = 0; i < SIZE; ++i){ std::cout << num::ran1(&seed) << ' '; } std::cout << std::endl;
  200.     
  201.     std::cout << "ran2: ";
  202.     for(i = 0; i < SIZE; ++i){ std::cout << num::ran2(&seed) << ' '; } std::cout << std::endl;
  203.     
  204.     std::cout << "ran3: ";
  205.     for(i = 0; i < SIZE; ++i){ std::cout << num::ran3(&seed) << ' '; } std::cout << std::endl;
  206.     
  207.     std::cout << "berdev(.5): ";
  208.     for(i = 0; i < SIZE; ++i){ std::cout << num::berdev(.5, &seed) << ' '; } std::cout << std::endl;
  209.     
  210.     std::cout << "betadev(1, 2): ";
  211.     for(i = 0; i < SIZE; ++i){ std::cout << num::betadev(1, 2, &seed) << ' '; } std::cout << std::endl;
  212.     
  213.     std::cout << "bnldev(10, .5): ";
  214.     for(i = 0; i < SIZE; ++i){ std::cout << num::bnldev(10, .5, &seed) << ' '; } std::cout << std::endl;
  215.     
  216.     std::cout << "caudev(1, 2): ";
  217.     for(i = 0; i < SIZE; ++i){ std::cout << num::caudev(1, 2, &seed) << ' '; } std::cout << std::endl;
  218.     
  219.     std::cout << "chidev(1): ";
  220.     for(i = 0; i < SIZE; ++i){ std::cout << num::chidev(1, &seed) << ' '; } std::cout << std::endl;
  221.     
  222.     std::cout << "dexpdev(1, 2): ";
  223.     for(i = 0; i < SIZE; ++i){ std::cout << num::dexpdev(1, 2, &seed) << ' '; } std::cout << std::endl;
  224.     
  225.     std::cout << "expdev(1): ";
  226.     for(i = 0; i < SIZE; ++i){ std::cout << num::expdev(1, &seed) << ' '; } std::cout << std::endl;
  227.     
  228.     std::cout << "fdev(1, 2): ";
  229.     for(i = 0; i < SIZE; ++i){ std::cout << num::fdev(1, 2, &seed) << ' '; } std::cout << std::endl;
  230.     
  231.     std::cout << "gamdev(1, 2): ";
  232.     for(i = 0; i < SIZE; ++i){ std::cout << num::gamdev(1, 2, &seed) << ' '; } std::cout << std::endl;
  233.     
  234.     std::cout << "gasdev(1, 2): ";
  235.     for(i = 0; i < SIZE; ++i){ std::cout << num::gasdev(1, 2, &seed) << ' '; } std::cout << std::endl;
  236.     
  237.     std::cout << "logdev(1, 2): ";
  238.     for(i = 0; i < SIZE; ++i){ std::cout << num::logdev(1, 2, &seed) << ' '; } std::cout << std::endl;
  239.     
  240.     std::cout << "nchidev(1, 2): ";
  241.     for(i = 0; i < SIZE; ++i){ std::cout << num::nchidev(1, 2, &seed) << ' '; } std::cout << std::endl;
  242.     
  243.     std::cout << "nfdev(1, 2, 3): ";
  244.     for(i = 0; i < SIZE; ++i){ std::cout << num::nfdev(1, 2, 3, &seed) << ' '; } std::cout << std::endl;
  245.     
  246.     std::cout << "poisdev(1): ";
  247.     for(i = 0; i < SIZE; ++i){ std::cout << num::poisdev(1, &seed) << ' '; } std::cout << std::endl;
  248.     
  249.     std::cout << "tdev(1): ";
  250.     for(i = 0; i < SIZE; ++i){ std::cout << num::tdev(1, &seed) << ' '; } std::cout << std::endl;
  251.     
  252.     std::cout << "unifdev(1, 2): ";
  253.     for(i = 0; i < SIZE; ++i){ std::cout << num::unifdev(1, 2, &seed) << ' '; } std::cout << std::endl;
  254.     
  255.     std::cout << "weibdev(1, 2): ";
  256.     for(i = 0; i < SIZE; ++i){ std::cout << num::weibdev(1, 2, &seed) << ' '; } std::cout << std::endl;
  257. }
  258.  
  259. void testDiscDist()
  260. {
  261.     double p;
  262.     int x1, x2;
  263.  
  264.     p = num::binomialp(5, 10, .5); num::binomialv(p, 10, .5, &x1, &x2);
  265.     std::cout << "binomialp(5, 10, .5): " << p <<
  266.         "\t binomialv(" << p << ", 10, .5): " << x1 << ", " << x2 << std::endl;
  267.  
  268.     p = num::geomp(2, .5); num::geomv(p, .5, &x1, &x2);
  269.     std::cout << "geomp(2, .5):     " << p <<
  270.         "\t geomv(" << p << ", .5):     " << x1 << ", " << x2 << std::endl;
  271.  
  272.     p = num::hyperp(5, 10, 20, 30); num::hyperv(p, 10, 20, 30, &x1, &x2);
  273.     std::cout << "hyperp(5, 10, 20, 30):    " << p <<
  274.         "\t hyperv(" << p << ", 10, 20, 30):    " << x1 << ", " << x2 << std::endl;
  275.  
  276.     p = num::negbnlp(5, 10, .5); num::negbnlv(p, 10, .5, &x1, &x2);
  277.     std::cout << "negbnlp(5, 10, .5):   " << p <<
  278.         "\t negbnlv(" << p << ", 10, .5):   " << x1 << ", " << x2 << std::endl;
  279.  
  280.     p = num::poissonp(5, 10); num::poissonv(p, 10, &x1, &x2);
  281.     std::cout << "poissonp(5, 10):  " << p <<
  282.         "\t poissonv(" << p << ", 10):  " << x1 << ", " << x2 << std::endl;
  283. }
  284.  
  285. void testEigenSys()
  286. {
  287.     num::Matrix<double> a(SIZE, SIZE);
  288.     num::Vector<double> d;
  289.     num::Vector<double> e;
  290.     int i, j;
  291.     double aij;
  292.  
  293.     for(i = 0; i < SIZE; ++i){
  294.         for(j = i + 1; j < SIZE; ++j){
  295.             aij = randomValue(SIZE);
  296.             a[i][j] = aij;
  297.             a[j][i] = aij;
  298.         }
  299.         a[i][i] = randomValue(SIZE);
  300.     }
  301.  
  302.     std::cout << "a: " << std::endl << a << std::endl;
  303.     num::tridred(a, d, e);
  304.     std::cout << "tridred(a):" << std::endl;
  305.     std::cout << "Diagonal: " << d;
  306.     std::cout << "Off Diag: " << e << std::endl;
  307.  
  308.     num::trideig(d, e);
  309.     std::cout << "trideig(d, e):" << std::endl;
  310.     std::cout << "Eigenvalues: " << d;
  311. }
  312.  
  313. void testLinAlg()
  314. {
  315.     num::Matrix<double> a(SIZE, SIZE);
  316.     num::Vector<double> p(SIZE);
  317.     num::Matrix<double> aa(SIZE, SIZE);
  318.     int aij;
  319.  
  320.     for(int i = 0; i < SIZE; ++i){
  321.         for(int j = i + 1; j < SIZE; ++j){
  322.             aij = randomValue(SIZE);
  323.             aa[i][j] = aij;
  324.             aa[j][i] = aij;
  325.         }
  326.         aa[i][i] = randomValue(SIZE);
  327.         p[i] = randomValue(SIZE);
  328.     }
  329.  
  330.     a = aa;
  331.     std::cout << "a: " << std::endl << a;
  332.     std::cout << "determinant(a): " << num::determinant(a) << std::endl;
  333.  
  334.     std::cout << "inverse(a):" << std::endl;
  335.     num::inverse(a);
  336.     std::cout << a;
  337.  
  338.     a = aa;
  339.     std::cout << "linsolve(a, b): " << std::endl;
  340.     std::cout << "b: " << p;
  341.     num::linsolve(a, p, p);
  342.     std::cout << "x: " << p;
  343. }
  344.  
  345. void testRanking()
  346. {
  347.     num::Vector<double> x(SIZE);
  348.     num::Vector<double> rank(SIZE);
  349.     num::Vector<double> cp(SIZE);
  350.  
  351.     initArray(cp.begin(), cp.end());
  352.  
  353.     x = cp;
  354.     std::cout << "x: " << x << std::endl;
  355.  
  356.     num::index(x.begin(), x.end(), rank.begin());
  357.     std::cout << "index(x): " << rank;
  358.  
  359.     x = cp;
  360.     num::rank(x.begin(), x.end(), rank.begin());
  361.     std::cout << "rank(x): " << rank;
  362.  
  363.     x = cp;
  364.     num::arank(x.begin(), x.end(), rank.begin());
  365.     std::cout << "arank(x): " << rank;
  366.  
  367.     x = cp;
  368.     num::abrank(x.begin(), x.end(), rank.begin());
  369.     std::cout << "abrank(x): " << rank;
  370.  
  371.     std::cout << "countGreater(x, 5): " << num::countGreater(x.begin(), x.end(), 5) << std::endl;
  372. }
  373.  
  374. void testResidual()
  375. {
  376.     num::Vector<double> x(SIZE);
  377.     num::Vector<double> cp(SIZE);
  378.  
  379.     initArray(cp.begin(), cp.end());
  380.  
  381.     x = cp;
  382.     std::cout << "x: " << x << std::endl;
  383.  
  384.     num::resid(x.begin(), x.end(), 3.0);
  385.     std::cout << "resid(x, 3): " << x;
  386.  
  387.     x = cp;
  388.     num::aresid(x.begin(), x.end(), 3.0);
  389.     std::cout << "aresid(x, 3): " << x;
  390.  
  391.     x = cp;
  392.     num::nresid(x.begin(), x.end(), 5, 1.5);
  393.     std::cout << "nresid(x, 5, 1.5): " << x;
  394. }
  395.  
  396. double function1(double x)
  397. {
  398.     return sin(x) * cos(x) / exp(x);
  399. }
  400.  
  401. double function2(double x)
  402. {
  403.     return sin(x) - cos(x) / exp(x);
  404. }
  405.  
  406. double function2Prime(double x)
  407. {
  408.     return cos(x) + (cos(x) + sin(x)) / exp(x);
  409. }
  410.  
  411. void testRootfind()
  412. {
  413.     double x = num::bisection(function1, 1.0, 2.0);
  414.     double y = function1(x);
  415.  
  416.     std::cout << "bisection(function1, 1, 2): " << x << std::endl;
  417.     std::cout << "function1(" << x << "): " << y << std::endl << std::endl;
  418.  
  419.     x = num::newton(function2, function2Prime, 3.0);
  420.     y = function2(x);
  421.     std::cout << "newton(function2, function2Prime, 3): " << x << std::endl;
  422.     std::cout << "function2(" << x << "): " << y << std::endl << std::endl;
  423.  
  424.     x = num::secant(function1, 3.0, 4.0);
  425.     std::cout << "secant(function1, 3, 4): " << x << std::endl;
  426.     std::cout << "function1(" << x << "): " << y << std::endl;
  427. }
  428.  
  429. void testSort()
  430. {
  431.     num::Vector<double> x(SIZE * 10);
  432.     num::Vector<double> cpx(SIZE * 10);
  433.     num::Vector<int> y(SIZE * 10);
  434.     num::Vector<int> cpy(SIZE * 10);
  435.  
  436.     initArray(cpx.begin(), cpx.end());
  437.     initArray(cpy.begin(), cpy.end());
  438.  
  439.     x = cpx;
  440.     y = cpy;
  441.     std::cout << "x: " << x << "y: " << y << std::endl;
  442.  
  443.     num::isort(x.begin(), x.end());
  444.     std::cout << "isort(x): " << std::endl << x << std::endl;
  445.  
  446.     x = cpx;
  447.     num::isort2(x.begin(), x.end(), y.begin());
  448.     std::cout << "isort2(x, y): " << std::endl << x << y << std::endl;
  449.  
  450.     x = cpx;
  451.     num::qsort(x.begin(), x.end());
  452.     std::cout << "qsort(x): " << std::endl << x << std::endl;
  453.  
  454.     x = cpx;
  455.     y = cpy;
  456.     num::qsort2(x.begin(), x.end(), y.begin());
  457.     std::cout << "qsort2(x, y): " << std::endl << x << y << std::endl;
  458.  
  459.     x = cpx;
  460.     y = cpy;
  461.     num::qsort2key(x.begin(), x.end(), y.begin());
  462.     std::cout << "qsort2key(x, y): " << std::endl << x << y << std::endl;
  463. }
  464.  
  465. int main()
  466. {
  467.     srand(10203040L);
  468.  
  469.     std::cout << "Testing Algorithms:" << std::endl;
  470.     testAlgorithm();
  471.  
  472.     std::cout << std::endl << std::endl << "Testing Continuous Distributions:" << std::endl;
  473.     testContDist();
  474.  
  475.     std::cout << std::endl << std::endl << "Testing Correlation:" << std::endl;
  476.     testCorrelation();
  477.  
  478.     std::cout << std::endl << std::endl << "Testing Descriptive Statistics:" << std::endl;
  479.     testDescript();
  480.  
  481.     std::cout << std::endl << std::endl << "Testing Deviate:" << std::endl;
  482.     testDeviate();
  483.  
  484.     std::cout << std::endl << std::endl << "Testing Discrete Distributions:" << std::endl;
  485.     testDiscDist();
  486.  
  487.     std::cout << std::endl << std::endl << "Testing Eigen Systems:" << std::endl;
  488.     testEigenSys();
  489.  
  490.     std::cout << std::endl << std::endl << "Testing Linear Algebra:" << std::endl;
  491.     testLinAlg();
  492.  
  493.     std::cout << std::endl << std::endl << "Testing Ranking:" << std::endl;
  494.     testRanking();
  495.  
  496.     std::cout << std::endl << std::endl << "Testing Residuals:" << std::endl;
  497.     testResidual();
  498.  
  499.     std::cout << std::endl << std::endl << "Testing Rootfind:" << std::endl;
  500.     testRootfind();
  501.  
  502.     std::cout << std::endl << std::endl << "Testing Sorting:" << std::endl;
  503.     testSort();
  504.  
  505.     return 0;
  506. }
  507.  
  508. //===================================================================
  509. // Revision History
  510. //
  511. // Version 1.0 - 05/17/1999 - New.
  512. //===================================================================
  513.