[Debian-med-packaging] MP patch for PhyML
Andreas Tille
andreas at an3as.eu
Mon Mar 14 18:02:28 UTC 2011
Hi,
I noticed that there is a new version of PhyML (20110304) available.
Does this regard the MP feature yet? (I have not checked the source
myself - that's why I'm asking.)
Kind regards
Andreas.
On Tue, Feb 15, 2011 at 11:38:16PM +0100, Diego Darriba López wrote:
> Dear Andreas,
>
>
> My OpenMP implementation perform the MP just using pragma directives, so if the -fopenmp flag is not set, the compiler will just ignore these lines and it will perform a single processor build. I checked both builds with gcc. We would be in a different situation if I had used OpenMP API calls. In fact, the #include<omp.h> I included in my patch it is not necessary at all (sorry). This way, the MP or sequential build just depends on the use of the -fopenmp flag.
>
>
> Keep in mind that this flag differs between compilares (e.g, -openmp for icc, -mp for openUH, etc.) if you want to bring support for different compilers in the autoconf file.
>
>
> Kind regards,
>
>
> Diego.
>
>
> De: "Andreas Tille" <andreas at an3as.eu>
> Para: "Diego Darriba" <ddarriba at udc.es>
> CC: "s guindon" <s.guindon at auckland.ac.nz>, "Stephane Guindon" <guindon at stat.auckland.ac.nz>, "Olivier Gascuel" <gascuel at lirmm.fr>, "Jean-Francois Dufayard" <jeanfrancois.dufayard at gmail.com>, "David Posada" <dposada at uvigo.es>, "Debian Med Packaging Team" <debian-med-packaging at lists.alioth.debian.org>
> Enviados: Martes, 15 de Febrero 2011 21:25:51
> Asunto: Re: MP patch for PhyML
>
> Dear Diego,
>
> I admit I have no experience with MP builds but I wonder whether you
> should not include the changes into #ifdef USEMP / #else / #endif
> statements to enable both - single processor and MP builds. Or does the
> patch work with both and it is just the build option which controls
> everything?
>
> Kind regards
>
> Andreas.
>
> On Tue, Feb 15, 2011 at 06:58:19PM +0100, Diego Darriba wrote:
> > Dear Stephan et al,
> >
> > I send you the patch for the current phyml trunk. I checked the results
> > with several input data, and I think I solved correctly some
> > dependencies in the parallelized loops from my computer science point of
> > view. In the configure.ac it is just necessary to include the *-fopenmp*
> > flag for compiling (if using gcc compiler). I would modify it, but I
> > prefer you to do so. That way you can create a conditional Makefile, or
> > include the flag as a commentary, or whatever.
> >
> > If the compilation flag is not set, the application will compile with no
> > error, but it will run as before (i.e, sequentially).
> >
> > Please, let me know if you have any doubt or any problem.
> >
> > Kind Regards,
> >
> > Diego.
> >
> > El 15/02/11 02:51, Stephane Guindon escribió:
> >> Dear Andreas,
> >>
> >> Thank you very much for contributing to PhyML.
> >> Yes, it makes perfect sense to try and include this functionality
> >> into the latest release of the program. Please feel free to take
> >> a look at the most recent version of the sources on our svn server
> >> (http://code.google.com/p/phyml/source/browse/#svn%2Ftrunk%2Fsrc)
> >> and let me know how to amend the configure.ac, *.c and *.h file
> >> in an appropriate manner.
> >>
> >> Regards,
> >>
> >> -Stephane-
> >>
> >> On 15/02/11 03:26, Andreas Tille wrote:
> >>> Hi PhyML authors,
> >>>
> >>> when discssing with Diego about a program which uses PhyML he told me
> >>> about his MP patch which would enable making use of hybrid memory when
> >>> running PhyML on clusters. I was considering to port this patch to the
> >>> Debian packaged version where I intended to support two versions of
> >>> the package: one conventional phylip and one phylip-mp.
> >>>
> >>> When inspecting the patch I noticed that it was done against an older
> >>> version of PhyML (20100123). I would consider it a reasonable idea if
> >>> this patch would be included into the upstream version by just wrapping
> >>> it into "#ifdef"s properly to enable switching between both versions by
> >>> simply providing different options to make. If the patch would be taken
> >>> over into the official PhyML code this would simplify our work as
> >>> distributor (as well as Diego's work) and it might be of extra value for
> >>> other users.
> >>>
> >>> What do you think about this idea?
> >>>
> >>> Kind regards
> >>>
> >>> Andreas.
> >>>
> >
>
> > Index: src/lk.c
> > ===================================================================
> > --- src/lk.c (revisi??n: 582)
> > +++ src/lk.c (copia de trabajo)
> > @@ -11,8 +11,8 @@
> > */
> >
> > #include "lk.h"
> > +#include <omp.h>
> >
> > -
> > /* int LIM_SCALE; */
> > /* phydbl LIM_SCALE_VAL; */
> > /* phydbl BIG; */
> > @@ -320,8 +320,9 @@
> >
> > phydbl Lk(t_tree *tree)
> > {
> > - int br;
> > + int br, site;
> > int n_patterns;
> > + phydbl c_lnL;
> >
> > tree->old_lnL = tree->c_lnL;
> >
> > @@ -349,10 +350,17 @@
> >
> > tree->c_lnL = .0;
> > tree->sum_min_sum_scale = .0;
> > - For(tree->curr_site,n_patterns)
> > +
> > + c_lnL = .0;
> > + #pragma omp parallel for reduction(+:c_lnL)
> > + For(site,n_patterns)
> > {
> > - if(tree->data->wght[tree->curr_site] > SMALL) Lk_Core(tree->noeud[0]->b[0],tree);
> > + if(tree->data->wght[site] > SMALL) {
> > + Lk_Core(tree->noeud[0]->b[0],tree,site);
> > + c_lnL += tree->c_lnL_sorted[site];
> > + }
> > }
> > + tree->c_lnL = c_lnL;
> >
> > /* Qksort(tree->c_lnL_sorted,NULL,0,n_patterns-1); */
> >
> > @@ -372,7 +380,8 @@
> >
> > phydbl Lk_At_Given_Edge(t_edge *b_fcus, t_tree *tree)
> > {
> > - int n_patterns;
> > + int n_patterns,site;
> > + phydbl c_lnL;
> >
> > n_patterns = tree->n_pattern;
> >
> > @@ -396,10 +405,17 @@
> >
> > tree->c_lnL = .0;
> > tree->sum_min_sum_scale = .0;
> > - For(tree->curr_site,n_patterns)
> > +
> > + c_lnL = .0;
> > + #pragma omp parallel for reduction(+:c_lnL)
> > + For(site,n_patterns)
> > {
> > - if(tree->data->wght[tree->curr_site] > SMALL) Lk_Core(b_fcus,tree);
> > + if(tree->data->wght[site] > SMALL) {
> > + Lk_Core(b_fcus,tree,site);
> > + c_lnL += tree->c_lnL_sorted[site];
> > + }
> > }
> > + tree->c_lnL = c_lnL;
> >
> > /* Qksort(tree->c_lnL_sorted,NULL,0,n_patterns-1); */
> >
> > @@ -422,7 +438,7 @@
> >
> > #ifndef USE_OLD_LK
> >
> > -phydbl Lk_Core(t_edge *b, t_tree *tree)
> > +phydbl Lk_Core(t_edge *b, t_tree *tree, int site)
> > {
> > phydbl log_site_lk;
> > phydbl site_lk_cat, site_lk;
> > @@ -431,7 +447,7 @@
> > phydbl max_sum_scale;
> > phydbl sum;
> > int ambiguity_check,state;
> > - int catg,ns,k,l,site;
> > + int catg,ns,k,l;
> > int dim1,dim2,dim3;
> > int *sum_scale_left_cat,*sum_scale_rght_cat;
> > phydbl multiplier;
> > @@ -439,6 +455,7 @@
> > phydbl tmp;
> > phydbl logbig;
> > phydbl inv_site_lk;
> > + phydbl tree_site_lk_cat[tree->mod->n_catg];
> >
> > logbig = LOG((phydbl)BIG);
> >
> > @@ -451,7 +468,6 @@
> > site_lk_cat = .0;
> > ambiguity_check = -1;
> > state = -1;
> > - site = tree->curr_site;
> > ns = tree->mod->ns;
> >
> > if((b->rght->tax) && (!tree->mod->s_opt->greedy))
> > @@ -537,7 +553,7 @@
> > }
> > }
> > }
> > - tree->site_lk_cat[catg] = site_lk_cat;
> > + tree_site_lk_cat[catg] = site_lk_cat;
> > }
> >
> > max_sum_scale = (phydbl)BIG;
> > @@ -562,7 +578,7 @@
> > Warn_And_Exit("\n");
> > }
> >
> > - tmp = sum + (logbig - LOG(tree->site_lk_cat[catg]))/(phydbl)LOG2;
> > + tmp = sum + (logbig - LOG(tree_site_lk_cat[catg]))/(phydbl)LOG2;
> > if(tmp < max_sum_scale) max_sum_scale = tmp; /* min of the maxs */
> > }
> >
> > @@ -575,7 +591,7 @@
> > {
> > exponent = -(sum_scale_left_cat[catg]+sum_scale_rght_cat[catg])+fact_sum_scale;
> >
> > - site_lk_cat = tree->site_lk_cat[catg];
> > + site_lk_cat = tree_site_lk_cat[catg];
> >
> > if(exponent > 0)
> > {
> > @@ -605,12 +621,12 @@
> > PhyML_Printf("\n+ site=%4d cat=%4d site_lk_cat=%G sum_scale=%d max=%G fact=%d expo=%d dbl=%G",
> > tree->curr_site,
> > catg,
> > - tree->site_lk_cat[catg],
> > + tree_site_lk_cat[catg],
> > sum_scale_left_cat[catg]+sum_scale_rght_cat[catg],
> > max_sum_scale,
> > fact_sum_scale,
> > -(sum_scale_left_cat[catg]+sum_scale_rght_cat[catg])+fact_sum_scale,
> > - (double)tree->site_lk_cat[catg] * pow(2.,-(sum_scale_left_cat[catg]+sum_scale_rght_cat[catg])+fact_sum_scale));
> > + (double)tree_site_lk_cat[catg] * pow(2.,-(sum_scale_left_cat[catg]+sum_scale_rght_cat[catg])+fact_sum_scale));
> >
> > site_lk_cat = BIG / 10;
> > }
> > @@ -623,12 +639,12 @@
> > PhyML_Printf("\n- site=%4d cat=%4d site_lk_cat=%G sum_scale=%d max=%G fact=%d expo=%d dbl=%G",
> > tree->curr_site,
> > catg,
> > - tree->site_lk_cat[catg],
> > + tree_site_lk_cat[catg],
> > sum_scale_left_cat[catg]+sum_scale_rght_cat[catg],
> > max_sum_scale,
> > fact_sum_scale,
> > -(sum_scale_left_cat[catg]+sum_scale_rght_cat[catg])+fact_sum_scale,
> > - (double)tree->site_lk_cat[catg] * pow(2.,-(sum_scale_left_cat[catg]+sum_scale_rght_cat[catg])+fact_sum_scale));
> > + (double)tree_site_lk_cat[catg] * pow(2.,-(sum_scale_left_cat[catg]+sum_scale_rght_cat[catg])+fact_sum_scale));
> >
> > Exit("\n");
> > }
> > @@ -636,13 +652,13 @@
> > site_lk_cat = .0;
> > }
> >
> > - tree->site_lk_cat[catg] = site_lk_cat;
> > + tree_site_lk_cat[catg] = site_lk_cat;
> > }
> >
> > site_lk = .0;
> > For(catg,tree->mod->n_catg)
> > {
> > - site_lk += tree->site_lk_cat[catg] * tree->mod->gamma_r_proba[catg];
> > + site_lk += tree_site_lk_cat[catg] * tree->mod->gamma_r_proba[catg];
> > }
> >
> > /* The substitution model does include invariable sites */
> > @@ -677,7 +693,7 @@
> >
> > log_site_lk = LOG(site_lk) - (phydbl)LOG2 * fact_sum_scale;
> >
> > - For(catg,tree->mod->n_catg) tree->log_site_lk_cat[catg][site] = LOG(tree->site_lk_cat[catg]) - (phydbl)LOG2 * fact_sum_scale;
> > + For(catg,tree->mod->n_catg) tree->log_site_lk_cat[catg][site] = LOG(tree_site_lk_cat[catg]) - (phydbl)LOG2 * fact_sum_scale;
> >
> > if(isinf(log_site_lk) || isnan(log_site_lk))
> > {
> > @@ -694,7 +710,6 @@
> > /* Multiply log likelihood by the number of times this site pattern is found in the data */
> > tree->c_lnL_sorted[site] = tree->data->wght[site]*log_site_lk;
> >
> > - tree->c_lnL += tree->data->wght[site]*log_site_lk;
> > /* tree->sum_min_sum_scale += (int)tree->data->wght[site]*min_sum_scale; */
> >
> > return log_site_lk;
> > @@ -1118,6 +1133,7 @@
> > Pij2 = d->b[dir2]->Pij_rr;
> >
> > /* For every site in the alignment */
> > + #pragma omp parallel for private(state_v1,state_v2,ambiguity_check_v1,ambiguity_check_v2,catg,i,p1_lk1,p2_lk2)
> > For(site,n_patterns)
> > {
> > state_v1 = state_v2 = -1;
> > @@ -1399,6 +1415,7 @@
> > Pij2 = d->b[dir2]->Pij_rr;
> >
> > /* For every site in the alignment */
> > + #pragma omp parallel for private(state_v1,state_v2,ambiguity_check_v1,ambiguity_check_v2,catg,i,p1_lk1,p2_lk2)
> > For(site,n_patterns)
> > {
> > state_v1 = state_v2 = -1;
> > @@ -2234,13 +2251,13 @@
> > sides of t_edge *b are up-to-date. Calculate the log-likelihood at one site.
> > */
> > #ifdef USE_OLD_LK
> > -phydbl Lk_Core(t_edge *b, t_tree *tree)
> > +phydbl Lk_Core(t_edge *b, t_tree *tree, int site)
> > {
> > phydbl log_site_lk, site_lk, site_lk_cat;
> > phydbl scale_left, scale_rght;
> > phydbl sum;
> > int ambiguity_check,state;
> > - int catg,ns,k,l,site;
> > + int catg,ns,k,l;
> > int dim1,dim2,dim3;
> >
> >
> > @@ -2253,7 +2270,6 @@
> > site_lk_cat = .0;
> > ambiguity_check = -1;
> > state = -1;
> > - site = tree->curr_site;
> > ns = tree->mod->ns;
> >
> > tree->old_site_lk[site] = tree->cur_site_lk[site];
> > @@ -2383,7 +2399,7 @@
> > tree->cur_site_lk[site] = log_site_lk;
> > /* Multiply log likelihood by the number of times this site pattern is found oin the data */
> > tree->c_lnL_sorted[site] = tree->data->wght[site]*log_site_lk;
> > - tree->c_lnL += tree->c_lnL_sorted[site];
> > +
> > return log_site_lk;
> > }
> >
> > @@ -2492,6 +2508,7 @@
> > Pij2 = d->b[dir2]->Pij_rr;
> >
> > /* For every site in the alignment */
> > + #pragma omp parallel for private(scale_v1,scale_v2,state_v1,state_v2,ambiguity_check_v1,ambiguity_check_v2,catg,i,j,p1_lk1,p2_lk2,max_p_lk)
> > For(site,n_patterns)
> > {
> > /* Current scaling values at that site */
> > Index: src/lk.h
> > ===================================================================
> > --- src/lk.h (revisi??n: 582)
> > +++ src/lk.h (copia de trabajo)
> > @@ -51,7 +51,7 @@
> > phydbl Update_Lk_At_Given_Edge(t_edge *b_fcus, t_tree *tree);
> > void Update_P_Lk_Greedy(t_tree *tree, t_edge *b_fcus, t_node *n);
> > void Get_All_Partial_Lk_Scale_Greedy(t_tree *tree, t_edge *b_fcus, t_node *a, t_node *d);
> > -phydbl Lk_Core(t_edge *b, t_tree *tree);
> > +phydbl Lk_Core(t_edge *b, t_tree *tree, int site);
> > phydbl Lk_Triplet(t_node *a, t_node *d, t_tree *tree);
> > void Print_Lk_Given_Edge_Recurr(t_node *a, t_node *d, t_edge *b, t_tree *tree);
> > phydbl *Post_Prob_Rates_At_Given_Edge(t_edge *b, phydbl *post_prob, t_tree *tree);
>
>
> --
> http://fam-tille.de
--
http://fam-tille.de
More information about the Debian-med-packaging
mailing list