[med-svn] [Git][med-team/phyml][upstream] 2 commits: New upstream version 3.3.20211118
Étienne Mollier (@emollier)
gitlab at salsa.debian.org
Fri Feb 11 15:50:46 GMT 2022
Étienne Mollier pushed to branch upstream at Debian Med / phyml
Commits:
37244647 by Étienne Mollier at 2021-11-22T23:29:53+01:00
New upstream version 3.3.20211118
- - - - -
55932608 by Étienne Mollier at 2022-02-11T15:39:13+01:00
New upstream version 3.3.20211231
- - - - -
24 changed files:
- .travis.yml
- configure.ac
- src/Makefile.am
- src/alrt.c
- src/date.c
- src/free.c
- src/free.h
- src/init.c
- src/io.c
- src/lk.c
- src/main.c
- src/make.c
- src/mcmc.c
- src/mixt.c
- src/mixt.h
- src/models.c
- src/mpi_boot.c
- src/optimiz.c
- src/phyrex.c
- src/slfv.c
- src/spr.c
- src/utilities.c
- src/utilities.h
- src/xml.c
Changes:
=====================================
.travis.yml
=====================================
@@ -1,34 +1,34 @@
-sudo: required
-cache: apt
-language: c
+# sudo: required
+# cache: apt
+# language: c
-os:
- - linux
- - osx
+# os:
+# - linux
+# - osx
-compiler:
- - gcc
+# compiler:
+# - gcc
-script:
-# build and install software into $HOME/local
- - ./autogen.sh
- - ./configure
- - make
- - sudo make install
- - ./configure --enable-phyrex
- - make clean
- - make
- - ./configure --enable-phytime
- - make clean
- - make
+# script:
+# # build and install software into $HOME/local
+# - ./autogen.sh
+# - ./configure
+# - make
+# - sudo make install
+# - ./configure --enable-phyrex
+# - make clean
+# - make
+# - ./configure --enable-phytime
+# - make clean
+# - make
-# put the installed binary into the $PATH
- - export PATH=$HOME/local/bin:$PATH
+# # put the installed binary into the $PATH
+# - export PATH=$HOME/local/bin:$PATH
-# install testiphy
-# - cd
-# - git clone https://gitlab.com/testiphy/testiphy.git
+# # install testiphy
+# # - cd
+# # - git clone https://gitlab.com/testiphy/testiphy.git
-# run testiphy
-# - cd testiphy
-# - ./testiphy phyml
\ No newline at end of file
+# # run testiphy
+# # - cd testiphy
+# # - ./testiphy phyml
\ No newline at end of file
=====================================
configure.ac
=====================================
@@ -1,16 +1,14 @@
# -*- Autoconf -*-
# Process this file with autoconf to produce a configure script.
-AC_INIT([PhyML],esyscmd([sh -c "echo "3.3." | tr -d '\n' ; echo "20211021" | tr -d '\n'"]),[guindon at lirmm.fr])
+AC_INIT([PhyML],esyscmd([sh -c "echo "3.3." | tr -d '\n' ; echo "20211231" | tr -d '\n'"]),[guindon at lirmm.fr])
AM_INIT_AUTOMAKE([foreign])
-
-AC_CONFIG_SRCDIR([src/simu.c],[doc/phyml-manual.tex])
+AC_CONFIG_SRCDIR([src/main.c],[doc/phyml-manual.tex])
AC_CONFIG_HEADERS([config.h])
AC_DEFINE([UNIX],[1],[Unix tag on])
AC_DEFINE([DEBUG],[1],[Debug tag on])
-AM_INIT_AUTOMAKE
AC_CANONICAL_HOST
AC_PROG_CC([gcc])
@@ -109,16 +107,6 @@ AS_IF([test "x$enable_win" = "xyes"],[CC="i686-w64-mingw32-gcc"])
AS_IF([test "x$enable_win" = "xyes"],AC_DEFINE([WIN32],[1],[WIN32 tag on]))
AM_CONDITIONAL([WANT_WIN], [test "$enable_win" = yes])
-
-dnl AC_ARG_ENABLE([beagle], [AS_HELP_STRING([--enable-beagle], [Compute likelihoods using BEAGLE library.])], [beagle=yes],[beagle=no])
-dnl AS_IF([test "x$enable_beagle" = "xyes"],[PKG_CHECK_MODULES([BEAGLE],[hmsbeagle-1])])
-dnl AS_IF([test "x$enable_beagle" = "xyes"],
-dnl [CFLAGS="${CFLAGS} ${BEAGLE_CFLAGS} ${BEAGLE_LIBS} -lstdc++"],
-dnl [CFLAGS=${CFLAGS}])
-AS_IF([test "x$enable_beagle" = "xyes"], AC_DEFINE([BEAGLE],[1],[BEAGLE tag on]))
-AM_CONDITIONAL([WANT_BEAGLE], [test "$enable_beagle" = yes])
-
-
AC_ARG_ENABLE([phytime],[AS_HELP_STRING([--enable-phytime],[Compile PhyTime])],[phytime=yes],[phytime=no])
AM_CONDITIONAL([WANT_PHYTIME], [test "$phytime" = yes])
if test "$phytime" = yes; then
@@ -132,30 +120,12 @@ if test "$phyml" = yes; then
fi
-AC_ARG_ENABLE([tiporder],[AS_HELP_STRING([--enable-tiporder],[Compile tiporder])],[tiporder=yes],[tiporder=no])
-AM_CONDITIONAL([WANT_TIPORDER], [test "$tiporder" = yes])
-if test "$tiporder" = yes; then
- AC_DEFINE([TIPORDER],[1],[TIPORDER tag on])
-fi
-
-AC_ARG_ENABLE([part],[AS_HELP_STRING([--enable-part],[Compile Part])],[part=yes],[part=no])
-AM_CONDITIONAL([WANT_PART], [test "$part" = yes])
-if test "$part" = yes; then
- AC_DEFINE([PART],[1],[PART tag on])
-fi
-
AC_ARG_ENABLE([rwrap],[AS_HELP_STRING([--enable-rwrap],[Compile Rwrap])],[rwrap=yes],[rwrap=no])
AM_CONDITIONAL([WANT_RWRAP], [test "$rwrap" = yes])
if test "$rwrap" = yes; then
AC_DEFINE([RWRAP],[1],[RWRAP tag on])
fi
-AC_ARG_ENABLE([phycont],[AS_HELP_STRING([--enable-phycont],[Compile PhyCont])],[phycont=yes],[phycont=no])
-AM_CONDITIONAL([WANT_PHYCONT], [test "$phycont" = yes])
-if test "$phycont" = yes; then
- AC_DEFINE([PHYCONT],[1],[PHYCONT tag on])
-fi
-
AC_ARG_ENABLE([rf],[AS_HELP_STRING([--enable-rf],[Compile RF])],[rf=yes],[rf=no])
AM_CONDITIONAL([WANT_RF], [test "$rf" = yes])
if test "$rf" = yes; then
@@ -174,25 +144,6 @@ if test "$evolve" = yes; then
AC_DEFINE([EVOLVE],[1],[EVOLVE tag on])
fi
-AC_ARG_ENABLE([invitee],[AS_HELP_STRING([--enable-invitee],[Compile invitee])],[invitee=yes],[invitee=no])
-AM_CONDITIONAL([WANT_INVITEE], [test "$invitee" = yes])
-if test "$invitee" = yes; then
- AC_DEFINE([INVITEE],[1],[INVITEE tag on])
-fi
-
-AC_ARG_ENABLE([geo],[AS_HELP_STRING([--enable-geo],[Compile geo])],[geo=yes],[geo=no])
-AM_CONDITIONAL([WANT_GEO], [test "$geo" = yes])
-if test "$geo" = yes; then
- AC_DEFINE([GEO],[1],[GEO tag on])
-fi
-
-AC_ARG_ENABLE([checkpoint],[AS_HELP_STRING([--enable-checkpoint],[Compile checkpoint])],[checkpoint=yes],[checkpoint=no])
-AM_CONDITIONAL([WANT_CHECKPOINT], [test "$checkpoint" = yes])
-if test "$checkpoint" = yes; then
- AC_DEFINE([CHECKPOINT],[1],[CHECKPOINT tag on])
-fi
-
-
AC_ARG_ENABLE([phyrex],[AS_HELP_STRING([--enable-phyrex],[Compile phyrex])],[phyrex=yes],[phyrex=no])
AM_CONDITIONAL([WANT_PHYREX], [test "$phyrex" = yes])
if test "$phyrex" = yes; then
@@ -219,16 +170,10 @@ if test "$date" = yes; then
fi
if test "$phytime" = no; then
-if test "$tiporder" = no; then
-if test "$part" = no; then
if test "$rwrap" = no; then
-if test "$phycont" = no; then
if test "$test" = no; then
if test "$rf" = no; then
-if test "$invitee" = no; then
-if test "$geo" = no; then
if test "$evolve" = no; then
-if test "$checkpoint" = no; then
if test "$phyrex" = no; then
if test "$phyrexsim" = no; then
if test "$date" = no; then
@@ -241,12 +186,6 @@ fi
fi
fi
fi
-fi
-fi
-fi
-fi
-fi
-fi
AC_CONFIG_FILES([Makefile src/Makefile doc/Makefile])
AC_OUTPUT
=====================================
src/Makefile.am
=====================================
@@ -11,21 +11,6 @@ if WANT_PHYTIME
bin_PROGRAMS = phytime
PROG = PHYTIME
else
-if WANT_PHYCONT
-bin_PROGRAMS = phycont
-PROG = PHYCONT
-else
-if WANT_PART
-bin_PROGRAMS = part
-PROG = PART
-else
-if WANT_RWRAP
-PROG = RWRAP
-else
-if WANT_TIPORDER
-bin_PROGRAMS = tiporder
-PROG = TIPORDER
-else
if WANT_RF
bin_PROGRAMS = rf
PROG = RF
@@ -42,26 +27,10 @@ if WANT_TEST
bin_PROGRAMS = test
PROG = TEST
else
-if WANT_INVITEE
-bin_PROGRAMS = invitee
-PROG = INVITEE
-else
-if WANT_GEO
-bin_PROGRAMS = phylogeo
-PROG = GEO
-else
if WANT_EVOLVE
bin_PROGRAMS = evolve
PROG = EVOLVE
else
-if WANT_CHECKPOINT
-bin_PROGRAMS = checkpoint
-PROG = CHECKPOINT
-else
-if WANT_BEAGLE
-bin_PROGRAMS = phyml-beagle
-PROG = PHYML
-else
if WANT_PHYREX
bin_PROGRAMS = phyrex
PROG = PHYREX
@@ -86,14 +55,6 @@ endif
endif
endif
endif
-endif
-endif
-endif
-endif
-endif
-endif
-endif
-endif
if WANT_PHYTIME
@@ -118,7 +79,6 @@ rates.c rates.h\
mcmc.c mcmc.h\
stats.c stats.h\
mg.c mg.h\
-tiporder.c tiporder.h\
io.c io.h\
make.c make.h\
mixt.c mixt.h\
@@ -132,33 +92,6 @@ ancestral.c ancestral.h\
xml.c xml.h
phytime_LDADD = -lm
else
-if WANT_PHYCONT
-phycont_SOURCES = main.c \
-utilities.c utilities.h\
-optimiz.c optimiz.h\
-lk.c lk.h\
-bionj.c bionj.h\
-models.c models.h\
-free.c free.h\
-help.c help.h\
-simu.c simu.h\
-eigen.c eigen.h\
-pars.c pars.h\
-alrt.c alrt.h\
-interface.c interface.h\
-cl.c cl.h\
-spr.c spr.h\
-times.c times.h\
-draw.c draw.h\
-rates.c rates.h\
-mcmc.c mcmc.h\
-tbe.c tbe.h \
-stats.c stats.h\
-mg.c mg.h\
-tiporder.c tiporder.h
-# continuous.c continuous.h
-phycont_LDADD = -lm
-else
# if WANT_RWRAP
# lib_LTLIBRARIES = librwrap.la
# librwrap_la_SOURCES = main.c \
@@ -194,55 +127,6 @@ else
# librwrap_la_LDFLAGS = -I/usr/share/R/include -shared -module -flat_namespace
# librwrap_la_CFLAGS=-std=gnu99 -fPIC -Wl,-z,defs
# else
-if WANT_PART
-part_SOURCES = main.c \
-utilities.c utilities.h\
-optimiz.c optimiz.h\
-lk.c lk.h\
-bionj.c bionj.h\
-models.c models.h\
-free.c free.h\
-help.c help.h\
-simu.c simu.h\
-eigen.c eigen.h\
-pars.c pars.h\
-alrt.c alrt.h\
-interface.c interface.h\
-cl.c cl.h\
-mg.c mg.h\
-spr.c spr.h\
-tbe.c tbe.h \
-draw.c draw.h\
-stats.c stats.h\
-tiporder.c tiporder.h
-part_LDADD = -lm
-else
-if WANT_TIPORDER
-tiporder_SOURCES = main.c \
-utilities.c utilities.h\
-optimiz.c optimiz.h\
-lk.c lk.h\
-bionj.c bionj.h\
-models.c models.h\
-free.c free.h\
-help.c help.h\
-simu.c simu.h\
-eigen.c eigen.h\
-pars.c pars.h\
-alrt.c alrt.h\
-interface.c interface.h\
-cl.c cl.h\
-tbe.c tbe.h \
-mg.c mg.h\
-times.c times.h\
-mcmc.c mcmc.h\
-rates.c rates.h\
-spr.c spr.h\
-draw.c draw.h\
-stats.c stats.h\
-tiporder.c tiporder.h
-tiporder_LDADD = -lm
-else
if WANT_RF
rf_SOURCES = main.c \
utilities.c utilities.h\
@@ -261,10 +145,6 @@ cl.c cl.h\
spr.c spr.h\
draw.c draw.h\
stats.c stats.h\
-rates.c rates.h\
-mcmc.c mcmc.h\
-times.c times.h\
-tiporder.c tiporder.h\
mg.c mg.h\
io.c io.h\
make.c make.h\
@@ -273,7 +153,6 @@ init.c init.h\
xml.c xml.h\
mixt.c mixt.h\
tbe.c tbe.h \
-date.c date.h\
ancestral.c ancestral.h\
sse.c sse.h\
avx.c avx.h
@@ -332,7 +211,6 @@ stats.c stats.h\
rates.c rates.h\
mcmc.c mcmc.h\
times.c times.h\
-tiporder.c tiporder.h\
mg.c mg.h\
io.c io.h\
make.c make.h\
@@ -368,7 +246,6 @@ stats.c stats.h\
rates.c rates.h\
mcmc.c mcmc.h\
times.c times.h\
-tiporder.c tiporder.h\
mg.c mg.h\
io.c io.h\
make.c make.h\
@@ -382,39 +259,6 @@ ancestral.c ancestral.h\
date.c date.h
test_LDADD = -lm
else
-if WANT_INVITEE
-invitee_SOURCES = main.c \
-utilities.c utilities.h\
-optimiz.c optimiz.h\
-lk.c lk.h\
-bionj.c bionj.h\
-models.c models.h\
-free.c free.h\
-help.c help.h\
-simu.c simu.h\
-eigen.c eigen.h\
-pars.c pars.h\
-alrt.c alrt.h\
-interface.c interface.h\
-cl.c cl.h\
-spr.c spr.h\
-draw.c draw.h\
-stats.c stats.h\
-rates.c rates.h\
-mcmc.c mcmc.h\
-times.c times.h\
-tiporder.c tiporder.h\
-mg.c mg.h\
-io.c io.h\
-make.c make.h\
-nexus.c nexus.h\
-tbe.c tbe.h \
-init.c init.h\
-xml.c xml.h\
-invitee.c invitee.h\
-mixt.c mixt.h
-invitee_LDADD = -lm
-else
if WANT_EVOLVE
evolve_SOURCES = main.c\
utilities.c utilities.h\
@@ -437,7 +281,6 @@ stats.c stats.h\
rates.c rates.h\
mcmc.c mcmc.h\
times.c times.h\
-tiporder.c tiporder.h\
mg.c mg.h\
io.c io.h\
make.c make.h\
@@ -447,72 +290,6 @@ xml.c xml.h\
mixt.c mixt.h
evolve_LDADD = -lm
else
-if WANT_GEO
-phylogeo_SOURCES = main.c\
-utilities.c utilities.h\
-optimiz.c optimiz.h\
-lk.c lk.h\
-bionj.c bionj.h\
-models.c models.h\
-free.c free.h\
-help.c help.h\
-simu.c simu.h\
-eigen.c eigen.h\
-pars.c pars.h\
-alrt.c alrt.h\
-interface.c interface.h\
-cl.c cl.h\
-spr.c spr.h\
-draw.c draw.h\
-stats.c stats.h\
-rates.c rates.h\
-mcmc.c mcmc.h\
-times.c times.h\
-tiporder.c tiporder.h\
-mg.c mg.h\
-io.c io.h\
-make.c make.h\
-tbe.c tbe.h \
-nexus.c nexus.h\
-init.c init.h\
-xml.c xml.h\
-mixt.c mixt.h\
-geo.c geo.h
-phylogeo_LDADD = -lm
-else
-if WANT_CHECKPOINT
-checkpoint_SOURCES = main.c\
-utilities.c utilities.h\
-optimiz.c optimiz.h\
-lk.c lk.h\
-bionj.c bionj.h\
-models.c models.h\
-free.c free.h\
-help.c help.h\
-simu.c simu.h\
-eigen.c eigen.h\
-pars.c pars.h\
-alrt.c alrt.h\
-interface.c interface.h\
-cl.c cl.h\
-spr.c spr.h\
-draw.c draw.h\
-stats.c stats.h\
-rates.c rates.h\
-mcmc.c mcmc.h\
-times.c times.h\
-tiporder.c tiporder.h\
-mg.c mg.h\
-io.c io.h\
-make.c make.h\
-tbe.c tbe.h \
-nexus.c nexus.h\
-init.c init.h\
-xml.c xml.h\
-mixt.c mixt.h\
-checkpoint.c checkpoint.h
-checkpoint_LDADD = -lm
-else
if WANT_PHYREX
phyrex_SOURCES = main.c\
utilities.c utilities.h\
@@ -573,7 +350,6 @@ stats.c stats.h\
rates.c rates.h\
mcmc.c mcmc.h\
times.c times.h\
-tiporder.c tiporder.h\
mg.c mg.h\
io.c io.h\
make.c make.h\
@@ -589,39 +365,6 @@ ancestral.c ancestral.h\
phyrex.c phyrex.h
phyrexsim_LDADD = -lm
else
-if WANT_BEAGLE
-phyml_beagle_SOURCES = main.c \
-utilities.c utilities.h\
-optimiz.c optimiz.h\
-lk.c lk.h\
-bionj.c bionj.h\
-models.c models.h\
-free.c free.h\
-help.c help.h\
-simu.c simu.h\
-eigen.c eigen.h\
-pars.c pars.h\
-alrt.c alrt.h\
-interface.c interface.h\
-cl.c cl.h\
-spr.c spr.h\
-draw.c draw.h\
-stats.c stats.h\
-rates.c rates.h\
-mcmc.c mcmc.h\
-times.c times.h\
-tiporder.c tiporder.h\
-mg.c mg.h\
-tbe.c tbe.h \
-io.c io.h\
-make.c make.h\
-nexus.c nexus.h\
-init.c init.h\
-xml.c xml.h\
-mixt.c mixt.h\
-beagle_utils.c beagle_utils.h
-phyml_beagle_LDADD = -lm -lhmsbeagle
-else
if WANT_DATE
date_SOURCES = main.c \
utilities.c utilities.h\
@@ -643,7 +386,6 @@ stats.c stats.h\
rates.c rates.h\
mcmc.c mcmc.h\
times.c times.h\
-tiporder.c tiporder.h\
mg.c mg.h\
io.c io.h\
make.c make.h\
@@ -692,13 +434,6 @@ endif
endif
endif
endif
-endif
-endif
-endif
-endif
-endif
-endif
-endif
all-am: intro $(bin_PROGRAMS)
=====================================
src/alrt.c
=====================================
@@ -61,15 +61,14 @@ int Check_NNI_Five_Branches(t_tree *tree)
init_lnL = Lk(NULL,tree);
MIXT_Set_Alias_Subpatt(NO,tree);
-
- //For every branch
+ // For every branch
for(i=0;i<2*tree->n_otu-3;++i)
- {
+ {
//if this branch is not terminal
if((!tree->a_edges[i]->left->tax) && (!tree->a_edges[i]->rght->tax))
{
result = -1;
-
+
//Try every NNIs for the tested branch
result = NNI_Neigh_BL(tree->a_edges[i],tree);
@@ -389,7 +388,7 @@ int NNI_Neigh_BL(t_edge *b_fcus, t_tree *tree)
PhyML_Fprintf(stderr,"\n. Err. in file %s at line %d\n\n",__FILE__,__LINE__);
Exit("");
}
-
+
/*! edges */
e1 = b_fcus->left->b[b_fcus->l_v1];
e2 = b_fcus->left->b[b_fcus->l_v2];
@@ -406,6 +405,8 @@ int NNI_Neigh_BL(t_edge *b_fcus, t_tree *tree)
var_e3 = Duplicate_Scalar_Dbl(e3->l_var);
var_e4 = Duplicate_Scalar_Dbl(e4->l_var);
+
+
/*! Optimize branch lengths and update likelihoods for
the original configuration.
*/
@@ -430,6 +431,7 @@ int NNI_Neigh_BL(t_edge *b_fcus, t_tree *tree)
lk_temp = Br_Len_Opt(&(b_fcus->rght->b[i]->l->v),b_fcus->rght->b[i],tree);
}
Update_Partial_Lk(tree,b_fcus,b_fcus->rght);
+ lk_temp = Br_Len_Opt(&(b_fcus->l->v),b_fcus,tree);
if(lk_temp < lk0 - tree->mod->s_opt->min_diff_lk_local)
{
@@ -440,6 +442,7 @@ int NNI_Neigh_BL(t_edge *b_fcus, t_tree *tree)
}
while(fabs(lk_temp-lk0) > tree->mod->s_opt->min_diff_lk_global);
//until no significative improvement is detected
+
lk0 = tree->c_lnL;
@@ -485,6 +488,7 @@ int NNI_Neigh_BL(t_edge *b_fcus, t_tree *tree)
}
+
/*! Do first possible swap */
Swap(v2,b_fcus->left,b_fcus->rght,v3,tree);
MIXT_Set_Alias_Subpatt(YES,tree);
@@ -494,6 +498,7 @@ int NNI_Neigh_BL(t_edge *b_fcus, t_tree *tree)
lk1 = Lk(b_fcus,tree);
MIXT_Set_Alias_Subpatt(NO,tree);
+
for(i=0;i<3;i++)
if(b_fcus->left->v[i] != b_fcus->rght)
Update_Partial_Lk(tree,b_fcus->left->b[i],b_fcus->left);
@@ -527,7 +532,8 @@ int NNI_Neigh_BL(t_edge *b_fcus, t_tree *tree)
}
Update_Partial_Lk(tree,b_fcus,b_fcus->rght);
-
+ lk_temp = Br_Len_Opt(&(b_fcus->l->v),b_fcus,tree);
+
if(lk_temp < lk1 - tree->mod->s_opt->min_diff_lk_local)
{
PhyML_Fprintf(stderr,"\n. lk_temp = %f lk1 = %f\n",lk_temp,lk1);
@@ -576,6 +582,15 @@ int NNI_Neigh_BL(t_edge *b_fcus, t_tree *tree)
MIXT_Set_Alias_Subpatt(YES,tree);
Update_Partial_Lk(tree,b_fcus,b_fcus->rght);
Update_Partial_Lk(tree,b_fcus,b_fcus->left);
+
+ for(i=0;i<3;i++)
+ if(b_fcus->left->v[i] != b_fcus->rght)
+ Update_Partial_Lk(tree,b_fcus->left->b[i],b_fcus->left);
+
+ for(i=0;i<3;i++)
+ if(b_fcus->rght->v[i] != b_fcus->left)
+ Update_Partial_Lk(tree,b_fcus->rght->b[i],b_fcus->rght);
+
lk_temp = Lk(b_fcus,tree);
MIXT_Set_Alias_Subpatt(NO,tree);
if((lk_temp > lk_init + tree->mod->s_opt->min_diff_lk_local) || (lk_temp < lk_init - tree->mod->s_opt->min_diff_lk_local))
@@ -584,7 +599,7 @@ int NNI_Neigh_BL(t_edge *b_fcus, t_tree *tree)
PhyML_Fprintf(stderr,"\n== Err. in file %s at line %d\n",__FILE__,__LINE__);
Warn_And_Exit("");
}
-
+
/*! do the second possible swap */
Swap(v2,b_fcus->left,b_fcus->rght,v4,tree);
Set_Both_Sides(YES,tree);
@@ -652,7 +667,7 @@ int NNI_Neigh_BL(t_edge *b_fcus, t_tree *tree)
}
Update_Partial_Lk(tree,b_fcus,b_fcus->rght);
-
+ lk_temp = Br_Len_Opt(&(b_fcus->l->v),b_fcus,tree);
}
while(FABS(lk_temp-lk2) > tree->mod->s_opt->min_diff_lk_global);
//until no significative improvement is detected
@@ -694,6 +709,7 @@ int NNI_Neigh_BL(t_edge *b_fcus, t_tree *tree)
Update_PMat_At_Given_Edge(e4,tree);
Update_PMat_At_Given_Edge(b_fcus,tree);
+
MIXT_Set_Alias_Subpatt(YES,tree);
Update_Partial_Lk(tree,b_fcus,b_fcus->left);
Update_Partial_Lk(tree,b_fcus,b_fcus->rght);
@@ -717,7 +733,7 @@ int NNI_Neigh_BL(t_edge *b_fcus, t_tree *tree)
Warn_And_Exit("");
}
-
+
//save likelihoods in NNI structures
b_fcus->nni->lk0 = lk0;
b_fcus->nni->lk1 = lk1;
@@ -755,7 +771,7 @@ int NNI_Neigh_BL(t_edge *b_fcus, t_tree *tree)
Free_Scalar_Dbl(var_e3);
Free_Scalar_Dbl(var_e4);
Free_Scalar_Dbl(v_init);
-
+
return result;
}
=====================================
src/date.c
=====================================
@@ -306,17 +306,20 @@ void DATE_XML(char *xml_filename)
MIXT_Check_Model_Validity(mixt_tree);
- MIXT_Init_Model(mixt_tree);
+ MIXT_Chain_Models(mixt_tree);
+ Init_Model(mixt_tree->mod->io->cdata,mixt_tree->mod,mixt_tree->mod->io);
+ Set_Model_Parameters(mixt_tree->mod);
Print_Data_Structure(NO,stdout,mixt_tree);
tree = MIXT_Starting_Tree(mixt_tree);
Add_Root(tree->a_edges[0],tree);
Copy_Tree(tree,mixt_tree);
Free_Tree(tree);
- MIXT_Connect_Cseqs_To_Nodes(mixt_tree);
- MIXT_Init_T_Beg(mixt_tree);
- MIXT_Make_Tree_For_Lk(mixt_tree);
- MIXT_Make_Tree_For_Pars(mixt_tree);
- MIXT_Make_Spr(mixt_tree);
+ Copy_Tree(mixt_tree,mixt_tree->next);
+ Connect_CSeqs_To_Nodes(mixt_tree->mod->io->cdata,mixt_tree->mod->io,mixt_tree);
+ Init_T_Beg(mixt_tree);
+ Make_Tree_For_Lk(mixt_tree);
+ Make_Tree_For_Pars(mixt_tree);
+ Make_Spr(mixt_tree);
MIXT_Chain_All(mixt_tree);
MIXT_Check_Edge_Lens_In_All_Elem(mixt_tree);
MIXT_Turn_Branches_OnOff_In_All_Elem(ON,mixt_tree);
=====================================
src/free.c
=====================================
@@ -24,52 +24,40 @@ void Free_Clade(t_clad *this)
//////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////
-void Free_All_Nodes_Light(t_tree *mixt_tree)
+void Free_All_Nodes_Light(t_tree *tree)
{
- int i;
- t_tree *tree;
-
- tree = mixt_tree;
- do
+ if(tree->a_nodes != NULL)
{
- for(i=0;i<2*tree->n_otu-1;++i) if(tree->a_nodes[i] != NULL) Free_Node(tree->a_nodes[i]);
+ for(int i=0;i<2*tree->n_otu-1;++i) if(tree->a_nodes[i] != NULL) Free_Node(tree->a_nodes[i]);
Free(tree->a_nodes);
- tree = tree->next;
}
- while(tree);
-
}
//////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////
-void Free_All_Edges_Light(t_tree *mixt_tree)
+void Free_All_Edges_Light(t_tree *tree)
{
- int i;
- t_tree *tree;
-
- tree = mixt_tree;
-
- for(i=0;i<2*tree->n_otu-1;++i)
+ if(tree->a_edges != NULL)
{
- if(tree->a_edges[i] != NULL)
- {
- Free_Scalar_Dbl(tree->a_edges[i]->l);
- Free_Scalar_Dbl(tree->a_edges[i]->l_old);
- Free_Scalar_Dbl(tree->a_edges[i]->l_var);
- Free_Scalar_Dbl(tree->a_edges[i]->l_var_old);
- }
+ for(int i=0;i<2*tree->n_otu-1;++i)
+ if(tree->a_edges[i] != NULL)
+ Free_Edge(tree->a_edges[i]);
+ Free(tree->a_edges);
}
+}
- do
+//////////////////////////////////////////////////////////////
+//////////////////////////////////////////////////////////////
+
+void Free_All_Edges_Lens(t_tree *tree)
+{
+ if(tree->a_edges != NULL)
{
- for(i=0;i<2*tree->n_otu-1;++i)
+ for(int i=0;i<2*tree->n_otu-1;++i)
if(tree->a_edges[i] != NULL)
- Free_Edge(tree->a_edges[i]);
- Free(tree->a_edges);
- tree = tree->next;
+ Free_Edge_Length(tree->a_edges[i]);
}
- while(tree);
}
//////////////////////////////////////////////////////////////
@@ -77,10 +65,10 @@ void Free_All_Edges_Light(t_tree *mixt_tree)
void Free_Edge_Length(t_edge *b)
{
- if(b->l) Free_Scalar_Dbl(b->l);
- if(b->l_old) Free_Scalar_Dbl(b->l_old);
- if(b->l_var) Free_Scalar_Dbl(b->l_var);
- if(b->l_var_old) Free_Scalar_Dbl(b->l_var_old);
+ if(b->l != NULL) Free_Scalar_Dbl(b->l);
+ if(b->l_old != NULL) Free_Scalar_Dbl(b->l_old);
+ if(b->l_var != NULL) Free_Scalar_Dbl(b->l_var);
+ if(b->l_var_old != NULL) Free_Scalar_Dbl(b->l_var_old);
}
//////////////////////////////////////////////////////////////
@@ -120,11 +108,7 @@ void Free_Node(t_node *n)
}
if(n->ori_name) { Free(n->ori_name); n->ori_name = NULL; }
-
- /* if(n->name) { Free(n->name); n->name = NULL; } */
- /* Don't do that: see Copy_Tax_Names_To_Tip_Labels
- tree->a_nodes[i]->ori_name = tree->a_nodes[i]->name; */
-
+
Free(n);
}
@@ -160,58 +144,35 @@ void Free_Mat(matrix *mat)
void Free_Partial_Lk(phydbl *p_lk, int len, int n_catg)
{
Free(p_lk);
-
-/* int i,j; */
-/* for(i=0;i<len;i++) */
-/* { */
-/* for(j=0;j<n_catg;j++) Free((*p_lk)[i][j]); */
-/* Free((*p_lk)[i]); */
-/* } */
-/* Free((*p_lk)); */
-/* (*p_lk) = NULL; */
}
//////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////
-void Free_Tree(t_tree *mixt_tree)
+void Free_Tree(t_tree *tree)
{
- t_tree *tree;
- t_tree *next;
-
- tree = mixt_tree;
- do
+ if(tree->is_mixt_tree == YES)
+ {
+ MIXT_Free_Tree(tree);
+ }
+ else
{
if(tree->mat) Free_Mat(tree->mat);
- Free(tree->t_dir);
+ if(tree->t_dir) Free(tree->t_dir);
if(tree->short_l) Free(tree->short_l);
- if(tree->mutmap) Free(tree->mutmap);
+ if(tree->mutmap) Free(tree->mutmap);
Free_Bip(tree);
- Free(tree->curr_path);
- tree = tree->next;
- }
- while(tree);
-
- Free_All_Edges_Light(mixt_tree);
- Free_All_Nodes_Light(mixt_tree);
-
- tree = mixt_tree;
- next = mixt_tree->next;
- do
- {
+ Free(tree->curr_path);
+ Free_All_Edges_Lens(tree);
+ Free_All_Edges_Light(tree);
+ Free_All_Nodes_Light(tree);
Free(tree);
- tree = next;
- if(!tree) break;
- next = next->next;
}
- while(tree);
-
}
//////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////
-
void Free_Bip(t_tree *tree)
{
int i,j;
@@ -219,11 +180,11 @@ void Free_Bip(t_tree *tree)
if(tree->has_bip)
{
For(i,2*tree->n_otu-2)
- {
- Free(tree->a_nodes[i]->bip_size);
- for(j=0;j<3;j++) Free(tree->a_nodes[i]->bip_node[j]);
- Free(tree->a_nodes[i]->bip_node);
- }
+ {
+ Free(tree->a_nodes[i]->bip_size);
+ for(j=0;j<3;j++) Free(tree->a_nodes[i]->bip_node[j]);
+ Free(tree->a_nodes[i]->bip_node);
+ }
}
tree->has_bip = NO;
}
@@ -334,36 +295,27 @@ void Free_Tree_Ins_Tar(t_tree *tree)
//////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////
-void Free_Tree_Pars(t_tree *mixt_tree)
+void Free_Tree_Pars(t_tree *tree)
{
int i;
- t_tree *tree;
- tree = mixt_tree;
- do
+ Free(tree->step_mat);
+ Free(tree->site_pars);
+
+ for(i=0;i<2*tree->n_otu-3;++i) Free_Edge_Pars(tree->a_edges[i]);
+
+ if(tree->n_root)
{
- Free(tree->step_mat);
- Free(tree->site_pars);
-
- for(i=0;i<2*tree->n_otu-3;++i)
- {
- Free_Edge_Pars(tree->a_edges[i]);
+ Free_Edge_Pars_Left(tree->n_root->b[1]);
+ Free_Edge_Pars_Left(tree->n_root->b[2]);
}
-
- if(tree->n_root)
- {
- Free_Edge_Pars_Left(tree->n_root->b[1]);
- Free_Edge_Pars_Left(tree->n_root->b[2]);
- }
- else
- {
- Free_Edge_Pars(tree->a_edges[2*tree->n_otu-3]);
- Free_Edge_Pars(tree->a_edges[2*tree->n_otu-2]);
- }
-
- tree = tree->next;
+ else
+ {
+ Free_Edge_Pars(tree->a_edges[2*tree->n_otu-3]);
+ Free_Edge_Pars(tree->a_edges[2*tree->n_otu-2]);
}
- while(tree);
+
+ if(tree->is_mixt_tree == YES) MIXT_Repeat_Task(Free_Tree_Pars,tree);
}
//////////////////////////////////////////////////////////////
@@ -400,81 +352,76 @@ void Free_Edge_Pars(t_edge *b)
//////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////
-void Free_Tree_Lk(t_tree *mixt_tree)
+void Free_Tree_Lk(t_tree *tree)
{
int i;
- t_tree *tree;
- tree = mixt_tree;
- do
+ Free(tree->big_lk_array);
+ Free(tree->c_lnL_sorted);
+ Free(tree->cur_site_lk);
+ Free(tree->old_site_lk);
+ Free(tree->site_lk_cat);
+ Free(tree->fact_sum_scale);
+ Free(tree->unscaled_site_lk_cat);
+ Free(tree->expl);
+
+ for(i=0;i<3;i++) Free(tree->log_lks_aLRT[i]);
+ Free(tree->log_lks_aLRT);
+
+ for(i=0;i<2*tree->n_otu-1;++i) Free_NNI(tree->a_edges[i]->nni);
+ for(i=0;i<2*tree->n_otu-3;++i) Free_Edge_Lk(tree->a_edges[i]);
+ for(i=0;i<2*tree->n_otu-3;++i) Free_Edge_Loc(tree->a_edges[i]);
+
+ if(tree->is_mixt_tree == NO)
{
- Free(tree->c_lnL_sorted);
- Free(tree->cur_site_lk);
- Free(tree->old_site_lk);
- Free(tree->site_lk_cat);
- Free(tree->fact_sum_scale);
- Free(tree->unscaled_site_lk_cat);
- Free(tree->expl);
-
- if(tree->is_mixt_tree == NO)
- {
- Free(tree->dot_prod);
- Free(tree->p_lk_left_pi);
- Free(tree->l_ev);
-
+ Free(tree->dot_prod);
+ Free(tree->p_lk_left_pi);
+ Free(tree->l_ev);
+
#if (defined(__AVX__) || defined(__AVX2__) || defined(__SSE__) || defined(__SSE2__) || defined(__SSE3__))
- Free(tree->_tPij1);
- Free(tree->_tPij2);
- Free(tree->_pmat1plk1);
- Free(tree->_pmat2plk2);
- Free(tree->_plk0);
+ Free(tree->_tPij1);
+ Free(tree->_tPij2);
+ Free(tree->_pmat1plk1);
+ Free(tree->_pmat2plk2);
+ Free(tree->_plk0);
+ Free(tree->_l_ev);
+ Free(tree->_r_ev);
+ Free(tree->_prod_left);
+ Free(tree->_prod_rght);
#endif
-
- for(i=0;i<3;i++) Free(tree->log_lks_aLRT[i]);
- Free(tree->log_lks_aLRT);
-
- Free(tree->div_post_pred_extra_0);
- Free(tree->sum_scale_cat_extra_0);
- Free(tree->sum_scale_extra_0);
- Free(tree->patt_id_extra_0);
- Free(tree->p_lk_extra_0);
- Free(tree->p_lk_tip_extra_0);
-
- Free(tree->div_post_pred_extra_1);
- Free(tree->sum_scale_cat_extra_1);
- Free(tree->sum_scale_extra_1);
- Free(tree->patt_id_extra_1);
- Free(tree->p_lk_extra_1);
- Free(tree->p_lk_tip_extra_1);
-
- for(i=0;i<2*tree->n_otu-1;++i) Free_NNI(tree->a_edges[i]->nni);
- for(i=0;i<2*tree->n_otu-3;++i) Free_Edge_Lk(tree->a_edges[i]);
- for(i=0;i<2*tree->n_otu-3;++i) Free_Edge_Loc(tree->a_edges[i]);
-
- if(tree->n_root != NULL)
- {
- Free(tree->n_root->b[1]->Pij_rr);
- Free(tree->n_root->b[2]->Pij_rr);
- Free(tree->n_root->b[1]->tPij_rr);
- Free(tree->n_root->b[2]->tPij_rr);
- Free_Edge_Lk_Left(tree->n_root->b[1]);
- Free_Edge_Lk_Left(tree->n_root->b[2]);
- Free_Edge_Loc_Left(tree->n_root->b[1]);
- Free_Edge_Loc_Left(tree->n_root->b[2]);
- }
- else
- {
- Free_Edge_Lk(tree->a_edges[2*tree->n_otu-3]);
- Free_Edge_Lk(tree->a_edges[2*tree->n_otu-2]);
- Free_Edge_Loc(tree->a_edges[2*tree->n_otu-3]);
- Free_Edge_Loc(tree->a_edges[2*tree->n_otu-2]);
- }
- }
+
+ Free(tree->div_post_pred_extra_0);
+ Free(tree->sum_scale_cat_extra_0);
+ Free(tree->sum_scale_extra_0);
+ Free(tree->patt_id_extra_0);
+ Free(tree->p_lk_extra_0);
+ Free(tree->p_lk_tip_extra_0);
+
+ Free(tree->div_post_pred_extra_1);
+ Free(tree->sum_scale_cat_extra_1);
+ Free(tree->sum_scale_extra_1);
+ Free(tree->patt_id_extra_1);
+ Free(tree->p_lk_extra_1);
+ Free(tree->p_lk_tip_extra_1);
- tree = tree->next;
+ if(tree->n_root != NULL)
+ {
+ Free_Edge_Lk_Left(tree->n_root->b[1]);
+ Free_Edge_Lk_Left(tree->n_root->b[2]);
+ Free_Edge_Loc_Left(tree->n_root->b[1]);
+ Free_Edge_Loc_Left(tree->n_root->b[2]);
+ }
+ else
+ {
+ Free_Edge_Lk(tree->a_edges[2*tree->n_otu-3]);
+ Free_Edge_Lk(tree->a_edges[2*tree->n_otu-2]);
+ Free_Edge_Loc(tree->a_edges[2*tree->n_otu-3]);
+ Free_Edge_Loc(tree->a_edges[2*tree->n_otu-2]);
+ }
}
- while(tree);
+
+ if(tree->is_mixt_tree == YES) MIXT_Repeat_Task(Free_Tree_Lk,tree);
}
@@ -518,7 +465,6 @@ void Free_Edge_Lk_Rght(t_edge *b)
if(b->p_lk_rght)
{
- Free(b->p_lk_rght);
if(b->sum_scale_rght) Free(b->sum_scale_rght);
}
@@ -537,7 +483,6 @@ void Free_Edge_Lk_Left(t_edge *b)
if(b->p_lk_left)
{
- Free(b->p_lk_left);
if(b->sum_scale_left) Free(b->sum_scale_left);
}
@@ -550,10 +495,7 @@ void Free_Edge_Lk_Left(t_edge *b)
//////////////////////////////////////////////////////////////
void Free_Edge_Lk(t_edge *b)
-{
- Free(b->tPij_rr);
- Free(b->Pij_rr);
-
+{
Free_Edge_Lk_Left(b);
Free_Edge_Lk_Rght(b);
}
@@ -988,7 +930,6 @@ void Free_Tree_List(t_treelist *list)
//////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////
-
void Free_St(supert_tree *st)
{
int i;
@@ -1049,44 +990,40 @@ void Free_One_Spr(t_spr *this_spr)
//////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////
-void Free_Spr_List_One_Edge(t_tree *mixt_tree)
+void Free_Spr_List_One_Edge(t_tree *tree)
{
int i;
- t_tree *tree;
-
- tree = mixt_tree;
- do
- {
- for(i=0;i<tree->size_spr_list_one_edge+1;++i) Free_One_Spr(tree->spr_list_one_edge[i]);
- Free(tree->spr_list_one_edge);
- tree = tree->next;
- }
- while(tree);
+ for(i=0;i<tree->size_spr_list_one_edge+1;++i) Free_One_Spr(tree->spr_list_one_edge[i]);
+ if(tree->is_mixt_tree == YES) MIXT_Repeat_Task(Free_Spr_List_One_Edge,tree);
+ Free(tree->spr_list_one_edge);
}
//////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////
-void Free_Spr_List_All_Edge(t_tree *mixt_tree)
+void Free_Spr_List_All_Edge(t_tree *tree)
{
int i;
- t_tree *tree;
- tree = mixt_tree;
- do
- {
- for(i=0;i<tree->size_spr_list_all_edge+1;++i) Free_One_Spr(tree->spr_list_all_edge[i]);
- Free(tree->spr_list_all_edge);
- tree = tree->next;
- }
- while(tree);
+ for(i=0;i<tree->size_spr_list_all_edge+1;++i) Free_One_Spr(tree->spr_list_all_edge[i]);
+ if(tree->is_mixt_tree == YES) MIXT_Repeat_Task(Free_Spr_List_All_Edge,tree);
+ Free(tree->spr_list_all_edge);
+}
+//////////////////////////////////////////////////////////////
+//////////////////////////////////////////////////////////////
+
+void Free_Best_Spr(t_tree *tree)
+{
+ Free_One_Spr(tree->best_spr);
+ if(tree->is_mixt_tree == YES) MIXT_Repeat_Task(Free_Best_Spr,tree);
}
//////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////
+
void Free_Actual_CSeq(calign *data)
{
int i;
=====================================
src/free.h
=====================================
@@ -99,5 +99,8 @@ void Free_Clade(t_clad *this);
void Free_Label(t_label *lab);
void Free_Contrasts(t_ctrst *ctrst);
void TIMES_Free_Times(t_time *times);
+void Free_Best_Spr(t_tree *tree);
+void Free_All_Edges_Light(t_tree *tree);
+void Free_All_Edges_Lens(t_tree *tree);
#endif
=====================================
src/init.c
=====================================
@@ -149,7 +149,7 @@ void Init_Tree(t_tree *tree, int n_otu)
tree->eval_alnL = YES;
tree->eval_rlnL = YES;
tree->eval_glnL = YES;
- tree->scaling_method = SCALE_FAST;
+ tree->scaling_method = NO;
tree->perform_spr_right_away = YES;
tree->tip_root = 0;
tree->n_edges_traversed = 0;
@@ -1466,7 +1466,7 @@ void Init_Model(calign *data, t_mod *mod, option *io)
Init_Eigen_Struct(mod->eigen);
- Set_Model_Parameters(mod);
+ if(mod->is_mixt_mod == YES) MIXT_Init_Model(mod);
free(dr);
free(di);
=====================================
src/io.c
=====================================
@@ -4425,7 +4425,10 @@ option *PhyML_XML(char *xml_filename)
for(num_rand_tree=0;num_rand_tree<io->mod->s_opt->n_rand_starts;num_rand_tree++)
{
MIXT_Check_Model_Validity(mixt_tree);
- MIXT_Init_Model(mixt_tree);
+ MIXT_Chain_Models(mixt_tree);
+ Init_Model(mixt_tree->mod->io->cdata,mixt_tree->mod,mixt_tree->mod->io);
+ Set_Model_Parameters(mixt_tree->mod);
+
Print_Data_Structure(NO,stdout,mixt_tree);
tree = MIXT_Starting_Tree(mixt_tree);
Copy_Tree(tree,mixt_tree);
@@ -4437,12 +4440,13 @@ option *PhyML_XML(char *xml_filename)
Random_Tree(mixt_tree);
}
- MIXT_Connect_Cseqs_To_Nodes(mixt_tree);
- MIXT_Init_T_Beg(mixt_tree);
+ Copy_Tree(mixt_tree,mixt_tree->next);
+ Connect_CSeqs_To_Nodes(mixt_tree->mod->io->cdata,mixt_tree->mod->io,mixt_tree);
+ Init_T_Beg(mixt_tree);
- MIXT_Make_Tree_For_Pars(mixt_tree);
- MIXT_Make_Tree_For_Lk(mixt_tree);
- MIXT_Make_Spr(mixt_tree);
+ Make_Tree_For_Pars(mixt_tree);
+ Make_Tree_For_Lk(mixt_tree);
+ Make_Spr(mixt_tree);
MIXT_Chain_All(mixt_tree);
MIXT_Check_Edge_Lens_In_All_Elem(mixt_tree);
@@ -4506,7 +4510,9 @@ option *PhyML_XML(char *xml_filename)
MIXT_Init_T_End(mixt_tree);
Print_Data_Structure(YES,mixt_tree->io->fp_out_stats,mixt_tree);
+ Free_Best_Spr(mixt_tree);
Free_Spr_List_One_Edge(mixt_tree);
+ Free_Spr_List_All_Edge(mixt_tree);
Free_Tree_Pars(mixt_tree);
Free_Tree_Lk(mixt_tree);
}
=====================================
src/lk.c
=====================================
@@ -1169,7 +1169,7 @@ void Lk_dLk_Core_One_Class_Eigen_Lr(phydbl *dot_prod, phydbl *expl, unsigned int
//////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////
-phydbl Lk_Core_One_Class_No_Eigen_Lr(phydbl *p_lk_left, phydbl *p_lk_rght, phydbl *Pij,phydbl *pi, int ns, int ambiguity_check, int state)
+phydbl Lk_Core_One_Class_No_Eigen_Lr(phydbl *p_lk_left, phydbl *p_lk_rght, phydbl *Pij, phydbl *pi, int ns, int ambiguity_check, int state)
{
unsigned int l,k;
phydbl lk = 0.0;
@@ -1280,7 +1280,6 @@ void Update_Partial_Lk(t_tree *tree, t_edge *b, t_node *d)
if((tree->io->do_alias_subpatt == YES) &&
(tree->update_alias_subpatt == YES))
Alias_One_Subpatt((d==b->left)?(b->rght):(b->left),d,tree);
-
if(d->tax) return;
=====================================
src/main.c
=====================================
@@ -33,8 +33,6 @@ the GNU public licence. See http://www.opensource.org for details.
#include "beagle_utils.h"
#endif
-
-
#if (defined PHYML || EVOLVE)
int main(int argc, char **argv)
@@ -141,7 +139,8 @@ int main(int argc, char **argv)
if(!io->quiet) PhyML_Printf("\n\n. [Random start %3d/%3d]",num_rand_tree+1,io->mod->s_opt->n_rand_starts);
Init_Model(cdata,mod,io);
-
+ Set_Model_Parameters(mod);
+
#ifdef M4
if(io->mod->use_m4mod) M4_Init_Model(mod->m4mod,cdata,mod);
#endif
@@ -323,7 +322,7 @@ int main(int argc, char **argv)
#ifdef BEAGLE
finalize_beagle_instance(tree);
#endif
- Free_One_Spr(tree->best_spr);
+ Free_Best_Spr(tree);
Free_Spr_List_One_Edge(tree);
Free_Spr_List_All_Edge(tree);
Free_Tree_Pars(tree);
=====================================
src/make.c
=====================================
@@ -19,6 +19,13 @@ void Make_Tree_For_Lk(t_tree *tree)
int i;
calign *cdata;
+
+ if(tree->is_mixt_tree == YES)
+ {
+ MIXT_Make_Tree_For_Lk(tree);
+ return;
+ }
+
const unsigned int ns = tree->mod->ns;
const unsigned int nsns = ns * ns;
@@ -29,7 +36,6 @@ void Make_Tree_For_Lk(t_tree *tree)
cdata = tree->data;
assert(cdata);
-
tree->c_lnL_sorted = (phydbl *)mCalloc(tree->n_pattern,sizeof(phydbl));
tree->cur_site_lk = (phydbl *)mCalloc(tree->n_pattern,sizeof(phydbl));
@@ -37,7 +43,6 @@ void Make_Tree_For_Lk(t_tree *tree)
tree->site_lk_cat = (phydbl *)mCalloc(MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes),sizeof(phydbl));
tree->unscaled_site_lk_cat = (phydbl *)mCalloc(MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->n_pattern,sizeof(phydbl));
tree->fact_sum_scale = (int *)mCalloc(tree->n_pattern,sizeof(int));
-
#if (defined(__AVX__) || defined(__AVX2__))
#ifndef WIN32
@@ -54,6 +59,7 @@ void Make_Tree_For_Lk(t_tree *tree)
if(posix_memalign((void **)&tree->_r_ev,BYTE_ALIGN,(size_t) nsns / sz * sizeof(__m256d))) Generic_Exit(__FILE__,__LINE__,__FUNCTION__);
if(posix_memalign((void **)&tree->_prod_left,BYTE_ALIGN,(size_t) ns / sz * sizeof(__m256d))) Generic_Exit(__FILE__,__LINE__,__FUNCTION__);
if(posix_memalign((void **)&tree->_prod_rght,BYTE_ALIGN,(size_t) ns / sz * sizeof(__m256d))) Generic_Exit(__FILE__,__LINE__,__FUNCTION__);
+ if(posix_memalign((void **)&tree->big_lk_array,BYTE_ALIGN,(size_t) ((3*tree->n_otu-2)*tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns + 2*(2*tree->n_otu-1)*tree->mod->ras->n_catg*tree->mod->ns*tree->mod->ns) * sizeof(phydbl))) Generic_Exit(__FILE__,__LINE__,__FUNCTION__);
#else
tree->dot_prod = _aligned_malloc(tree->n_pattern*tree->mod->ns*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*sizeof(phydbl),BYTE_ALIGN);
tree->expl = _aligned_malloc(3*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns*sizeof(phydbl),BYTE_ALIGN);
@@ -68,6 +74,7 @@ void Make_Tree_For_Lk(t_tree *tree)
tree->_r_ev = _aligned_malloc(ns * ns / sz * sizeof(__m256d),BYTE_ALIGN);
tree->_prod_left = _aligned_malloc(ns / sz * sizeof(__m256d),BYTE_ALIGN);
tree->_prod_rght = _aligned_malloc(ns / sz * sizeof(__m256d),BYTE_ALIGN);
+ tree->big_lk_array = _aligned_malloc(((3*tree->n_otu-2)*tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns + 2*(2*tree->n_otu-1)*tree->mod->ras->n_catg*tree->mod->ns*tree->mod->ns) * sizeof(phydbl),BYTE_ALIGN);
#endif
#elif (defined(__SSE__) || defined(__SSE2__) || defined(__SSE3__))
#ifndef WIN32
@@ -84,6 +91,7 @@ void Make_Tree_For_Lk(t_tree *tree)
if(posix_memalign((void **)&tree->_r_ev,BYTE_ALIGN,(size_t) nsns / sz * sizeof(__m128d))) Generic_Exit(__FILE__,__LINE__,__FUNCTION__);
if(posix_memalign((void **)&tree->_prod_left,BYTE_ALIGN,(size_t) ns / sz * sizeof(__m128d))) Generic_Exit(__FILE__,__LINE__,__FUNCTION__);
if(posix_memalign((void **)&tree->_prod_rght,BYTE_ALIGN,(size_t) ns / sz * sizeof(__m128d))) Generic_Exit(__FILE__,__LINE__,__FUNCTION__);
+ if(posix_memalign((void **)&tree->big_lk_array,BYTE_ALIGN,(size_t) ((3*tree->n_otu-2)*tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns + 2*(2*tree->n_otu-1)*tree->mod->ras->n_catg*tree->mod->ns*tree->mod->ns) * sizeof(phydbl))) Generic_Exit(__FILE__,__LINE__,__FUNCTION__);
#else
tree->dot_prod = _aligned_malloc(tree->n_pattern*tree->mod->ns*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*sizeof(phydbl),BYTE_ALIGN);
tree->expl = _aligned_malloc(3*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns*sizeof(phydbl),BYTE_ALIGN);
@@ -98,15 +106,18 @@ void Make_Tree_For_Lk(t_tree *tree)
tree->_r_ev = _aligned_malloc(ns * ns / sz * sizeof(__m128d),BYTE_ALIGN);
tree->_prod_left = _aligned_malloc(ns / sz * sizeof(__m128d),BYTE_ALIGN);
tree->_prod_rght = _aligned_malloc(ns / sz * sizeof(__m128d),BYTE_ALIGN);
+ tree->big_lk_array = _aligned_malloc(((3*tree->n_otu-2)*tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns + 2*(2*tree->n_otu-1)*tree->mod->ras->n_catg*tree->mod->ns*tree->mod->ns) * sizeof(phydbl),BYTE_ALIGN);
#endif
#elif (!(defined(__AVX__) || defined(__AVX2__) || defined(__SSE__) || defined(__SSE2__) || defined(__SSE3__)))
tree->dot_prod = (phydbl *)mCalloc(tree->n_pattern*tree->mod->ns*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes),sizeof(phydbl));
tree->expl = (phydbl *)mCalloc(3*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns,sizeof(phydbl));
tree->p_lk_left_pi = (phydbl *)mCalloc(ns,sizeof(phydbl));
tree->l_ev = (phydbl *)mCalloc(nsns,sizeof(phydbl));
+ tree->big_lk_array = (phydbl *)mCalloc(((3*tree->n_otu-2)*tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns + 2*(2*tree->n_otu-1)*tree->mod->ras->n_catg*tree->mod->ns*tree->mod->ns),sizeof(phydbl));
#endif
-
-
+
+ tree->big_lk_array_pos = 0;
+
tree->log_lks_aLRT = (phydbl **)mCalloc(3,sizeof(phydbl *));
for(i=0;i<3;i++) tree->log_lks_aLRT[i] = (phydbl *)mCalloc(tree->data->init_len,sizeof(phydbl));
@@ -125,8 +136,10 @@ void Make_Tree_For_Lk(t_tree *tree)
if(tree->n_root != NULL)
{
- Free_Edge_Lk_Rght(tree->n_root->b[1]);
- Free_Edge_Lk_Rght(tree->n_root->b[2]);
+ /* !!!!!!!!!!!!!!!!!!!!!!!!!!!!!1 */
+ /* Free_Edge_Lk_Rght(tree->n_root->b[1]); */
+ /* Free_Edge_Lk_Rght(tree->n_root->b[2]); */
+
Free_Edge_Loc_Rght(tree->n_root->b[1]);
Free_Edge_Loc_Rght(tree->n_root->b[2]);
@@ -159,6 +172,12 @@ void Make_Tree_For_Pars(t_tree *tree)
int i;
calign *cdata;
+ if(tree->is_mixt_tree == YES)
+ {
+ MIXT_Make_Tree_For_Pars(tree);
+ return;
+ }
+
cdata = tree->data;
assert(cdata);
@@ -297,18 +316,26 @@ void Make_Edge_Lk(t_edge *b, t_tree *tree)
b->l_old->v = b->l->v;
-#if (defined(__AVX__) || defined(__AVX2__) || defined(__SSE__) || defined(__SSE2__) || defined(__SSE3__))
-#ifndef WIN32
- if(posix_memalign((void *)&b->Pij_rr,BYTE_ALIGN,(size_t)tree->mod->ras->n_catg*tree->mod->ns*tree->mod->ns*sizeof(phydbl))) Generic_Exit(__FILE__,__LINE__,__FUNCTION__);
- if(posix_memalign((void *)&b->tPij_rr,BYTE_ALIGN,(size_t)tree->mod->ras->n_catg*tree->mod->ns*tree->mod->ns*sizeof(phydbl))) Generic_Exit(__FILE__,__LINE__,__FUNCTION__);
-#else
- b->Pij_rr = _aligned_malloc(tree->mod->ras->n_catg*tree->mod->ns*tree->mod->ns*sizeof(phydbl),BYTE_ALIGN);
- b->tPij_rr = _aligned_malloc(tree->mod->ras->n_catg*tree->mod->ns*tree->mod->ns*sizeof(phydbl),BYTE_ALIGN);
-#endif
-#else
- b->Pij_rr = (phydbl *)mCalloc(tree->mod->ras->n_catg*tree->mod->ns*tree->mod->ns,sizeof(phydbl));
- b->tPij_rr = (phydbl *)mCalloc(tree->mod->ras->n_catg*tree->mod->ns*tree->mod->ns,sizeof(phydbl));
-#endif
+/* #if (defined(__AVX__) || defined(__AVX2__) || defined(__SSE__) || defined(__SSE2__) || defined(__SSE3__)) */
+/* #ifndef WIN32 */
+/* if(posix_memalign((void *)&b->Pij_rr,BYTE_ALIGN,(size_t)tree->mod->ras->n_catg*tree->mod->ns*tree->mod->ns*sizeof(phydbl))) Generic_Exit(__FILE__,__LINE__,__FUNCTION__); */
+/* if(posix_memalign((void *)&b->tPij_rr,BYTE_ALIGN,(size_t)tree->mod->ras->n_catg*tree->mod->ns*tree->mod->ns*sizeof(phydbl))) Generic_Exit(__FILE__,__LINE__,__FUNCTION__); */
+/* #else */
+/* b->Pij_rr = _aligned_malloc(tree->mod->ras->n_catg*tree->mod->ns*tree->mod->ns*sizeof(phydbl),BYTE_ALIGN); */
+/* b->tPij_rr = _aligned_malloc(tree->mod->ras->n_catg*tree->mod->ns*tree->mod->ns*sizeof(phydbl),BYTE_ALIGN); */
+/* #endif */
+/* #else */
+/* b->Pij_rr = (phydbl *)mCalloc(tree->mod->ras->n_catg*tree->mod->ns*tree->mod->ns,sizeof(phydbl)); */
+/* b->tPij_rr = (phydbl *)mCalloc(tree->mod->ras->n_catg*tree->mod->ns*tree->mod->ns,sizeof(phydbl)); */
+/* #endif */
+
+
+ b->Pij_rr = tree->big_lk_array + tree->big_lk_array_pos;
+ tree->big_lk_array_pos += tree->mod->ras->n_catg*tree->mod->ns*tree->mod->ns;
+
+ b->tPij_rr = tree->big_lk_array + tree->big_lk_array_pos;
+ tree->big_lk_array_pos += tree->mod->ras->n_catg*tree->mod->ns*tree->mod->ns;
+
Make_Edge_Lk_Left(b,tree);
Make_Edge_Lk_Rght(b,tree);
@@ -334,16 +361,21 @@ void Make_Edge_Lk_Left(t_edge *b, t_tree *tree)
{
if((!b->left->tax) || (tree->mod->s_opt->greedy))
{
-#if (defined(__AVX__) || defined(__AVX2__) || defined(__SSE__) || defined(__SSE2__) || defined(__SSE3__))
-#ifndef WIN32
- if(posix_memalign((void **)&b->p_lk_left,BYTE_ALIGN,(size_t)tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns*sizeof(phydbl))) Generic_Exit(__FILE__,__LINE__,__FUNCTION__);
-#else
- b->p_lk_left = _aligned_malloc(tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns*sizeof(phydbl),BYTE_ALIGN);
-#endif
-#else
- b->p_lk_left = (phydbl *)mCalloc(tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns,sizeof(phydbl));
-#endif
+
+/* #if (defined(__AVX__) || defined(__AVX2__) || defined(__SSE__) || defined(__SSE2__) || defined(__SSE3__)) */
+/* #ifndef WIN32 */
+/* if(posix_memalign((void **)&b->p_lk_left,BYTE_ALIGN,(size_t)tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns*sizeof(phydbl))) Generic_Exit(__FILE__,__LINE__,__FUNCTION__); */
+/* #else */
+/* b->p_lk_left = _aligned_malloc(tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns*sizeof(phydbl),BYTE_ALIGN); */
+/* #endif */
+/* #else */
+/* b->p_lk_left = (phydbl *)mCalloc(tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns,sizeof(phydbl)); */
+/* #endif */
+
b->p_lk_tip_l = NULL;
+
+ b->p_lk_left = tree->big_lk_array + tree->big_lk_array_pos;
+ tree->big_lk_array_pos += tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns;
}
else if(b->left->tax)
{
@@ -370,15 +402,18 @@ void Make_Edge_Lk_Left(t_edge *b, t_tree *tree)
{
b->sum_scale_left = (int *)mCalloc(tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes),sizeof(int));
-#if (defined(__AVX__) || defined(__AVX2__) || defined(__SSE__) || defined(__SSE2__) || defined(__SSE3__))
-#ifndef WIN32
- if(posix_memalign((void **)&b->p_lk_left,BYTE_ALIGN,(size_t)tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns*sizeof(phydbl))) Generic_Exit(__FILE__,__LINE__,__FUNCTION__);
-#else
- b->p_lk_left = _aligned_malloc(tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns*sizeof(phydbl),BYTE_ALIGN);
-#endif
-#else
- b->p_lk_left = (phydbl *)mCalloc(tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns,sizeof(phydbl));
-#endif
+/* #if (defined(__AVX__) || defined(__AVX2__) || defined(__SSE__) || defined(__SSE2__) || defined(__SSE3__)) */
+/* #ifndef WIN32 */
+/* if(posix_memalign((void **)&b->p_lk_left,BYTE_ALIGN,(size_t)tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns*sizeof(phydbl))) Generic_Exit(__FILE__,__LINE__,__FUNCTION__); */
+/* #else */
+/* b->p_lk_left = _aligned_malloc(tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns*sizeof(phydbl),BYTE_ALIGN); */
+/* #endif */
+/* #else */
+/* b->p_lk_left = (phydbl *)mCalloc(tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns,sizeof(phydbl)); */
+/* #endif */
+
+ b->p_lk_left = tree->big_lk_array + tree->big_lk_array_pos;
+ tree->big_lk_array_pos += tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns;
}
b->patt_id_left = (int *)mCalloc(tree->data->crunch_len,sizeof(int));
@@ -405,16 +440,21 @@ void Make_Edge_Lk_Rght(t_edge *b, t_tree *tree)
{
if((!b->rght->tax) || (tree->mod->s_opt->greedy))
{
-#if (defined(__AVX__) || defined(__AVX2__) || defined(__SSE__) || defined(__SSE2__) || defined(__SSE3__))
-#ifndef WIN32
- if(posix_memalign((void **)&b->p_lk_rght,BYTE_ALIGN,(size_t)tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns*sizeof(phydbl))) Generic_Exit(__FILE__,__LINE__,__FUNCTION__);
-#else
- b->p_lk_rght = _aligned_malloc(tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns*sizeof(phydbl),BYTE_ALIGN);
-#endif
-#else
- b->p_lk_rght = (phydbl *)mCalloc(tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns,sizeof(phydbl));
-#endif
+/* #if (defined(__AVX__) || defined(__AVX2__) || defined(__SSE__) || defined(__SSE2__) || defined(__SSE3__)) */
+/* #ifndef WIN32 */
+/* if(posix_memalign((void **)&b->p_lk_rght,BYTE_ALIGN,(size_t)tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns*sizeof(phydbl))) Generic_Exit(__FILE__,__LINE__,__FUNCTION__); */
+/* #else */
+/* b->p_lk_rght = _aligned_malloc(tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns*sizeof(phydbl),BYTE_ALIGN); */
+/* #endif */
+/* #else */
+/* b->p_lk_rght = (phydbl *)mCalloc(tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns,sizeof(phydbl)); */
+/* #endif */
+
b->p_lk_tip_r = NULL;
+
+
+ b->p_lk_rght = tree->big_lk_array + tree->big_lk_array_pos;
+ tree->big_lk_array_pos += tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns;
}
else if(b->rght->tax)
{
@@ -439,17 +479,20 @@ void Make_Edge_Lk_Rght(t_edge *b, t_tree *tree)
if(b->num >= 2*tree->n_otu-3)
{
b->sum_scale_rght = (int *)mCalloc(tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes),sizeof(int));
-#if (defined(__AVX__) || defined(__AVX2__) || defined(__SSE__) || defined(__SSE2__) || defined(__SSE3__))
-#ifndef WIN32
- if(posix_memalign((void **)&b->p_lk_rght,BYTE_ALIGN,(size_t)tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns*sizeof(phydbl))) Generic_Exit(__FILE__,__LINE__,__FUNCTION__);
-#else
- b->p_lk_rght = _aligned_malloc(tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns*sizeof(phydbl),BYTE_ALIGN);
-#endif
+/* #if (defined(__AVX__) || defined(__AVX2__) || defined(__SSE__) || defined(__SSE2__) || defined(__SSE3__)) */
+/* #ifndef WIN32 */
+/* if(posix_memalign((void **)&b->p_lk_rght,BYTE_ALIGN,(size_t)tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns*sizeof(phydbl))) Generic_Exit(__FILE__,__LINE__,__FUNCTION__); */
+/* #else */
+/* b->p_lk_rght = _aligned_malloc(tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns*sizeof(phydbl),BYTE_ALIGN); */
+/* #endif */
-#else
- b->p_lk_rght = (phydbl *)mCalloc(tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns,sizeof(phydbl));
-#endif
+/* #else */
+/* b->p_lk_rght = (phydbl *)mCalloc(tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns,sizeof(phydbl)); */
+/* #endif */
+
+ b->p_lk_rght = tree->big_lk_array + tree->big_lk_array_pos;
+ tree->big_lk_array_pos += tree->data->crunch_len*MAX(tree->mod->ras->n_catg,tree->mod->n_mixt_classes)*tree->mod->ns;
}
b->patt_id_rght = (int *)mCalloc(tree->data->crunch_len,sizeof(int));
@@ -714,7 +757,6 @@ t_tree *Make_Tree(int n_otu)
//////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////
-
void Make_Tree_Path(t_tree *tree)
{
tree->curr_path = (t_node **)mCalloc(tree->n_otu,sizeof(t_node *));
@@ -1244,6 +1286,7 @@ void Make_Spr(t_tree *tree)
Make_Spr_List_One_Edge(tree);
Make_Spr_List_All_Edge(tree);
Make_Best_Spr(tree);
+ if(tree->is_mixt_tree == YES) MIXT_Repeat_Task(Make_Spr,tree);
}
//////////////////////////////////////////////////////////////
=====================================
src/mcmc.c
=====================================
@@ -1502,7 +1502,9 @@ void MCMC_Root_Time(t_tree *tree, int print)
Set_Both_Sides(NO,tree);
RATES_Update_Edge_Lengths(tree);
+#ifdef PHYREX
PHYREX_Update_Ldsk_Rates_Given_Edges(tree);
+#endif
if(tree->eval_alnL == YES) Lk(NULL,tree);
}
@@ -6124,7 +6126,7 @@ void MCMC_Complete_MCMC(t_mcmc *mcmc, t_tree *tree)
mcmc->num_move_phyrex_ldscape_lim = mcmc->n_moves; mcmc->n_moves += 1;
mcmc->num_move_phyrex_sigsq = mcmc->n_moves; mcmc->n_moves += 1;
mcmc->num_move_time_neff = mcmc->n_moves; mcmc->n_moves += 1;
- mcmc->num_move_time_neff_growth = mcmc->n_moves; mcmc->n_moves += 1;
+ mcmc->num_move_time_neff_growth = mcmc->n_moves; mcmc->n_moves += 1;
mcmc->num_move_phyrex_sim = mcmc->n_moves; mcmc->n_moves += 1;
mcmc->num_move_phyrex_traj = mcmc->n_moves; mcmc->n_moves += 1;
mcmc->num_move_phyrex_indel_disk_serial = mcmc->n_moves; mcmc->n_moves += 1;
@@ -8265,7 +8267,7 @@ void MCMC_PHYREX_Swap_Disk(t_tree *tree, int print)
/* Uniform selection of a valid disk */
target_disk = valid_disks[permut[j]];
- ori_time[j] = target_disk->time;
+ ori_time[j] = target_disk->time;
t_min = (target_disk->ldsk->prev != NULL) ? target_disk->ldsk->prev->disk->time : -TWO_TO_THE_LARGE;
t_max = +INFINITY;
=====================================
src/mixt.c
=====================================
@@ -51,6 +51,7 @@ void MIXT_Chain_All(t_tree *mixt_tree)
}
while(next);
+ MIXT_Chain_Models(mixt_tree);
MIXT_Chain_Edges(mixt_tree);
MIXT_Chain_Nodes(mixt_tree);
MIXT_Chain_Sprs(mixt_tree);
@@ -98,6 +99,50 @@ void MIXT_Chain_Edges(t_tree *mixt_tree)
//////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////
+void MIXT_Chain_Models(t_tree *mixt_tree)
+{
+ t_tree *tree;
+
+ tree = mixt_tree;
+ do
+ {
+ if(tree->next != NULL)
+ {
+ tree->mod->next = tree->next->mod;
+ tree->next->mod->prev = tree->mod;
+ }
+ else
+ {
+ tree->mod->next = NULL;
+ }
+ tree = tree->next;
+ }
+ while(tree);
+
+
+ tree = mixt_tree;
+ do
+ {
+ if(tree->next_mixt != NULL)
+ {
+ tree->mod->next_mixt = tree->next_mixt->mod;
+ tree->next_mixt->mod->prev_mixt = tree->mod;
+ }
+ else
+ {
+ tree->mod->next_mixt = NULL;
+ }
+
+ tree = tree->next_mixt;
+ }
+ while(tree);
+
+}
+
+//////////////////////////////////////////////////////////////
+//////////////////////////////////////////////////////////////
+
+
void MIXT_Chain_Nodes(t_tree *mixt_tree)
{
int i;
@@ -1595,7 +1640,6 @@ void MIXT_Multiply_Scalar_Dbl(scalar_dbl *this, phydbl scalar)
void MIXT_Br_Len_Opt(t_edge *mixt_b, t_tree *mixt_tree)
{
t_edge *b;
- t_tree *tree;
scalar_dbl *l;
MIXT_Turn_Branches_OnOff_In_All_Elem(ON,mixt_tree);
@@ -1603,33 +1647,14 @@ void MIXT_Br_Len_Opt(t_edge *mixt_b, t_tree *mixt_tree)
mixt_tree->ignore_mixt_info = YES;
b = mixt_b;
- tree = mixt_tree;
- l = NULL;
+ l = b->l;
do
{
- while(tree->is_mixt_tree == YES)
- {
- tree = tree->next;
- b = b->next;
- }
-
- l = b->l;
- do
- {
- if(l->onoff == ON)
- {
- Br_Len_Opt(&(l->v),mixt_b,mixt_tree);
- l->onoff = OFF;
- }
- l = l->next;
- }
- while(l);
-
- tree = tree->next;
- b = b->next;
+ if(l->onoff == ON) Br_Len_Opt(&(l->v),mixt_b,mixt_tree);
+ l = l->next;
}
- while(tree);
+ while(l);
mixt_tree->ignore_mixt_info = NO;
}
@@ -1750,13 +1775,17 @@ void MIXT_Set_Ignore_Root(int yesno, t_tree *mixt_tree)
{
t_tree *tree;
- tree = mixt_tree;
- do
+ tree = mixt_tree->next;
+
+ if(tree != NULL)
{
- tree->ignore_root = yesno;
- tree = tree->next;
+ do
+ {
+ Set_Ignore_Root(yesno,tree);
+ tree = tree->next;
+ }
+ while(tree);
}
- while(tree);
}
//////////////////////////////////////////////////////////////
@@ -2232,7 +2261,7 @@ void MIXT_RATES_Update_Edge_Lengths(t_tree *mixt_tree)
#if (defined PHYREX || PHYTIME)
RATES_Update_Edge_Lengths(tree);
#endif
-
+ if(tree->is_mixt_tree == YES) return;
tree = tree->next;
}
while(tree);
@@ -2278,36 +2307,20 @@ void MIXT_Update_Br_Len_Multipliers(t_mod *mod)
//////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////
-void MIXT_Init_Model(t_tree *mixt_tree)
+void MIXT_Init_Model(t_mod *mod)
{
- t_mod *mod,*mixt_mod;
- option *io;
- t_tree *tree;
-
- assert(mixt_tree);
-
- mixt_mod = mixt_tree->mod;
- io = mixt_tree->io;
-
- mod = mixt_mod;
- do
- {
- Init_Model(mod->io->cdata,mod,io);
- mod = mod->next;
- }
- while(mod);
+ mod = mod->next;
- tree = mixt_tree;
- do
+ if(mod != NULL)
{
- if(tree->next_mixt != NULL)
+ do
{
- tree->mod->next_mixt = tree->next_mixt->mod;
- tree->next_mixt->mod->prev_mixt = tree->mod;
+ Init_Model(mod->io->cdata,mod,mod->io);
+ if(mod->is_mixt_mod == YES) return;
+ mod = mod->next;
}
- tree = tree->next_mixt;
+ while(mod);
}
- while(tree);
}
//////////////////////////////////////////////////////////////
@@ -2410,17 +2423,18 @@ t_tree *MIXT_Starting_Tree(t_tree *mixt_tree)
void MIXT_Connect_Cseqs_To_Nodes(t_tree *mixt_tree)
{
t_tree *tree;
-
- Copy_Tree(mixt_tree,mixt_tree->next);
- tree = mixt_tree;
- do
+ tree = mixt_tree->next;
+
+ if(tree != NULL)
{
- Connect_CSeqs_To_Nodes(tree->data,mixt_tree->io,tree);
- tree = tree->next;
+ do
+ {
+ Connect_CSeqs_To_Nodes(tree->data,mixt_tree->io,tree);
+ tree = tree->next;
+ }
+ while(tree);
}
- while(tree);
-
}
//////////////////////////////////////////////////////////////
@@ -2429,15 +2443,10 @@ void MIXT_Connect_Cseqs_To_Nodes(t_tree *mixt_tree)
void MIXT_Init_T_Beg(t_tree *mixt_tree)
{
t_tree *tree;
-
+
/*! Initialize t_beg in each mixture tree */
- tree = mixt_tree;
- do
- {
- time(&(tree->t_beg));
- tree = tree->next_mixt;
- }
- while(tree);
+ tree = mixt_tree->next_mixt;
+ if(tree != NULL) Init_T_Beg(tree);
}
//////////////////////////////////////////////////////////////
@@ -2465,7 +2474,8 @@ void MIXT_Prepare_All(int num_rand_tree, t_tree *mixt_tree)
t_tree *tree;
MIXT_Check_Model_Validity(mixt_tree);
- MIXT_Init_Model(mixt_tree);
+ Init_Model(mixt_tree->mod->io->cdata,mixt_tree->mod,mixt_tree->mod->io);
+ Set_Model_Parameters(mixt_tree->mod);
Print_Data_Structure(NO,stdout,mixt_tree);
tree = MIXT_Starting_Tree(mixt_tree);
Copy_Tree(tree,mixt_tree);
@@ -2481,8 +2491,8 @@ void MIXT_Prepare_All(int num_rand_tree, t_tree *mixt_tree)
MIXT_Init_T_Beg(mixt_tree);
MIXT_Make_Tree_For_Pars(mixt_tree);
- MIXT_Make_Tree_For_Lk(mixt_tree);
- MIXT_Make_Spr(mixt_tree);
+ Make_Tree_For_Lk(mixt_tree);
+ Make_Spr(mixt_tree);
MIXT_Chain_All(mixt_tree);
MIXT_Check_Edge_Lens_In_All_Elem(mixt_tree);
@@ -2494,22 +2504,6 @@ void MIXT_Prepare_All(int num_rand_tree, t_tree *mixt_tree)
//////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////
-void MIXT_Make_Spr(t_tree *mixt_tree)
-{
- t_tree *tree;
-
- tree = mixt_tree;
- do
- {
- Make_Spr(tree);
- tree = tree->next;
- }
- while(tree);
-}
-
-//////////////////////////////////////////////////////////////
-//////////////////////////////////////////////////////////////
-
void MIXT_Make_Tree_For_Pars(t_tree *mixt_tree)
{
t_tree *tree;
@@ -3250,16 +3244,18 @@ void MIXT_Set_Both_Sides(int yesno, t_tree *mixt_tree)
assert(mixt_tree->is_mixt_tree == YES);
- tree = mixt_tree;
-
- do
+ tree = mixt_tree->next;
+
+ if(tree != NULL)
{
- if(tree->is_mixt_tree == YES) tree = tree->next;
- Set_Both_Sides(yesno,tree);
- tree = tree->next;
+ do
+ {
+ Set_Both_Sides(yesno,tree);
+ if(tree->is_mixt_tree == YES) return;
+ tree = tree->next;
+ }
+ while(tree);
}
- while(tree);
-
}
//////////////////////////////////////////////////////////////
@@ -3274,6 +3270,7 @@ void MIXT_Set_Model_Parameters(t_mod *mixt_mod)
do
{
Set_Model_Parameters(mod);
+ if(mod->is_mixt_mod == YES) return;
mod = mod->next;
}
while(mod);
@@ -3308,23 +3305,21 @@ void MIXT_Update_Eigen_Lr(t_edge *mixt_b, t_tree *mixt_tree)
t_edge *b;
t_tree *tree;
- tree = mixt_tree;
- b = mixt_b;
+ tree = mixt_tree->next;
+ b = mixt_b->next;
- do
+ if(tree != NULL)
{
- if(tree->is_mixt_tree == YES)
- {
+ do
+ {
+ Update_Eigen_Lr(b,tree);
+ if(tree->is_mixt_tree == YES) return;
+
tree = tree->next;
b = b->next;
}
-
- Update_Eigen_Lr(b,tree);
-
- tree = tree->next;
- b = b->next;
+ while(tree);
}
- while(tree);
}
@@ -3535,14 +3530,18 @@ void MIXT_Set_Bl_From_Rt(int yn, t_tree *mixt_tree)
{
t_tree *tree;
- tree = mixt_tree;
- do
+ tree = mixt_tree->next;
+
+ if(tree != NULL)
{
- assert(tree->rates);
- tree->rates->bl_from_rt = yn;
- tree = tree->next;
+ do
+ {
+ Set_Bl_From_Rt(yn,tree);
+ if(tree->is_mixt_tree == YES) return;
+ tree = tree->next;
+ }
+ while(tree);
}
- while(tree);
}
//////////////////////////////////////////////////////////////
@@ -3611,13 +3610,8 @@ void MIXT_Init_NNI_Score(phydbl val, t_edge *mixt_b, t_tree *mixt_tree)
b = mixt_b->next;
do
{
- if(tree->is_mixt_tree)
- {
- tree = tree->next;
- b = b->next;
- }
-
Init_NNI_Score(val,b,tree);
+ if(tree->is_mixt_tree == YES) return;
tree = tree->next;
b = b->next;
}
@@ -3630,39 +3624,87 @@ void MIXT_Init_NNI_Score(phydbl val, t_edge *mixt_b, t_tree *mixt_tree)
t_tree *MIXT_Duplicate_Tree(t_tree *ori)
{
- t_tree *cpy,*first,*prev;
- int dum;
+ t_tree *tree,*cpy;
+ int i;
+
+ assert(ori->is_mixt_tree == YES);
ori->is_mixt_tree = NO;
- first = Duplicate_Tree(ori);
- first->prev = NULL;
- first->is_mixt_tree = YES;
+ cpy = Duplicate_Tree(ori);
ori->is_mixt_tree = YES;
-
- ori = ori->next;
- prev = cpy = first;
+ cpy->is_mixt_tree = YES;
+ cpy->prev = NULL;
+ tree = ori;
+
do
{
- dum = ori->is_mixt_tree;
- ori->is_mixt_tree = NO;
-
- prev = cpy;
+ if(tree->next != NULL)
+ {
+ i = tree->next->is_mixt_tree;
+ tree->next->is_mixt_tree = NO;
+ cpy->next = Duplicate_Tree(tree->next);
+ tree->next->is_mixt_tree = i;
- cpy = Duplicate_Tree(ori);
+ cpy->next->is_mixt_tree = tree->next->is_mixt_tree;
+ cpy->next->prev = cpy;
+ }
+ else break;
+
+ tree = tree->next;
+ cpy = cpy->next;
+ }
+ while(tree != NULL);
- ori->is_mixt_tree = dum;
- cpy->is_mixt_tree = ori->is_mixt_tree;
+ do cpy = cpy->prev; while(cpy->prev);
- cpy->prev = prev;
- prev->next = cpy;
+ do
+ {
+ for(i=0;i<2*tree->n_otu-1;++i)
+ {
+ Free_Scalar_Dbl(cpy->a_edges[i]->l);
+ Free_Scalar_Dbl(cpy->a_edges[i]->l_old);
+ Free_Scalar_Dbl(cpy->a_edges[i]->l_var);
+ Free_Scalar_Dbl(cpy->a_edges[i]->l_var_old);
+ }
- ori = ori->next;
+ if(cpy->next == NULL) break;
+
+ cpy = cpy->next;
}
- while(ori);
+ while(1);
+
+ do cpy = cpy->prev; while(cpy->prev);
+
+ for(i=0;i<2*tree->n_otu-1;++i)
+ {
+ cpy->a_edges[i]->l = Duplicate_Scalar_Dbl(ori->a_edges[i]->l);
+ cpy->a_edges[i]->l_old = Duplicate_Scalar_Dbl(ori->a_edges[i]->l_old);
+ cpy->a_edges[i]->l_var = Duplicate_Scalar_Dbl(ori->a_edges[i]->l_var);
+ cpy->a_edges[i]->l_var_old = Duplicate_Scalar_Dbl(ori->a_edges[i]->l_var_old);
+ }
+
- cpy->next = NULL;
+ tree = cpy;
- return(first);
+ do
+ {
+ for(i=0;i<2*tree->n_otu-1;++i)
+ {
+ cpy->a_edges[i]->l = tree->a_edges[i]->l;
+ cpy->a_edges[i]->l_old = tree->a_edges[i]->l_old;
+ cpy->a_edges[i]->l_var = tree->a_edges[i]->l_var;
+ cpy->a_edges[i]->l_var_old = tree->a_edges[i]->l_var_old;
+ }
+
+ if(cpy->next == NULL) break;
+
+ cpy = cpy->next;
+ }
+ while(1);
+
+ do cpy = cpy->prev; while(cpy->prev);
+
+ return(cpy);
}
//////////////////////////////////////////////////////////////
@@ -3828,6 +3870,7 @@ void MIXT_Exchange_Nodes(t_node *a, t_node *d, t_node *w, t_node *v, t_tree *mix
tree->a_nodes[w->num],
tree->a_nodes[v->num],
tree);
+ if(tree->is_mixt_tree == YES) return;
i++;
tree = tree->next;
}
@@ -3836,8 +3879,69 @@ void MIXT_Exchange_Nodes(t_node *a, t_node *d, t_node *w, t_node *v, t_tree *mix
//////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////
+
+void MIXT_Free_Tree(t_tree *mixt_tree)
+{
+ t_tree *tree;
+ t_tree *next;
+
+ Free_All_Edges_Lens(mixt_tree);
+
+ tree = mixt_tree;
+ do
+ {
+ if(tree->mat) Free_Mat(tree->mat);
+ if(tree->t_dir) Free(tree->t_dir);
+ if(tree->short_l) Free(tree->short_l);
+ if(tree->mutmap) Free(tree->mutmap);
+ Free_Bip(tree);
+ Free(tree->curr_path);
+ Free_All_Nodes_Light(tree);
+ Free_All_Edges_Light(tree);
+ tree = tree->next;
+ }
+ while(tree);
+
+
+ tree = mixt_tree;
+ next = mixt_tree->next;
+ do
+ {
+ Free(tree);
+ tree = next;
+ if(!tree) break;
+ next = next->next;
+ }
+ while(tree);
+}
+
//////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////
+void MIXT_Repeat_Task(void (*Task_Function)(), t_tree *mixt_tree)
+{
+ t_tree *tree,*next;
+
+ tree = mixt_tree->next;
+
+ if(tree != NULL)
+ {
+ do
+ {
+ next = tree->next;
+ (*Task_Function)(tree);
+ if(tree->is_mixt_tree == YES) return;
+ tree = next;
+ }
+ while(tree != NULL);
+ }
+}
+
+//////////////////////////////////////////////////////////////
+//////////////////////////////////////////////////////////////
+//////////////////////////////////////////////////////////////
+//////////////////////////////////////////////////////////////
+//////////////////////////////////////////////////////////////
+//////////////////////////////////////////////////////////////
=====================================
src/mixt.h
=====================================
@@ -71,7 +71,7 @@ phydbl MIXT_Rescale_Free_Rate_Tree(t_tree *mixt_tree);
void MIXT_Set_Br_Len_Var(t_edge *mixt_b, t_tree *mixt_tree);
void MIXT_Optimize_Br_Len_Multiplier(t_tree *mixt_tree);
void MIXT_Update_Br_Len_Multipliers(t_mod *mod);
-void MIXT_Init_Model(t_tree *mixt_tree);
+void MIXT_Init_Model(t_mod *mod);
t_tree *MIXT_Starting_Tree(t_tree *mixt_tree);
void MIXT_Connect_Cseqs_To_Nodes(t_tree *mixt_tree);
void MIXT_Check_Edge_Lens_In_One_Elem(t_tree *mixt_tree);
@@ -118,5 +118,8 @@ t_tree *MIXT_Duplicate_Tree(t_tree *ori);
void MIXT_Set_Model_Parameters(t_mod *mixt_mod);
void MIXT_Print_Site_Lk(t_tree *mixt_tree, FILE *fp);
void MIXT_Exchange_Nodes(t_node *a, t_node *d, t_node *w, t_node *v, t_tree *mixt_tree);
+void MIXT_Chain_Models(t_tree *mixt_tree);
+void MIXT_Repeat_Task(void (*Task_Function)(),t_tree *mixt_tree);
+void MIXT_Free_Tree(t_tree *mixt_tree);
#endif
=====================================
src/models.c
=====================================
@@ -780,7 +780,7 @@ int Update_Efrq(t_mod *mod)
//////////////////////////////////////////////////////////////
int Set_Model_Parameters(t_mod *mod)
-{
+{
if(!Update_Boundaries(mod)) return 0;
if(!Update_RAS(mod)) return 0;
if(!Update_Efrq(mod)) return 0;
@@ -829,6 +829,7 @@ int Update_Boundaries(t_mod *mod)
for(i=0;i<6;i++) if(mod->r_mat->rr->v[i] > RR_MAX) mod->r_mat->rr->v[i] = RR_MAX;
}
+
for(i=0;i<mod->ns;i++)
{
if(mod->e_frq->pi_unscaled->v[i] < UNSCALED_E_FRQ_MIN)
@@ -836,7 +837,7 @@ int Update_Boundaries(t_mod *mod)
if(mod->e_frq->pi_unscaled->v[i] > UNSCALED_E_FRQ_MAX)
mod->e_frq->pi_unscaled->v[i] = UNSCALED_E_FRQ_MAX;
-
+
if(mod->e_frq->pi->v[i] < E_FRQ_MIN)
mod->e_frq->pi->v[i] = E_FRQ_MIN;
=====================================
src/mpi_boot.c
=====================================
@@ -187,8 +187,10 @@ void Bootstrap_MPI(t_tree *tree)
requires to leave the value of s_opt unchanged during the boostrap. */
boot_mod->io = tree->io; /* WARNING: re-using the same address here instead of creating a copying
requires to leave the value of io unchanged during the boostrap. */
- Init_Model(boot_data,boot_mod,tree->io);
+ Init_Model(boot_data,boot_mod,tree->io);
+ Set_Model_Parameters(boot_mod);
+
if(tree->io->in_tree == 2)
{
switch(tree->io->tree_file_format)
=====================================
src/optimiz.c
=====================================
@@ -612,7 +612,6 @@ phydbl Br_Len_Opt(phydbl *l, t_edge *b, t_tree *tree)
{
phydbl lk_begin, lk_end;
-
if(tree->is_mixt_tree == YES && tree->ignore_mixt_info == NO)
{
MIXT_Br_Len_Opt(b,tree);
@@ -718,16 +717,6 @@ void Optimize_Br_Len_Serie(int n_max_iter, t_tree *tree)
{
phydbl lk_init,lk_end;
int iter;
-
-
- /* if(tree->n_tot_bl_opt == 0) */
- /* { */
- /* for(int i=0;i<2*tree->n_otu-3;++i) Set_Scalar_Dbl(1.E-6,tree->a_edges[i]->l); */
- /* for(int i=0;i<2*tree->n_otu-3;++i) Set_Scalar_Dbl_Min_Thresh(tree->mod->l_min,tree->a_edges[i]->l); */
- /* for(int i=0;i<2*tree->n_otu-3;++i) Set_Scalar_Dbl_Max_Thresh(tree->mod->l_max,tree->a_edges[i]->l); */
- /* } */
-
-
Set_Both_Sides(NO,tree);
Lk(NULL,tree);
@@ -795,12 +784,7 @@ void Optimize_Br_Len_Serie(int n_max_iter, t_tree *tree)
}
iter++;
-
-
- /* PhyML_Printf("\n. lnL: %f %f %f %d",tree->c_lnL,Rgamma((phydbl)(iter+1),(phydbl)(1./(iter+1))),(phydbl)(iter+1)*(phydbl)(1./(iter+1))*(phydbl)(1./(iter+1)),iter); */
-
}
- /* while(lk_end - lk_init > tree->mod->s_opt->min_diff_lk_local && iter < n_max_iter); */
while(iter < 1);
}
=====================================
src/phyrex.c
=====================================
@@ -455,17 +455,21 @@ void PHYREX_XML(char *xml_filename)
mixt_tree->io->r_seed = seed;
MIXT_Check_Model_Validity(mixt_tree);
- MIXT_Init_Model(mixt_tree);
+ MIXT_Chain_Models(mixt_tree);
+ Init_Model(mixt_tree->mod->io->cdata,mixt_tree->mod,mixt_tree->mod->io);
+ Set_Model_Parameters(mixt_tree->mod);
Print_Data_Structure(NO,stdout,mixt_tree);
tree = MIXT_Starting_Tree(mixt_tree);
if(mixt_tree->io->in_tree < 2) Add_Root(tree->a_edges[0],tree);
Copy_Tree(tree,mixt_tree);
Free_Tree(tree);
- MIXT_Connect_Cseqs_To_Nodes(mixt_tree);
- MIXT_Init_T_Beg(mixt_tree);
- MIXT_Make_Tree_For_Lk(mixt_tree);
- MIXT_Make_Tree_For_Pars(mixt_tree);
- MIXT_Make_Spr(mixt_tree);
+ Copy_Tree(mixt_tree,mixt_tree->next);
+ Connect_CSeqs_To_Nodes(mixt_tree->mod->io->cdata,mixt_tree->mod->io,mixt_tree);
+ Init_T_Beg(mixt_tree);
+ Make_Tree_For_Lk(mixt_tree);
+ Make_Tree_For_Pars(mixt_tree);
+ Make_Spr(mixt_tree);
+
MIXT_Chain_All(mixt_tree);
MIXT_Check_Edge_Lens_In_All_Elem(mixt_tree);
MIXT_Turn_Branches_OnOff_In_All_Elem(ON,mixt_tree);
@@ -614,8 +618,8 @@ void PHYREX_XML(char *xml_filename)
MCMC_Randomize_Clock_Rate(mixt_tree);
MCMC_Randomize_Sigsq_Scale(mixt_tree);
- MIXT_Set_Ignore_Root(YES,mixt_tree);
- MIXT_Set_Bl_From_Rt(YES,mixt_tree);
+ Set_Ignore_Root(YES,mixt_tree);
+ Set_Bl_From_Rt(YES,mixt_tree);
PHYREX_Oldest_Sampled_Disk(mixt_tree);
for(int i=0;i<3;++i) if(mixt_tree->aux_tree && mixt_tree->aux_tree[i]) PHYREX_Oldest_Sampled_Disk(mixt_tree->aux_tree[i]);
=====================================
src/slfv.c
=====================================
@@ -1771,7 +1771,8 @@ t_tree *SLFV_Simulate(int n_otu, int n_sites, phydbl w, phydbl h, phydbl lbda,
RATES_Update_Edge_Lengths(tree);
Init_Model(cdata,mod,io);
-
+ Set_Model_Parameters(mod);
+
tree->mod->ras->n_catg = 1;
tree->mod->whichmodel = HKY85;
tree->mod->kappa->v = 4.0;
=====================================
src/spr.c
=====================================
@@ -128,8 +128,6 @@ int Spr(phydbl init_lnL, phydbl prop_spr, t_tree *tree)
Set_Both_Sides(YES,tree);
Lk(NULL,tree);
tree->best_lnL = tree->c_lnL;
-
- /* PhyML_Printf("\n. init: %f",tree->c_lnL); */
for(i=0;i<MAX(1,(int)((2*tree->n_otu-3)*prop_spr));++i)
{
@@ -145,7 +143,6 @@ int Spr(phydbl init_lnL, phydbl prop_spr, t_tree *tree)
Spr_Subtree(b,b->rght,tree);
}
}
-
Free(br_idx);
@@ -750,7 +747,6 @@ void Global_Spr_Search(t_tree *tree)
if(tree->io->print_json_trace == YES) JSON_Tree_Io(tree,tree->io->fp_out_json_trace);
if(tree->verbose > VL0 && tree->io->quiet == NO) PhyML_Printf("\n\n. Score of initial tree: %.2f",tree->c_lnL);
-
tree->mod->s_opt->min_diff_lk_move = 1.E-1;
tree->mod->s_opt->min_diff_lk_local = 1.E-1;
@@ -781,7 +777,6 @@ void Global_Spr_Search(t_tree *tree)
iter = 0;
do
{
-
if(!(iter%freq))
{
for(int i=0;i<2*tree->n_otu-3;++i) MIXT_Multiply_Scalar_Dbl(tree->a_edges[i]->l,Rgamma((phydbl)(tune_l_mult*iter+1),(phydbl)(1./(tune_l_mult*iter+1))));
@@ -864,8 +859,6 @@ void Global_Spr_Search(t_tree *tree)
tree->mod->s_opt->spr_pars = NO;
tree->mod->s_opt->min_diff_lk_move = 1.E-2;
tree->mod->s_opt->min_diff_lk_local = 1.E-2;
- tree->mod->s_opt->apply_spr_right_away = YES;
- tree->mod->s_opt->apply_spr = YES;
tree->mod->s_opt->eval_list_regraft = NO;
tree->mod->s_opt->max_delta_lnL_spr_current = 0.0;
tree->mod->s_opt->min_n_triple_moves = 1;
@@ -881,7 +874,7 @@ void Global_Spr_Search(t_tree *tree)
{
for(int i=0;i<2*tree->n_otu-3;++i) MIXT_Multiply_Scalar_Dbl(tree->a_edges[i]->l,Rgamma((phydbl)(tune_l_mult*iter+1),(phydbl)(1./(tune_l_mult*iter+1))));
for(int i=0;i<2*tree->n_otu-3;++i) Set_Scalar_Dbl_Min_Thresh(tree->mod->l_min,tree->a_edges[i]->l);
- for(int i=0;i<2*tree->n_otu-3;++i) Set_Scalar_Dbl_Max_Thresh(tree->mod->l_max,tree->a_edges[i]->l);
+ for(int i=0;i<2*tree->n_otu-3;++i) Set_Scalar_Dbl_Max_Thresh(tree->mod->l_max,tree->a_edges[i]->l);
freq--;
if(freq < 1) freq=1;
}
@@ -890,7 +883,7 @@ void Global_Spr_Search(t_tree *tree)
Optimize_Br_Len_Serie(2,tree);
if(!(iter%round_freq) && iter > 0) Round_Optimize(tree,1000);
-
+
if(tree->verbose > VL0 && tree->io->quiet == NO)
{
time(&t_cur);
@@ -952,8 +945,6 @@ void Global_Spr_Search(t_tree *tree)
tree->mod->s_opt->l_min_spr = 1.E-4;
tree->mod->s_opt->min_diff_lk_move = 1.E-2;
tree->mod->s_opt->min_diff_lk_local = 1.E-2;
- tree->mod->s_opt->apply_spr_right_away = YES;
- tree->mod->s_opt->apply_spr = YES;
tree->mod->s_opt->eval_list_regraft = YES;
tree->mod->s_opt->max_delta_lnL_spr_current = 0.0;
tree->mod->s_opt->min_n_triple_moves = 5;
@@ -964,7 +955,7 @@ void Global_Spr_Search(t_tree *tree)
{
Spr(tree->c_lnL,1.0,tree);
Optimize_Br_Len_Serie(2,tree);
-
+
if(!(iter%round_freq)) Round_Optimize(tree,1000);
if(tree->verbose > VL0 && tree->io->quiet == NO)
=====================================
src/utilities.c
=====================================
@@ -3181,7 +3181,8 @@ void Bootstrap(t_tree *tree)
requires to leave the value of io unchanged during the boostrap. */
Init_Model(boot_data,boot_mod,tree->io);
-
+ Set_Model_Parameters(boot_mod);
+
if(tree->io->in_tree == 2)
{
rewind(tree->io->fp_in_tree);
@@ -4947,102 +4948,106 @@ void Copy_Tree(t_tree *ori, t_tree *cpy)
MIXT_Copy_Tree(ori,cpy);
return;
}
-
- for(i=0;i<2*ori->n_otu-1;++i)
- {
- if(ori->a_nodes[i] != NULL)
- {
- cpy->a_nodes[i]->anc =
- (ori->a_nodes[i]->anc != NULL) ?
- cpy->a_nodes[ori->a_nodes[i]->anc->num] :
- NULL;
-
- for(j=0;j<3;++j)
+ else
+ {
+ for(i=0;i<2*ori->n_otu-1;++i)
+ {
+ if(ori->a_nodes[i] != NULL)
{
- if(ori->a_nodes[i]->v[j] != NULL)
- {
- cpy->a_nodes[i]->v[j] = cpy->a_nodes[ori->a_nodes[i]->v[j]->num];
- cpy->a_nodes[i]->b[j] = cpy->a_edges[ori->a_nodes[i]->b[j]->num];
- }
- else
+ cpy->a_nodes[i]->anc =
+ (ori->a_nodes[i]->anc != NULL) ?
+ cpy->a_nodes[ori->a_nodes[i]->anc->num] :
+ NULL;
+
+ for(j=0;j<3;++j)
{
- cpy->a_nodes[i]->v[j] = NULL;
- cpy->a_nodes[i]->b[j] = NULL;
+ if(ori->a_nodes[i]->v[j] != NULL)
+ {
+ cpy->a_nodes[i]->v[j] = cpy->a_nodes[ori->a_nodes[i]->v[j]->num];
+ cpy->a_nodes[i]->b[j] = cpy->a_edges[ori->a_nodes[i]->b[j]->num];
+ }
+ else
+ {
+ cpy->a_nodes[i]->v[j] = NULL;
+ cpy->a_nodes[i]->b[j] = NULL;
+ }
}
}
+ cpy->a_nodes[i]->c_seq = ori->a_nodes[i]->c_seq;
}
- cpy->a_nodes[i]->c_seq = ori->a_nodes[i]->c_seq;
- }
-
- for(i=0;i<2*ori->n_otu-1;++i)
- {
- if(ori->a_edges[i] != NULL)
+
+ for(i=0;i<2*ori->n_otu-1;++i)
{
- cpy->a_edges[i]->l->v = ori->a_edges[i]->l->v;
- cpy->a_edges[i]->l_old->v = ori->a_edges[i]->l_old->v;
- cpy->a_edges[i]->l_var->v = ori->a_edges[i]->l_var->v;
- cpy->a_edges[i]->l_var_old->v = ori->a_edges[i]->l_var_old->v;
- cpy->a_edges[i]->left = ori->a_edges[i]->left ? cpy->a_nodes[ori->a_edges[i]->left->num] : NULL;
- cpy->a_edges[i]->rght = ori->a_edges[i]->rght ? cpy->a_nodes[ori->a_edges[i]->rght->num] : NULL;
- cpy->a_edges[i]->l_v1 = ori->a_edges[i]->l_v1;
- cpy->a_edges[i]->l_v2 = ori->a_edges[i]->l_v2;
- cpy->a_edges[i]->r_v1 = ori->a_edges[i]->r_v1;
- cpy->a_edges[i]->r_v2 = ori->a_edges[i]->r_v2;
- cpy->a_edges[i]->l_r = ori->a_edges[i]->l_r;
- cpy->a_edges[i]->r_l = ori->a_edges[i]->r_l;
- cpy->a_edges[i]->does_exist = ori->a_edges[i]->does_exist;
- cpy->a_edges[i]->support_val = ori->a_edges[i]->support_val;
-
+ if(ori->a_edges[i] != NULL)
+ {
+ cpy->a_edges[i]->l->v = ori->a_edges[i]->l->v;
+ cpy->a_edges[i]->l_old->v = ori->a_edges[i]->l_old->v;
+ cpy->a_edges[i]->l_var->v = ori->a_edges[i]->l_var->v;
+ cpy->a_edges[i]->l_var_old->v = ori->a_edges[i]->l_var_old->v;
+ cpy->a_edges[i]->left = ori->a_edges[i]->left ? cpy->a_nodes[ori->a_edges[i]->left->num] : NULL;
+ cpy->a_edges[i]->rght = ori->a_edges[i]->rght ? cpy->a_nodes[ori->a_edges[i]->rght->num] : NULL;
+ cpy->a_edges[i]->l_v1 = ori->a_edges[i]->l_v1;
+ cpy->a_edges[i]->l_v2 = ori->a_edges[i]->l_v2;
+ cpy->a_edges[i]->r_v1 = ori->a_edges[i]->r_v1;
+ cpy->a_edges[i]->r_v2 = ori->a_edges[i]->r_v2;
+ cpy->a_edges[i]->l_r = ori->a_edges[i]->l_r;
+ cpy->a_edges[i]->r_l = ori->a_edges[i]->r_l;
+ cpy->a_edges[i]->does_exist = ori->a_edges[i]->does_exist;
+ cpy->a_edges[i]->support_val = ori->a_edges[i]->support_val;
+
#ifdef BEAGLE
- cpy->a_edges[i]->p_lk_left_idx = ori->a_edges[i]->p_lk_left_idx;
- cpy->a_edges[i]->p_lk_rght_idx = ori->a_edges[i]->p_lk_rght_idx;
- cpy->a_edges[i]->p_lk_tip_idx = ori->a_edges[i]->p_lk_tip_idx;
+ cpy->a_edges[i]->p_lk_left_idx = ori->a_edges[i]->p_lk_left_idx;
+ cpy->a_edges[i]->p_lk_rght_idx = ori->a_edges[i]->p_lk_rght_idx;
+ cpy->a_edges[i]->p_lk_tip_idx = ori->a_edges[i]->p_lk_tip_idx;
#endif
+ }
}
- }
- for(i=0;i<ori->n_otu;++i)
- {
- cpy->a_nodes[i]->tax = YES;
-
- Free(cpy->a_nodes[i]->name);
-
- cpy->a_nodes[i]->name = (char *)mCalloc(strlen(ori->a_nodes[i]->name)+1,sizeof(char));
- cpy->a_nodes[i]->ori_name = cpy->a_nodes[i]->name ;
-
- strcpy(cpy->a_nodes[i]->name,ori->a_nodes[i]->name);
- }
-
-
- if(ori->n_root)
- {
- cpy->e_root = cpy->a_edges[ori->e_root->num];
- cpy->n_root = cpy->a_nodes[ori->n_root->num];
- cpy->n_root_pos = ori->n_root_pos;
+ for(i=0;i<ori->n_otu;++i)
+ {
+ cpy->a_nodes[i]->tax = YES;
+
+ Free(cpy->a_nodes[i]->name);
+
+ cpy->a_nodes[i]->name = (char *)mCalloc(strlen(ori->a_nodes[i]->name)+1,sizeof(char));
+ cpy->a_nodes[i]->ori_name = cpy->a_nodes[i]->name ;
+
+ strcpy(cpy->a_nodes[i]->name,ori->a_nodes[i]->name);
+ }
+
+
+ if(ori->n_root)
+ {
+ cpy->e_root = cpy->a_edges[ori->e_root->num];
+ cpy->n_root = cpy->a_nodes[ori->n_root->num];
+ cpy->n_root_pos = ori->n_root_pos;
+
+ cpy->n_root->b[1] = cpy->a_edges[ori->n_root->b[1]->num];
+ cpy->n_root->b[2] = cpy->a_edges[ori->n_root->b[2]->num];
+ }
+
+ cpy->num_curr_branch_available = 0;
+ cpy->t_beg = ori->t_beg;
+ cpy->verbose = ori->verbose;
- cpy->n_root->b[1] = cpy->a_edges[ori->n_root->b[1]->num];
- cpy->n_root->b[2] = cpy->a_edges[ori->n_root->b[2]->num];
- }
-
- cpy->num_curr_branch_available = 0;
- cpy->t_beg = ori->t_beg;
- cpy->verbose = ori->verbose;
-
#ifdef BEAGLE
- cpy->b_inst = ori->b_inst;
+ cpy->b_inst = ori->b_inst;
#endif
+ }
}
-
//////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////
t_tree *Duplicate_Tree(t_tree *ori)
{
if(ori->is_mixt_tree == YES) return MIXT_Duplicate_Tree(ori);
- t_tree *cpy = Make_Tree_From_Scratch(ori->n_otu,ori->data);
- Copy_Tree(ori,cpy);
- return(cpy);
+ else
+ {
+ t_tree *cpy = Make_Tree_From_Scratch(ori->n_otu,ori->data);
+ Copy_Tree(ori,cpy);
+ return(cpy);
+ }
}
//////////////////////////////////////////////////////////////
@@ -10273,6 +10278,34 @@ char *To_Upper_String(char *in)
//////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////
+void Set_Ignore_Root(int yesno, t_tree *tree)
+{
+ tree->ignore_root = yesno;
+ if(tree->is_mixt_tree == YES) MIXT_Set_Ignore_Root(yesno,tree);
+}
+
+//////////////////////////////////////////////////////////////
+//////////////////////////////////////////////////////////////
+
+void Set_Bl_From_Rt(int yesno, t_tree *tree)
+{
+ assert(tree->rates);
+ tree->rates->bl_from_rt = yesno;
+ if(tree->is_mixt_tree == YES) MIXT_Set_Bl_From_Rt(yesno,tree);
+}
+
+//////////////////////////////////////////////////////////////
+//////////////////////////////////////////////////////////////
+
+void Init_T_Beg(t_tree *tree)
+{
+ time(&(tree->t_beg));
+ if(tree->is_mixt_tree == YES) MIXT_Init_T_Beg(tree);
+}
+
+//////////////////////////////////////////////////////////////
+//////////////////////////////////////////////////////////////
+
void Connect_CSeqs_To_Nodes(calign *cdata, option *io, t_tree *tree)
{
int i,j,n_otu_tree,n_otu_cdata;
@@ -10302,6 +10335,8 @@ void Connect_CSeqs_To_Nodes(calign *cdata, option *io, t_tree *tree)
}
tree->a_nodes[i]->c_seq = cdata->c_seq[j];
}
+
+ if(tree->is_mixt_tree == YES) MIXT_Connect_Cseqs_To_Nodes(tree);
}
//////////////////////////////////////////////////////////////
@@ -11217,7 +11252,7 @@ void Set_Scalar_Dbl_Min_Thresh(phydbl thresh, scalar_dbl *from)
f = from;
do
{
- if(f->v < thresh) f->v = thresh;;
+ if(f->v < thresh) f->v = thresh;
f = f->next;
}
while(f);
@@ -11233,7 +11268,7 @@ void Set_Scalar_Dbl_Max_Thresh(phydbl thresh, scalar_dbl *from)
f = from;
do
{
- if(f->v > thresh) f->v = thresh;;
+ if(f->v > thresh) f->v = thresh;
f = f->next;
}
while(f);
=====================================
src/utilities.h
=====================================
@@ -757,6 +757,8 @@ typedef struct __Tree{
#endif
phydbl *p_lk_left_pi,*l_ev;
+ phydbl *big_lk_array;
+ int big_lk_array_pos;
short int eval_alnL; /*! Evaluate likelihood for genetic data */
short int eval_rlnL; /*! Evaluate likelihood for rates along the tree */
@@ -2477,6 +2479,9 @@ void Reset_Prior(t_tree *tree);
void Set_Prior(t_tree *tree);
t_edge *Get_Edge(t_node *a, t_node *d, t_tree *tree);
void Exchange_Nodes(t_node *a, t_node *d, t_node *w, t_node *v, t_tree *tree);
+void Init_T_Beg(t_tree *tree);
+void Set_Ignore_Root(int yesno, t_tree *tree);
+void Set_Bl_From_Rt(int yesno, t_tree *tree);
#include "xml.h"
=====================================
src/xml.c
=====================================
@@ -462,15 +462,7 @@ t_tree *XML_Process_Base(char *xml_filename)
tree->next = mixt_tree;
mixt_tree->prev = tree;
}
-
- /*! Do the same for the model
- */
- if(mod)
- {
- mod->next = iomod;
- iomod->prev = mod;
- }
-
+
if(!root_tree) root_tree = mixt_tree;
/*! Tree size scaling factor */
@@ -1036,6 +1028,7 @@ t_tree *XML_Process_Base(char *xml_filename)
}
else
{
+ /* Free existing edge lengths, which will be replaced with those read from <branchlengths> !*/
for(i=0;i<2*tree->n_otu-1;++i)
{
Free_Scalar_Dbl(tree->a_edges[i]->l);
@@ -1068,7 +1061,7 @@ t_tree *XML_Process_Base(char *xml_filename)
Exit("\n");
}
- For(i,2*tree->n_otu-1)
+ for(i=0;i<2*tree->n_otu-1;++i)
{
tree->a_edges[i]->l = lens[i];
mixt_tree->a_edges[i]->l = lens[i];
@@ -1089,7 +1082,7 @@ t_tree *XML_Process_Base(char *xml_filename)
if(first_m_elem > 1) // Done with this component, move to the next tree and model
{
if(tree->next) tree = tree->next;
- if(mod->next) mod = mod->next;
+ if(mod->next) mod = mod->next;
}
}
else if(list[j] != ' ')
@@ -1118,13 +1111,13 @@ t_tree *XML_Process_Base(char *xml_filename)
/*! Finish making the models */
- mod = mixt_tree->mod;
+ tree = mixt_tree;
do
{
- Make_Model_Complete(mod);
- mod = mod->next;
+ Make_Model_Complete(tree->mod);
+ tree = tree->next;
}
- while(mod);
+ while(tree);
Check_Taxa_Sets(mixt_tree);
View it on GitLab: https://salsa.debian.org/med-team/phyml/-/compare/2dbaad0f1e21dfda9016e20d0af8dabd2247a43e...55932608de89ccef431f9d5e2bd5adf20c8822c6
--
View it on GitLab: https://salsa.debian.org/med-team/phyml/-/compare/2dbaad0f1e21dfda9016e20d0af8dabd2247a43e...55932608de89ccef431f9d5e2bd5adf20c8822c6
You're receiving this email because of your account on salsa.debian.org.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://alioth-lists.debian.net/pipermail/debian-med-commit/attachments/20220211/232f288f/attachment-0001.htm>
More information about the debian-med-commit
mailing list