pll.h
Go to the documentation of this file.
1 
50 #ifndef __pll__
51 #define __pll__
52 
53 #include <stdint.h>
54 #include <stdio.h>
55 #include <errno.h>
56 
57 #ifdef __cplusplus
58 extern "C" {
59 #endif
60 
61 #ifdef __AVX
62 
63 #include <xmmintrin.h>
64 #include <immintrin.h>
65 #include <pmmintrin.h>
66 
67 #define PLL_BYTE_ALIGNMENT 32
68 
69 #else
70 
71 #ifdef __SSE3
72 
73 #include <xmmintrin.h>
74 #include <pmmintrin.h>
75 
76 #define PLL_BYTE_ALIGNMENT 16
77 
78 #else
79 #define PLL_BYTE_ALIGNMENT 1
80 #endif
81 #endif
82 
83 
84 #include "genericParallelization.h"
85 #include "errcodes.h"
86 #include "stack.h"
87 #include "queue.h"
88 #include "hash.h"
89 #include "newick.h"
90 #include "lexer.h"
91 #include "parsePartition.h"
92 #include "mem_alloc.h"
93 
94 #define PLL_MAX_TIP_EV 0.999999999 /* max tip vector value, sum of EVs needs to be smaller than 1.0, otherwise the numerics break down */
95 #define PLL_MAX_LOCAL_SMOOTHING_ITERATIONS 32
96 #define PLL_ITERATIONS 10 /* maximum iterations of iterations per insert */
97 #define PLL_NEWZPERCYCLE 1 /* iterations of makenewz per tree traversal */
98 #define PLL_NMLNGTH 256 /* number of characters in species name */
99 #define PLL_DELTAZ 0.00001 /* test of net branch length change in update */
100 #define PLL_DEFAULTZ 0.9 /* value of z assigned as starting point */
101 #define PLL_UNLIKELY -1.0E300 /* low likelihood for initialization */
102 
103 
104 #define PLL_SUMMARIZE_LENGTH -3
105 #define PLL_SUMMARIZE_LH -2
106 #define PLL_NO_BRANCHES -1
107 
108 #define PLL_MASK_LENGTH 32
109 #define GET_BITVECTOR_LENGTH(x) ((x % PLL_MASK_LENGTH) ? (x / PLL_MASK_LENGTH + 1) : (x / PLL_MASK_LENGTH))
110 
111 #define PLL_ZMIN 1.0E-15 /* max branch prop. to -log(PLL_ZMIN) (= 34) */
112 #define PLL_ZMAX (1.0 - 1.0E-6) /* min branch prop. to 1.0-zmax (= 1.0E-6) */
113 
114 #define PLL_TWOTOTHE256 \
115  115792089237316195423570985008687907853269984665640564039457584007913129639936.0
116  /* 2**256 (exactly) */
117 
118 #define PLL_MINLIKELIHOOD (1.0/PLL_TWOTOTHE256)
119 #define PLL_MINUSMINLIKELIHOOD -PLL_MINLIKELIHOOD
120 
121 
122 #define PLL_DEEP_COPY 1 << 0
123 #define PLL_SHALLOW_COPY 1 << 1
124 
125 
126 #define PLL_FORMAT_PHYLIP 1
127 #define PLL_FORMAT_FASTA 2
128 #define PLL_FORMAT_NEWICK 3
129 
130 #define PLL_NNI_P_NEXT 1
131 #define PLL_NNI_P_NEXTNEXT 2
133 /* 18446744073709551616.0 */
134 
135 /*4294967296.0*/
136 
137 /* 18446744073709551616.0 */
138 
139 /* 2**64 (exactly) */
140 /* 4294967296 2**32 */
141 
142 #define PLL_BADREAR -1
143 
144 #define PLL_NUM_BRANCHES 16
145 
146 #define PLL_TRUE 1
147 #define PLL_FALSE 0
148 
149 #define PLL_REARRANGE_SPR 0
150 #define PLL_REARRANGE_TBR 1
151 #define PLL_REARRANGE_NNI 2
152 
153 
154 #define PLL_AA_SCALE 10.0
155 #define PLL_AA_SCALE_PLUS_EPSILON 10.001
156 
157 /* ALPHA_MIN is critical -> numerical instability, eg for 4 discrete rate cats */
158 /* and alpha = 0.01 the lowest rate r_0 is */
159 /* 0.00000000000000000000000000000000000000000000000000000000000034878079110511010487 */
160 /* which leads to numerical problems Table for alpha settings below: */
161 /* */
162 /* 0.010000 0.00000000000000000000000000000000000000000000000000000000000034878079110511010487 */
163 /* 0.010000 yielded nasty numerical bugs in at least one case ! */
164 /* 0.020000 0.00000000000000000000000000000044136090435925743185910935350715027016962154188875 */
165 /* 0.030000 0.00000000000000000000476844846859006690412039180149775802624789852441798419292220 */
166 /* 0.040000 0.00000000000000049522423236954066431210260930029681736928018820007024736185030633 */
167 /* 0.050000 0.00000000000050625351310359203371872643495343928538368616365517027588794007897377 */
168 /* 0.060000 0.00000000005134625283884191118711474021861409372524676086868566926568746566772461 */
169 /* 0.070000 0.00000000139080650074206434685544624965062437960128249869740102440118789672851562 */
170 /* 0.080000 0.00000001650681201563587066858709818343436959153791576682124286890029907226562500 */
171 /* 0.090000 0.00000011301977332931251259273962858978301859735893231118097901344299316406250000 */
172 /* 0.100000 0.00000052651925834844387815526344648331402709118265192955732345581054687500000000 */
173 
174 
175 #define PLL_ALPHA_MIN 0.02
176 #define PLL_ALPHA_MAX 1000.0
177 
178 #define PLL_RATE_MIN 0.0000001
179 #define PLL_RATE_MAX 1000000.0
180 
181 
182 
183 /*
184  previous values between 0.001 and 0.000001
185 
186  TO AVOID NUMERICAL PROBLEMS WHEN FREQ == 0 IN PARTITIONED MODELS, ESPECIALLY WITH AA
187  previous value of FREQ_MIN was: 0.000001, but this seemed to cause problems with some
188  of the 7-state secondary structure models with some rather exotic small toy test datasets,
189  on the other hand 0.001 caused problems with some of the 16-state secondary structure models
190 
191  For some reason the frequency settings seem to be repeatedly causing numerical problems
192 
193 */
194 
195 #define PLL_ITMAX 100 /* max number of iterations in brent's algorithm */
196 
197 
198 
199 #define PLL_SHFT(a,b,c,d) (a)=(b);(b)=(c);(c)=(d);
200 #define PLL_SIGN(a,b) ((b) > 0.0 ? fabs(a) : -fabs(a))
201 
202 #define PLL_ABS(x) (((x)<0) ? (-(x)) : (x))
203 #define PLL_MIN(x,y) (((x)<(y)) ? (x) : (y))
204 #define PLL_MAX(x,y) (((x)>(y)) ? (x) : (y))
205 #define PLL_FABS(x) fabs(x)
206 
207 #ifdef _USE_FPGA_LOG
208 extern double log_approx (double input);
209 #define PLL_LOG(x) log_approx(x)
210 #else
211 #define PLL_LOG(x) log(x)
212 #endif
213 
214 
215 #ifdef _USE_FPGA_EXP
216 extern double exp_approx (double x);
217 #define EXP(x) exp_approx(x)
218 #else
219 #define EXP(x) exp(x)
220 #endif
221 
222 
223 
224 #define PLL_SWAP(x,y) do{ __typeof__ (x) _t = x; x = y; y = _t; } while(0)
225 
226 
227 #define PointGamma(prob,alpha,beta) PointChi2(prob,2.0*(alpha))/(2.0*(beta))
228 
229 #define PLL_LIB_NAME "PLL"
230 #define PLL_LIB_VERSION "1.0.0"
231 #define PLL_LIB_DATE "September 2013"
232 
233 
234 #define PLL_DAYHOFF 0
235 #define PLL_DCMUT 1
236 #define PLL_JTT 2
237 #define PLL_MTREV 3
238 #define PLL_WAG 4
239 #define PLL_RTREV 5
240 #define PLL_CPREV 6
241 #define PLL_VT 7
242 #define PLL_BLOSUM62 8
243 #define PLL_MTMAM 9
244 #define PLL_LG 10
245 #define PLL_MTART 11
246 #define PLL_MTZOA 12
247 #define PLL_PMB 13
248 #define PLL_HIVB 14
249 #define PLL_HIVW 15
250 #define PLL_JTTDCMUT 16
251 #define PLL_FLU 17
252 #define PLL_AUTO 18
253 #define PLL_LG4 19
254 #define PLL_GTR 20 /* GTR always needs to be the last one */
255 
256 #define PLL_NUM_PROT_MODELS 21
257 
258 /* bipartition stuff */
259 
260 #define PLL_BIPARTITIONS_RF 4
261 
262 
263 #define PLL_TIP_TIP 0
264 #define PLL_TIP_INNER 1
265 #define PLL_INNER_INNER 2
266 
267 #define PLL_MIN_MODEL -1
268 #define PLL_BINARY_DATA 0
269 #define PLL_DNA_DATA 1
270 #define PLL_AA_DATA 2
271 #define PLL_SECONDARY_DATA 3
272 #define PLL_SECONDARY_DATA_6 4
273 #define PLL_SECONDARY_DATA_7 5
274 #define PLL_GENERIC_32 6
275 #define PLL_GENERIC_64 7
276 #define PLL_MAX_MODEL 8
277 
278 #define PLL_SEC_6_A 0
279 #define PLL_SEC_6_B 1
280 #define PLL_SEC_6_C 2
281 #define PLL_SEC_6_D 3
282 #define PLL_SEC_6_E 4
283 
284 #define PLL_SEC_7_A 5
285 #define PLL_SEC_7_B 6
286 #define PLL_SEC_7_C 7
287 #define PLL_SEC_7_D 8
288 #define PLL_SEC_7_E 9
289 #define PLL_SEC_7_F 10
290 
291 #define PLL_SEC_16 11
292 #define PLL_SEC_16_A 12
293 #define PLL_SEC_16_B 13
294 #define PLL_SEC_16_C 14
295 #define PLL_SEC_16_D 15
296 #define PLL_SEC_16_E 16
297 #define PLL_SEC_16_F 17
298 #define PLL_SEC_16_I 18
299 #define PLL_SEC_16_J 19
300 #define PLL_SEC_16_K 20
301 
302 #define PLL_ORDERED_MULTI_STATE 0
303 #define PLL_MK_MULTI_STATE 1
304 #define PLL_GTR_MULTI_STATE 2
305 
306 #define PLL_CAT 0
307 #define PLL_GAMMA 1
308 
309 /* recomp */
310 #define PLL_SLOT_UNUSED -2 /* value to mark an available vector */
311 #define PLL_NODE_UNPINNED -3 /* marks an inner node as not available in RAM */
312 #define PLL_INNER_NODE_INIT_STLEN -1 /* initialization */
313 
314 #define PLL_MIN_RECOM_FRACTION 0.1 /* at least this % of inner nodes will be allocated in RAM */
315 #define PLL_MAX_RECOM_FRACTION 1.0 /* always 1, just there for boundary checks */
316 
317 
318 typedef int boolean;
319 
320 /* @brief PLL instance attribute structure */
321 typedef struct
322 {
323  int rateHetModel;
324  int fastScaling;
325  int saveMemory;
326  int useRecom;
327  long randomNumberSeed;
328  int numberOfThreads;
330 
332 typedef struct
333 {
335  int *iVector;
336  int *iNode;
337  int *stlen;
338  int *unpinnable;
339  int maxVectorsUsed;
340  boolean allSlotsBusy;
341 } recompVectors;
342 /* E recomp */
343 
346 struct ent
347 {
348  unsigned int *bitVector;
349  unsigned int *treeVector;
350  unsigned int amountTips;
351  int *supportVector;
352  unsigned int bipNumber;
353  unsigned int bipNumber2;
354  unsigned int supportFromTreeset[2];
355  struct ent *next;
356 };
357 
358 typedef struct ent entry;
359 
360 typedef unsigned int hashNumberType;
361 
362 
363 
364 /*typedef uint_fast32_t parsimonyNumber;*/
365 
366 #define PLL_PCF 32
367 
370 typedef struct
371 {
372  hashNumberType tableSize;
373  entry **table;
374  hashNumberType entryCount;
375 }
376  hashtable;
377 struct stringEnt
378 {
379  int nodeNumber;
380  char *word;
381  struct stringEnt *next;
382 };
383 
384 typedef struct stringEnt stringEntry;
385 typedef struct
386 {
387  hashNumberType tableSize;
388  stringEntry **table;
389 }
391 
392 
393 
394 
398 typedef struct ratec
399 {
400  double accumulatedSiteLikelihood;
401  double rate;
403 
418 typedef struct
419 {
420  int tipCase;
421  int pNumber;
422  int qNumber;
423  int rNumber;
424  double qz[PLL_NUM_BRANCHES];
425  double rz[PLL_NUM_BRANCHES];
426  /* recom */
427  int slot_p;
428  int slot_q;
429  int slot_r;
430  /* E recom */
431 } traversalInfo;
432 
437 typedef struct
438 {
440  int count;
441  int functionType;
442  boolean traversalHasChanged;
443  boolean *executeModel;
444  double *parameterValues;
445 } traversalData;
446 
456 struct noderec;
457 
462 typedef struct
463 {
464  unsigned int *vector;
465  int support;
466  struct noderec *oP;
467  struct noderec *oQ;
468 } branchInfo;
469 
470 
471 
472 
473 
478 typedef struct
479 {
480  boolean valid;
481  int partitions;
482  int *partitionList;
483 }
484  linkageData;
485 typedef struct
486 {
487  int entries;
488  linkageData* ld;
489 }
490  linkageList;
491 
492 
493 
702 typedef struct noderec
703 {
704 
705  branchInfo *bInf;
706  double z[PLL_NUM_BRANCHES];
707  struct noderec *next;
708  struct noderec *back;
709  hashNumberType hash;
710  int support;
711  int number;
712  char x;
713  char xPars;
714  char xBips;
715 }
716  node, *nodeptr;
717 
731 typedef struct
732  {
733  double lh;
734  int number;
735  }
736  info;
737 
738 typedef struct bInf {
739  double likelihood;
740  nodeptr node;
741 } bestInfo;
742 
743 typedef struct iL {
744  bestInfo *list;
745  int n;
746  int valid;
747 } infoList;
748 
749 
750 typedef unsigned int parsimonyNumber;
751 
752 /* @brief Alignment, transition model, model of rate heterogenety and likelihood vectors for one partition.
753  *
754  * @todo De-couple into smaller data structures
755  *
756  * ALIGNMENT DATA
757  * This depends only on the type of data in this partition of the alignment
758  *
759  * MODEL OF RATE HETEROGENETY, We use either GAMMA or PSR
760  * Rate heterogenety: Per Site Categories (PSR) model aka CAT,
761  * Rate of site i is given by perSiteRates[rateCategory[i]]
762  *
763  * TRANSITION MODEL: We always assume General Time Reversibility
764  * Transistion probability matrix: P(t) = exp(Qt)
765  * Branch length t is the expected number of substitutions per site
766  * Pij(t) is the probability of going from state i to state j in a branch of length t
767  * Relative substitution rates (Entries in the Q matrix)
768  * In GTR we can write Q = S * D, where S is a symmetrical matrix and D a diagonal with the state frequencies
769 
770  @var protModels
771  @brief Protein models
772 
773  @detail Detailed protein models descriptiopn
774 
775  @var autoProtModels
776  @brief Auto prot models
777  @detail Detailed auto prot models
778  */
779 
780 
781 
924 typedef struct {
925  int dataType;
926  int states;
929  int lower;
930  int upper;
931  int width;
932  int *wgt;
934 
935 
936  /* MODEL OF RATE HETEROGENETY, We use either GAMMA or PSR */
937  /* Rate heterogenety: Per Site Categories (PSR) model aka CAT, see updatePerSiteRates() */
938  /* Rate of site i is given by perSiteRates[rateCategory[i]] */
939  double *perSiteRates;
942  /* Rate heterogenety: GAMMA model of rate heterogenety */
943  double alpha;
944  double *gammaRates;
945 
946 
947  /* TRANSITION MODEL: We always assume General Time Reversibility */
948  /* Transistion probability matrix: P(t) = exp(Qt)*/
949  /* Branch length t is the expected number of substitutions per site */
950  /* Pij(t) is the probability of going from state i to state j in a branch of length t */
951  /* Relative substitution rates (Entries in the Q matrix) */
952  /* In GTR we can write Q = S * D, where S is a symmetrical matrix and D a diagonal with the state frequencies */
953  double *substRates;
954  double *frequencies;
955  double *freqExponents;
956  /* Matrix decomposition: @todo map this syntax to Explanation of the mathematical background */
957  double *EIGN;
958  double *EV;
959  double *EI;
960  double *left;
961  double *right;
962  double *tipVector;
963 
964  /* LG4 */
965 
966  double *EIGN_LG4[4];
967  double *EV_LG4[4];
968  double *EI_LG4[4];
969 
970  double *frequencies_LG4[4];
971  double *tipVector_LG4[4];
972  double *substRates_LG4[4];
973 
974  /* LG4 */
975 
976  /* Protein specific ?? */
979  int protFreqs;
980  /* specific for secondary structures ?? */
981  boolean nonGTR;
982  boolean optimizeBaseFrequencies;
983  boolean optimizeAlphaParameter;
984  boolean optimizeSubstitutionRates;
985  int *symmetryVector;
987 
988 
989  /* LIKELIHOOD VECTORS */
990 
991  /* partial LH Inner vectors ancestral vectors, we have 2*tips - 3 inner nodes */
992  double **xVector;
993  unsigned char **yVector;
994  unsigned int *globalScaler;
996  /* data structures for conducting per-site likelihood scaling.
997  this allows to compute the per-site log likelihood scores
998  needed for RELL-based bootstrapping and all sorts of statistical
999  tests for comparing trees ! */
1000  int **expVector;
1001  size_t *expSpaceVector;
1003  /* These are for the saveMemory option (tracking gaps to skip computations and memory) */
1004  size_t *xSpaceVector; /* Size of conditional likelihood vectors per inner node */
1005  int gapVectorLength;
1006  unsigned int *gapVector;
1007  double *gapColumn;
1008 
1009  /* Parsimony vectors at each node */
1010  size_t parsimonyLength;
1011  parsimonyNumber *parsVect;
1012 
1013  /* This buffer of size width is used to store intermediate values for the branch length optimization under
1014  newton-raphson. The data in here can be re-used for all iterations irrespective of the branch length.
1015  */
1016  double *sumBuffer;
1017 
1018  /* Buffer to store the per-site log likelihoods */
1019  double *perSiteLikelihoods;
1020 
1021  /* This buffer of size width is used to store the ancestral state at a node of the tree. */
1022  double *ancestralBuffer;
1023 
1024  /* From tree */
1025  boolean executeModel;
1026  double fracchange;
1027  double partitionContribution;
1028  double partitionLH;
1029 
1030 // #if (defined(_USE_PTHREADS) || defined(_FINE_GRAIN_MPI))
1031  int partitionAssignment;
1032 // #endif
1033 
1034 } pInfo;
1035 
1036 typedef struct
1037  {
1038  pInfo **partitionData;
1039  int numberOfPartitions;
1040  boolean perGeneBranchLengths;
1041  boolean dirty;
1042  linkageList *alphaList;
1043  linkageList *rateList;
1044  linkageList *freqList;
1045  } partitionList;
1046 
1047 
1048 
1049 #define PLL_REARR_SETTING 1
1050 #define PLL_FAST_SPRS 2
1051 #define PLL_SLOW_SPRS 3
1052 
1053 
1058 typedef struct {
1059 
1060  int state;
1061 
1062  /*unsigned int vLength;*/
1063  double accumulatedTime;
1064  int rearrangementsMax;
1065  int rearrangementsMin;
1066  int thoroughIterations;
1067  int fastIterations;
1068  int mintrav;
1069  int maxtrav;
1070  int bestTrav;
1071  double startLH;
1072  double lh;
1073  double previousLh;
1074  double difference;
1075  double epsilon;
1076  boolean impr;
1077  boolean cutoff;
1078 
1079  double tr_startLH;
1080  double tr_endLH;
1081  double tr_likelihood;
1082  double tr_bestOfNode;
1083  double tr_lhCutoff;
1084  double tr_lhAVG;
1085  double tr_lhDEC;
1086  int tr_NumberOfCategories;
1087  int tr_itCount;
1088  int tr_doCutoff;
1089  int tr_thoroughInsertion;
1090  int tr_optimizeRateCategoryInvocations;
1091 
1092 
1093  /* prevent users from doing stupid things */
1094 
1095 
1096  int searchConvergenceCriterion;
1097  int rateHetModel;
1098  int maxCategories;
1099  int NumberOfModels;
1100  int numBranches;
1101  int originalCrunchedLength;
1102  int mxtips;
1103  char seq_file[1024];
1104 } checkPointState;
1105 
1106 
1107 
1108 /* recomp */
1109 #ifdef _DEBUG_RECOMPUTATION
1110 typedef struct {
1111  unsigned long int numTraversals;
1112  unsigned long int tt;
1113  unsigned long int ti;
1114  unsigned long int ii;
1115  unsigned int *travlenFreq;
1117 #endif
1118 /* E recomp */
1119 
1120 
1125 typedef struct {
1126 
1127  int *ti;
1128 
1129  /* recomp */
1133  boolean useRecom;
1134 #ifdef _DEBUG_RECOMPUTATION
1135  traversalCounter *travCounter;
1136  double stlenTime;
1137 #endif
1138  /* E recomp */
1139 
1140 
1141  boolean fastScaling;
1142  boolean saveMemory;
1143  int startingTree;
1144  long randomNumberSeed;
1145 
1146  double *lhs;
1147  double *patrat;
1148  double *patratStored;
1149  int *rateCategory;
1150  int *aliaswgt;
1151  boolean manyPartitions;
1152 
1153  boolean grouped;
1154  boolean constrained;
1155  int threadID;
1156  volatile int numberOfThreads;
1157 
1158 //#if (defined(_USE_PTHREADS) || defined(_FINE_GRAIN_MPI))
1159 
1160  unsigned char *y_ptr;
1161 
1162  double lower_spacing;
1163  double upper_spacing;
1164 
1165  double *ancestralVector;
1166 //#endif
1167 
1168 
1169 
1170 
1171 
1172 
1173 
1174  stringHashtable *nameHash;
1175 
1176  char *secondaryStructureInput;
1177 
1178  traversalData td[1];
1179 
1180  int maxCategories;
1181  int categories;
1182 
1183  double coreLZ[PLL_NUM_BRANCHES];
1184 
1185 
1186  branchInfo *bInf;
1187 
1188  int multiStateModel;
1189 
1190 
1191  boolean curvatOK[PLL_NUM_BRANCHES];
1192  /* the stuff below is shared among DNA and AA, span does
1193  not change depending on datatype */
1194 
1195 
1196 
1197  /* model stuff end */
1198  int bDeep;
1200  unsigned char **yVector;
1202  int secondaryStructureModel;
1206  int *secondaryStructurePairs;
1207 
1208 
1209  double fracchange;
1210  double lhCutoff;
1211  double lhAVG;
1212  unsigned long lhDEC;
1213  unsigned long itCount;
1214  int numberOfInvariableColumns;
1215  int weightOfInvariableColumns;
1216  int rateHetModel;
1217 
1218  double startLH;
1219  double endLH;
1220  double likelihood;
1223  nodeptr nodeBaseAddress;
1225  int mxtips;
1228  int numberOfSecondaryColumns;
1229  boolean searchConvergenceCriterion;
1230  int ntips;
1231  int nextnode;
1232 
1233 
1234  boolean bigCutoff;
1235  boolean partitionSmoothed[PLL_NUM_BRANCHES];
1236  boolean partitionConverged[PLL_NUM_BRANCHES];
1237  boolean rooted;
1238  boolean doCutoff;
1239 
1240  double gapyness;
1241 
1242  char **nameList;
1243  char *tree_string;
1244  char *tree0;
1245  char *tree1;
1246  int treeStringLength;
1247 
1248  unsigned int bestParsimony;
1249  unsigned int *parsimonyScore;
1250 
1251  double bestOfNode;
1252  nodeptr removeNode;
1253  nodeptr insertNode;
1255  double zqr[PLL_NUM_BRANCHES];
1256  double currentZQR[PLL_NUM_BRANCHES];
1257 
1258  double currentLZR[PLL_NUM_BRANCHES];
1259  double currentLZQ[PLL_NUM_BRANCHES];
1260  double currentLZS[PLL_NUM_BRANCHES];
1261  double currentLZI[PLL_NUM_BRANCHES];
1262  double lzs[PLL_NUM_BRANCHES];
1263  double lzq[PLL_NUM_BRANCHES];
1264  double lzr[PLL_NUM_BRANCHES];
1265  double lzi[PLL_NUM_BRANCHES];
1266 
1267 
1268  unsigned int **bitVectors;
1269 
1270  unsigned int vLength;
1271 
1274  int optimizeRateCategoryInvocations;
1275 
1276  checkPointState ckp;
1278  boolean useMedian;
1279 
1280  pllStack * rearrangeHistory;
1281 
1282 
1283  /* analdef defines */
1284  /* TODO: Do some initialization */
1285  int bestTrav;
1288  int initial;
1289  boolean initialSet;
1290  int mode;
1291  boolean perGeneBranchLengths;
1293  boolean compressPatterns;
1294  double likelihoodEpsilon;
1295  boolean useCheckpoint;
1296 
1297 } pllInstance;
1298 
1300 typedef struct {
1301  pllInstance * tr;
1302  nodeptr p;
1303  int nniType;
1304  double z[PLL_NUM_BRANCHES]; // optimize branch lengths
1305  double z0[PLL_NUM_BRANCHES]; // unoptimized branch lengths
1306  double likelihood;
1307  double deltaLH;
1308 } nniMove;
1309 
1310 /***************************************************************/
1311 
1312 typedef struct {
1313  int partitionNumber;
1314  int partitionLength;
1315 } partitionType;
1316 
1317 typedef struct
1318 {
1319  double z[PLL_NUM_BRANCHES];
1320  nodeptr p, q;
1321  int cp, cq;
1322 }
1324 
1325 typedef struct
1326 {
1327  connectRELL *connect;
1328  int start;
1329  double likelihood;
1330 }
1331  topolRELL;
1332 
1333 
1334 typedef struct
1335 {
1336  int max;
1337  topolRELL **t;
1338 }
1340 
1341 
1342 
1343 
1344 
1345 /**************************************************************/
1346 
1347 
1348 
1351 typedef struct conntyp {
1352  double z[PLL_NUM_BRANCHES];
1353  node *p, *q;
1354  void *valptr;
1355  int descend;
1356  int sibling;
1357  } connect, *connptr;
1358 
1361 typedef struct {
1362  double likelihood;
1363  int initialTreeNumber;
1365  node *start;
1366  int nextlink;
1368  int ntips;
1369  int nextnode;
1370  int scrNum;
1371  int tplNum;
1372  } topol;
1373 
1379 typedef struct {
1380  double *probs;
1381  char c;
1382  int states;
1383 } ancestralState;
1384 
1388 typedef struct {
1389  double best;
1390  double worst;
1392  topol **byScore;
1393  topol **byTopol;
1394  int nkeep;
1395  int nvalid;
1396  int ninit;
1397  int numtrees;
1398  boolean improved;
1399  } bestlist;
1400 
1403 typedef struct
1404 {
1408  int evLength;
1409  int eiLength;
1412  int tipVectorLength; /* ??? */
1413  int symmetryVectorLength;
1414  int frequencyGroupingLength;
1415 
1416  boolean nonGTR;
1417  boolean optimizeBaseFrequencies;
1418 
1419  int undetermined;
1420 
1421  const char *inverseMeaning;
1422 
1423  int states; /* s */
1424 
1425  boolean smoothFrequencies;
1426 
1427  const unsigned int *bitVector;
1428 
1430 
1431 typedef struct
1432 {
1433  int rearrangeType;
1434  double likelihood;
1435 
1436  union {
1437  struct {
1438  double * zp;
1439  double * zpn;
1440  double * zpnn;
1441  double * zqr;
1442  nodeptr pn;
1443  nodeptr pnn;
1444  nodeptr r;
1445  nodeptr p;
1446  nodeptr q;
1447  } SPR;
1448  struct {
1449  nodeptr origin;
1450  int swapType;
1451  double z[PLL_NUM_BRANCHES];
1452  } NNI;
1453  };
1454 } pllRollbackInfo;
1455 
1456 
1475 typedef struct
1476  {
1477  nodeptr p;
1478  int mintrav;
1479  int maxtrav;
1480  } pllRearrangeAttr;
1481 
1506 typedef struct
1507  {
1509  double likelihood;
1510  union {
1511  struct {
1512  nodeptr removeNode;
1513  nodeptr insertNode;
1514  double zqr[PLL_NUM_BRANCHES];
1515  } SPR;
1516  struct {
1517  nodeptr originNode;
1518  int swapType;
1519  } NNI;
1520  };
1521  } pllRearrangeInfo;
1522 
1523 
1524 typedef struct
1525  {
1526  int max_entries;
1527  int entries;
1528  pllRearrangeInfo * rearr;
1529  } pllRearrangeList;
1530 
1532 typedef struct
1533  {
1536  char ** sequenceLabels;
1537  unsigned char ** sequenceData;
1538  int * siteWeights;
1539  } pllAlignmentData;
1540 
1541 
1542 /****************************** FUNCTIONS ****************************************************/
1543 
1544 
1545 #if (defined(_USE_PTHREADS) || defined(_FINE_GRAIN_MPI))
1546 boolean isThisMyPartition(partitionList *pr, int tid, int model);
1547 void printParallelTimePerRegion(void);
1548 #endif
1549 
1550 #ifdef _FINE_GRAIN_MPI
1551 extern void pllFinalizeMPI (void);
1552 #endif
1553 
1554 
1555 extern void pllStartPthreads (pllInstance *tr, partitionList *pr);
1556 extern void pllStopPthreads (pllInstance * tr);
1557 extern int lookupWord(char *s, stringHashtable *h);
1558 extern void pllLockMPI (pllInstance * tr);
1559 extern void pllInitMPI(int * argc, char **argv[]);
1560 
1561 extern void getDataTypeString(pllInstance *tr, pInfo *partitionInfo, char typeOfData[1024]);
1562 
1563 extern int countTips(nodeptr p, int numsp);
1564 extern entry *initEntry(void);
1565 extern unsigned int precomputed16_bitcount(unsigned int n, char *bits_in_16bits);
1566 
1567 
1568 /* Handling branch lengths*/
1569 extern double pllGetBranchLength (pllInstance *tr, nodeptr p, int partition_id);
1570 extern void pllSetBranchLength (pllInstance *tr, nodeptr p, int partition_id, double bl);
1571 
1572 
1573 extern size_t discreteRateCategories(int rateHetModel);
1574 
1575 extern const partitionLengths * getPartitionLengths(pInfo *p);
1576 extern boolean getSmoothFreqs(int dataType);
1577 extern const unsigned int *getBitVector(int dataType);
1578 extern int getUndetermined(int dataType);
1579 extern int getStates(int dataType);
1580 extern char getInverseMeaning(int dataType, unsigned char state);
1581 extern double gettime ( void );
1582 extern int gettimeSrand ( void );
1583 extern double randum ( long *seed );
1584 
1585 extern void getxnode ( nodeptr p );
1586 extern void hookup ( nodeptr p, nodeptr q, double *z, int numBranches);
1587 extern void hookupFull ( nodeptr p, nodeptr q, double *z);
1588 extern void hookupDefault ( nodeptr p, nodeptr q);
1589 extern boolean whitechar ( int ch );
1590 extern void printLog ( pllInstance *tr);
1591 extern void initReversibleGTR( pllInstance *tr, partitionList *pr, int model);
1592 extern double LnGamma ( double alpha );
1593 extern double IncompleteGamma ( double x, double alpha, double ln_gamma_alpha );
1594 extern double PointNormal ( double prob );
1595 extern double PointChi2 ( double prob, double v );
1596 extern void makeGammaCats (double alpha, double *gammaRates, int K, boolean useMedian);
1597 extern void initModel ( pllInstance *tr, double **empiricalFrequencies, partitionList * partitions);
1598 
1599 extern void resetBranches ( pllInstance *tr );
1600 extern void modOpt ( pllInstance *tr, partitionList *pr, double likelihoodEpsilon);
1601 
1602 extern void initializePartitionData(pllInstance *localTree, partitionList * localPartitions);
1603 extern void initMemorySavingAndRecom(pllInstance *tr, partitionList *pr);
1604 
1605 
1606 
1607 extern void makeRandomTree ( pllInstance *tr);
1608 extern void nodeRectifier ( pllInstance *tr );
1609 extern void makeParsimonyTreeFast(pllInstance *tr, partitionList *pr);
1610 extern void allocateParsimonyDataStructures(pllInstance *tr, partitionList *pr);
1611 extern void freeParsimonyDataStructures(pllInstance *tr, partitionList *pr);
1612 extern void parsimonySPR(nodeptr p, partitionList *pr, pllInstance *tr);
1613 
1614 extern FILE *myfopen(const char *path, const char *mode);
1615 
1616 
1617 extern boolean initrav ( pllInstance *tr, partitionList *pr, nodeptr p );
1618 extern void initravPartition ( pllInstance *tr, nodeptr p, int model );
1619 extern void update ( pllInstance *tr, partitionList *pr, nodeptr p );
1620 extern void smooth ( pllInstance *tr, partitionList *pr, nodeptr p );
1621 extern void smoothTree ( pllInstance *tr, partitionList *pr, int maxtimes );
1622 extern void localSmooth ( pllInstance *tr, partitionList *pr, nodeptr p, int maxtimes );
1623 extern boolean localSmoothMulti(pllInstance *tr, nodeptr p, int maxtimes, int model);
1624 extern int pllNniSearch(pllInstance * tr, partitionList *pr, int estimateModel);
1625 extern int NNI(pllInstance * tr, nodeptr p, int swap);
1626 
1627 extern void smoothRegion ( pllInstance *tr, partitionList *pr, nodeptr p, int region );
1628 extern void regionalSmooth ( pllInstance *tr, partitionList *pr, nodeptr p, int maxtimes, int region );
1629 extern nodeptr removeNodeBIG ( pllInstance *tr, partitionList *pr, nodeptr p, int numBranches);
1630 extern nodeptr removeNodeRestoreBIG ( pllInstance *tr, partitionList *pr, nodeptr p );
1631 extern boolean insertBIG ( pllInstance *tr, partitionList *pr, nodeptr p, nodeptr q);
1632 extern boolean insertRestoreBIG ( pllInstance *tr, partitionList *pr, nodeptr p, nodeptr q );
1633 extern boolean testInsertBIG ( pllInstance *tr, partitionList *pr, nodeptr p, nodeptr q );
1634 extern void addTraverseBIG ( pllInstance *tr, partitionList *pr, nodeptr p, nodeptr q, int mintrav, int maxtrav );
1635 extern int rearrangeBIG ( pllInstance *tr, partitionList *pr, nodeptr p, int mintrav, int maxtrav );
1636 extern void traversalOrder ( nodeptr p, int *count, nodeptr *nodeArray );
1637 extern boolean testInsertRestoreBIG ( pllInstance *tr, partitionList *pr, nodeptr p, nodeptr q );
1638 extern void restoreTreeFast ( pllInstance *tr, partitionList *pr );
1639 
1640 extern void pllTreeEvaluate ( pllInstance *tr, partitionList *pr, int maxSmoothIterations );
1641 
1642 extern void initTL ( topolRELL_LIST *rl, pllInstance *tr, int n );
1643 extern void freeTL ( topolRELL_LIST *rl);
1644 extern void restoreTL ( topolRELL_LIST *rl, pllInstance *tr, int n, int numBranches );
1645 extern void resetTL ( topolRELL_LIST *rl );
1646 extern void saveTL ( topolRELL_LIST *rl, pllInstance *tr, int index );
1647 
1648 extern int saveBestTree (bestlist *bt, pllInstance *tr, int numBranches);
1649 extern int recallBestTree (bestlist *bt, int rank, pllInstance *tr, partitionList *pr);
1650 extern int initBestTree ( bestlist *bt, int newkeep, int numsp );
1651 extern void resetBestTree ( bestlist *bt );
1652 extern boolean freeBestTree ( bestlist *bt );
1653 
1654 
1655 extern char *Tree2String ( char *treestr, pllInstance *tr, partitionList *pr, nodeptr p, boolean printBranchLengths, boolean printNames, boolean printLikelihood,
1656  boolean rellTree, boolean finalPrint, int perGene, boolean branchLabelSupport, boolean printSHSupport);
1657 void printTopology(pllInstance *tr, partitionList *pr, boolean printInner);
1658 
1659 
1660 
1661 extern int treeReadLen (FILE *fp, pllInstance *tr, boolean readBranches, boolean readNodeLabels, boolean topologyOnly);
1662 extern void treeReadTopologyString(char *treeString, pllInstance *tr);
1663 
1664 extern void getStartingTree (pllInstance *tr);
1665 extern double treeLength (pllInstance *tr, int model);
1666 
1667 extern double evaluatePartialGeneric (pllInstance *, partitionList *pr, int i, double ki, int _model);
1668 extern void pllEvaluateGeneric (pllInstance *tr, partitionList *pr, nodeptr p, boolean fullTraversal, boolean getPerSiteLikelihoods);
1669 extern void pllNewviewGeneric (pllInstance *tr, partitionList *pr, nodeptr p, boolean masked);
1670 
1671 extern void pllNewviewGenericAncestral(pllInstance *tr, partitionList *pr, nodeptr p);
1673 extern void printAncestralState(nodeptr p, boolean printStates, boolean printProbs, pllInstance *tr, partitionList *pr);
1674 
1675 extern void makenewzGeneric(pllInstance *tr, partitionList * pr, nodeptr p, nodeptr q, double *z0, int maxiter, double *result, boolean mask);
1676 extern void makenewzGenericDistance(pllInstance *tr, int maxiter, double *z0, double *result, int taxon1, int taxon2);
1677 extern double evaluatePartitionGeneric (pllInstance *tr, nodeptr p, int model);
1678 extern void newviewPartitionGeneric (pllInstance *tr, nodeptr p, int model);
1679 extern double evaluateGenericVector (pllInstance *tr, nodeptr p);
1680 extern void categorizeGeneric (pllInstance *tr, nodeptr p);
1681 extern double makenewzPartitionGeneric(pllInstance *tr, nodeptr p, nodeptr q, double z0, int maxiter, int model);
1682 extern boolean isTip(int number, int maxTips);
1683 
1684 /* recom functions */
1685 extern void computeTraversal(pllInstance *tr, nodeptr p, boolean partialTraversal, int numBranches);
1686 extern void allocRecompVectorsInfo(pllInstance *tr);
1687 extern void allocTraversalCounter(pllInstance *tr);
1688 extern boolean getxVector(recompVectors *rvec, int nodenum, int *slot, int mxtips);
1689 extern boolean needsRecomp(boolean recompute, recompVectors *rvec, nodeptr p, int mxtips);
1690 extern void unpinNode(recompVectors *v, int nodenum, int mxtips);
1691 extern void protectNode(recompVectors *rvec, int nodenum, int mxtips);
1692 
1693 extern void computeTraversalInfoStlen(nodeptr p, int maxTips, recompVectors *rvec, int *count);
1694 extern void computeFullTraversalInfoStlen(nodeptr p, int maxTips, recompVectors *rvec);
1695 extern void printTraversalInfo(pllInstance *tr);
1696 extern void countTraversal(pllInstance *tr);
1697 
1698 
1699 extern void pllNewviewIterative(pllInstance *tr, partitionList *pr, int startIndex);
1700 extern void pllEvaluateIterative(pllInstance *tr, partitionList *pr, boolean getPerSiteLikelihoods);
1701 extern void storeExecuteMaskInTraversalDescriptor(pllInstance *tr, partitionList *pr);
1702 extern void storeValuesInTraversalDescriptor(pllInstance *tr, partitionList *pr, double *value);
1703 extern void makenewzIterative(pllInstance *, partitionList *pr);
1704 extern void execCore(pllInstance *, partitionList *pr, volatile double *dlnLdlz, volatile double *d2lnLdlz2);
1705 
1706 
1707 
1708 
1709 
1710 extern nodeptr findAnyTip(nodeptr p, int numsp);
1711 extern void putWAG(double *ext_initialRates);
1712 extern unsigned int **initBitVector(int mxtips, unsigned int *vectorLength);
1713 extern hashtable *initHashTable(unsigned int n);
1714 extern void cleanupHashTable(hashtable *h, int state);
1715 extern double convergenceCriterion(hashtable *h, int mxtips);
1716 extern void freeBitVectors(unsigned int **v, int n);
1717 extern void freeHashTable(hashtable *h);
1718 extern stringHashtable *initStringHashTable(hashNumberType n);
1719 extern void addword(char *s, stringHashtable *h, int nodeNumber);
1720 
1721 extern void printBothOpen(const char* format, ... );
1722 extern void initRateMatrix(pllInstance *tr, partitionList *pr);
1723 
1724 extern void bitVectorInitravSpecial(unsigned int **bitVectors, nodeptr p, int numsp, unsigned int vectorLength, hashtable *h, int treeNumber, int function, branchInfo *bInf,
1725  int *countBranches, int treeVectorLength, boolean traverseOnly, boolean computeWRF, int processID);
1726 
1727 
1728 extern unsigned int bitcount_32_bit(unsigned int i);
1729 extern inline unsigned int bitcount_64_bit(unsigned long i);
1730 
1731 
1732 extern void perSiteLogLikelihoods(pllInstance *tr, partitionList *pr, double *logLikelihoods);
1733 
1734 extern void updatePerSiteRates(pllInstance *tr, partitionList *pr, boolean scaleRates);
1735 
1736 extern void restart(pllInstance *tr, partitionList *pr);
1737 
1738 extern double getBranchLength(pllInstance *tr, partitionList *pr, int perGene, nodeptr p);
1739 
1740 inline boolean isGap(unsigned int *x, int pos);
1741 inline boolean noGap(unsigned int *x, int pos);
1742 
1743 /* newick parser declarations */
1744 extern pllNewickTree * pllNewickParseString (const char * newick);
1745 extern pllNewickTree * pllNewickParseFile (const char * filename);
1746 extern int pllValidateNewick (pllNewickTree *);
1747 extern void pllNewickParseDestroy (pllNewickTree **);
1748 extern int pllNewickUnroot (pllNewickTree * t);
1749 
1750 /* partition parser declarations */
1751 extern void pllQueuePartitionsDestroy (pllQueue ** partitions);
1752 extern pllQueue * pllPartitionParse (const char * filename);
1753 extern void pllPartitionDump (pllQueue * partitions);
1754 
1755 /* alignment data declarations */
1757 extern int pllAlignmentDataDumpFile (pllAlignmentData *, int, const char *);
1758 extern void pllAlignmentDataDumpConsole (pllAlignmentData * alignmentData);
1759 extern pllAlignmentData * pllInitAlignmentData (int, int);
1760 extern pllAlignmentData * pllParseAlignmentFile (int fileType, const char *);
1761 
1762 extern char * pllReadFile (const char *, long *);
1763 extern int * pllssort1main (char ** x, int n);
1764 /* from utils.h */
1765 linkageList* initLinkageList(int *linkList, partitionList *pr);
1766 
1767 int pllLinkAlphaParameters(char *string, partitionList *pr);
1768 int pllLinkFrequencies(char *string, partitionList *pr);
1769 int pllLinkRates(char *string, partitionList *pr);
1770 int pllSetSubstitutionRateMatrixSymmetries(char *string, partitionList * pr, int model);
1771 
1772 void pllSetFixedAlpha(double alpha, int model, partitionList * pr, pllInstance *tr);
1773 void pllSetFixedBaseFrequencies(double *f, int length, int model, partitionList * pr, pllInstance *tr);
1774 int pllSetOptimizeBaseFrequencies(int model, partitionList * pr, pllInstance *tr);
1775 void pllSetFixedSubstitutionMatrix(double *q, int length, int model, partitionList * pr, pllInstance *tr);
1776 
1777 nodeptr pllGetRandomSubtree(pllInstance *);
1778 void makeParsimonyTree(pllInstance *tr);
1780 int pllPartitionsValidate (pllQueue * parts, pllAlignmentData * alignmentData);
1781 partitionList * pllPartitionsCommit (pllQueue * parts, pllAlignmentData * alignmentData);
1782 extern void pllAlignmentRemoveDups (pllAlignmentData * alignmentData, partitionList * pl);
1783 double ** pllBaseFrequenciesGTR (partitionList * pl, pllAlignmentData * alignmentData);
1785 int pllLoadAlignment (pllInstance * tr, pllAlignmentData * alignmentData, partitionList *, int);
1786 void pllEmpiricalFrequenciesDestroy (double *** empiricalFrequencies, int models);
1787 void pllTreeInitTopologyRandom (pllInstance * tr, int tips, char ** nameList);
1788 void pllTreeInitTopologyForAlignment (pllInstance * tr, pllAlignmentData * alignmentData);
1789 void pllBaseSubstitute (pllAlignmentData * alignmentData, partitionList * partitions);
1794 int pllOptimizeModelParameters(pllInstance *tr, partitionList *pr, double likelihoodEpsilon);
1795 
1797 void pllDestroyRearrangeList (pllRearrangeList ** bestList);
1798 void pllRearrangeSearch (pllInstance * tr, partitionList * pr, int rearrangeType, nodeptr p, int mintrav, int maxtrav, pllRearrangeList * bestList);
1799 void pllRearrangeCommit (pllInstance * tr, partitionList * pr, pllRearrangeInfo * rearr, int saveRollbackInfo);
1802 
1803 #if (defined(_FINE_GRAIN_MPI) || defined(_USE_PTHREADS) )
1804 /* work tags for parallel regions */
1805 
1806 #define PLL_THREAD_NEWVIEW 0
1807 #define PLL_THREAD_EVALUATE 1
1808 #define PLL_THREAD_MAKENEWZ 2
1809 #define PLL_THREAD_MAKENEWZ_FIRST 3
1810 #define PLL_THREAD_RATE_CATS 4
1811 #define PLL_THREAD_COPY_RATE_CATS 5
1812 #define PLL_THREAD_COPY_INIT_MODEL 6
1813 #define PLL_THREAD_INIT_PARTITION 7
1814 #define PLL_THREAD_OPT_ALPHA 8
1815 #define PLL_THREAD_OPT_RATE 9
1816 #define PLL_THREAD_COPY_ALPHA 10
1817 #define PLL_THREAD_COPY_RATES 11
1818 #define PLL_THREAD_PER_SITE_LIKELIHOODS 12
1819 #define PLL_THREAD_NEWVIEW_ANCESTRAL 13
1820 #define PLL_THREAD_GATHER_ANCESTRAL 14
1821 #define PLL_THREAD_EXIT_GRACEFULLY 15
1822 #define PLL_THREAD_EVALUATE_PER_SITE_LIKES 16
1823 
1824 
1825 typedef struct
1826 {
1827  pllInstance *tr;
1828 
1829  partitionList *pr;
1830  int threadNumber;
1831 }
1832  threadData;
1833 extern void optRateCatPthreads(pllInstance *tr, partitionList *pr, double lower_spacing, double upper_spacing, double *lhs, int n, int tid);
1834 extern void pllMasterBarrier(pllInstance *, partitionList *, int);
1835 #endif
1836 
1837 
1838 #ifdef __AVX
1839 
1840 extern void newviewGTRGAMMAPROT_AVX_LG4(int tipCase,
1841  double *x1, double *x2, double *x3, double *extEV[4], double *tipVector[4],
1842  int *ex3, unsigned char *tipX1, unsigned char *tipX2, int n,
1843  double *left, double *right, int *wgt, int *scalerIncrement, const boolean useFastScaling);
1844 
1845 
1846 extern void newviewGTRCAT_AVX_GAPPED_SAVE(int tipCase, double *EV, int *cptr,
1847  double *x1_start, double *x2_start, double *x3_start, double *tipVector,
1848  int *ex3, unsigned char *tipX1, unsigned char *tipX2,
1849  int n, double *left, double *right, int *wgt, int *scalerIncrement, const boolean useFastScaling,
1850  unsigned int *x1_gap, unsigned int *x2_gap, unsigned int *x3_gap,
1851  double *x1_gapColumn, double *x2_gapColumn, double *x3_gapColumn, const int maxCats);
1852 
1853 extern void newviewGTRCATPROT_AVX_GAPPED_SAVE(int tipCase, double *extEV,
1854  int *cptr,
1855  double *x1, double *x2, double *x3, double *tipVector,
1856  int *ex3, unsigned char *tipX1, unsigned char *tipX2,
1857  int n, double *left, double *right, int *wgt, int *scalerIncrement, const boolean useFastScaling,
1858  unsigned int *x1_gap, unsigned int *x2_gap, unsigned int *x3_gap,
1859  double *x1_gapColumn, double *x2_gapColumn, double *x3_gapColumn, const int maxCats);
1860 
1861 extern void newviewGTRGAMMA_AVX_GAPPED_SAVE(int tipCase,
1862  double *x1_start, double *x2_start, double *x3_start,
1863  double *extEV, double *tipVector,
1864  int *ex3, unsigned char *tipX1, unsigned char *tipX2,
1865  const int n, double *left, double *right, int *wgt, int *scalerIncrement, const boolean useFastScaling,
1866  unsigned int *x1_gap, unsigned int *x2_gap, unsigned int *x3_gap,
1867  double *x1_gapColumn, double *x2_gapColumn, double *x3_gapColumn
1868  );
1869 
1870 extern void newviewGTRGAMMAPROT_AVX_GAPPED_SAVE(int tipCase,
1871  double *x1_start, double *x2_start, double *x3_start, double *extEV, double *tipVector,
1872  int *ex3, unsigned char *tipX1, unsigned char *tipX2, int n,
1873  double *left, double *right, int *wgt, int *scalerIncrement, const boolean useFastScaling,
1874  unsigned int *x1_gap, unsigned int *x2_gap, unsigned int *x3_gap,
1875  double *x1_gapColumn, double *x2_gapColumn, double *x3_gapColumn);
1876 
1877 extern void newviewGTRCAT_AVX(int tipCase, double *EV, int *cptr,
1878  double *x1_start, double *x2_start, double *x3_start, double *tipVector,
1879  int *ex3, unsigned char *tipX1, unsigned char *tipX2,
1880  int n, double *left, double *right, int *wgt, int *scalerIncrement, const boolean useFastScaling);
1881 
1882 
1883 extern void newviewGenericCATPROT_AVX(int tipCase, double *extEV,
1884  int *cptr,
1885  double *x1, double *x2, double *x3, double *tipVector,
1886  int *ex3, unsigned char *tipX1, unsigned char *tipX2,
1887  int n, double *left, double *right, int *wgt, int *scalerIncrement, const boolean useFastScaling);
1888 
1889 
1890 extern void newviewGTRGAMMA_AVX(int tipCase,
1891  double *x1_start, double *x2_start, double *x3_start,
1892  double *EV, double *tipVector,
1893  int *ex3, unsigned char *tipX1, unsigned char *tipX2,
1894  const int n, double *left, double *right, int *wgt, int *scalerIncrement, const boolean useFastScaling);
1895 
1896 extern void newviewGTRGAMMAPROT_AVX(int tipCase,
1897  double *x1, double *x2, double *x3, double *extEV, double *tipVector,
1898  int *ex3, unsigned char *tipX1, unsigned char *tipX2, int n,
1899  double *left, double *right, int *wgt, int *scalerIncrement, const boolean useFastScaling);
1900 
1901 extern void newviewGTRCATPROT_AVX(int tipCase, double *extEV,
1902  int *cptr,
1903  double *x1, double *x2, double *x3, double *tipVector,
1904  int *ex3, unsigned char *tipX1, unsigned char *tipX2,
1905  int n, double *left, double *right, int *wgt, int *scalerIncrement, const boolean useFastScaling);
1906 
1907 #endif
1908 
1909 extern int virtual_width( int n );
1910 extern void computeAllAncestralVectors(nodeptr p, pllInstance *tr, partitionList *pr);
1911 
1912 #ifdef __cplusplus
1913 } /* extern "C" */
1914 #endif
1915 
1916 #endif
void pllDestroyRearrangeList(pllRearrangeList **bestList)
Deallocator for topology rearrangements list.
Definition: searchAlgo.c:2214
node * start
Definition: pll.h:1224
void pllEvaluateGeneric(pllInstance *tr, partitionList *pr, nodeptr p, boolean fullTraversal, boolean getPerSiteLikelihoods)
Evaluate the log likelihood of the tree topology.
Definition: evaluateGenericSpecial.c:1412
int pllOptimizeModelParameters(pllInstance *tr, partitionList *pr, double likelihoodEpsilon)
Optimize all free model parameters of the likelihood model.
Definition: utils-merge.c:3185
void pllEvaluateIterative(pllInstance *tr, partitionList *pr, boolean getPerSiteLikelihoods)
Evaluate the log likelihood of a specific branch of the topology.
Definition: evaluateGenericSpecial.c:967
double * gapColumn
Definition: pll.h:1007
boolean noGap(unsigned int *x, int pos)
Check whether the position pos in bitvector x is NOT a gap.
Definition: newviewGenericSpecial.c:4675
int pllNniSearch(pllInstance *tr, partitionList *pr, int estimateModel)
Perform an NNI search.
Definition: searchAlgo.c:2058
void pllQueuePartitionsDestroy(pllQueue **partitions)
Destroy queue structure that contains parsed information from a partition file.
Definition: parsePartition.c:54
void pllStartPthreads(pllInstance *tr, partitionList *pr)
Definition: genericParallelization.c:327
void makenewzIterative(pllInstance *tr, partitionList *pr)
Precompute values (sumtable) from the 2 likelihood vectors of a given branch.
Definition: makenewzGenericSpecial.c:638
boolean isThisMyPartition(partitionList *localPr, int tid, int model)
Check if a partition is assign to a thread/process.
Definition: genericParallelization.c:486
void pllSetFixedAlpha(double alpha, int model, partitionList *pr, pllInstance *tr)
Set the alpha parameter of the Gamma model to a fixed value for a partition.
Definition: utils-merge.c:2600
void pllAlignmentRemoveDups(pllAlignmentData *alignmentData, partitionList *pl)
Remove duplicate sites from alignment and update weights vector.
Definition: utils-merge.c:1331
topol * start
Definition: pll.h:1391
void computeTraversalInfoStlen(nodeptr p, int maxTips, recompVectors *rvec, int *count)
Annotes unoriented tree nodes tr with their subtree size.
Definition: recom.c:411
struct noderec * back
Outer node.
Definition: pll.h:708
void makenewzGeneric(pllInstance *tr, partitionList *pr, nodeptr p, nodeptr q, double *z0, int maxiter, double *result, boolean mask)
Optimize branch length value(s) of a given branch with the Newton-Raphtson procedure.
Definition: makenewzGenericSpecial.c:1151
int * aliaswgt
Definition: pll.h:1150
int slot_p
Definition: pll.h:427
void pllSetFixedSubstitutionMatrix(double *q, int length, int model, partitionList *pr, pllInstance *tr)
Set all substitution rates to a fixed value for a specific partition.
Definition: utils-merge.c:2794
void protectNode(recompVectors *rvec, int nodenum, int mxtips)
Locks node nodenum to force it remains availably in memory.
Definition: recom.c:34
int pllLinkRates(char *string, partitionList *pr)
Link Substitution matrices across partitions.
Definition: utils-merge.c:3054
int descend
Definition: pll.h:1355
void pllAlignmentDataDumpConsole(pllAlignmentData *alignmentData)
Prints the alignment to the console.
Definition: alignment.c:171
unsigned int * globalScaler
Definition: pll.h:994
double z[PLL_NUM_BRANCHES]
Definition: pll.h:1352
boolean permuteTreeoptimize
Definition: pll.h:1292
void pllSetFixedBaseFrequencies(double *f, int length, int model, partitionList *pr, pllInstance *tr)
Set all base freuqncies to a fixed value for a partition.
Definition: utils-merge.c:2652
pllInstance * pllCreateInstance(pllInstanceAttr *)
Create the main instance of PLL.
Definition: utils-merge.c:1749
int nextnode
Definition: pll.h:1369
int nvalid
Definition: pll.h:1395
void addTraverseBIG(pllInstance *tr, partitionList *pr, nodeptr p, nodeptr q, int mintrav, int maxtrav)
Recursively traverse tree and test insertion.
Definition: searchAlgo.c:964
boolean insertBIG(pllInstance *tr, partitionList *pr, nodeptr p, nodeptr q)
Connect two disconnected tree components.
Definition: searchAlgo.c:651
Definition: pll.h:321
Definition: pll.h:738
Definition: pll.h:1110
partitionList * pllPartitionsCommit(pllQueue *parts, pllAlignmentData *alignmentData)
Constructs the proposed partition scheme.
Definition: utils.c:1225
Stores the recomputation-state of likelihood vectors.
Definition: pll.h:332
void pllNewviewGeneric(pllInstance *tr, partitionList *pr, nodeptr p, boolean masked)
Computes the conditional likelihood vectors of all nodes in the subtree rooted at p...
Definition: newviewGenericSpecial.c:2073
int autoProtModels
Best fitted protein model for the PLL_AUTO partitions.
Definition: pll.h:978
void computeTraversal(pllInstance *tr, nodeptr p, boolean partialTraversal, int numBranches)
Compute the traversal descriptor of the subtree rooted at p.
Definition: newviewGenericSpecial.c:2037
void smooth(pllInstance *tr, partitionList *pr, nodeptr p)
Branch length optimization of subtree.
Definition: searchAlgo.c:189
void pllFinalizeMPI(void)
Definition: genericParallelization.c:166
double likelihood
Definition: pll.h:1509
void pllTreeEvaluate(pllInstance *tr, partitionList *pr, int maxSmoothIterations)
Optimize branch lenghts and evaluate likelihood of topology.
Definition: searchAlgo.c:1811
Definition: pll.h:1506
pllNewickTree * pllNewickParseString(const char *newick)
Parse a newick tree string.
Definition: newick.c:470
double likelihood
Definition: pll.h:1220
void pllBaseSubstitute(pllAlignmentData *alignmentData, partitionList *partitions)
Encode the alignment data to the PLL numerical representation.
Definition: utils-merge.c:2102
boolean getxVector(recompVectors *rvec, int nodenum, int *slot, int mxtips)
Get a pinned slot slot that holds the likelihood vector for inner node nodenum.
Definition: recom.c:361
struct noderec * next
Next node in the roundabout.
Definition: pll.h:707
void pllComputeRandomizedStepwiseAdditionParsimonyTree(pllInstance *tr, partitionList *partitions)
Compute a randomized stepwise addition oder parsimony tree.
Definition: utils-merge.c:2083
A brief line.
Definition: pll.h:731
pllRearrangeList * pllCreateRearrangeList(int max)
Create a list for storing topology rearrangements.
Definition: searchAlgo.c:2192
node * q
Definition: pll.h:1353
int numberOfCategories
CAT size of the set of possible categories.
Definition: pll.h:941
Checkpointing states.
Definition: pll.h:1058
void pllNewviewIterative(pllInstance *tr, partitionList *pr, int startIndex)
Compute the conditional likelihood for each entry (node) of the traversal descriptor.
Definition: newviewGenericSpecial.c:1447
void pllPartitionDump(pllQueue *partitions)
Dump a parsed partition file in the console.
Definition: parsePartition.c:268
double * left
P matrix for the left term of the conditional likelihood equation.
Definition: pll.h:960
int pllSetSubstitutionRateMatrixSymmetries(char *string, partitionList *pr, int model)
Set symmetries among parameters in the Q matrix.
Definition: utils-merge.c:2570
void pllStopPthreads(pllInstance *tr)
Definition: genericParallelization.c:379
Definition: pll.h:743
int rearrangeType
Definition: pll.h:1508
void pllDestroyInstance(pllInstance *)
Deallocate the PLL instance.
Definition: utils-merge.c:2249
unsigned int * gapVector
Definition: pll.h:1006
void saveTL(topolRELL_LIST *rl, pllInstance *tr, int index)
Save.
Definition: topologies.c:219
int ** expVector
An entry per inner node. Each element is an array of size the number of sites in the current partitio...
Definition: pll.h:1000
boolean testInsertBIG(pllInstance *tr, partitionList *pr, nodeptr p, nodeptr q)
Test the.
Definition: searchAlgo.c:872
double * empiricalFrequencies
Empirical frequency of each state in the current partition.
Definition: pll.h:933
nodeptr removeNode
Definition: pll.h:1252
nodeptr removeNodeRestoreBIG(pllInstance *tr, partitionList *pr, nodeptr p)
Split the tree into two components and recompute likelihood.
Definition: searchAlgo.c:617
boolean useRecom
Definition: pll.h:1133
char ** sequenceLabels
An array of where the i-th element is the name of the i-th sequence.
Definition: pll.h:1536
int rNumber
Definition: pll.h:423
int pNumber
Definition: pll.h:421
int rightLength
Definition: pll.h:1406
Tree node record.
Definition: pll.h:702
void regionalSmooth(pllInstance *tr, partitionList *pr, nodeptr p, int maxtimes, int region)
Wrapper function for optimizing the branch length of a region maxtimes times.
Definition: searchAlgo.c:510
void pllNewickParseDestroy(pllNewickTree **t)
Deallocate newick parser stack structure.
Definition: newick.c:509
int mxtips
Definition: pll.h:1225
void pllClearRearrangeHistory(pllInstance *tr)
Definition: utils-merge.c:2231
Tree topology.
Definition: pll.h:1125
int pllPartitionsValidate(pllQueue *parts, pllAlignmentData *alignmentData)
Correspondance check between partitions and alignment.
Definition: utils.c:1032
boolean constrained
Definition: pll.h:1154
double lh
this is lh
Definition: pll.h:733
int slot_r
Definition: pll.h:429
void initTL(topolRELL_LIST *rl, pllInstance *tr, int n)
Initializes space as large as the tree.
Definition: topologies.c:148
void initModel(pllInstance *tr, double **empiricalFrequencies, partitionList *partitions)
Initialize the model parameters.
Definition: models.c:4243
int states
Definition: pll.h:1382
int qNumber
Definition: pll.h:422
char ** nameList
Definition: pll.h:1242
???
Definition: pll.h:346
double * EV
Eigenvectors of Q matrix.
Definition: pll.h:958
boolean allSlotsBusy
Definition: pll.h:340
void makeGammaCats(double alpha, double *gammaRates, int K, boolean useMedian)
Compute the gamma rates.
Definition: models.c:3799
int originalCrunchedLength
Definition: pll.h:1203
void resetTL(topolRELL_LIST *rl)
Reset this structure.
Definition: topologies.c:202
int * stlen
Definition: pll.h:337
pllAlignmentData * pllInitAlignmentData(int sequenceCount, int sequenceLength)
Initialize alignment structure fields.
Definition: alignment.c:112
int * frequencyGrouping
Definition: pll.h:986
int pllLoadAlignment(pllInstance *tr, pllAlignmentData *alignmentData, partitionList *, int)
Load alignment to the PLL instance.
Definition: utils.c:1661
int ninit
Definition: pll.h:1396
int initial
Definition: pll.h:1288
void update(pllInstance *tr, partitionList *pr, nodeptr p)
Optimize the length of a specific branch.
Definition: searchAlgo.c:122
void pllRearrangeSearch(pllInstance *tr, partitionList *pr, int rearrangeType, nodeptr p, int mintrav, int maxtrav, pllRearrangeList *bestList)
Search for rearrangement topologies.
Definition: searchAlgo.c:3071
Definition: queue.h:10
Single Topology.
Definition: pll.h:1361
int bDeep
Definition: pll.h:1198
void pllSetBranchLength(pllInstance *tr, nodeptr p, int partition_id, double bl)
Set the length of a specific branch.
Definition: utils-merge.c:858
int pllAlignmentDataDumpFile(pllAlignmentData *alignmentData, int fileFormat, const char *filename)
Dump the alignment to a file of format fileFormat.
Definition: alignment.c:202
int mode
Definition: pll.h:1290
Definition: pll.h:1825
int protModels
Protein model for current partition.
Definition: pll.h:977
int stepwidth
Definition: pll.h:1287
int width
Number of sites in the partition (i.e. upper - lower)
Definition: pll.h:931
double accumulatedTime
Definition: globalVariables.h:34
unsigned char ** yVector
Definition: pll.h:993
int * siteWeights
An array where the i-th element indicates how many times site i appeared (prior to duplicates removal...
Definition: pll.h:1538
void computeFullTraversalInfoStlen(nodeptr p, int maxTips, recompVectors *rvec)
Annotes all tree nodes tr with their subtree size.
Definition: recom.c:491
void unpinNode(recompVectors *v, int nodenum, int mxtips)
Marks node nodenum as unpinnable.
Definition: recom.c:325
int count
Definition: pll.h:440
int nkeep
Definition: pll.h:1394
void newviewAncestralIterative(pllInstance *tr, partitionList *pr)
A very simple iterative function, we only access the conditional likelihood vector at node p...
Definition: newviewGenericSpecial.c:2371
void pllRearrangeCommit(pllInstance *tr, partitionList *pr, pllRearrangeInfo *rearr, int saveRollbackInfo)
Commit a rearrangement move.
Definition: searchAlgo.c:3002
double * substRates
Entries of substitution matrix, e.g. 6 free parameters in DNA.
Definition: pll.h:953
void smoothRegion(pllInstance *tr, partitionList *pr, nodeptr p, int region)
Optimize branch lengths of region.
Definition: searchAlgo.c:466
unsigned char ** yVector
Definition: pll.h:1200
boolean needsRecomp(boolean recompute, recompVectors *rvec, nodeptr p, int mxtips)
Checks if the likelihood entries at node p should be updated.
Definition: recom.c:86
int frequenciesLength
Definition: pll.h:1411
Traversal descriptor entry.
Definition: pll.h:418
Stores data related to a NNI move.
Definition: pll.h:1300
char * tree_string
Definition: pll.h:1243
void smoothTree(pllInstance *tr, partitionList *pr, int maxtimes)
Optimize all branch lenghts of a tree.
Definition: searchAlgo.c:254
int tplNum
Definition: pll.h:1371
int pllRearrangeRollback(pllInstance *tr, partitionList *pr)
Rollback the last committed rearrangement move.
Definition: searchAlgo.c:2954
Definition: stack.h:4
void * valptr
Definition: pll.h:1354
int eignLength
Definition: pll.h:1407
int nextlink
Definition: pll.h:1366
int rearrangeBIG(pllInstance *tr, partitionList *pr, nodeptr p, int mintrav, int maxtrav)
Compute the best SPR movement.
Definition: searchAlgo.c:1008
Intermediate structure for storing a newick tree.
Definition: newick.h:8
float maxMegabytesMemory
Definition: pll.h:1131
double * perSiteRates
Per Site Categories (PSR) or aka CAT values for each rate.
Definition: pll.h:939
char * pllReadFile(const char *, long *)
Read the contents of a file.
Definition: utils-merge.c:3216
pllQueue * pllPartitionParse(const char *filename)
Parse a partition (model) file.
Definition: parsePartition.c:313
boolean insertRestoreBIG(pllInstance *tr, partitionList *pr, nodeptr p, nodeptr q)
Connect two disconnected tree components without optimizing branch lengths.
Definition: searchAlgo.c:757
int saveBestTree(bestlist *bt, pllInstance *tr, int numBranches)
Save the current tree in the bestlist structure.
Definition: topologies.c:693
void pllAlignmentDataDestroy(pllAlignmentData *alignmentData)
Deallocates the memory associated with the alignment data structure.
Definition: alignment.c:148
double fracchange
Definition: pll.h:1209
int pllLinkFrequencies(char *string, partitionList *pr)
Link base frequency parameters across partitions.
Definition: utils-merge.c:3025
void pllTreeInitTopologyForAlignment(pllInstance *tr, pllAlignmentData *alignmentData)
Initialize a tree that corresponds to a given (already parsed) alignment.
Definition: utils-merge.c:2050
double * patrat
Definition: pll.h:1147
Definition: pll.h:377
hashtable * h
Definition: pll.h:1272
Linkage of partitions.
Definition: pll.h:478
char seq_file[1024]
Definition: globalVariables.h:38
nodeptr insertNode
Definition: pll.h:1253
int sequenceCount
Number of sequences.
Definition: pll.h:1534
int max_rearrange
Definition: pll.h:1286
node ** nodep
Definition: pll.h:1222
void putWAG(double *ext_initialRates)
Hardcoded values for the WAG model.
Definition: models.c:89
double * EI
Inverse eigenvectors of Q matrix.
Definition: pll.h:959
char * partitionName
Name of partition.
Definition: pll.h:928
boolean thoroughInsertion
Definition: pll.h:1277
Definition: pll.h:485
void pllPartitionsDestroy(pllInstance *, partitionList **)
free all data structures associated to a partition
Definition: utils-merge.c:933
double worst
Definition: pll.h:1390
double * frequencies
State frequencies, entries in D are initialized as empiricalFrequencies.
Definition: pll.h:954
double ** xVector
Definition: pll.h:992
int * iVector
Definition: pll.h:335
int scrNum
Definition: pll.h:1370
boolean initialSet
Definition: pll.h:1289
int lower
Position of the first site in the alignment that is part of this partition [1, tr-&gt;originalCrunchedLe...
Definition: pll.h:929
Generic structure for storing a multiple sequence alignment.
Definition: pll.h:1532
int sibling
Definition: pll.h:1356
int bestTrav
Definition: pll.h:1285
void printAncestralState(nodeptr p, boolean printStates, boolean printProbs, pllInstance *tr, partitionList *pr)
Prints the ancestral state information for a node p to the terminal.
Definition: newviewGenericSpecial.c:2616
traversalInfo * ti
Definition: pll.h:439
void pllNewviewGenericAncestral(pllInstance *tr, partitionList *pr, nodeptr p)
Computes the conditional likelihood vectors of all nodes in the subtree rooted at p and the marginal ...
Definition: newviewGenericSpecial.c:2493
Definition: pll.h:1334
boolean isTip(int number, int maxTips)
Check whether a node is a tip.
Definition: utils-merge.c:370
double * tipVector
Precomputed (based on current P matrix) conditional likelihood vectors for every possible base...
Definition: pll.h:962
void execCore(pllInstance *tr, partitionList *pr, volatile double *_dlnLdlz, volatile double *_d2lnLdlz2)
Compute first and second derivatives of the likelihood with respect to a given branch length...
Definition: makenewzGenericSpecial.c:770
int * wgt
Weight of site.
Definition: pll.h:932
recompVectors * rvec
Definition: pll.h:1130
double * right
P matrix for the right term of the conditional likelihood equation.
Definition: pll.h:961
int maxTipStates
Number of undetermined states (possible states at the tips)
Definition: pll.h:927
pllNewickTree * pllNewickParseFile(const char *filename)
Parse a newick tree file.
Definition: newick.c:539
Connection within a topology.
Definition: pll.h:1351
unsigned char ** sequenceData
The actual sequence data.
Definition: pll.h:1537
connect * links
Definition: pll.h:1364
int pllSetOptimizeBaseFrequencies(int model, partitionList *pr, pllInstance *tr)
Set that the base freuqencies are optimized under ML.
Definition: utils-merge.c:2710
Structure holding attributes for searching possible tree rearrangements.
Definition: pll.h:1475
double * EIGN
Eigenvalues of Q matrix.
Definition: pll.h:957
small helper data structure for printing out/downstream use of marginal ancestral probability vectors...
Definition: pll.h:1379
int initBestTree(bestlist *bt, int newkeep, int numsp)
Initialize a list of best trees.
Definition: topologies.c:492
void resetBranches(pllInstance *tr)
Reset all branch lengths to default values.
Definition: optimizeModel.c:2335
struct conntyp connect
Connection within a topology.
void initReversibleGTR(pllInstance *tr, partitionList *pr, int model)
Initialize GTR.
Definition: models.c:3438
int pllValidateNewick(pllNewickTree *t)
Validate if a newick tree is a valid phylogenetic tree.
Definition: newick.c:367
int * constraintVector
Definition: pll.h:1227
int NNI(pllInstance *tr, nodeptr p, int swap)
Perform an NNI move.
Definition: searchAlgo.c:1845
int substRatesLength
Definition: pll.h:1410
int number
This is a number Detailed description of number.
Definition: pll.h:734
Definition: pll.h:1312
boolean nonGTR
Definition: pll.h:981
double pllGetBranchLength(pllInstance *tr, nodeptr p, int partition_id)
Get the length of a specific branch.
Definition: utils-merge.c:828
Definition: pll.h:1317
struct ratec rateCategorize
Per-site Rate category entry: likelihood per-site and CAT rate applied ???
int pllLinkAlphaParameters(char *string, partitionList *pr)
Link alpha parameters across partitions.
Definition: utils-merge.c:2993
void localSmooth(pllInstance *tr, partitionList *pr, nodeptr p, int maxtimes)
Optimize the branch length of edges around a specific node.
Definition: searchAlgo.c:303
pllAlignmentData * pllParseAlignmentFile(int fileType, const char *filename)
Parse a file that contains a multiple sequence alignment.
Definition: alignment.c:664
float vectorRecomFraction
Definition: pll.h:1132
int numtrees
Definition: pll.h:1397
int leftLength
Definition: pll.h:1405
int * iNode
Definition: pll.h:336
int * unpinnable
Definition: pll.h:338
Branch length information.
Definition: pll.h:462
int states
Number of states.
Definition: pll.h:926
int recallBestTree(bestlist *bt, int rank, pllInstance *tr, partitionList *pr)
Restore the best tree from bestlist structure.
Definition: topologies.c:776
This is used to look up some hard-coded data for each data type.
Definition: pll.h:1403
int pllInitModel(pllInstance *, partitionList *, pllAlignmentData *)
Initialize partitions according to model parameters.
Definition: utils-merge.c:3088
int slot_q
Definition: pll.h:428
Definition: pll.h:385
int tipCase
Definition: pll.h:420
Per-site Rate category entry: likelihood per-site and CAT rate applied ???
Definition: pll.h:398
int upper
Position of the last site that is part of this partition plus one (i.e. position of the first site th...
Definition: pll.h:930
size_t * expSpaceVector
Each entry represents an inner node and states the size of the corresponding element in expVector...
Definition: pll.h:1001
int number
Node identifier.
Definition: pll.h:711
double alpha
Gamma parameter to be optimized.
Definition: pll.h:943
void pllLockMPI(pllInstance *tr)
Lock the MPI slave processes prior allocating partitions.
Definition: genericParallelization.c:135
int pllNewickUnroot(pllNewickTree *t)
Convert a binary rooted trree to a binary unrooted tree.
Definition: newick.c:423
void hookup(nodeptr p, nodeptr q, double *z, int numBranches)
Connect two nodes and assign branch lengths.
Definition: utils-merge.c:419
Definition: pll.h:1036
boolean grouped
Definition: pll.h:1153
List of topologies.
Definition: pll.h:1388
Definition: pll.h:1431
double * lhs
Definition: pll.h:1146
double * gammaRates
Values of the 4 gamma categories (rates) computed given an alpha.
Definition: pll.h:944
int numVectors
Definition: pll.h:334
double best
Definition: pll.h:1389
int sequenceLength
Length of sequences.
Definition: pll.h:1535
void allocRecompVectorsInfo(pllInstance *tr)
Allocates memory for recomputation structure.
Definition: recom.c:103
int * rateCategory
CAT category index for each site.
Definition: pll.h:940
???Hash tables
Definition: pll.h:370
void pllMasterBarrier(pllInstance *tr, partitionList *pr, int jobType)
A generic master barrier for executing parallel parts of the code.
Definition: genericParallelization.c:1782
nodeptr removeNodeBIG(pllInstance *tr, partitionList *pr, nodeptr p, int numBranches)
Split the tree into two components and optimize new branch length.
Definition: searchAlgo.c:572
double ** pllBaseFrequenciesGTR(partitionList *pl, pllAlignmentData *alignmentData)
Compute the empirical base frequencies for all partitions.
Definition: utils-merge.c:1582
double * probs
Definition: pll.h:1380
Partition information structure.
Definition: pll.h:924
boolean isGap(unsigned int *x, int pos)
Check whether the position pos in bitvector x is a gap.
Definition: newviewGenericSpecial.c:4656
void freeTL(topolRELL_LIST *rl)
Deallocate the space associated with this structure.
Definition: topologies.c:171
void pllTreeInitTopologyRandom(pllInstance *tr, int tips, char **nameList)
Initialize PLL tree with a random topology.
Definition: utils-merge.c:2019
void pllInitMPI(int *argc, char **argv[])
Sets up the MPI environment.
Definition: genericParallelization.c:186
Traversal descriptor.
Definition: pll.h:437
Definition: pll.h:1524
void pllTreeInitTopologyNewick(pllInstance *, pllNewickTree *, int)
Initializes the PLL tree topology according to a parsed newick tree.
Definition: utils.c:2064
Definition: pll.h:1325
void initRateMatrix(pllInstance *tr, partitionList *pr)
Initialize the substitution rates matrix.
Definition: models.c:3894
int dataType
Type of data this partition contains.
Definition: pll.h:925
nodeptr pllGetRandomSubtree(pllInstance *)
Get a random subtree.
Definition: utils-merge.c:540
char c
Definition: pll.h:1381
void getxnode(nodeptr p)
Set the orientation of a node.
Definition: utils-merge.c:391