1 /******************************************************************************
2 * Author: Laurent Kneip *
3 * Contact: kneip.laurent@gmail.com *
4 * License: Copyright (c) 2013 Laurent Kneip, ANU. All rights reserved. *
5 * *
6 * Redistribution and use in source and binary forms, with or without *
7 * modification, are permitted provided that the following conditions *
8 * are met: *
9 * * Redistributions of source code must retain the above copyright *
10 * notice, this list of conditions and the following disclaimer. *
11 * * Redistributions in binary form must reproduce the above copyright *
12 * notice, this list of conditions and the following disclaimer in the *
13 * documentation and/or other materials provided with the distribution. *
14 * * Neither the name of ANU nor the names of its contributors may be *
15 * used to endorse or promote products derived from this software without *
16 * specific prior written permission. *
17 * *
18 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"*
19 * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE *
20 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE *
21 * ARE DISCLAIMED. IN NO EVENT SHALL ANU OR THE CONTRIBUTORS BE LIABLE *
22 * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL *
23 * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR *
24 * SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER *
25 * CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT *
26 * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY *
27 * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF *
28 * SUCH DAMAGE. *
29 ******************************************************************************/
30
31 #include <stdlib.h>
32 #include <stdio.h>
33 #include <iostream>
34 #include <opengv/math/Sturm.hpp>
35
36 #include "time_measurement.hpp"
37
38
39 using namespace std;
40 using namespace Eigen;
41 using namespace opengv;
42
main(int argc,char ** argv)43 int main( int argc, char** argv )
44 {
45 std::vector<double> coeffs;
46 coeffs.push_back(1.0);
47 coeffs.push_back(4.0);
48 coeffs.push_back(1.0);
49 coeffs.push_back(-6.0);
50
51 //timer
52 struct timeval tic;
53 struct timeval toc;
54 size_t iterations = 50;
55
56 //for now just construct the problem
57 gettimeofday( &tic, 0 );
58 for(size_t i = 0; i < iterations; i++)
59 math::Sturm sturm(coeffs);
60 gettimeofday( &toc, 0 );
61 double time = TIMETODOUBLE(timeval_minus(toc,tic)) / iterations;
62
63 std::cout << "the initialization takes " << time << " seconds" << std::endl;
64
65 math::Sturm sturm(coeffs);
66
67 //test the lagrangian bounds
68 gettimeofday( &tic, 0 );
69 double bound;
70 for(size_t i = 0; i < iterations; i++)
71 bound = sturm.computeLagrangianBound();
72 gettimeofday( &toc, 0 );
73 time = TIMETODOUBLE(timeval_minus(toc,tic)) / iterations;
74
75 std::cout << "the initial bound computation takes " << time << " seconds" << std::endl;
76 std::cout << "the bound evaluates to " << bound << std::endl;
77
78 //now evaluate the chain to verify that the number of roots is 3
79 gettimeofday( &tic, 0 );
80 size_t numberRoots;
81 for(size_t i = 0; i < iterations; i++)
82 numberRoots = sturm.evaluateChain(-bound) - sturm.evaluateChain(bound);
83 gettimeofday( &toc, 0 );
84 time = TIMETODOUBLE(timeval_minus(toc,tic)) / iterations;
85
86 std::cout << "the evaluation of two bounds takes " << time << " seconds" << std::endl;
87 std::cout << "the obtained number of roots is " << numberRoots << std::endl;
88
89 //now test the bracketing mechanism
90 gettimeofday( &tic, 0 );
91 std::vector<double> roots;
92 for(size_t i = 0; i < iterations; i++)
93 {
94 roots.clear();
95 sturm.bracketRoots(roots,0.5);
96 }
97 gettimeofday( &toc, 0 );
98 time = TIMETODOUBLE(timeval_minus(toc,tic)) / iterations;
99
100 std::cout << "the bracketing of the roots took " << time << " seconds" << std::endl;
101 std::cout << "the obtained brackets are:" << std::endl;
102 std::vector<double>::iterator it = roots.begin();
103 while( it != roots.end() )
104 {
105 std::cout << "root:" << (*it) << std::endl;
106 it++;
107 }
108
109 //now test the entire root-finding mechanism
110 gettimeofday( &tic, 0 );
111 for(size_t i = 0; i < iterations; i++)
112 roots = sturm.findRoots();
113 gettimeofday( &toc, 0 );
114 time = TIMETODOUBLE(timeval_minus(toc,tic)) / iterations;
115
116 std::cout << "the entire finding of the roots took " << time << " seconds" << std::endl;
117 std::cout << "the obtained roots are:" << std::endl;
118 std::vector<double>::iterator it2 = roots.begin();
119 while( it2 != roots.end() )
120 {
121 std::cout << (*it2) << std::endl;
122 it2++;
123 }
124
125 //now test the new root-finding with inbuild gradient descent
126 gettimeofday( &tic, 0 );
127 for(size_t i = 0; i < iterations; i++)
128 {
129 roots.clear();
130 sturm.findRoots2(roots);
131 }
132 gettimeofday( &toc, 0 );
133 time = TIMETODOUBLE(timeval_minus(toc,tic)) / iterations;
134
135 std::cout << "the entire finding of the roots took " << time << " seconds" << std::endl;
136 std::cout << "the obtained roots are:" << std::endl;
137 it2 = roots.begin();
138 while( it2 != roots.end() )
139 {
140 std::cout << (*it2) << std::endl;
141 it2++;
142 }
143
144 return 0;
145 }
146