praktikum-holons/src/holeg/algorithm/topologie/PsoAlgorithm.java

package holeg.algorithm.topologie;

import holeg.algorithm.objective_function.TopologieObjectiveFunction;
import holeg.api.TopologieAlgorithmFramework;
import holeg.model.Model;
import holeg.utility.math.Maths;
import holeg.utility.math.Random;
import java.util.ArrayList;
import java.util.Arrays;
import java.util.LinkedList;
import java.util.List;
import java.util.TreeSet;
import java.util.stream.Collectors;

public class PsoAlgorithm extends TopologieAlgorithmFramework {

  private int popsize = 20;
  private int maxGenerations = 500;
  private double dependency = 2.07;
  private double c1, c2, w;
  private double maxVelocity = 10.0;
  private double deviation = 0.5;
  //mutation
  private int mutationInterval = 1;
  private boolean useIntervalMutation = false;
  private double mutationRate = 0.005;
  private double mutateProbabilityInterval = 0.01;
  private double maxMutationPercent = 0.01;

  private boolean moreInformation = false;

  public PsoAlgorithm() {
    addIntParameter("Swarmsize", popsize, intValue -> popsize = intValue, () -> popsize, 1);
    addIntParameter("Iterations", maxGenerations, intValue -> maxGenerations = intValue,
        () -> maxGenerations, 1);
    addSeperator();
    addDoubleParameter("Deviation", deviation, doubleValue -> deviation = doubleValue,
        () -> deviation, 0);
    addDoubleParameter("Dependency", dependency, doubleValue -> dependency = doubleValue,
        () -> dependency, true, 2.001, 2.4);
    addDoubleParameter("Particle Max-Velocity", maxVelocity,
        doubleValue -> maxVelocity = doubleValue, () -> maxVelocity, 0.0);
    addSeperator();
    addIntParameter("Mutation Frequency (Iteration)", mutationInterval,
        intValue -> mutationInterval = intValue, () -> mutationInterval, 0);
    addBooleanParameter("Use Interval Mutation", useIntervalMutation,
        booleanValue -> useIntervalMutation = booleanValue,
        Arrays.asList("Probability for Frequency Mutation", "Scope of Mutation"),
        Arrays.asList("Probability for Bit-wise Mutation"));
    addDoubleParameter("Probability for Frequency Mutation", mutateProbabilityInterval,
        doubleValue -> mutateProbabilityInterval = doubleValue, () -> mutateProbabilityInterval,
        useIntervalMutation, 0.0, 1.0);
    addDoubleParameter("Probability for Bit-wise Mutation", mutationRate,
        doubleValue -> mutationRate = doubleValue, () -> mutationRate, !useIntervalMutation, 0.0,
        1.0);
    addDoubleParameter("Scope of Mutation", maxMutationPercent,
        doubleValue -> maxMutationPercent = doubleValue, () -> maxMutationPercent,
        useIntervalMutation, 0.0, 1.0);
    addSeperator();
    addBooleanParameter("Print Auxiliary Information", moreInformation,
        booleanValue -> moreInformation = booleanValue, new LinkedList<String>(),
        new LinkedList<String>());
  }


  @Override
  protected double evaluateState(Model model, int amountOfAddedSwitch, double addedCableMeters,
      boolean moreInformation) {
    return TopologieObjectiveFunction.getFitnessValueForState(model, amountOfAddedSwitch,
        addedCableMeters, moreInformation);
  }

