src/annealing.cpp

Go to the documentation of this file.

/*
  Copyright 2015 Nicolas Melot

  This file is part of Crown.

  Crown 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.

  Crown 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 Crown. If not, see <http://www.gnu.org/licenses/>.

*/


#include <fstream>
#include <set>

#include <cstdlib>

#include <crown/allocation.h>
#include <crown/mapping.h>
#include <crown/scaling.h>
#include <crown/crown.h>

#include <pelib/AmplInput.hpp>
#include <pelib/AmplOutput.hpp>
#include <pelib/Vector.hpp>
#include <pelib/Matrix.hpp>
#include <pelib/Set.hpp>

#include <crown/allocation.h>
#include <crown/CrownScheduler.hpp>

#ifdef debug
#undef debug
#endif

#if 01
#define debug(var) cout << "[" << __FILE__ << ":" << __FUNCTION__ << ":" << __LINE__ << "] " << #var << " = \"" << var << "\"" << endl;
#else
#define debug(var)
#endif

using namespace pelib;
using namespace std;

namespace pelib
{
    namespace crown
    {

        Algebra
        neighbor_allocation(const Algebra &solution, float tasks, float cores)
        {
            Algebra schedule = solution;
            // distance == 1 -> change all neighbors by +- Wi processors
            // distance == 0 -> change only one neighbor by 1 processor

            map<int, float> vector_wi = schedule.find<Vector<int, float> >("wi")->getValues();
            int b = (int)schedule.find<Scalar<float> >("b")->getValue();
            int p = (int)schedule.find<Scalar<float> >("p")->getValue();
            const Vector<int, float> *Wi = schedule.find<Vector<int, float> >("Wi");
            map<int, float> vector_Wi = Wi->getValues();
            int n = vector_wi.size();

            int num_tasks = ((int)((float)rand() / (float)RAND_MAX * n * tasks));
            while (num_tasks > 0 && vector_Wi.size() > 0)
            {
                int size = vector_Wi.size();
                map<int, float>::iterator itask = vector_Wi.begin();
                int index = rand() % size;

                for(int i = 0; i < index; i++)
                {
                    itask++;
                }

                // If a task is sequential, there is no neighbor to find there, remove that task from list and continue
                if(itask->second < 2)
                {
                    vector_Wi.erase(itask);
                    continue;
                }

                // Task maximal width > 1 (Wit > 0)
                int task = itask->first;
                int Wit = (int)vector_Wi[task];
                Wit = Wit > p ? p : Wit;
                Wit = (int)(log(Wit) / log(b));

                // task current allocated width
                int wi = (int)(log(vector_wi.find(task)->second) / log(b));

                int amplitude = Wit * (int)cores; // + 1 to count the number of possibilities, -1 to remove the current allocation from possible new allocation
                amplitude = amplitude < 1 ? 1 : amplitude; // Make sure at least one transformation is performed

                // Count how many spots are available below and above current wi
                int avail_lower = wi;
                int avail_higher = Wit - wi;
        
                // Cound how many spots we need below and above current wi
                int wi_lower, wi_higher;
                if(rand() % 2 == 0)
                {
                    wi_lower = (int)floor((float)amplitude / 2);
                    wi_higher = (int)ceil((float)amplitude / 2);
                }
                else
                {
                    wi_lower = (int)ceil((float)amplitude / 2);
                    wi_higher = (int)floor((float)amplitude / 2);
                }
                
                int len_lower, len_higher;
                if(avail_lower >= wi_lower && avail_higher >=wi_higher)
                {
                    len_lower = wi_lower;
                    len_higher = wi_higher;
                }
                else if(avail_lower < wi_lower)
                {
                    len_lower = avail_lower;
                    len_higher = wi_higher + (wi_lower - avail_lower);
                }
                else if(avail_higher < wi_higher)
                {
                    len_higher = avail_higher;
                    len_lower = wi_lower + (wi_higher - avail_higher);
                }
                else
                {
                    // Not enough space to fit the amplitude, aborting
                    cerr << "# Not enough allocation opportunities to fit the amplitude" << endl;
                    cerr << "# Task " << task << "; p = " << p << "; log2(Wi) = " << Wit << "; log2(wi) = " << wi << "; amplitude = " << amplitude << "; avail_lower = " << avail_lower << "; avail_higher = " << avail_higher << "; wi_lower = " << wi_lower << "; wi_higher = " << wi_higher << endl; 
                    abort();
                }

                int alternative = ((int)rand() % (amplitude)) + 1; // Start counting the new allocation number with 1
                int new_wi = 0;
                if(alternative <= len_lower)
                {
                    new_wi = wi - (len_lower - alternative + 1);
                }
                else if(alternative <= len_lower + len_higher)
                {
                    new_wi = wi + alternative - len_lower;
                }
                else
                {
                    cerr << "# Chosen alternative is beyond the spots available" << endl;
                    cerr << "# Task " << task << "; p = " << p << "; log2(Wi) = " << Wit << "; log2(wi) = " << wi << "; amplitude = " << amplitude << "; avail_lower = " << avail_lower << "; avail_higher = " << avail_higher << "; wi_lower = " << wi_lower << "; wi_higher = " << wi_higher << "; alternative = " << alternative << endl; 
                    len_lower = avail_lower;
                    len_higher = wi_higher + (wi_lower - avail_lower);
                }
                        
                if(wi == new_wi || new_wi > Wit || new_wi < 0)
                {
                    cerr << "# Task " << task << "; Wi " << Wit << "; wi " << wi << "; amplitude " << amplitude << "; new_wi " << new_wi << endl;
                    abort();
                }

                new_wi = (int)pow((double)b, (double)new_wi);
                if(new_wi == 64)
                {
                    cerr << "# Task " << task << "; Wi " << Wit << "; wi " << wi << "; amplitude " << amplitude << "; new_wi " << new_wi << endl;
                }
                vector_wi[task] = new_wi;

                // Done with this task, continue with the next
                vector_Wi.erase(itask);
                num_tasks--;

            }

            Vector<int, float> wi("wi", vector_wi);
            schedule.insert(&wi);

            return schedule;
        }

