diff --git a/src/FollowUp/Nomad/FStatNomad.cpp b/src/FollowUp/Nomad/FStatNomad.cpp index 6490b6624de34473271ab763f223aa8fc989a182..8e6d7cf0434a76798832e58f088756c092600bcd 100644 --- a/src/FollowUp/Nomad/FStatNomad.cpp +++ b/src/FollowUp/Nomad/FStatNomad.cpp @@ -37,10 +37,12 @@ #include <omp.h> #endif -#include "FStatNomad.h" +#include <nomad.hpp> +#include "FStatNomad.h" using namespace std; +using namespace NOMAD; #define EARTHEPHEMERIS "earth00-19-DE405.dat" #define SUNEPHEMERIS "sun00-19-DE405.dat" @@ -63,10 +65,10 @@ using namespace std; FCSearch *MFCSearch; static UserInput_t empty_UserInput; -bool SCFStatNomadEvaluator::eval_x( NOMAD::Eval_Point & x, const NOMAD::Double & h_max, bool & count_eval) { +bool SCFStatNomadEvaluator::eval_x( Eval_Point & x, const Double & h_max, bool & count_eval) { int n = 0; ExPoint p; - NOMAD::Point y(x.size()); + Point y(x.size()); #ifdef DEBUG MFCSearch->PrintNomadPoint(x,"TRY TO EVAL_X at: "); #endif @@ -260,11 +262,11 @@ bool SCFStatNomadEvaluator::eval_x( NOMAD::Eval_Point & x, const NOMAD::Double return TRUE; // the evaluation succeeded } -bool FCFStatNomadEvaluator::eval_x( NOMAD::Eval_Point & x, const NOMAD::Double & h_max, bool & count_eval) { +bool FCFStatNomadEvaluator::eval_x( Eval_Point & x, const Double & h_max, bool & count_eval) { //printf("%d XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX eval_x %d\n",omp_get_thread_num(),s); int n = 0; ExPoint p; - NOMAD::Point y(x.size()); + Point y(x.size()); if (MFCSearch->FixedDirection.at(PID_ALPHA)) { p.alpha = MFCSearch->FixedPoint.alpha; @@ -472,7 +474,7 @@ INT4 FCSearch::write_pfs_file(REAL8 snr2, INT4 nseg) { }; -INT4 FCSearch::write_cache(NOMAD::Cache *cache, INT4 cnt, INT4 type) { +INT4 FCSearch::write_cache(Cache *cache, INT4 cnt, INT4 type) { FILE *fpCACHE; if ( (fpCACHE = fopen (CacheFile(cnt,type), "wb")) == NULL) { @@ -596,8 +598,8 @@ int main(int argc, char *argv[]) { } - NOMAD::Display out ( std::cout ); - out.precision ( NOMAD::DISPLAY_PRECISION_STD ); + Display out ( std::cout ); + out.precision ( DISPLAY_PRECISION_STD ); try { vector<INT4> searchtype; @@ -605,7 +607,7 @@ int main(int argc, char *argv[]) { searchtype.push_back(MFCSearch->Type()); if (MFCSearch->LTMADS() > 0) { - searchtype.push_back(NOMAD::LT_NP1); + searchtype.push_back(LT_NP1); } @@ -687,8 +689,8 @@ int main(int argc, char *argv[]) { cerr << "\nNOMAD has been interrupted (" << e.what() << ")\n\n"; } - NOMAD::Slave::stop_slaves ( out ); - NOMAD::end(); + Slave::stop_slaves ( out ); + end(); } // Close loop SCStages } // Close SCSEARCH @@ -730,8 +732,8 @@ int main(int argc, char *argv[]) { } - NOMAD::Display out ( std::cout ); - out.precision ( NOMAD::DISPLAY_PRECISION_STD ); + Display out ( std::cout ); + out.precision ( DISPLAY_PRECISION_STD ); try { vector<INT4> searchtype; @@ -739,7 +741,7 @@ int main(int argc, char *argv[]) { searchtype.push_back(MFCSearch->Type()); if (MFCSearch->LTMADS() > 0) { - searchtype.push_back(NOMAD::LT_NP1); + searchtype.push_back(LT_NP1); } MFCSearch->init_SP(FCSEARCH,MFCSearch->singleSegment()); MFCSearch->ComputeSteps(FCSEARCH); @@ -821,8 +823,8 @@ int main(int argc, char *argv[]) { cerr << "\nNOMAD has been interrupted (" << e.what() << ")\n\n"; } - NOMAD::Slave::stop_slaves ( out ); - NOMAD::end(); + Slave::stop_slaves ( out ); + end(); } break; @@ -884,8 +886,8 @@ int main(int argc, char *argv[]) { - NOMAD::Display out ( std::cout ); - out.precision ( NOMAD::DISPLAY_PRECISION_STD ); + Display out ( std::cout ); + out.precision ( DISPLAY_PRECISION_STD ); try { vector<INT4> searchtype; @@ -893,7 +895,7 @@ int main(int argc, char *argv[]) { searchtype.push_back(MFCSearch->Type()); if (MFCSearch->LTMADS() > 0) { - searchtype.push_back(NOMAD::LT_NP1); + searchtype.push_back(LT_NP1); } if (MFCSearch->SignalIsKnown && MFCSearch->uvar->ComputeDTS) { @@ -1022,7 +1024,7 @@ int main(int argc, char *argv[]) { REAL8 tmp_dist = MFCSearch->compute_metric_distance(MFCSearch->SType(),MFCSearch->Semicohmetric(),&sigp,&p); MFCSearch->PrintPoint(MFCSearch->Nseg(),&p,"dst point",MFCSearch->BinarySearch()); printf("%.16e vs %.16e\n",stage_sc_distance_to_the_signal_after_the_search,tmp_dist); - NOMAD::Point y(2); + Point y(2); y[0] = p.freq; y[1] = p.f1dot; MFCSearch->PrintNomadPoint(y,"nmd point"); @@ -1065,8 +1067,8 @@ int main(int argc, char *argv[]) { cerr << "\nNOMAD has been interrupted (" << e.what() << ")\n\n"; } - NOMAD::Slave::stop_slaves ( out ); - NOMAD::end(); + Slave::stop_slaves ( out ); + end(); } @@ -1089,8 +1091,8 @@ int main(int argc, char *argv[]) { MFCSearch->reset_trials(); - NOMAD::Display out ( std::cout ); - out.precision ( NOMAD::DISPLAY_PRECISION_STD ); + Display out ( std::cout ); + out.precision ( DISPLAY_PRECISION_STD ); try { vector<INT4> searchtype; @@ -1098,7 +1100,7 @@ int main(int argc, char *argv[]) { searchtype.push_back(MFCSearch->uvar->FinalType); if (MFCSearch->FinalLTMADS() > 0) { - searchtype.push_back(NOMAD::LT_NP1); + searchtype.push_back(LT_NP1); } if (MFCSearch->SignalIsKnown) { @@ -1217,8 +1219,8 @@ int main(int argc, char *argv[]) { cerr << "\nNOMAD has been interrupted (" << e.what() << ")\n\n"; } - NOMAD::Slave::stop_slaves ( out ); - NOMAD::end(); + Slave::stop_slaves ( out ); + end(); } else { LogPrintf(LOG_NORMAL,"Simulated semicoherent search. Exit now.\n"); @@ -1273,8 +1275,8 @@ INT4 FCSearch::XLALInitUserVars ( int argc, char *argv[] ) uvar->avg2F = 0; - uvar->Type = NOMAD::ORTHO_NP1_QUAD; - uvar->FinalType = NOMAD::ORTHO_NP1_QUAD; + uvar->Type = ORTHO_NP1_QUAD; + uvar->FinalType = ORTHO_NP1_QUAD; uvar->DisplayDegree = 1; uvar->LTMADS = -1; uvar->FinalLTMADS = -1; diff --git a/src/FollowUp/Nomad/libFStatNomad.cpp b/src/FollowUp/Nomad/libFStatNomad.cpp index 3ca41c9b2f4a9d85d6a0880af429c3af152d2573..c47211bdca5af284aa9e73ed8f0dd2ded70a659b 100644 --- a/src/FollowUp/Nomad/libFStatNomad.cpp +++ b/src/FollowUp/Nomad/libFStatNomad.cpp @@ -35,6 +35,7 @@ static FstatCandidate empty_FstatCandidate; static ConfigVariables empty_ConfigVariables; using namespace std; +using namespace NOMAD; REAL8 FStatNomad::DoPredictFStat(INT4 fstatsearchtype, LALSegList *segList) { @@ -2480,7 +2481,7 @@ REAL8 FStatNomad::Randomize(REAL8 s, REAL8 u, REAL8 l, INT4 p) { return l + 0.5*(u-l); } -BOOLEAN FStatNomad::check_point_inside(INT4 t,gsl_matrix *g,NOMAD::Point p) { +BOOLEAN FStatNomad::check_point_inside(INT4 t,gsl_matrix *g,Point p) { gsl_vector *c = gsl_vector_calloc(g->size1); gsl_vector *s = gsl_vector_calloc(g->size2); @@ -2611,7 +2612,7 @@ BOOLEAN FStatNomad::check_point_inside(INT4 t,gsl_matrix *g,NOMAD::Point p) { return point_inside(t,g,c,s); } -gsl_vector* FStatNomad::RandomPointInside(gsl_matrix *g, NOMAD::Point start_point,NOMAD::Point upper_bound, NOMAD::Point lower_bound, NOMAD::Point scale_point, INT4 p) { +gsl_vector* FStatNomad::RandomPointInside(gsl_matrix *g, Point start_point,Point upper_bound, Point lower_bound, Point scale_point, INT4 p) { /* This is an example mapping for the directed case. [0] nu <---> f <---> f [0] @@ -2812,7 +2813,7 @@ gsl_vector* FStatNomad::RandomPointInside(gsl_matrix *g, NOMAD::Point start_poin return c; } -INT4 FStatNomad::DoSearch(NOMAD::Display out,vector<ExPoint> &epv,INT4 fstatsearchtype,vector<INT4> searchtype, INT4 &stcounter, BOOLEAN randomstartpoint) { +INT4 FStatNomad::DoSearch(Display out,vector<ExPoint> &epv,INT4 fstatsearchtype,vector<INT4> searchtype, INT4 &stcounter, BOOLEAN randomstartpoint) { for (INT4 s=0; s<numsteps; s++) { bFirstPoint.push_back(true); @@ -2831,18 +2832,18 @@ INT4 FStatNomad::DoSearch(NOMAD::Display out,vector<ExPoint> &epv,INT4 fstatsear #pragma omp parallel for ordered for (INT4 s=0; s<numsteps; s++) { - NOMAD::Point lower_bound(directions.size()); - NOMAD::Point upper_bound(directions.size()); - NOMAD::Point start_point(directions.size()); - NOMAD::Point scale_point(directions.size()); + Point lower_bound(directions.size()); + Point upper_bound(directions.size()); + Point start_point(directions.size()); + Point scale_point(directions.size()); - NOMAD::Parameters p ( out ); + Parameters p ( out ); p.set_DIMENSION (directions.size()); // number of variables - vector<NOMAD::bb_output_type> bbot (1); // definition of - bbot[0] = NOMAD::OBJ; // output types + vector<bb_output_type> bbot (1); // definition of + bbot[0] = OBJ; // output types p.set_BB_OUTPUT_TYPE ( bbot ); p.set_DISPLAY_STATS ( "time - bbe ( %.16gsol ) obj" ); @@ -2855,7 +2856,7 @@ INT4 FStatNomad::DoSearch(NOMAD::Display out,vector<ExPoint> &epv,INT4 fstatsear p.set_MAX_BB_EVAL(MaxBBEval()); // the algorithm terminates after - p.set_DIRECTION_TYPE((NOMAD::direction_type)searchtype.at(stcounter)); + p.set_DIRECTION_TYPE((direction_type)searchtype.at(stcounter)); p.set_MAX_ITERATIONS(MaxIterations()); @@ -2882,7 +2883,7 @@ INT4 FStatNomad::DoSearch(NOMAD::Display out,vector<ExPoint> &epv,INT4 fstatsear p.set_OPPORTUNISTIC_EVAL(OpportunisticEval()); p.set_OPPORTUNISTIC_LUCKY_EVAL(OpportunisticLuckyEval()); - NOMAD::Point init_mesh(directions.size()); + Point init_mesh(directions.size()); for (INT4 i = 0; i < directions.size(); i++) { init_mesh[i] = directions.at(i).startstep * directions.at(i).scale; } @@ -2953,25 +2954,35 @@ INT4 FStatNomad::DoSearch(NOMAD::Display out,vector<ExPoint> &epv,INT4 fstatsear if (fstatsearchtype == SCSEARCH) { SCFStatNomadEvaluator scev( p, s ); - NOMAD::Mads mads ( p , &scev ); - NOMAD::Evaluator_Control *ec; + Mads mads ( p , &scev ); + Evaluator_Control *ec; ec = &mads.get_evaluator_control(); if (ClearCache()) { - NOMAD::Cache *cache = &ec->get_cache(); + Cache *cache = &ec->get_cache(); cache->clear(); } - NOMAD::Evaluator *pev; + Evaluator *pev; pev = ec->get_evaluator(); FirstPoint(s,true); SkipSearch(s,false); trials_init(); + + // user search: + if (HybridSearch()){ + SCHybridSearch *hybrid_search = new SCHybridSearch ( p ); + hybrid_search->Tiling(tiling); + hybrid_search->NumberPoints(var_HybridSearchPoints); + printf("SET HYBRID SEARCH\n"); + mads.set_user_search ( hybrid_search ); + } + mads.run(); CacheIn(&ec->get_cache()); mads.reset(true,true); // algorithm creation and execution: countstep++; - const NOMAD::Eval_Point *bfp = mads.get_best_feasible(); + const Eval_Point *bfp = mads.get_best_feasible(); ExPoint tExPoint; if (bfp != NULL) { @@ -2990,25 +3001,31 @@ INT4 FStatNomad::DoSearch(NOMAD::Display out,vector<ExPoint> &epv,INT4 fstatsear } else { FCFStatNomadEvaluator fcev( p, s ); - NOMAD::Mads mads ( p , &fcev ); - NOMAD::Evaluator_Control *ec; + Mads mads ( p , &fcev ); + Evaluator_Control *ec; ec = &mads.get_evaluator_control(); if (ClearCache()) { - NOMAD::Cache *cache = &ec->get_cache(); + Cache *cache = &ec->get_cache(); cache->clear(); } - NOMAD::Evaluator *pev; + Evaluator *pev; pev = ec->get_evaluator(); FirstPoint(s,true); SkipSearch(s,false); trials_init(); + + if (HybridSearch()){ + FCHybridSearch hybrid_search ( p ); + mads.set_user_search ( &hybrid_search ); + } + mads.run(); CacheIn(&ec->get_cache()); mads.reset(true,true); // algorithm creation and execution: countstep++; - const NOMAD::Eval_Point *bfp = mads.get_best_feasible(); + const Eval_Point *bfp = mads.get_best_feasible(); ExPoint tExPoint; @@ -3115,12 +3132,12 @@ INT4 DumpSegmentsToFile(const char *fname, LALSegList *segments) { } -void FStatNomad::CacheIn(NOMAD::Cache *cache) { +void FStatNomad::CacheIn(Cache *cache) { LogPrintf(LOG_NORMAL,"CacheIn %d points.\n",cache->size()); LogPrintf(LOG_NORMAL,"Cache size before CacheIn: %d\n",eCache.size()); INT4 n; for (INT4 i = 0; i < cache->size(); i++) { - const NOMAD::Eval_Point *x; + const Eval_Point *x; n = 0; if (i == 0) { x = cache->begin(); @@ -3240,7 +3257,7 @@ FStatNomad::~FStatNomad() {}; void FStatNomad::HybridSearchCreateLattice(gsl_matrix *g, REAL4 m,CHAR *lattice) { - LatticeTiling *tiling = NULL; + tiling = NULL; tiling = XLALCreateLatticeTiling(DimMetric(SType(),USETYPE_SEARCH)); @@ -3313,11 +3330,11 @@ XLALGetDetectorIDs( } /* XLALGetDetectorIDs() */ -void SCHybridSearch::search ( NOMAD::Mads &mads, int &nb_search_pts, bool &stop, NOMAD::stop_type &stop_reason, NOMAD::success_type &success, bool &count_search , const NOMAD::Eval_Point *& new_feas_inc, const NOMAD::Eval_Point *& new_infeas_inc ) { - +void SCHybridSearch::search( Mads & mads, int &nb_search_pts, bool &stop, stop_type &stop_reason, success_type &success, bool &count_search, const Eval_Point *& new_feas_inc, const Eval_Point *& new_infeas_inc){ + printf("ENTER HYBRID SEARCH!\n"); new_feas_inc = new_infeas_inc = NULL; nb_search_pts = 0; - success = NOMAD::UNSUCCESSFUL; + success = UNSUCCESSFUL; count_search = !stop; if ( stop ) { @@ -3325,25 +3342,23 @@ void SCHybridSearch::search ( NOMAD::Mads &mads, int &nb_search_pts, bool &stop, } // initial display: - const NOMAD::Display &out = _p.out(); - NOMAD::dd_type display_degree = out.get_search_dd(); - if ( display_degree == NOMAD::FULL_DISPLAY ){ + const Display &out = _p.out(); + dd_type display_degree = out.get_search_dd(); + if ( display_degree == FULL_DISPLAY ){ std::ostringstream oss; - oss << NOMAD::VNS_SEARCH; - out << std::endl << NOMAD::open_block ( oss.str() ) << std::endl; + oss << USER_SEARCH; + out << std::endl << open_block ( oss.str() ) << std::endl; } - bool opt_only_sgte = _p.get_opt_only_sgte(); - // the barriers: - NOMAD::Barrier & true_barrier = mads.get_true_barrier(); - NOMAD::Barrier & sgte_barrier = mads.get_sgte_barrier(); - const NOMAD::Barrier & active_barrier = mads.get_active_barrier(); + Barrier & true_barrier = mads.get_true_barrier(); + Barrier & sgte_barrier = mads.get_sgte_barrier(); + const Barrier & active_barrier = mads.get_active_barrier(); // point x: - NOMAD::Double best_f; + Double best_f; bool x_feas = true; - const NOMAD::Eval_Point * x = active_barrier.get_best_feasible(); + const Eval_Point * x = active_barrier.get_best_feasible(); if ( x ) { best_f = x->get_f(); } @@ -3353,33 +3368,86 @@ void SCHybridSearch::search ( NOMAD::Mads &mads, int &nb_search_pts, bool &stop, } if ( !x ) { - if ( display_degree == NOMAD::FULL_DISPLAY ) { + if ( display_degree == FULL_DISPLAY ) { out.close_block ( "end of HybridSearch (no incumbent)" ); } return; } // get the signature: - NOMAD::Signature * signature = x->get_signature(); + Signature * signature = x->get_signature(); if ( !signature ) { - if ( display_degree == NOMAD::FULL_DISPLAY ) + if ( display_degree == FULL_DISPLAY ) out.close_block ( "end of HybridSearch (no signature)" ); return; } int n = signature->get_n(); if ( n != x->size() ){ - if ( display_degree == NOMAD::FULL_DISPLAY ) + if ( display_degree == FULL_DISPLAY ) out.close_block ( "end of HybridSearch (incompatible signature)" ); return; } - // build template bank + // template bank + + LatticeTilingLocator *tloc = XLALCreateLatticeTilingLocator(tiling); + + gsl_matrix *points = gsl_matrix_calloc(n,1); + for (INT4 i = 0; i < n; i++){ + gsl_matrix_set(points,i,1,(*x)[i].value()); + } + gsl_matrix **nearest_points = NULL; + UINT8VectorSequence **nearest_indexes = NULL; + + XLALNearestLatticeTilingPoints (tloc, points, nearest_points, nearest_indexes); + printf("EXIT HYBRID SEARCH\n!"); + +Eval_Point * tk = new Eval_Point; + tk->set ( 1 , 1 ); + tk->set_signature ( signature ); + // xk: + double xk = (*new_feas_inc)[0].value(); + if ( xk < 0 ) + return; + + (*tk)[0] = static_cast<int>(ceil(1.0/xk)) * xk - 1.0; + + + // Evaluator_Control: + Evaluator_Control & ev_control = mads.get_evaluator_control(); + + // add the new point to the ordered list of search trial points: + ev_control.add_eval_point ( tk , + _p.out().get_search_dd() , + false , + Double() , + Double() , + Double() , + Double() ); + + nb_search_pts = 1; + + // evaluation: + new_feas_inc = new_infeas_inc = NULL; + ev_control.eval_list_of_points ( _type , + mads.get_true_barrier() , + mads.get_sgte_barrier() , + mads.get_pareto_front() , + stop , + stop_reason , + new_feas_inc , + new_infeas_inc , + success ); + + + } -void FCHybridSearch::search ( NOMAD::Mads &mads, int &nb_search_pts, bool &stop, NOMAD::stop_type &stop_reason, NOMAD::success_type &success, bool &count_search , const NOMAD::Eval_Point *& new_feas_inc, const NOMAD::Eval_Point *& new_infeas_inc ) { +void FCHybridSearch::search ( Mads &mads, int &nb_search_pts, bool &stop, stop_type &stop_reason, success_type &success, bool &count_search , const Eval_Point *& new_feas_inc, const Eval_Point *& new_infeas_inc ) { } + diff --git a/src/FollowUp/Nomad/libFStatNomad.h b/src/FollowUp/Nomad/libFStatNomad.h index 327f43d0a6e0cca0bf1ba392760d03e7c995bc16..dc98e95631afc6e69fb98a7b6852bcd7307a59a5 100644 --- a/src/FollowUp/Nomad/libFStatNomad.h +++ b/src/FollowUp/Nomad/libFStatNomad.h @@ -84,6 +84,8 @@ #include <gsl/gsl_blas.h> #include <Mads.hpp> +using namespace NOMAD; + /* ---------- Macros -------------------- */ #define HSMAX(x,y) ( (x) > (y) ? (x) : (y) ) #define HSMIN(x,y) ( (x) < (y) ? (x) : (y) ) @@ -384,11 +386,22 @@ private: int s; }; -class SCHybridSearch : public NOMAD::Search { +class SCHybridSearch : public Search { public: - SCHybridSearch ( NOMAD::Parameters & p ):NOMAD::Search ( p , NOMAD::USER_SEARCH ) {} + SCHybridSearch ( NOMAD::Parameters & p ):Search ( p , USER_SEARCH ) {} ~SCHybridSearch ( void ) {} - void search ( NOMAD::Mads & mads, int &nb_search_pts, bool &stop, NOMAD::stop_type &stop_reason, NOMAD::success_type &success, bool &count_search, const NOMAD::Eval_Point *& new_feas_inc, const NOMAD::Eval_Point *& new_infeas_inc); + void search ( Mads & mads, int &nb_search_pts, bool &stop, stop_type &stop_reason, success_type &success, bool &count_search, const Eval_Point *& new_feas_inc, const Eval_Point *& new_infeas_inc); + + void Tiling(LatticeTiling *v){ + tiling = v; + } + void NumberPoints(INT4 v){ + np = v; + } +private: + LatticeTiling *tiling; + INT4 np; + }; class FCHybridSearch : public NOMAD::Search { @@ -396,6 +409,11 @@ public: FCHybridSearch ( NOMAD::Parameters & p ):NOMAD::Search ( p , NOMAD::USER_SEARCH ) {} ~FCHybridSearch ( void ) {} void search ( NOMAD::Mads & mads, int &nb_search_pts, bool &stop, NOMAD::stop_type &stop_reason, NOMAD::success_type &success, bool &count_search, const NOMAD::Eval_Point *& new_feas_inc, const NOMAD::Eval_Point *& new_infeas_inc); + void Tiling(LatticeTiling *v){ + tiling = v; + }; +private: + LatticeTiling *tiling; }; class FStatNomad { @@ -2310,6 +2328,8 @@ private: /* template and grid variables */ INT8 *freq_index; + + LatticeTiling *tiling; #ifdef OPENMP