New openMP development for test cases

git-svn-id: svn://scm.gforge.inria.fr/svnroot/paradiseo@1489 331e1502-861f-0410-8da2-ba01fb791d7f
This commit is contained in:
wcancino 2009-03-09 12:52:05 +00:00
commit 7a8ab5fb50
6 changed files with 118 additions and 253 deletions

View file

@ -25,6 +25,7 @@
#include <likelihoodcalculator.h>
//#include <peo>
#include <iostream>
#include <utils.h>
extern ProbMatrixContainer *probmatrixs_ptr;
@ -58,4 +59,32 @@ private:
ParsimonyCalculator &parseval;
};
class PhyloLikelihoodTimeEval : public moeoEvalFunc < PhyloMOEO >
{
public:
PhyloLikelihoodTimeEval( LikelihoodCalculator &y): likeval(y) { }
void operator () (PhyloMOEO & _sol)
{
ObjectiveVector objVec;
likeval.set_tree(_sol.get_tree());
struct timeval start;
gettimeofday(&start,NULL);
objVec[0] = -likeval.calculate_likelihood();
struct timeval end;
gettimeofday(&end,NULL);
print_elapsed_time_short(&start,&end);
_sol.objectiveVector(objVec);
cout << endl;
}
private:
LikelihoodCalculator &likeval;
};
#endif

View file

