source: src/mg.cpp@ 4571da

/**
 * @file   mg.cpp
 * @author Julian Iseringhausen <isering@ins.uni-bonn.de>
 * @date   Sat Jun 12 20:36:24 2010
 *
 * @brief  A multigrid solver
 *
 * This file contains the implementation of the main class for
 * a multigrid solver.
 *
 */

#ifdef HAVE_CONFIG_H
#include <config.h>
#endif

#ifdef DEBUG_MEASURE_TIME
#ifdef HAVE_MPI
#include <mpi.h>
#ifdef HAVE_MARMOT
#include <enhancempicalls.h>
#include <sourceinfompicalls.h>
#endif
#endif
#endif

#include <cmath>
#include <cstdio>
#include <cfloat>
#include <ctime>
#include <string>
#include <fstream>
#include <sstream>

#include "base/command_list.hpp"
#include "base/discretization.hpp"
#include "base/factory.hpp"
#include "base/timer.hpp"
#include "comm/comm.hpp"
#include "grid/grid.hpp"
#include "grid/tempgrid.hpp"
#include "interface/interface.hpp"
#include "interface/interface_particles_cf.hpp"
#include "level/level_operator.hpp"
#include "smoother/smoother.hpp"
#include "solver/solver.hpp"
#include "mg.hpp"

using namespace VMG;

#define REGISTER_COMMAND(a) extern void Initialize##a();Initialize##a();

VMG::CommandFactory MG::command_factory;

static void VMGRegisterBuiltinCommands()
{
  REGISTER_COMMAND(VMGCommandCheckConsistency);
  REGISTER_COMMAND(VMGCommandCheckIterationCounter);
  REGISTER_COMMAND(VMGCommandCheckRelativeResidual);
  REGISTER_COMMAND(VMGCommandCheckResidual);
  REGISTER_COMMAND(VMGCommandClearCoarseLevels);
  REGISTER_COMMAND(VMGCommandClearGrid);
  REGISTER_COMMAND(VMGCommandComputeResidualNorm);
  REGISTER_COMMAND(VMGCommandCopyBoundary);
  REGISTER_COMMAND(VMGCommandExecuteCycle);
  REGISTER_COMMAND(VMGCommandExecuteCycleLoop);
  REGISTER_COMMAND(VMGCommandExecuteFullCycle);
  REGISTER_COMMAND(VMGCommandExecuteFullCycleLoop);
  REGISTER_COMMAND(VMGCommandExportSolution);
  REGISTER_COMMAND(VMGCommandForceDiscreteCompatibility);
  REGISTER_COMMAND(VMGCommandImportRightHandSide);
  REGISTER_COMMAND(VMGCommandInterpolateFMG);
  REGISTER_COMMAND(VMGCommandInitializeIterationCounter);
  REGISTER_COMMAND(VMGCommandInitializeResidualNorm);
  REGISTER_COMMAND(VMGCommandNOP);
  REGISTER_COMMAND(VMGCommandPrintAllSettings);
  REGISTER_COMMAND(VMGCommandPrintDefect);
  REGISTER_COMMAND(VMGCommandPrintGridStructure);
  REGISTER_COMMAND(VMGCommandPrintGrid);
  REGISTER_COMMAND(VMGCommandPrintResidualNorm);
  REGISTER_COMMAND(VMGCommandPrintRunningTime);
  REGISTER_COMMAND(VMGCommandProlongate);
  REGISTER_COMMAND(VMGCommandRestrict);
  REGISTER_COMMAND(VMGCommandSetAverageToZero);
  REGISTER_COMMAND(VMGCommandSetCoarserDirichletValues);
  REGISTER_COMMAND(VMGCommandSetLevel);
  REGISTER_COMMAND(VMGCommandSmooth);
  REGISTER_COMMAND(VMGCommandSolve);
}

MG::MG()
{
  state = 0;
  VMGRegisterBuiltinCommands();
}

MG::~MG()
{
  MG::Destroy();
}

// Brings Multigrid back to starting state.
void MG::Destroy()
{
  MG::Instance()->factories.clear();
  MG::Instance()->state = 0;
  Timer::Clear();
}

void MG::SetInterface(VMG::Interface* interface, VMG::Comm* comm)
{
  interface->Register("INTERFACE");

  Multigrid* sol = new Multigrid(comm, interface);
  Multigrid* rhs = new Multigrid(comm, interface);
  TempGrid* temp = new TempGrid();

  sol->Register("SOL");
  rhs->Register("RHS");
  temp->Register("TEMPGRID");

  new ObjectStorage<int>("GLOBAL_MAXLEVEL", sol->GlobalMaxLevel());
  new ObjectStorage<int>("MINLEVEL", sol->MinLevel());
  new ObjectStorage<int>("MAXLEVEL", sol->MaxLevel());
}

void MG::InitializeFromParticleConfigFile(const char* filename)
{
  InterfaceParticlesCF interface_particles(filename);
}