        Algebra
        neighbor_allocation_full(const Algebra &solution, float distance)
        {
            return neighbor_allocation(solution, distance, distance);
        }

        Algebra
        neighbor_allocation_tasks(const Algebra &solution, float distance)
        {
            return neighbor_allocation(solution, distance, 0);
        }

        Algebra
        neighbor_allocation_cores(const Algebra &solution, float distance)
        {
            return neighbor_allocation(solution, 0, distance);
        }

        Algebra
        neighbor_efficiency(const Algebra &solution, float distance)
        {
            Algebra schedule = solution;
        
            const Scalar<float> *gemin = schedule.find<Scalar<float> >("gemin");
            if(gemin == NULL)
            {
                gemin = new Scalar<float>("gemin", 0);
                schedule.insert(gemin);
                delete gemin;
                gemin = schedule.find<Scalar<float> >("gemin");
            }
            float gemin_val = gemin->getValue();

            const Scalar<float> *gemax = schedule.find<Scalar<float> >("gemax");
            if(gemax == NULL)
            {
                gemax = new Scalar<float>("gemax", 1);
                schedule.insert(gemax);
                delete gemax;
                gemax = schedule.find<Scalar<float> >("gemax");
            }
            float gemax_val = gemax->getValue();

            // The whole universe we can vary the global minimal efficiency in 
            float universe = gemax_val - gemin_val;

            const Scalar<float> *geuse = schedule.find<Scalar<float> >("geuse");
            float geuse_val = 0;
            if(geuse == NULL)
            {
                // The initial value is somewhere between gemin and gemax ( = gemin + rand() * universe )
                geuse = new Scalar<float>("geuse", gemin_val + (universe * ((float)rand() / RAND_MAX)));
                schedule.insert(geuse);
                delete geuse;
                geuse = schedule.find<Scalar<float> >("geuse");
            }
            else
            {
                geuse_val = geuse->getValue();
                float random = (float)rand() / RAND_MAX;
                geuse_val = geuse_val + (random * universe * distance);

                if(geuse_val < gemin_val)
                {
                    geuse_val += gemin_val - geuse_val;
                }
                else if(geuse_val > gemax_val)
                {
                    geuse_val -= geuse_val - gemax_val;
                }

                if (geuse_val < gemin_val || geuse_val > gemax_val)
                {
                    // Something went terribly wrong
                    cerr << "Geuse val went off the chart despite (or because of) the out-of-universe problem fix" << endl;
                    cerr << "geuse_val = " << geuse_val << "; gemin = " << gemin_val << "; gemax = " << gemax_val << "; universe = " << universe << "; random = " << random << ";distance = " << distance << endl; 
                }
            }

            // Reallocate
            return allocation_gemin(schedule, geuse_val);
        }

