HMSBEAGLE  1.0.0
libhmsbeagle/BeagleImpl.h
00001 /*
00002  *  BeagleImpl.h
00003  *  BEAGLE
00004  *
00005  * Copyright 2009 Phylogenetic Likelihood Working Group
00006  *
00007  * This file is part of BEAGLE.
00008  *
00009  * BEAGLE is free software: you can redistribute it and/or modify
00010  * it under the terms of the GNU Lesser General Public License as
00011  * published by the Free Software Foundation, either version 3 of
00012  * the License, or (at your option) any later version.
00013  *
00014  * BEAGLE is distributed in the hope that it will be useful,
00015  * but WITHOUT ANY WARRANTY; without even the implied warranty of
00016  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00017  * GNU Lesser General Public License for more details.
00018  *
00019  * You should have received a copy of the GNU Lesser General Public
00020  * License along with BEAGLE.  If not, see
00021  * <http://www.gnu.org/licenses/>.
00022  * 
00023  * @author Andrew Rambaut
00024  * @author Marc Suchard
00025  * @author Daniel Ayres
00026  */
00027 
00028 #ifndef __beagle_impl__
00029 #define __beagle_impl__
00030 
00031 #include "libhmsbeagle/beagle.h"
00032 
00033 #ifdef DOUBLE_PRECISION
00034 #define REAL    double
00035 #else
00036 #define REAL    float
00037 #endif
00038 
00039 namespace beagle {
00040 
00041 class BeagleImpl
00042 {
00043 public:
00044     virtual ~BeagleImpl(){}
00045     
00046     virtual int createInstance(int tipCount,
00047                                int partialsBufferCount,
00048                                int compactBufferCount,
00049                                int stateCount,
00050                                int patternCount,
00051                                int eigenBufferCount,
00052                                int matrixBufferCount,                               
00053                                int categoryCount,
00054                                int scaleBufferCount,
00055                                int resourceNumber,
00056                                long preferenceFlags,
00057                                long requirementFlags) = 0;
00058     
00059     virtual int getInstanceDetails(BeagleInstanceDetails* returnInfo) = 0;
00060     
00061     virtual int setTipStates(int tipIndex,
00062                              const int* inStates) = 0;
00063 
00064     virtual int setTipPartials(int tipIndex,
00065                                const double* inPartials) = 0;
00066     
00067     virtual int setPartials(int bufferIndex,
00068                             const double* inPartials) = 0;
00069     
00070     virtual int getPartials(int bufferIndex,
00071                                                         int scaleIndex,
00072                             double* outPartials) = 0;
00073     
00074     virtual int setEigenDecomposition(int eigenIndex,
00075                                       const double* inEigenVectors,
00076                                       const double* inInverseEigenVectors,
00077                                       const double* inEigenValues) = 0;
00078     
00079     virtual int setStateFrequencies(int stateFrequenciesIndex,
00080                                   const double* inStateFrequencies) = 0;    
00081     
00082     virtual int setCategoryWeights(int categoryWeightsIndex,
00083                                  const double* inCategoryWeights) = 0;
00084     
00085     virtual int setPatternWeights(const double* inPatternWeights) = 0;
00086     
00087     virtual int setCategoryRates(const double* inCategoryRates) = 0;
00088     
00089     virtual int setTransitionMatrix(int matrixIndex,
00090                                     const double* inMatrix,
00091                                     double paddedValue) = 0;
00092 
00093     virtual int setTransitionMatrices(const int* matrixIndices,
00094                                       const double* inMatrices,
00095                                       const double* paddedValues,
00096                                       int count) = 0;    
00097     
00098     virtual int getTransitionMatrix(int matrixIndex,
00099                                     double* outMatrix) = 0;
00100 
00101     virtual int updateTransitionMatrices(int eigenIndex,
00102                                          const int* probabilityIndices,
00103                                          const int* firstDerivativeIndices,
00104                                          const int* secondDerivativeIndices,
00105                                          const double* edgeLengths,
00106                                          int count) = 0;
00107     
00108     virtual int updatePartials(const int* operations,
00109                                int operationCount,
00110                                int cumulativeScalingIndex) = 0;
00111     
00112     virtual int waitForPartials(const int* destinationPartials,
00113                                 int destinationPartialsCount) = 0;
00114     
00115     virtual int accumulateScaleFactors(const int* scalingIndices,
00116                                                                            int count,
00117                                                                            int cumulativeScalingIndex) = 0;
00118     
00119     virtual int removeScaleFactors(const int* scalingIndices,
00120                                      int count,
00121                                      int cumulativeScalingIndex) = 0;   
00122     
00123     virtual int resetScaleFactors(int cumulativeScalingIndex) = 0;   
00124     
00125     virtual int copyScaleFactors(int destScalingIndex,
00126                                  int srcScalingIndex) = 0; 
00127     
00128     virtual int calculateRootLogLikelihoods(const int* bufferIndices,
00129                                             const int* categoryWeightsIndices,
00130                                             const int* stateFrequenciesIndices,
00131                                             const int* scalingFactorsIndices,
00132                                             int count,
00133                                             double* outSumLogLikelihood) = 0;
00134     
00135     virtual int calculateEdgeLogLikelihoods(const int* parentBufferIndices,
00136                                             const int* childBufferIndices,
00137                                             const int* probabilityIndices,
00138                                             const int* firstDerivativeIndices,
00139                                             const int* secondDerivativeIndices,
00140                                             const int* categoryWeightsIndices,
00141                                             const int* stateFrequenciesIndices,
00142                                             const int* scalingFactorsIndices,
00143                                             int count,
00144                                             double* outSumLogLikelihood,
00145                                             double* outSumFirstDerivative,
00146                                             double* outSumSecondDerivative) = 0;
00147     
00148     virtual int getSiteLogLikelihoods(double* outLogLikelihoods) = 0;
00149     
00150     virtual int getSiteDerivatives(double* outFirstDerivatives,
00151                                    double* outSecondDerivatives) = 0;
00152 //protected:
00153     int resourceNumber;
00154 };
00155 
00156 class BeagleImplFactory {
00157 public:
00158     virtual BeagleImpl* createImpl(int tipCount,
00159                                    int partialsBufferCount,
00160                                    int compactBufferCount,
00161                                    int stateCount,
00162                                    int patternCount,
00163                                    int eigenBufferCount,
00164                                    int matrixBufferCount,
00165                                    int categoryCount,
00166                                    int scaleBufferCount,
00167                                    int resourceNumber,
00168                                    long preferenceFlags,
00169                                    long requirementFlags,
00170                                    int* errorCode) = 0; // pure virtual
00171     
00172     virtual const char* getName() = 0; // pure virtual
00173     
00174     virtual const long getFlags() = 0; // pure virtual
00175 };
00176 
00177 } // end namespace beagle
00178 
00179 #endif // __beagle_impl__