503 return {
"srvn",
"srvn.ph",
"srvn.cs",
"flat",
"flat.cs",
"flat.ph",
"moment3",
"default"};
515 opt.layer_ssa.verbose =
false;
517 for (
char& ch : opt.interlock_refpath_scope)
518 ch =
char(std::tolower(
static_cast<unsigned char>(ch)));
519 if (opt.interlock_refpath_scope !=
"merging" && opt.interlock_refpath_scope !=
"all")
520 throw InputError(
"Unknown config.interlock_refpath_scope '" +
521 opt.interlock_refpath_scope +
"', use 'merging' or 'all'.");
522 if (std::isnan(opt.interlock_maxpaths) || opt.interlock_maxpaths < 0.0)
523 throw InputError(
"config.interlock_maxpaths must be a non-negative number.");
533 if (opt.method ==
"mw.upper" || opt.method ==
"mw.lower")
return box_bounds();
537 return is_ph_encoding() ? aggregate_ph() : aggregate();
553 if (lqn.nentries == 0)
return entrycdfrespt;
554 if (entrycdfrespt.size() <= lqn.nentries || entrycdfrespt[1].empty()) {
561 if (is_ph_encoding())
563 "getCdfRespT needs the routing encoding of the activity graph, which "
564 "method='srvn.ph' does not build. Rebuild the solver with "
565 "method='srvn.cs' or method='moment3'.");
572 const std::string saved = opt.method;
573 const std::string saved_resolved = lnmethod;
574 opt.method =
"moment3";
575 lnmethod =
"moment3";
579 lnmethod = saved_resolved;
581 return entrycdfrespt;
592 if (opt.ln_transient ==
"decoupled")
return tran_avg_decoupled();
593 if (opt.ln_transient ==
"coupled")
return tran_avg_coupled();
594 throw InputError(
"SolverLN: unknown ln_transient mode '" + opt.ln_transient +
595 "' (use 'coupled' or 'decoupled')");
611 if (results.empty()) iterate();
615 for (std::size_t e = 0; e < ensemble.size(); ++e) {
618 const bool exact_available =
619 opt.layer_solver ==
"mva" || opt.layer_solver ==
"nc";
621 ensemble[e], so, exact_available && !fj_tr[e].active(),
622 [
this, e]() {
return this->solve_layer(e); });
626 row.layer = ensemble[e].name;
633 out.
rows.push_back(row);
640 if (m.empty())
continue;
656 const bool upper = opt.method !=
"mw.lower";
661 const std::size_t N = lqn.nidx;
662 auto blank = [&](std::vector<T>& v, std::vector<bool>& d) {
663 v.assign(N + 1, Tzero());
664 d.assign(N + 1,
false);
672 for (std::size_t i = 1; i <= N; ++i) {
673 if (b.defined_T[i]) {
674 s.
TN[i] = upper ? b.TN_up[i] : b.TN_lo[i];
677 if (b.defined_U[i]) {
678 s.
UN[i] = upper ? b.UN_up[i] : b.UN_lo[i];
688 std::size_t
nlayers()
const {
return ensemble.size(); }
689 const std::vector<qn::Layer<T>>&
layers()
const {
return ensemble; }
701 const std::size_t E = ensemble.size();
706 for (std::size_t e = 0; e < E; ++e) {
709 b.
msz[e] = ensemble[e].nstations;
710 b.
ksz[e] = ensemble[e].nclasses;
740 const std::size_t E = ensemble.size();
741 layer_tran_init.assign(E, std::vector<double>());
742 if (n.
empty())
return;
746 "SolverLN::init_from_marginal: the marginal is " + std::to_string(n.
rows()) +
"x" +
747 std::to_string(n.
cols()) +
" where the block-diagonal union of the layers is " +
748 std::to_string(b.
M) +
"x" + std::to_string(b.
K) +
749 "; the aggregate view is layerBlocks', not the LQN element count");
750 for (std::size_t e = 0; e < E; ++e) {
753 std::vector<double> y(lay.
nstates, 0.0);
755 for (std::size_t i = 0; i < L.
nstations && i < b.
msz[e]; ++i)
756 for (std::size_t r = 0; r < L.
nclasses && r < b.
ksz[e]; ++r) {
757 if (!lay.
enabled[i][r])
continue;
758 const double q = std::max(0.0, n(b.
roff[e] + i, b.
coff[e] + r));
759 y[lay.
qidx[i][r]] = q;
760 if (q > 0.0) any =
true;
766 if (any) layer_tran_init[e] = y;
813 std::vector<qn::Layer<T>> ensemble;
814 std::vector<long> idxhash;
815 std::vector<bool> ignore;
819 std::size_t idx, aidx, node, cls;
821 std::vector<UpdRow> servt_map, thinkt_map, call_map, actthinkt_map;
831 std::vector<UpdRow> arv_call_map;
834 std::size_t idx, tidx_caller, eidx, nodefrom, nodeto, cfrom, cto;
836 std::vector<RouteRow> route_map;
837 std::vector<std::size_t> unique_route_idx;
838 std::vector<std::size_t> route_reset, svc_reset;
841 std::string interlockMethod =
"ilrate";
849 std::size_t idx, node, cls;
853 std::size_t fromentry;
855 std::vector<RefPathRow> refpath_map;
857 std::vector<std::array<std::size_t, 3>> refpath_keys;
858 std::vector<std::size_t> refpath_group;
860 std::vector<double> refpath_stage_prev;
861 std::vector<T> refpath_stage_prev_v;
867 std::map<std::size_t, std::vector<std::size_t>> layer_normclass;
869 std::vector<T> entryvisits;
871 Matrix<double> njobs;
874 std::vector<T> servt, residt, tput, util, thinkt;
875 std::vector<T> callservt, callresidt;
882 bool has_phase2 =
false;
883 std::vector<T> servt_ph1, servt_ph2, prOvertake;
885 std::vector<T> entry_servt_last, entry_servt_ph2;
887 std::set<std::size_t> single_replica_tasks;
888 std::vector<Distrib<T>> servtproc, thinkproc, thinktproc, tputproc, callservtproc;
894 std::vector<fluid::FluidPassage> servtcdf, callservtcdf;
895 std::vector<LnCdf> entrycdfrespt;
896 std::vector<mam::AphPair<T>> entryproc;
898 std::vector<std::vector<std::vector<fluid::FluidPassage>>> cdf_repo;
899 std::vector<double> servt_prev, residt_prev, tput_prev, thinkt_prev;
900 std::vector<double> callservt_prev, callresidt_prev;
901 std::vector<T> servt_prev_v, residt_prev_v, tput_prev_v, thinkt_prev_v, callservt_prev_v;
902 Matrix<T> servtmatrix;
905 Matrix<T> il_all, il_ph1;
906 std::vector<std::vector<std::size_t>> il_common_entries, il_src_all, il_src_ph2;
912 std::vector<std::vector<std::vector<double>>> layer_interlock;
913 std::vector<double> il_num_sources;
915 std::vector<std::vector<LayerResult<T>>> results;
923 std::vector<LayerResult<T>> results_raw_last;
924 bool results_raw_valid =
false;
925 std::vector<Matrix<T>> layer_init_sol;
930 std::vector<std::vector<double>> layer_tran_init;
942 std::vector<mva::FjMmt<T>> fj_tr;
943 std::vector<std::vector<T>> fj_lambda;
944 double relax_omega = 1.0;
946 std::shared_ptr<LnStochController<T>> stoch_ctl;
947 long averagingstart = -1;
948 bool hasconverged =
false;
950 bool moment_pass_done =
false;
951 std::vector<double> maxitererr;
952 int iterations_done = 0;
953 bool did_converge =
false;
955 std::size_t NT()
const {
return lqn.tshift + lqn.ntasks; }
957 static T Tzero() {
return num_traits<T>::from_int(0); }
958 static T Tone() {
return num_traits<T>::from_int(1); }
959 static double dbl(
const T& x) {
return num_traits<T>::to_double(x); }
971 reject_unsupported();
974 const std::size_t N = lqn.nidx;
975 ignore.assign(N + 1,
false);
978 std::vector<long> comp(N + 1, -1);
980 for (std::size_t v = 1; v <= N; ++v) {
981 if (comp[v] >= 0)
continue;
982 std::vector<std::size_t> stack{v};
984 while (!stack.empty()) {
985 const std::size_t u = stack.back();
987 for (std::size_t w : lqn.graph.succ(u))
992 for (std::size_t w : lqn.graph.pred(u))
1001 std::vector<bool> has_ref(ncomp,
false);
1002 for (std::size_t t = 1; t <= lqn.ntasks; ++t) {
1003 const std::size_t tidx = lqn.tshift + t;
1004 if (lqn.sched[tidx] == SchedStrategy::REF) has_ref[comp[tidx]] =
true;
1006 for (std::size_t v = 1; v <= N; ++v)
1007 if (!has_ref[comp[v]]) ignore[v] =
true;
1014 for (std::size_t i = 1; i <= N; ++i) {
1015 servtproc[i] = lqn.hostdem[i];
1016 thinkproc[i] = lqn.think[i];
1019 for (std::size_t c = 1; c <= lqn.ncalls; ++c)
1020 callservtproc[c] = lqn.hostdem[lqn.callpair_dst[c]];
1024 entryvisits.assign(N + 1, Tone());
1025 refpath_map.clear();
1026 layer_normclass.clear();
1028 njobs = Matrix<double>(NT() + 1, NT() + 1, 0.0);
1032 std::vector<bool> rr(NT() + 1,
false), sr(NT() + 1,
false);
1033 for (
const RouteRow& r : route_map) rr[r.idx] =
true;
1034 for (
const UpdRow& r : thinkt_map) sr[r.idx] =
true;
1035 for (
const UpdRow& r : call_map) sr[r.idx] =
true;
1041 for (
const RefPathRow& r : refpath_map) sr[r.idx] =
true;
1044 for (
const UpdRow& r : arv_call_map) rr[r.idx] =
true;
1045 for (std::size_t i = 1; i <= NT(); ++i) {
1046 if (rr[i] && idxhash[i] >= 0) route_reset.push_back(std::size_t(idxhash[i]));
1047 if (sr[i] && idxhash[i] >= 0) svc_reset.push_back(std::size_t(idxhash[i]));
1049 std::sort(route_reset.begin(), route_reset.end());
1050 route_reset.erase(std::unique(route_reset.begin(), route_reset.end()), route_reset.end());
1051 std::sort(svc_reset.begin(), svc_reset.end());
1052 svc_reset.erase(std::unique(svc_reset.begin(), svc_reset.end()), svc_reset.end());
1055 refpath_keys.clear();
1056 refpath_group.clear();
1057 for (
const RefPathRow& r : refpath_map) refpath_keys.push_back({r.idx, r.node, r.cls});
1058 std::sort(refpath_keys.begin(), refpath_keys.end());
1059 refpath_keys.erase(std::unique(refpath_keys.begin(), refpath_keys.end()),
1060 refpath_keys.end());
1061 for (
const RefPathRow& r : refpath_map) {
1062 const std::array<std::size_t, 3> key{r.idx, r.node, r.cls};
1063 refpath_group.push_back(std::size_t(
1064 std::lower_bound(refpath_keys.begin(), refpath_keys.end(), key) -
1065 refpath_keys.begin()));
1067 if (refpath_map.empty() && interlockMethod ==
"refpath") {
1074 interlockMethod =
"ilrate";
1076 refpath_stage_prev.assign(refpath_keys.size(), std::numeric_limits<double>::quiet_NaN());
1077 refpath_stage_prev_v.assign(refpath_keys.size(), Tzero());
1085 void detect_phase2() {
1087 for (std::size_t a = 1; a <= lqn.nacts; ++a)
1088 if (lqn.actphase[a] > 1) has_phase2 =
true;
1100 bool is_fwd_target(std::size_t eidx)
const {
1101 for (std::size_t c = 1; c <= lqn.ncalls; ++c)
1102 if (lqn.calltype[c] == CallType::FWD && lqn.callpair_dst[c] == eidx)
return true;
1119 double open_arrival_rate_of(std::size_t tidx)
const {
1120 if (lqn.isref[tidx])
return 0.0;
1121 for (std::size_t e : lqn.entriesof[tidx])
1122 if (lqn.issynccaller.any_col(e) || lqn.isasynccaller.any_col(e) || is_fwd_target(e))
1125 for (std::size_t e : lqn.entriesof[tidx]) {
1126 if (!lqn.has_arrival[e])
continue;
1127 const double m = dbl(lqn.arrival[e].mean);
1134 bool has_cache_node(std::size_t e)
const {
1135 for (std::size_t n = 0; n < ensemble[e].nodes.size(); ++n)
1136 if (ensemble[e].nodes[n].nodetype == NodeType::Cache)
return true;
1141 bool async_only_activity(std::size_t aidx)
const {
1142 if (aidx <= lqn.ashift || aidx > lqn.ashift + lqn.nacts)
return false;
1143 for (std::size_t e = 1; e <= lqn.nentries; ++e) {
1144 const std::size_t eidx = lqn.eshift + e;
1145 if (lqn.graph.get(eidx, aidx) == Tzero())
continue;
1146 return lqn.isasynccaller.any_col(eidx) && !lqn.issynccaller.any_col(eidx);
1161 void reject_unsupported()
const {}
1168 void resolve_interlock_method() {
1169 const LnInterlockChoice ch =
1171 interlockMethod = ch.method;
1172 if (!ch.why.empty()) std::cerr <<
"[LINE] Warning: " << ch.why <<
"\n";
1178 void build_layers() {
1185 const std::string requested = ln_requested_method(opt.method);
1186 assert_call_groups(requested ==
"flat.cs");
1187 if (requested ==
"flat.cs") {
1188 lnmethod =
"flat.cs";
1189 resolve_interlock_method();
1193 if (requested ==
"flat.ph") {
1198 ph_laws_ready =
false;
1199 lnmethod =
"flat.ph";
1200 build_layers_ph(
true);
1203 if (requested ==
"srvn.ph" || requested ==
"srvn") {
1204 ph_laws_ready =
false;
1205 if (requested ==
"srvn.ph" || probe_srvn_ph()) {
1206 lnmethod =
"srvn.ph";
1211 lnmethod = (requested ==
"moment3") ?
"moment3" :
"srvn.cs";
1212 resolve_interlock_method();
1213 std::vector<qn::Layer<T>> raw(NT() + 1);
1214 std::vector<bool> present(NT() + 1,
false);
1216 for (std::size_t hidx = 1; hidx <= lqn.nhosts; ++hidx) {
1217 if (ignore[hidx])
continue;
1218 build_layer(raw[hidx], {hidx}, lqn.tasksof[hidx],
true,
false);
1219 present[hidx] =
true;
1221 for (std::size_t t = 1; t <= lqn.ntasks; ++t) {
1222 const std::size_t tidx = lqn.tshift + t;
1223 if (ignore[tidx] || lqn.isref[tidx])
continue;
1224 bool any_caller = lqn.iscaller.any_row(tidx) || lqn.iscaller.any_col(tidx);
1225 if (!any_caller)
continue;
1227 std::vector<std::size_t> callers;
1228 for (std::size_t ct = 1; ct <= lqn.ntasks; ++ct) {
1229 const std::size_t c = lqn.tshift + ct;
1231 for (std::size_t e : lqn.entriesof[tidx])
1232 if (lqn.iscaller.get(c, e)) calls =
true;
1233 if (calls) callers.push_back(c);
1235 if (callers.empty())
continue;
1236 build_layer(raw[tidx], {tidx}, callers,
false,
false);
1237 present[tidx] =
true;
1240 idxhash.assign(lqn.nidx + 1, -1);
1242 for (std::size_t i = 1; i <= NT(); ++i)
1244 idxhash[i] = next++;
1245 ensemble.push_back(std::move(raw[i]));
1247 layer_init_sol.assign(ensemble.size(), Matrix<T>());
1274 void assert_call_groups(
bool flat)
const {
1275 if (lqn.callgroups.empty())
return;
1278 "Call groups routed by a routing strategy require the squashed layering; use "
1279 "method='flat'. Under srvn the targets never share a submodel, so the dispatch "
1280 "order cannot be represented.");
1281 if (opt.layer_solver !=
"ssa")
1283 "Routed call groups need a layer solver that resolves the strategy from the "
1284 "state; set layer_solver='ssa'. MVA, NC and FLD read the routing matrix, into "
1285 "which refresh_routing has expanded the strategy as a uniform split, and would "
1286 "return that split under a round-robin or JSQ label.");
1289 std::vector<std::size_t> flat_server_set()
const {
1290 std::vector<std::size_t> servers;
1291 for (std::size_t i = 1; i <= NT(); ++i) {
1292 if (lqn.repl[i] > 1.0)
1294 "Flat layering does not support replicated processors or tasks, use the "
1295 "default 'srvn' layering.");
1298 "Flat layering does not support cache tasks, use the default 'srvn' "
1300 if (lqn.hassetup[i])
1302 "Flat layering does not support setup tasks, use the default 'srvn' "
1305 for (std::size_t hidx = 1; hidx <= lqn.nhosts; ++hidx)
1306 if (!ignore[hidx] && !lqn.tasksof[hidx].empty()) servers.push_back(hidx);
1307 for (std::size_t t = 1; t <= lqn.ntasks; ++t) {
1308 const std::size_t tidx = lqn.tshift + t;
1309 if (ignore[tidx] || lqn.isref[tidx])
continue;
1310 if (!lqn.iscaller.any_row(tidx) && !lqn.iscaller.any_col(tidx))
continue;
1311 bool has_task_caller =
false;
1312 for (std::size_t eidx : lqn.entriesof[tidx])
1313 for (std::size_t c = 1; c <= lqn.ntasks; ++c)
1314 if (lqn.iscaller.get(lqn.tshift + c, eidx)) has_task_caller =
true;
1315 if (has_task_caller) servers.push_back(tidx);
1317 if (servers.empty())
1319 "Flat layering found no server: the model has no processor with tasks.");
1331 void build_flat_layer() {
1332 const std::vector<std::size_t> flat_servers = flat_server_set();
1333 std::vector<std::size_t> flat_callers;
1334 for (std::size_t t = 1; t <= lqn.ntasks; ++t) {
1335 const std::size_t tidx = lqn.tshift + t;
1336 if (!ignore[tidx]) flat_callers.push_back(tidx);
1339 build_layer(layer, flat_servers, flat_callers,
false,
true);
1341 ensemble.push_back(std::move(layer));
1342 idxhash.assign(lqn.nidx + 1, -1);
1343 for (std::size_t s : flat_servers) idxhash[s] = 0;
1344 layer_init_sol.assign(ensemble.size(), Matrix<T>());
1352 void build_fork_views() {
1353 fj_tr.assign(ensemble.size(), mva::FjMmt<T>());
1354 fj_lambda.assign(ensemble.size(), {});
1355 for (std::size_t e = 0; e < ensemble.size(); ++e) {
1356 if (!ensemble[e].has_fork())
continue;
1363 if (ensemble[e].sourceIdx != 0)
1365 "SolverLN: layer '" + ensemble[e].name +
1366 "' carries both an AND fork and an open stream (an async call or an entry "
1367 "arrival); the fork-join transform needs a Source of its own");
1369 fj_lambda[e].assign(fj_tr[e].V.classes.size() + 1,
1390 const std::shared_ptr<std::vector<std::vector<std::size_t>>>& chaincols, std::size_t R) {
1391 return [f, cols, chaincols, R](
const std::vector<T>& n) -> std::vector<T> {
1392 const std::size_t L = n.size();
1393 const std::vector<std::vector<std::size_t>>& use =
1394 (L == R || chaincols->empty()) ? cols : *chaincols;
1395 std::vector<T> nop(use.size(), num_traits<T>::from_int(0));
1396 for (std::size_t j = 0; j < use.size(); ++j)
1397 for (std::size_t k = 0; k < use[j].size(); ++k)
1398 if (use[j][k] < L) nop[j] = T(nop[j] + n[use[j][k]]);
1399 const std::vector<T> w = f(nop);
1400 std::vector<T> v(L, num_traits<T>::from_int(1));
1401 if (w.empty())
return v;
1402 for (std::size_t j = 0; j < use.size(); ++j) {
1403 const T wj = w[j < w.size() ? j : w.size() - 1];
1404 for (std::size_t k = 0; k < use[j].size(); ++k)
1405 if (use[j][k] < L) v[use[j][k]] = wj;
1412 static std::vector<T> layer_peak(
const std::vector<T>& peakPerOperand,
1413 const std::vector<std::vector<std::size_t>>& cols,
1415 std::vector<T> peak(R, num_traits<T>::from_int(1));
1416 if (peakPerOperand.empty())
return peak;
1417 for (std::size_t j = 0; j < cols.size(); ++j) {
1418 const T pj = peakPerOperand[j < peakPerOperand.size() ? j : peakPerOperand.size() - 1];
1419 for (std::size_t k = 0; k < cols[j].size(); ++k)
1420 if (cols[j][k] < R) peak[cols[j][k]] = pj;
1434 void apply_host_task_priorities(qn::Layer<T>& m,
const std::vector<std::size_t>& idxSet)
const {
1437 return s == SchedStrategy::HOL || s == SchedStrategy::FCFSPRPRIO ||
1438 s == SchedStrategy::FCFSPIPRIO || s == SchedStrategy::LCFSPRIO ||
1439 s == SchedStrategy::LCFSPRPRIO || s == SchedStrategy::LCFSPIPRIO ||
1440 s == SchedStrategy::PSPRIO || s == SchedStrategy::DPSPRIO ||
1441 s == SchedStrategy::GPSPRIO;
1443 const std::size_t R = m.classes.size();
1444 std::vector<std::size_t> owner(R + 1, 0);
1445 auto own = [&](std::size_t k, std::size_t t) {
1446 if (k >= 1 && k <= R && owner[k] == 0) owner[k] = t;
1448 for (
const auto& a : m.attr_tasks) own(a.first, a.second);
1449 for (
const auto& a : m.attr_entries) own(a.first, lqn.parent[a.second]);
1450 for (
const auto& a : m.attr_activities) own(a.first, lqn.parent[a.second]);
1451 for (
const auto& a : m.attr_calls) own(a[0], lqn.parent[a[2]]);
1452 for (std::size_t s : idxSet) {
1453 if (s < 1 || s > lqn.nhosts || !is_prio(lqn.sched[s]))
continue;
1455 for (std::size_t t = lqn.tshift + 1; t <= lqn.tshift + lqn.ntasks; ++t)
1456 if (lqn.parent[t] == s) pmax = std::max(pmax, lqn.prio[t]);
1457 for (std::size_t k = 1; k <= R; ++k)
1458 if (owner[k] != 0 && lqn.parent[owner[k]] == s)
1459 m.classes[k - 1].prio = pmax - lqn.prio[owner[k]];
1472 void build_layer(qn::Layer<T>& m,
const std::vector<std::size_t>& idxSet,
1473 const std::vector<std::size_t>& callers,
bool ishostlayer,
bool flat) {
1474 const T one = Tone();
1478 const std::size_t idx = idxSet[0];
1479 m.name = flat ? lqn.hashnames[idx] +
".Flat" : lqn.hashnames[idx];
1482 const double rawrepl = lqn.repl[idx];
1483 std::size_t nreplicas = 1;
1484 if (!flat && rawrepl > 1.0 && !callers.empty()) {
1493 for (std::size_t c : callers)
1494 if (lqn.repl[c] != rawrepl) reduce =
false;
1496 for (std::size_t c : callers)
1497 if (lqn.fanout_at(c, idx) < rawrepl) reduce =
false;
1499 nreplicas = reduce ? 1 :
static_cast<std::size_t
>(std::llround(rawrepl));
1500 if (reduce && !ishostlayer) single_replica_tasks.insert(idx);
1502 const bool reduce_fanout = (nreplicas == 1 && rawrepl > 1.0 && !callers.empty());
1503 const std::vector<double>& mult = lqn.maxmult;
1507 m.clientIdx = m.add_station(qn::Station<T>{
"Clients", NodeType::Delay, SchedStrategy::INF,
1508 std::numeric_limits<double>::infinity(),
false, 0});
1509 const std::size_t clientNode = m.node_of_station(m.clientIdx);
1510 m.serverIdx = m.clientIdx + 1;
1514 std::vector<std::vector<std::size_t>> srv(lqn.nidx + 1);
1515 std::vector<std::vector<std::size_t>> srvnode(lqn.nidx + 1);
1516 m.server_idx_of.assign(lqn.nidx + 1, 0);
1517 for (std::size_t sidx : idxSet) {
1518 const bool sishost = sidx <= lqn.nhosts;
1519 srv[sidx].resize(nreplicas);
1520 for (std::size_t r = 0; r < nreplicas; ++r) {
1522 st.name = r == 0 ? lqn.hashnames[sidx]
1523 : lqn.hashnames[sidx] +
"." + std::to_string(r + 1);
1526 lqn.sched[sidx] == SchedStrategy::INF ? NodeType::Delay : NodeType::Queue;
1527 st.sched = lqn.sched[sidx];
1529 st.nservers = lqn.sched[sidx] == SchedStrategy::INF
1530 ? std::numeric_limits<double>::infinity()
1532 st.attr_ishost = flat ? sishost : ishostlayer;
1534 srv[sidx][r] = m.add_station(st);
1535 srvnode[sidx].push_back(m.node_of_station(srv[sidx][r]));
1537 m.server_idx_of[sidx] = srv[sidx][0];
1539 m.host_stations.push_back(srv[sidx][0]);
1541 m.task_stations.push_back(srv[sidx][0]);
1544 const std::vector<std::size_t>& server = srv[idx];
1545 const std::vector<std::size_t>& serverNode = srvnode[idx];
1548 auto servers_for = [&](std::size_t elem) ->
const std::vector<std::size_t>& {
1549 static const std::vector<std::size_t> none;
1550 if (elem >= 1 && elem <= lqn.nidx && !srv[elem].empty())
return srv[elem];
1553 auto server_nodes_for = [&](std::size_t elem) ->
const std::vector<std::size_t>& {
1554 static const std::vector<std::size_t> none;
1555 if (elem >= 1 && elem <= lqn.nidx && !srvnode[elem].empty())
return srvnode[elem];
1559 auto host_is_server = [&](std::size_t tidx_) {
1560 return !servers_for(lqn.parent[tidx_]).empty();
1563 auto is_layer_client = [&](std::size_t tidx_) {
1564 for (std::size_t s : idxSet)
1565 for (std::size_t e : lqn.entriesof[s])
1566 if (lqn.issynccaller.get(tidx_, e))
return true;
1570 auto set_all_servers = [&](std::size_t cl,
const Distrib<T>& d) {
1571 for (std::size_t s : idxSet)
1572 for (std::size_t st : srv[s]) m.set_service(st, cl, d);
1580 std::vector<std::size_t> group_of_call(lqn.ncalls + 1, 0);
1581 std::vector<std::vector<std::size_t>> group_members(lqn.callgroups.size() + 1);
1582 for (std::size_t g = 0; g < lqn.callgroups.size(); ++g) {
1583 const LqnCallGroup& grp = lqn.callgroups[g];
1584 for (std::size_t tgt : grp.
targets)
1585 for (std::size_t cidx : lqn.callsof[grp.
caller])
1586 if (lqn.callpair_dst[cidx] == tgt && lqn.calltype[cidx] == CallType::SYNC &&
1587 group_of_call[cidx] == 0) {
1588 group_of_call[cidx] = g + 1;
1589 group_members[g + 1].push_back(cidx);
1594 std::vector<std::size_t> grp_router(lqn.callgroups.size() + 1, 0);
1595 std::vector<std::size_t> grp_dispatch(lqn.callgroups.size() + 1, 0);
1596 std::vector<std::size_t> grp_class(lqn.callgroups.size() + 1, 0);
1599 std::vector<std::array<std::size_t, 3>> routed_group_sites;
1602 std::vector<std::size_t> acts_in_caller;
1603 for (std::size_t c : callers)
1604 for (std::size_t a : lqn.actsof[c]) acts_in_caller.push_back(a);
1605 bool hasfork =
false, hasjoin =
false;
1606 std::size_t maxfanout = 1;
1607 for (std::size_t a : acts_in_caller) {
1608 if (lqn.actposttype[a] == PrecedenceType::POST_AND) hasfork =
true;
1609 if (lqn.actpretype[a] == PrecedenceType::PRE_AND) hasjoin =
true;
1610 std::size_t nand = 0;
1611 for (std::size_t sx : lqn.graph.succ(a))
1612 if (lqn.actposttype[sx] == PrecedenceType::POST_AND) ++nand;
1613 if (nand > maxfanout) maxfanout = nand;
1615 std::size_t forkNode = 0, joinNode = 0, joinStation = 0;
1616 std::vector<std::size_t> forkRouter;
1618 forkNode = m.add_node(
"Fork_PostAnd", NodeType::Fork,
false);
1619 for (std::size_t f = 1; f <= maxfanout; ++f)
1620 forkRouter.push_back(
1621 m.add_node(
"Fork_PostAnd_" + std::to_string(f), NodeType::Router,
true));
1625 js.name =
"Join_PreAnd";
1626 js.nodetype = NodeType::Join;
1627 js.sched = SchedStrategy::INF;
1628 js.nservers = std::numeric_limits<double>::infinity();
1629 js.attr_ishost =
false;
1631 joinStation = m.add_station(js);
1632 joinNode = m.node_of_station(joinStation);
1633 if (forkNode) m.fj.emplace_back(forkNode, joinNode);
1640 bool iscachelayer = !flat && ishostlayer && !callers.empty();
1641 for (std::size_t c : callers)
1642 if (!lqn.iscache[c]) iscachelayer =
false;
1651 std::size_t cacheNode = 0;
1652 qn::CacheParam<T> cachepar;
1654 const std::size_t ct = callers[0];
1655 cachepar.nitems = lqn.nitems[ct];
1656 cachepar.itemcap = lqn.itemcap[ct];
1657 cachepar.replacestrat = lqn.replacestrat[ct];
1658 cacheNode = m.add_node(lqn.hashnames[ct], NodeType::Cache,
true);
1668 std::vector<std::size_t> async_here;
1669 for (std::size_t c : callers)
1670 for (std::size_t aidx : lqn.actsof[c])
1671 for (std::size_t cidx : lqn.callsof[aidx])
1672 if (lqn.calltype[cidx] == CallType::ASYNC &&
1673 !servers_for(lqn.parent[lqn.callpair_dst[cidx]]).empty())
1674 async_here.push_back(cidx);
1675 std::vector<std::size_t> open_entries;
1676 for (std::size_t c : callers) {
1684 for (std::size_t eidx : lqn.entriesof[c])
1685 if (lqn.has_arrival[eidx]) open_entries.push_back(eidx);
1688 std::size_t sourceStation = 0, sourceNode = 0, sinkNode = 0;
1689 if (!async_here.empty() || !open_entries.empty()) {
1691 src.name =
"Source";
1692 src.nodetype = NodeType::Source;
1694 src.sched = SchedStrategy::EXT;
1696 src.attr_ishost =
false;
1698 sourceStation = m.add_station(src);
1699 sourceNode = m.node_of_station(sourceStation);
1700 sinkNode = m.add_node(
"Sink", NodeType::Sink,
false);
1701 m.sourceIdx = sourceStation;
1702 m.sinkNode = sinkNode;
1706 std::vector<std::size_t> forkClassStack;
1709 std::vector<std::size_t> cls(lqn.nidx + 1, 0);
1710 std::vector<std::size_t> callcls(lqn.ncalls + 1, 0);
1713 std::vector<std::size_t> auxcallcls(lqn.ncalls + 1, 0);
1721 auto caller_needs_class = [&](std::size_t tidx_caller) {
1722 if (host_is_server(tidx_caller)) {
1723 if (lqn.isref[tidx_caller])
return true;
1724 for (std::size_t e : lqn.entriesof[tidx_caller])
1725 if (lqn.issynccaller.any_col(e) || lqn.isasynccaller.any_col(e) ||
1726 is_fwd_target(e) || lqn.has_arrival[e])
1729 return is_layer_client(tidx_caller);
1731 auto caller_acts_visible = [&](std::size_t tidx_caller) {
1732 return host_is_server(tidx_caller) || is_layer_client(tidx_caller);
1740 auto caller_pool = [&](std::size_t tidx_) ->
double {
1741 if (njobs(tidx_, idx) != 0.0)
return njobs(tidx_, idx);
1743 const bool single_replica = reduce_fanout || single_replica_tasks.count(tidx_) > 0;
1744 double nj = single_replica ? mult[tidx_] : mult[tidx_] * lqn.repl[tidx_];
1745 if (std::isinf(nj)) {
1747 for (std::size_t c = 1; c <= NT(); ++c)
1748 if (lqn.taskgraph.get(c, tidx_) != Tzero()) s += mult[c];
1750 if (std::isinf(nj)) {
1752 for (std::size_t c = 1; c <= NT(); ++c)
1753 if (std::isfinite(mult[c])) s2 += mult[c] * lqn.repl[c];
1754 nj = std::min(s2, 1000.0);
1770 std::vector<bool> rp_member(lqn.nidx + 1,
false);
1771 std::vector<bool> rp_spliced(lqn.nidx + 1,
false);
1772 std::vector<bool> rp_hop(lqn.nidx + 1,
false);
1773 std::vector<bool> rp_expand(lqn.ncalls + 1,
false);
1774 std::vector<std::size_t> rp_child_task(lqn.ncalls + 1, 0);
1775 std::vector<std::size_t> rp_child_entry(lqn.ncalls + 1, 0);
1776 std::vector<T> rp_ret_share(lqn.ncalls + 1, Tzero());
1777 std::vector<T> rp_call_mean(lqn.ncalls + 1, Tzero());
1778 std::vector<std::vector<std::size_t>> rp_hop_children(lqn.nidx + 1);
1779 std::vector<double> rp_hop_seed(lqn.nidx + 1, 0.0);
1780 std::vector<std::vector<std::size_t>> rp_member_inbound(lqn.nidx + 1);
1781 std::vector<std::vector<std::size_t>> rp_hop_inbound(lqn.nidx + 1);
1782 std::vector<std::vector<std::size_t>> rp_group_roots;
1783 std::vector<std::size_t> rp_group_reftask;
1784 std::vector<bool> rp_group_head_is_caller;
1786 std::vector<std::size_t> rp_ref_stage;
1787 std::vector<std::size_t> rp_hop_cls(lqn.nidx + 1, 0), rp_hop_ret(lqn.nidx + 1, 0),
1788 rp_ret_cls(lqn.nidx + 1, 0);
1789 std::vector<std::size_t> rp_gate(lqn.ncalls + 1, 0), rp_aux(lqn.ncalls + 1, 0),
1790 rp_resume(lqn.ncalls + 1, 0);
1793 auto host_demand_of_entry = [&](std::size_t eidx_) {
1795 for (std::size_t a_ : lqn.actsof[eidx_]) {
1796 if (lqn.hostdem[a_].disabled)
continue;
1797 const double m_ = dbl(lqn.hostdem[a_].mean);
1798 if (std::isfinite(m_)) d += m_;
1802 const auto in_idxset = [&](std::size_t t_) {
1803 return std::find(idxSet.begin(), idxSet.end(), t_) != idxSet.end();
1810 if (interlockMethod !=
"refpath" || flat || iscachelayer || hasfork || hasjoin)
return;
1813 if (reduce_fanout || nreplicas > 1)
return;
1814 for (std::size_t t_ : callers)
1815 if (lqn.repl[t_] > 1.0 || single_replica_tasks.count(t_) > 0)
return;
1816 const lqn::LqnRefRoutes<T> R =
1818 if (!R.why.empty()) {
1820 std::cerr <<
"[LINE] Warning: layer '" << lqn.hashnames[idx]
1821 <<
"' falls back to the 'ilrate' interlock method: " << R.why
1825 std::vector<std::size_t> keep;
1826 for (std::size_t g_ = 0; g_ < R.groups.size(); ++g_) {
1827 const lqn::LqnRefGroup<T>& G = R.groups[g_];
1828 bool wanted = G.members.size() >= 2;
1831 if (!wanted && opt.interlock_refpath_scope ==
"all" && !G.head_is_caller &&
1834 if (!wanted)
continue;
1839 for (std::size_t t_ : G.members) pool += caller_pool(t_);
1840 const double chainPop = mult[G.reftask] * lqn.repl[G.reftask];
1841 if (!(pool > chainPop))
continue;
1844 if (keep.empty())
return;
1847 std::vector<bool> sMember(lqn.nidx + 1,
false), sSpliced(lqn.nidx + 1,
false),
1848 sHop(lqn.nidx + 1,
false), sExpand(lqn.ncalls + 1,
false);
1849 std::vector<std::size_t> sChildTask(lqn.ncalls + 1, 0),
1850 sChildEntry(lqn.ncalls + 1, 0);
1851 std::vector<T> sRetShare(lqn.ncalls + 1, Tzero()), sCallMean(lqn.ncalls + 1, Tzero());
1852 std::vector<std::vector<std::size_t>> sHopChildren(lqn.nidx + 1),
1853 sMemberInbound(lqn.nidx + 1), sHopInbound(lqn.nidx + 1);
1854 std::vector<double> sHopSeed(lqn.nidx + 1, 0.0);
1855 std::vector<std::vector<std::size_t>> sRoots;
1856 std::vector<std::size_t> sRefTask;
1857 std::vector<bool> sHeadIsCaller;
1858 for (std::size_t g_ : keep) {
1859 const lqn::LqnRefGroup<T>& G = R.groups[g_];
1860 const std::size_t nE = G.entries.size();
1862 std::vector<bool> leads(G.ismember.begin(), G.ismember.end());
1863 for (std::size_t i_ = nE; i_-- > 0;)
1864 for (
const lqn::LqnRefCall<T>& c_ : G.calls)
1865 if (c_.from == i_ && leads[c_.to]) leads[i_] =
true;
1866 for (std::size_t t_ : G.members)
1867 if (!caller_needs_class(t_))
return;
1868 if (G.head_is_caller && !caller_needs_class(G.reftask))
return;
1870 if (!G.head_is_caller && !std::isfinite(lqn.maxmult[G.reftask]))
return;
1871 for (std::size_t i_ = 0; i_ < nE; ++i_) {
1872 if (!leads[i_] || G.ismember[i_])
continue;
1874 if (in_idxset(G.etask[i_]) || lqn.repl[G.etask[i_]] > 1.0)
return;
1876 for (
const lqn::LqnRefCall<T>& c_ : G.calls) {
1877 if (!leads[c_.to])
continue;
1878 if (in_idxset(G.etask[c_.to]))
return;
1881 const std::size_t head = G.reftask;
1882 for (std::size_t t_ : G.members) {
1884 sSpliced[t_] = !(G.head_is_caller && t_ == head);
1886 std::vector<T> inTotMember(lqn.nidx + 1, Tzero()), inTotHop(lqn.nidx + 1, Tzero());
1887 for (
const lqn::LqnRefCall<T>& c_ : G.calls) {
1888 if (!leads[c_.to])
continue;
1889 const std::size_t cidx_ = c_.cidx;
1890 const std::size_t fromEidx = G.entries[c_.from];
1891 sExpand[cidx_] =
true;
1892 sCallMean[cidx_] = lqn.callproc_mean[cidx_];
1893 sHopChildren[fromEidx].push_back(cidx_);
1894 if (G.ismember[c_.to]) {
1895 const std::size_t tv = G.etask[c_.to];
1896 sChildTask[cidx_] = tv;
1897 sMemberInbound[tv].push_back(cidx_);
1898 inTotMember[tv] = T(inTotMember[tv] + c_.vcall);
1900 const std::size_t ev = G.entries[c_.to];
1901 sChildEntry[cidx_] = ev;
1902 sHopInbound[ev].push_back(cidx_);
1903 inTotHop[ev] = T(inTotHop[ev] + c_.vcall);
1905 sRetShare[cidx_] = c_.vcall;
1907 for (
const lqn::LqnRefCall<T>& c_ : G.calls) {
1908 const std::size_t cidx_ = c_.cidx;
1909 if (!sExpand[cidx_])
continue;
1910 const T tot = sChildTask[cidx_] > 0 ? inTotMember[sChildTask[cidx_]]
1911 : inTotHop[sChildEntry[cidx_]];
1913 ? T(sRetShare[cidx_] / tot)
1916 for (std::size_t i_ = 0; i_ < nE; ++i_)
1917 if (leads[i_] && !G.ismember[i_]) {
1918 const std::size_t ev = G.entries[i_];
1920 sHopSeed[ev] = host_demand_of_entry(ev);
1922 if (!G.head_is_caller) {
1926 for (std::size_t ev : lqn.entriesof[G.reftask])
1929 sHopSeed[ev] = host_demand_of_entry(ev);
1932 sRoots.push_back(lqn.entriesof[head]);
1933 sRefTask.push_back(head);
1934 sHeadIsCaller.push_back(G.head_is_caller);
1936 rp_member = sMember;
1937 rp_spliced = sSpliced;
1939 rp_expand = sExpand;
1940 rp_child_task = sChildTask;
1941 rp_child_entry = sChildEntry;
1942 rp_ret_share = sRetShare;
1943 rp_call_mean = sCallMean;
1944 rp_hop_children = sHopChildren;
1945 rp_hop_seed = sHopSeed;
1946 rp_member_inbound = sMemberInbound;
1947 rp_hop_inbound = sHopInbound;
1948 rp_group_roots = sRoots;
1949 rp_group_reftask = sRefTask;
1950 rp_group_head_is_caller = sHeadIsCaller;
1955 for (std::size_t tidx_caller : callers) {
1956 if (caller_needs_class(tidx_caller)) {
1957 double nj = njobs(tidx_caller, idx);
1959 nj = caller_pool(tidx_caller);
1960 njobs(tidx_caller, idx) = nj;
1965 const bool rp_spliced_here = rp_on && rp_spliced[tidx_caller];
1967 jc.name = lqn.hashnames[tidx_caller];
1968 jc.type = JobClassType::CLOSED;
1969 jc.population = rp_spliced_here ? 0.0 : nj;
1970 jc.refstat = m.clientIdx;
1971 jc.completes =
false;
1974 jc.is_ref_class = !rp_spliced_here;
1975 jc.attr_kind = int(LqnElement::TASK);
1976 jc.attr_idx = tidx_caller;
1977 cls[tidx_caller] = m.add_class(jc);
1978 m.attr_tasks.emplace_back(cls[tidx_caller], tidx_caller);
1979 if (lqn.isref[tidx_caller]) {
1980 m.set_service(m.clientIdx, cls[tidx_caller], thinkproc[tidx_caller]);
1981 }
else if (rp_spliced_here) {
1990 thinkt_map.push_back({idx, tidx_caller, m.clientIdx, cls[tidx_caller]});
1993 for (std::size_t eidx : lqn.entriesof[tidx_caller]) {
1995 ec.name = lqn.hashnames[eidx];
1996 ec.type = JobClassType::CLOSED;
1997 ec.population = 0.0;
1998 ec.refstat = m.clientIdx;
1999 ec.completes =
false;
2000 ec.attr_kind = int(LqnElement::ENTRY);
2002 cls[eidx] = m.add_class(ec);
2003 m.attr_entries.emplace_back(cls[eidx], eidx);
2008 for (std::size_t aidx : lqn.actsof[tidx_caller]) {
2009 if (caller_acts_visible(tidx_caller)) {
2011 ac.name = lqn.hashnames[aidx];
2012 ac.type = JobClassType::CLOSED;
2013 ac.population = 0.0;
2014 ac.refstat = m.clientIdx;
2015 ac.completes =
false;
2016 ac.attr_kind = int(LqnElement::ACTIVITY);
2018 cls[aidx] = m.add_class(ac);
2019 m.attr_activities.emplace_back(cls[aidx], aidx);
2023 const std::size_t hidx = lqn.parent[lqn.parent[aidx]];
2024 if (servers_for(hidx).empty())
2025 m.set_service(m.clientIdx, cls[aidx], servtproc[aidx]);
2027 for (std::size_t cidx : lqn.callsof[aidx]) {
2028 if (lqn.calltype[cidx] == CallType::ASYNC) {
2032 const std::size_t adst = lqn.parent[lqn.callpair_dst[cidx]];
2033 if (servers_for(adst).empty())
continue;
2035 oc.name = lqn.callhashnames[cidx];
2036 oc.type = JobClassType::OPEN;
2037 oc.population = std::numeric_limits<double>::infinity();
2038 oc.refstat = sourceStation;
2039 oc.completes =
false;
2040 oc.is_ref_class =
false;
2041 oc.attr_kind = int(LqnElement::CALL);
2043 callcls[cidx] = m.add_class(oc);
2044 m.attr_calls.push_back({callcls[cidx], cidx, lqn.callpair_src[cidx],
2045 lqn.callpair_dst[cidx]});
2053 m.set_service(sourceStation, callcls[cidx],
2056 T minRespTA = Tzero();
2057 for (std::size_t ta : lqn.actsof[adst]) minRespTA += lqn.hostdem[ta].mean;
2058 for (std::size_t st : servers_for(adst)) {
2060 call_map.push_back({idx, cidx, st, callcls[cidx]});
2062 arv_call_map.push_back({idx, cidx, sourceStation, callcls[cidx]});
2065 if (lqn.calltype[cidx] != CallType::SYNC)
continue;
2066 const std::size_t gid = group_of_call[cidx];
2075 if (grp_dispatch[gid] == 0) {
2081 const std::string tag =
2082 lqn.hashnames[aidx] +
".Dispatch" + std::to_string(gid);
2083 grp_router[gid] = m.add_node(tag +
".Router", NodeType::Router,
true);
2086 dc.type = JobClassType::CLOSED;
2087 dc.population = 0.0;
2088 dc.refstat = m.clientIdx;
2089 dc.completes =
false;
2090 dc.attr_kind = int(LqnElement::CALL);
2092 grp_dispatch[gid] = m.add_class(dc);
2095 gc.name = lqn.callhashnames[cidx] +
".Group" + std::to_string(gid);
2096 gc.type = JobClassType::CLOSED;
2097 gc.population = 0.0;
2098 gc.refstat = m.clientIdx;
2099 gc.completes =
false;
2100 gc.attr_kind = int(LqnElement::CALL);
2102 grp_class[gid] = m.add_class(gc);
2104 routed_group_sites.push_back(
2105 {grp_router[gid], grp_dispatch[gid], gid});
2107 callcls[cidx] = grp_dispatch[gid];
2108 m.attr_calls.push_back({callcls[cidx], cidx, lqn.callpair_src[cidx],
2109 lqn.callpair_dst[cidx]});
2110 for (std::size_t st2 : servers_for(lqn.parent[lqn.callpair_dst[cidx]]))
2111 m.set_service(st2, callcls[cidx], callservtproc[cidx]);
2115 cc.name = lqn.callhashnames[cidx];
2116 cc.type = JobClassType::CLOSED;
2117 cc.population = 0.0;
2118 cc.refstat = m.clientIdx;
2119 cc.completes =
false;
2120 cc.attr_kind = int(LqnElement::CALL);
2122 callcls[cidx] = m.add_class(cc);
2123 m.attr_calls.push_back({callcls[cidx], cidx, lqn.callpair_src[cidx],
2124 lqn.callpair_dst[cidx]});
2132 const std::size_t seedidx = flat ? lqn.parent[lqn.callpair_dst[cidx]] : idx;
2133 T minRespT = Tzero();
2134 for (std::size_t ta : lqn.actsof[seedidx]) minRespT += lqn.hostdem[ta].mean;
2135 for (std::size_t st : servers_for(seedidx))
2138 if (lqn.callproc_mean[cidx] !=
2139 num_traits<T>::from_double(
double(nreplicas))) {
2141 xc.name = lqn.callhashnames[cidx] +
".Aux";
2142 xc.type = JobClassType::CLOSED;
2143 xc.population = 0.0;
2144 xc.refstat = m.clientIdx;
2145 xc.completes =
false;
2146 xc.attr_kind = int(LqnElement::CALL);
2148 auxcallcls[cidx] = m.add_class(xc);
2160 auto rp_class = [&](
const std::string& nm,
double pop, std::size_t attr) {
2163 rc.type = JobClassType::CLOSED;
2164 rc.population = pop;
2165 rc.refstat = m.clientIdx;
2166 rc.completes =
false;
2167 rc.is_ref_class =
false;
2170 const std::size_t k_ = m.add_class(rc);
2175 for (std::size_t g_ = 0; g_ < rp_group_reftask.size(); ++g_) {
2176 if (rp_group_head_is_caller[g_]) {
2177 rp_ref_stage.push_back(0);
2180 const std::size_t rt = rp_group_reftask[g_];
2181 const std::size_t k_ =
2182 rp_class(lqn.hashnames[rt] +
".RefPath", mult[rt] * lqn.repl[rt], rt);
2183 m.classes[k_ - 1].is_ref_class =
true;
2184 m.set_service(m.clientIdx, k_, thinkproc[rt]);
2185 rp_ref_stage.push_back(k_);
2187 for (std::size_t ev = 1; ev <= lqn.nidx; ++ev) {
2188 if (!rp_hop[ev])
continue;
2189 const std::size_t k_ = rp_class(lqn.hashnames[ev] +
".RefHop", 0.0, ev);
2191 m.set_service(m.clientIdx, k_,
2193 rp_hop_cls[ev] = k_;
2196 refpath_map.push_back({idx, m.clientIdx, k_, 1, ev, 1.0, 0});
2197 for (std::size_t c_ : rp_hop_children[ev])
2198 refpath_map.push_back({idx, m.clientIdx, k_, 2, c_, -1.0, ev});
2199 refpath_map.push_back({idx, m.clientIdx, k_, 3, lqn.parent[ev], 1.0, 0});
2200 rp_hop_ret[ev] = rp_class(lqn.hashnames[ev] +
".RefHopRet", 0.0, ev);
2202 for (std::size_t tv = 1; tv <= lqn.nidx; ++tv)
2203 if (rp_spliced[tv]) rp_ret_cls[tv] = rp_class(lqn.hashnames[tv] +
".RefRet", 0.0, tv);
2204 for (std::size_t cidx_ = 1; cidx_ <= lqn.ncalls; ++cidx_) {
2205 if (!rp_expand[cidx_])
continue;
2206 if (callcls[cidx_] != 0) {
2207 rp_gate[cidx_] = callcls[cidx_];
2209 if (auxcallcls[cidx_] != 0) rp_aux[cidx_] = auxcallcls[cidx_];
2211 rp_gate[cidx_] = rp_class(lqn.callhashnames[cidx_] +
".RefGate", 0.0, cidx_);
2212 if (rp_call_mean[cidx_] != Tone())
2214 rp_class(lqn.callhashnames[cidx_] +
".RefGate.Aux", 0.0, cidx_);
2219 if (rp_child_task[cidx_] > 0)
2220 refpath_map.push_back(
2221 {idx, m.clientIdx, rp_gate[cidx_], 3, rp_child_task[cidx_], 1.0, 0});
2222 rp_resume[cidx_] = rp_class(lqn.callhashnames[cidx_] +
".RefResume", 0.0, cidx_);
2228 std::size_t curclass;
2234 std::vector<std::size_t> curnodes;
2236 const int atClient = 1, atServer = 2, atCache = 3;
2237 std::vector<int> jobposkey(lqn.nidx + 1, atClient);
2238 std::vector<std::size_t> curclasskey(lqn.nidx + 1, 0);
2239 std::vector<std::vector<std::size_t>> curnodeskey(lqn.nidx + 1);
2246 auto rp_descent = [&](std::size_t cidx_, Ctx st) -> Ctx {
2247 const std::size_t gate = rp_gate[cidx_], res = rp_resume[cidx_], aux = rp_aux[cidx_];
2248 const T m_ = rp_call_mean[cidx_];
2249 const std::size_t enter = rp_child_task[cidx_] > 0 ? cls[rp_child_task[cidx_]]
2250 : rp_hop_cls[rp_child_entry[cidx_]];
2251 if (st.jobpos == atClient)
2252 m.set_route(st.curclass, gate, clientNode, clientNode, Tone());
2254 for (std::size_t nd : st.curnodes) m.set_route(st.curclass, gate, nd, clientNode, Tone());
2257 m.set_route(gate, enter, clientNode, clientNode, Tone());
2259 }
else if (m_ > Tone()) {
2260 m.set_route(gate, enter, clientNode, clientNode, Tone());
2261 m.set_route(res, enter, clientNode, clientNode, T(Tone() - Tone() / m_));
2262 m.set_route(res, aux, clientNode, clientNode, T(Tone() / m_));
2265 m.set_route(gate, enter, clientNode, clientNode, m_);
2266 m.set_route(gate, aux, clientNode, clientNode, T(Tone() - m_));
2267 m.set_route(res, aux, clientNode, clientNode, Tone());
2270 st.jobpos = atClient;
2271 st.curnodes.clear();
2276 std::function<Ctx(std::size_t, std::size_t, Ctx)> recur =
2277 [&](std::size_t tidx_caller, std::size_t aidx, Ctx st) -> Ctx {
2278 jobposkey[aidx] = st.jobpos;
2279 curclasskey[aidx] = st.curclass;
2280 curnodeskey[aidx] = st.curnodes;
2281 const std::vector<std::size_t> nexts = lqn.graph.succ(aidx);
2282 std::size_t lastEntryClass = st.curclass;
2284 bool next_is_fork =
false;
2285 for (std::size_t sx : nexts)
2286 if (lqn.actposttype[sx] == PrecedenceType::POST_AND) next_is_fork =
true;
2287 const Ctx preFork = st;
2288 std::vector<std::size_t> andSuccs;
2289 for (std::size_t sx : nexts)
2290 if (lqn.actposttype[sx] == PrecedenceType::POST_AND) andSuccs.push_back(sx);
2292 for (std::size_t k = 0; k < nexts.size(); ++k) {
2293 const std::size_t nextaidx = nexts[k];
2294 if (next_is_fork) st = preFork;
2295 bool isLoop = lqn.graph.get(aidx, nextaidx) != lqn.dag.get(aidx, nextaidx);
2296 if (lqn.parent[aidx] != lqn.parent[nextaidx]) {
2298 std::size_t cidx = 0;
2299 for (std::size_t c : lqn.callsof[aidx])
2300 if (lqn.callpair_dst[c] == nextaidx) cidx = c;
2301 if (cidx == 0)
continue;
2303 if (lqn.calltype[cidx] != CallType::SYNC)
continue;
2304 const std::size_t gid = group_of_call[cidx];
2308 if (group_members[gid].empty() || group_members[gid][0] != cidx)
continue;
2315 std::vector<std::size_t> tnode2, tstat2, tcall2;
2316 for (std::size_t mc : group_members[gid]) {
2317 const std::size_t tt = lqn.parent[lqn.callpair_dst[mc]];
2318 if (servers_for(tt).empty())
continue;
2319 tnode2.push_back(server_nodes_for(tt)[0]);
2320 tstat2.push_back(servers_for(tt)[0]);
2321 tcall2.push_back(mc);
2323 if (tnode2.size() < 2)
continue;
2324 const std::size_t fromNode2 =
2325 st.jobpos == atClient ? clientNode : st.curnodes[0];
2326 const std::size_t dc = grp_dispatch[gid], gc = grp_class[gid];
2327 m.set_route(st.curclass, dc, fromNode2, grp_router[gid], Tone());
2329 T(Tone() / num_traits<T>::from_int(
int(tnode2.size())));
2330 for (std::size_t d = 0; d < tnode2.size(); ++d) {
2331 m.set_route(dc, dc, grp_router[gid], tnode2[d], share2);
2332 m.set_route(dc, gc, tnode2[d], clientNode, Tone());
2333 m.set_service(tstat2[d], dc, callservtproc[tcall2[d]]);
2334 call_map.push_back({idx, tcall2[d], tstat2[d], dc});
2337 st.jobpos = atClient;
2338 st.curnodes.clear();
2341 if (rp_on && rp_expand[cidx]) {
2344 st = rp_descent(cidx, st);
2347 const std::size_t ctgt = lqn.parent[lqn.callpair_dst[cidx]];
2348 st = route_sync_call(m, idx, cidx, st, server_nodes_for(ctgt),
2349 servers_for(ctgt), callcls, auxcallcls, clientNode,
2350 atClient, atServer, flat);
2354 bool any_entry_succ =
false;
2355 for (std::size_t sx : nexts)
2356 if (sx > lqn.eshift && sx <= lqn.eshift + lqn.nentries) any_entry_succ =
true;
2357 if (!any_entry_succ) {
2358 st.jobpos = jobposkey[aidx];
2359 st.curclass = curclasskey[aidx];
2360 st.curnodes = curnodeskey[aidx];
2362 if (k > 0 && nexts[k - 1] > lqn.eshift && nexts[k - 1] <= lqn.eshift + lqn.nentries)
2363 lastEntryClass = st.curclass;
2364 st.jobpos = atClient;
2365 st.curclass = lastEntryClass;
2366 st.curnodes.clear();
2368 const T w = lqn.graph.get(aidx, nextaidx);
2376 if (iscachelayer && lqn.nitems[aidx] > 0 && cls[nextaidx] != 0) {
2377 const std::size_t readcls = cls[nextaidx];
2378 m.set_route(st.curclass, readcls, clientNode, cacheNode, w);
2379 if (cachepar.pread.size() < m.classes.size())
2380 cachepar.pread.resize(m.classes.size());
2381 cachepar.pread[readcls - 1] = lqn.itemproc[aidx];
2382 const std::vector<std::size_t> hm = lqn.graph.succ(nextaidx);
2384 throw InputError(
"SolverLN: the cache read '" + lqn.names[nextaidx] +
2385 "' needs exactly one hit and one miss successor");
2386 if (cachepar.hitclass.size() < m.classes.size()) {
2387 cachepar.hitclass.resize(m.classes.size(), 0);
2388 cachepar.missclass.resize(m.classes.size(), 0);
2390 cachepar.hitclass[readcls - 1] = cls[hm[0]];
2391 cachepar.missclass[readcls - 1] = cls[hm[1]];
2392 st.jobpos = atCache;
2393 st.curclass = readcls;
2394 st.curnodes.clear();
2395 st = recur(tidx_caller, nextaidx, st);
2399 const bool is_and_join_tail = lqn.actpretype[aidx] == PrecedenceType::PRE_AND;
2401 std::size_t fbranch = 0;
2403 for (std::size_t q = 0; q < andSuccs.size(); ++q)
2404 if (andSuccs[q] == nextaidx) fbranch = q + 1;
2411 const std::size_t hidxOf = lqn.parent[lqn.parent[nextaidx]];
2412 const std::vector<std::size_t>& actStations = servers_for(hidxOf);
2413 const std::vector<std::size_t>& actNodes = server_nodes_for(hidxOf);
2414 const bool actAtServer = !actStations.empty();
2415 const std::size_t from = st.jobpos == atClient ? clientNode
2416 : st.jobpos == atCache ? cacheNode
2419 for (std::size_t r = 0; r < nreplicas; ++r) {
2420 const std::size_t fromNode =
2421 st.jobpos == atClient ? clientNode
2422 : st.jobpos == atCache ? cacheNode
2423 : st.curnodes[std::min(r, st.curnodes.size() - 1)];
2424 const std::size_t toNode = actAtServer ? actNodes[r] : clientNode;
2425 if (next_is_fork && fbranch > 0) {
2426 m.set_route(st.curclass, st.curclass, fromNode, forkNode, Tone());
2427 if (r == 0) forkClassStack.push_back(st.curclass);
2428 m.set_route(st.curclass, st.curclass, forkNode, forkRouter[fbranch - 1],
2430 m.set_route(st.curclass, cls[nextaidx], forkRouter[fbranch - 1], toNode,
2432 }
else if (is_and_join_tail) {
2435 if (forkClassStack.empty())
2436 throw InputError(
"SolverLN: an AND join has no matching fork in '" +
2437 lqn.names[aidx] +
"'");
2438 const std::size_t forkClass = forkClassStack.back();
2439 if (r + 1 == nreplicas) forkClassStack.pop_back();
2440 m.set_route(st.curclass, forkClass, fromNode, joinNode, Tone());
2441 m.set_route(forkClass, cls[nextaidx], joinNode, toNode, Tone());
2443 m.set_route(st.curclass, cls[nextaidx], fromNode, toNode, w);
2446 m.set_service(actStations[r], cls[nextaidx], lqn.hostdem[nextaidx]);
2450 st.jobpos = atServer;
2451 st.curclass = cls[nextaidx];
2452 st.curnodes = actNodes;
2453 servt_map.push_back({idx, nextaidx, actStations[0], cls[nextaidx]});
2455 st.jobpos = atClient;
2456 st.curclass = cls[nextaidx];
2457 st.curnodes.clear();
2458 m.set_service(m.clientIdx, cls[nextaidx], servtproc[nextaidx]);
2459 thinkt_map.push_back({idx, nextaidx, m.clientIdx, cls[nextaidx]});
2461 if (aidx != nextaidx && !isLoop) {
2462 st = recur(tidx_caller, nextaidx, st);
2466 const std::size_t reply =
2467 (rp_on && rp_spliced[tidx_caller]) ? rp_ret_cls[tidx_caller]
2469 if (st.jobpos == atClient) {
2470 m.set_route(st.curclass, reply, clientNode, clientNode, Tone());
2472 for (std::size_t nd : st.curnodes)
2473 m.set_route(st.curclass, reply, nd, clientNode, Tone());
2476 if (!is_aux_class(m.classes[st.curclass - 1].name))
2477 m.classes[st.curclass - 1].completes =
true;
2483 for (std::size_t tidx_caller : callers) {
2484 if (!caller_needs_class(tidx_caller))
continue;
2485 const std::vector<std::size_t>& ents = lqn.entriesof[tidx_caller];
2486 const T share = T(Tone() / num_traits<T>::from_int(
int(ents.size())));
2487 for (std::size_t eidx : ents) {
2488 m.set_route(cls[tidx_caller], cls[eidx], clientNode, clientNode, share);
2489 if (ents.size() > 1)
2490 route_map.push_back({idx, tidx_caller, eidx, m.clientIdx, m.clientIdx,
2491 cls[tidx_caller], cls[eidx]});
2492 Ctx st{cls[eidx], atClient, {}};
2493 recur(tidx_caller, eidx, st);
2501 for (std::size_t g_ = 0; g_ < rp_group_reftask.size(); ++g_) {
2502 if (rp_group_head_is_caller[g_])
continue;
2503 std::vector<std::size_t> roots;
2504 for (std::size_t r_ : rp_group_roots[g_])
2505 if (rp_hop[r_]) roots.push_back(r_);
2506 if (roots.empty())
continue;
2508 const T share = T(Tone() / num_traits<T>::from_int(
int(roots.size())));
2509 for (std::size_t r_ : roots) {
2510 m.set_route(rp_ref_stage[g_], rp_hop_cls[r_], clientNode, clientNode, share);
2511 m.set_route(rp_hop_ret[r_], rp_ref_stage[g_], clientNode, clientNode, Tone());
2515 for (std::size_t ev = 1; ev <= lqn.nidx; ++ev) {
2516 if (!rp_hop[ev])
continue;
2517 Ctx cur{rp_hop_cls[ev], atClient, {}};
2518 for (std::size_t c_ : rp_hop_children[ev]) cur = rp_descent(c_, cur);
2519 m.set_route(cur.curclass, rp_hop_ret[ev], clientNode, clientNode, Tone());
2523 for (std::size_t tv = 1; tv <= lqn.nidx; ++tv) {
2524 if (!rp_spliced[tv])
continue;
2525 for (std::size_t c_ : rp_member_inbound[tv])
2526 if (rp_ret_share[c_] > Tzero())
2527 m.set_route(rp_ret_cls[tv], rp_resume[c_], clientNode, clientNode,
2530 for (std::size_t ev = 1; ev <= lqn.nidx; ++ev) {
2531 if (!rp_hop[ev])
continue;
2532 for (std::size_t c_ : rp_hop_inbound[ev])
2533 if (rp_ret_share[c_] > Tzero())
2534 m.set_route(rp_hop_ret[ev], rp_resume[c_], clientNode, clientNode,
2545 if (sourceStation != 0) {
2546 const T nrep = num_traits<T>::from_int(
int(nreplicas));
2547 for (std::size_t cidx : async_here) {
2548 const std::size_t oc = callcls[cidx];
2549 if (oc == 0)
continue;
2552 const std::vector<std::size_t>& anode =
2553 server_nodes_for(lqn.parent[lqn.callpair_dst[cidx]]);
2554 const T callmean = lqn.callproc_mean[cidx];
2555 if (callmean < Tone()) {
2557 m.set_route(oc, oc, sourceNode, sinkNode, T(Tone() - callmean));
2558 for (std::size_t r = 0; r < anode.size(); ++r) {
2559 m.set_route(oc, oc, sourceNode, anode[r], T(callmean / nrep));
2560 m.set_route(oc, oc, anode[r], sinkNode, Tone());
2564 const T p = T(Tone() / callmean);
2565 for (std::size_t r = 0; r < anode.size(); ++r) {
2566 m.set_route(oc, oc, sourceNode, anode[r], T(Tone() / nrep));
2567 for (std::size_t q = 0; q < anode.size(); ++q)
2568 m.set_route(oc, oc, anode[r], anode[q], T((Tone() - p) / nrep));
2569 m.set_route(oc, oc, anode[r], sinkNode, p);
2573 for (std::size_t eidx : open_entries) {
2575 eo.name = lqn.hashnames[eidx] +
"_Open";
2576 eo.type = JobClassType::OPEN;
2577 eo.population = std::numeric_limits<double>::infinity();
2578 eo.refstat = sourceStation;
2579 eo.completes =
false;
2580 eo.is_ref_class =
false;
2581 eo.attr_kind = int(LqnElement::ENTRY);
2583 const std::size_t ec = m.add_class(eo);
2585 m.set_service(sourceStation, ec, lqn.arrival[eidx]);
2587 std::size_t bound = 0;
2588 for (std::size_t sx : lqn.graph.succ(eidx))
2589 if (bound == 0 && sx > lqn.ashift) bound = sx;
2590 const Distrib<T>& svc = bound != 0 ? servtproc[bound] : servtproc[eidx];
2593 const std::vector<std::size_t>& ostat =
2594 flat ? servers_for(lqn.parent[eidx]) : server;
2595 const std::vector<std::size_t>& onode =
2596 flat ? server_nodes_for(lqn.parent[eidx]) : serverNode;
2597 for (std::size_t r = 0; r < ostat.size(); ++r) {
2598 m.set_service(ostat[r], ec, svc);
2599 m.set_route(ec, ec, sourceNode, onode[r], T(Tone() / nrep));
2600 m.set_route(ec, ec, onode[r], sinkNode, Tone());
2614 std::vector<std::vector<T>> Arows;
2615 std::vector<T> brows;
2616 if (idx < lqn.lincon_A.size() && lqn.lincon_A[idx].rows() > 0) {
2617 const Matrix<T>& Aelem = lqn.lincon_A[idx];
2618 Matrix<T> Alayer(Aelem.rows(), m.classes.size(), Tzero());
2619 const std::vector<std::size_t>& constrained =
2620 ishostlayer ? lqn.tasksof[idx] : lqn.entriesof[idx];
2621 for (std::size_t j = 0; j < constrained.size() && j < Aelem.cols(); ++j) {
2622 std::vector<std::size_t> layerClasses;
2624 for (std::size_t a : lqn.actsof[constrained[j]])
2625 if (cls[a] != 0) layerClasses.push_back(cls[a]);
2627 for (std::size_t c = 1; c <= lqn.ncalls; ++c)
2628 if (lqn.callpair_dst[c] == constrained[j] && callcls[c] != 0)
2629 layerClasses.push_back(callcls[c]);
2631 for (std::size_t k = 0; k < layerClasses.size(); ++k)
2632 for (std::size_t rr = 0; rr < Aelem.rows(); ++rr)
2633 Alayer(rr, layerClasses[k] - 1) =
2634 T(Alayer(rr, layerClasses[k] - 1) + Aelem(rr, j));
2636 for (std::size_t rr = 0; rr < Alayer.rows(); ++rr) {
2637 std::vector<T> row(m.classes.size(), Tzero());
2638 for (std::size_t k = 0; k < m.classes.size(); ++k) row[k] = Alayer(rr, k);
2639 Arows.push_back(row);
2640 brows.push_back(lqn.lincon_b[idx][rr]);
2649 if (rp_on && layer_solver_declares_region()) {
2650 std::vector<std::size_t> srvEntries;
2651 for (std::size_t sidx : idxSet)
2652 for (std::size_t e : lqn.entriesof[sidx]) srvEntries.push_back(e);
2653 for (std::size_t tv = 1; tv <= lqn.nidx; ++tv) {
2654 if (!rp_spliced[tv])
continue;
2655 const double cap = lqn.maxmult[tv];
2656 if (!std::isfinite(cap) || cap <= 0.0)
continue;
2657 std::vector<T> row(m.classes.size(), Tzero());
2661 for (std::size_t a_ : lqn.actsof[tv])
2663 row[cls[a_] - 1] = Tone();
2668 for (std::size_t a_ : lqn.actsof[tv])
2669 for (std::size_t c_ : lqn.callsof[a_])
2670 if (std::find(srvEntries.begin(), srvEntries.end(),
2671 lqn.callpair_dst[c_]) != srvEntries.end() &&
2673 row[callcls[c_] - 1] = Tone();
2678 Arows.push_back(row);
2679 brows.push_back(num_traits<T>::from_double(cap));
2685 for (
const std::vector<T>& row : Arows)
2686 for (
const T& v : row)
2687 if (v != Tzero()) any =
true;
2689 Matrix<T> Alayer(Arows.size(), m.classes.size(), Tzero());
2690 for (std::size_t rr = 0; rr < Arows.size(); ++rr)
2691 for (std::size_t k = 0; k < m.classes.size(); ++k) Alayer(rr, k) = Arows[rr][k];
2695 typename qn::NetworkStruct<T>::Region rg;
2696 const std::size_t M = m.stations.size(), K = m.classes.size();
2697 rg.cap.assign(M, std::vector<double>(K + 1, -1.0));
2698 rg.maxmem.assign(M, -1.0);
2699 rg.members.assign(M,
false);
2701 rg.weight.assign(K, Tone());
2702 rg.size.assign(K, Tone());
2703 for (std::size_t r = 0; r < nreplicas; ++r) rg.members[server[r] - 1] =
true;
2704 rg.lincon_A = Alayer;
2705 rg.lincon_b = brows;
2706 m.regions.push_back(rg);
2713 if (cacheNode != 0) {
2714 cachepar.pread.resize(m.classes.size());
2715 cachepar.hitclass.resize(m.classes.size(), 0);
2716 cachepar.missclass.resize(m.classes.size(), 0);
2717 m.nodeparam[cacheNode] = cachepar;
2725 for (
const std::array<std::size_t, 3>& site : routed_group_sites) {
2726 qn::NodeDef& nd = m.nodes[site[0] - 1];
2727 if (nd.routing.size() < m.classes.size())
2728 nd.routing.resize(m.classes.size(), RoutingStrategy::PROB);
2729 nd.routing[site[1] - 1] = lqn.callgroups[site[2] - 1].strategy;
2741 std::vector<std::pair<std::shared_ptr<std::vector<std::vector<std::size_t>>>,
2742 std::vector<std::vector<std::size_t>>>> deferred_chaincols;
2743 for (std::size_t sidx : idxSet) {
2744 const bool hasld = sidx < lqn.lldscaling.size() && !lqn.lldscaling[sidx].empty();
2745 const bool hascd = sidx < lqn.cdscaling.size() && bool(lqn.cdscaling[sidx]);
2746 const bool hasjd = sidx < lqn.jdscaling.size() && bool(lqn.jdscaling[sidx]);
2747 const bool haspools = sidx < lqn.pools.size() && !lqn.pools[sidx].empty();
2748 if (!(hasld || hascd || hasjd || haspools))
continue;
2749 const bool sishost = sidx <= lqn.nhosts;
2750 const std::vector<std::size_t>& operandIdx =
2751 sishost ? lqn.tasksof[sidx] : lqn.entriesof[sidx];
2752 std::vector<std::vector<std::size_t>> cols(operandIdx.size());
2753 for (std::size_t j = 0; j < operandIdx.size(); ++j) {
2755 for (std::size_t a : lqn.actsof[operandIdx[j]])
2756 if (cls[a] != 0) cols[j].push_back(cls[a] - 1);
2758 for (std::size_t c = 1; c <= lqn.ncalls; ++c)
2759 if (lqn.callpair_dst[c] == operandIdx[j] && callcls[c] != 0)
2760 cols[j].push_back(callcls[c] - 1);
2763 const std::size_t R = m.classes.size();
2764 std::shared_ptr<std::vector<std::vector<std::size_t>>> chaincols =
2765 std::make_shared<std::vector<std::vector<std::size_t>>>();
2766 deferred_chaincols.push_back(std::make_pair(chaincols, cols));
2767 for (std::size_t r = 0; r < nreplicas; ++r) {
2772 qn::Station<T>& stn = m.stations[srv[sidx][r] - 1];
2773 if (hasld) stn.lldscaling = lqn.lldscaling[sidx];
2779 bool one_class_each =
true;
2780 for (std::size_t j = 0; j < cols.size(); ++j)
2781 if (cols[j].size() > 1) one_class_each =
false;
2783 layer_dep_handle(lqn.cdscaling[sidx], cols, chaincols, R);
2784 const std::vector<T> pk = layer_peak(lqn.cdscalingpeak[sidx], cols, R);
2785 if (one_class_each) {
2787 stn.cdscalingpeak = pk;
2790 stn.jdscalingpeak = pk;
2794 stn.jdscaling = layer_dep_handle(lqn.jdscaling[sidx], cols, chaincols, R);
2795 stn.jdscalingpeak = layer_peak(lqn.jdscalingpeak[sidx], cols, R);
2806 const lqn::ServerPools<T>& pl = lqn.pools[sidx];
2807 const Matrix<T> compat = pl.compat;
2808 const std::vector<double> counts = pl.counts;
2809 const std::vector<T> rates = pl.rates;
2811 rates](
const std::vector<T>& nop) {
2812 return std::vector<T>(
2815 stn.jdscaling = layer_dep_handle(etaPool, cols, chaincols, R);
2816 stn.jdscalingpeak = layer_peak(
2817 std::vector<T>(cols.size(), num_traits<T>::from_int(1)), cols, R);
2827 std::vector<std::size_t> nc(m.classes.size(), 0);
2828 for (std::size_t tv = 1; tv <= lqn.nidx; ++tv) {
2829 if (!rp_member[tv] || cls[tv] == 0)
continue;
2830 const std::size_t own = cls[tv];
2831 std::vector<std::size_t> owned{cls[tv]};
2832 for (std::size_t e_ : lqn.entriesof[tv]) owned.push_back(cls[e_]);
2833 for (std::size_t a_ : lqn.actsof[tv]) {
2834 owned.push_back(cls[a_]);
2835 for (std::size_t c_ : lqn.callsof[a_]) {
2836 owned.push_back(callcls[c_]);
2837 owned.push_back(auxcallcls[c_]);
2841 for (std::size_t q_ : owned)
2842 if (q_ != 0 && m.classes[q_ - 1].type == JobClassType::CLOSED) nc[q_ - 1] = own;
2844 layer_normclass[idx] = nc;
2845 ref_path_audit(m, idx, callers, caller_needs_class, rp_spliced, rp_member, rp_hop,
2846 rp_expand, rp_group_head_is_caller, rp_ref_stage, rp_hop_cls,
2847 rp_hop_ret, rp_ret_cls, rp_gate, rp_aux, rp_resume, cls);
2851 apply_host_task_priorities(m, idxSet);
2856 for (std::size_t d = 0; d < deferred_chaincols.size(); ++d) {
2857 const std::vector<std::vector<std::size_t>>& cols = deferred_chaincols[d].second;
2858 std::vector<std::vector<std::size_t>>& out = *deferred_chaincols[d].first;
2859 out.assign(cols.size(), std::vector<std::size_t>());
2860 for (std::size_t j = 0; j < cols.size(); ++j) {
2861 for (std::size_t k = 0; k < cols[j].size(); ++k)
2862 for (std::size_t cc = 0; cc < m.chains.size(); ++cc)
2863 if (cols[j][k] < m.chains[cc].size() && m.chains[cc][cols[j][k]]) {
2865 for (std::size_t q = 0; q < out[j].size(); ++q)
2866 if (out[j][q] == cc) seen =
true;
2867 if (!seen) out[j].push_back(cc);
2888 bool layer_solver_declares_region()
const {
2889 return opt.layer_solver ==
"nc" || opt.layer_solver ==
"fluid" || opt.layer_solver ==
"ssa";
2898 template <
class NeedsClass>
2899 void ref_path_audit(
const qn::Layer<T>& m, std::size_t idx,
2900 const std::vector<std::size_t>& callers,
const NeedsClass& needs_class,
2901 const std::vector<bool>& spliced,
const std::vector<bool>& member,
2902 const std::vector<bool>& hop,
const std::vector<bool>& expand,
2903 const std::vector<bool>& head_is_caller,
2904 const std::vector<std::size_t>& ref_stage,
2905 const std::vector<std::size_t>& hop_cls,
2906 const std::vector<std::size_t>& hop_ret,
2907 const std::vector<std::size_t>& ret_cls,
2908 const std::vector<std::size_t>& gate,
const std::vector<std::size_t>& aux,
2909 const std::vector<std::size_t>& resume,
2910 const std::vector<std::size_t>& cls)
const {
2911 std::set<std::size_t> watch;
2912 for (std::size_t k : ref_stage)
2913 if (k) watch.insert(k);
2914 for (std::size_t e = 1; e < hop.size(); ++e)
2916 if (hop_cls[e]) watch.insert(hop_cls[e]);
2917 if (hop_ret[e]) watch.insert(hop_ret[e]);
2919 for (std::size_t t = 1; t < member.size(); ++t)
2921 if (ret_cls[t]) watch.insert(ret_cls[t]);
2922 if (cls[t]) watch.insert(cls[t]);
2924 for (std::size_t c = 1; c < expand.size(); ++c)
2926 if (gate[c]) watch.insert(gate[c]);
2927 if (aux[c]) watch.insert(aux[c]);
2928 if (resume[c]) watch.insert(resume[c]);
2930 const std::size_t I = m.nodes.size();
2931 for (std::size_t r : watch) {
2932 for (std::size_t n = 1; n <= I; ++n) {
2935 for (
const auto& kv : m.P) {
2936 if (kv.first.first != r || n > kv.second.rows())
continue;
2937 for (std::size_t j = 0; j < kv.second.cols(); ++j) {
2938 const double v = dbl(kv.second(n - 1, j));
2944 std::cerr <<
"[LINE] Warning: refpath layer '" << lqn.hashnames[idx]
2945 <<
"': class '" << m.classes[r - 1].name <<
"' routes out of node "
2946 << n <<
" with a negative probability.\n";
2948 std::cerr <<
"[LINE] Warning: refpath layer '" << lqn.hashnames[idx]
2949 <<
"': class '" << m.classes[r - 1].name <<
"' leaves node " << n
2950 <<
" with total probability " << tot <<
".\n";
2953 std::size_t nref = 0;
2954 for (
const qn::JobClass& jc : m.classes)
2955 if (jc.is_ref_class) ++nref;
2956 std::size_t expect = 0;
2957 for (
bool h : head_is_caller)
2959 for (std::size_t t : callers)
2960 if (!spliced[t] && needs_class(t)) ++expect;
2962 std::cerr <<
"[LINE] Warning: refpath layer '" << lqn.hashnames[idx] <<
"' declares "
2963 << nref <<
" reference classes where " << expect
2964 <<
" were built; a chain holding two makes refresh_chains throw.\n";
2996 std::size_t curclass;
3001 static bool is_aux_class(
const std::string& name) {
3002 return name.size() >= 4 && name.compare(name.size() - 4, 4,
".Aux") == 0;
3004 template <
class Ctx>
3005 Ctx route_sync_call(qn::Layer<T>& m, std::size_t idx, std::size_t cidx, Ctx st,
3006 const std::vector<std::size_t>& tnode,
3007 const std::vector<std::size_t>& tstat,
3008 const std::vector<std::size_t>& callcls,
3009 const std::vector<std::size_t>& auxcallcls, std::size_t clientNode,
3010 int atClient,
int atServer,
bool flat) {
3011 const T one = Tone();
3015 const std::size_t ntgt = tnode.size();
3016 const bool to_this_server = ntgt > 0;
3017 const T nrep = num_traits<T>::from_int(
int(ntgt > 0 ? ntgt : 1));
3018 const T share = T(one / nrep);
3019 const T callmean = lqn.callproc_mean[cidx];
3020 const std::size_t cc = callcls[cidx];
3021 const std::size_t ax = auxcallcls[cidx];
3022 const bool below = callmean < one;
3023 const bool above = callmean > one;
3024 if (st.jobpos == atClient) {
3025 if (to_this_server) {
3027 m.set_route(st.curclass, ax, clientNode, clientNode, T(one - callmean));
3028 for (std::size_t r = 0; r < ntgt; ++r) {
3029 m.set_route(st.curclass, cc, clientNode, tnode[r], T(callmean / nrep));
3030 m.set_route(cc, cc, tnode[r], clientNode, one);
3033 m.set_route(ax, cc, clientNode, clientNode, one);
3035 for (std::size_t r = 0; r < ntgt; ++r) {
3036 m.set_route(st.curclass, cc, clientNode, tnode[r], share);
3037 m.set_route(cc, ax, tnode[r], clientNode, one);
3038 m.set_route(ax, cc, clientNode, tnode[r], T((one - one / callmean) / nrep));
3040 m.set_route(ax, cc, clientNode, clientNode, T(one / callmean));
3042 for (std::size_t r = 0; r < ntgt; ++r) {
3043 m.set_route(st.curclass, cc, clientNode, tnode[r], share);
3044 m.set_route(cc, cc, tnode[r], clientNode, one);
3047 for (std::size_t r = 0; r < ntgt; ++r) {
3048 m.set_service(tstat[r], cc, callservtproc[cidx]);
3049 call_map.push_back({idx, cidx, tstat[r], cc});
3052 st.jobpos = atClient;
3053 st.curnodes.clear();
3056 m.set_route(st.curclass, cc, clientNode, clientNode, one);
3057 if (below || above) {
3058 m.set_route(cc, ax, clientNode, clientNode, one);
3063 st.jobpos = atClient;
3064 st.curnodes.clear();
3065 m.set_service(m.clientIdx, cc, callservtproc[cidx]);
3066 call_map.push_back({idx, cidx, m.clientIdx, cc});
3070 auto fromNode = [&](std::size_t r) {
3071 return st.curnodes[std::min(r, st.curnodes.size() - 1)];
3073 if (to_this_server) {
3082 for (std::size_t r = 0; r < ntgt; ++r) {
3083 m.set_route(st.curclass, ax, fromNode(r), clientNode, T(one - callmean));
3084 m.set_route(st.curclass, cc, fromNode(r), tnode[r], T(callmean / nrep));
3085 m.set_route(cc, cc, tnode[r], clientNode, one);
3087 m.set_route(ax, cc, clientNode, clientNode, one);
3089 st.jobpos = atClient;
3090 st.curnodes.clear();
3096 for (std::size_t r = 0; r < ntgt; ++r) {
3097 m.set_route(st.curclass, cc, fromNode(r), tnode[r], one);
3098 m.set_route(cc, ax, tnode[r], clientNode, one);
3099 m.set_route(ax, cc, clientNode, tnode[r],
3100 T((one - one / callmean) / nrep));
3102 m.set_route(ax, cc, clientNode, clientNode, T(one / callmean));
3106 for (std::size_t r = 0; r < ntgt; ++r) {
3107 m.set_route(st.curclass, cc, fromNode(r), tnode[r], one);
3108 m.set_route(cc, cc, tnode[r], tnode[r], T(one - one / callmean));
3109 m.set_route(cc, ax, tnode[r], clientNode, T(one / callmean));
3113 st.jobpos = atClient;
3114 st.curnodes.clear();
3116 for (std::size_t r = 0; r < ntgt; ++r)
3117 m.set_route(st.curclass, cc, fromNode(r), tnode[r], one);
3121 for (std::size_t r = 0; r < ntgt; ++r)
3122 m.set_route(cc, cc, tnode[r], clientNode, one);
3124 st.jobpos = atClient;
3125 st.curnodes.clear();
3127 st.jobpos = atServer;
3128 st.curnodes = tnode;
3132 for (std::size_t r = 0; r < ntgt; ++r) {
3133 m.set_service(tstat[r], cc, callservtproc[cidx]);
3134 call_map.push_back({idx, cidx, tstat[r], cc});
3137 for (std::size_t nd : st.curnodes)
3138 m.set_route(st.curclass, cc, nd, clientNode, one);
3139 if (below || above) {
3140 m.set_route(cc, ax, clientNode, clientNode, one);
3145 st.jobpos = atClient;
3146 st.curnodes.clear();
3147 m.set_service(m.clientIdx, cc, callservtproc[cidx]);
3148 call_map.push_back({idx, cidx, m.clientIdx, cc});
3163 void build_entry_service_matrix() {
3164 const std::size_t dim = lqn.nidx + lqn.ncalls;
3165 servtmatrix =
Matrix<T>(dim + 1, dim + 1, Tzero());
3166 std::function<void(std::size_t, std::size_t)> rec = [&](std::size_t aidx, std::size_t eidx) {
3167 for (std::size_t nextaidx : lqn.graph.succ(aidx)) {
3168 const bool isLoop = lqn.graph.get(aidx, nextaidx) != lqn.dag.get(aidx, nextaidx);
3169 if (lqn.parent[aidx] != lqn.parent[nextaidx]) {
3170 for (std::size_t cidx : lqn.callsof[aidx])
3171 if (lqn.calltype[cidx] == CallType::SYNC)
3172 servtmatrix(eidx, lqn.nidx + cidx) = Tone();
3173 }
else if (nextaidx != aidx && !isLoop) {
3174 servtmatrix(eidx, nextaidx) = Tone();
3175 rec(nextaidx, eidx);
3179 for (std::size_t e = 1; e <= lqn.nentries; ++e) {
3180 const std::size_t eidx = lqn.eshift + e;
3190 void init_interlock() {
3191 const std::size_t NE = lqn.nentries;
3192 il_all =
Matrix<T>(NE + 1, NE + 1, Tzero());
3193 il_ph1 =
Matrix<T>(NE + 1, NE + 1, Tzero());
3195 std::function<void(std::size_t, std::size_t, T, T, std::vector<bool>&,
int)> trace =
3196 [&](std::size_t eidx, std::size_t root_e, T pall, T pph1, std::vector<bool>& visited,
3198 if (eidx <= lqn.eshift || eidx > lqn.eshift + NE)
return;
3199 const std::size_t e = eidx - lqn.eshift;
3200 if (visited[e])
return;
3202 il_all(root_e, e) = T(il_all(root_e, e) + pall);
3203 il_ph1(root_e, e) = T(il_ph1(root_e, e) + pph1);
3204 for (std::size_t aidx : lqn.actsof[eidx]) {
3205 if (aidx <= lqn.ashift || aidx > lqn.ashift + lqn.nacts)
continue;
3206 const std::size_t a = aidx - lqn.ashift;
3207 if (depth > 0 && lqn.actphase[a] > 1)
continue;
3208 const bool is_ph1 = lqn.actphase[a] <= 1;
3209 for (std::size_t cidx : lqn.callsof[aidx]) {
3210 if (lqn.calltype[cidx] != CallType::SYNC)
continue;
3211 if (!(lqn.callproc_mean[cidx] > Tzero()))
continue;
3212 const std::size_t dst = lqn.callpair_dst[cidx];
3213 if (dst <= lqn.eshift || dst > lqn.eshift + NE)
continue;
3214 trace(dst, root_e, T(pall * lqn.callproc_mean[cidx]),
3215 is_ph1 ? T(pph1 * lqn.callproc_mean[cidx]) : Tzero(), visited,
3221 for (std::size_t e = 1; e <= NE; ++e) {
3222 std::vector<bool> visited(NE + 1,
false);
3223 trace(lqn.eshift + e, e, Tone(), Tone(), visited, 0);
3226 il_common_entries.assign(NT() + 1, {});
3227 il_src_all.assign(NT() + 1, {});
3228 il_src_ph2.assign(NT() + 1, {});
3229 il_num_sources.assign(NT() + 1, 0.0);
3231 for (std::size_t t = 1; t <= lqn.ntasks; ++t) {
3232 const std::size_t tidx = lqn.tshift + t;
3233 if (lqn.isref[tidx] || lqn.sched[tidx] == SchedStrategy::INF)
continue;
3234 interlock_for_server(tidx);
3236 for (std::size_t h = 1; h <= lqn.nhosts; ++h) {
3237 if (lqn.sched[h] == SchedStrategy::INF)
continue;
3238 interlock_for_server(h);
3242 std::vector<std::size_t> server_entry_nums(std::size_t serverIdx)
const {
3243 std::vector<std::size_t> out;
3244 if (serverIdx <= lqn.nhosts) {
3245 for (std::size_t tidx : lqn.tasksof[serverIdx])
3246 for (std::size_t se : lqn.entriesof[tidx]) out.push_back(se - lqn.eshift);
3248 for (std::size_t se : lqn.entriesof[serverIdx]) out.push_back(se - lqn.eshift);
3253 std::vector<std::size_t> client_tasks(std::size_t serverIdx)
const {
3254 std::vector<std::size_t> out;
3255 if (serverIdx <= lqn.nhosts)
return lqn.tasksof[serverIdx];
3256 for (std::size_t se : lqn.entriesof[serverIdx])
3257 for (std::size_t ci : lqn.iscaller.col(se))
3258 if (ci > lqn.tshift && ci <= NT()) out.push_back(ci);
3259 std::sort(out.begin(), out.end());
3260 out.erase(std::unique(out.begin(), out.end()), out.end());
3264 std::vector<std::size_t> call_dst_tasks(std::size_t src_eidx, std::size_t target_e)
const {
3265 std::vector<std::size_t> out;
3266 for (std::size_t aidx : lqn.actsof[src_eidx]) {
3267 if (aidx <= lqn.ashift || aidx > lqn.ashift + lqn.nacts)
continue;
3268 for (std::size_t cidx : lqn.callsof[aidx]) {
3269 if (lqn.calltype[cidx] != CallType::SYNC)
continue;
3270 const std::size_t dst = lqn.callpair_dst[cidx];
3271 const std::size_t de = dst - lqn.eshift;
3272 if (de >= 1 && de <= lqn.nentries && il_all(de, target_e) > Tzero())
3273 out.push_back(lqn.parent[dst]);
3276 std::sort(out.begin(), out.end());
3277 out.erase(std::unique(out.begin(), out.end()), out.end());
3281 bool is_branch_point(std::size_t srcX, std::size_t entryA, std::size_t srcY,
3282 std::size_t entryB)
const {
3283 const std::size_t taskA = lqn.parent[entryA], taskB = lqn.parent[entryB];
3284 const std::size_t taskX = lqn.parent[srcX];
3285 if (taskX == taskA && taskX == taskB)
return false;
3286 if (srcX == entryA || srcY == entryB)
return true;
3287 const std::vector<std::size_t> dx = call_dst_tasks(srcX, entryA - lqn.eshift);
3288 const std::vector<std::size_t> dy = call_dst_tasks(srcY, entryB - lqn.eshift);
3289 for (std::size_t a : dx)
3290 for (std::size_t b : dy)
3291 if (a != b)
return true;
3295 void trace_to_server(std::size_t eidx, std::size_t serverIdx, std::vector<bool>& visited,
3296 std::vector<std::size_t>& itasks,
bool isHead)
const {
3297 if (eidx <= lqn.eshift || eidx > lqn.eshift + lqn.nentries)
return;
3298 const std::size_t e = eidx - lqn.eshift;
3299 if (visited[e])
return;
3300 const std::size_t owner = lqn.parent[eidx];
3301 if (owner == serverIdx)
return;
3302 if (serverIdx <= lqn.nhosts && lqn.parent[owner] == serverIdx)
return;
3305 const std::vector<std::size_t> sen = server_entry_nums(serverIdx);
3306 for (std::size_t aidx : lqn.actsof[eidx]) {
3307 if (aidx <= lqn.ashift || aidx > lqn.ashift + lqn.nacts)
continue;
3308 for (std::size_t cidx : lqn.callsof[aidx]) {
3309 if (lqn.calltype[cidx] != CallType::SYNC)
continue;
3310 const std::size_t dst = lqn.callpair_dst[cidx];
3311 const std::size_t dtask = lqn.parent[dst];
3312 bool reaches = dtask == serverIdx ||
3313 (serverIdx <= lqn.nhosts && lqn.parent[dtask] == serverIdx);
3315 const std::size_t de = dst - lqn.eshift;
3316 for (std::size_t sn : sen)
3317 if (il_all(de, sn) > Tzero()) reaches =
true;
3320 trace_to_server(dst, serverIdx, visited, itasks,
false);
3325 if (found && !isHead) {
3326 itasks.push_back(owner);
3327 std::sort(itasks.begin(), itasks.end());
3328 itasks.erase(std::unique(itasks.begin(), itasks.end()), itasks.end());
3333 void interlock_for_server(std::size_t serverIdx) {
3334 const std::vector<std::size_t> sen = server_entry_nums(serverIdx);
3335 if (sen.empty())
return;
3336 const std::vector<std::size_t> cts = client_tasks(serverIdx);
3337 if (cts.empty())
return;
3339 std::vector<std::pair<std::size_t, std::size_t>> pairs;
3340 for (std::size_t ct : cts)
3341 for (std::size_t ce : lqn.entriesof[ct]) {
3342 const std::size_t cen = ce - lqn.eshift;
3343 if (cen < 1 || cen > lqn.nentries)
continue;
3344 for (std::size_t se : sen)
3345 if (il_all(cen, se) > Tzero()) {
3346 pairs.emplace_back(ct, cen);
3350 if (pairs.size() < 2)
return;
3352 std::vector<std::size_t> common;
3353 for (std::size_t i = 0; i < pairs.size(); ++i)
3354 for (std::size_t j = i + 1; j < pairs.size(); ++j) {
3355 if (pairs[i].first == pairs[j].first)
continue;
3356 const std::size_t eA = pairs[i].second, eC = pairs[j].second;
3357 for (std::size_t t = 1; t <= lqn.ntasks; ++t) {
3358 const std::size_t tidx = lqn.tshift + t;
3359 for (std::size_t ex : lqn.entriesof[tidx])
3360 for (std::size_t ey : lqn.entriesof[tidx]) {
3361 const std::size_t xn = ex - lqn.eshift, yn = ey - lqn.eshift;
3362 if (xn < 1 || yn < 1 || xn > lqn.nentries || yn > lqn.nentries)
continue;
3363 if (il_all(xn, eA) > Tzero() && il_all(yn, eC) > Tzero() &&
3364 is_branch_point(ex, eA + lqn.eshift, ey, eC + lqn.eshift))
3365 common.push_back(ex);
3369 std::sort(common.begin(), common.end());
3370 common.erase(std::unique(common.begin(), common.end()), common.end());
3371 if (common.empty())
return;
3373 std::vector<std::size_t> interlocked;
3374 for (std::size_t ce : common) {
3375 std::vector<bool> visited(lqn.nentries + 1,
false);
3376 std::vector<std::size_t> it;
3377 trace_to_server(ce, serverIdx, visited, it,
true);
3378 for (std::size_t x : it) interlocked.push_back(x);
3380 std::sort(interlocked.begin(), interlocked.end());
3381 interlocked.erase(std::unique(interlocked.begin(), interlocked.end()), interlocked.end());
3383 std::vector<std::size_t> src_all;
3384 for (std::size_t ce : common) src_all.push_back(lqn.parent[ce]);
3385 std::sort(src_all.begin(), src_all.end());
3386 src_all.erase(std::unique(src_all.begin(), src_all.end()), src_all.end());
3388 std::vector<std::size_t> diff;
3389 for (std::size_t x : src_all)
3390 if (!std::binary_search(interlocked.begin(), interlocked.end(), x)) diff.push_back(x);
3398 std::vector<std::size_t> src_ph2;
3399 for (std::size_t it : interlocked)
3400 for (std::size_t ie : lqn.entriesof[it])
3401 for (std::size_t ci : lqn.iscaller.col(ie))
3402 if (ci > lqn.tshift && ci <= NT() &&
3403 !std::binary_search(interlocked.begin(), interlocked.end(), ci))
3404 src_all.push_back(ci);
3405 std::sort(src_all.begin(), src_all.end());
3406 src_all.erase(std::unique(src_all.begin(), src_all.end()), src_all.end());
3409 for (std::size_t st : src_all) nsrc += lqn.mult[st];
3411 il_common_entries[serverIdx] = common;
3412 il_src_all[serverIdx] = src_all;
3413 il_src_ph2[serverIdx] = src_ph2;
3414 il_num_sources[serverIdx] = nsrc;
3421 const std::size_t N = lqn.nidx;
3422 tput.assign(N + 1, Tzero());
3423 util.assign(N + 1, Tzero());
3424 servt.assign(N + 1, Tzero());
3425 residt.assign(N + 1, Tzero());
3426 thinkt.assign(N + 1, Tzero());
3427 callservt.assign(lqn.ncalls + 1, Tzero());
3428 callresidt.assign(lqn.ncalls + 1, Tzero());
3430 servt_ph1.assign(N + 1, Tzero());
3431 servt_ph2.assign(N + 1, Tzero());
3432 prOvertake.assign(lqn.nentries + 1, Tzero());
3433 build_entry_service_matrix();
3435 relax_omega = (opt.relax ==
"fixed" || opt.relax ==
"adaptive") ? opt.relax_factor : 1.0;
3437 servt_prev.assign(N + 1, std::numeric_limits<double>::quiet_NaN());
3438 residt_prev.assign(N + 1, std::numeric_limits<double>::quiet_NaN());
3439 tput_prev.assign(N + 1, std::numeric_limits<double>::quiet_NaN());
3440 thinkt_prev.assign(N + 1, std::numeric_limits<double>::quiet_NaN());
3441 callservt_prev.assign(lqn.ncalls + 1, std::numeric_limits<double>::quiet_NaN());
3442 callresidt_prev.assign(lqn.ncalls + 1, std::numeric_limits<double>::quiet_NaN());
3443 servt_prev_v.assign(N + 1, Tzero());
3444 residt_prev_v.assign(N + 1, Tzero());
3445 tput_prev_v.assign(N + 1, Tzero());
3446 thinkt_prev_v.assign(N + 1, Tzero());
3447 callservt_prev_v.assign(lqn.ncalls + 1, Tzero());
3449 unique_route_idx.clear();
3450 for (
const RouteRow& r : route_map) unique_route_idx.push_back(r.idx);
3451 std::sort(unique_route_idx.begin(), unique_route_idx.end());
3452 unique_route_idx.erase(std::unique(unique_route_idx.begin(), unique_route_idx.end()),
3453 unique_route_idx.end());
3456 if (interlockMethod ==
"ilrate") init_interlock();
3458 maxitererr.assign(opt.iter_max + 2, 0.0);
3459 results_raw_last.clear();
3460 results_raw_valid =
false;
3461 averagingstart = -1;
3462 hasconverged =
false;
3463 moment_pass_done =
false;
3466 servtcdf.assign(N + 1, fluid::FluidPassage());
3467 callservtcdf.assign(lqn.ncalls + 1, fluid::FluidPassage());
3468 entrycdfrespt.assign(lqn.nentries + 1, LnCdf());
3469 entryproc.assign(lqn.nentries + 1, mam::AphPair<T>());
3470 cdf_repo.assign(ensemble.size(), std::vector<std::vector<fluid::FluidPassage>>());
3484 const bool stoch = opt.layer_solver ==
"ssa";
3487 cfg.iter_tol = opt.iter_tol;
3488 cfg.relax_burnin = relax_omega;
3489 stoch_ctl = std::make_shared<LnStochController<T>>(cfg);
3494 bool done = !stoch && converged(it);
3495 while (!done && it < opt.iter_max) {
3496 if (results_raw_valid) {
3497 results.back() = results_raw_last;
3498 results_raw_valid =
false;
3501 results.emplace_back(ensemble.size());
3502 for (std::size_t e = 0; e < ensemble.size(); ++e) analyze(it, e);
3505 std::vector<double> jobs(ensemble.size(), 0.0);
3506 for (std::size_t e = 0; e < ensemble.size(); ++e)
3507 jobs[e] = ensemble[e].total_jobs();
3508 const bool stop = stoch_ctl->update(it, results.back(), jobs, servt, residt);
3511 relax_omega = stoch_ctl->relax_omega();
3513 did_converge =
true;
3514 hasconverged =
true;
3518 done = converged(it);
3521 iterations_done = it;
3524 if (stoch && stoch_ctl && stoch_ctl->averaging_count() > 0) {
3525 const std::vector<LayerResult<T>>& avg = stoch_ctl->averaged_results();
3526 for (std::size_t e = 0; e < results.back().size() && e < avg.size(); ++e)
3527 results.back()[e] = avg[e];
3528 servt = stoch_ctl->averaged_servt();
3529 residt = stoch_ctl->averaged_residt();
3553 Matrix<T> filter_metric(
const qn::Layer<T>& L,
const Matrix<T>& metric,
3554 const std::vector<std::vector<bool>>* zero_mask)
const {
3555 const std::size_t M = L.nstations, K = L.nclasses;
3556 Matrix<T> out(M, K, Tzero());
3557 for (std::size_t i = 0; i < M; ++i)
3558 for (std::size_t k = 0; k < K; ++k)
3559 if (!L.disabled[i][k]) out(i, k) = metric(i, k);
3561 for (std::size_t i = 0; i < M; ++i)
3562 for (std::size_t k = 0; k < K; ++k)
3563 if ((*zero_mask)[i][k]) out(i, k) = Tzero();
3564 for (std::size_t i = 0; i < M; ++i)
3565 for (std::size_t k = 0; k < K; ++k)
3567 for (std::size_t k = 0; k < K; ++k) {
3568 std::size_t c = L.nchains;
3569 for (std::size_t cc = 0; cc < L.nchains; ++cc)
3570 if (L.chains[cc][k]) c = cc;
3571 if (c == L.nchains)
continue;
3572 for (std::size_t i = 0; i < M; ++i)
3573 if (L.visits[c](L.stateful_of_station(i + 1) - 1, k) == Tzero())
3574 out(i, k) = Tzero();
3621 mva::MvaSolution<T> solve_layer_ctmc(std::size_t e) {
3622 qn::Layer<T>& L = ensemble[e];
3627 ctmc::CtmcOptions co;
3632 mva::MvaSolution<T> out;
3633 out.method =
"ctmc";
3652 mva::MvaSolution<T> solve_layer_mam(std::size_t e) {
3653 qn::Layer<T>& L = ensemble[e];
3655 L.refresh_capacity();
3657 mo.method =
"dec.poisson";
3673 mva::MvaSolution<T> solve_layer_nc(std::size_t e) {
3674 if (fj_tr[e].active())
3676 "SolverLN: layer '" + ensemble[e].name +
3677 "' carries a fork, whose fixed point is driven by MVA in this port; solve this "
3678 "model with layer_solver 'mva'");
3679 qn::Layer<T>& L = ensemble[e];
3709 mva::MvaSolution<T> solve_layer_ssa(std::size_t e) {
3710 if (fj_tr[e].active())
3712 "SolverLN: layer '" + ensemble[e].name +
3713 "' carries a fork, whose fixed point is driven by MVA in this port; solve this "
3714 "model with layer_solver 'mva'");
3715 qn::Layer<T>& L = ensemble[e];
3718 mva::MvaSolution<T> out;
3721 out.Q =
Matrix<T>(L.nstations, L.nclasses, Tzero());
3722 out.U =
Matrix<T>(L.nstations, L.nclasses, Tzero());
3723 out.R =
Matrix<T>(L.nstations, L.nclasses, Tzero());
3724 out.Tp =
Matrix<T>(L.nstations, L.nclasses, Tzero());
3725 for (std::size_t i = 0; i < L.nstations; ++i)
3726 for (std::size_t r = 0; r < L.nclasses; ++r) {
3727 out.Q(i, r) = num_traits<T>::from_double(s.QN(i, r));
3728 out.U(i, r) = num_traits<T>::from_double(s.UN(i, r));
3729 out.R(i, r) = num_traits<T>::from_double(s.RN(i, r));
3730 out.Tp(i, r) = num_traits<T>::from_double(s.TN(i, r));
3732 out.C.assign(L.nclasses, Tzero());
3733 out.X.assign(L.nclasses, Tzero());
3734 for (std::size_t r = 0; r < L.nclasses && r < s.XN.size(); ++r) {
3735 out.X[r] = num_traits<T>::from_double(s.XN[r]);
3736 out.C[r] = num_traits<T>::from_double(s.CN[r]);
3741 mva::MvaSolution<T> solve_layer_fluid(std::size_t e) {
3742 if (fj_tr[e].active())
3744 "SolverLN: layer '" + std::to_string(e) +
3745 "' carries a fork, whose fixed point is driven by MVA in this port; solve this "
3746 "model with layer_solver 'mva'");
3747 qn::Layer<T>& L = ensemble[e];
3757 mva::MvaSolution<T> out;
3758 out.method =
"fluid";
3760 out.Q =
Matrix<T>(L.nstations, L.nclasses, Tzero());
3761 out.U =
Matrix<T>(L.nstations, L.nclasses, Tzero());
3762 out.R =
Matrix<T>(L.nstations, L.nclasses, Tzero());
3763 out.Tp =
Matrix<T>(L.nstations, L.nclasses, Tzero());
3764 out.C.assign(L.nclasses, Tzero());
3765 out.X.assign(L.nclasses, Tzero());
3766 detail::ln_fluid_solve(L, opt.layer_fluid, out);
3770 mva::MvaSolution<T> solve_layer(std::size_t e) {
3780 ensemble[e].refresh_capacity();
3788 if (!ensemble[e].regions.empty())
return solve_layer_ctmc(e);
3801 if (has_cache_node(e)) ensemble[e].refresh_rt();
3802 if (opt.layer_solver ==
"fluid")
return solve_layer_fluid(e);
3803 if (opt.layer_solver ==
"nc")
return solve_layer_nc(e);
3804 if (opt.layer_solver ==
"ssa")
return solve_layer_ssa(e);
3805 if (opt.layer_solver !=
"mva")
3807 "' is not available; use 'mva', 'nc', 'fluid' or 'ssa'");
3814 mva::MvaOptions eopt = opt.layer;
3815 if (e < layer_interlock.size()) eopt.interlock = layer_interlock[e];
3816 if (!fj_tr[e].active())
3818 mva::MvaOptions lopt = eopt;
3819 if (lopt.method ==
"default") lopt.method =
"amva";
3821 ensemble[e], fj_tr[e], fj_lambda[e], eopt,
3822 [&lopt](qn::NetworkStruct<T>& V) {
3827 void analyze(
int it, std::size_t e) {
3828 qn::Layer<T>& L = ensemble[e];
3829 const mva::MvaSolution<T> s = solve_layer(e);
3830 LayerResult<T>& r = results[it - 1][e];
3831 r.RN = filter_metric(L, s.R,
nullptr);
3832 std::vector<std::vector<bool>> zmask(L.nstations, std::vector<bool>(L.nclasses,
false));
3833 for (std::size_t i = 0; i < L.nstations; ++i)
3834 for (std::size_t k = 0; k < L.nclasses; ++k)
3836 r.QN = filter_metric(L, s.Q, &zmask);
3837 r.UN = filter_metric(L, s.U, &zmask);
3838 r.TN = filter_metric(L, s.Tp,
nullptr);
3839 r.WN = residence_from_response(L, r.RN);
3849 if (!fj_tr[e].active() && opt.layer_solver ==
"mva") {
3850 Matrix<T> Qch(L.nstations, L.nchains, Tzero());
3851 for (std::size_t c = 0; c < L.nchains; ++c)
3852 for (std::size_t i = 0; i < L.nstations; ++i) {
3854 for (std::size_t k : L.inchain[c]) s2 += s.Q(i, k - 1);
3857 layer_init_sol[e] = Qch;
3862 Matrix<T> residence_from_response(
const qn::Layer<T>& L,
const Matrix<T>& RN)
const {
3863 Matrix<T> V(L.nstations, L.nclasses, Tzero());
3864 for (std::size_t c = 0; c < L.nchains; ++c)
3865 for (std::size_t i = 0; i < L.nstations; ++i) {
3866 const std::size_t sf = L.stateful_of_station(i + 1) - 1;
3867 for (std::size_t k = 0; k < L.nclasses; ++k)
3868 V(i, k) = T(V(i, k) + L.visits[c](sf, k));
3870 Matrix<T> WN(L.nstations, L.nclasses, Tzero());
3871 for (std::size_t i = 0; i < L.nstations; ++i)
3872 for (std::size_t k = 0; k < L.nclasses; ++k) {
3873 if (L.disabled[i][k])
continue;
3874 if (!(RN(i, k) > Tzero()))
continue;
3876 WN(i, k) = RN(i, k);
3880 for (std::size_t cc = 0; cc < L.nchains; ++cc)
3881 if (L.chains[cc][k]) c = cc;
3882 const std::size_t rstat = L.classes[k].refstat;
3884 if (L.refclass[c] > 0) {
3885 den = V(rstat - 1, L.refclass[c] - 1);
3887 for (std::size_t kk : L.inchain[c]) den += V(rstat - 1, kk - 1);
3889 if (den == Tzero())
continue;
3890 WN(i, k) = T(RN(i, k) * V(i, k) / den);
3892 for (std::size_t i = 0; i < L.nstations; ++i)
3893 for (std::size_t k = 0; k < L.nclasses; ++k)
3902 if (interlockMethod ==
"ilrate") update_populations(it);
3903 update_think_times(it);
3905 update_routing_probabilities(it);
3906 for (std::size_t e : route_reset) {
3907 ensemble[e].refresh_chains();
3913 for (std::size_t e : svc_reset) ensemble[e].refresh_rates();
3917 bool converged(
int it) {
3918 const std::size_t E = ensemble.size();
3919 const int iter_min = std::max<int>(2 *
int(E),
int(std::ceil(opt.iter_max / 4.0)));
3920 const int wnd_size = std::max(5,
int(std::ceil(iter_min / 5.0)));
3922 if (it >= iter_min &&
int(results.size()) >= wnd_size) {
3923 results_raw_last = results.back();
3924 results_raw_valid =
true;
3925 const T w = T(Tone() / num_traits<T>::from_int(wnd_size));
3926 for (std::size_t e = 0; e < E; ++e) {
3927 LayerResult<T>& cur = results[results.size() - 1][e];
3928 auto scale = [&](Matrix<T>& A) {
3929 for (std::size_t i = 0; i < A.rows(); ++i)
3930 for (std::size_t j = 0; j < A.cols(); ++j) A(i, j) = T(A(i, j) * w);
3932 Matrix<T> Q = cur.QN, U = cur.UN, R = cur.RN, Tp = cur.TN, W = cur.WN;
3933 scale(Q); scale(U); scale(R); scale(Tp); scale(W);
3934 for (
int k = 1; k < wnd_size; ++k) {
3935 const LayerResult<T>& old = results[results.size() - 1 - k][e];
3936 auto add = [&](Matrix<T>& A,
const Matrix<T>& B) {
3937 for (std::size_t i = 0; i < A.rows(); ++i)
3938 for (std::size_t j = 0; j < A.cols(); ++j) A(i, j) = T(A(i, j) + B(i, j) * w);
3940 add(Q, old.QN); add(U, old.UN); add(R, old.RN); add(Tp, old.TN); add(W, old.WN);
3942 cur.QN = Q; cur.UN = U; cur.RN = R; cur.TN = Tp; cur.WN = W;
3948 for (std::size_t e = 0; e < E; ++e) {
3949 const Matrix<T>& Q = results[results.size() - 1][e].QN;
3950 const Matrix<T>& Q1 = results[results.size() - 2][e].QN;
3951 const double Njobs = ensemble[e].total_jobs();
3952 if (!(Njobs > 0.0))
continue;
3954 for (std::size_t i = 0; i < Q.rows(); ++i)
3955 for (std::size_t j = 0; j < Q.cols(); ++j)
3956 mx = std::max(mx, std::fabs(dbl(Q(i, j)) - dbl(Q1(i, j))));
3959 maxitererr[it] = err;
3961 static_cast<long>(it),
3962 "layer iteration %zu: max queue-length change %.3e (tolerance %.3e)",
3963 static_cast<std::size_t
>(it), err, opt.iter_tol);
3964 if (it == iter_min) {
3966 averagingstart = it;
3970 if (it > iter_min && maxitererr[it] < opt.iter_tol && maxitererr[it - 1] < opt.iter_tol &&
3971 maxitererr[it - 2] < opt.iter_tol) {
3972 if (!hasconverged) {
3973 hasconverged =
true;
3975 did_converge =
true;
3979 hasconverged =
false;
4014 std::string lnmethod;
4016 bool ph_laws_ready =
false;
4033 static std::string ln_requested_method(
const std::string& method) {
4035 for (
char c : method) m +=
static_cast<char>(std::tolower(
static_cast<unsigned char>(c)));
4036 if (m.empty() || m ==
"srvn" || m ==
"default" || m ==
"auto")
return "srvn";
4037 if (m ==
"srvn.ph" || m ==
"ph")
return "srvn.ph";
4038 if (m ==
"srvn.cs" || m ==
"srvncs" || m ==
"cs")
return "srvn.cs";
4039 if (m ==
"flat.cs" || m ==
"flatcs" || m ==
"flat" || m ==
"squashed")
return "flat.cs";
4040 if (m ==
"flat.ph" || m ==
"flatph" || m ==
"squashed.ph")
return "flat.ph";
4041 if (m ==
"moment3")
return "moment3";
4048 bool is_srvn_ph()
const {
return lnmethod ==
"srvn.ph"; }
4056 bool is_ph_encoding()
const {
return lnmethod ==
"srvn.ph" || lnmethod ==
"flat.ph"; }
4064 static std::size_t station_idx_of(
const qn::Layer<T>& L, std::size_t elem) {
4065 if (elem >= 1 && elem < L.server_idx_of.size() && L.server_idx_of[elem] != 0)
4066 return L.server_idx_of[elem];
4074 std::size_t station_idx_of_class(
const qn::Layer<T>& L, std::size_t k)
const {
4075 const int kind = L.classes[k].attr_kind;
4076 const std::size_t a = L.classes[k].attr_idx;
4077 if (kind ==
int(LqnElement::ACTIVITY))
4078 return station_idx_of(L, lqn.parent[lqn.parent[a]]);
4079 if (kind ==
int(LqnElement::CALL))
4080 return station_idx_of(L, lqn.parent[lqn.callpair_dst[a]]);
4092 bool probe_srvn_ph() {
4094 ph_assert_supported();
4096 ph_laws_ready =
true;
4098 }
catch (
const std::exception&) {
4099 ph_laws_ready =
false;
4106 std::size_t idx = 0;
4107 bool ishost =
false;
4108 std::vector<std::size_t> callers;
4110 std::vector<std::size_t> class_of_caller;
4111 std::size_t nreplicas = 1;
4113 std::vector<std::size_t> qstations;
4115 std::vector<T> svcmean_by_class;
4117 std::vector<std::pair<std::size_t, long>> open_arrivals;
4127 std::vector<workflow::Workflow<T>> ph_wf, ph_wfhost;
4128 std::vector<std::unordered_map<std::size_t, T>> ph_execs;
4129 std::vector<bool> ph_has_wf;
4130 std::vector<workflow::PhLaw<T>> ph_hostlaw, ph_entrylaw;
4131 std::vector<T> ph_hostmean, ph_entrymean, ph_entryscv;
4132 std::vector<T> ph_share, ph_overlap, ph_setupshare, ph_xdemand;
4133 Matrix<T> ph_ncalls, ph_calltime;
4134 std::vector<T> ph_procresid, ph_actthinkt, ph_calltotal;
4135 std::vector<PHLayer> ph_layers;
4136 std::vector<bool> ph_has_layer;
4139 T act_think_time(std::size_t aidx)
const {
4140 if (aidx >= lqn.actthink.size() || lqn.actthink[aidx].disabled)
return Tzero();
4141 const double m = dbl(lqn.actthink[aidx].mean);
4143 return lqn.actthink[aidx].mean;
4164 double setup_charge(std::size_t tidx)
const {
4165 if (tidx >= lqn.hassetup.size() || !lqn.hassetup[tidx])
return 0.0;
4166 const double s = tidx < lqn.setuptime.size() && !lqn.setuptime[tidx].disabled
4167 ? dbl(lqn.setuptime[tidx].mean) : 0.0;
4168 const double d = tidx < lqn.delayofftime.size() && !lqn.delayofftime[tidx].disabled
4169 ? dbl(lqn.delayofftime[tidx].mean) : 0.0;
4171 const double mult = lqn.mult[tidx];
4172 if (!std::isfinite(mult) || mult <= 0.0)
return 0.0;
4173 if (tidx >= tput.size() || tidx >= util.size())
return s;
4174 const double X = dbl(tput[tidx]);
4176 double rho = dbl(util[tidx]);
4177 if (!std::isfinite(rho) || rho < 0.0) rho = 0.0;
4179 const double b = rho * mult;
4180 const double EI = (std::max(1.0, b) - b) / X;
4181 return s * EI / (EI + d);
4190 double ph_setup_prob(std::size_t eidx)
const {
4191 const std::size_t tidx = lqn.parent[eidx];
4192 if (tidx >= lqn.hassetup.size() || !lqn.hassetup[tidx])
return 0.0;
4193 const double s = tidx < lqn.setuptime.size() && !lqn.setuptime[tidx].disabled
4194 ? dbl(lqn.setuptime[tidx].mean) : 0.0;
4195 const double d = tidx < lqn.delayofftime.size() && !lqn.delayofftime[tidx].disabled
4196 ? dbl(lqn.delayofftime[tidx].mean) : 0.0;
4198 return std::min(1.0, std::max(0.0, setup_charge(tidx) / s));
4202 double ph_host_servers(std::size_t hidx)
const {
4203 if (lqn.sched[hidx] == SchedStrategy::INF)
return 1.0;
4204 const double m = lqn.maxmult[hidx];
4205 return (std::isfinite(m) && m > 0.0) ? m : 1.0;
4209 bool ph_any_caller_of(std::size_t eidx)
const {
4210 return lqn.issynccaller.any_col(eidx) || lqn.isasynccaller.any_col(eidx);
4218 bool ph_open_arrival_only(std::size_t tidx)
const {
4219 if (lqn.isref[tidx])
return false;
4220 for (std::size_t eidx : lqn.entriesof[tidx])
4221 if (ph_any_caller_of(eidx))
return false;
4222 for (std::size_t eidx : lqn.entriesof[tidx])
4223 if (lqn.has_arrival[eidx])
return true;
4228 std::vector<std::size_t> ph_async_calls_into(std::size_t tidx)
const {
4229 std::vector<std::size_t> out;
4230 for (std::size_t cidx = 1; cidx <= lqn.ncalls; ++cidx) {
4231 if (lqn.calltype[cidx] != CallType::ASYNC)
continue;
4232 for (std::size_t e : lqn.entriesof[tidx])
4233 if (lqn.callpair_dst[cidx] == e) { out.push_back(cidx);
break; }
4242 void ph_assert_supported(
bool flat =
false)
const {
4245 const std::string mname = flat ?
"flat.ph" :
"srvn.ph";
4251 "method='" + mname +
"' does not support second-phase activities: the composed "
4252 "entry law has no reply point. Use method='default'.");
4258 const std::string why = ph_method_refusal(mname);
4277 std::string ph_method_refusal(
const std::string& method)
const {
4278 if (method !=
"srvn.ph" && method !=
"flat.ph")
return std::string();
4281 if (method ==
"flat.ph") {
4282 for (std::size_t i = 1; i <= NT(); ++i) {
4283 if (lqn.repl[i] > 1.0)
4284 return "method='flat.ph' does not support replicated processors or tasks, "
4285 "whose replicas need a submodel each. Use method='srvn.ph'.";
4286 if (lqn.hassetup[i])
4287 return "method='flat.ph' does not support setup tasks, whose powered-down "
4288 "threads are per-layer state. Use method='srvn.ph'.";
4291 for (std::size_t cidx = 1; cidx <= lqn.ncalls; ++cidx)
4292 if (lqn.calltype[cidx] == CallType::FWD)
4293 return "method='" + method +
4294 "' does not support forwarding calls, whose target is not part of the "
4295 "caller's activity graph. Use method='default'.";
4296 for (std::size_t i = 1; i < lqn.iscache.size(); ++i)
4298 return "method='" + method +
4299 "' does not support cache tasks. Use method='default'.";
4300 for (std::size_t i = 1; i < lqn.hassetup.size(); ++i) {
4301 if (!lqn.hassetup[i])
continue;
4302 if (lqn.sched[i] == SchedStrategy::INF || !std::isfinite(lqn.mult[i]))
4303 return "method='" + method +
"': task '" + lqn.names[i] +
4304 "' declares a setup time on an infinite-server task, which holds no "
4305 "thread to power down; give it a finite multiplicity.";
4307 for (std::size_t i = 0; i < lqn.lincon_A.size(); ++i)
4308 if (lqn.lincon_A[i].rows() > 0)
4309 return "method='" + method +
4310 "' does not support admission constraints on a layer station. "
4311 "Use method='default'.";
4312 for (std::size_t i = 1; i < lqn.lldscaling.size(); ++i) {
4313 const char* fname =
nullptr;
4314 if (!lqn.lldscaling[i].empty()) fname =
"a load dependence";
4315 else if (lqn.cdscaling[i]) fname =
"a class dependence";
4316 else if (lqn.jdscaling[i]) fname =
"a joint dependence";
4317 else if (!lqn.pools[i].empty()) fname =
"server pools";
4318 if (fname !=
nullptr)
4319 return "method='" + method +
4320 "' does not support queue-dependent service rates on a layer station ('" +
4321 lqn.names[i] +
"' declares " + fname +
"). Use method='srvn.cs'.";
4323 return std::string();
4327 void ph_init_laws() {
4328 const std::size_t N = lqn.nidx;
4329 ph_wf.assign(N + 1, workflow::Workflow<T>(
"empty"));
4330 ph_wfhost.assign(N + 1, workflow::Workflow<T>(
"empty"));
4331 ph_execs.assign(N + 1, {});
4332 ph_has_wf.assign(N + 1,
false);
4333 ph_hostlaw.assign(N + 1, workflow::PhLaw<T>());
4334 ph_entrylaw.assign(N + 1, workflow::PhLaw<T>());
4335 ph_hostmean.assign(N + 1, Tzero());
4336 ph_entrymean.assign(N + 1, Tzero());
4337 ph_entryscv.assign(N + 1, Tone());
4338 ph_share.assign(N + 1, Tzero());
4339 ph_overlap.assign(N + 1, Tone());
4340 ph_setupshare.assign(N + 1, Tzero());
4341 ph_xdemand.assign(N + 1, Tzero());
4342 ph_ncalls =
Matrix<T>(N + 1, N + 1, Tzero());
4343 ph_calltime =
Matrix<T>(N + 1, N + 1, Tzero());
4344 ph_procresid.assign(N + 1, Tzero());
4345 ph_actthinkt.assign(N + 1, Tzero());
4346 ph_calltotal.assign(N + 1, Tzero());
4347 ph_layers.assign(NT() + 1, PHLayer());
4348 ph_has_layer.assign(NT() + 1,
false);
4350 for (std::size_t e = 1; e <= lqn.nentries; ++e) {
4351 const std::size_t eidx = lqn.eshift + e;
4352 const std::size_t tidx = lqn.parent[eidx];
4353 if (ignore[tidx])
continue;
4355 ph_wf[eidx] = std::move(ew.wf);
4356 ph_execs[eidx] = ew.execs;
4357 ph_has_wf[eidx] =
true;
4359 ph_wfhost[eidx] = std::move(eh.wf);
4369 for (std::size_t t = 1; t <= lqn.ntasks; ++t) {
4370 const std::size_t tidx = lqn.tshift + t;
4371 const std::size_t n = lqn.entriesof[tidx].size();
4372 if (n == 0)
continue;
4373 for (std::size_t eidx : lqn.entriesof[tidx])
4374 ph_share[eidx] = num_traits<T>::from_double(1.0 /
double(n));
4385 Distrib<T> ph_call_burst_law(std::size_t cidx)
const {
4386 const double m = dbl(lqn.callproc_mean[cidx]);
4387 const std::size_t eidx = lqn.callpair_dst[cidx];
4389 const T R = T(callservt[cidx] / lqn.callproc_mean[cidx]);
4390 double scvd = dbl(ph_entryscv[eidx]);
4396 const workflow::PhLaw<T> loop =
4410 Distrib<T> ph_station_law(
const workflow::PhLaw<T>& law)
const {
4418 static workflow::PhLaw<T> ph_immediate_law() {
4419 workflow::PhLaw<T> out;
4420 out.alpha.assign(1, Tone());
4434 void ph_compose_entry_laws() {
4435 const std::size_t N = lqn.nidx;
4436 std::vector<T> entry_setup_share(N + 1, Tzero());
4437 ph_overlap.assign(N + 1, Tone());
4439 for (std::size_t e = 1; e <= lqn.nentries; ++e) {
4440 const std::size_t eidx = lqn.eshift + e;
4441 if (!ph_has_wf[eidx])
continue;
4442 workflow::Workflow<T>& w = ph_wf[eidx];
4443 const std::unordered_map<std::size_t, T>& ex = ph_execs[eidx];
4444 T entrysum = Tzero(), procsum = Tzero();
4445 for (std::size_t aidx : lqn.actsof[eidx]) {
4446 T m = T(residt[aidx] + act_think_time(aidx));
4447 const T xa = ex.at(aidx);
4448 procsum += T(xa * m);
4449 const double md = dbl(m);
4450 w.set_activity_demand_mean(
4454 for (std::size_t cidx : lqn.callsof[aidx]) {
4455 if (lqn.calltype[cidx] != CallType::SYNC)
continue;
4456 w.set_activity_demand(lqn.callhashnames[cidx], ph_call_burst_law(cidx));
4457 m = T(m + callservt[cidx]);
4459 entrysum += T(xa * m);
4463 T m1 = mm.first, scv = mm.second;
4470 const T f = T(m1 / procsum);
4471 for (std::size_t i = 0; i < law.S.rows(); ++i)
4472 for (std::size_t j = 0; j < law.S.cols(); ++j) law.S(i, j) = T(law.S(i, j) * f);
4481 const double p = ph_setup_prob(eidx);
4483 const std::size_t tidx = lqn.parent[eidx];
4484 const double sm = dbl(lqn.setuptime[tidx].mean);
4485 double sscv = dbl(lqn.setuptime[tidx].scv);
4489 num_traits<T>::from_double(sm), num_traits<T>::from_double(sscv)));
4490 std::vector<workflow::PhLaw<T>> mix;
4493 std::vector<T> probs;
4494 probs.push_back(num_traits<T>::from_double(p));
4495 probs.push_back(num_traits<T>::from_double(1.0 - p));
4505 entry_setup_share[eidx] = num_traits<T>::from_double(p * sm / denom);
4508 ph_entrylaw[eidx] = law;
4509 ph_entrymean[eidx] = m1;
4510 ph_entryscv[eidx] = scv;
4512 const double r = std::min(1.0, dbl(m1) / dbl(entrysum));
4513 ph_overlap[eidx] = num_traits<T>::from_double(r);
4519 ph_setupshare.assign(N + 1, Tzero());
4520 for (std::size_t t = 1; t <= lqn.ntasks; ++t) {
4521 const std::size_t tidx = lqn.tshift + t;
4522 if (ignore[tidx])
continue;
4523 for (std::size_t eidx : lqn.entriesof[tidx])
4524 ph_setupshare[tidx] = T(ph_setupshare[tidx] + ph_share[eidx] * entry_setup_share[eidx]);
4528 ph_ncalls =
Matrix<T>(N + 1, N + 1, Tzero());
4529 ph_calltime =
Matrix<T>(N + 1, N + 1, Tzero());
4530 ph_procresid.assign(N + 1, Tzero());
4531 ph_actthinkt.assign(N + 1, Tzero());
4532 ph_calltotal.assign(N + 1, Tzero());
4533 for (std::size_t t = 1; t <= lqn.ntasks; ++t) {
4534 const std::size_t tidx = lqn.tshift + t;
4535 if (ignore[tidx])
continue;
4536 for (std::size_t eidx : lqn.entriesof[tidx]) {
4537 if (!ph_has_wf[eidx])
continue;
4538 const T w = ph_share[eidx];
4539 if (!(dbl(w) > 0.0))
continue;
4540 const std::unordered_map<std::size_t, T>& ex = ph_execs[eidx];
4541 const T r = ph_overlap[eidx];
4542 for (std::size_t aidx : lqn.actsof[eidx]) {
4543 const T xa = ex.at(aidx);
4544 ph_procresid[tidx] = T(ph_procresid[tidx] + w * r * xa * residt[aidx]);
4545 ph_actthinkt[tidx] = T(ph_actthinkt[tidx] + w * r * xa * act_think_time(aidx));
4546 for (std::size_t cidx : lqn.callsof[aidx]) {
4547 if (lqn.calltype[cidx] != CallType::SYNC)
continue;
4548 const std::size_t tgte = lqn.callpair_dst[cidx];
4549 const std::size_t tgtt = lqn.parent[tgte];
4552 ph_ncalls(tidx, tgte) =
4553 T(ph_ncalls(tidx, tgte) + w * xa * lqn.callproc_mean[cidx]);
4554 ph_calltime(tidx, tgtt) =
4555 T(ph_calltime(tidx, tgtt) + w * r * xa * callservt[cidx]);
4556 ph_calltotal[tidx] = T(ph_calltotal[tidx] + w * r * xa * callservt[cidx]);
4564 void build_layers_ph(
bool flat =
false) {
4565 if (!ph_laws_ready) ph_assert_supported(flat);
4569 opt.interlocking =
false;
4570 interlockMethod =
"none";
4574 if (!ph_laws_ready) ph_init_laws();
4577 const std::size_t N = lqn.nidx;
4578 residt.assign(N + 1, Tzero());
4579 servt.assign(N + 1, Tzero());
4580 callservt.assign(lqn.ncalls + 1, Tzero());
4581 callresidt.assign(lqn.ncalls + 1, Tzero());
4582 for (std::size_t aidx = lqn.ashift + 1; aidx <= lqn.ashift + lqn.nacts; ++aidx)
4583 residt[aidx] = lqn.hostdem[aidx].disabled ? Tzero() : lqn.hostdem[aidx].mean;
4584 for (std::size_t cidx = 1; cidx <= lqn.ncalls; ++cidx) {
4585 if (lqn.calltype[cidx] != CallType::SYNC && lqn.calltype[cidx] != CallType::ASYNC)
4587 const std::size_t eidx = lqn.callpair_dst[cidx];
4588 callservt[cidx] = T(lqn.callproc_mean[cidx] * ph_hostmean[eidx]);
4589 callresidt[cidx] = callservt[cidx];
4591 ph_compose_entry_laws();
4595 build_ph_flat_layer();
4596 tput.assign(N + 1, Tzero());
4597 util.assign(N + 1, Tzero());
4598 thinkt.assign(N + 1, Tzero());
4599 update_layers_ph(0);
4603 std::vector<qn::Layer<T>> raw(NT() + 1);
4604 std::vector<bool> present(NT() + 1,
false);
4606 for (std::size_t hidx = 1; hidx <= lqn.nhosts; ++hidx) {
4607 if (ignore[hidx])
continue;
4608 std::vector<std::size_t> callers;
4609 for (std::size_t tidx : lqn.tasksof[hidx]) {
4610 if (ignore[tidx])
continue;
4611 if (lqn.isref[tidx]) { callers.push_back(tidx);
continue; }
4612 for (std::size_t eidx : lqn.entriesof[tidx])
4613 if (ph_any_caller_of(eidx) || lqn.has_arrival[eidx]) {
4614 callers.push_back(tidx);
4618 if (callers.empty())
continue;
4619 build_layer_ph(raw[hidx], hidx, callers,
true);
4620 present[hidx] =
true;
4622 for (std::size_t t = 1; t <= lqn.ntasks; ++t) {
4623 const std::size_t tidx = lqn.tshift + t;
4624 if (ignore[tidx] || lqn.isref[tidx])
continue;
4625 std::vector<std::size_t> callers;
4626 for (std::size_t ct = 1; ct <= lqn.ntasks; ++ct) {
4627 const std::size_t c = lqn.tshift + ct;
4628 if (c == tidx || ignore[c])
continue;
4629 for (std::size_t e : lqn.entriesof[tidx])
4630 if (lqn.issynccaller.get(c, e)) { callers.push_back(c);
break; }
4632 if (callers.empty() && ph_async_calls_into(tidx).empty())
continue;
4633 build_layer_ph(raw[tidx], tidx, callers,
false);
4634 present[tidx] =
true;
4637 idxhash.assign(lqn.nidx + 1, -1);
4639 for (std::size_t i = 1; i <= NT(); ++i)
4641 idxhash[i] = next++;
4642 ensemble.push_back(std::move(raw[i]));
4644 layer_init_sol.assign(ensemble.size(), Matrix<T>());
4649 tput.assign(N + 1, Tzero());
4650 util.assign(N + 1, Tzero());
4651 thinkt.assign(N + 1, Tzero());
4652 update_layers_ph(0);
4668 void build_ph_flat_layer() {
4669 const std::vector<std::size_t> servers = ph_flat_server_set();
4670 const std::size_t nsrv = servers.size();
4674 m.clientIdx = m.add_station(qn::Station<T>{
4675 "Clients", NodeType::Delay, SchedStrategy::INF,
4676 std::numeric_limits<double>::infinity(),
false, 0});
4677 const std::size_t clientNode = m.node_of_station(m.clientIdx);
4678 m.server_idx_of.assign(lqn.nidx + 1, 0);
4680 std::vector<std::size_t> station_of(lqn.nidx + 1, 0);
4681 std::vector<std::size_t> serverNode(nsrv, 0);
4682 std::vector<std::size_t> srvStation(nsrv, 0);
4683 for (std::size_t si = 0; si < nsrv; ++si) {
4684 const std::size_t idx = servers[si];
4685 const bool ishost = idx <= lqn.nhosts;
4687 st.name = lqn.hashnames[idx];
4688 st.nodetype = lqn.sched[idx] == SchedStrategy::INF ? NodeType::Delay : NodeType::Queue;
4689 st.sched = lqn.sched[idx];
4690 st.nservers = lqn.sched[idx] == SchedStrategy::INF
4691 ? std::numeric_limits<double>::infinity()
4693 st.attr_ishost = ishost;
4695 const std::size_t s = m.add_station(st);
4697 station_of[idx] = s;
4698 serverNode[si] = m.node_of_station(s);
4699 m.server_idx_of[idx] = s;
4701 m.host_stations.push_back(s);
4703 m.task_stations.push_back(s);
4706 m.serverIdx = station_of[servers[0]];
4710 std::vector<std::vector<std::size_t>> callers_of(lqn.nidx + 1);
4711 std::vector<std::size_t> all_callers;
4712 for (std::size_t si = 0; si < nsrv; ++si) {
4713 const std::size_t idx = servers[si];
4714 std::vector<std::size_t> cs;
4715 if (idx <= lqn.nhosts) {
4716 for (std::size_t tidx : lqn.tasksof[idx]) {
4717 if (ignore[tidx])
continue;
4718 if (lqn.isref[tidx]) { cs.push_back(tidx);
continue; }
4719 for (std::size_t eidx : lqn.entriesof[tidx])
4720 if (ph_any_caller_of(eidx) || lqn.has_arrival[eidx]) {
4726 for (std::size_t ct = 1; ct <= lqn.ntasks; ++ct) {
4727 const std::size_t c = lqn.tshift + ct;
4728 if (c == idx || ignore[c])
continue;
4729 for (std::size_t e : lqn.entriesof[idx])
4730 if (lqn.issynccaller.get(c, e)) { cs.push_back(c);
break; }
4733 callers_of[idx] = cs;
4734 for (std::size_t c : cs)
4735 if (std::find(all_callers.begin(), all_callers.end(), c) == all_callers.end())
4736 all_callers.push_back(c);
4738 std::sort(all_callers.begin(), all_callers.end());
4741 std::vector<std::size_t> class_of_caller(lqn.nidx + 1, 0);
4743 for (std::size_t c : all_callers) {
4746 double nj = lqn.maxmult[c];
4747 if (std::isinf(nj)) {
4749 for (std::size_t k = 1; k <= NT(); ++k)
4750 if (lqn.taskgraph.get(k, c) != Tzero()) sacc += lqn.maxmult[k];
4752 if (std::isinf(nj) || nj == 0.0) {
4754 for (std::size_t k = 1; k <= NT(); ++k)
4755 if (std::isfinite(lqn.maxmult[k])) s2 += lqn.maxmult[k];
4756 nj = std::min(s2, 1000.0);
4760 jc.name = lqn.hashnames[c];
4761 jc.type = JobClassType::CLOSED;
4763 jc.refstat = m.clientIdx;
4764 jc.is_ref_class =
true;
4765 jc.attr_kind = int(LqnElement::TASK);
4767 const std::size_t k = m.add_class(jc);
4768 class_of_caller[c] = k;
4769 m.attr_tasks.emplace_back(k, c);
4771 const double zt = dbl(ref_think_time(c));
4772 m.set_service(m.clientIdx, k,
4781 for (std::size_t si = 0; si < nsrv; ++si)
4783 for (std::size_t si = 0; si < nsrv; ++si) {
4784 const std::size_t idx = servers[si];
4785 const std::vector<std::size_t>& cs = callers_of[idx];
4786 if (std::find(cs.begin(), cs.end(), c) == cs.end())
continue;
4787 m.set_service(srvStation[si], k,
4791 thinkt_map.push_back({idx, c, m.clientIdx, k});
4792 servt_map.push_back({idx, c, station_of[idx], k});
4797 std::vector<std::vector<std::pair<std::size_t, long>>> open_of(lqn.nidx + 1);
4798 std::size_t sourceStation = 0, sinkNode = 0;
4799 for (std::size_t si = 0; si < nsrv; ++si) {
4800 const std::size_t hidx = servers[si];
4801 if (hidx > lqn.nhosts)
continue;
4802 for (std::size_t c : callers_of[hidx]) {
4808 if (ph_open_arrival_only(c))
continue;
4809 for (std::size_t eidx : lqn.entriesof[c]) {
4810 if (!lqn.has_arrival[eidx])
continue;
4811 if (sourceStation == 0) {
4812 sourceStation = m.add_station(qn::Station<T>{
4813 "Source", NodeType::Source, SchedStrategy::EXT,
4814 std::numeric_limits<double>::infinity(),
false, 0});
4815 m.sourceIdx = sourceStation;
4816 sinkNode = m.add_node(
"Sink", NodeType::Sink,
false);
4817 m.sinkNode = sinkNode;
4820 oc.name = lqn.hashnames[eidx] +
".Open";
4821 oc.type = JobClassType::OPEN;
4822 oc.population = std::numeric_limits<double>::infinity();
4823 oc.refstat = sourceStation;
4824 oc.attr_kind = int(LqnElement::ENTRY);
4826 const std::size_t k = m.add_class(oc);
4827 m.set_service(sourceStation, k, lqn.arrival[eidx]);
4829 for (std::size_t s2 = 0; s2 < nsrv; ++s2)
4832 m.set_service(srvStation[si], k,
4834 open_of[hidx].emplace_back(k,
long(eidx));
4835 m.attr_entries.emplace_back(k, eidx);
4839 for (std::size_t si = 0; si < nsrv; ++si) {
4840 const std::size_t tidx = servers[si];
4841 if (tidx <= lqn.nhosts)
continue;
4842 for (std::size_t cidx : ph_async_calls_into(tidx)) {
4843 if (sourceStation == 0) {
4844 sourceStation = m.add_station(qn::Station<T>{
4845 "Source", NodeType::Source, SchedStrategy::EXT,
4846 std::numeric_limits<double>::infinity(),
false, 0});
4847 m.sourceIdx = sourceStation;
4848 sinkNode = m.add_node(
"Sink", NodeType::Sink,
false);
4849 m.sinkNode = sinkNode;
4851 const std::size_t eidx = lqn.callpair_dst[cidx];
4853 oc.name = lqn.callhashnames[cidx];
4854 oc.type = JobClassType::OPEN;
4855 oc.population = std::numeric_limits<double>::infinity();
4856 oc.refstat = sourceStation;
4857 oc.attr_kind = int(LqnElement::CALL);
4859 const std::size_t k = m.add_class(oc);
4862 for (std::size_t s2 = 0; s2 < nsrv; ++s2)
4865 m.set_service(srvStation[si], k,
4867 open_of[tidx].emplace_back(k, -
long(cidx));
4868 m.attr_calls.push_back({k, cidx, lqn.callpair_src[cidx], eidx});
4869 arv_call_map.push_back({tidx, cidx, sourceStation, k});
4870 call_map.push_back({tidx, cidx, station_of[tidx], k});
4877 for (std::size_t c : all_callers) {
4878 const std::size_t k = class_of_caller[c];
4879 std::size_t prev = clientNode;
4880 bool visited =
false;
4881 for (std::size_t si = 0; si < nsrv; ++si) {
4882 const std::vector<std::size_t>& cs = callers_of[servers[si]];
4883 if (std::find(cs.begin(), cs.end(), c) == cs.end())
continue;
4884 m.set_route(k, k, prev, serverNode[si], Tone());
4885 prev = serverNode[si];
4888 if (visited) m.set_route(k, k, prev, clientNode, Tone());
4890 if (sourceStation != 0) {
4891 const std::size_t srcNode = m.node_of_station(sourceStation);
4892 for (std::size_t si = 0; si < nsrv; ++si)
4893 for (
const auto& oa : open_of[servers[si]]) {
4894 m.set_route(oa.first, oa.first, srcNode, serverNode[si], Tone());
4895 m.set_route(oa.first, oa.first, serverNode[si], sinkNode, Tone());
4898 apply_host_task_priorities(m, servers);
4901 idxhash.assign(lqn.nidx + 1, -1);
4902 for (std::size_t idx : servers) idxhash[idx] = 0;
4904 ensemble.push_back(std::move(m));
4905 layer_init_sol.assign(ensemble.size(), Matrix<T>());
4908 const std::size_t nclasses = ensemble[0].classes.size();
4909 for (std::size_t idx : servers) {
4912 L.ishost = idx <= lqn.nhosts;
4913 L.callers = callers_of[idx];
4914 L.class_of_caller = class_of_caller;
4916 L.qstations.assign(1, station_of[idx]);
4917 L.svcmean_by_class.assign(nclasses + 1, Tzero());
4918 L.open_arrivals = open_of[idx];
4919 L.npop = npop < 1.0 ? 1.0 : npop;
4921 ph_has_layer[idx] =
true;
4934 std::vector<std::size_t> ph_flat_server_set()
const {
4935 for (std::size_t i = 1; i <= NT(); ++i) {
4936 if (lqn.repl[i] > 1.0)
4938 "method='flat.ph' does not support replicated processors or tasks, whose "
4939 "replicas need a submodel each. Use method='srvn.ph'.");
4942 "method='flat.ph' does not support cache tasks. Use method='default'.");
4943 if (lqn.hassetup[i])
4945 "method='flat.ph' does not support setup tasks, whose powered-down threads "
4946 "are per-layer state. Use method='srvn.ph'.");
4948 std::vector<std::size_t> servers;
4949 for (std::size_t hidx = 1; hidx <= lqn.nhosts; ++hidx) {
4950 if (ignore[hidx] || lqn.tasksof[hidx].empty())
continue;
4952 for (std::size_t tidx : lqn.tasksof[hidx]) {
4953 if (ignore[tidx])
continue;
4954 if (lqn.isref[tidx]) { any =
true;
break; }
4955 for (std::size_t eidx : lqn.entriesof[tidx])
4956 if (ph_any_caller_of(eidx) || lqn.has_arrival[eidx]) { any =
true;
break; }
4959 if (any) servers.push_back(hidx);
4961 for (std::size_t t = 1; t <= lqn.ntasks; ++t) {
4962 const std::size_t tidx = lqn.tshift + t;
4963 if (ignore[tidx] || lqn.isref[tidx])
continue;
4964 bool has_caller =
false;
4965 for (std::size_t ct = 1; ct <= lqn.ntasks && !has_caller; ++ct) {
4966 const std::size_t c = lqn.tshift + ct;
4967 if (c == tidx || ignore[c])
continue;
4968 for (std::size_t e : lqn.entriesof[tidx])
4969 if (lqn.issynccaller.get(c, e)) { has_caller =
true;
break; }
4971 if (!has_caller && ph_async_calls_into(tidx).empty())
continue;
4972 servers.push_back(tidx);
4974 if (servers.empty())
4976 "method='flat.ph' found no server: the model has no processor with tasks.");
4981 void build_layer_ph(qn::Layer<T>& m, std::size_t idx,
4982 const std::vector<std::size_t>& callers,
bool ishost) {
4983 m.name = lqn.hashnames[idx];
4988 const double rawrepl = lqn.repl[idx];
4989 std::size_t nreplicas = 1;
4990 if (rawrepl > 1.0 && !callers.empty()) {
4993 for (std::size_t c : callers)
4994 if (lqn.repl[c] != rawrepl) reduce =
false;
4996 for (std::size_t c : callers)
4997 if (lqn.fanout_at(c, idx) < rawrepl) reduce =
false;
4999 nreplicas = reduce ? 1 :
static_cast<std::size_t
>(std::llround(rawrepl));
5000 if (reduce && !ishost) single_replica_tasks.insert(idx);
5002 const bool reduce_fanout = (nreplicas == 1 && rawrepl > 1.0 && !callers.empty());
5004 m.clientIdx = m.add_station(qn::Station<T>{
5005 "Clients", NodeType::Delay, SchedStrategy::INF,
5006 std::numeric_limits<double>::infinity(),
false, 0});
5007 const std::size_t clientNode = m.node_of_station(m.clientIdx);
5008 m.serverIdx = m.clientIdx + 1;
5012 L.callers = callers;
5013 L.nreplicas = nreplicas;
5014 L.class_of_caller.assign(lqn.nidx + 1, 0);
5015 std::vector<std::size_t> serverNode(nreplicas);
5016 for (std::size_t r = 0; r < nreplicas; ++r) {
5018 st.name = r == 0 ? lqn.hashnames[idx] : lqn.hashnames[idx] +
"." + std::to_string(r + 1);
5019 st.nodetype = lqn.sched[idx] == SchedStrategy::INF ? NodeType::Delay : NodeType::Queue;
5020 st.sched = lqn.sched[idx];
5021 st.nservers = lqn.sched[idx] == SchedStrategy::INF
5022 ? std::numeric_limits<double>::infinity()
5024 st.attr_ishost = ishost;
5026 const std::size_t s = m.add_station(st);
5027 L.qstations.push_back(s);
5028 serverNode[r] = m.node_of_station(s);
5032 for (std::size_t c : callers) {
5033 double nj = njobs(c, idx);
5035 const bool caller_single_replica =
5036 reduce_fanout || single_replica_tasks.count(c) > 0;
5037 nj = caller_single_replica ? lqn.maxmult[c] : lqn.maxmult[c] * lqn.repl[c];
5038 if (std::isinf(nj)) {
5040 for (std::size_t k = 1; k <= NT(); ++k)
5041 if (lqn.taskgraph.get(k, c) != Tzero()) s += lqn.maxmult[k];
5043 if (std::isinf(nj) || nj == 0.0) {
5045 for (std::size_t k = 1; k <= NT(); ++k)
5046 if (std::isfinite(lqn.maxmult[k])) s2 += lqn.maxmult[k] * lqn.repl[k];
5047 nj = std::min(s2, 1000.0);
5053 jc.name = lqn.hashnames[c];
5054 jc.type = JobClassType::CLOSED;
5056 jc.refstat = m.clientIdx;
5057 jc.is_ref_class =
true;
5058 jc.attr_kind = int(LqnElement::TASK);
5060 const std::size_t k = m.add_class(jc);
5061 L.class_of_caller[c] = k;
5062 m.attr_tasks.emplace_back(k, c);
5063 const double zt = dbl(ref_think_time(c));
5064 m.set_service(m.clientIdx, k,
5067 for (std::size_t s : L.qstations)
5072 thinkt_map.push_back({idx, c, m.clientIdx, k});
5073 servt_map.push_back({idx, c, L.qstations[0], k});
5077 std::size_t sourceStation = 0, sinkNode = 0;
5079 for (std::size_t c : callers) {
5086 if (ph_open_arrival_only(c))
continue;
5087 for (std::size_t eidx : lqn.entriesof[c]) {
5088 if (!lqn.has_arrival[eidx])
continue;
5089 if (sourceStation == 0) {
5090 sourceStation = m.add_station(qn::Station<T>{
5091 "Source", NodeType::Source, SchedStrategy::EXT,
5092 std::numeric_limits<double>::infinity(),
false, 0});
5093 m.sourceIdx = sourceStation;
5094 sinkNode = m.add_node(
"Sink", NodeType::Sink,
false);
5095 m.sinkNode = sinkNode;
5098 oc.name = lqn.hashnames[eidx] +
".Open";
5099 oc.type = JobClassType::OPEN;
5100 oc.population = std::numeric_limits<double>::infinity();
5101 oc.refstat = sourceStation;
5102 oc.attr_kind = int(LqnElement::ENTRY);
5104 const std::size_t k = m.add_class(oc);
5105 m.set_service(sourceStation, k, lqn.arrival[eidx]);
5106 for (std::size_t s : L.qstations) {
5110 L.open_arrivals.emplace_back(k,
long(eidx));
5111 m.attr_entries.emplace_back(k, eidx);
5115 for (std::size_t cidx : ph_async_calls_into(idx)) {
5116 if (sourceStation == 0) {
5117 sourceStation = m.add_station(qn::Station<T>{
5118 "Source", NodeType::Source, SchedStrategy::EXT,
5119 std::numeric_limits<double>::infinity(),
false, 0});
5120 m.sourceIdx = sourceStation;
5121 sinkNode = m.add_node(
"Sink", NodeType::Sink,
false);
5122 m.sinkNode = sinkNode;
5124 const std::size_t eidx = lqn.callpair_dst[cidx];
5126 oc.name = lqn.callhashnames[cidx];
5127 oc.type = JobClassType::OPEN;
5128 oc.population = std::numeric_limits<double>::infinity();
5129 oc.refstat = sourceStation;
5130 oc.attr_kind = int(LqnElement::CALL);
5132 const std::size_t k = m.add_class(oc);
5134 for (std::size_t s : L.qstations) {
5138 L.open_arrivals.emplace_back(k, -
long(cidx));
5139 m.attr_calls.push_back({k, cidx, lqn.callpair_src[cidx], eidx});
5140 arv_call_map.push_back({idx, cidx, sourceStation, k});
5141 call_map.push_back({idx, cidx, L.qstations[0], k});
5147 const T share = num_traits<T>::from_double(1.0 /
double(nreplicas));
5148 for (std::size_t c : callers) {
5149 const std::size_t k = L.class_of_caller[c];
5150 for (std::size_t r = 0; r < nreplicas; ++r) {
5151 m.set_route(k, k, clientNode, serverNode[r], share);
5152 m.set_route(k, k, serverNode[r], clientNode, Tone());
5155 for (
const auto& oa : L.open_arrivals) {
5156 const std::size_t k = oa.first;
5157 const std::size_t srcNode = m.node_of_station(sourceStation);
5158 for (std::size_t r = 0; r < nreplicas; ++r) {
5159 m.set_route(k, k, srcNode, serverNode[r], share);
5160 m.set_route(k, k, serverNode[r], sinkNode, Tone());
5163 L.svcmean_by_class.assign(m.classes.size() + 1, Tzero());
5165 for (std::size_t c : callers) {
5166 const double v = njobs(c, idx);
5167 if (std::isfinite(v) && v > 0.0) np += v;
5169 L.npop = np < 1.0 ? 1.0 : np;
5170 apply_host_task_priorities(m, std::vector<std::size_t>{idx});
5173 ph_has_layer[idx] =
true;
5177 workflow::PhLaw<T> ph_service_law(std::size_t idx,
bool ishost, std::size_t c)
const {
5180 std::vector<workflow::PhLaw<T>> laws;
5181 std::vector<T> probs;
5183 for (std::size_t eidx : lqn.entriesof[c]) {
5184 if (ph_hostlaw[eidx].S.rows() == 0 || !(dbl(ph_share[eidx]) > 0.0))
continue;
5185 laws.push_back(ph_hostlaw[eidx]);
5186 probs.push_back(ph_share[eidx]);
5187 tot += dbl(ph_share[eidx]);
5189 if (laws.empty())
return ph_immediate_law();
5190 for (T& p : probs) p = T(p / num_traits<T>::from_double(tot));
5196 bool started =
false;
5197 workflow::PhLaw<T> out;
5198 for (std::size_t eidx : lqn.entriesof[idx]) {
5199 const T n = ph_ncalls(c, eidx);
5201 const workflow::PhLaw<T> lp =
5206 if (!started)
return ph_immediate_law();
5214 T ph_delay_mean(std::size_t idx, std::size_t c)
const {
5224 T z = T(thinkt[c] + ref_think_time(c));
5225 if (!std::isfinite(dbl(z)) || dbl(z) < 0.0) z = Tzero();
5226 z = T(z + ph_actthinkt[c]);
5227 const std::size_t hidx = lqn.parent[c];
5228 if (!ph_served_here(idx, hidx)) z = T(z + ph_procresid[c]);
5229 for (std::size_t t = 1; t <= lqn.ntasks; ++t) {
5230 const std::size_t tidx = lqn.tshift + t;
5231 if (!ph_served_here(idx, tidx)) z = T(z + ph_calltime(c, tidx));
5233 if (!std::isfinite(dbl(z)) || dbl(z) < 0.0)
5243 bool ph_served_here(std::size_t idx, std::size_t elem)
const {
5244 if (elem < 1 || elem >= idxhash.size() || idx < 1 || idx >= idxhash.size())
5246 return idxhash[elem] >= 0 && idxhash[elem] == idxhash[idx];
5258 void update_layers_ph(
int) {
5259 for (std::size_t idx = 1; idx <= NT(); ++idx) {
5260 if (idxhash[idx] < 0 || !ph_has_layer[idx])
continue;
5261 PHLayer& L = ph_layers[idx];
5262 qn::Layer<T>& m = ensemble[std::size_t(idxhash[idx])];
5263 for (std::size_t c : L.callers) {
5264 const std::size_t k = L.class_of_caller[c];
5265 const workflow::PhLaw<T> sl = ph_service_law(idx, L.ishost, c);
5267 const Distrib<T> law = ph_station_law(sl);
5268 for (std::size_t s : L.qstations) m.set_service(s, k, law);
5269 const double zd = std::max(dbl(ph_delay_mean(idx, c)),
5271 m.set_service(m.clientIdx, k,
5274 for (
const auto& oa : L.open_arrivals) {
5275 const std::size_t k = oa.first;
5276 if (oa.second > 0) {
5278 L.svcmean_by_class[k] = ph_hostmean[std::size_t(oa.second)];
5281 const std::size_t cidx = std::size_t(-oa.second);
5282 const std::size_t eidx = lqn.callpair_dst[cidx];
5283 L.svcmean_by_class[k] = ph_entrymean[eidx];
5284 const Distrib<T> law = ph_station_law(ph_entrylaw[eidx]);
5285 for (std::size_t s : L.qstations) m.set_service(s, k, law);
5286 const std::size_t aidx = lqn.callpair_src[cidx];
5287 double rate = dbl(tput[aidx]) * dbl(lqn.callproc_mean[cidx]);
5290 m.set_service(m.sourceIdx, k,
5303 static T ph_residence(
const T& Q,
const T& X,
const T& RN) {
5304 const double q = num_traits<T>::to_double(Q), x = num_traits<T>::to_double(X);
5315 static T ph_inflation_of(
const T& R,
const T& S,
double npop) {
5317 const double s = num_traits<T>::to_double(S), r = num_traits<T>::to_double(R);
5319 if (!std::isfinite(f) || f < 1.0) f = 1.0;
5320 if (std::isfinite(npop) && npop >= 1.0 && f > npop) f = npop;
5321 return num_traits<T>::from_double(f);
5325 double ph_layer_pop(
const PHLayer& L, std::size_t idx)
const {
5328 if (L.npop >= 1.0)
return L.npop;
5330 for (std::size_t c : L.callers) {
5331 const double v = njobs(c, idx);
5332 if (std::isfinite(v) && v > 0.0) n += v;
5334 return n < 1.0 ? 1.0 : n;
5350 void update_metrics_ph(
int it) {
5351 const std::size_t N = lqn.nidx;
5352 servt.assign(N + 1, Tzero());
5353 residt.assign(N + 1, Tzero());
5354 callservt.assign(lqn.ncalls + 1, Tzero());
5355 callresidt.assign(lqn.ncalls + 1, Tzero());
5357 std::vector<T> inflNum(N + 1, Tzero()), inflDen(N + 1, Tzero());
5358 std::vector<T> taskTput(N + 1, Tzero()), openTput(N + 1, Tzero());
5361 for (std::size_t hidx = 1; hidx <= lqn.nhosts; ++hidx) {
5362 if (idxhash[hidx] < 0 || !ph_has_layer[hidx])
continue;
5363 const PHLayer& L = ph_layers[hidx];
5364 const LayerResult<T>& res = results.back()[std::size_t(idxhash[hidx])];
5365 const double npop = ph_layer_pop(L, hidx);
5366 for (std::size_t c : L.callers) {
5367 const std::size_t k = L.class_of_caller[c];
5368 T X = Tzero(), Q = Tzero();
5369 for (std::size_t s : L.qstations) {
5370 X = T(X + res.TN(s - 1, k - 1));
5371 Q = T(Q + res.QN(s - 1, k - 1));
5373 const T R = ph_residence(Q, X, res.RN(L.qstations[0] - 1, k - 1));
5374 const T f = ph_inflation_of(R, L.svcmean_by_class[k], npop);
5375 if (!std::isfinite(dbl(X)) || dbl(X) < 0.0) X = Tzero();
5381 taskTput[c] = T(taskTput[c] + num_traits<T>::from_double(lqn.repl[c]) * X);
5382 for (std::size_t eidx : lqn.entriesof[c]) {
5383 const T sh = dbl(ph_share[eidx]) > 0.0 ? ph_share[eidx] : Tzero();
5384 const T w = T(sh * X);
5385 inflNum[eidx] = T(inflNum[eidx] + w * f);
5386 inflDen[eidx] = T(inflDen[eidx] + w);
5389 for (
const auto& oa : L.open_arrivals) {
5390 if (oa.second <= 0)
continue;
5391 const std::size_t k = oa.first;
5392 const std::size_t eidx = std::size_t(oa.second);
5393 T X = Tzero(), Q = Tzero();
5394 for (std::size_t s : L.qstations) {
5395 X = T(X + res.TN(s - 1, k - 1));
5396 Q = T(Q + res.QN(s - 1, k - 1));
5398 if (!std::isfinite(dbl(X)) || dbl(X) <= 0.0)
continue;
5399 const T f = ph_inflation_of(ph_residence(Q, X, res.RN(L.qstations[0] - 1, k - 1)),
5400 L.svcmean_by_class[k], npop);
5401 inflNum[eidx] = T(inflNum[eidx] + X * f);
5402 inflDen[eidx] = T(inflDen[eidx] + X);
5403 openTput[eidx] = T(openTput[eidx] + X);
5404 taskTput[lqn.parent[eidx]] = T(taskTput[lqn.parent[eidx]] + X);
5408 for (std::size_t e = 1; e <= lqn.nentries; ++e) {
5409 const std::size_t eidx = lqn.eshift + e;
5412 f = dbl(inflNum[eidx]) / dbl(inflDen[eidx]);
5413 if (!std::isfinite(f) || f < 1.0) f = 1.0;
5414 const T fv = num_traits<T>::from_double(f);
5415 for (std::size_t aidx : lqn.actsof[eidx])
5416 residt[aidx] = T(fv * (lqn.hostdem[aidx].disabled ? Tzero()
5417 : lqn.hostdem[aidx].mean));
5421 std::vector<T> relw(N + 1, Tzero());
5422 for (std::size_t t = 1; t <= lqn.ntasks; ++t) {
5423 const std::size_t tidx = lqn.tshift + t;
5424 if (idxhash[tidx] < 0 || !ph_has_layer[tidx])
continue;
5425 const PHLayer& L = ph_layers[tidx];
5426 const LayerResult<T>& res = results.back()[std::size_t(idxhash[tidx])];
5427 const double npop = ph_layer_pop(L, tidx);
5428 for (std::size_t c : L.callers) {
5429 const std::size_t k = L.class_of_caller[c];
5430 T X = Tzero(), Q = Tzero();
5431 for (std::size_t s : L.qstations) {
5432 X = T(X + res.TN(s - 1, k - 1));
5433 Q = T(Q + res.QN(s - 1, k - 1));
5435 if (!std::isfinite(dbl(X)) || dbl(X) < 0.0) X = Tzero();
5436 const T g = ph_inflation_of(ph_residence(Q, X, res.RN(L.qstations[0] - 1, k - 1)),
5437 L.svcmean_by_class[k], npop);
5438 for (std::size_t cidx = 1; cidx <= lqn.ncalls; ++cidx) {
5439 if (lqn.calltype[cidx] != CallType::SYNC)
continue;
5440 if (lqn.parent[lqn.callpair_src[cidx]] != c)
continue;
5441 if (lqn.parent[lqn.callpair_dst[cidx]] != tidx)
continue;
5442 const std::size_t eidx = lqn.callpair_dst[cidx];
5443 callservt[cidx] = T(lqn.callproc_mean[cidx] * g * ph_entrymean[eidx]);
5444 callresidt[cidx] = callservt[cidx];
5446 for (std::size_t eidx : lqn.entriesof[tidx])
5447 relw[eidx] = T(relw[eidx] + X * ph_ncalls(c, eidx));
5449 for (
const auto& oa : L.open_arrivals) {
5450 if (oa.second >= 0)
continue;
5451 const std::size_t k = oa.first;
5452 const std::size_t cidx = std::size_t(-oa.second);
5453 const std::size_t eidx = lqn.callpair_dst[cidx];
5455 for (std::size_t s : L.qstations) X = T(X + res.TN(s - 1, k - 1));
5456 const T R = res.RN(L.qstations[0] - 1, k - 1);
5457 if (std::isfinite(dbl(R)) && dbl(R) > 0.0) {
5458 callservt[cidx] = T(R * lqn.callproc_mean[cidx]);
5459 callresidt[cidx] = callservt[cidx];
5461 if (std::isfinite(dbl(X)) && dbl(X) > 0.0) relw[eidx] = T(relw[eidx] + X);
5466 for (std::size_t t = 1; t <= lqn.ntasks; ++t) {
5467 const std::size_t tidx = lqn.tshift + t;
5468 const std::vector<std::size_t>& entries = lqn.entriesof[tidx];
5469 if (entries.empty())
continue;
5476 for (std::size_t eidx : entries) tot = T(tot + relw[eidx] + openTput[eidx]);
5478 for (std::size_t eidx : entries)
5479 ph_share[eidx] = T((relw[eidx] + openTput[eidx]) / tot);
5481 for (std::size_t eidx : entries)
5482 ph_share[eidx] = num_traits<T>::from_double(1.0 /
double(entries.size()));
5490 for (std::size_t eidx : entries) tput[eidx] = T(tput[tidx] * ph_share[eidx]);
5496 const T nrep = num_traits<T>::from_double(std::max(1.0, lqn.repl[tidx]));
5498 : T(tput[tidx] / nrep);
5502 for (std::size_t aidx = lqn.ashift + 1; aidx <= lqn.ashift + lqn.nacts; ++aidx) {
5504 if (!std::isfinite(dbl(v)) && it > 1 && !std::isnan(residt_prev[aidx]))
5505 v = residt_prev_v[aidx];
5506 if (relax_omega < 1.0 && it > 1 && !std::isnan(residt_prev[aidx])) {
5507 const T om = num_traits<T>::from_double(relax_omega);
5508 const T om1 = num_traits<T>::from_double(1.0 - relax_omega);
5509 v = T(om * v + om1 * residt_prev_v[aidx]);
5512 residt_prev[aidx] = dbl(v);
5513 residt_prev_v[aidx] = v;
5515 for (std::size_t cidx = 1; cidx <= lqn.ncalls; ++cidx) {
5516 T v = callservt[cidx];
5517 if (!std::isfinite(dbl(v)))
5518 v = (it > 1 && std::isfinite(callservt_prev[cidx])) ? callservt_prev_v[cidx]
5520 if (relax_omega < 1.0 && it > 1 && !std::isnan(callservt_prev[cidx])) {
5521 const T om = num_traits<T>::from_double(relax_omega);
5522 const T om1 = num_traits<T>::from_double(1.0 - relax_omega);
5523 v = T(om * v + om1 * callservt_prev_v[cidx]);
5525 callservt[cidx] = v;
5526 callresidt[cidx] = v;
5527 callservt_prev[cidx] = dbl(v);
5528 callresidt_prev[cidx] = dbl(v);
5529 callservt_prev_v[cidx] = v;
5537 ph_compose_entry_laws();
5539 for (std::size_t e = 1; e <= lqn.nentries; ++e) {
5540 const std::size_t eidx = lqn.eshift + e;
5541 if (!ph_has_wf[eidx])
continue;
5542 const std::unordered_map<std::size_t, T>& ex = ph_execs[eidx];
5543 for (std::size_t aidx : lqn.actsof[eidx]) {
5544 T sa = T(residt[aidx] + act_think_time(aidx));
5545 for (std::size_t cidx : lqn.callsof[aidx])
5546 if (lqn.calltype[cidx] == CallType::SYNC) sa = T(sa + callservt[cidx]);
5548 servt_prev[aidx] = dbl(sa);
5549 servt_prev_v[aidx] = sa;
5550 tput[aidx] = T(tput[eidx] * ex.at(aidx));
5551 tput_prev[aidx] = dbl(tput[aidx]);
5552 tput_prev_v[aidx] = tput[aidx];
5555 : Distrib<T>::disabled_dist();
5558 servt[eidx] = ph_entrymean[eidx];
5559 residt[eidx] = ph_entrymean[eidx];
5575 void update_think_times_ph(
int it) {
5578 for (std::size_t t = 1; t <= lqn.ntasks; ++t) {
5579 const std::size_t tidx = lqn.tshift + t;
5580 if (ignore[tidx])
continue;
5581 const T ztask = ref_think_time(tidx);
5582 if (idxhash[tidx] < 0) {
5587 const double arvrate = open_arrival_rate_of(tidx);
5589 double nja = lqn.maxmult[tidx];
5590 if (!std::isfinite(nja) || nja <= 0.0) {
5592 for (std::size_t c = 1; c <= NT(); ++c) nja = std::max(nja, njobs(tidx, c));
5595 const std::size_t hidx = lqn.parent[tidx];
5596 if (hidx >= 1 && hidx < idxhash.size() && idxhash[hidx] >= 0 &&
5597 ph_has_layer[hidx] && ph_layers[hidx].class_of_caller[tidx] > 0) {
5598 const PHLayer& HL = ph_layers[hidx];
5599 const LayerResult<T>& hr = results.back()[std::size_t(idxhash[hidx])];
5600 const T rr = hr.RN(HL.qstations[0] - 1, HL.class_of_caller[tidx] - 1);
5601 if (!std::isnan(dbl(rr))) hres = rr;
5603 T za = T(num_traits<T>::from_double(nja / arvrate) - hres - ztask);
5604 if (za < floorv) za = floorv;
5605 if (relax_omega < 1.0 && it > 1 && !std::isnan(thinkt_prev[tidx])) {
5606 const T om = num_traits<T>::from_double(relax_omega);
5607 const T om1 = num_traits<T>::from_double(1.0 - relax_omega);
5608 za = T(om * za + om1 * thinkt_prev_v[tidx]);
5610 tput[tidx] = num_traits<T>::from_double(arvrate);
5612 thinkt_prev[tidx] = dbl(za);
5613 thinkt_prev_v[tidx] = za;
5623 const PHLayer& L = ph_layers[tidx];
5624 const qn::Layer<T>& m = ensemble[std::size_t(idxhash[tidx])];
5625 const LayerResult<T>& r = results.back()[std::size_t(idxhash[tidx])];
5627 for (std::size_t k = 0; k < m.nclasses; ++k) {
5628 const T v = r.UN(L.qstations[0] - 1, k);
5629 if (!std::isnan(dbl(v))) U = T(U + v);
5637 if (dbl(ph_setupshare[tidx]) > 0.0)
5638 U = T(U * (Tone() - ph_setupshare[tidx]));
5641 T X = ph_xdemand[tidx];
5644 double nj = lqn.maxmult[tidx];
5645 if (!std::isfinite(nj) || nj <= 0.0) {
5647 for (std::size_t c = 1; c <= NT(); ++c) nj = std::max(nj, njobs(tidx, c));
5651 if (lqn.sched[tidx] == SchedStrategy::INF) {
5653 z = T((num_traits<T>::from_double(nj) - U) / X - ztask);
5655 const T om = U > Tone() ? T(U - Tone()) : T(Tone() - U);
5656 z = T(num_traits<T>::from_double(nj) * om / X - ztask);
5661 if (z < floorv) z = floorv;
5662 if (it > 1 && !std::isnan(thinkt_prev[tidx]) && !std::isfinite(dbl(z)))
5663 z = thinkt_prev_v[tidx];
5664 if (relax_omega < 1.0 && it > 1 && !std::isnan(thinkt_prev[tidx])) {
5665 const T om = num_traits<T>::from_double(relax_omega);
5666 const T om1 = num_traits<T>::from_double(1.0 - relax_omega);
5667 z = T(om * z + om1 * thinkt_prev_v[tidx]);
5670 thinkt_prev[tidx] = dbl(z);
5671 thinkt_prev_v[tidx] = z;
5685 LnSolution<T> aggregate_ph() {
5686 const std::size_t N = lqn.nidx;
5688 out.QN.assign(N + 1, Tzero());
5689 out.UN.assign(N + 1, Tzero());
5690 out.RN.assign(N + 1, Tzero());
5691 out.TN.assign(N + 1, Tzero());
5692 out.AN.assign(N + 1, Tzero());
5693 out.WN.assign(N + 1, Tzero());
5694 out.defined_Q.assign(N + 1,
false);
5695 out.defined_U.assign(N + 1,
false);
5696 out.defined_R.assign(N + 1,
false);
5697 out.defined_T.assign(N + 1,
false);
5698 out.defined_A.assign(N + 1,
false);
5699 out.defined_W.assign(N + 1,
false);
5700 out.iterations = iterations_done;
5701 out.converged = did_converge;
5703 std::vector<T> PN(N + 1, Tzero()), UT(N + 1, Tzero());
5704 std::vector<bool> hasPN(N + 1,
false), hasUT(N + 1,
false);
5706 for (std::size_t a = 1; a <= lqn.nacts; ++a) {
5707 const std::size_t aidx = lqn.ashift + a;
5708 const std::size_t tidx = lqn.parent[aidx];
5709 if (ignore[tidx])
continue;
5710 const std::size_t hidx = lqn.parent[tidx];
5711 out.TN[aidx] = tput[aidx];
5712 out.defined_T[aidx] =
true;
5713 out.RN[aidx] = servt[aidx];
5714 out.defined_R[aidx] =
true;
5715 UT[aidx] = T(tput[aidx] * servt[aidx]);
5721 const T hd = lqn.hostdem[aidx].disabled ? Tzero() : lqn.hostdem[aidx].mean;
5722 PN[aidx] = T(tput[aidx] * hd / num_traits<T>::from_double(ph_host_servers(hidx)));
5724 PN[hidx] = T(PN[hidx] + PN[aidx]);
5728 for (std::size_t e = 1; e <= lqn.nentries; ++e) {
5729 const std::size_t eidx = lqn.eshift + e;
5730 const std::size_t tidx = lqn.parent[eidx];
5731 if (ignore[tidx])
continue;
5732 out.TN[eidx] = tput[eidx];
5733 out.defined_T[eidx] =
true;
5734 out.RN[eidx] = servt[eidx];
5735 out.defined_R[eidx] =
true;
5736 UT[eidx] = T(tput[eidx] * servt[eidx]);
5738 for (std::size_t aidx : lqn.actsof[eidx]) {
5739 PN[eidx] = T(PN[eidx] + PN[aidx]);
5746 if (ph_has_wf[eidx]) {
5747 const std::unordered_map<std::size_t, T>& ex = ph_execs[eidx];
5748 for (std::size_t aidx : lqn.actsof[eidx]) {
5749 out.WN[aidx] = T(ph_share[eidx] * ex.at(aidx) * residt[aidx]);
5750 out.defined_W[aidx] =
true;
5753 UT[tidx] = T(UT[tidx] + UT[eidx]);
5757 for (std::size_t t = 1; t <= lqn.ntasks; ++t) {
5758 const std::size_t tidx = lqn.tshift + t;
5759 if (ignore[tidx])
continue;
5760 out.TN[tidx] = tput[tidx];
5761 out.defined_T[tidx] =
true;
5764 for (std::size_t aidx : lqn.actsof[tidx]) {
5765 PN[tidx] = T(PN[tidx] + PN[aidx]);
5767 if (out.defined_W[aidx]) { w = T(w + out.WN[aidx]); anyw =
true; }
5769 if (anyw) { out.WN[tidx] = w; out.defined_W[tidx] =
true; }
5772 for (std::size_t hidx = 1; hidx <= lqn.nhosts; ++hidx)
5773 out.defined_T[hidx] =
false;
5775 for (std::size_t idx = 1; idx <= N; ++idx) {
5776 out.QN[idx] = UT[idx];
5777 out.defined_Q[idx] = hasUT[idx];
5778 out.UN[idx] = PN[idx];
5779 out.defined_U[idx] = hasPN[idx];
5785 out.UN[idx] = Tzero();
5786 out.defined_U[idx] =
true;
5787 out.defined_A[idx] =
false;
5788 const bool host = lqn.type[idx] == LqnElement::HOST;
5789 const bool task = lqn.type[idx] == LqnElement::TASK;
5790 const bool entry = lqn.type[idx] == LqnElement::ENTRY;
5791 out.QN[idx] = Tzero();
5792 out.defined_Q[idx] = !host;
5793 out.RN[idx] = Tzero();
5794 out.defined_R[idx] = !host && !task;
5795 out.WN[idx] = Tzero();
5796 out.defined_W[idx] = !host && !entry;
5797 out.TN[idx] = Tzero();
5798 out.defined_T[idx] = !host;
5804 void update_metrics(
int it) {
5805 if (is_ph_encoding()) {
5806 update_metrics_ph(it);
5809 if (lnmethod ==
"moment3") {
5810 update_metrics_moment_based(it);
5813 update_metrics_default(it);
5824 std::size_t layer_refclass(std::size_t idx,
const qn::Layer<T>& L, std::size_t k)
const {
5825 const auto it = layer_normclass.find(idx);
5826 if (it != layer_normclass.end() && k < it->second.size() && it->second[k] > 0)
5827 return it->second[k];
5829 for (std::size_t cc = 0; cc < L.nchains; ++cc)
5830 if (L.chains[cc][k]) c = cc;
5831 return L.refclass[c];
5834 void update_metrics_default(
int it) {
5835 const std::size_t N = lqn.nidx;
5836 servt.assign(N + 1, Tzero());
5837 residt.assign(N + 1, Tzero());
5838 const int iter_min = std::min(30,
int(std::ceil(opt.iter_max / 4.0)));
5839 const bool averaging = averagingstart >= 0 && it >= iter_min;
5840 const int wnd = averaging ? int(it - averagingstart + 1) : 1;
5842 for (
const UpdRow& row : servt_map) {
5843 const std::size_t e = std::size_t(idxhash[row.idx]);
5844 const qn::Layer<T>& L = ensemble[e];
5845 const std::size_t k = row.cls - 1;
5846 const std::size_t refclass_c = layer_refclass(row.idx, L, k);
5847 const std::size_t refstat_k = L.classes[k].refstat;
5849 T sv = Tzero(), rs = Tzero(), tp = Tzero();
5850 const T wT = T(Tone() / num_traits<T>::from_int(wnd));
5851 for (
int w = 0; w < wnd; ++w) {
5852 const LayerResult<T>& r = results[results.size() - 1 - w][e];
5853 sv += r.RN(row.node - 1, k) * wT;
5854 const T TN_ref = (refclass_c > 0 && refstat_k > 0) ? r.TN(refstat_k - 1, refclass_c - 1) : Tzero();
5856 rs += r.QN(row.node - 1, k) / TN_ref * wT;
5858 rs += r.WN(row.node - 1, k) * wT;
5859 tp += r.TN(row.node - 1, k) * wT;
5861 servt[row.aidx] = sv;
5862 residt[row.aidx] = rs;
5863 tput[row.aidx] = tp;
5866 const Distrib<T>& at = lqn.actthink[row.aidx];
5868 servt[row.aidx] = T(servt[row.aidx] + at.mean);
5869 residt[row.aidx] = T(residt[row.aidx] + at.mean);
5877 if (async_only_activity(row.aidx)) residt[row.aidx] = servt[row.aidx];
5879 if (relax_omega < 1.0 && it > 1) {
5880 const T om = num_traits<T>::from_double(relax_omega);
5881 const T om1 = num_traits<T>::from_double(1.0 - relax_omega);
5882 if (!std::isnan(servt_prev[row.aidx]))
5883 servt[row.aidx] = T(om * servt[row.aidx] + om1 * servt_prev_v[row.aidx]);
5884 if (!std::isnan(residt_prev[row.aidx]))
5885 residt[row.aidx] = T(om * residt[row.aidx] + om1 * residt_prev_v[row.aidx]);
5886 if (!std::isnan(tput_prev[row.aidx]))
5887 tput[row.aidx] = T(om * tput[row.aidx] + om1 * tput_prev_v[row.aidx]);
5889 servt_prev[row.aidx] = dbl(servt[row.aidx]);
5890 residt_prev[row.aidx] = dbl(residt[row.aidx]);
5891 tput_prev[row.aidx] = dbl(tput[row.aidx]);
5892 servt_prev_v[row.aidx] = servt[row.aidx];
5893 residt_prev_v[row.aidx] = residt[row.aidx];
5894 tput_prev_v[row.aidx] = tput[row.aidx];
5896 if (servt[row.aidx] > Tzero() && dbl(servt[row.aidx]) <= 1e10)
5906 servt_ph1.assign(N + 1, Tzero());
5907 servt_ph2.assign(N + 1, Tzero());
5908 for (std::size_t a = 1; a <= lqn.nacts; ++a) {
5909 const std::size_t aidx = lqn.ashift + a;
5910 if (lqn.actphase[a] == 1)
5911 servt_ph1[aidx] = servt[aidx];
5913 servt_ph2[aidx] = servt[aidx];
5920 for (
const UpdRow& row : thinkt_map) {
5921 if (!tputproc[row.aidx].disabled)
continue;
5922 const std::size_t e = std::size_t(idxhash[row.idx]);
5924 const T wT = T(Tone() / num_traits<T>::from_int(wnd));
5925 for (
int w = 0; w < wnd; ++w)
5926 tp += results[results.size() - 1 - w][e].TN(row.node - 1, row.cls - 1) * wT;
5927 tput[row.aidx] = tp;
5932 callservt.assign(lqn.ncalls + 1, Tzero());
5933 callresidt.assign(lqn.ncalls + 1, Tzero());
5934 for (
const UpdRow& row : call_map) {
5935 if (row.node <= 1)
continue;
5936 const std::size_t e = std::size_t(idxhash[row.idx]);
5937 const LayerResult<T>& r = results.back()[e];
5938 callservt[row.aidx] =
5939 T(r.RN(row.node - 1, row.cls - 1) * lqn.callproc_mean[row.aidx]);
5945 const qn::Layer<T>& L = ensemble[e];
5946 const std::size_t k = row.cls - 1;
5947 const std::size_t refclass_c = layer_refclass(row.idx, L, k);
5948 const std::size_t refstat_k = L.classes[k].refstat;
5949 const T TN_ref = (refclass_c > 0 && refstat_k > 0) ? r.TN(refstat_k - 1, refclass_c - 1) : Tzero();
5951 ? T(r.QN(row.node - 1, k) / TN_ref)
5952 : r.WN(row.node - 1, k);
5954 const T rw = region_wait(e, row.cls);
5956 callservt[row.aidx] = T(callservt[row.aidx] + rw);
5957 callresidt[row.aidx] = T(callresidt[row.aidx] + rw);
5959 if (relax_omega < 1.0 && it > 1 && !std::isnan(callservt_prev[row.aidx])) {
5960 const T om = num_traits<T>::from_double(relax_omega);
5961 const T om1 = num_traits<T>::from_double(1.0 - relax_omega);
5962 callservt[row.aidx] =
5963 T(om * callservt[row.aidx] + om1 * callservt_prev_v[row.aidx]);
5965 callservt_prev[row.aidx] = dbl(callservt[row.aidx]);
5966 callservt_prev_v[row.aidx] = callservt[row.aidx];
5967 callresidt_prev[row.aidx] = dbl(callresidt[row.aidx]);
5970 resolve_entry_service();
5982 for (std::size_t e = 1; e <= lqn.nentries; ++e) {
5983 const std::size_t eidx = lqn.eshift + e;
5984 servt_ph2[eidx] = (eidx < entry_servt_last.size() && entry_servt_last[eidx] > Tzero())
5985 ? T(entry_servt_ph2[eidx] * servt[eidx] / entry_servt_last[eidx])
5987 servt_ph1[eidx] = T(servt[eidx] - servt_ph2[eidx]);
5989 for (std::size_t e = 1; e <= lqn.nentries; ++e) {
5990 const std::size_t eidx = lqn.eshift + e;
5991 const std::size_t tidx = lqn.parent[eidx];
5993 if (lqn.isref[tidx] || !lqn.issynccaller.any_col(eidx)) {
5994 residt[eidx] = servt[eidx];
5997 T entry_tput = Tzero();
5999 entry_tput = tput[eidx];
6001 entry_tput = tput[tidx];
6007 residt[eidx] = T(servt_ph1[eidx] + prOvertake[e] * servt_ph2[eidx]);
6011 for (
const UpdRow& row : call_map) {
6012 if (row.node <= 1)
continue;
6013 const std::size_t eidx = lqn.callpair_dst[row.aidx];
6016 for (
const UpdRow& row : call_map) {
6017 if (row.node <= 1)
continue;
6018 const std::size_t eidx = lqn.callpair_dst[row.aidx];
6020 callservt[row.aidx] = servt[eidx];
6021 callservtproc[row.aidx] = servtproc[eidx];
6022 }
else if (callservt[row.aidx] > Tzero()) {
6032 for (std::size_t cidx = 1; cidx <= lqn.ncalls; ++cidx) {
6033 if (lqn.calltype[cidx] != CallType::SYNC)
continue;
6034 const std::size_t target = lqn.callpair_dst[cidx];
6035 if (target <= lqn.eshift || target > lqn.eshift + lqn.nentries)
continue;
6037 const T eff = residt[target];
6038 if (!(eff > Tzero()))
continue;
6039 const T w = T(eff * lqn.callproc_mean[cidx]);
6040 callservt[cidx] = w;
6041 callresidt[cidx] = w;
6061 static void cdf_moments(
const fluid::FluidPassage& c,
double& m1,
double& m2,
double& m3) {
6063 for (std::size_t i = 0; i + 1 < c.t.size(); ++i) {
6064 const double x = 0.5 * (c.t[i + 1] + c.t[i]);
6065 const double w = c.cdf[i + 1] - c.cdf[i];
6068 m3 += x * x * x * w;
6073 const std::vector<std::vector<fluid::FluidPassage>>& layer_cdf(std::size_t e) {
6074 if (cdf_repo[e].empty()) {
6081 ensemble[e].refresh_rt();
6083 cdf_repo[e] = detail::ln_fluid_cdf_respt(ensemble[e], opt.layer_fluid);
6084 }
catch (
const std::exception&) {
6091 cdf_repo[e] = layer_cdf_exp_fallback(e);
6099 std::vector<std::vector<fluid::FluidPassage>> layer_cdf_exp_fallback(std::size_t e) {
6100 const qn::NetworkStruct<T>& L = ensemble[e];
6101 std::vector<std::vector<fluid::FluidPassage>> out(
6102 L.nstations, std::vector<fluid::FluidPassage>(L.nclasses));
6103 if (results.empty() || e >= results.back().size())
return out;
6104 const Matrix<T>& RN = results.back()[e].RN;
6105 const std::size_t npts = 100;
6106 for (std::size_t i = 0; i < L.nstations; ++i) {
6107 if (L.stations[i].nodetype == qn::NodeType::Source)
continue;
6108 for (std::size_t r = 0; r < L.nclasses; ++r) {
6109 if (L.disabled[i][r])
continue;
6110 const double rn = (i < static_cast<std::size_t>(RN.rows()) &&
6111 r < static_cast<std::size_t>(RN.cols()))
6114 if (!(std::isfinite(rn) && rn > 0.0))
continue;
6115 fluid::FluidPassage& cell = out[i][r];
6116 cell.t.reserve(npts);
6117 cell.cdf.reserve(npts);
6118 for (std::size_t j = 0; j < npts; ++j) {
6120 0.001 + (0.999 - 0.001) *
static_cast<double>(j) / (npts - 1);
6121 cell.cdf.push_back(q);
6122 cell.t.push_back(-std::log(1.0 - q) * rn);
6143 int entry_visit_ratio(std::size_t eidx, T& ratio)
const {
6144 const std::size_t tidx = lqn.parent[eidx];
6145 const std::size_t hidx = lqn.parent[tidx];
6146 if (ignore[tidx] || ignore[hidx] || idxhash[hidx] < 0)
return 0;
6147 if (!lqn.issynccaller.any_col(eidx))
return 1;
6148 const std::size_t hl = std::size_t(idxhash[hidx]);
6149 const qn::Layer<T>& L = ensemble[hl];
6150 const LayerResult<T>& r = results.back()[hl];
6151 T task_tput = Tzero(), entry_tput = Tzero();
6152 for (
const auto& kv : L.attr_tasks)
6153 if (kv.second == tidx) task_tput += r.TN(L.clientIdx - 1, kv.first - 1);
6154 for (
const auto& kv : L.attr_entries)
6155 if (kv.second == eidx) entry_tput += r.TN(L.clientIdx - 1, kv.first - 1);
6157 ratio = T(task_tput / (entry_tput > floor ? entry_tput : floor));
6162 Matrix<T> entry_service_resolvent()
const {
6163 const std::size_t dim = lqn.nidx + lqn.ncalls;
6164 Matrix<T> A(dim + 1, dim + 1, Tzero());
6165 for (std::size_t i = 0; i <= dim; ++i) {
6167 for (std::size_t j = 0; j <= dim; ++j) A(i, j) = T(A(i, j) - servtmatrix(i, j));
6169 return ::line::inverse(A);
6195 void update_metrics_moment_based(
int it) {
6196 const std::size_t N = lqn.nidx;
6197 servt.assign(N + 1, Tzero());
6198 residt.assign(N + 1, Tzero());
6199 callservt.assign(lqn.ncalls + 1, Tzero());
6200 callresidt.assign(lqn.ncalls + 1, Tzero());
6203 for (
const UpdRow& row : servt_map) {
6204 const std::size_t e = std::size_t(idxhash[row.idx]);
6205 const qn::Layer<T>& L = ensemble[e];
6206 const LayerResult<T>& r = results.back()[e];
6207 const std::size_t k = row.cls - 1;
6209 for (std::size_t cc = 0; cc < L.nchains; ++cc)
6210 if (L.chains[cc][k]) c = cc;
6211 const std::size_t refclass_c = L.refclass[c];
6212 const std::size_t refstat_k = L.classes[k].refstat;
6213 const T TN_ref = (refclass_c > 0 && refstat_k > 0) ? r.TN(refstat_k - 1, refclass_c - 1) : Tzero();
6215 tput[row.aidx] = r.TN(row.node - 1, k);
6217 ? T(r.QN(row.node - 1, k) / TN_ref)
6218 : r.WN(row.node - 1, k);
6219 if (!hasconverged) {
6220 servt[row.aidx] = r.RN(row.node - 1, k);
6222 const Distrib<T>& at = lqn.actthink[row.aidx];
6224 servt[row.aidx] = T(servt[row.aidx] + at.mean);
6225 residt[row.aidx] = T(residt[row.aidx] + at.mean);
6230 if (async_only_activity(row.aidx)) residt[row.aidx] = servt[row.aidx];
6232 servtcdf[row.aidx] = layer_cdf(e)[row.node - 1][k];
6236 for (
const UpdRow& row : call_map) {
6237 if (row.node <= 1)
continue;
6238 const std::size_t e = std::size_t(idxhash[row.idx]);
6239 const LayerResult<T>& r = results.back()[e];
6240 callresidt[row.aidx] = r.WN(row.node - 1, row.cls - 1);
6242 callservt[row.aidx] =
6243 T(r.RN(row.node - 1, row.cls - 1) * lqn.callproc_mean[row.aidx]);
6245 callservtcdf[row.aidx] = layer_cdf(e)[row.node - 1][row.cls - 1];
6249 moment3_means_pass(it);
6251 moment3_distribution_pass(it);
6255 void moment3_means_pass(
int it) {
6256 const std::size_t dim = lqn.nidx + lqn.ncalls;
6257 std::vector<T> x(dim + 1, Tzero());
6258 for (std::size_t i = 1; i <= lqn.nidx; ++i) x[i] = residt[i];
6259 for (std::size_t c = 1; c <= lqn.ncalls; ++c) x[lqn.nidx + c] = callresidt[c];
6262 const Matrix<T> Rinv = entry_service_resolvent();
6263 std::vector<T> entry_servt(dim + 1, Tzero());
6264 for (std::size_t i = 1; i <= dim; ++i) {
6266 for (std::size_t j = 1; j <= dim; ++j)
6267 if (Rinv(i, j) != Tzero()) s += Rinv(i, j) * x[j];
6270 for (std::size_t i = 1; i <= lqn.eshift; ++i) entry_servt[i] = Tzero();
6281 for (std::size_t e = 1; e <= lqn.nentries; ++e)
6282 servt[lqn.eshift + e] = entry_servt[lqn.eshift + e];
6285 std::vector<T> entry_residt(dim + 1, Tzero());
6286 for (std::size_t i = lqn.eshift + 1; i <= lqn.eshift + lqn.nentries; ++i) {
6288 for (std::size_t j = 1; j <= dim; ++j)
6289 if (servtmatrix(i, j) != Tzero()) s += servtmatrix(i, j) * x[j];
6290 entry_residt[i] = s;
6293 for (std::size_t e = 1; e <= lqn.nentries; ++e) {
6294 const std::size_t eidx = lqn.eshift + e;
6296 const int state = entry_visit_ratio(eidx, ratio);
6297 if (state == 0)
continue;
6299 residt[eidx] = entry_residt[eidx];
6302 servt[eidx] = T(entry_servt[eidx] * ratio);
6303 residt[eidx] = T(entry_residt[eidx] * ratio);
6306 for (
const UpdRow& row : call_map) {
6307 if (row.node <= 1)
continue;
6308 const std::size_t eidx = lqn.callpair_dst[row.aidx];
6311 for (
const UpdRow& row : call_map) {
6312 if (row.node <= 1)
continue;
6313 const std::size_t eidx = lqn.callpair_dst[row.aidx];
6316 callservt[row.aidx] = servt[eidx];
6317 callservtproc[row.aidx] = servtproc[eidx];
6318 }
else if (callservt[row.aidx] > Tzero()) {
6334 void moment3_distribution_pass(
int it) {
6335 if constexpr (!num_traits<T>::has_transcendental) {
6337 "SolverLN: the 'moment3' method fits an APH to a response-time CDF, which needs "
6338 "square roots and a matrix exponential; rerun with --arith double or real");
6340 const std::size_t dim = lqn.nidx + lqn.ncalls;
6341 const Matrix<T> Rinv = entry_service_resolvent();
6344 mam::AphPair<T> zero_law;
6345 zero_law.alpha.push_back(Tone());
6348 for (std::size_t en = 1; en <= lqn.nentries; ++en) {
6349 const std::size_t eidx = lqn.eshift + en;
6350 std::vector<mam::AphPair<T>> seq;
6351 for (std::size_t fitidx = 1; fitidx <= dim; ++fitidx) {
6352 if (!(Rinv(eidx, fitidx) > Tzero()))
continue;
6356 if (fitidx <= lqn.eshift + lqn.nentries)
continue;
6357 const bool is_call = fitidx > lqn.nidx;
6358 const fluid::FluidPassage& curve =
6359 is_call ? callservtcdf[fitidx - lqn.nidx] : servtcdf[fitidx];
6360 double m1 = 0.0, m2 = 0.0, m3 = 0.0;
6361 cdf_moments(curve, m1, m2, m3);
6365 if (!is_call && !lqn.actthink[fitidx].disabled &&
6367 const Distrib<T>& zt = lqn.actthink[fitidx];
6371 m3 = m3 + 3.0 * m2 * t1 + 3.0 * m1 * t2 + t3;
6372 m2 = m2 + 2.0 * m1 * t1 + t2;
6383 num_traits<T>::from_double(m1), num_traits<T>::from_double(m2),
6384 num_traits<T>::from_double(m3), 10u,
6386 const mam::AphPair<T> law = aph_pair_of_map(fit.aph);
6388 double reps = dbl(Rinv(eidx, fitidx));
6391 if (is_call) reps *= dbl(lqn.callproc_mean[fitidx - lqn.nidx]);
6392 const long whole =
static_cast<long>(std::floor(reps));
6393 const double frac = reps -
static_cast<double>(whole);
6394 for (
long q = 0; q < whole; ++q) seq.push_back(law);
6397 num_traits<T>::from_double(frac),
6398 num_traits<T>::from_double(1.0 - frac),
6402 const std::size_t cidx = fitidx - lqn.nidx;
6403 callservt[cidx] = num_traits<T>::from_double(m1);
6406 servt[fitidx] = num_traits<T>::from_double(m1);
6412 servt[eidx] = Tzero();
6416 entryproc[en] = entry_law;
6419 entrycdfrespt[en] = aph_eval_cdf(entry_law);
6427 std::vector<T> x(dim + 1, Tzero());
6428 for (std::size_t i = 1; i <= lqn.nidx; ++i) x[i] = residt[i];
6429 for (std::size_t c = 1; c <= lqn.ncalls; ++c) x[lqn.nidx + c] = callresidt[c];
6430 for (std::size_t en = 1; en <= lqn.nentries; ++en) {
6431 const std::size_t eidx = lqn.eshift + en;
6433 for (std::size_t j = 1; j <= dim; ++j)
6434 if (servtmatrix(eidx, j) != Tzero()) s += servtmatrix(eidx, j) * x[j];
6436 const int state = entry_visit_ratio(eidx, ratio);
6437 if (state == 0)
continue;
6438 residt[eidx] = state == 2 ? T(s * ratio) : s;
6442 for (
const UpdRow& row : call_map) {
6443 if (row.node <= 1)
continue;
6444 const std::size_t eidx = lqn.callpair_dst[row.aidx];
6445 callservt[row.aidx] = servt[eidx];
6452 moment_pass_done =
true;
6464 static mam::AphPair<T> aph_pair_of_map(
const mam::Map<T>& m) {
6465 mam::AphPair<T> out;
6466 const std::size_t n = m.D0.rows();
6468 out.alpha.assign(n, Tzero());
6469 for (std::size_t i = 0; i < n; ++i) {
6471 for (std::size_t j = 0; j < n; ++j) rowsum += m.D1(i, j);
6473 for (std::size_t j = 0; j < n; ++j) out.alpha[j] = T(m.D1(i, j) / rowsum);
6478 out.alpha[0] = Tone();
6486 static LnCdf aph_eval_cdf(
const mam::AphPair<T>& law) {
6490 const double var = m2 - m1 * m1;
6491 const double sigma = var > 0.0 ? std::sqrt(var) : 0.0;
6492 const double tmax = m1 + 10.0 * sigma;
6493 const std::size_t P = 500;
6496 for (std::size_t k = 0; k < P; ++k) {
6497 const double t = tmax *
static_cast<double>(k) /
static_cast<double>(P - 1);
6499 Matrix<T> St(law.S.rows(), law.S.cols(), Tzero());
6500 for (std::size_t i = 0; i < law.S.rows(); ++i)
6501 for (std::size_t j = 0; j < law.S.cols(); ++j)
6502 St(i, j) = T(law.S(i, j) * num_traits<T>::from_double(t));
6505 for (std::size_t i = 0; i < E.rows(); ++i)
6506 for (std::size_t j = 0; j < E.cols(); ++j) surv += law.alpha[i] * E(i, j);
6507 out.cdf[k] = 1.0 - dbl(surv);
6526 std::vector<T> join_excess()
const {
6527 std::vector<T> excess(lqn.nidx + 1, Tzero());
6528 std::vector<std::size_t> joined;
6529 for (std::size_t tail = 1; tail <= lqn.nidx; ++tail) {
6530 if (lqn.actpretype[tail] != PrecedenceType::PRE_AND)
continue;
6531 for (std::size_t sx : lqn.graph.succ(tail))
6532 if (sx > lqn.ashift && sx <= lqn.ashift + lqn.nacts) joined.push_back(sx);
6534 std::sort(joined.begin(), joined.end());
6535 joined.erase(std::unique(joined.begin(), joined.end()), joined.end());
6536 if (joined.empty())
return excess;
6538 fj::LqnBranchView<T> view;
6539 view.graph =
Matrix<T>(lqn.nidx, lqn.nidx, Tzero());
6540 for (std::size_t i = 1; i <= lqn.nidx; ++i)
6541 for (
const auto& e : lqn.graph.row[i]) view.graph(i - 1, e.first - 1) = e.second;
6542 view.ashift = lqn.ashift;
6543 view.nacts = lqn.nacts;
6544 view.actposttype.assign(lqn.nidx, 0);
6545 for (std::size_t i = 1; i <= lqn.nidx; ++i)
6546 view.actposttype[i - 1] =
static_cast<int>(lqn.actposttype[i]);
6548 for (std::size_t aidx : joined) {
6549 const std::vector<std::vector<std::size_t>> branches =
6551 if (branches.empty())
continue;
6552 std::vector<T> means;
6553 for (
const auto& br : branches) {
6555 for (std::size_t a : br) {
6559 if (a < lqn.callsof.size())
6560 for (std::size_t cidx : lqn.callsof[a])
6561 if (cidx < lqn.calltype.size() && lqn.calltype[cidx] == CallType::SYNC
6562 && cidx < callresidt.size())
6563 s += callresidt[cidx];
6567 if (means.size() == 1)
continue;
6568 std::size_t quorum = means.size();
6569 if (lqn.actquorum[aidx] >= 1 && lqn.actquorum[aidx] <= means.size())
6570 quorum = lqn.actquorum[aidx];
6571 std::vector<T> vars;
6572 for (
const T& m2 : means) vars.push_back(T(m2 * m2));
6574 for (
const T& m2 : means) serial += m2;
6575 if constexpr (num_traits<T>::has_transcendental) {
6577 excess[aidx] = T(q.m - serial);
6580 "SolverLN: the AND-join completion time is a k-th order statistic fitted "
6581 "through a three-point distribution, which needs a square root; use the "
6582 "double or real backend for a model with an AND join");
6589 void resolve_entry_service() {
6590 const std::size_t dim = lqn.nidx + lqn.ncalls;
6591 std::vector<T> x(dim + 1, Tzero());
6592 for (std::size_t i = 1; i <= lqn.nidx; ++i) x[i] = residt[i];
6593 for (std::size_t c = 1; c <= lqn.ncalls; ++c) x[lqn.nidx + c] = callresidt[c];
6594 std::vector<T> entry_servt(dim + 1, Tzero());
6595 for (std::size_t i = 1; i <= dim; ++i) {
6597 for (std::size_t j = 1; j <= dim; ++j)
6598 if (servtmatrix(i, j) != Tzero()) s += servtmatrix(i, j) * x[j];
6604 const std::vector<T> excess = join_excess();
6605 for (std::size_t e = 1; e <= lqn.nentries; ++e) {
6606 const std::size_t eidx = lqn.eshift + e;
6607 for (std::size_t aidx = 1; aidx <= lqn.nidx; ++aidx) {
6608 if (excess[aidx] == Tzero())
continue;
6609 if (servtmatrix(eidx, aidx) > Tzero())
6610 entry_servt[eidx] = T(entry_servt[eidx] + excess[aidx]);
6612 if (entry_servt[eidx] < Tzero()) entry_servt[eidx] = Tzero();
6622 for (std::size_t e = 1; e <= lqn.nentries; ++e) {
6623 const std::size_t eidx = lqn.eshift + e;
6624 const double c = setup_charge(lqn.parent[eidx]);
6626 entry_servt[eidx] = T(entry_servt[eidx] + num_traits<T>::from_double(c));
6632 std::vector<T> x2(dim + 1, Tzero());
6633 for (std::size_t a = 1; a <= lqn.nacts; ++a)
6634 if (lqn.actphase[a] > 1) x2[lqn.ashift + a] = x[lqn.ashift + a];
6635 for (std::size_t c = 1; c <= lqn.ncalls; ++c) {
6636 const std::size_t src = lqn.callpair_src[c];
6637 if (src > lqn.ashift && src <= lqn.ashift + lqn.nacts && lqn.actphase[src - lqn.ashift] > 1)
6638 x2[lqn.nidx + c] = x[lqn.nidx + c];
6640 entry_servt_ph2.assign(dim + 1, Tzero());
6641 for (std::size_t i = lqn.eshift + 1; i <= dim; ++i) {
6643 for (std::size_t j = 1; j <= dim; ++j)
6644 if (servtmatrix(i, j) != Tzero()) s2 += servtmatrix(i, j) * x2[j];
6645 entry_servt_ph2[i] = s2;
6647 entry_servt_last = entry_servt;
6652 for (std::size_t e = 1; e <= lqn.nentries; ++e) {
6653 const std::size_t eidx = lqn.eshift + e;
6654 const std::size_t tidx = lqn.parent[eidx];
6655 const std::size_t hidx = lqn.parent[tidx];
6656 if (ignore[tidx] || ignore[hidx])
continue;
6657 if (idxhash[hidx] < 0)
continue;
6658 const bool has_sync = lqn.issynccaller.any_col(eidx);
6660 servt[eidx] = entry_servt[eidx];
6661 residt[eidx] = entry_servt[eidx];
6664 const std::size_t hl = std::size_t(idxhash[hidx]);
6665 const qn::Layer<T>& L = ensemble[hl];
6666 const LayerResult<T>& r = results.back()[hl];
6667 T task_tput = Tzero(), entry_tput = Tzero();
6668 for (
const auto& kv : L.attr_tasks)
6669 if (kv.second == tidx) task_tput += r.TN(L.clientIdx - 1, kv.first - 1);
6670 for (
const auto& kv : L.attr_entries)
6671 if (kv.second == eidx) entry_tput += r.TN(L.clientIdx - 1, kv.first - 1);
6673 servt[eidx] = T(entry_servt[eidx] * task_tput / entry_tput);
6674 residt[eidx] = servt[eidx];
6678 entryvisits[eidx] = T(entry_tput / task_tput);
6680 servt[eidx] = entry_servt[eidx];
6681 residt[eidx] = entry_servt[eidx];
6697 bool layer_takes_interlock(std::size_t e)
const {
6698 return opt.layer_solver ==
"mva" && e < ensemble.size() &&
6714 std::vector<std::vector<double>> build_layer_interlock(std::size_t e,
6715 const std::vector<std::size_t>& hts,
6716 const std::vector<T>& prIL,
6717 const std::vector<T>& PrIL)
const {
6718 const qn::Layer<T>& L = ensemble[e];
6719 const std::size_t R = L.nclasses;
6720 std::vector<double> cls_ir(R, 0.0);
6721 std::vector<double> cls_pr(R, 0.0);
6722 auto stamp = [&](std::size_t classIdx, std::size_t tidx) {
6723 if (classIdx < 1 || classIdx > R)
return;
6724 for (std::size_t i = 0; i < hts.size(); ++i)
6725 if (hts[i] == tidx) {
6726 cls_ir[classIdx - 1] = dbl(prIL[i]);
6727 cls_pr[classIdx - 1] = dbl(PrIL[i]);
6730 for (
const auto& a : L.attr_tasks) stamp(a.first, a.second);
6731 for (
const auto& a : L.attr_entries) stamp(a.first, lqn.parent[a.second]);
6732 for (
const auto& a : L.attr_activities) stamp(a.first, lqn.parent[a.second]);
6733 for (
const auto& a : L.attr_calls) stamp(a[0], lqn.parent[a[2]]);
6735 std::vector<std::vector<double>> IL(R, std::vector<double>(R, 0.0));
6737 for (std::size_t r = 0; r < R; ++r) {
6739 for (std::size_t sIl = 0; sIl < R; ++sIl) {
6741 IL[r][sIl] = cls_pr[sIl] * cls_ir[sIl] * cls_ir[r];
6745 return any ? IL : std::vector<std::vector<double>>();
6758 std::pair<T, T> interlock_prob(std::size_t client_tidx, std::size_t server_idx,
6759 bool isProcessorHost)
const {
6760 const std::vector<std::size_t>& common = il_common_entries[server_idx];
6761 const double nsrc = il_num_sources[server_idx];
6762 if (nsrc == 0.0 || common.empty())
return std::make_pair(Tzero(), Tzero());
6763 const std::vector<std::size_t>& allsrc = il_src_all[server_idx];
6764 T sum_flow = Tzero();
6765 T sum_pril = Tzero();
6766 for (std::size_t ce : common) {
6767 const std::size_t srcTask = lqn.parent[ce];
6768 const std::size_t cen = ce - lqn.eshift;
6770 double m_src = lqn.mult[srcTask];
6771 if (!std::isfinite(m_src) || m_src < 1.0) m_src = 1.0;
6772 double m_eff = m_src;
6773 if (isProcessorHost) m_eff = (m_src > 3.0) ? (m_src + m_src) : (m_src * m_src);
6774 for (std::size_t dst : lqn.entriesof[client_tidx]) {
6775 const std::size_t dn = dst - lqn.eshift;
6776 if (dn < 1 || dn > lqn.nentries)
continue;
6777 if (!(il_all(cen, dn) > Tzero()))
continue;
6778 T ce_tput = tput[ce];
6780 ce_tput = tput[lqn.actsof[ce][0]];
6784 if (std::find(allsrc.begin(), allsrc.end(), srcTask) != allsrc.end()) {
6785 const T contrib = T(ce_tput * il_all(cen, dn));
6786 sum_flow += contrib;
6787 sum_pril += T(contrib / num_traits<T>::from_double(m_eff));
6791 T client_tput = tput[client_tidx];
6793 for (std::size_t e : lqn.entriesof[client_tidx]) {
6796 et = tput[lqn.actsof[e][0]];
6797 client_tput = T(client_tput + et);
6801 return std::make_pair(Tzero(), Tzero());
6802 const T flow = sum_flow < client_tput ? sum_flow : client_tput;
6803 T IR = T(flow / client_tput);
6804 if (IR > Tone()) IR = Tone();
6805 if (dbl(IR) < 0.0) IR = Tzero();
6808 if (pr > Tone()) pr = Tone();
6809 if (dbl(pr) < 0.0) pr = Tzero();
6810 return std::make_pair(IR, pr);
6813 void update_populations(
int) {
6814 const std::vector<T> residt_orig = residt;
6815 const std::vector<T> callresidt_orig = callresidt;
6816 bool adjusted =
false;
6818 for (std::size_t cidx = 1; cidx <= lqn.ncalls; ++cidx) {
6819 if (lqn.calltype[cidx] != CallType::SYNC)
continue;
6820 const std::size_t dst = lqn.callpair_dst[cidx];
6821 const std::size_t server_tidx = lqn.parent[dst];
6822 std::size_t server = 0;
6823 if (server_tidx <= NT() && !il_common_entries[server_tidx].empty()) {
6824 server = server_tidx;
6825 }
else if (server_tidx > lqn.tshift) {
6826 const std::size_t h = lqn.parent[server_tidx];
6827 if (h >= 1 && h <= NT() && !il_common_entries[h].empty()) server = h;
6829 if (server == 0)
continue;
6830 const std::size_t client_tidx = lqn.parent[lqn.callpair_src[cidx]];
6834 const std::pair<T, T> ilp = interlock_prob(client_tidx, server,
false);
6835 const T prIL = T(ilp.first * ilp.second);
6837 const T S = servt[dst];
6838 const T cm = lqn.callproc_mean[cidx];
6839 if (!(cm > Tzero()) || !(callservt[cidx] > Tzero()))
continue;
6840 const T RN = T(callservt[cidx] / cm);
6841 const T W = RN > S ? T(RN - S) : Tzero();
6843 const T RN_adj = T(S + (Tone() - prIL) * W);
6844 const T scale = T(RN_adj / RN);
6845 callservt[cidx] = T(callservt[cidx] * scale);
6846 callresidt[cidx] = T(callresidt[cidx] * scale);
6847 if (callservt[cidx] > Tzero())
6854 layer_interlock.assign(ensemble.size(), {});
6855 for (std::size_t h = 1; h <= lqn.nhosts; ++h) {
6856 if (il_common_entries[h].empty())
continue;
6857 const std::vector<std::size_t>& hts = lqn.tasksof[h];
6858 std::vector<T> prIL(hts.size(), Tzero()), PrIL(hts.size(), Tzero()),
6859 tutil(hts.size(), Tzero());
6860 for (std::size_t i = 0; i < hts.size(); ++i) {
6862 const std::pair<T, T> ilp = interlock_prob(hts[i], h,
true);
6863 prIL[i] = ilp.first;
6864 PrIL[i] = ilp.second;
6865 for (std::size_t e : lqn.entriesof[hts[i]])
6866 for (std::size_t a : lqn.actsof[e]) tutil[i] += tput[a] * lqn.hostdem[a].mean;
6868 T Utot = Tzero(), Uil = Tzero();
6869 for (std::size_t i = 0; i < hts.size(); ++i) {
6875 const T frac = T(Uil / Utot);
6883 if (idxhash[h] >= 0 && layer_takes_interlock(std::size_t(idxhash[h]))) {
6884 const std::size_t e = std::size_t(idxhash[h]);
6885 layer_interlock[e] = build_layer_interlock(e, hts, prIL, PrIL);
6889 for (std::size_t i = 0; i < hts.size(); ++i) {
6897 const T eff = T(prIL[i] * PrIL[i] * frac);
6898 for (std::size_t e : lqn.entriesof[hts[i]])
6899 for (std::size_t a : lqn.actsof[e]) {
6900 const T D = lqn.hostdem[a].mean;
6902 residt[a] = T(D + (Tone() - eff) * (residt[a] - D));
6909 if (!adjusted)
return;
6910 const std::size_t dim = lqn.nidx + lqn.ncalls;
6911 auto esum = [&](
const std::vector<T>& rs,
const std::vector<T>& cr, std::size_t i) {
6913 for (std::size_t j = 1; j <= dim; ++j) {
6914 if (servtmatrix(i, j) == Tzero())
continue;
6915 s += servtmatrix(i, j) * (j <= lqn.nidx ? rs[j] : cr[j - lqn.nidx]);
6919 for (std::size_t e = 1; e <= lqn.nentries; ++e) {
6920 const std::size_t eidx = lqn.eshift + e;
6921 const T oldv = esum(residt_orig, callresidt_orig, eidx);
6923 const T ratio = T(esum(residt, callresidt, eidx) / oldv);
6932 if (!moment_pass_done) {
6933 servt[eidx] = T(servt[eidx] * ratio);
6936 residt[eidx] = T(residt[eidx] * ratio);
6952 T ref_think_time(std::size_t tidx)
const {
6953 if (tidx >= lqn.isref.size() || !lqn.isref[tidx])
return Tzero();
6954 if (lqn.think[tidx].disabled)
return Tzero();
6955 return lqn.think[tidx].mean;
6964 T phase2_tail(std::size_t tidx)
const {
6966 if (!has_phase2)
return tail;
6969 for (std::size_t eidx : lqn.entriesof[tidx]) {
6970 if (servt_ph2[eidx] <= fine)
continue;
6972 if (w <= fine) w = Tone();
6973 const T d = T(servt[eidx] - residt[eidx]);
6974 tail += T(w * (d > Tzero() ? d : Tzero()));
6977 if (wtot > fine) tail = T(tail / wtot);
6981 void update_think_times(
int it) {
6984 if (is_ph_encoding()) {
6985 update_think_times_ph(it);
6989 for (std::size_t t = 1; t <= lqn.ntasks; ++t) {
6990 const std::size_t tidx = lqn.tshift + t;
6992 const T ztask = ref_think_time(tidx);
6993 if (idxhash[tidx] < 0) {
6997 const double arvrate = open_arrival_rate_of(tidx);
7000 for (std::size_t c = 1; c <= NT(); ++c) nja = std::max(nja, njobs(tidx, c));
7001 if (!(nja > 0.0)) nja = lqn.maxmult[tidx];
7003 for (std::size_t eidx : lqn.entriesof[tidx])
7004 for (std::size_t aidx : lqn.actsof[eidx])
7005 if (!std::isnan(dbl(residt[aidx]))) hres += residt[aidx];
7007 T za = T(num_traits<T>::from_double(nja / arvrate) - hres - ztask);
7008 if (za < floora) za = floora;
7009 if (relax_omega < 1.0 && it > 1 && !std::isnan(thinkt_prev[tidx])) {
7010 const T om = num_traits<T>::from_double(relax_omega);
7011 const T om1 = num_traits<T>::from_double(1.0 - relax_omega);
7012 za = T(om * za + om1 * thinkt_prev_v[tidx]);
7014 tput[tidx] = num_traits<T>::from_double(arvrate);
7016 thinkt_prev[tidx] = dbl(za);
7017 thinkt_prev_v[tidx] = za;
7025 const std::size_t e = std::size_t(idxhash[tidx]);
7026 const qn::Layer<T>& L = ensemble[e];
7027 const LayerResult<T>& r = results.back()[e];
7029 for (std::size_t c = 1; c <= NT(); ++c) nj = std::max(nj, njobs(tidx, c));
7031 const std::size_t ts = station_idx_of(L, tidx);
7032 T tp = Tzero(), ut = Tzero();
7033 for (std::size_t k = 0; k < L.nclasses; ++k) {
7034 tp += r.TN(ts - 1, k);
7035 ut += r.UN(ts - 1, k);
7037 tput[tidx] = T(num_traits<T>::from_double(lqn.repl[tidx]) * tp);
7042 const double nrep = std::max(1.0, lqn.repl[tidx]);
7043 const T xcyc = T(tput[tidx] / num_traits<T>::from_double(nrep));
7044 if (nrep > 1.0 && std::isfinite(lqn.maxmult[tidx]) && lqn.maxmult[tidx] > 0.0) nj = lqn.maxmult[tidx];
7046 if (lqn.sched[tidx] == SchedStrategy::INF) {
7048 raw = xcyc == Tzero()
7050 : T((num_traits<T>::from_double(nj) - util[tidx]) / xcyc - ztask);
7052 const T om = util[tidx] > Tone() ? T(util[tidx] - Tone()) : T(Tone() - util[tidx]);
7053 raw = xcyc == Tzero()
7055 : T(num_traits<T>::from_double(nj) * om / xcyc - ztask);
7058 thinkt[tidx] = raw < floorv ? floorv : raw;
7062 if (lnmethod !=
"flat.cs") {
7063 const T lessz = T(thinkt[tidx] - phase2_tail(tidx));
7064 thinkt[tidx] = lessz < floorv ? floorv : lessz;
7076 const double c = setup_charge(tidx);
7078 const T withc = T(thinkt[tidx] + num_traits<T>::from_double(c));
7079 thinkt[tidx] = withc < floorv ? floorv : withc;
7083 if (relax_omega < 1.0 && it > 1 && !std::isnan(thinkt_prev[tidx])) {
7084 const double rawd = dbl(thinkt[tidx]);
7086 thinkt_prev[tidx] = rawd;
7087 thinkt_prev_v[tidx] = thinkt[tidx];
7089 const T om = num_traits<T>::from_double(relax_omega);
7090 const T om1 = num_traits<T>::from_double(1.0 - relax_omega);
7091 thinkt[tidx] = T(om * thinkt[tidx] + om1 * thinkt_prev_v[tidx]);
7093 thinkt_prev[tidx] = dbl(thinkt[tidx]);
7094 thinkt_prev_v[tidx] = thinkt[tidx];
7102 void update_layers(
int it) {
7105 if (is_ph_encoding()) {
7106 update_layers_ph(it);
7109 const bool elevator = (it % 2) == 1;
7110 const std::size_t nt = thinkt_map.size();
7111 for (std::size_t r = 0; r < nt; ++r) {
7112 const UpdRow& row = thinkt_map[elevator ? nt - 1 - r : r];
7113 qn::Layer<T>& L = ensemble[std::size_t(idxhash[row.idx])];
7114 if (row.aidx <= NT() && opt.interlocking &&
7115 L.classes[row.cls - 1].type == JobClassType::CLOSED)
7116 L.classes[row.cls - 1].population = njobs(row.aidx, row.idx);
7117 if (row.node == L.clientIdx) {
7118 if (lqn.type[row.aidx] == LqnElement::TASK) {
7119 if (lqn.sched[row.aidx] != SchedStrategy::REF) {
7120 if (!thinktproc[row.aidx].disabled)
7121 L.set_service(row.node, row.cls, thinktproc[row.aidx]);
7123 L.set_service(row.node, row.cls, servtproc[row.aidx]);
7126 L.set_service(row.node, row.cls, servtproc[row.aidx]);
7129 L.set_service(row.node, row.cls, servtproc[row.aidx]);
7132 const std::size_t nc = call_map.size();
7133 for (std::size_t r = 0; r < nc; ++r) {
7134 const UpdRow& row = call_map[elevator ? nc - 1 - r : r];
7135 qn::Layer<T>& L = ensemble[std::size_t(idxhash[row.idx])];
7136 if (row.node == L.clientIdx) {
7137 L.set_service(row.node, row.cls, callservtproc[row.aidx]);
7139 const std::size_t eidx = lqn.callpair_dst[row.aidx];
7146 L.set_service(row.node, row.cls, servtproc[eidx]);
7157 if (interlockMethod ==
"refpath") update_ref_path_stages(it);
7158 for (
const UpdRow& row : arv_call_map) {
7159 qn::Layer<T>& L = ensemble[std::size_t(idxhash[row.idx])];
7160 const Distrib<T>& d = tputproc[lqn.callpair_src[row.aidx]];
7161 if (!d.disabled) L.set_service(row.node, row.cls, d);
7175 void update_ref_path_stages(
int it) {
7176 if (refpath_map.empty())
return;
7177 const std::size_t nkeys = refpath_keys.size();
7178 std::vector<T> stage(nkeys, Tzero());
7179 for (std::size_t r = 0; r < refpath_map.size(); ++r) {
7180 const RefPathRow& row = refpath_map[r];
7184 term = residt[row.elem];
7189 v = entryvisits[row.fromentry];
7190 term = T(callresidt[row.elem] / v);
7194 term = thread_wait(row.elem);
7197 if (!std::isfinite(dbl(term))) term = Tzero();
7198 stage[refpath_group[r]] = T(stage[refpath_group[r]] + term * num_traits<T>::from_double(row.sign));
7200 const double omega = relax_omega;
7201 for (std::size_t k = 0; k < nkeys; ++k) {
7202 if (stage[k] < Tzero()) stage[k] = Tzero();
7204 if (omega < 1.0 && it > 1 && !std::isnan(refpath_stage_prev[k])) {
7205 T prevS = refpath_stage_prev_v[k];
7208 stage[k] = T(num_traits<T>::from_double(omega) * stage[k] +
7209 num_traits<T>::from_double(1.0 - omega) * prevS);
7211 refpath_stage_prev[k] = dbl(stage[k]);
7212 refpath_stage_prev_v[k] = stage[k];
7213 qn::Layer<T>& L = ensemble[std::size_t(idxhash[refpath_keys[k][0]])];
7228 T thread_wait(std::size_t tidx)
const {
7229 if (tidx < 1 || tidx > lqn.nidx)
return Tzero();
7230 if (lqn.sched[tidx] == SchedStrategy::INF || !std::isfinite(lqn.maxmult[tidx]))
return Tzero();
7231 if (idxhash[tidx] < 0 || results.empty())
return Tzero();
7232 const std::size_t e = std::size_t(idxhash[tidx]);
7233 const qn::Layer<T>& L = ensemble[e];
7234 const LayerResult<T>& r = results.back()[e];
7235 const std::size_t ks = station_idx_of(L, tidx);
7237 for (std::size_t c = 0; c < L.nclasses; ++c) Qtot += r.QN(ks - 1, c);
7238 const T Btot = T(num_traits<T>::from_double(lqn.maxmult[tidx]) * util[tidx]);
7239 const T Xtot = T(tput[tidx] / num_traits<T>::from_double(std::max(1.0, lqn.repl[tidx])));
7241 T d = T(Qtot - Btot);
7260 T region_wait(std::size_t e, std::size_t cls)
const {
7261 const qn::Layer<T>& L = ensemble[e];
7262 if (L.regions.empty() || results.empty())
return Tzero();
7263 const LayerResult<T>& r = results.back()[e];
7264 const std::size_t k = cls - 1;
7265 std::size_t c = L.nchains;
7266 for (std::size_t cc = 0; cc < L.nchains; ++cc)
7267 if (L.chains[cc][k]) c = cc;
7268 if (c == L.nchains)
return Tzero();
7272 std::size_t rs = L.serverIdx;
7273 for (std::size_t i = 0; i < L.regions[0].members.size(); ++i)
7274 if (L.regions[0].members[i]) { rs = i + 1;
break; }
7277 T inside = Tzero(), tput_srv = Tzero();
7278 for (std::size_t j = 0; j < L.nclasses; ++j) {
7279 if (!L.chains[c][j])
continue;
7280 const double p = L.classes[j].population;
7281 if (!std::isfinite(p))
return Tzero();
7283 for (std::size_t i = 0; i < L.nstations; ++i) inside += r.QN(i, j);
7284 tput_srv += r.TN(rs - 1, j);
7286 const T deficit = T(num_traits<T>::from_double(pop) - inside);
7287 if (!(deficit > Tzero()) || !(tput_srv > Tzero()))
return Tzero();
7288 return T(deficit / tput_srv);
7292 void update_routing_probabilities(
int) {
7293 for (std::size_t u = 0; u < unique_route_idx.size(); ++u) {
7296 const std::size_t idx = unique_route_idx[unique_route_idx.size() - 1 - u];
7297 qn::Layer<T>& L = ensemble[std::size_t(idxhash[idx])];
7298 bool updated =
false;
7299 for (
const RouteRow& r : route_map) {
7300 if (r.idx != idx)
continue;
7301 if (idxhash[r.tidx_caller] < 0)
continue;
7302 const std::size_t cl = std::size_t(idxhash[r.tidx_caller]);
7303 const qn::Layer<T>& CL = ensemble[cl];
7304 const LayerResult<T>& cr = results.back()[cl];
7307 const std::size_t cs = station_idx_of(CL, r.tidx_caller);
7309 for (std::size_t k = 0; k < CL.nclasses; ++k) Xtot += cr.TN(cs - 1, k);
7310 if (!(Xtot > Tzero()))
continue;
7311 T entry_tput = Tzero();
7312 for (
const auto& a : CL.attr_calls)
7314 entry_tput += cr.TN(station_idx_of(CL, lqn.parent[a[3]]) - 1, a[0] - 1);
7315 L.set_route(r.cfrom, r.cto, r.nodefrom, r.nodeto, T(entry_tput / Xtot));
7318 if (updated) L.refresh_chains();
7337 void require_transient_ready()
const {
7338 if (opt.layer_solver !=
"fluid")
7340 "SolverLN: the layered transient is the transient OF EACH LAYER, which only the "
7341 "fluid layer solver has; rebuild the ensemble with layer_solver 'fluid'");
7345 bool has_finite_horizon()
const {
7346 return std::isfinite(opt.timespan_end) && opt.timespan_end > 0.0;
7351 std::vector<std::vector<std::vector<double>>> Q, U, Tp, R;
7358 void run_layer_transient(std::size_t e,
const std::vector<FluidRateSched>& sched,
7359 LnTranLayer& block, TranTraj& traj) {
7360 qn::Layer<T>& L = ensemble[e];
7362 fluid::FluidOptions fo = opt.layer_fluid;
7363 fo.rate_sched = sched;
7369 if (e < layer_tran_init.size() && !layer_tran_init[e].empty())
7370 fo.init_sol = layer_tran_init[e];
7375 const double t_end =
7377 const std::vector<fluid::FluidTranPoint> pts =
7378 detail::ln_fluid_transient(L, fo, t_end, opt.tran_points, opt.tran_grid);
7379 const std::size_t M = L.nstations, K = L.nclasses, P = pts.size();
7381 auto alloc = [&](std::vector<std::vector<std::vector<double>>>& A) {
7382 A.assign(M, std::vector<std::vector<double>>(K, std::vector<double>(P, 0.0)));
7391 for (std::size_t p = 0; p < P; ++p) {
7392 block.t[p] = pts[p].t;
7393 for (std::size_t i = 0; i < M; ++i)
7394 for (std::size_t r = 0; r < K; ++r) {
7395 const double q = pts[p].QN(i, r), u = pts[p].UN(i, r), x = pts[p].TN(i, r);
7396 block.QN[i][r][p] = q;
7397 block.UN[i][r][p] = u;
7398 block.TN[i][r][p] = x;
7399 traj.Q[i][r][p] = q;
7400 traj.U[i][r][p] = u;
7401 traj.Tp[i][r][p] = x;
7413 LnTranSolution tran_avg_decoupled() {
7414 require_transient_ready();
7415 if (results.empty()) iterate();
7417 out.mode =
"decoupled";
7418 out.layers.resize(ensemble.size());
7419 for (std::size_t e = 0; e < ensemble.size(); ++e) {
7421 run_layer_transient(e, std::vector<FluidRateSched>(), out.layers[e], tj);
7443 LnTranSolution tran_avg_coupled() {
7444 require_transient_ready();
7449 if (!has_finite_horizon())
return tran_avg_decoupled();
7450 if (results.empty()) iterate();
7451 const std::size_t E = ensemble.size();
7454 out.mode =
"coupled";
7455 out.layers.resize(E);
7456 std::vector<TranTraj> traj(E), prev(E);
7457 for (std::size_t e = 0; e < E; ++e)
7458 run_layer_transient(e, std::vector<FluidRateSched>(), out.layers[e], traj[e]);
7459 const std::vector<double> tgrid = out.layers.empty() ? std::vector<double>() : out.
layers[0].t;
7461 for (
long iter = 1; iter <= opt.ln_transient_iter_max; ++iter) {
7463 const std::vector<std::vector<FluidRateSched>> sched =
7464 build_rate_sched(recompute_demand(traj, tgrid), tgrid);
7465 for (std::size_t e = 0; e < E; ++e)
7466 run_layer_transient(e, sched[e], out.layers[e], traj[e]);
7468 for (std::size_t e = 0; e < E; ++e)
7469 for (std::size_t i = 0; i < traj[e].Q.size(); ++i)
7470 for (std::size_t r = 0; r < traj[e].Q[i].size(); ++r)
7471 for (std::size_t p = 0; p < traj[e].Q[i][r].size(); ++p)
7472 gap = std::max(gap, std::fabs(traj[e].Q[i][r][p] - prev[e].Q[i][r][p]));
7473 out.iterations = iter;
7475 if (gap < opt.ln_transient_tol)
break;
7482 std::map<std::size_t, std::vector<double>> thinkt;
7483 std::map<std::size_t, std::vector<double>> callservt;
7494 TranDemand recompute_demand(
const std::vector<TranTraj>& traj,
7495 const std::vector<double>& tgrid)
const {
7497 const std::size_t ng = tgrid.size();
7499 for (std::size_t t = 1; t <= lqn.ntasks; ++t) {
7500 const std::size_t tidx = lqn.tshift + t;
7501 if (idxhash[tidx] < 0 || lqn.isref[tidx])
continue;
7502 const std::size_t e = std::size_t(idxhash[tidx]);
7503 const qn::Layer<T>& L = ensemble[e];
7504 const std::size_t s = station_idx_of(L, tidx) - 1;
7506 for (std::size_t c = 1; c <= NT(); ++c) nj = std::max(nj, njobs(tidx, c));
7508 const double userthink = dbl(ref_think_time(tidx));
7509 std::vector<double> tk(ng, 0.0);
7510 for (std::size_t p = 0; p < ng; ++p) {
7511 double U = 0.0, X = 0.0;
7512 for (std::size_t r = 0; r < L.nclasses; ++r) {
7513 U += traj[e].U[s][r][p];
7514 X += traj[e].Tp[s][r][p];
7517 double v = lqn.sched[tidx] == SchedStrategy::INF
7518 ? (nj - U) / Xs - userthink
7519 : nj * std::fabs(1.0 - U) / Xs - userthink;
7521 tk[p] = v + userthink;
7523 out.thinkt[tidx] = tk;
7526 for (std::size_t cidx = 1; cidx <= lqn.ncalls; ++cidx) {
7527 if (lqn.calltype[cidx] != CallType::SYNC)
continue;
7528 const std::size_t eidx = lqn.callpair_dst[cidx];
7529 const std::size_t tidx = lqn.parent[eidx];
7530 if (tidx > NT() || idxhash[tidx] < 0)
continue;
7531 const std::size_t e = std::size_t(idxhash[tidx]);
7532 const qn::Layer<T>& L = ensemble[e];
7533 const std::size_t s = station_idx_of(L, tidx) - 1;
7534 std::vector<double> Rc(ng, 0.0);
7536 for (std::size_t r = 0; r < L.nclasses; ++r) {
7537 if (L.classes[r].attr_kind !=
int(LqnElement::ENTRY))
continue;
7538 if (L.classes[r].attr_idx != eidx)
continue;
7540 for (std::size_t p = 0; p < ng; ++p) Rc[p] += traj[e].R[s][r][p];
7545 for (std::size_t r = 0; r < L.nclasses; ++r)
7546 for (std::size_t p = 0; p < ng; ++p) Rc[p] += traj[e].R[s][r][p];
7548 const double cm = dbl(lqn.callproc_mean[cidx]);
7549 for (std::size_t p = 0; p < ng; ++p) Rc[p] *= cm;
7550 out.callservt[cidx] = Rc;
7566 std::vector<std::vector<FluidRateSched>> build_rate_sched(
7567 const TranDemand& demand,
const std::vector<double>& tgrid)
const {
7568 std::vector<std::vector<FluidRateSched>> out(ensemble.size());
7569 const bool want_think =
7570 opt.ln_transient_channels ==
"both" || opt.ln_transient_channels ==
"thinkt";
7571 const bool want_call =
7572 opt.ln_transient_channels ==
"both" || opt.ln_transient_channels ==
"callservt";
7574 auto add = [&](std::size_t e, std::size_t station, std::size_t cls,
7575 const std::vector<double>& d) {
7576 if (d.empty())
return;
7577 const double dend = d.back();
7582 const double cap = 20.0;
7587 s.rates.resize(d.size());
7588 for (std::size_t p = 0; p < d.size(); ++p) {
7589 const double v = std::min(std::max(d[p], dend / cap), dend * cap);
7590 s.rates[p] = 1.0 / v;
7592 s.nominal = 1.0 / dend;
7593 out[e].push_back(s);
7597 for (
const UpdRow& row : thinkt_map) {
7598 if (idxhash[row.idx] < 0)
continue;
7599 const std::size_t e = std::size_t(idxhash[row.idx]);
7600 if (row.node != ensemble[e].clientIdx)
continue;
7601 if (lqn.type[row.aidx] != LqnElement::TASK)
continue;
7602 if (lqn.sched[row.aidx] == SchedStrategy::REF)
continue;
7603 const auto it = demand.thinkt.find(row.aidx);
7604 if (it == demand.thinkt.end())
continue;
7605 add(e, row.node, row.cls, it->second);
7608 for (
const UpdRow& row : call_map) {
7609 if (idxhash[row.idx] < 0)
continue;
7610 const std::size_t e = std::size_t(idxhash[row.idx]);
7611 if (row.node != ensemble[e].clientIdx)
continue;
7612 const auto it = demand.callservt.find(row.aidx);
7613 if (it == demand.callservt.end())
continue;
7614 add(e, row.node, row.cls, it->second);
7622 LnSolution<T> aggregate() {
7623 const std::size_t N = lqn.nidx;
7625 auto mk = [&](std::vector<T>& v, std::vector<bool>& d) {
7626 v.assign(N + 1, Tzero());
7627 d.assign(N + 1,
false);
7629 std::vector<T>
QN,
UN, RN, TN, PN, SN, WN, AN;
7630 std::vector<bool> dQ, dU, dR, dT, dP, dS, dW, dA;
7631 mk(QN, dQ); mk(UN, dU); mk(RN, dR); mk(TN, dT);
7632 mk(PN, dP); mk(SN, dS); mk(WN, dW); mk(AN, dA);
7633 std::vector<bool> wn_done(N + 1,
false);
7635 for (std::size_t e = 0; e < ensemble.size(); ++e) {
7636 const qn::Layer<T>& L = ensemble[e];
7637 const LayerResult<T>& r = results.back()[e];
7638 const std::size_t clientIdx = L.clientIdx;
7643 const bool has_host_server = !L.host_stations.empty();
7644 for (std::size_t hs : L.host_stations) {
7645 const std::size_t hidx = L.stations[hs - 1].attr_idx;
7650 for (std::size_t k = 0; k < L.nclasses; ++k) {
7651 if (L.classes[k].completes) {
7652 T t = clientIdx > 0 ? r.TN(clientIdx - 1, k) : Tzero();
7653 const T ts = r.TN(hs - 1, k);
7654 TN[hidx] = T(TN[hidx] + (t > ts ? t : ts));
7656 if (L.classes[k].attr_kind ==
int(LqnElement::ACTIVITY)) {
7658 if (station_idx_of_class(L, k) != hs)
continue;
7659 const std::size_t aidx = L.classes[k].attr_idx;
7660 const std::size_t tidx = lqn.parent[aidx];
7663 PN[aidx] = T(PN[aidx] + r.UN(hs - 1, k));
7664 PN[tidx] = T(PN[tidx] + r.UN(hs - 1, k));
7665 PN[hidx] = T(PN[hidx] + r.UN(hs - 1, k));
7671 for (std::size_t k = 0; k < L.nclasses; ++k) {
7672 const int kind = L.classes[k].attr_kind;
7676 const std::size_t serverIdx = station_idx_of_class(L, k);
7677 if (kind ==
int(LqnElement::TASK)) {
7678 const std::size_t tidx = L.classes[k].attr_idx;
7679 if (has_host_server && !dT[tidx]) {
7681 TN[tidx] = r.TN(clientIdx - 1, k);
7683 }
else if (kind ==
int(LqnElement::ENTRY)) {
7684 const std::size_t eidx = L.classes[k].attr_idx;
7691 if (has_host_server && !dT[eidx]) {
7693 TN[eidx] = r.TN(clientIdx - 1, k);
7695 }
else if (kind ==
int(LqnElement::CALL)) {
7696 const std::size_t cidx = L.classes[k].attr_idx;
7697 const std::size_t aidx = lqn.callpair_src[cidx];
7698 if (lqn.calltype[cidx] == CallType::SYNC) {
7700 SN[aidx] = T(SN[aidx] + r.RN(serverIdx - 1, k) * lqn.callproc_mean[cidx]);
7703 QN[aidx] = T(QN[aidx] + r.QN(serverIdx - 1, k));
7704 }
else if (kind ==
int(LqnElement::ACTIVITY)) {
7705 const std::size_t aidx = L.classes[k].attr_idx;
7706 const std::size_t tidx = lqn.parent[aidx];
7708 QN[tidx] = T(QN[tidx] + r.QN(serverIdx - 1, k));
7711 TN[aidx] = T(TN[aidx] + r.TN(serverIdx - 1, k));
7713 SN[aidx] = T(SN[aidx] + r.RN(serverIdx - 1, k));
7715 RN[aidx] = T(RN[aidx] + r.RN(serverIdx - 1, k));
7718 WN[aidx] = residt[aidx];
7719 if (!wn_done[aidx]) {
7720 WN[tidx] = T(WN[tidx] + residt[aidx]);
7721 wn_done[aidx] =
true;
7723 QN[aidx] = T(QN[aidx] + r.QN(serverIdx - 1, k));
7731 for (std::size_t t = 1; t <= lqn.ntasks; ++t) {
7732 const std::size_t tidx = lqn.tshift + t;
7733 const double nrep = lqn.repl[tidx];
7734 if (!(nrep > 1.0))
continue;
7735 const T scale = num_traits<T>::from_double(nrep);
7736 std::vector<std::size_t> scaleIdx{tidx};
7737 scaleIdx.insert(scaleIdx.end(), lqn.entriesof[tidx].begin(), lqn.entriesof[tidx].end());
7738 scaleIdx.insert(scaleIdx.end(), lqn.actsof[tidx].begin(), lqn.actsof[tidx].end());
7739 for (std::size_t idx : scaleIdx) {
7740 if (dT[idx]) TN[idx] = T(scale * TN[idx]);
7741 if (dP[idx]) PN[idx] = T(scale * PN[idx]);
7744 for (std::size_t hidx = 1; hidx <= lqn.nhosts; ++hidx)
7745 if (lqn.repl[hidx] > 1.0 && dP[hidx]) PN[hidx] = T(num_traits<T>::from_double(lqn.repl[hidx]) * PN[hidx]);
7747 for (std::size_t e = 1; e <= lqn.nentries; ++e) {
7748 const std::size_t eidx = lqn.eshift + e;
7749 const std::size_t tidx = lqn.parent[eidx];
7756 ? T(TN[eidx] * (servt_ph1[eidx] + servt_ph2[eidx]))
7757 : T(TN[eidx] * SN[eidx]);
7760 for (std::size_t a : lqn.actsof[eidx])
7765 if (!lqn.actsof[eidx].empty()) {
7767 PN[eidx] = any ? ps : Tzero();
7769 for (std::size_t a : lqn.actsof[tidx]) {
7771 UN[a] = T(TN[a] * SN[a]);
7773 UN[tidx] = T(UN[tidx] + UN[eidx]);
7785 for (std::size_t i = 1; i <= N; ++i)
7790 const bool host = lqn.type[i] == LqnElement::HOST;
7791 const bool task = lqn.type[i] == LqnElement::TASK;
7792 const bool entry = lqn.type[i] == LqnElement::ENTRY;
7796 dS[i] = !host && !task;
7798 dW[i] = !host && !entry;
7804 s.QN =
UN; s.defined_Q = dU;
7805 s.UN = PN; s.defined_U = dP;
7806 s.RN = SN; s.defined_R = dS;
7807 s.TN = TN; s.defined_T = dT;
7808 s.AN = AN; s.defined_A = dA;
7809 s.WN = WN; s.defined_W = dW;
7810 s.iterations = iterations_done;
7811 s.converged = did_converge;