123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169 |
- #include <iostream>
- #include <iomanip>
- #include <random>
- #include "ra/bench.hh"
- using std::cout, std::endl, std::flush, ra::TestRecorder;
- using real = double;
- int nx = 100;
- int ny = 100;
- int nz = 100;
- int ts = 10;
- auto I = ra::iota(nx-2, 1);
- auto J = ra::iota(ny-2, 1);
- auto K = ra::iota(nz-2, 1);
- constexpr ra::Small<real, 3, 3, 3> mask = { 0, 0, 0, 0, 1, 0, 0, 0, 0,
- 0, 1, 0, 1, -6, 1, 0, 1, 0,
- 0, 0, 0, 0, 1, 0, 0, 0, 0 };
- #define THEOP template <class A_, class Anext_, class Astencil_> __attribute__((noinline)) \
- auto operator()(A_ & A, Anext_ & Anext, Astencil_ & Astencil)
- struct f_raw
- {
- THEOP
- {
- for (int i=1; i+1<nx; ++i) {
- for (int j=1; j+1<ny; ++j) {
- for (int k=1; k+1<nz; ++k) {
- Anext(i, j, k) = -6*A(i, j, k)
- + A(i+1, j, k) + A(i, j+1, k) + A(i, j, k+1)
- + A(i-1, j, k) + A(i, j-1, k) + A(i, j, k-1);
- }
- }
- }
- std::swap(A.cp, Anext.cp);
- };
- };
- struct f_slices
- {
- THEOP
- {
- Anext(I, J, K) = -6*A(I, J, K)
- + A(I+1, J, K) + A(I, J+1, K) + A(I, J, K+1)
- + A(I-1, J, K) + A(I, J-1, K) + A(I, J, K-1);
- std::swap(A.cp, Anext.cp);
- };
- };
- struct f_stencil_explicit
- {
- THEOP
- {
- Astencil.cp = A.data();
- Anext(I, J, K) = map([](auto && A) { return -6*A(1, 1, 1)
- + A(2, 1, 1) + A(1, 2, 1) + A(1, 1, 2)
- + A(0, 1, 1) + A(1, 0, 1) + A(1, 1, 0); },
- iter<3>(Astencil));
- std::swap(A.cp, Anext.cp);
- };
- };
- struct f_stencil_arrayop
- {
- THEOP
- {
- Astencil.cp = A.data();
- Anext(I, J, K) = map([](auto && s) { return sum(s*mask); }, iter<3>(Astencil));
- std::swap(A.cp, Anext.cp);
- };
- };
- struct f_sumprod
- {
- THEOP
- {
- Astencil.cp = A.data();
- Anext(I, J, K) = 0;
- Anext(I, J, K) += map(ra::wrank<3, 3>(std::multiplies<>()), Astencil, mask);
- std::swap(A.cp, Anext.cp);
- };
- };
- struct f_sumprod2
- {
- THEOP
- {
- Astencil.cp = A.data();
- Anext(I, J, K) = 0;
- ply_fixed(map(ra::wrank<0, 3, 3>([](auto && A, auto && B, auto && C) { A += B*C; }), Anext(I, J, K), Astencil, mask));
- std::swap(A.cp, Anext.cp);
- };
- };
- int main()
- {
- TestRecorder tr(std::cout);
- std::random_device rand;
- real value = rand();
- auto bench = [&](auto & A, auto & Anext, auto & Astencil, auto && ref, auto && tag, auto && f)
- {
- auto bv = Benchmark().repeats(ts).runs(3)
- .once_f([&](auto && repeat)
- {
- Anext = 0.;
- A = value;
- repeat([&]() { f(A, Anext, Astencil); });
- });
- tr.info(std::setw(5), std::fixed, Benchmark::avg(bv)/A.size()/1e-9, " ns [",
- Benchmark::stddev(bv)/A.size()/1e-9 ,"] ", tag)
- .test_rel(ref, A, 1e-11);
- };
- ra::Big<real, 3> Aref;
- tr.section("static rank");
- {
- ra::Big<real, 3> A({nx, ny, nz}, 1.);
- ra::Big<real, 3> Anext({nx, ny, nz}, 0.);
- auto Astencil = stencil(A, 1, 1);
- cout << "Astencil " << format_array(Astencil(0, 0, 0, ra::dots<3>), "|", " ") << endl;
- #define BENCH(ref, op) bench(A, Anext, Astencil, ref, STRINGIZE(op), op {});
- BENCH(A, f_raw);
- Aref = ra::Big<real, 3>(A);
- BENCH(Aref, f_slices);
- BENCH(Aref, f_stencil_explicit);
- BENCH(Aref, f_stencil_arrayop);
- BENCH(Aref, f_sumprod);
- BENCH(Aref, f_sumprod2);
- #undef BENCH
- }
- tr.section("dynamic rank");
- {
- ra::Big<real> B({nx, ny, nz}, 1.);
- ra::Big<real> Bnext({nx, ny, nz}, 0.);
- auto Bstencil = stencil(B, 1, 1);
- cout << "Bstencil " << format_array(Bstencil(0, 0, 0, ra::dots<3>), "|", " ") << endl;
- #define BENCH(ref, op) bench(B, Bnext, Bstencil, ref, STRINGIZE(op), op {});
-
- BENCH(Aref, f_slices);
- BENCH(Aref, f_stencil_explicit);
- BENCH(Aref, f_stencil_arrayop);
- #undef BENCH
- }
- return tr.summary();
- }
|