        float
        quality_simple(const Algebra &schedule)
        {
            //int n = schedule.find<Scalar<float> >("n")->getValue();
            int p = (int)schedule.find<Scalar<float> >("p")->getValue();
            float alpha = schedule.find<Scalar<float> >("alpha")->getValue();
            float zeta = schedule.find<Scalar<float> >("zeta")->getValue();
            float M = schedule.find<Scalar<float> >("M")->getValue();
            Vector<int, float> tau = schedule.find<Vector<int, float> >("Tau");
            Vector<int, float> Wi = schedule.find<Vector<int, float> >("Wi");
            Vector<int, float> wi = schedule.find<Vector<int, float> >("wi");
            Vector<int, float> frequency = schedule.find<Vector<int, float> >("frequency");
            Vector<int, float> group = schedule.find<Vector<int, float> >("group");
            Matrix<int, int, float> e = schedule.find<Matrix<int, int, float> >("e");
        
            double static_energy = M * p;
            double dynamic_energy = 0;

            // Now evaluate the quality of this schedule
            for(map<int, float>::const_iterator i = group.getValues().begin(); i != group.getValues().end(); i++)
            {
                int task = i->first;

                float tau_t = tau.find(task);
                int width_t = (int)wi.find(task);
                width_t = width_t > p ? p : width_t;
                int f_t = (int)frequency.find(task);

                double area = (double)tau_t / (double)(e.find(width_t, task)) / (double)frequency.find(task);
                static_energy -= area;
                dynamic_energy += area * pow(f_t, alpha);
            }

            return dynamic_energy + zeta * static_energy;
        }

        float
        energy_static_dynamic_busy_idle_active(const Algebra &solution)
        {
            Algebra rec = CrownScheduler::crownToSchedule(solution);
            int n = rec.find<Scalar<float> >("n")->getValue();
            int p = rec.find<Scalar<float> >("p")->getValue();
            float M = rec.find<Scalar<float> >("M")->getValue();
            float kappa = rec.find<Scalar<float> >("kappa")->getValue();
            float eta = rec.find<Scalar<float> >("eta")->getValue();
            float zeta = rec.find<Scalar<float> >("zeta")->getValue();
            float alpha = rec.find<Scalar<float> >("alpha")->getValue();
            Vector<int, float> tau = rec.find<Vector<int, float> >("Tau");
            const Vector<int, string> *nameptr = rec.find<Vector<int, string> >("name");
            Vector<int, string> name("name", nameptr != NULL ? nameptr->getValues() : map<int, string>());
        
            Vector<int, float> Wi = rec.find<Vector<int, float> >("Wi");
            Vector<int, float> frequency = rec.find<Vector<int, float> >("frequency");
            Set<float> F = rec.find<Set<float> >("F");
            Matrix<int, int, float> e = rec.find<Matrix<int, int, float> >("e");
            Matrix<int, int, float> schedule = rec.find<Matrix<int, int, float> >("schedule");
            Vector<int, float> wi = rec.find<Vector<int, float> >("wi");
            vector<double> busy_task(n, 0);
            vector<double> area_task(n, 0);
            vector<double> met_task(n, 0);
        
            float lowest_freq = *(F.getValues().begin());
            kappa = kappa * lowest_freq;
            eta = eta * (zeta * pow(lowest_freq, alpha) + (1 - zeta) * (lowest_freq + kappa)) / zeta;
            double idle_energy = p * M * pow(lowest_freq, alpha);
            double busy_energy = 0;

            map<float, int> core_active;
            for(size_t i = 0; i < (size_t)p; i++)
            {
                core_active[i + 1] = 0;
            }
            int total_core_active = 0;

            // Iterate through all task names until we find the dummy task;
            vector<size_t> dummy_index;
            size_t index = 0;
            for(map<int, string>::const_iterator i = name.getValues().begin(); i != name.getValues().end(); i++, index++)
            {
                stringstream ss;
                ss << DUMMY_TASK << " " << index + 1;
                if(string(i->second).compare(ss.str()) == 0)
                {
                    dummy_index.push_back(index + 1);   
                }
            }

            // Now evaluate the quality of this schedule
            float total_area = 0;
            for(map<int, map<int, float> >::const_iterator i = schedule.getValues().begin(); i != schedule.getValues().end(); i++)
            {
                //int core = i->first;
                for(map<int, float>::const_iterator j = i->second.begin(); j != i->second.end(); j++)
                {
                    //int index = j->first;
                    int task = j->second;
                    //cout << "Core " << core << ", task " << task << endl;

                    // Check if the task index is a dummy task
                    bool dummy = false;
                    for(size_t k = 0; k < dummy_index.size(); k++)
                    {
                        if((size_t)task == dummy_index[k])
                        {
                            dummy = true;
                            break;
                        }
                    }

                    // We don't take dummy task power consumption into account
                    if(task > 0 && !dummy)
                    {
                        float width_t = wi.getValues().find(task)->second;
                        double time = (double)tau.find(task) / (width_t * e.find(width_t, task) * (double)frequency.find(task));
                        total_area += time;
                        float power = (zeta * pow(frequency.find(task), alpha) + (1 - zeta) * (frequency.find(task) + kappa)) / zeta;
                        busy_energy += time * power;

                        // Mark this core as active
                        if(core_active[i->first] == 0)
                        {
                            core_active[i->first] = 1;
                            total_core_active++;
                        }
                    }
                }
            }

            idle_energy = (M * (total_core_active) - total_area) * eta;

            return zeta * busy_energy + (1 - zeta) * idle_energy;
        }

