123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778 |
- #include <iomanip>
- #include <chrono>
- #include "ra/test.hh"
- using std::cout, std::endl, std::flush, ra::TestRecorder;
- auto now() { return std::chrono::high_resolution_clock::now(); }
- template <class DT> auto ms(DT && dt) { return std::chrono::duration_cast<std::chrono::milliseconds>(dt).count(); }
- int main()
- {
- int numIters = 301;
- int N = 64;
- ra::Big<float, 3> A({N, N, N}, ra::none);
- ra::Big<float, 3> B({N, N, N}, ra::none);
- auto interior = ra::iota(N/2, N/4);
- A = 22.;
- A(interior, interior, interior) = 30.;
- auto I = ra::iota(N-2, 1);
- auto J = ra::iota(N-2, 1);
- auto K = ra::iota(N-2, 1);
- {
- std::chrono::duration<float> dt(0);
- double c = 1/6.5;
- for (int i=0; i<numIters; ++i) {
- auto t0 = now();
- B(I, J, K) = c * (.5 * A(I, J, K) + A(I+1, J, K) + A(I-1, J, K)
- + A(I, J+1, K) + A(I, J-1, K) + A(I, J, K+1) + A(I, J, K-1));
- A(I, J, K) = c * (.5 * B(I, J, K) + B(I+1, J, K) + B(I-1, J, K)
- + B(I, J+1, K) + B(I, J-1, K) + B(I, J, K+1) + B(I, J, K-1));
- dt += (now()-t0);
- cout << std::setprecision(4) << std::fixed << ra::noshape << A(N/2, N/2, ra::iota(8, 0, N/8)) << endl;
- }
- cout << std::setw(10) << std::fixed << (ms(dt)/double(numIters)) << " ms / iter " << endl;
- }
- ra::Big<float, 3> first_A = A;
- A = 22.;
- A(interior, interior, interior) = 30.;
- auto Im = ra::iota(N-2, 0);
- auto Ip = ra::iota(N-2, 2);
- {
- std::chrono::duration<float> dt(0);
- double c = 1/6.5;
- for (int i=0; i<numIters; ++i) {
- auto t0 = now();
- B(I, I, I) = c * (.5 * A(I, I, I) + A(Ip, I, I) + A(Im, I, I)
- + A(I, Ip, I) + A(I, Im, I) + A(I, I, Ip) + A(I, I, Im));
- A(I, I, I) = c * (.5 * B(I, I, I) + B(Ip, I, I) + B(Im, I, I)
- + B(I, Ip, I) + B(I, Im, I) + B(I, I, Ip) + B(I, I, Im));
- dt += (now()-t0);
- cout << std::setprecision(4) << std::fixed << ra::noshape << A(N/2, N/2, ra::iota(8, 0, N/8)) << endl;
- }
- cout << std::setw(10) << std::fixed << (ms(dt)/double(numIters)) << " ms / iter " << endl;
- }
- TestRecorder tr(std::cout);
- tr.quiet().test_rel(first_A, A, 0.);
- return tr.summary();
- }
|