  @Override
  /**
   *  <p>Algo from Paper:</p><font size="3"><pre>
   *
   *  Begin
   *  	t = 0; {t: generation index}
   *  	initialize particles x<sub>p,i,j</sub>(t);
   *  	evaluation x<sub>p,i,j</sub>(t);
   *  	while (termination condition &ne; true) do
   *  		v<sub>i,j</sub>(t) = update v<sub>i,j</sub>(t); {by Eq. (6)}
   *  		x<sub>g,i,j</sub>(t) = update x<sub>g,i,j</sub>(t); {by Eq. (7)}
   *  		x<sub>g,i,j</sub>(t) = mutation x<sub>g,i,j</sub>(t); {by Eq. (11)}
   *  		x<sub>p,i,j</sub>(t) = decode x<sub>g,i,j</sub>(t); {by Eqs. (8) and (9)}
   *  		evaluate x<sub>p,i,j</sub>(t);
   *  		t = t + 1;
   *  	end while
   *  End</pre></font>
   *  <p>with:</p><font size="3">
   *
   *  x<sub>g,i,j</sub>: genotype ->genetic information -> in continuous space<br>
   *  x<sub>p,i,j</sub>: phenotype -> observable characteristics-> in binary space<br>
   *  X<sub>g,max</sub>: is the Maximum here set to 4.<br>
   *  Eq. (6):v<sub>i,j</sub>(t + 1) = wv<sub>i,j</sub>+c<sub>1</sub>R<sub>1</sub>(P<sub>best,i,j</sub>-x<sub>p,i,j</sub>(t))+c<sub>2</sub>R<sub>2</sub>(g<sub>best,i,j</sub>-x<sub>p,i,j</sub>(t))<br>
   *  Eq. (7):x<sub>g,i,j</sub>(t + 1) = x<sub>g,i,j</sub>(t) + v<sub>i,j</sub>(t + 1)<br>
   *  Eq. (11):<b>if(</b>rand()&lt;r<sub>mu</sub><b>)then</b> x<sub>g,i,j</sub>(t + 1) = -x<sub>g,i,j</sub>(t + 1)<br>
   *  Eq. (8):x<sub>p,i,j</sub>(t + 1) = <b>(</b>rand() &lt; S(x<sub>g,i,j</sub>(t + 1))<b>) ?</b> 1 <b>:</b> 0<br>
   *  Eq. (9) Sigmoid:S(x<sub>g,i,j</sub>(t + 1)) := 1/(1 + e<sup>-x<sub>g,i,j</sub>(t + 1)</sup>)<br></font>
   *  <p>Parameter:</p>
   *  w inertia, calculated from phi(Variable:{@link #dependency})<br>
   *  c1:	influence, calculated from phi(Variable:{@link #dependency}) <br>
   *  c2:	influence, calculated from phi(Variable:{@link #dependency})<br>
   *  r<sub>mu</sub>: probability that the proposed operation is conducted defined by limit(Variable:{@link #limit})<br>
   *
   *
   */
  protected Individual executeAlgo() {
    resetWildcards();
    initDependentParameter();
    Individual globalBest = new Individual();
    globalBest.position = extractPositionAndAccess();
    globalBest.fitness = evaluatePosition(globalBest.position);
    console.println("Start Value:" + globalBest.fitness);
    int dimensions = globalBest.position.size();
    List<Particle> swarm = initializeParticles(dimensions);
    List<Double> runList = new ArrayList<Double>();
    runList.add(globalBest.fitness);
    evaluation(globalBest, swarm);
    runList.add(globalBest.fitness);
    for (int iteration = 0; iteration < this.maxGenerations; iteration++) {
      int mutationAllowed = iteration % mutationInterval;
      double bitsFlipped = 0;
      for (int particleNumber = 0; particleNumber < this.popsize; particleNumber++) {
        Particle particle = swarm.get(particleNumber);

        if (this.useIntervalMutation) {
          boolean allowMutation = (Random.nextDouble() < this.mutateProbabilityInterval);
          TreeSet<Integer> mutationLocationSet = null;
          if (allowMutation) {
            mutationLocationSet = locationsToMutate(dimensions);
          }
          for (int index = 0; index < dimensions; index++) {
            updateVelocity(particle, index, globalBest);
            updateGenotype(particle, index);
            if (allowMutation && mutationAllowed == 0 && iteration != 0
                && mutationLocationSet.contains(index)) {
              mutation(particle, index);
            }
            decode(particle, index);
          }
        } else {
          int count = 0;
          for (int index = 0; index < dimensions; index++) {
            updateVelocity(particle, index, globalBest);
            updateGenotype(particle, index);
            if (mutationAllowed == 0 && iteration != 0 && Random.nextDouble() < mutationRate) {
              count++;
              mutation(particle, index);
            }
            decode(particle, index);
          }
          bitsFlipped += count;
        }
      }
      if (moreInformation) {
        console.println("\t\t\t\t\t\tAvgBitsMutate: " + (bitsFlipped / (double) popsize));
      }
      if (cancel) {
        return null;
      }
      evaluation(globalBest, swarm);
      runList.add(globalBest.fitness);
      if (moreInformation) {
        console.println("------------------------");
      }
    }
    console.println(" End Value:" + globalBest.fitness);
    this.runList = runList;
    return globalBest;
  }