/*
 * Post init communication class
 */
void MG::PostInit()
{
  GetComm()->PostInit(*GetSol(), *GetRhs());
}

/**
 * Solves a given system with a multigrid method
 *
 */
void MG::Solve()
{
#ifdef DEBUG_MEASURE_TIME
#ifdef HAVE_MPI
  MPI_Barrier(MPI_COMM_WORLD);
#endif
  Timer::Start("CompleteRunningTime");
#endif

  CommandList* cl_init = MG::GetFactory().Get("COMMANDLIST_INIT")->Cast<CommandList>();
  CommandList* cl_loop = MG::GetFactory().Get("COMMANDLIST_LOOP")->Cast<CommandList>();
  CommandList* cl_finalize = MG::GetFactory().Get("COMMANDLIST_FINALIZE")->Cast<CommandList>();

  cl_init->ExecuteList();

  while (cl_loop->ExecuteList() == Continue);

  cl_finalize->ExecuteList();

#ifdef DEBUG_MEASURE_TIME
#ifdef HAVE_MPI
  MPI_Barrier(MPI_COMM_WORLD);
#endif
  Timer::Stop("CompleteRunningTime");
#ifdef DEBUG_MEASURE_TIME_OUTPUT
#ifdef HAVE_MPI
  Timer::PrintGlobal();
#else
  Timer::Print();
#endif
#endif
#endif
}

void MG::SetState(const int& state_)
{
  MG::Instance()->state = state_;
}

VMG::Factory& MG::GetFactory()
{
  std::map<int, VMG::Factory>::iterator iter = MG::Instance()->factories.find(MG::Instance()->state);

  if (iter == MG::Instance()->factories.end())
    iter = MG::Instance()->factories.insert(std::make_pair(MG::Instance()->state, Factory())).first;

  assert(iter != MG::Instance()->factories.end());

  return iter->second;
}

VMG::CommandFactory& MG::GetCommands()
{
  return MG::command_factory;
}

Comm* MG::GetComm()
{
  return MG::GetFactory().Get("COMM")->Cast<VMG::Comm>();
}

Discretization* MG::GetDiscretization()
{
  return MG::GetFactory().Get("DISCRETIZATION")->Cast<VMG::Discretization>();
}

LevelOperator* MG::GetLevelOperator()
{
  return MG::GetFactory().Get("LEVEL_OPERATOR")->Cast<VMG::LevelOperator>();
}

Multigrid* MG::GetRhs()
{
  return MG::GetFactory().Get("RHS")->Cast<VMG::Multigrid>();
}

Multigrid* MG::GetSol()
{
  return MG::GetFactory().Get("SOL")->Cast<VMG::Multigrid>();
}

VMG::Grid& MG::GetRhsMaxLevel()
{
  return (*MG::GetRhs())(MG::GetRhs()->MaxLevel());
}

VMG::Grid& MG::GetSolMaxLevel()
{
  return (*MG::GetSol())(MG::GetSol()->MaxLevel());
}

Smoother* MG::GetSmoother()
{
  return MG::GetFactory().Get("SMOOTHER")->Cast<VMG::Smoother>();
}

Solver* MG::GetSolver()
{
  return MG::GetFactory().Get("SOLVER")->Cast<VMG::Solver>();
}

TempGrid* MG::GetTempGrid()
{
  return MG::GetFactory().Get("TEMPGRID")->Cast<VMG::TempGrid>();
}

Interface* MG::GetInterface()
{
  return MG::GetFactory().Get("INTERFACE")->Cast<VMG::Interface>();
}

static bool CheckObject(std::string id)
{
  Object *obj = MG::GetFactory().Get(id);

#ifdef DEBUG_OUTPUT
  if (obj == NULL)
    printf("\nMultigrid: CLASS %s NOT INITIALIZED\n\n", id.c_str());
#endif

  return obj != NULL;
}

bool MG::IsInitialized()
{
  bool init = true;

  init &= CheckObject("COMM");
  init &= CheckObject("LEVEL_OPERATOR");
  init &= CheckObject("RHS");
  init &= CheckObject("SOL");
  init &= CheckObject("SOLVER");
  init &= CheckObject("SMOOTHER");
  init &= CheckObject("DISCRETIZATION");
  init &= CheckObject("MAX_ITERATION");
  init &= CheckObject("PRECISION");
  init &= CheckObject("PRESMOOTHSTEPS");
  init &= CheckObject("POSTSMOOTHSTEPS");
  init &= CheckObject("COMMANDLIST_INIT");
  init &= CheckObject("COMMANDLIST_LOOP");
  init &= CheckObject("COMMANDLIST_FINALIZE");
  init &= CheckObject("MINLEVEL");
  init &= CheckObject("MAXLEVEL");
  init &= CheckObject("GLOBAL_MAXLEVEL");
  init &= CheckObject("INTERFACE");

  return init;
}