        static
        float
        energy_static_idle(const Algebra &schedule, float epsilon)
        {
            int p = (int)schedule.find<Scalar<float> >("p")->getValue();
            float alpha = schedule.find<Scalar<float> >("alpha")->getValue();
            float zeta = schedule.find<Scalar<float> >("zeta")->getValue();
            float minF = *(schedule.find<Set<float> >("F")->getValues().begin());
        
            float M;
            const Scalar<float> *float_M = schedule.find<Scalar<float> >("M");
            if(float_M != NULL)
            {
                M = float_M->getValue();
            }
            else
            {
                const Scalar<float> *int_M = schedule.find<Scalar<float> >("M");
                M = int_M->getValue();
            }

            Vector<int, float> tau = schedule.find<Vector<int, float> >("Tau");
            Vector<int, float> Wi = schedule.find<Vector<int, float> >("Wi");
            Vector<int, float> wi = schedule.find<Vector<int, float> >("wi");
            Vector<int, float> frequency = schedule.find<Vector<int, float> >("frequency");
            Vector<int, float> group = schedule.find<Vector<int, float> >("group");
            Matrix<int, int, float> e = schedule.find<Matrix<int, int, float> >("e");
        
            double idle_energy = 0;
            double busy_energy = 0;
            double total_area = 0;

            // Now evaluate the quality of this schedule
            for(map<int, float>::const_iterator i = group.getValues().begin(); i != group.getValues().end(); i++)
            {
                int task = i->first;

                float tau_t = tau.find(task);
                float width_t = wi.find(task);
                width_t = width_t > p ? p : width_t;
                //double e_wt = e.find(width_t, task);
                float f_t = frequency.find(task);

                double area = (double)tau_t / (double)(e.find((int)width_t, task)) / (double)frequency.find(task);

                total_area += area;
                busy_energy += area * pow(f_t, alpha) + zeta * f_t;
            }
            idle_energy = (p * M - total_area) * pow(minF, alpha) + zeta * minF;

            return busy_energy + epsilon * idle_energy;
        }

