FIrst intents of OpenMP implementation of the likelihood function

git-svn-id: svn://scm.gforge.inria.fr/svnroot/paradiseo@1471 331e1502-861f-0410-8da2-ba01fb791d7f
This commit is contained in:
wcancino 2009-03-04 14:21:19 +00:00
commit 790594f18c
2 changed files with 31 additions and 11 deletions

View file

@ -32,8 +32,26 @@ SET( PHYLOMOEA_SOURCES eigensolver.cpp
utils.cpp
PhyloMOEO_packunpack.cpp)
SET( TESTOMP_SOURCES eigensolver.cpp
likelihoodcalculator.cpp
matrixutils.cpp
probmatrixcontainer.cpp
ProbMatrix.cpp
RandomNr.cpp
Sequences.cpp
SubsModel.cpp
phylotreeIND.cpp
treeIterator.cpp
utils.cpp
likoptimizer.cpp
testomp.cpp
)
ADD_EXECUTABLE( PhyloMOEA ${PHYLOMOEA_SOURCES} )
ADD_EXECUTABLE( testomp ${TESTOMP_SOURCES} )
TARGET_LINK_LIBRARIES(PhyloMOEA gsl gslcblas GTL eo eoutils ga moeo cma peo rmc_mpi xml2)
TARGET_LINK_LIBRARIES(testomp gsl gslcblas GTL eo eoutils ga moeo cma peo rmc_mpi xml2 gomp)
INSTALL( TARGETS PhyloMOEA RUNTIME DESTINATION bin)

View file

@ -449,23 +449,25 @@ double LikelihoodCalculator::sum_partials_a_to_b_notaxon( int pos)
void LikelihoodCalculator::calculate_node_partial( node father, node son, edge edgeaux)
{
register double sum;
double corr_factor;
long index;
double len;
int r,i,j;
//unsigned char l;
ProbMatrix *p;
register int seqlen = tree_ptr->number_of_positions();
#pragma omp parallel for private(p) schedule(dynamic) num_threads(2)
for(int k=0; k<seqlen;k++)
{
index = k*nrates*4;
long index = k*nrates*4;
// accumulatre
Factors[father][k]+=Factors[son][k];
corr_factor = MDBL_MIN;
double corr_factor = MDBL_MIN;
for(int r=0; r<nrates; r++)
{
len = tree_ptr->get_branch_length(edgeaux) * rates_prob[r];
double len = tree_ptr->get_branch_length(edgeaux) * rates_prob[r];
len = (len < BL_MIN ? BL_MIN : len);
ProbMatrix &p = (*probmatrixs)[ len ];
double sum = 0;
#pragma omp critical
p = &((*probmatrixs)[ len ]);
for(int i=0; i < 4; i++)
{
@ -473,12 +475,12 @@ void LikelihoodCalculator::calculate_node_partial( node father, node son, edge e
if(tree_ptr->istaxon( son))
{
unsigned char l=SeqData->pattern_pos( k, tree_ptr->taxon_id( son));
if(SeqData->is_defined( l) ) sum = p.p_ij_t( i, l );
if(SeqData->is_defined( l) ) sum = p->p_ij_t( i, l );
else if(SeqData->is_ambiguous( l))
{
unsigned char *meaning = SeqData->ambiguos_meaning( l);
for(int j=0; j < 4; j++)
sum +=meaning[j]* p.p_ij_t( i, j );
sum +=meaning[j]* p->p_ij_t( i, j );
//sum +=Partials[son][k*4+j]* p.p_ij_t( i, j );
}
else sum = 1;
@ -486,7 +488,7 @@ void LikelihoodCalculator::calculate_node_partial( node father, node son, edge e
else{
for(int j=0; j < 4; j++)
{
sum +=Partials[son][index+ r*4 +j]* p.p_ij_t( i, j );
sum +=Partials[son][index+ r*4 +j]* p->p_ij_t( i, j );
}
}