  @Override
  protected int getProgressBarMaxCount() {
    return rounds * maxGenerations * popsize + 1;
  }

  @Override
  protected String algoInformationToPrint() {
    return "PsoAlgo" + " Rounds:" + rounds
        + " maxIterations:" + maxGenerations
        + " swarmSize:" + popsize
        + " dependency:" + dependency
        + " mutationInterval:" + mutationInterval
        + " maxVelocity: " + maxVelocity
        + " deviation: " + deviation
        + (useIntervalMutation ?
        (" mutateProbabilityInterval:" + mutateProbabilityInterval
            + " maxMutationPercent:" + maxMutationPercent) : " mutationRate:" + mutationRate);
  }

  @Override
  protected String plottFileName() {
    return "pso-topologie.txt";
  }


  /**
   * @param problemSize maximum index of position in the particle
   * @return
   */
  private List<Particle> initializeParticles(int problemSize) {
    List<Particle> swarm = new ArrayList<Particle>();
    //Create The Particle
    for (int particleNumber = 0; particleNumber < popsize; particleNumber++) {
      //Create a Random position
      List<Integer> aRandomPosition = new ArrayList<Integer>();
      for (int index = 0; index < problemSize; index++) {
        aRandomPosition.add(Random.nextIntegerInRange(0, this.getMaximumIndexObjects(index) + 1));
      }
      swarm.add(new Particle(aRandomPosition));
    }
    return swarm;
  }

  /**
   * Calculate w, c1, c2
   */
  private void initDependentParameter() {
    w = 1.0 / (dependency - 1 + Math.sqrt(dependency * dependency - 2 * dependency));
    c1 = c2 = dependency * w;
  }

  /**
   * Evaluate each particle and update the global Best position;
   *
   * @param globalBest
   * @param swarm
   */
  private void evaluation(Individual globalBest, List<Particle> swarm) {
    for (Particle p : swarm) {
      double localEvaluationValue = evaluatePosition(p.xPhenotype);
      if (moreInformation) {
        console.println("Fitness " + localEvaluationValue);
      }
      p.checkNewEvaluationValue(localEvaluationValue);
      if (localEvaluationValue < globalBest.fitness) {
        globalBest.fitness = localEvaluationValue;
        globalBest.position = p.localBest.position;
      }
    }
  }

  /**
   * Eq. (6):v<sub>i,j</sub>(t + 1) = wv<sub>i,j</sub>+c<sub>1</sub>R<sub>1</sub>(P<sub>best,i,j</sub>-x<sub>p,i,j</sub>(t))+c<sub>2</sub>R<sub>2</sub>(g<sub>best,i,j</sub>-x<sub>p,i,j</sub>(t))<br>
   *
   * @param particle
   * @param index
   * @param globalBest
   */
  private void updateVelocity(Particle particle, int index, Individual globalBest) {
    double r1 = Random.nextDouble();
    double r2 = Random.nextDouble();
    double posValue = particle.xPhenotype.get(index);
    particle.velocity.set(index, clamp(
        w * particle.velocity.get(index) + c1 * r1 * ((particle.localBest.position.get(index))
            - posValue) + c2 * r2 * ((globalBest.position.get(index)) - posValue)));
  }

