Definition:
class rsc
{
public:
// build objects
int nbX;
int nbY;
int nbParams;
bool outsideOption;
double* dist_pars;
arma::mat zeta;
arma::mat aux_ord;
arma::cube aux_Influence_lhs;
arma::cube aux_Influence_rhs;
arma::cube aux_DinvPsigma;
arma::cube aux_Psigma;
double (*aux_cdf_eps)(double x, double* dist_pars);
double (*aux_quant_eps)(double x, double* dist_pars);
double (*aux_pdf_eps)(double x, double* dist_pars);
double (*aux_pot_eps)(double x, double* dist_pars);
arma::vec (*aux_cdf_eps_vec)(arma::vec x, double* dist_pars);
arma::vec (*aux_quant_eps_vec)(arma::vec x, double* dist_pars);
arma::vec (*aux_pdf_eps_vec)(arma::vec x, double* dist_pars);
arma::vec (*aux_pot_eps_vec)(arma::vec x, double* dist_pars);
// input objects
arma::mat U;
arma::mat mu;
// equilibrium objects
arma::mat U_sol;
arma::mat mu_sol;
// member functions
~rsc(){};
rsc(){};
explicit rsc(int nbX_inp, int nbY_inp);
explicit rsc(arma::mat zeta_inp, bool outsideOption_inp);
explicit rsc(arma::mat zeta_inp, double alpha, double beta);
void build(int nbX_inp, int nbY_inp);
void build(arma::mat zeta_inp, bool outsideOption_inp);
void build_beta(arma::mat zeta_inp, double alpha, double beta);
double G(const arma::vec& n);
double G(const arma::vec& n, const arma::mat& U_inp, arma::mat& mu_out);
double Gx(const arma::mat& U_x_inp, arma::mat& mu_x_out, int x);
double Gstar(const arma::vec& n);
double Gstar(const arma::vec& n, const arma::mat& mu_inp, arma::mat& U_out);
double Gstarx(const arma::mat& mu_x_inp, arma::mat &U_x_out, int x);
double Gbar(const arma::mat& Ubar, const arma::mat& mubar, const arma::vec& n, arma::mat& U_out, arma::mat& mu_out);
double Gbarx(const arma::vec& Ubar_x, const arma::vec& mubar_x, arma::mat& U_x_out, arma::mat& mu_x_out, int x);
arma::mat D2Gstar(const arma::vec& n, bool xFirst);
void D2Gstar(arma::mat &H, const arma::vec& n, bool xFirst);
void D2Gstar(arma::mat &H, const arma::vec& n, const arma::mat& mu_inp, bool xFirst);
void dtheta_NablaGstar(arma::mat& ret, const arma::vec& n, arma::mat* dtheta, bool x_first);
empirical simul();
empirical simul(int* nbDraws, int* seed);
void simul(empirical& obj_out);
void simul(empirical& obj_out, int* nbDraws, int* seed);
private:
double cdf (double x);
arma::vec cdf (arma::vec x);
double pdf (double x);
arma::vec pdf (arma::vec x);
double quantile (double x);
arma::vec quantile (arma::vec x);
double pot (double x);
arma::vec pot (arma::vec x);
static double Gbar_opt_objfn(const arma::vec& vals_inp, arma::vec* grad, void* opt_data);
static double Gbar_opt_constr(const arma::vec& vals_inp, arma::vec* grad, void* constr_data);
static double Gstarx(arma::vec& U_x, const arma::vec& mu_x_inp, const arma::mat& zeta,
const arma::mat& aux_DinvPsigma, const arma::mat& aux_Psigma,
const arma::mat& aux_Influence_lhs, const arma::mat& aux_Influence_rhs,
arma::vec (*pot_eps_vec)(arma::vec pot_inp, double* dist_pars),
arma::vec (*quantile_eps_vec)(arma::vec quant_inp, double* dist_pars),
double* dist_pars, int nbY, int x);
};
Example:
arma::mat U(2,3);
U << 1.6 << 3.2 << 1.1 << arma::endr
<< 2.9 << 1.0 << 3.1 << arma::endr;
arma::mat mu(2,3);
mu << 1.0 << 3.0 << 1.0 << arma::endr
<< 2.0 << 1.0 << 3.0 << arma::endr;
//
// setup RSC class object
int nbX = U.n_rows;
int nbY = U.n_cols;
arma::vec n = arma::sum(mu,1) + 1.0;
arma::mat z_temp(1,nbY+1);
z_temp(0,0) = 0.1; z_temp(0,1) = 0.2; z_temp(0,2) = 0.3; z_temp(0,3) = 0.0;
arma::mat zeta = arma::ones(nbX,1) * z_temp;
arma::cout << "zeta: \n" << zeta << arma::endl;
//
// RSC object
trame::rsc rsc_obj;
rsc_obj.U = U;
rsc_obj.mu = mu;
rsc_obj.build_beta(zeta,2.0,2.0);
//
// empirical object:
int sim_seed = 1777;
int n_draws = 1000;
trame::empirical rsc_sim;
rsc_obj.simul(rsc_sim, &n_draws, &sim_seed);
rsc_sim.U = U;
rsc_sim.mu = mu;
//
// first compute optimal assignment (mu)
double G_val = rsc_obj.G(n);
double G_sim_val = rsc_sim.G(n);
std::cout << "G(U) and G-sim(U): \n" << G_val << " and " << G_sim_val << std::endl;
arma::cout << "\nG -> mu: \n" << rsc_obj.mu_sol << arma::endl;
arma::cout << "G-sim -> mu: \n" << rsc_sim.mu_sol << arma::endl;
//
// solution to dual problem U*
double Gstar_val = rsc_obj.Gstar(n);
double Gstar_sim_val = rsc_sim.Gstar(n);
std::cout << "G*(mu) and G*-sim(mu): \n" << Gstar_val << " and " << Gstar_sim_val << std::endl;
arma::cout << "\n\\nabla G*(\\nabla G(U)): \n" << rsc_obj.U_sol << arma::endl;
arma::cout << "\\nabla G-sim*(\\nabla G-sim(U)): \n" << rsc_sim.U_sol << arma::endl;
//
// Gbar
arma::mat mu_bar(2,3);
mu_bar.fill(2);
arma::mat U_bar_temp, mu_bar_temp;
arma::mat U_bar_sim_temp, mu_bar_sim_temp;
double val_Gbar = rsc_obj.Gbar(U,mu_bar,n,U_bar_temp,mu_bar_temp);
double val_Gbar_sim = rsc_sim.Gbar(U,mu_bar,n,U_bar_sim_temp,mu_bar_sim_temp);
std::cout << "Gbar val: \n" << val_Gbar << std::endl;
std::cout << "Gbar-sim val: \n" << val_Gbar_sim << std::endl;