Fixed compilation bugs

src/init.c

@@ -3489,6 +3489,9 @@ void PHYREX_Init_Migrep_Mod(t_phyrex_mod *t, int n_dim, phydbl min_lat, phydbl m
 void PHYREX_Set_Default_Migrep_Mod(int n_otu, t_phyrex_mod *t)
 {
   for(int i=0;i<2*n_otu-2;++i) t->sigsq_scale[i] = 1.0;
+
+  t->sigsq_scale_min = 0.01;
+  t->sigsq_scale_max = 100.;
 
   t->id          = -1;
   t->safe_phyrex = NO;

src/lk.c

@@ -602,8 +602,7 @@ phydbl Lk(t_edge *b, t_tree *tree)
 
   p_lk_left = b->p_lk_left;
   p_lk_rght = b->rght->tax ? b->p_lk_tip_r : b->p_lk_rght;
-
-
+
   for(site=0;site<npatterns;++site)
     {
       ambiguity_check = -1;

src/mcmc.c

@@ -6619,8 +6619,8 @@ void MCMC_PHYREX_Neff(t_tree *tree)
                                 0.0,
                                 1.E+6,
                                 tree->mcmc->num_move_phyrex_neff,
-                                NULL,&(tree->times->c_lnL_times),
-                                NULL,(tree->eval_glnL == YES) ? TIMES_Wrap_Lk_Coalescent : NULL,
+                                NULL,&(tree->mmod->c_lnL),
+                                NULL,(tree->eval_glnL == YES) ? PHYREX_Wrap_Lk : NULL,
                                 tree->mcmc->move_type[tree->mcmc->num_move_phyrex_neff],
                                 NO,NULL,tree,NULL);
     }
@@ -9032,7 +9032,8 @@ void MCMC_PHYREX_Ldsk_Given_Disk(t_tree *tree)
   all_disks   = NULL;
   n_all_disks = 0;
   c           = -1.;
-  K           = tree->mcmc->tune_move[tree->mcmc->num_move_phyrex_ldsk_given_disk];
+  /* K           = tree->mcmc->tune_move[tree->mcmc->num_move_phyrex_ldsk_given_disk]; */
+  K           = 1.0;
   rad         = tree->mmod->rad;
 
 
@@ -9103,13 +9104,15 @@ void MCMC_PHYREX_Ldsk_Given_Disk(t_tree *tree)
               c = disk->ldsk->coord->lonlat[j];
             }
 
+
           err = NO;
           disk->ldsk->coord->lonlat[j] =
             Rnorm_Trunc(c,
                         K*rad,
                         tree->mmod->lim_do->lonlat[j],
                         tree->mmod->lim_up->lonlat[j],&err);
 
+
           if(err == YES) Generic_Exit(__FILE__,__LINE__,__FUNCTION__);
 
           hr -= Log_Dnorm_Trunc(disk->ldsk->coord->lonlat[j],
@@ -9143,7 +9146,6 @@ void MCMC_PHYREX_Ldsk_Given_Disk(t_tree *tree)
       ratio = exp(ratio);
       alpha = MIN(1.,ratio);
 
-      if(disk->prev == NULL) PhyML_Printf("\n. new: %f cur:%f hr: %f",new_glnL,cur_glnL,hr);
 
       u = Uni();
 
@@ -10361,35 +10363,39 @@ void MCMC_PHYREX_Sigsq_Scale(t_tree *tree)
 
       m = exp(K*(Uni()-.5));
       new_scale = m * cur_scale;
-      hr += log(m);
-
-      tree->mmod->sigsq_scale[permut[i]] = new_scale;
-
-      if(tree->eval_glnL == YES) new_lnL_geo = PHYREX_Lk(tree);
 
-      ratio += (new_lnL_geo - cur_lnL_geo);
-
-      ratio += hr;
-
-      ratio = exp(ratio);
-      alpha = MIN(1.,ratio);
-
-        /* Always accept move */
-      if(tree->mcmc->always_yes == YES && new_lnL_geo > UNLIKELY) alpha = 1.0;
-
-      u = Uni();
-
-      assert(isnan(u) == NO && isinf(fabs(u)) == NO);
-
-      if(u > alpha) /* Reject */
+      if(new_scale > tree->mmod->sigsq_scale_min && new_scale < tree->mmod->sigsq_scale_max)
         {
-          tree->mmod->c_lnL = cur_lnL_geo;
-          tree->mmod->sigsq_scale[permut[i]] = cur_scale;
+          hr += log(m);
+
+          tree->mmod->sigsq_scale[permut[i]] = new_scale;
+
+          if(tree->eval_glnL == YES) new_lnL_geo = PHYREX_Lk(tree);
+
+          ratio += (new_lnL_geo - cur_lnL_geo);
+
+          ratio += hr;
+
+          ratio = exp(ratio);
+          alpha = MIN(1.,ratio);
+
+          /* Always accept move */
+          if(tree->mcmc->always_yes == YES && new_lnL_geo > UNLIKELY) alpha = 1.0;
+
+          u = Uni();
+
+          assert(isnan(u) == NO && isinf(fabs(u)) == NO);
+
+          if(u > alpha) /* Reject */
+            {
+              tree->mmod->c_lnL = cur_lnL_geo;
+              tree->mmod->sigsq_scale[permut[i]] = cur_scale;
+            }
+          else
+            {
+              tree->mcmc->acc_move[tree->mcmc->num_move_phyrex_sigsq_scale]++;
+            }
         }
-      else
-        {
-          tree->mcmc->acc_move[tree->mcmc->num_move_phyrex_sigsq_scale]++;
-        }    
     }
 
   Free(permut);