@ -241,55 +241,6 @@ double LikelihoodCalculator::calculate_likelihood_exp(edge focus)
}
// double LikelihoodCalculator::calculate_all_likelihoods()
// {
// oldroot = invalid_node;
// graph::edge_iterator it = tree_ptr->TREE.edges_begin();
// graph::edge_iterator it_end = tree_ptr->TREE.edges_end();
// while(it != it_end)
// {
//
// // select a new root
// bfocus = *it;
// a = tree_ptr->istaxon(bfocus.target()) ? bfocus.source() : bfocus.target();
// b = bfocus.opposite(a);
//
// for(int i=0; i< nrates; i++)
// {
// double len = tree_ptr->get_branch_length(bfocus)*rates_prob[i];
// len = len < BL_MIN ? BL_MIN : len;
// ProbMatrix &p = (*probmatrixs)[ len];
// prob[i] = p.p;
// }
//
//
// // make the new tree
// tree_ptr->convert_graph_to_tree(a, NULL);
//
// if( oldroot!=invalid_node)
// {
// update_partials(); //cj( oldroot, a, *it);
// cout << "likelihood ..." << sum_site_liks() <<endl;
// }
// else { //first iteration
// init_partials();
// calculate_partials( a, &b);
// calculate_partials( b, &a);
// cout << "likelihood ..." << sum_site_liks() << endl;
// }
//
// for(int i=0; i< nrates; i++)
// {
// double factor = rates_prob[i];
// double len = tree_ptr->get_branch_length(bfocus);
// len = len < BL_MIN ? BL_MIN : len;
// probmatrixs->change_matrix( len*factor, tree_ptr->get_branch_length(bfocus)*factor);
// }
//
// oldroot = a;
// ++it;
// }
// }
void LikelihoodCalculator::update_partials() //node *oldroot, node newroot, edge newedge)
{
@ -334,7 +285,6 @@ void LikelihoodCalculator::recalculate_cj(node n, node father)
double LikelihoodCalculator::calculate_likelihood()
{
double lik;
// pthread_t threads[2]; /* holds thread info */
@ -360,29 +310,16 @@ double LikelihoodCalculator::calculate_likelihood()
for(int i=0; i<2; i++)
pthread_join(threads[i],NULL);*/
//cout << "calculando partials..." << endl;
struct timeval tempo1, tempo2, result;
gettimeofday(&tempo1, NULL);
struct timeval start;
gettimeofday(&start,NULL);
calculate_partials( a, &b);
calculate_partials( b, &a);
//cout << "somando..." << endl;
// sum all partials
lik = sum_site_liks();
gettimeofday(&tempo2, NULL);
timeval_subtract(&result,&tempo2,&tempo1);
long remainder = result.tv_sec % 3600;
long hours = (result.tv_sec - remainder)/3600;
long seconds = remainder % 60;
long minutes = (remainder - seconds) / 60;
cout << "Execution time : ";
cout.width(3);
cout.fill(' ');
cout << hours << ":";
cout.width(2);
cout.fill('0');
cout << minutes << ":";
cout.width(2);
cout.fill('0');
cout << seconds << "." << result.tv_usec << "(" << result.tv_sec << ")" << endl;
struct timeval end;
gettimeofday(&end,NULL);
print_elapsed_time_short(&start, &end);
return lik;
}
@ -402,12 +339,10 @@ double LikelihoodCalculator::sum_site_liks( )
prob[i] = p->p;
}
//#pragma omp parallel for private(factor_correct) schedule(dynamic) num_threads(2) reduction(+:lik)
for(int i=0; i < seqlen; i++)
{
factor_correct = Factors[a][i] + Factors[b][i] ;
site_liks[i] = sum_partials(i);
//#pragma omp critical
lik += ( log(site_liks[i]) + factor_correct)* SeqData->pattern_count(i);
}
return lik;
@ -487,7 +422,7 @@ void LikelihoodCalculator::calculate_node_partial( node father, node son, edge e
int r,i,j;
//unsigned char l;
register int seqlen = tree_ptr->number_of_positions();
#pragma omp parallel for
#pragma omp parallel for
for(int k=0; k<seqlen;k++)
{
long index = k*nrates*4;
@ -498,7 +433,6 @@ void LikelihoodCalculator::calculate_node_partial( node father, node son, edge e
for(int r=0; r<nrates; r++)
{
double sum = 0;
//#pragma omp critical
ProbMatrix *p = edgeprobmatrix[edgeaux][r];
for(int i=0; i < 4; i++)
@ -594,7 +528,6 @@ double LikelihoodCalculator::calculate_likelihood_omp()
init_partials();
struct timeval tempo1, tempo2, result;
int seqlen = tree_ptr->number_of_positions();
@ -608,7 +541,8 @@ double LikelihoodCalculator::calculate_likelihood_omp()
}
gettimeofday(&tempo1, NULL);
struct timeval start;
gettimeofday(&start,NULL);
#pragma omp parallel for reduction(+:lik)
for(int i=0; i< seqlen; i++)
@ -619,22 +553,10 @@ double LikelihoodCalculator::calculate_likelihood_omp()
// sum all partials
lik += sum_site_liks_omp(i);
}
gettimeofday(&tempo2, NULL);
timeval_subtract(&result,&tempo2,&tempo1);
long remainder = result.tv_sec % 3600;
long hours = (result.tv_sec - remainder)/3600;
long seconds = remainder % 60;
long minutes = (remainder - seconds) / 60;
cout << "Execution time : ";
cout.width(3);
cout.fill(' ');
cout << hours << ":";
cout.width(2);
cout.fill('0');
cout << minutes << ":";
cout.width(2);
cout.fill('0');
cout << seconds << "." << result.tv_usec << "(" << result.tv_sec << ")" << endl;
struct timeval end;
gettimeofday(&end,NULL);
cout << " !! ";
print_elapsed_time_short(&start,&end);
return lik;
}
@ -645,7 +567,6 @@ double LikelihoodCalculator::sum_site_liks_omp( int pos )
register double lik = 0;
register double factor_correct;
//#pragma omp parallel for private(factor_correct) schedule(dynamic) num_threads(2) reduction(+:lik)
factor_correct = Factors[a][pos] + Factors[b][pos] ;
site_liks[pos] = sum_partials(pos);
lik = ( log(site_liks[pos]) + factor_correct)* SeqData->pattern_count(pos);
@ -668,7 +589,6 @@ void LikelihoodCalculator::calculate_node_partial_omp( node father, node son, ed
{
register double sum;
//unsigned char l;
//#pragma omp parallel for
long index = pos*nrates*4;
// accumulatre
Factors[father][pos]+=Factors[son][pos];
@ -677,7 +597,6 @@ void LikelihoodCalculator::calculate_node_partial_omp( node father, node son, ed
for(int r=0; r<nrates; r++)
{
double sum = 0;
//#pragma omp critical
ProbMatrix *p = edgeprobmatrix[edgeaux][r];
for(int i=0; i < 4; i++)

View file

@ -89,8 +89,12 @@ private:
long total_pos = SeqData->pattern_count();
int ntaxons = SeqData->num_seqs();
//int nedges = tree_ptr->TREE.number_of_edges();
part_memory_internal = new double[ (ntaxons-2) * nrates * total_pos * 4 ];
part_memory_factors = new double [ (2*ntaxons-2) * total_pos];
unsigned long tamanho = (ntaxons-2) * nrates * total_pos * 4;
cout << "tamanho de alocacao:" << tamanho << endl;
part_memory_internal = new double[ tamanho ];
tamanho = (2*ntaxons-2) * total_pos;
cout << "tamanho de alocacao:" << tamanho << endl;
part_memory_factors = new double [ tamanho ];
part_memory_probmatrix_ptr = new ProbMatrix* [ (2*ntaxons-3) * nrates ] ;
site_liks = new double[total_pos];
cout << "allocating done..." << endl;

View file

@ -47,180 +47,55 @@ int main(int argc, char *argv[])
{
// measures execution time
struct timeval tempo1, tempo2, result;
eoParser parser(argc, argv);
datafile = parser.createParam(string(), "data", "Datafile", 'd',"Param").value();
int nthreads = parser.createParam(omp_get_max_threads(), "nthreads", "Numthreads", 't',"Param").value();
int ntrees = parser.createParam(20, "ntrees", "NumTrees", 'n',"Param").value();
int nexp = parser.createParam(1, "nexps", "NumExps", 'e',"Param").value();
cout << "\n\nReading Sequence Datafile...";
// Sequences seq("/home/wcancino/experimentos/PhyloMOEA_0.2/omp_tests/datasets/TEST.50_5000");
// Sequences seq2("/home/wcancino/experimentos/PhyloMOEA_0.2/omp_tests/datasets/TEST.50_50000");
// Sequences seq3("/home/wcancino/experimentos/PhyloMOEA_0.2/omp_tests/datasets/TEST.50_500000");
/* Sequences seq4("/home/wcancino/experimentos/PhyloMOEA_0.2/omp_tests/datasets/TEST.250_5000");
Sequences seq5("/home/wcancino/experimentos/PhyloMOEA_0.2/omp_tests/datasets/TEST.250_50000");*/
Sequences seq6("/home/wcancino/experimentos/PhyloMOEA_0.2/omp_tests/datasets/TEST.250_500000");
Sequences seq6(datafile.c_str());
// Sequences seq7("/home/wcancino/experimentos/PhyloMOEA_0.2/omp_tests/datasets/TEST.500_5000");
cout << " done.\n";
// calculate datafile
cout << "calculating pattersn..." << endl;
/* seq.calculate_patterns();
seq.calculate_frequences();
seq2.calculate_patterns();
seq2.calculate_frequences();
seq3.calculate_patterns();
seq3.calculate_frequences();*/
/* seq4.calculate_patterns();
seq4.calculate_frequences();
seq5.calculate_patterns();
seq5.calculate_frequences();*/
seq6.calculate_patterns();
seq6.calculate_frequences();
// seq7.calculate_patterns();
// seq7.calculate_frequences();
ostringstream os;
os << datafile << "_results_" << nthreads << ".txt";
ofstream of(os.str().c_str());
gsl_rng *rn2 = gsl_rng_alloc(gsl_rng_default);
RandomNr *rn = new RandomNr(time(NULL));
/* phylotreeIND templatetree( rn, seq, rn2);
phylotreeIND templatetree2( rn, seq2, rn2);
phylotreeIND templatetree3( rn, seq3, rn2);*/
/* phylotreeIND templatetree4( rn, seq4, rn2);
phylotreeIND templatetree5( rn, seq5, rn2);*/
phylotreeIND templatetree6( rn, seq6, rn2);
// phylotreeIND templatetree7( rn, seq7, rn2);
/* SubstModel modelHKY( seq, SubstModel::HKY85);
SubstModel modelHKY2( seq2, SubstModel::HKY85);
SubstModel modelHKY3( seq3, SubstModel::HKY85);*/
/* SubstModel modelHKY4( seq4, SubstModel::HKY85);
SubstModel modelHKY5( seq5, SubstModel::HKY85);*/
SubstModel modelHKY6( seq6, SubstModel::HKY85);
// SubstModel modelHKY7( seq7, SubstModel::HKY85);
/* modelHKY.init();
modelHKY2.init();
modelHKY3.init();*/
/* modelHKY4.init();
modelHKY5.init();*/
modelHKY6.init();
// modelHKY7.init();
//modelHKY.set_kappa(3.890);
/* ProbMatrixContainer probmatrixs(modelHKY);
ProbMatrixContainer probmatrixs2(modelHKY2);
ProbMatrixContainer probmatrixs3(modelHKY3);*/
/* ProbMatrixContainer probmatrixs4(modelHKY4);
ProbMatrixContainer probmatrixs5(modelHKY5);*/
ProbMatrixContainer probmatrixs6(modelHKY6);
// ProbMatrixContainer probmatrixs7(modelHKY7);
//probmatrixs_ptr = &probmatrixs;
/* LikelihoodCalculator lik_calc(templatetree, modelHKY, probmatrixs, 4);
LikelihoodCalculator lik_calc2(templatetree2, modelHKY2, probmatrixs2, 4);
LikelihoodCalculator lik_calc3(templatetree3, modelHKY3, probmatrixs3, 4);*/
/* LikelihoodCalculator lik_calc4(templatetree4, modelHKY4, probmatrixs4, 4);
LikelihoodCalculator lik_calc5(templatetree5, modelHKY5, probmatrixs5, 4); */
LikelihoodCalculator lik_calc6(templatetree6, modelHKY6, probmatrixs6, 4);
// LikelihoodCalculator lik_calc7(templatetree7, modelHKY7, probmatrixs7, 4);
//lik_calc.set_alpha(0.950);
/* modelHKY.init();
modelHKY2.init();
modelHKY3.init();*/
/* modelHKY4.init();
modelHKY5.init();*/
modelHKY6.init();
// modelHKY7.init();
/* Phyloraninit initializer(templatetree);
Phyloraninit initializer2(templatetree2);
Phyloraninit initializer3(templatetree3);*/
/* Phyloraninit initializer4(templatetree4);
Phyloraninit initializer5(templatetree5);*/
Phyloraninit initializer6(templatetree6);
// Phyloraninit initializer7(templatetree7);
/* eoPop<PhyloMOEO> population(1, initializer);
eoPop<PhyloMOEO> population2(1, initializer2);
eoPop<PhyloMOEO> population3(1, initializer3);*/
/* eoPop<PhyloMOEO> population4(1, initializer4);
eoPop<PhyloMOEO> population5(1, initializer5);*/
eoPop<PhyloMOEO> population6(1, initializer6);
// eoPop<PhyloMOEO> population7(1, initializer7);
cout << "Reading trees..." << endl;
/* readtrees("/home/wcancino/experimentos/PhyloMOEA_0.2/omp_tests/datasets/start.TEST.50_5000", population);
readtrees("/home/wcancino/experimentos/PhyloMOEA_0.2/omp_tests/datasets/start.TEST.50_50000", population2);
readtrees("/home/wcancino/experimentos/PhyloMOEA_0.2/omp_tests/datasets/start.TEST.50_500000", population3);*/
/* readtrees("/home/wcancino/experimentos/PhyloMOEA_0.2/omp_tests/datasets/start.TEST.250_5000", population4);
readtrees("/home/wcancino/experimentos/PhyloMOEA_0.2/omp_tests/datasets/start.TEST.250_50000", population5);*/
readtrees("/home/wcancino/experimentos/PhyloMOEA_0.2/omp_tests/datasets/start.TEST.250_500000", population6);
// readtrees("/home/wcancino/experimentos/PhyloMOEA_0.2/omp_tests/datasets/start.TEST.500_5000", population7);
gettimeofday(&tempo1, NULL);
cout << " Number of processors available:" << omp_get_num_procs() << " MAX Number of threads " << omp_get_max_threads() << endl;
cout.precision(15);
// lik_calc.set_tree(population[0].get_tree());
// cout << lik_calc.calculate_likelihood() << endl;
// cout << lik_calc.calculate_likelihood_omp() << endl;
//
// lik_calc2.set_tree(population2[0].get_tree());
// cout << lik_calc2.calculate_likelihood() << endl;
// cout << lik_calc2.calculate_likelihood_omp() << endl;
//
// lik_calc3.set_tree(population3[0].get_tree());
// cout << lik_calc3.calculate_likelihood() << endl;
// cout << lik_calc3.calculate_likelihood_omp() << endl;
// lik_calc4.set_tree(population4[0].get_tree());
// cout << lik_calc4.calculate_likelihood() << endl;
// cout << lik_calc4.calculate_likelihood_omp() << endl;
//
// lik_calc5.set_tree(population5[0].get_tree());
// cout << lik_calc5.calculate_likelihood() << endl;
// cout << lik_calc5.calculate_likelihood_omp() << endl;
lik_calc6.set_tree(population6[0].get_tree());
cout << lik_calc6.calculate_likelihood() << endl;
cout << lik_calc6.calculate_likelihood_omp() << endl;
//
// lik_calc7.set_tree(population7[0].get_tree());
// cout << lik_calc7.calculate_likelihood() << endl;
// cout << lik_calc7.calculate_likelihood_omp() << endl;
/*
for(int i=0; i < population.size(); i++)
omp_set_num_threads(nthreads);
for(int i=0; i < nexp; i++)
{
phylotreeIND templatetree6( rn, seq6, rn2);
SubstModel modelHKY6( seq6, SubstModel::HKY85);
modelHKY6.init();
ProbMatrixContainer probmatrixs6(modelHKY6);
LikelihoodCalculator lik_calc6(templatetree6, modelHKY6, probmatrixs6, 4);
modelHKY6.init();
Phyloraninit initializer6(templatetree6);
eoPop<PhyloMOEO> population6(ntrees, initializer6);
cout.precision(15);
PhyloLikelihoodTimeEval eval( lik_calc6 );
cout << " Number of processors available:" << omp_get_num_procs() << " MAX Number of threads " << omp_get_max_threads() << endl;
gettimeofday(&tempo1, NULL);
lik_calc.set_tree(population[i].get_tree());
// gettimeofday(&tempo1, NULL);
cout << lik_calc.calculate_likelihood_omp() << endl;
/* gettimeofday(&tempo2, NULL);
timeval_subtract(&result,&tempo2,&tempo1);
long remainder = result.tv_sec % 3600;
long hours = (result.tv_sec - remainder)/3600;
long seconds = remainder % 60;
long minutes = (remainder - seconds) / 60;
cout << "Execution time : ";
cout.width(3);
cout.fill(' ');
cout << hours << ":";
cout.width(2);
cout.fill('0');
cout << minutes << ":";
cout.width(2);
cout.fill('0');
cout << seconds << "." << result.tv_usec << "(" << result.tv_sec << ")" << endl;*/
// }
apply<PhyloMOEO> (eval, population6);
gettimeofday(&tempo2, NULL);
cout << "\n"; print_elapsed_time(&tempo1,&tempo2);
print_elapsed_time_short(&tempo1,&tempo2,of);
of << endl;
}
of.close();
gsl_rng_free(rn2);
// delete probmatrixs;
delete rn;

View file

@ -61,6 +61,41 @@ int timeval_subtract (struct timeval *result, struct timeval *x, struct timeva
return x->tv_sec < y->tv_sec;
}
void print_elapsed_time(struct timeval *x, struct timeval *y)
{
struct timeval result;
timeval_subtract(&result,y,x);
long remainder = result.tv_sec % 3600;
long hours = (result.tv_sec - remainder)/3600;
long seconds = remainder % 60;
long minutes = (remainder - seconds) / 60;
cout << "Execution time : ";
cout.width(3);
cout.fill(' ');
cout << hours << ":";
cout.width(2);
cout.fill('0');
cout << minutes << ":";
cout.width(2);
cout.fill('0');
cout << seconds << "." << result.tv_usec << "(" << result.tv_sec << ")" << endl;
}
void print_elapsed_time_short(struct timeval *x, struct timeval *y, ostream &os)
{
struct timeval result;
timeval_subtract(&result,y,x);
os << " " << result.tv_sec << "." << result.tv_usec;
}
void print_cpu_time(clock_t start, clock_t end)
{
double cpu_time_used = ((double) (end - start)) / CLOCKS_PER_SEC;
cout << " " << cpu_time_used << " ";
}
void welcome_message()
{
cout << "\nPhyloMOEA, a program for multi-criteria phylogenetic inference\n";

View file

@ -32,5 +32,8 @@ void optimize_solutions( eoPop<PhyloMOEO> &);
void optimize_solution( PhyloMOEO &);
void readtrees(const char *, eoPop<PhyloMOEO> &);
int timeval_subtract (struct timeval *, struct timeval *, struct timeval *);
void print_cpu_time(clock_t ,clock_t);
void print_elapsed_time(struct timeval *, struct timeval *);
void print_elapsed_time_short(struct timeval *, struct timeval *, ostream &os=cout);
#endif