source: branches/BottomUpTreeEvaluation/BakedTreeEvaluator.cs @ 328

#region License Information
/* HeuristicLab
 * Copyright (C) 2002-2008 Heuristic and Evolutionary Algorithms Laboratory (HEAL)
 *
 * This file is part of HeuristicLab.
 *
 * HeuristicLab is free software: you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation, either version 3 of the License, or
 * (at your option) any later version.
 *
 * HeuristicLab is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with HeuristicLab. If not, see <http://www.gnu.org/licenses/>.
 */
#endregion

using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using HeuristicLab.DataAnalysis;
using HeuristicLab.Core;
using System.Xml;

namespace HeuristicLab.Functions {
  internal static class BakedTreeEvaluator {
    private const int MAX_TREE_SIZE = 4096;
    private const int MAX_TREE_DEPTH = 20;
    private class Instr {
      public double result;
      public double d_arg0;
      public int i_arg0;
      public int i_arg1;
      public int arity;
      public int symbol;
    }

    private static int[] nInstr;
    private static Instr[,] evaluationTable;
    private static Dataset dataset;
    private static int sampleIndex;


    static BakedTreeEvaluator() {
      evaluationTable = new Instr[MAX_TREE_SIZE, MAX_TREE_DEPTH];
      nInstr = new int[MAX_TREE_DEPTH];
      for(int j = 0; j < MAX_TREE_DEPTH; j++) {
        for(int i = 0; i < MAX_TREE_SIZE; i++) {
          evaluationTable[i, j] = new Instr();
        }
      }
    }

    public static void ResetEvaluator(List<LightWeightFunction> linearRepresentation) {
      int length;
      for(int i = 0; i < MAX_TREE_DEPTH; i++) nInstr[i] = 0;
      //TranslateToInstr(0, linearRepresentation, out length);
      int[] heights = new int[linearRepresentation.Count];
      CalcHeights(linearRepresentation, heights, 0, out length);
      TranslateToTable(0, linearRepresentation, heights);
    }

    private static int CalcHeights(List<LightWeightFunction> linearRepresentation, int[] heights, int p, out int branchLength) {
      if(linearRepresentation[p].arity == 0) {
        heights[p] = 1;
        branchLength = 1;
        return 1;
      }
      int height = 0;
      int length = 1;
      for(int i = 0; i < linearRepresentation[p].arity; i++) {
        int curBranchLength;
        int curHeight = CalcHeights(linearRepresentation, heights, p + length, out curBranchLength);
        if(curHeight > height) {
          height = curHeight;
        }
        length += curBranchLength;
      }
      heights[p] = height+1;
      branchLength = length;
      return height+1;
    }

    private static int TranslateToTable(int pos, List<LightWeightFunction> list, int[] heights) {
      LightWeightFunction f = list[pos];
      if(f.arity == 0) {
        Instr instr = evaluationTable[nInstr[0], 0];
        instr.symbol = EvaluatorSymbolTable.MapFunction(f.functionType);
        switch(instr.symbol) {
          case EvaluatorSymbolTable.VARIABLE: {
              instr.i_arg0 = (int)f.data[0]; // var
              instr.d_arg0 = f.data[1]; // weight
              instr.i_arg1 = (int)f.data[2]; // sample-offset
              break;
            }
          case EvaluatorSymbolTable.CONSTANT: {
              instr.result = f.data[0]; // value
              break;
            }
        }
        nInstr[0]++;
        return 1;
      } else {
        int length = 1;
        int height = heights[pos];
        for(int i = 0; i < f.arity; i++) {
          int curBranchHeight = heights[pos + length];
          if(curBranchHeight < height - 1) {
            for(int j = curBranchHeight; j < height - 1; j++) {
              evaluationTable[nInstr[j], j].symbol = EvaluatorSymbolTable.IDENTITY;
              nInstr[j]++;
            }
          }
          int curBranchLength = TranslateToTable(pos + length, list, heights);
          length += curBranchLength;
        }

        Instr cell = evaluationTable[nInstr[height-1], height-1];
        nInstr[height-1]++;
        cell.arity = f.arity;
        cell.symbol = EvaluatorSymbolTable.MapFunction(f.functionType);
        return length;
      }
    }


    //private static int TranslateToInstr(int pos, List<LightWeightFunction> linearRepresentation, out int branchLength) {
    //  int height = 0;
    //  int length = 1;
    //  LightWeightFunction f = linearRepresentation[pos];
    //  for(int i = 0; i < f.arity; i++) {
    //    int curBranchLength;
    //    int curBranchHeight = TranslateToInstr(pos + length, linearRepresentation, out curBranchLength);
    //    if(curBranchHeight > height) height = curBranchHeight;
    //    length += curBranchLength;
    //  }
    //  Instr instr = evaluationTable[nInstr[height], height];
    //  instr.arity = f.arity;
    //  instr.symbol = EvaluatorSymbolTable.MapFunction(f.functionType);
    //  switch(instr.symbol) {
    //    case EvaluatorSymbolTable.VARIABLE: {
    //        instr.i_arg0 = (int)f.data[0]; // var
    //        instr.d_arg0 = f.data[1]; // weight
    //        instr.i_arg1 = (int)f.data[2]; // sample-offset
    //        break;
    //      }
    //    case EvaluatorSymbolTable.CONSTANT: {
    //        instr.result = f.data[0]; // value
    //        break;
    //      }
    //  }
    //  nInstr[height]++;
    //  branchLength = length;
    //  return height + 1;
    //}