src/mpi_boot.h

@@ -28,7 +28,7 @@ the GNU public licence. See http://www.opensource.org for details.
 #define BootTreeTag 0
 #define BootStatTag 1
 
-int Global_numTask, Global_myRank;
+extern int Global_numTask, Global_myRank;
 
 
 void Bootstrap_MPI(t_tree *tree);

src/phyrex.c

@@ -400,8 +400,8 @@ void PHYREX_XML(char *xml_filename)
   if(mixt_tree->mmod->id == RW || mixt_tree->mmod->id == RRW)
     {
       Make_Contrasts(mixt_tree);
-      mixt_tree->mmod->max_num_of_intervals = 3*tree->n_otu-2;
-  }
+      mixt_tree->mmod->max_num_of_intervals = 3*mixt_tree->n_otu-2;
+    }
 
   assert(PHYREX_Check_Struct(mixt_tree));
   PHYREX_Lk(mixt_tree);
@@ -807,7 +807,10 @@ phydbl *PHYREX_MCMC(t_tree *tree)
       PhyML_Fprintf(fp_stats,"%s\t","tuneLdskGivenDisk");
       PhyML_Fprintf(fp_stats,"%s\t","tuneIndelDiskSerial");
       PhyML_Fprintf(fp_stats,"%s\t","tuneIndelHitSerial");
-      for(int i=0;i<2*tree->n_otu-2;++i) PhyML_Fprintf(fp_stats,"s%d\t",i);
+      for(int i=0;i<2*tree->n_otu-2;++i) PhyML_Fprintf(fp_stats,"s%d%c\t",
+                                                       i,
+                                                       (tree->a_nodes[i]==tree->n_root->v[1] ||
+                                                        tree->a_nodes[i]==tree->n_root->v[2]) ? '*' : ' ');
     }
 
 

src/rrw.c

@@ -39,9 +39,11 @@ phydbl RRW_Lk(t_tree *tree)
   d_sigsq_scale = RRW_Prior_Sigsq_Scale(tree);
 
   // Make sure node times are set back to their original values
-  assert(fabs(t_dum - tree->times->nd_t[idx_dum]) < 1.E-5);
-
-  tree->mmod->c_lnL = d_fwd - d_ic + d_coal + d_sigsq_scale;
+  assert(fabs(t_dum - tree->times->nd_t[idx_dum]) < 1.E-4);
+
+  tree->times->c_lnL_times = d_coal;
+  /* tree->mmod->c_lnL = d_fwd - d_ic + d_coal + d_sigsq_scale; */
+  tree->mmod->c_lnL = d_fwd + d_coal + d_sigsq_scale;
 
   return(tree->mmod->c_lnL);
 }

src/utilities.h

@@ -1940,6 +1940,9 @@ typedef struct __Migrep_Model{
   struct __Geo_Coord         *lim_up; // max longitude and lattitude                        
   struct __Geo_Coord         *lim_do; // min longitude and lattitude                        
   phydbl                *sigsq_scale; // Scaling factors for the variance parameter in RRW model
+
+  phydbl             sigsq_scale_min;
+  phydbl             sigsq_scale_max;
 
   short int                       id;
   int                          n_dim;