1 // Copyright (c) 2014 Hartmut Kaiser
2 // Copyright (c) 2014 Patricia Grubel
3 //
4 // Distributed under the Boost Software License, Version 1.0. (See accompanying
5 // file LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
6
7 // This is the first in a series of examples demonstrating the development of a
8 // fully distributed solver for a simple 1D heat distribution problem.
9 //
10 // This example provides a serial base line implementation. No parallelization
11 // is performed.
12
13 #include <hpx/hpx_init.hpp>
14 #include <hpx/hpx.hpp>
15
16 #include <cstddef>
17 #include <cstdint>
18 #include <iostream>
19 #include <vector>
20
21 #include "print_time_results.hpp"
22
23 ///////////////////////////////////////////////////////////////////////////////
24 // Command-line variables
25 bool header = true; // print csv heading
26 double k = 0.5; // heat transfer coefficient
27 double dt = 1.; // time step
28 double dx = 1.; // grid spacing
29
30 ///////////////////////////////////////////////////////////////////////////////
31 //[stepper_1
32 struct stepper
33 {
34 // Our partition type
35 typedef double partition;
36
37 // Our data for one time step
38 typedef std::vector<partition> space;
39
40 // Our operator
heatstepper41 static double heat(double left, double middle, double right)
42 {
43 return middle + (k*dt/(dx*dx)) * (left - 2*middle + right);
44 }
45
46 // do all the work on 'nx' data points for 'nt' time steps
do_workstepper47 space do_work(std::size_t nx, std::size_t nt)
48 {
49 // U[t][i] is the state of position i at time t.
50 std::vector<space> U(2);
51 for (space& s : U)
52 s.resize(nx);
53
54 // Initial conditions: f(0, i) = i
55 for (std::size_t i = 0; i != nx; ++i)
56 U[0][i] = double(i);
57
58 // Actual time step loop
59 for (std::size_t t = 0; t != nt; ++t)
60 {
61 space const& current = U[t % 2];
62 space& next = U[(t + 1) % 2];
63
64 next[0] = heat(current[nx-1], current[0], current[1]);
65
66 for (std::size_t i = 1; i != nx-1; ++i)
67 next[i] = heat(current[i-1], current[i], current[i+1]);
68
69 next[nx-1] = heat(current[nx-2], current[nx-1], current[0]);
70 }
71
72 // Return the solution at time-step 'nt'.
73 return U[nt % 2];
74 }
75 };
76 //]
77 ///////////////////////////////////////////////////////////////////////////////
hpx_main(boost::program_options::variables_map & vm)78 int hpx_main(boost::program_options::variables_map& vm)
79 {
80 std::uint64_t nx = vm["nx"].as<std::uint64_t>(); // Number of grid points.
81 std::uint64_t nt = vm["nt"].as<std::uint64_t>(); // Number of steps.
82
83 if (vm.count("no-header"))
84 header = false;
85
86 // Create the stepper object
87 stepper step;
88
89 // Measure execution time.
90 std::uint64_t t = hpx::util::high_resolution_clock::now();
91
92 // Execute nt time steps on nx grid points.
93 stepper::space solution = step.do_work(nx, nt);
94
95 // Print the final solution
96 if (vm.count("results"))
97 {
98 for (std::size_t i = 0; i != nx; ++i)
99 std::cout << "U[" << i << "] = " << solution[i] << std::endl;
100 }
101
102 std::uint64_t elapsed = hpx::util::high_resolution_clock::now() - t;
103
104 std::uint64_t const os_thread_count = hpx::get_os_thread_count();
105 print_time_results(os_thread_count, elapsed, nx, nt, header);
106
107 return hpx::finalize();
108 }
109
main(int argc,char * argv[])110 int main(int argc, char* argv[])
111 {
112 namespace po = boost::program_options;
113
114 po::options_description desc_commandline;
115 desc_commandline.add_options()
116 ("results", "print generated results (default: false)")
117 ("nx", po::value<std::uint64_t>()->default_value(100),
118 "Local x dimension")
119 ("nt", po::value<std::uint64_t>()->default_value(45),
120 "Number of time steps")
121 ("k", po::value<double>(&k)->default_value(0.5),
122 "Heat transfer coefficient (default: 0.5)")
123 ("dt", po::value<double>(&dt)->default_value(1.0),
124 "Timestep unit (default: 1.0[s])")
125 ("dx", po::value<double>(&dx)->default_value(1.0),
126 "Local x dimension")
127 ( "no-header", "do not print out the csv header row")
128 ;
129
130 // Initialize and run HPX
131 return hpx::init(desc_commandline, argc, argv);
132 }
133