solver_mam_sync.h

Go to the documentation of this file.

/*
 * Copyright (c) 2012-2026, QORE Lab, Imperial College London
 * All rights reserved.
 */
#ifndef LINE_SOLVERS_MAM_SOLVER_MAM_SYNC_H
#define LINE_SOLVERS_MAM_SOLVER_MAM_SYNC_H
 
/**
 * @file
 * @ingroup line_solvers
 * The 'sync.faster' and 'sync.slower' methods of SolverMAM: the synchronous
 * approximation of an M/G/c/K station with dual bounds, qsys_mgck_sync.
 *
 * Port of matlab/src/solvers/MAM/mam_sync_applicable.m and solver_mam_sync.m.
 * 'sync.faster' returns lower bounds on the mean queue length, response time
 * and loss probability, 'sync.slower' upper ones; options.config.sync_phases
 * (MamOptions::sync_phases) selects the phase expansion.
 *
 * THE RULE IS A SHAPE: the synchronous approximation is a model of one M/G/c/K
 * station, so the network must be a single open class with one Source (Poisson)
 * and one FCFS Queue with a finite capacity and the drop rule, every completion
 * leaving the network. The service law must be one the embedded chain can take
 * exactly: deterministic, or phase type. The phase expansion needs it
 * deterministic.
 *
 * ONE PREDICATE, THREE CALLERS: `runner_detail::check_model_method`, the
 * `mam_dispatch` arm and `solver_mam_sync` all ask `mam_sync_refusal`.
 */
 
#include <cmath>
#include <cstddef>
#include <string>
#include <vector>
 
#include "line/api/mam/map_moment.h"
#include "line/api/qsys/qsys_mgck_sync.h"
#include "line/lang/distribution.h"
#include "line/lang/qn/network_struct.h"
#include "line/solvers/mam/mam_types.h"
#include "line/util/error.h"
 
