1 /* ----------------------------------------------------------------------
2 LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
3 https://www.lammps.org/, Sandia National Laboratories
4 Steve Plimpton, sjplimp@sandia.gov
5
6 Copyright (2003) Sandia Corporation. Under the terms of Contract
7 DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
8 certain rights in this software. This software is distributed under
9 the GNU General Public License.
10
11 See the README file in the top-level LAMMPS directory.
12 ------------------------------------------------------------------------- */
13
14 #include "imbalance_time.h"
15
16 #include "atom.h"
17 #include "error.h"
18 #include "timer.h"
19
20 using namespace LAMMPS_NS;
21
22 #define BIG 1.0e20
23
24 /* -------------------------------------------------------------------- */
25
ImbalanceTime(LAMMPS * lmp)26 ImbalanceTime::ImbalanceTime(LAMMPS *lmp) : Imbalance(lmp) {}
27
28 /* -------------------------------------------------------------------- */
29
options(int narg,char ** arg)30 int ImbalanceTime::options(int narg, char **arg)
31 {
32 if (narg < 1) error->all(FLERR, "Illegal balance weight command");
33 factor = utils::numeric(FLERR, arg[0], false, lmp);
34 if (factor <= 0.0) error->all(FLERR, "Illegal balance weight command");
35 return 1;
36 }
37
38 /* ----------------------------------------------------------------------
39 reset last and timers if necessary
40 ------------------------------------------------------------------------- */
41
init(int flag)42 void ImbalanceTime::init(int flag)
43 {
44 last = 0.0;
45
46 // flag = 1 if called from FixBalance at start of run
47 // init Timer, so accumulated time not carried over from previous run
48 // should NOT init Timer if called from Balance, it uses time from last run
49
50 if (flag) timer->init();
51 }
52
53 /* -------------------------------------------------------------------- */
54
compute(double * weight)55 void ImbalanceTime::compute(double *weight)
56 {
57 if (!timer->has_normal()) return;
58
59 // cost = CPU time for relevant timers since last invocation
60 // localwt = weight assigned to each owned atom
61 // just return if no time yet tallied
62 // we 0.1 seconds as a minimum time to avoid computation of bogus
63 // load balancing weights due to limited timer resolution/precision
64
65 double cost = -last;
66 cost += timer->get_wall(Timer::PAIR);
67 cost += timer->get_wall(Timer::NEIGH);
68 cost += timer->get_wall(Timer::BOND);
69 cost += timer->get_wall(Timer::KSPACE);
70 cost += 0.1;
71
72 double maxcost;
73 MPI_Allreduce(&cost, &maxcost, 1, MPI_DOUBLE, MPI_MAX, world);
74 if (maxcost <= 0.1) return;
75
76 int nlocal = atom->nlocal;
77 double localwt = 0.0;
78 if (nlocal) localwt = cost / nlocal;
79
80 if (nlocal && localwt <= 0.0) error->one(FLERR, "Balance weight <= 0.0");
81
82 // apply factor if specified != 1.0
83 // wtlo,wthi = lo/hi values excluding 0.0 due to no atoms on this proc
84 // lo value does not change
85 // newhi = new hi value to give hi/lo ratio factor times larger/smaller
86 // expand/contract all localwt values from lo->hi to lo->newhi
87
88 if (factor != 1.0) {
89 double wtlo, wthi;
90 if (localwt == 0.0) localwt = BIG;
91 MPI_Allreduce(&localwt, &wtlo, 1, MPI_DOUBLE, MPI_MIN, world);
92 if (localwt == BIG) localwt = 0.0;
93 MPI_Allreduce(&localwt, &wthi, 1, MPI_DOUBLE, MPI_MAX, world);
94 if (wtlo == wthi) return;
95
96 double newhi = wthi * factor;
97 localwt = wtlo + ((localwt - wtlo) / (wthi - wtlo)) * (newhi - wtlo);
98 }
99
100 for (int i = 0; i < nlocal; i++) weight[i] *= localwt;
101
102 // record time up to this point
103
104 last += cost;
105 }
106
107 /* -------------------------------------------------------------------- */
108
info()109 std::string ImbalanceTime::info()
110 {
111 return fmt::format(" time weight factor: {}\n", factor);
112 }
113