  /**
   * Eq. (7):x<sub>g,i,j</sub>(t + 1) = x<sub>g,i,j</sub>(t) + v<sub>i,j</sub>(t + 1)<br>
   *
   * @param particle
   * @param index
   */
  private void updateGenotype(Particle particle, int index) {
    particle.xGenotype.set(index,
        clamp(particle.xGenotype.get(index) + particle.velocity.get(index)));
  }

  /**
   * Eq. (11):<b>if(</b>rand()&lt;r<sub>mu</sub><b>)then</b> x<sub>g,i,j</sub>(t + 1) =
   * -x<sub>g,i,j</sub>(t + 1)<br>
   *
   * @param particle
   * @param index
   */
  private void mutation(Particle particle, int index) {
    //if(Random.nextDouble() < limit)
    particle.xGenotype.set(index, -particle.xGenotype.get(index));
  }

  /**
   * Eq. (8):x<sub>p,i,j</sub>(t + 1) = <b>(</b>rand() &lt; S(x<sub>g,i,j</sub>(t + 1))<b>) ?</b> 1
   * <b>:</b> 0<br>
   *
   * @param particle
   * @param index
   */
  private void decode(Particle particle, int index) {
    double value = clamp(Random.nextGaussian(particle.xGenotype.get(index), deviation));
    double alpha = Maths.invLerp(-maxVelocity, +maxVelocity, value);
    double result = Maths.lerp(0, this.getMaximumIndexObjects(index), alpha);
    particle.xPhenotype.set(index, (int) Math.round(result));
  }

  /**
   * To clamp X<sub>g,j,i</sub> and v<sub>i,j</sub> in Range [-X<sub>g,max</sub>|+X<sub>g,max</sub>]
   * with {X<sub>g,max</sub>= 4}
   *
   * @param value
   * @return
   */
  private double clamp(double value) {
    return Math.max(-maxVelocity, Math.min(maxVelocity, value));
  }

  private TreeSet<Integer> locationsToMutate(int dimensions) {
    TreeSet<Integer> mutationLocation = new TreeSet<Integer>(); //sortedSet
    int maximumAmountOfMutatedBits = Math.max(1,
        (int) Math.round(((double) dimensions) * this.maxMutationPercent));
    int randomUniformAmountOfMutatedValues = Random.nextIntegerInRange(1,
        maximumAmountOfMutatedBits + 1);
    for (int i = 0; i < randomUniformAmountOfMutatedValues; i++) {
      boolean success = mutationLocation.add(Random.nextIntegerInRange(0, dimensions));
      if (!success) {
        i--; //can be add up to some series long loops if maximumAmountOfMutatedBits get closed to problemsize.
      }
    }
    return mutationLocation;
  }


  /**
   * Class to represent a Particle.
   */
  private class Particle {

    /**
     * The velocity of a particle.
     */
    public List<Double> velocity;
    /**
     * The positions genotype.
     */
    public List<Double> xGenotype;
    /**
     * The positions phenotype. Alias the current position.
     */
    public List<Integer> xPhenotype;

    public Individual localBest;

    Particle(List<Integer> position) {
      this.xPhenotype = position;
      //Init velocity, xGenotype with 0.0 values.
      this.velocity = position.stream().map(bool -> 0.0).collect(Collectors.toList());
      this.xGenotype = position.stream().map(bool -> 0.0).collect(Collectors.toList());
      localBest = new Individual();
      localBest.fitness = Double.MAX_VALUE;
    }

    public void checkNewEvaluationValue(double newEvaluationValue) {
      if (newEvaluationValue < localBest.fitness) {
        localBest.fitness = newEvaluationValue;
        localBest.position = xPhenotype.stream().collect(Collectors.toList());
      }
    }

    public String toString() {
      return "Particle with xPhenotype(Position), xGenotype, velocity:["
          + listToString(xPhenotype) + "],[" + listToString(xGenotype) + "],["
          + listToString(velocity) + "]";
    }

    private <Type> String listToString(List<Type> list) {
      return list.stream().map(Object::toString).collect(Collectors.joining(", "));
    }

  }

}