namespace line {
namespace mam {
 
/** What qsys_mgck_sync needs, read off the model by `mam_sync_refusal`. */
template <class T>
struct MamSyncParams {
    T lambda = num_traits<T>::from_int(0);
    std::size_t c = 1;
    std::size_t K = 1;
    qsys::ServiceLaw<T> svc;
    std::size_t nphases = 0;
    std::size_t srcIdx = 0;  ///< 1-based Source station
    std::size_t qIdx = 0;    ///< 1-based Queue station
};
 
namespace sync_detail {
 
inline const char* drop_to_text(lang::DropStrategy d) {
    switch (d) {
        case lang::DropStrategy::WAITQ: return "WAITQ";
        case lang::DropStrategy::DROP: return "DROP";
        case lang::DropStrategy::BAS: return "BAS";
        case lang::DropStrategy::BBS: return "BBS";
        case lang::DropStrategy::RSRD: return "RSRD";
        case lang::DropStrategy::RETRIAL: return "RETRIAL";
        case lang::DropStrategy::RETRIAL_WITH_LIMIT: return "RETRIAL_WITH_LIMIT";
    }
    return "unknown";
}
 
}  // namespace sync_detail
 
/**
 * `mam_sync_applicable.m`: why the sync methods cannot answer this model, or ""
 * when they can, in which case `par` (when given) is filled.
 */
template <class T>
std::string mam_sync_refusal(const qn::NetworkStruct<T>& L, std::size_t nph,
                             MamSyncParams<T>* par = nullptr) {
    using lang::ProcessType;
    using lang::SchedStrategy;
    const T zero = num_traits<T>::from_int(0), one = num_traits<T>::from_int(1);
    if (L.nclasses != 1 || !std::isinf(L.classes[0].population))
        return "The sync methods require a single open class.";
    std::size_t nsrc = 0, nq = 0, src = 0, q = 0;
    for (std::size_t i = 0; i < L.nstations; ++i) {
        if (L.stations[i].sched == SchedStrategy::EXT) {
            ++nsrc;
            src = i + 1;
        } else {
            ++nq;
            q = i + 1;
        }
    }
    if (L.nstations != 2 || nsrc != 1 || nq != 1)
        return "The sync methods require exactly one Source and one Queue station.";
    if (L.stations[q - 1].sched != SchedStrategy::FCFS)
        return std::string("The sync methods require an FCFS queue; the station uses ") +
               lang::sched_to_text(L.stations[q - 1].sched) +
               ". The synchronous approximation is a model of the M/G/c/K FIFO queue.";
    if (L.procid(src, 1) != ProcessType::EXP || L.phases_of(src, 1) > 1)
        return "The sync methods require Poisson (exponential) arrivals.";
    // every completion at the queue leaves the network: no feedback to the queue
    const std::size_t qrow = (q - 1) * L.nclasses;
    if (L.rt.rows() > qrow && L.rt(qrow, qrow) > zero)
        return "The sync methods do not support feedback: every completion must leave the network.";
    bool ld = false;
    for (const T& v : L.stations[q - 1].lldscaling)
        if (v != one) ld = true;
    if (L.setupparam.find(q) != L.setupparam.end() || ld)
        return "The sync methods do not support setup times or load dependence.";
    const double c = L.stations[q - 1].nservers;
    const double K = L.stations[q - 1].cap;
    if (!std::isfinite(c) || c < 1.0) return "The sync methods require a finite number of servers.";
    if (!std::isfinite(K))
        return "The sync methods require a finite capacity: the synchronous approximation is an "
               "M/G/c/K model. Set Queue.setCapacity(K).";
    if (K < c)
        return "The sync methods require a capacity no smaller than the number of servers.";
    if (!L.droprule.empty() && L.droprule[q - 1][0] != lang::DropStrategy::DROP)
        return std::string("The sync methods require the drop rule at the finite buffer; the "
                           "station uses ") +
               sync_detail::drop_to_text(L.droprule[q - 1][0]) + ".";
    // the declared law: sn_nonmarkov_toph retags a non-Markovian fit, only DET survives
    const lang::Distrib<T>& dq = L.service[q - 1][0];
    const lang::Distrib<T>& decl = dq.declared ? *dq.declared : dq;
    const ProcessType pid = decl.type;
    qsys::ServiceLaw<T> svc;
    if (pid == ProcessType::DET) {
        svc = qsys::ServiceLaw<T>::deterministic(T(one / L.rates(q - 1, 0)));
    } else if (pid == ProcessType::EXP || pid == ProcessType::ERLANG ||
               pid == ProcessType::HYPEREXP || pid == ProcessType::PH ||
               pid == ProcessType::APH || pid == ProcessType::COXIAN ||
               pid == ProcessType::COX2) {
        if (nph > 0)
            return "The synchronous phase expansion (options.config.sync_phases >= 1) requires a "
                   "deterministic service time.";
        const Map<T> m = lang::dist_to_map(dq);
        svc = qsys::ServiceLaw<T>::phase_type(map_pie(m), m.D0);
    } else {
        return std::string("The sync methods take a deterministic or phase-type service time; "
                           "the station declares ") +
               lang::process_to_text(pid) + ".";
    }
    if (par) {
        par->lambda = L.rates(src - 1, 0);
        par->c = static_cast<std::size_t>(c);
        par->K = static_cast<std::size_t>(K);
        par->svc = svc;
        par->nphases = nph;
        par->srcIdx = src;
        par->qIdx = q;
    }
    return std::string();
}
 
/** Port of `solver_mam_sync.m`: fills Q, U, R, T, C, X from qsys_mgck_sync. */
template <class T>
mva::MvaSolution<T> solver_mam_sync(const qn::NetworkStruct<T>& L, const MamOptions& opt) {
    MamSyncParams<T> par;
    const std::string why = mam_sync_refusal(L, opt.sync_phases, &par);
    if (!why.empty()) throw UnsupportedError(ErrorKind{"Structure"}, "SolverMAM: " + why);
    const std::string prefix = "sync.";
    if (opt.method.compare(0, prefix.size(), prefix) != 0)
        throw UnsupportedError(ErrorKind{"Method"},
                               "SolverMAM: solver_mam_sync answers 'sync.faster' and 'sync.slower' "
                               "only, not '" + opt.method + "'");
    const std::string variant = opt.method.substr(prefix.size());
    const qsys::MgckSyncResult<T> r =
        qsys::qsys_mgck_sync(par.lambda, par.c, par.K, par.svc, variant, par.nphases);
 
    const T zero = num_traits<T>::from_int(0);
    const std::size_t M = L.nstations, K = L.nclasses;
    mva::MvaSolution<T> s;
    s.Q = Matrix<T>(M, K, zero);
    s.U = Matrix<T>(M, K, zero);
    s.R = Matrix<T>(M, K, zero);
    s.Tp = Matrix<T>(M, K, zero);
    s.C.assign(K, zero);
    s.X.assign(K, zero);
    const std::size_t si = par.srcIdx - 1, qi = par.qIdx - 1;
    s.Tp(si, 0) = r.throughput;
    s.Q(qi, 0) = r.meanQueueLength;
    s.U(qi, 0) = r.utilization;
    s.R(qi, 0) = r.meanSojournTime;
    s.Tp(qi, 0) = r.throughput;
    s.X[0] = r.throughput;
    s.C[0] = r.meanSojournTime;
    s.iter = 1;
    s.method = opt.method;
    return s;
}
 
}  // namespace mam
}  // namespace line
 
#endif  // LINE_SOLVERS_MAM_SOLVER_MAM_SYNC_H