Sacado Package Browser (Single Doxygen Collection) Version of the Day
Loading...
Searching...
No Matches
fad_fe_jac_fill_range.cpp
Go to the documentation of this file.
1// @HEADER
2// ***********************************************************************
3//
4// Sacado Package
5// Copyright (2006) Sandia Corporation
6//
7// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
8// the U.S. Government retains certain rights in this software.
9//
10// This library is free software; you can redistribute it and/or modify
11// it under the terms of the GNU Lesser General Public License as
12// published by the Free Software Foundation; either version 2.1 of the
13// License, or (at your option) any later version.
14//
15// This library is distributed in the hope that it will be useful, but
16// WITHOUT ANY WARRANTY; without even the implied warranty of
17// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
18// Lesser General Public License for more details.
19//
20// You should have received a copy of the GNU Lesser General Public
21// License along with this library; if not, write to the Free Software
22// Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301
23// USA
24// Questions? Contact David M. Gay (dmgay@sandia.gov) or Eric T. Phipps
25// (etphipp@sandia.gov).
26//
27// ***********************************************************************
28// @HEADER
29
30#include "fe_jac_fill_funcs.hpp"
31
32std::vector<double>
33do_times(int work_count, int num_eqns_begin, int num_eqns_end,
34 int num_eqns_delta,
35 double (*func)(unsigned int,unsigned int,double)) {
36 std::vector<double> times;
37 for (int num_eqns = num_eqns_begin; num_eqns <= num_eqns_end;
38 num_eqns += num_eqns_delta) {
39 int num_nodes = work_count / num_eqns;
40 double mesh_spacing = 1.0 / (num_nodes - 1);
41 times.push_back(func(num_nodes, num_eqns, mesh_spacing));
42 }
43 return times;
44}
45
46template <template <typename> class FadType>
47std::vector<double>
48do_times_fad(int work_count, int num_eqns_begin, int num_eqns_end,
49 int num_eqns_delta) {
50 std::vector<double> times;
51 for (int num_eqns = num_eqns_begin; num_eqns <= num_eqns_end;
52 num_eqns += num_eqns_delta) {
53 int num_nodes = work_count / num_eqns;
54 double mesh_spacing = 1.0 / (num_nodes - 1);
55 times.push_back(fad_jac_fill<FadType<double> >(num_nodes, num_eqns, mesh_spacing));
56 }
57 return times;
58}
59
60template <template <typename,int> class FadType>
61std::vector<double>
62do_times_slfad(int work_count, int num_eqns_begin, int num_eqns_end,
63 int num_eqns_delta) {
64 const int slfad_max = 130; // Maximum number of derivative components for SLFad
65 std::vector<double> times;
66 for (int num_eqns = num_eqns_begin; num_eqns <= num_eqns_end;
67 num_eqns += num_eqns_delta) {
68 int num_nodes = work_count / num_eqns;
69 double mesh_spacing = 1.0 / (num_nodes - 1);
70 if (num_eqns*2 < slfad_max)
71 times.push_back(fad_jac_fill<FadType<double,slfad_max> >(num_nodes, num_eqns, mesh_spacing));
72 }
73 return times;
74}
75
76template <template <typename,int> class FadType>
77std::vector<double>
78do_times_sfad(int work_count, int num_eqns_begin, int num_eqns_end,
79 int num_eqns_delta) {
80 std::vector<double> times;
81 for (int num_eqns = num_eqns_begin; num_eqns <= num_eqns_end;
82 num_eqns += num_eqns_delta) {
83 int num_nodes = work_count / num_eqns;
84 double mesh_spacing = 1.0 / (num_nodes - 1);
85 if (num_eqns*2 == 10)
86 times.push_back(fad_jac_fill<FadType<double,10> >(num_nodes, num_eqns, mesh_spacing));
87 else if (num_eqns*2 == 30)
88 times.push_back(fad_jac_fill<FadType<double,30> >(num_nodes, num_eqns, mesh_spacing));
89 else if (num_eqns*2 == 50)
90 times.push_back(fad_jac_fill<FadType<double,50> >(num_nodes, num_eqns, mesh_spacing));
91 else if (num_eqns*2 == 70)
92 times.push_back(fad_jac_fill<FadType<double,70> >(num_nodes, num_eqns, mesh_spacing));
93 else if (num_eqns*2 == 90)
94 times.push_back(fad_jac_fill<FadType<double,90> >(num_nodes, num_eqns, mesh_spacing));
95 else if (num_eqns*2 == 110)
96 times.push_back(fad_jac_fill<FadType<double,110> >(num_nodes, num_eqns, mesh_spacing));
97 else if (num_eqns*2 == 130)
98 times.push_back(fad_jac_fill<FadType<double,130> >(num_nodes, num_eqns, mesh_spacing));
99 }
100 return times;
101}
102
103void print_times(const std::vector<double>& times,
104 const std::vector<double>& base,
105 const std::string& name, int p, int w, int w_name) {
106 std::cout.setf(std::ios::scientific);
107 std::cout.precision(p);
108 std::cout.setf(std::ios::right);
109 std::cout << std::setw(w_name) << name << " ";
110 std::cout.setf(std::ios::right);
111 for (unsigned int i=0; i<times.size(); i++)
112 std::cout << std::setw(w) << times[i]/base[i] << " ";
113 std::cout << std::endl;
114}
115
116int main(int argc, char* argv[]) {
117 int ierr = 0;
118 int p = 1;
119 int w = p+7;
120 int w_name = 13;
121
122 try {
123
124 // Set up command line options
125 Teuchos::CommandLineProcessor clp;
126 clp.setDocString("This program tests the speed of various forward mode AD implementations for a finite-element-like Jacobian fill");
127 int work_count = 200000;
128 int num_eqns_begin = 5;
129 int num_eqns_end = 65;
130 int num_eqns_delta = 10;
131 int rt = 0;
132 clp.setOption("wc", &work_count, "Work count = num_nodes*num_eqns");
133 clp.setOption("p_begin", &num_eqns_begin, "Intitial number of equations");
134 clp.setOption("p_end", &num_eqns_end, "Final number of equations");
135 clp.setOption("p_delta", &num_eqns_delta, "Step in number of equations");
136 clp.setOption("rt", &rt, "Include ADOL-C retaping test");
137
138 // Parse options
139 Teuchos::CommandLineProcessor::EParseCommandLineReturn
140 parseReturn= clp.parse(argc, argv);
141 if(parseReturn != Teuchos::CommandLineProcessor::PARSE_SUCCESSFUL)
142 return 1;
143
144 // Print header
145 std::cout.setf(std::ios::right);
146 std::cout << std::setw(w_name) << "Name" << " ";
147 for (int num_eqns = num_eqns_begin; num_eqns <= num_eqns_end;
148 num_eqns += num_eqns_delta)
149 std::cout << std::setw(w) << num_eqns << " ";
150 std::cout << std::endl;
151 for (int j=0; j<w_name; j++)
152 std::cout << '=';
153 std::cout << " ";
154 for (int num_eqns = num_eqns_begin; num_eqns <= num_eqns_end;
155 num_eqns += num_eqns_delta) {
156 for (int j=0; j<w; j++)
157 std::cout << '=';
158 std::cout << " ";
159 }
160 std::cout << std::endl;
161
162 // Analytic
163 std::vector<double> times_analytic =
164 do_times(work_count, num_eqns_begin, num_eqns_end, num_eqns_delta,
166 print_times(times_analytic, times_analytic, "Analytic", p, w, w_name);
167
168#ifdef HAVE_ADIC
169 // Note there seems to be a bug in ADIC where doing more than one num_eqns
170 // value results in incorrect timings after the first. Doing one value
171 // at a time seems to give correct values though.
172 std::vector<double> times_adic =
173 do_times(work_count, num_eqns_begin, num_eqns_end, num_eqns_delta,
174 adic_jac_fill);
175 print_times(times_adic, times_analytic, "ADIC", p, w, w_name);
176#endif
177
178 // Original Fad
179 std::vector<double> times_sfad =
180 do_times_sfad<Sacado::Fad::SFad>(
181 work_count, num_eqns_begin, num_eqns_end, num_eqns_delta);
182 print_times(times_sfad, times_analytic, "SFAD", p, w, w_name);
183
184 std::vector<double> times_slfad =
185 do_times_sfad<Sacado::Fad::SLFad>(
186 work_count, num_eqns_begin, num_eqns_end, num_eqns_delta);
187 print_times(times_slfad, times_analytic, "SLFAD", p, w, w_name);
188
189 std::vector<double> times_dfad =
190 do_times_fad<Sacado::Fad::DFad>(
191 work_count, num_eqns_begin, num_eqns_end, num_eqns_delta);
192 print_times(times_dfad, times_analytic, "DFAD", p, w, w_name);
193
194
195 // ELR Fad
196 std::vector<double> times_elr_sfad =
197 do_times_sfad<Sacado::ELRFad::SFad>(
198 work_count, num_eqns_begin, num_eqns_end, num_eqns_delta);
199 print_times(times_elr_sfad, times_analytic, "ELRSFAD", p, w, w_name);
200
201 std::vector<double> times_elr_slfad =
202 do_times_sfad<Sacado::ELRFad::SLFad>(
203 work_count, num_eqns_begin, num_eqns_end, num_eqns_delta);
204 print_times(times_elr_slfad, times_analytic, "ELRSLFAD", p, w, w_name);
205
206 std::vector<double> times_elr_dfad =
207 do_times_fad<Sacado::ELRFad::DFad>(
208 work_count, num_eqns_begin, num_eqns_end, num_eqns_delta);
209 print_times(times_elr_dfad, times_analytic, "ELRDFAD", p, w, w_name);
210
211
212 // Cache Fad
213 std::vector<double> times_cache_sfad =
214 do_times_sfad<Sacado::CacheFad::SFad>(
215 work_count, num_eqns_begin, num_eqns_end, num_eqns_delta);
216 print_times(times_cache_sfad, times_analytic, "CacheSFAD", p, w, w_name);
217
218 std::vector<double> times_cache_slfad =
219 do_times_sfad<Sacado::CacheFad::SLFad>(
220 work_count, num_eqns_begin, num_eqns_end, num_eqns_delta);
221 print_times(times_cache_slfad, times_analytic, "CacheSLFAD", p, w, w_name);
222
223 std::vector<double> times_cache_dfad =
224 do_times_fad<Sacado::CacheFad::DFad>(
225 work_count, num_eqns_begin, num_eqns_end, num_eqns_delta);
226 print_times(times_cache_dfad, times_analytic, "CacheDFAD", p, w, w_name);
227
228 // ELR Cache Fad
229 std::vector<double> times_cache_elr_sfad =
230 do_times_sfad<Sacado::ELRCacheFad::SFad>(
231 work_count, num_eqns_begin, num_eqns_end, num_eqns_delta);
232 print_times(times_cache_elr_sfad, times_analytic, "ELRCacheSFAD", p, w, w_name);
233
234 std::vector<double> times_cache_elr_slfad =
235 do_times_sfad<Sacado::ELRCacheFad::SLFad>(
236 work_count, num_eqns_begin, num_eqns_end, num_eqns_delta);
237 print_times(times_cache_elr_slfad, times_analytic, "ELRCacheSLFAD", p, w, w_name);
238
239 std::vector<double> times_cache_elr_dfad =
240 do_times_fad<Sacado::ELRCacheFad::DFad>(
241 work_count, num_eqns_begin, num_eqns_end, num_eqns_delta);
242 print_times(times_cache_elr_dfad, times_analytic, "ELRCacheDFAD", p, w, w_name);
243
244 }
245 catch (std::exception& e) {
246 std::cout << e.what() << std::endl;
247 ierr = 1;
248 }
249 catch (const char *s) {
250 std::cout << s << std::endl;
251 ierr = 1;
252 }
253 catch (...) {
254 std::cout << "Caught unknown exception!" << std::endl;
255 ierr = 1;
256 }
257
258 return ierr;
259}
int main()
Definition: ad_example.cpp:191
const T func(int n, T *x)
Definition: ad_example.cpp:49
std::vector< double > do_times_slfad(int work_count, int num_eqns_begin, int num_eqns_end, int num_eqns_delta)
std::vector< double > do_times(int work_count, int num_eqns_begin, int num_eqns_end, int num_eqns_delta, double(*func)(unsigned int, unsigned int, double))
void print_times(const std::vector< double > &times, const std::vector< double > &base, const std::string &name, int p, int w, int w_name)
std::vector< double > do_times_sfad(int work_count, int num_eqns_begin, int num_eqns_end, int num_eqns_delta)
std::vector< double > do_times_fad(int work_count, int num_eqns_begin, int num_eqns_end, int num_eqns_delta)
double analytic_jac_fill(unsigned int num_nodes, unsigned int num_eqns, double mesh_spacing)
double fad_jac_fill(unsigned int num_nodes, unsigned int num_eqns, double mesh_spacing)
const char * p