    internal static double Evaluate(Dataset dataset, int sampleIndex) {
      BakedTreeEvaluator.sampleIndex = sampleIndex;
      BakedTreeEvaluator.dataset = dataset;
      return EvaluateTable();
    }

    private static double EvaluateTable() {
      int terminalP = 0;
      // process remaining instr first
      for(int i = 0; i < nInstr[0] % 4; i++) {
        Instr curInstr = evaluationTable[terminalP++, 0];
        if(curInstr.symbol == EvaluatorSymbolTable.VARIABLE) {
          int row = sampleIndex + curInstr.i_arg1;
          if(row < 0 || row >= dataset.Rows) curInstr.result = double.NaN;
          else curInstr.result = curInstr.d_arg0 * dataset.GetValue(row, curInstr.i_arg0);
        }
      }
      // unrolled loop
      for(; terminalP < nInstr[0] - 4; terminalP += 4) {
        Instr curInstr0 = evaluationTable[terminalP, 0];
        Instr curInstr1 = evaluationTable[terminalP + 1, 0];
        Instr curInstr2 = evaluationTable[terminalP + 2, 0];
        Instr curInstr3 = evaluationTable[terminalP + 3, 0];
        if(curInstr0.symbol == EvaluatorSymbolTable.VARIABLE) {
          int row = sampleIndex + curInstr0.i_arg1;
          if(row < 0 || row >= dataset.Rows) curInstr0.result = double.NaN;
          else curInstr0.result = curInstr0.d_arg0 * dataset.GetValue(row, curInstr0.i_arg0);
        }
        if(curInstr1.symbol == EvaluatorSymbolTable.VARIABLE) {
          int row = sampleIndex + curInstr1.i_arg1;
          if(row < 0 || row >= dataset.Rows) curInstr1.result = double.NaN;
          else curInstr1.result = curInstr1.d_arg0 * dataset.GetValue(row, curInstr1.i_arg0);
        }
        if(curInstr2.symbol == EvaluatorSymbolTable.VARIABLE) {
          int row = sampleIndex + curInstr2.i_arg1;
          if(row < 0 || row >= dataset.Rows) curInstr2.result = double.NaN;
          else curInstr2.result = curInstr2.d_arg0 * dataset.GetValue(row, curInstr2.i_arg0);
        }
        if(curInstr3.symbol == EvaluatorSymbolTable.VARIABLE) {
          int row = sampleIndex + curInstr3.i_arg1;
          if(row < 0 || row >= dataset.Rows) curInstr3.result = double.NaN;
          else curInstr3.result = curInstr3.d_arg0 * dataset.GetValue(row, curInstr3.i_arg0);
        }
      }

      int curLevel = 1;
      while(nInstr[curLevel] > 0) {
        int lastLayerInstrP = 0;
        for(int curLayerInstrP = 0; curLayerInstrP < nInstr[curLevel]; curLayerInstrP++) {
          Instr curInstr = evaluationTable[curLayerInstrP, curLevel];
          switch(curInstr.symbol) {
            case EvaluatorSymbolTable.MULTIPLICATION: {
                curInstr.result = evaluationTable[lastLayerInstrP++, curLevel - 1].result;
                for(int i = 1; i < curInstr.arity; i++) {
                  curInstr.result *= evaluationTable[lastLayerInstrP++, curLevel - 1].result;
                }
                break;
              }
            case EvaluatorSymbolTable.ADDITION: {
                curInstr.result = evaluationTable[lastLayerInstrP++, curLevel - 1].result;
                for(int i = 1; i < curInstr.arity; i++) {
                  curInstr.result += evaluationTable[lastLayerInstrP++, curLevel - 1].result;
                }
                break;
              }
            case EvaluatorSymbolTable.SUBTRACTION: {
                if(curInstr.arity == 1) {
                  curInstr.result = -evaluationTable[lastLayerInstrP++, curLevel - 1].result;
                } else {
                  curInstr.result = evaluationTable[lastLayerInstrP++, curLevel - 1].result;
                  for(int i = 1; i < curInstr.arity; i++) {
                    curInstr.result -= evaluationTable[lastLayerInstrP++, curLevel - 1].result;
                  }
                }
                break;
              }
            case EvaluatorSymbolTable.DIVISION: {
                if(curInstr.arity == 1) {
                  curInstr.result = 1.0 / evaluationTable[lastLayerInstrP++, curLevel - 1].result;
                } else {
                  curInstr.result = evaluationTable[lastLayerInstrP++, curLevel - 1].result;
                  for(int i = 1; i < curInstr.arity; i++) {
                    curInstr.result /= evaluationTable[lastLayerInstrP++, curLevel - 1].result;
                  }
                }
                if(double.IsInfinity(curInstr.result)) curInstr.result = 0.0;
                break;
              }
            case EvaluatorSymbolTable.AVERAGE: {
                curInstr.result = evaluationTable[lastLayerInstrP++, curLevel - 1].result;
                for(int i = 1; i < curInstr.arity; i++) {
                  curInstr.result += evaluationTable[lastLayerInstrP++, curLevel - 1].result;
                }
                curInstr.result /= curInstr.arity;
                break;
              }
            case EvaluatorSymbolTable.COSINUS: {
                curInstr.result = Math.Cos(evaluationTable[lastLayerInstrP++, curLevel - 1].result);
                break;
              }
            case EvaluatorSymbolTable.SINUS: {
                curInstr.result = Math.Sin(evaluationTable[lastLayerInstrP++, curLevel - 1].result);
                break;
              }
            case EvaluatorSymbolTable.EXP: {
                curInstr.result = Math.Exp(evaluationTable[lastLayerInstrP++, curLevel - 1].result);
                break;
              }
            case EvaluatorSymbolTable.LOG: {
                curInstr.result = Math.Log(evaluationTable[lastLayerInstrP++, curLevel - 1].result);
                break;
              }
            case EvaluatorSymbolTable.POWER: {
                double x = evaluationTable[lastLayerInstrP++, curLevel - 1].result;
                double p = evaluationTable[lastLayerInstrP++, curLevel - 1].result;
                curInstr.result = Math.Pow(x, p);
                break;
              }
            case EvaluatorSymbolTable.SIGNUM: {
                double value = evaluationTable[lastLayerInstrP++, curLevel - 1].result;
                if(double.IsNaN(value)) curInstr.result = double.NaN;
                else curInstr.result = Math.Sign(value);
                break;
              }
            case EvaluatorSymbolTable.SQRT: {
                curInstr.result = Math.Sqrt(evaluationTable[lastLayerInstrP++, curLevel - 1].result);
                break;
              }
            case EvaluatorSymbolTable.TANGENS: {
                curInstr.result = Math.Tan(evaluationTable[lastLayerInstrP++, curLevel - 1].result);
                break;
              }
            //case EvaluatorSymbolTable.AND: {
            //    double result = 1.0;
            //    // have to evaluate all sub-trees, skipping would probably not lead to a big gain because 
            //    // we have to iterate over the linear structure anyway
            //    for(int i = 0; i < currInstr.arity; i++) {
            //      double x = Math.Round(EvaluateBakedCode());
            //      if(x == 0 || x == 1.0) result *= x;
            //      else result = double.NaN;
            //    }
            //    return result;
            //  }
            //case EvaluatorSymbolTable.EQU: {
            //    double x = EvaluateBakedCode();
            //    double y = EvaluateBakedCode();
            //    if(x == y) return 1.0; else return 0.0;
            //  }
            //case EvaluatorSymbolTable.GT: {
            //    double x = EvaluateBakedCode();
            //    double y = EvaluateBakedCode();
            //    if(x > y) return 1.0;
            //    else return 0.0;
            //  }
            //case EvaluatorSymbolTable.IFTE: {
            //    double condition = Math.Round(EvaluateBakedCode());
            //    double x = EvaluateBakedCode();
            //    double y = EvaluateBakedCode();
            //    if(condition < .5) return x;
            //    else if(condition >= .5) return y;
            //    else return double.NaN;
            //  }
            //case EvaluatorSymbolTable.LT: {
            //    double x = EvaluateBakedCode();
            //    double y = EvaluateBakedCode();
            //    if(x < y) return 1.0;
            //    else return 0.0;
            //  }
            //case EvaluatorSymbolTable.NOT: {
            //    double result = Math.Round(EvaluateBakedCode());
            //    if(result == 0.0) return 1.0;
            //    else if(result == 1.0) return 0.0;
            //    else return double.NaN;
            //  }
            //case EvaluatorSymbolTable.OR: {
            //    double result = 0.0; // default is false
            //    for(int i = 0; i < currInstr.arity; i++) {
            //      double x = Math.Round(EvaluateBakedCode());
            //      if(x == 1.0 && result == 0.0) result = 1.0; // found first true (1.0) => set to true
            //      else if(x != 0.0) result = double.NaN; // if it was not true it can only be false (0.0) all other cases are undefined => (NaN)
            //    }
            //    return result;
            //  }
            //case EvaluatorSymbolTable.XOR: {
            //    double x = Math.Round(EvaluateBakedCode());
            //    double y = Math.Round(EvaluateBakedCode());
            //    if(x == 0.0 && y == 0.0) return 0.0;
            //    if(x == 1.0 && y == 0.0) return 1.0;
            //    if(x == 0.0 && y == 1.0) return 1.0;
            //    if(x == 1.0 && y == 1.0) return 0.0;
            //    return double.NaN;
            //  }
            case EvaluatorSymbolTable.IDENTITY: {
                curInstr.result = evaluationTable[lastLayerInstrP++, curLevel - 1].result;
                break;
              }
            default: {
                throw new NotImplementedException();
              }
          }
        }
        curLevel++;
      }
      return evaluationTable[0, curLevel - 1].result;
    }
  }
}