        float
        quality_generic(const Algebra &schedule)
        {
            int p = (int)schedule.find<Scalar<float> >("p")->getValue();
            float alpha = schedule.find<Scalar<float> >("alpha")->getValue();
            float alpha_static = schedule.find<Scalar<float> >("alpha-static")->getValue();
            float zeta = schedule.find<Scalar<float> >("zeta")->getValue();
            float epsilon = schedule.find<Scalar<float> >("epsilon")->getValue();
            float minF = *(schedule.find<Set<float> >("F")->getValues().begin());
        
            float M;
            const Scalar<float> *float_M = schedule.find<Scalar<float> >("M");
            if(float_M != NULL)
            {
                M = float_M->getValue();
            }
            else
            {
                const Scalar<float> *int_M = schedule.find<Scalar<float> >("M");
                M = int_M->getValue();
            }

            Vector<int, float> tau = schedule.find<Vector<int, float> >("Tau");
            Vector<int, float> Wi = schedule.find<Vector<int, float> >("Wi");
            Vector<int, float> wi = schedule.find<Vector<int, float> >("wi");
            Vector<int, float> frequency = schedule.find<Vector<int, float> >("frequency");
            Vector<int, float> group = schedule.find<Vector<int, float> >("group");
            Matrix<int, int, float> e = schedule.find<Matrix<int, int, float> >("e");
        
            double idle_energy = 0;
            double busy_energy = 0;
            double total_area = 0;

            // Now evaluate the quality of this schedule
            for(map<int, float>::const_iterator i = group.getValues().begin(); i != group.getValues().end(); i++)
            {
                int task = i->first;

                float tau_t = tau.find(task);
                float width_t = wi.find(task);
                width_t = width_t > p ? p : width_t;
                float f_t = frequency.find(task);

                double area = (double)tau_t / (double)(e.find((int)width_t, task)) / (double)frequency.find(task);

                total_area += area;
                busy_energy += area * pow(f_t, alpha) + zeta * f_t;
            }
            idle_energy = (p * M - total_area) * pow(minF, alpha_static) + zeta * minF;

            return busy_energy + epsilon * idle_energy;
        }

        float
        quality_idle(const Algebra &schedule)
        {
            return energy_static_idle(schedule, 0.5);
        }

        float
        quality_off(const Algebra &schedule)
        {
            return energy_static_idle(schedule, 0);
        }

// This is not possible now that we allow energy saving
// by switchingu nused cores off
#if RETURN_EARLY
#define return_now_if_best(solution) if(solution.find<Scalar<float> >("all_sequential")->getValue() == 1 && solution.find<Scalar<float> >("frequency_optimal")->getValue() == 1) return solution
#else
#define return_now_if_best(solution)
#endif

        Algebra
        annealing(const Algebra &initial, const Algebra &initial_best, float (*quality)(const Algebra &solution), Algebra (*neighbor)(const Algebra &solution, float distance), float temperature, float cooling, float final, int max_transform, int max_new_state, float distance)
        {
            Algebra schedule = initial;
            schedule = neighbor(schedule, distance);

            float m;
            int counter = 0;
            float T = temperature;

            Algebra best = initial_best;
            float best_q = quality(best);
            Vector<int, float> *best_solution = initial_best.find<Vector<int, float> >("wi")->clone();
            float q = quality(schedule);

            while (T > final)
            {
                int new_state = 0;

                for(int i = 0; i < max_transform; i++)
                {
                    // Keep track of the solution of last iteration
                    Vector<int, float> *old_solution = schedule.find<Vector<int, float> >("wi")->clone();
                    float old_q = q;

                    // Transform the solution
                    schedule = neighbor(schedule, distance);
                    schedule = mapping_ltlg(schedule);
                    schedule = frequency_height(schedule);

                    m = schedule.find<Scalar<float> >("m")->getValue();
                    q = quality(schedule);
                        
                    // If the new energy is lower than the previous one, then we accept this solution
                    if (q - old_q <= 0 && m >= 0)
                    {
                        delete old_solution;

                        // If this solution is better than the best found so far, then keep this best solution
                        if (q < best_q)
                        {
                            best_solution = schedule.find<Vector<int, float> >("wi")->clone();
                            best_q = q;
                        }

                        // One new state validated
                        new_state++;
                    }
                    else
                    {
                        // Randomly accept a degrading solution, more likely when the temperature is higher
                        if (rand() / RAND_MAX <= exp(-(q - old_q) / T) && m >= 0)
                        {
                            // Update current solution
                            delete old_solution;

                            // One new state validated
                            new_state++;
                        }
                        else
                        {
                            // Revert old solution 
                            schedule.insert(old_solution);
                            q = old_q;
                        }
                    }
                        
                    // Too much successful transformations
                    if(new_state == max_new_state) 
                    {
                        break;
                    }
                }
        
                T = T * cooling;
                counter++;
            }

            schedule.insert(best_solution);
            schedule = mapping_ltlg(schedule);
            schedule = frequency_height(schedule);

            return schedule;
        }

        float
        annealing_complexity(float temperature, float cooling, float final, int max_transform, int max_new_state, float distance)
        {
            return max_transform * ceil(log(temperature / final) / log(1 / cooling));
        }

    }
}