Commit b64128cd authored by Damien Leroux's avatar Damien Leroux
Browse files

There might still be a crash when exiting but it's mostly stable otherwise.

parent 44494e56
This diff is collapsed.
......@@ -33,8 +33,9 @@
struct chrono_trace {
const std::string& name;
chrono_trace(const std::string& n) : name(n) { chrono::increment(name); chrono::start(name); }
~chrono_trace() { chrono::stop(name); }
double t0;
chrono_trace(const std::string& n) : name(n), t0(chrono::get_t0()) { chrono::increment(name); }
~chrono_trace() { chrono::stop(name, t0); }
};
#include "error.h"
......
......@@ -108,6 +108,7 @@ template <typename Ret, typename... Args>
*/
std::string func_name = get_func_name(func);
chrono_trace _(func_name);
/*
if (0) {
// msg_handler_t::cout() << "INVOKING " << func_name << "(…)" << std::endl;
std::stringstream ss;
......@@ -119,17 +120,18 @@ template <typename Ret, typename... Args>
}
msg_handler_t::cout() << ss.str() << ')';
}
std::thread::id this_id = std::this_thread::get_id();
// std::thread::id this_id = std::this_thread::get_id();
// active_settings->thread_stacks[this_id].push_back(func_name);
MSG_DEBUG("[" << this_id << ',' << (this_id == active_settings->main_thread) << "] ENTER " << func_name);
// MSG_DEBUG("[" << this_id << ',' << (this_id == active_settings->main_thread) << "] ENTER " << func_name);
msg_handler_t::run_hooks();
*/
auto ret = func(*args...);
// active_settings->thread_stacks[this_id].pop_back();
msg_handler_t::run_hooks();
unregister_task_in_progress(func, {*args}...);
TaskPool::release_slot();
// TaskPool::remove_task(m_thread->get_id());
MSG_DEBUG("[" << this_id << ',' << (this_id == active_settings->main_thread) << "] LEAVE " << func_name);
// MSG_DEBUG("[" << this_id << ',' << (this_id == active_settings->main_thread) << "] LEAVE " << func_name);
m_started_cv.notify_all();
return ret;
})
......@@ -171,6 +173,7 @@ template <typename Ret, typename... Args>
*/
std::string func_name = get_func_name(func);
chrono_trace _(func_name);
/*
if (0) {
// msg_handler_t::cout() << "INVOKING " << func_name << "(" << std::endl;
std::stringstream ss;
......@@ -182,13 +185,14 @@ template <typename Ret, typename... Args>
}
msg_handler_t::cout() << ss.str() << ')';
}
std::thread::id this_id = std::this_thread::get_id();
// std::thread::id this_id = std::this_thread::get_id();
// active_settings->thread_stacks[this_id].push_back(func_name);
MSG_DEBUG("[" << this_id << ',' << (this_id == active_settings->main_thread) << "] p ENTER " << func_name);
// MSG_DEBUG("[" << this_id << ',' << (this_id == active_settings->main_thread) << "] p ENTER " << func_name);
// msg_handler_t::run_hooks();
// m_promise.set_value(proxy(*args...));
*/
auto ret = proxy(*args...);
MSG_DEBUG("[" << this_id << ',' << (this_id == active_settings->main_thread) << "] p LEAVE " << func_name);
// MSG_DEBUG("[" << this_id << ',' << (this_id == active_settings->main_thread) << "] p LEAVE " << func_name);
// active_settings->thread_stacks[this_id].pop_back();
// msg_handler_t::run_hooks();
unregister_task_in_progress(func, {*args}...);
......
......@@ -72,6 +72,8 @@ struct chrono {
double t0;
double accum;
static std::mutex& mutex() { static std::mutex _; return _; }
static double time_()
{
struct timeval tv;
......@@ -84,10 +86,14 @@ struct chrono {
void start() { t0 = time_(); }
void stop() { accum += time_() - t0; }
static void start(const std::string& s) { chronos()[s].start(); }
static void stop(const std::string& s) { chronos()[s].stop(); }
static void start(const std::string& s) { std::unique_lock<std::mutex> _(mutex()); chronos()[s].start(); }
static void stop(const std::string& s) { std::unique_lock<std::mutex> _(mutex()); chronos()[s].stop(); }
static double get_t0() { return time_(); }
static void stop(const std::string& s, double t0) { std::unique_lock<std::mutex> _(mutex()); chronos()[s].accum += time_() - t0; }
static void increment(const std::string& s) { counters()[s]++; }
static void increment(const std::string& s) { std::unique_lock<std::mutex> _(mutex()); counters()[s]++; }
};
......
......@@ -2064,6 +2064,7 @@ struct model_manager {
auto sel = get_selection();
MSG_DEBUG("Computing LOD values for QTLs in selection " << sel);
MSG_INFO("Keys before " << vMcurrent->keys());
for (const auto& chr_lk: sel) {
const chromosome* chr = chr_lk.first;
for (double loc: chr_lk.second) {
......@@ -2075,6 +2076,7 @@ struct model_manager {
ret.emplace_back(chr->name, loc, sub.last_test_positions[chr], *sub.last_computation[chr]);
}
}
MSG_INFO("Keys after " << vMcurrent->keys());
#if 0
for (const auto& ul: uniq_loci) {
/*chrom_under_study = ul.selection.begin()->first;*/
......
......@@ -410,7 +410,7 @@ compute_along_chromosome(computation_along_chromosome& ret,
const std::vector<double>& loci,
const collection<parental_origin_per_locus_type>& popl)
{
MSG_DEBUG("INIT COMPUTATION ct=" << ct << " cr=" << cr);
// MSG_DEBUG("INIT COMPUTATION ct=" << ct << " cr=" << cr);
if (cr & RSS) {
ret.rss.resize(Mcurrent->Y().outerSize(), popl.size());
}
......@@ -430,14 +430,14 @@ compute_along_chromosome(computation_along_chromosome& ret,
ret.ftest_lod.resize(Mcurrent->Y().outerSize(), popl.size());
}
if (ct & Chi2) {
MSG_DEBUG("INIT CHI2");
MSG_DEBUG("Y.cols = " << Mcurrent->Y().cols() << " y_block_size = " << y_block_size << " popl.size = " << popl.size());
// MSG_DEBUG("INIT CHI2");
// MSG_DEBUG("Y.cols = " << Mcurrent->Y().cols() << " y_block_size = " << y_block_size << " popl.size = " << popl.size());
MSG_QUEUE_FLUSH();
ret.chi2.resize(Mcurrent->Y().cols() / y_block_size, popl.size());
}
if (ct & Chi2LOD) {
MSG_DEBUG("INIT CHI2 LOD");
MSG_DEBUG("Y.cols = " << Mcurrent->Y().cols() << " y_block_size = " << y_block_size << " popl.size = " << popl.size());
// MSG_DEBUG("INIT CHI2 LOD");
// MSG_DEBUG("Y.cols = " << Mcurrent->Y().cols() << " y_block_size = " << y_block_size << " popl.size = " << popl.size());
MSG_QUEUE_FLUSH();
ret.chi2_lod.resize(Mcurrent->Y().cols() / y_block_size, popl.size());
}
......
......@@ -881,7 +881,7 @@ struct model {
model(const model& mo)
: m_Y(mo.m_Y)
, m_blocks(mo.m_blocks), m_X(mo.m_X)
, m_blocks(), m_X(mo.m_X)
, m_rank(mo.m_rank), m_rss(mo.m_rss), m_coefficients(mo.m_coefficients)
, m_residuals(mo.m_residuals)
, m_solver_type(mo.m_solver_type)
......@@ -893,11 +893,16 @@ struct model {
, m_max_order(mo.m_max_order)
, m_threshold(mo.m_threshold)
/*{ MSG_DEBUG("new model " << __LINE__ << " with Y(" << m_Y->innerSize() << ',' << m_Y->outerSize() << ')'); }*/
{}
{
for (const auto& kv: mo.m_blocks) {
model_block_key mbk = std::make_shared<model_block_key_struc>(*kv.first);
m_blocks.emplace(mbk, kv.second);
}
}
model(const value<MatrixXd>& y, const model& mo)
: m_Y(y)
, m_blocks(mo.m_blocks), m_X(mo.m_X)
, m_blocks(), m_X(mo.m_X)
, m_rank(), m_rss(), m_coefficients()
, m_residuals()
, m_solver_type(mo.m_solver_type)
......@@ -909,7 +914,12 @@ struct model {
, m_max_order(mo.m_max_order)
, m_threshold(mo.m_threshold)
/*{ MSG_DEBUG("new model " << __LINE__ << " with Y(" << m_Y->innerSize() << ',' << m_Y->outerSize() << ')'); }*/
{}
{
for (const auto& kv: mo.m_blocks) {
model_block_key mbk = std::make_shared<model_block_key_struc>(*kv.first);
m_blocks.emplace(mbk, kv.second);
}
}
# if 0
model&
......@@ -1100,7 +1110,7 @@ struct model {
bool
add_block_if_test_is_good(const model_block_key& key, const value<model_block_type>& x, const model& ref)
{
MatrixXd pvalue(Y().cols(), 1);
MatrixXd pvalue(m_Y->cols(), 1);
model tmp(*this);
tmp.add_block(key, x);
tmp.compute();
......@@ -1243,7 +1253,7 @@ struct model {
void compute()
{
MSG_DEBUG("Recomputing model with selection " << keys());
// MSG_DEBUG("Recomputing model with selection " << keys());
if (m_computed) {
return;
}
......@@ -1260,7 +1270,7 @@ struct model {
}
}
}
MSG_DEBUG("Filtered selection " << keys());
// MSG_DEBUG("Filtered selection " << keys());
m_computed = true;
m_X = new immediate_value<MatrixXd>(MatrixXd());
......@@ -1293,7 +1303,7 @@ struct model {
}
/**m_X = MatrixXd(n_rows, n_cols);*/
*m_X = MatrixXd::Zero(n_rows, n_cols);
MSG_DEBUG("n_rows=" << n_rows << " n_cols=" << n_cols);
// MSG_DEBUG("n_rows=" << n_rows << " n_cols=" << n_cols);
/*MSG_DEBUG("DEBUG X" << std::endl << (*m_X));*/
/* now fill the matrix */
n_rows = 0;
......@@ -1309,7 +1319,7 @@ struct model {
for (const auto& cmap: constraints) {
n_cols = 0;
size_t row_incr = 0;
MSG_DEBUG("row=" << n_rows << " col=" << n_cols);
// MSG_DEBUG("row=" << n_rows << " col=" << n_cols);
MSG_QUEUE_FLUSH();
for (const auto& b: m_blocks) {
auto it = cmap.find(b.first);
......@@ -1325,7 +1335,7 @@ struct model {
n_rows += row_incr;
/*MSG_DEBUG("DEBUG X" << std::endl << (*m_X));*/
}
(void) Y();
// (void) Y();
// m_Y->conservativeResize(n_rows, Eigen::NoChange_t());
// int dr = n_rows - n_pop_rows;
// int c = m_Y->outerSize();
......@@ -1362,16 +1372,17 @@ struct model {
solver.setThreshold(sqrt(COMPONENT_EPSILON));
m_rank = solver.rank();
/*MSG_DEBUG("Rank: " << m_rank);*/
m_coefficients->resize(m_columns, Y().outerSize());
MatrixXd fullY = Y();
m_coefficients->resize(m_columns, fullY.outerSize());
/*MSG_DEBUG(MATRIX_SIZE((*m_coefficients)));*/
for (int i = 0; i < Y().outerSize(); ++i) {
for (int i = 0; i < fullY.outerSize(); ++i) {
/*MSG_DEBUG(MATRIX_SIZE(Y().col(i)) << std::endl << Y().col(i));*/
MatrixXd tmp = solver.solve(Y().col(i));
MatrixXd tmp = solver.solve(fullY.col(i));
/*MSG_DEBUG(MATRIX_SIZE(tmp) << std::endl << tmp);*/
m_coefficients->col(i) = tmp;
}
/*MSG_DEBUG("Coefficients:" << std::endl << m_coefficients);*/
*m_residuals = Y() - X() * (*m_coefficients);
*m_residuals = fullY - X() * (*m_coefficients);
/*MSG_DEBUG("Residuals:" << std::endl << m_residuals);*/
*m_rss = m_residuals->array().square().colwise().sum();
/*MSG_DEBUG("RSS:" << std::endl << m_rss);*/
......@@ -1416,22 +1427,26 @@ struct model {
return m_rank;
}
const MatrixXd& Y() const
MatrixXd Y() const
{
if (!m_Y) {
MSG_ERROR("NULL Y!", "Fix the code");
throw 0;
}
if (m_X) {
if (m_Y->rows() > m_X->rows()) {
model* unconst = const_cast<model*>(this);
unconst->m_Y->conservativeResize(m_X->rows(), Eigen::NoChange_t{});
} else if (m_Y->rows() < m_X->rows()) {
model* unconst = const_cast<model*>(this);
auto yrows = m_Y->rows();
unconst->m_Y->conservativeResize(m_X->rows(), Eigen::NoChange_t{});
unconst->m_Y->bottomRows(m_X->rows() - yrows).setZero();
}
MatrixXd ret;
ret = MatrixXd::Zero(m_X->rows(), m_Y->cols());
ret.topRows(std::min(m_X->rows(), m_Y->rows())) = m_Y->topRows(std::min(m_X->rows(), m_Y->rows()));
return ret;
// if (m_Y->rows() > m_X->rows()) {
// model* unconst = const_cast<model*>(this);
// unconst->m_Y->conservativeResize(m_X->rows(), Eigen::NoChange_t{});
// } else if (m_Y->rows() < m_X->rows()) {
// model* unconst = const_cast<model*>(this);
// auto yrows = m_Y->rows();
// unconst->m_Y->conservativeResize(m_X->rows(), Eigen::NoChange_t{});
// unconst->m_Y->bottomRows(m_X->rows() - yrows).setZero();
// }
}
return *m_Y;
}
......@@ -1884,14 +1899,14 @@ public:
{
static std::string traitstr("Trait");
std::vector<std::vector<char>> traitlabels;
for (int i = 0; i < m.Y().cols(); ++i) {
for (int i = 0; i < m.m_Y->cols(); ++i) {
std::string str(MESSAGE((i + 1)));
traitlabels.emplace_back(str.begin(), str.end());
}
/*os << "<model Y(" << m.m_Y->innerSize() << ',' << m.m_Y->outerSize() << "), " << m.m_blocks.size() << " blocks>";*/
MatrixXd big(m.X().rows(), m.Y().cols() + m.X().cols());
big.topLeftCorner(m.Y().rows(), m.Y().cols()) = m.Y();
big.bottomLeftCorner(big.rows() - m.Y().rows(), m.Y().cols()).setConstant(0);
MatrixXd big(m.X().rows(), m.m_Y->cols() + m.X().cols());
big.topLeftCorner(m.m_Y->rows(), m.m_Y->cols()) = *m.m_Y;
big.bottomLeftCorner(big.rows() - m.m_Y->rows(), m.m_Y->cols()).setConstant(0);
big.rightCols(m.X().cols()) = m.X();
model_print::matrix_with_sections<std::string, void, std::string, std::vector<char>> mws(big);
/*model_print::matrix_with_sections<std::string, void, model_block_key, void> mws(m.X());*/
......
......@@ -156,7 +156,7 @@ struct TaskPool {
m_slots_cv.wait(lock, [this]() { return m_free_slots > 0 || m_quit; });
}
m_free_slots--;
MSG_DEBUG("*** TOOK ONE SLOT. " << m_free_slots.load() << " LEFT.")
// MSG_DEBUG("*** TOOK ONE SLOT. " << m_free_slots.load() << " LEFT.")
}
void
......@@ -168,7 +168,7 @@ struct TaskPool {
{
std::unique_lock<std::mutex> lock(m_slots_mutex);
m_free_slots++;
MSG_DEBUG("*** FREED ONE SLOT. " << m_free_slots.load() << " LEFT.")
// MSG_DEBUG("*** FREED ONE SLOT. " << m_free_slots.load() << " LEFT.")
}
m_slots_cv.notify_all();
}
......
......@@ -132,6 +132,7 @@ int main(int argc, const char** argv)
atexit(delete_settings);
chrono::display() = true;
chrono::start("*** initialization");
active_settings = settings_t::from_args(argc, argv);
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment