/*
 *    vmg - a versatile multigrid solver
 *    Copyright (C) 2012 Institute for Numerical Simulation, University of Bonn
 *
 *  vmg 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.
 *
 *  vmg 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 this program.  If not, see <http://www.gnu.org/licenses/>.
 */

/**
 * @file   techniques.hpp
 * @author Julian Iseringhausen <isering@ins.uni-bonn.de>
 * @date   Mon Apr 18 13:06:42 2011
 *
 * @brief  Some examples for multigrid methods this library
 *         might be used for.
 *
 */

#ifndef TECHNIQUES_HPP_
#define TECHNIQUES_HPP_

#include <sstream>
#include <string>

#include "base/command_list.hpp"

namespace VMG
{

class Techniques
{
public:
  static void SetCorrectionSchemeDirichlet(int minLevel, int maxLevel, int gamma)
    {
      CommandList* init = new CommandList();
      CommandList* loop = new CommandList();
      CommandList* finalize = new CommandList();

      init->AddCommand("ClearGrid", "RHS");
      init->AddCommand("ClearGrid", "SOL");
      init->AddCommand("ImportRightHandSide");
      init->AddCommand("ForceDiscreteCompatibility");
      init->AddCommand("CheckConsistency", "RHS");
      init->AddCommand("CopyBoundary", "RHS:SOL");
      init->AddCommand("InitializeIterationCounter");
      init->AddCommand("InitializeResidualNorm", "INITIAL_RESIDUAL");

      loop->AddCommand("ClearCoarseLevels", "RHS");
      loop->AddCommand("ClearCoarseLevels", "SOL");

      AddCycleGamma(*loop, maxLevel-minLevel+1, gamma);

      loop->AddCommand("ComputeResidualNorm", "RESIDUAL");
      loop->AddCommand("CheckResidual", "RESIDUAL");
      loop->AddCommand("CheckRelativeResidual", "RESIDUAL:INITIAL_RESIDUAL");
      loop->AddCommand("CheckIterationCounter");

      finalize->AddCommand("ExportSolution");

      init->Register("COMMANDLIST_INIT");
      loop->Register("COMMANDLIST_LOOP");
      finalize->Register("COMMANDLIST_FINALIZE");
    }


  static void SetCorrectionSchemeDirichletDebug(int minLevel, int maxLevel, int gamma)
  {
    CommandList* init = new CommandList();
    CommandList* loop = new CommandList();
    CommandList* finalize = new CommandList();

    init->AddCommand("ClearGrid", "RHS");
    init->AddCommand("ClearGrid", "SOL");
    init->AddCommand("ImportRightHandSide");
    init->AddCommand("ForceDiscreteCompatibility");
    init->AddCommand("CheckConsistency", "RHS");
    init->AddCommand("CopyBoundary", "RHS:SOL");
    init->AddCommand("InitializeIterationCounter");
    init->AddCommand("InitializeResidualNorm", "INITIAL_RESIDUAL");
    init->AddCommand("PrintAllSettings");

    loop->AddCommand("ClearCoarseLevels", "RHS");
    loop->AddCommand("ClearCoarseLevels", "SOL");

    AddCycleGammaDebug(*loop, maxLevel-minLevel+1, gamma);

    loop->AddCommand("ComputeResidualNorm", "RESIDUAL");
    loop->AddCommand("CheckResidual", "RESIDUAL");
    loop->AddCommand("CheckRelativeResidual", "RESIDUAL:INITIAL_RESIDUAL");
    loop->AddCommand("CheckIterationCounter");

    finalize->AddCommand("ExportSolution");

    init->Register("COMMANDLIST_INIT");
    loop->Register("COMMANDLIST_LOOP");
    finalize->Register("COMMANDLIST_FINALIZE");
  }

  static void SetCorrectionSchemePeriodic(int minLevel, int maxLevel, int gamma)
  {
    CommandList* init = new CommandList();
    CommandList* loop = new CommandList();
    CommandList* finalize = new CommandList();

    init->AddCommand("ClearGrid", "RHS");
    init->AddCommand("ClearGrid", "SOL");
    init->AddCommand("ImportRightHandSide");
    init->AddCommand("ForceDiscreteCompatibility");
    init->AddCommand("CheckConsistency", "RHS");
    init->AddCommand("InitializeIterationCounter");
    init->AddCommand("InitializeResidualNorm", "INITIAL_RESIDUAL");

    loop->AddCommand("ClearCoarseLevels", "SOL");
    loop->AddCommand("ClearCoarseLevels", "RHS");

    AddCycleGamma(*loop, maxLevel-minLevel+1, gamma);

    loop->AddCommand("ComputeResidualNorm", "RESIDUAL");
    loop->AddCommand("CheckResidual", "RESIDUAL");
    loop->AddCommand("CheckRelativeResidual", "RESIDUAL:INITIAL_RESIDUAL");
    loop->AddCommand("CheckIterationCounter");

    finalize->AddCommand("SetAverageToZero", "SOL");
    finalize->AddCommand("ExportSolution");

    init->Register("COMMANDLIST_INIT");
    loop->Register("COMMANDLIST_LOOP");
    finalize->Register("COMMANDLIST_FINALIZE");
  }

  static void SetCorrectionSchemePeriodicDebug(int minLevel, int maxLevel, int gamma)
  {
    CommandList* init = new CommandList();
    CommandList* loop = new CommandList();
    CommandList* finalize = new CommandList();

    init->AddCommand("ClearGrid", "RHS");
    init->AddCommand("ClearGrid", "SOL");
    init->AddCommand("ImportRightHandSide");
    init->AddCommand("ForceDiscreteCompatibility");
    init->AddCommand("CheckConsistency", "RHS");
    init->AddCommand("InitializeIterationCounter");
    init->AddCommand("InitializeResidualNorm", "INITIAL_RESIDUAL");
    init->AddCommand("PrintAllSettings");

    loop->AddCommand("ClearCoarseLevels", "SOL");
    loop->AddCommand("ClearCoarseLevels", "RHS");

    AddCycleGammaDebug(*loop, maxLevel-minLevel+1, gamma);

    loop->AddCommand("SetAverageToZero", "SOL");
    loop->AddCommand("ComputeResidualNorm", "RESIDUAL");
    loop->AddCommand("CheckResidual", "RESIDUAL");
    loop->AddCommand("CheckRelativeResidual", "RESIDUAL:INITIAL_RESIDUAL");
    loop->AddCommand("CheckIterationCounter");

    finalize->AddCommand("ExportSolution");

    init->Register("COMMANDLIST_INIT");
    loop->Register("COMMANDLIST_LOOP");
    finalize->Register("COMMANDLIST_FINALIZE");
  }

  static void SetCorrectionSchemePeriodicParticle(int minLevel, int maxLevel, int gamma)
  {
    CommandList* init = new CommandList();
    CommandList* loop = new CommandList();
    CommandList* finalize = new CommandList();

    init->AddCommand("ClearGrid", "RHS");
    init->AddCommand("ClearGrid", "SOL");
    init->AddCommand("ImportRightHandSide");
    init->AddCommand("CheckConsistency", "RHS");
    init->AddCommand("InitializeIterationCounter");
    init->AddCommand("InitializeResidualNorm", "INITIAL_RESIDUAL");

    loop->AddCommand("ClearCoarseLevels", "SOL");
    loop->AddCommand("ClearCoarseLevels", "RHS");

    AddCycleGamma(*loop, maxLevel-minLevel+1, gamma);

    loop->AddCommand("SetAverageToZero", "SOL");
    loop->AddCommand("ComputeResidualNorm", "RESIDUAL");
    loop->AddCommand("CheckResidual", "RESIDUAL");
    loop->AddCommand("CheckRelativeResidual", "RESIDUAL:INITIAL_RESIDUAL");
    loop->AddCommand("CheckIterationCounter");

    finalize->AddCommand("ExportSolution");

    init->Register("COMMANDLIST_INIT");
    loop->Register("COMMANDLIST_LOOP");
    finalize->Register("COMMANDLIST_FINALIZE");
  }

  static void SetCorrectionSchemePeriodicParticleDebug(int minLevel, int maxLevel, int gamma)
  {
    CommandList* init = new CommandList();
    CommandList* loop = new CommandList();
    CommandList* finalize = new CommandList();

    init->AddCommand("ClearGrid", "RHS");
    init->AddCommand("ClearGrid", "SOL");
    init->AddCommand("ImportRightHandSide");
    init->AddCommand("CheckConsistency", "RHS");
    init->AddCommand("InitializeIterationCounter");
    init->AddCommand("InitializeResidualNorm", "INITIAL_RESIDUAL");
    init->AddCommand("PrintAllSettings");

    loop->AddCommand("ClearCoarseLevels", "SOL");
    loop->AddCommand("ClearCoarseLevels", "RHS");

    AddCycleGammaDebug(*loop, maxLevel-minLevel+1, gamma);

    loop->AddCommand("SetAverageToZero", "SOL");
    loop->AddCommand("ComputeResidualNorm", "RESIDUAL");
    loop->AddCommand("CheckResidual", "RESIDUAL");
    loop->AddCommand("CheckRelativeResidual", "RESIDUAL:INITIAL_RESIDUAL");
    loop->AddCommand("CheckIterationCounter");

    finalize->AddCommand("ExportSolution");

    init->Register("COMMANDLIST_INIT");
    loop->Register("COMMANDLIST_LOOP");
    finalize->Register("COMMANDLIST_FINALIZE");
  }

  static void SetFullApproximationSchemeDirichlet(int minLevel, int maxLevel, int gamma)
  {
    CommandList* init = new CommandList();
    CommandList* loop = new CommandList();
    CommandList* finalize = new CommandList();

    init->AddCommand("ClearGrid", "RHS");
    init->AddCommand("ClearGrid", "SOL");
    init->AddCommand("ImportRightHandSide");
    init->AddCommand("ForceDiscreteCompatibility");
    init->AddCommand("CheckConsistency", "RHS");
    init->AddCommand("CopyBoundary", "RHS:SOL");
    init->AddCommand("InitializeIterationCounter");
    init->AddCommand("InitializeResidualNorm", "INITIAL_RESIDUAL");

    AddCycleGamma(*loop, maxLevel-minLevel+1, gamma);

    loop->AddCommand("ComputeResidualNorm", "RESIDUAL");
    loop->AddCommand("CheckResidual", "RESIDUAL");
    loop->AddCommand("CheckRelativeResidual", "RESIDUAL:INITIAL_RESIDUAL");
    loop->AddCommand("CheckIterationCounter");

    finalize->AddCommand("ExportSolution");

    init->Register("COMMANDLIST_INIT");
    loop->Register("COMMANDLIST_LOOP");
    finalize->Register("COMMANDLIST_FINALIZE");
  }

  static void SetFullApproximationSchemeDirichletDebug(int minLevel, int maxLevel, int gamma)
  {
    CommandList* init = new CommandList();
    CommandList* loop = new CommandList();
    CommandList* finalize = new CommandList();

    init->AddCommand("ClearGrid", "RHS");
    init->AddCommand("ClearGrid", "SOL");
    init->AddCommand("ImportRightHandSide");
    init->AddCommand("ForceDiscreteCompatibility");
    init->AddCommand("CheckConsistency", "RHS");
    init->AddCommand("CopyBoundary", "RHS:SOL");
    init->AddCommand("InitializeIterationCounter");
    init->AddCommand("InitializeResidualNorm", "INITIAL_RESIDUAL");
    init->AddCommand("PrintAllSettings");

    AddCycleGammaDebug(*loop, maxLevel-minLevel+1, gamma);

    loop->AddCommand("ComputeResidualNorm", "RESIDUAL");
    loop->AddCommand("CheckResidual", "RESIDUAL");
    loop->AddCommand("CheckRelativeResidual", "RESIDUAL:INITIAL_RESIDUAL");
    loop->AddCommand("CheckIterationCounter");

    finalize->AddCommand("ExportSolution");

    init->Register("COMMANDLIST_INIT");
    loop->Register("COMMANDLIST_LOOP");
    finalize->Register("COMMANDLIST_FINALIZE");
  }

  static void SetFullApproximationSchemePeriodic(int minLevel, int maxLevel, int gamma)
  {
    CommandList* init = new CommandList();
    CommandList* loop = new CommandList();
    CommandList* finalize = new CommandList();

    init->AddCommand("ClearGrid", "RHS");
    init->AddCommand("ClearGrid", "SOL");
    init->AddCommand("ImportRightHandSide");
    init->AddCommand("ForceDiscreteCompatibility");
    init->AddCommand("CheckConsistency", "RHS");
    init->AddCommand("InitializeIterationCounter");
    init->AddCommand("InitializeResidualNorm", "INITIAL_RESIDUAL");

    AddCycleGamma(*loop, maxLevel-minLevel+1, gamma);

    loop->AddCommand("ComputeResidualNorm", "RESIDUAL");
    loop->AddCommand("CheckResidual", "RESIDUAL");
    loop->AddCommand("CheckRelativeResidual", "RESIDUAL:INITIAL_RESIDUAL");
    loop->AddCommand("CheckIterationCounter");

    finalize->AddCommand("ExportSolution");

    init->Register("COMMANDLIST_INIT");
    loop->Register("COMMANDLIST_LOOP");
    finalize->Register("COMMANDLIST_FINALIZE");
  }

  static void SetFullApproximationSchemePeriodicDebug(int minLevel, int maxLevel, int gamma)
  {
    CommandList* init = new CommandList();
    CommandList* loop = new CommandList();
    CommandList* finalize = new CommandList();

    init->AddCommand("ClearGrid", "RHS");
    init->AddCommand("ClearGrid", "SOL");
    init->AddCommand("ImportRightHandSide");
    init->AddCommand("ForceDiscreteCompatibility");
    init->AddCommand("CheckConsistency", "RHS");
    init->AddCommand("InitializeIterationCounter");
    init->AddCommand("InitializeResidualNorm", "INITIAL_RESIDUAL");
    init->AddCommand("PrintAllSettings");

    AddCycleGammaDebug(*loop, maxLevel-minLevel+1, gamma);

    loop->AddCommand("ComputeResidualNorm", "RESIDUAL");
    loop->AddCommand("CheckResidual", "RESIDUAL");
    loop->AddCommand("CheckRelativeResidual", "RESIDUAL:INITIAL_RESIDUAL");
    loop->AddCommand("CheckIterationCounter");

    finalize->AddCommand("ExportSolution");

    init->Register("COMMANDLIST_INIT");
    loop->Register("COMMANDLIST_LOOP");
    finalize->Register("COMMANDLIST_FINALIZE");
  }

private:
  static void AddCycleGamma(CommandList& list, int numLevels, int gamma)
  {
    for (int i=0; i<numLevels-1; i++) {
      list.AddCommand("Smooth", "PRESMOOTHSTEPS");
      list.AddCommand("Restrict");
    }

    list.AddCommand("Solve");

    for (int i=1; i<gamma; i++) {

      for (int j=0; j<numLevels-2; j++) {

        for (int k=0; k<=j; k++) {
          list.AddCommand("Prolongate");
          list.AddCommand("Smooth", "PRESMOOTHSTEPS");
        }

        for (int k=0; k<=j; k++) {
          list.AddCommand("Smooth", "PRESMOOTHSTEPS");
          list.AddCommand("Restrict");
        }

        list.AddCommand("Solve");

      }

      for (int j=numLevels-4; j>=0; j--) {

        for (int k=0; k<=j; k++) {
          list.AddCommand("Prolongate");
          list.AddCommand("Smooth", "POSTSMOOTHSTEPS");
        }

        for (int k=0; k<=j; k++) {
          list.AddCommand("Smooth", "PRESMOOTHSTEPS");
          list.AddCommand("Restrict");
        }

        list.AddCommand("Solve");
      }

    }

    for (int i=0; i<numLevels-1; i++) {
      list.AddCommand("Prolongate");
      list.AddCommand("Smooth", "POSTSMOOTHSTEPS");
    }
  }

  static void AddCycleGammaDebug(CommandList& list, int numLevels, int gamma)
  {
    for (int i=0; i<numLevels-1; i++) {
      list.AddCommand("PrintGrid", "RHS");
      list.AddCommand("Smooth", "PRESMOOTHSTEPS");
      list.AddCommand("PrintGrid", "SOL");
      list.AddCommand("PrintDefect");
      list.AddCommand("Restrict");
    }

    list.AddCommand("PrintGrid", "RHS");
    list.AddCommand("Solve");
    list.AddCommand("PrintGrid", "SOL");
    list.AddCommand("PrintDefect");

    for (int i=1; i<gamma; i++) {

      for (int j=0; j<numLevels-2; j++) {

        for (int k=0; k<=j; k++) {
          list.AddCommand("Prolongate");
          list.AddCommand("PrintGrid", "RHS");
          list.AddCommand("Smooth", "PRESMOOTHSTEPS");
          list.AddCommand("PrintGrid", "SOL");
          list.AddCommand("PrintDefect");
        }

        for (int k=0; k<=j; k++) {
          list.AddCommand("PrintGrid", "RHS");
          list.AddCommand("Smooth", "PRESMOOTHSTEPS");
          list.AddCommand("PrintGrid", "SOL");
          list.AddCommand("PrintDefect");
          list.AddCommand("Restrict");
        }

        list.AddCommand("PrintGrid", "RHS");
        list.AddCommand("Solve");
        list.AddCommand("PrintGrid", "SOL");
        list.AddCommand("PrintDefect");

      }

      for (int j=numLevels-4; j>=0; j--) {

        for (int k=0; k<=j; k++) {
          list.AddCommand("Prolongate");
          list.AddCommand("PrintGrid", "RHS");
          list.AddCommand("Smooth", "POSTSMOOTHSTEPS");
          list.AddCommand("PrintGrid", "SOL");
          list.AddCommand("PrintDefect");
        }

        for (int k=0; k<=j; k++) {
          list.AddCommand("PrintGrid", "RHS");
          list.AddCommand("Smooth", "PRESMOOTHSTEPS");
          list.AddCommand("PrintGrid", "SOL");
          list.AddCommand("PrintDefect");
          list.AddCommand("Restrict");
        }

        list.AddCommand("PrintGrid", "RHS");
        list.AddCommand("Solve");
        list.AddCommand("PrintGrid", "SOL");
        list.AddCommand("PrintDefect");
      }

    }

    for (int i=0; i<numLevels-1; i++) {
      list.AddCommand("Prolongate");
      list.AddCommand("PrintGrid", "RHS");
      list.AddCommand("Smooth", "POSTSMOOTHSTEPS");
      list.AddCommand("PrintGrid", "SOL");
      list.AddCommand("PrintDefect");
    }
  }
};

}

#endif /* TECHNIQUES_HPP_ */