procharv.h

Go to the documentation of this file.

/////////////////////////////////////////////////////////////////////////////
///
/// @file procharv.h
///
/// A plug-in for post-processing of Morse decompositions
/// for the harvesting model.
/// This is a replacement for "procmdec.h".
///
/// @author Pawel Pilarczyk
///
/////////////////////////////////////////////////////////////////////////////
 
// Copyright (C) 1997-2014 by Pawel Pilarczyk.
//
// This file is part of my research software package. This 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.
//
// This software 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 software; see the file "license.txt". If not,
// please, see <https://www.gnu.org/licenses/>.
 
// Started on June 27, 2011. Last revision: June 28, 2011.
 
 
#ifndef _CMGRAPHS_PROCHARV_H_
#define _CMGRAPHS_PROCHARV_H_
 
// include some standard C++ header files
#include <sstream>
#include <fstream>
#include <string>
 
// include selected header files from the CHomP library
#include "chomp/system/textfile.h"
#include "chomp/system/timeused.h"
#include "chomp/struct/digraph.h"
#include "chomp/struct/flatmatr.h"
#include "chomp/struct/multitab.h"
#include "chomp/cubes/neighbor.h"
 
// include local header files
#include "config.h"
#include "typedefs.h"
#include "conindex.h"
#include "morsedec.h"
#include "typedyns.h"
#include "utils.h"
 
 
namespace custom {
 
/// Post-processing of Morse decompositions for the harvesting model.
namespace procMorseDecHarv {
 
 
// --------------------------------------------------
// ------- post-process a Morse decomposition -------
// --------------------------------------------------
 
/// Post-processes a Morse decomposition by computing a decomposition
/// of each Morse set into cycle sets and analyzing this decomposition.
template <class CubSetType, class MorseDecType, class GraphType,
        class CubMapType>
inline void processMorseDec (const MorseDecType &morseDec,
        const CubSetType &allCubes, const GraphType &g,
        const CubMapType &theCubMap, const CubMapType &,
        const std::string &cubesFilePrefix,
        const std::string &procFilePrefix)
{
        using namespace chomp::homology;
 
        sbug << "===== Post-processing the Morse sets. =====\n";
        timeused timePostProc;
        timePostProc = 0;
 
        // prepare the size of attractor population size and attractor index
        double attrPopSize = 0;
        double attrPopSizes [spaceDim];
        for (int i = 0; i < spaceDim; ++ i)
                attrPopSizes [i] = 0;
        int attractor = -1;
 
        // analyze each Morse set
        int nSets = morseDec. count ();
        for (int n = 0; n < nSets; ++ n)
        {
                const CubSetType &morseCubes = morseDec [n];
                int_t morseCount = morseCubes. size ();
 
                // restrict the analysis to the given Morse set
                sbug << "Morse set no. " << n << ": " <<
                        morseCount << " cubes, ";
 
                // compute the real exit set for the Morse set
                sbug << "exit set: ";
                CubSetType exitSet;
                for (int i = 0; i < morseCount; ++ i)
                        theCubMap (morseCubes [i], &exitSet, 0, 0, 0, false);
                exitSet. remove (morseCubes);
                sbug << exitSet. size () << " cubes.\n";
 
                // save the exit set to a file if requested to
                if (!cubesFilePrefix. empty ())
                {
                        // prepare a file name for this set
                        std::ostringstream exitFileNameStr;
                        exitFileNameStr << cubesFilePrefix << n << "e.cub";
                        std::string exitFileName = exitFileNameStr. str ();
 
                        // save the exit set to a file
                        sbug << "Saving the exit set to '" <<
                                exitFileName << "'... ";
                        std::ofstream f (exitFileName. c_str ());
                        f << "; The exit set for the corresponding Morse "
                                "set (" << exitSet. size () << " cubes).\n";
                        f << exitSet;
                        f. close ();
                        sbug << "Done.\n";
                }
 
                // compute the average population size
                double popSize = 0;
                double popSizes [spaceDim];
                double boxSize [spaceDim];
                for (int i = 0; i < spaceDim; ++ i)
                {
                        popSizes [i] = 0;
                        boxSize [i] = spaceWidth [i] / (1 << finalDepth);
                }
                for (int k = 0; k < morseCount; ++ k)
                {
                        spcCoord c [spaceDim];
                        morseCubes [k]. coord (c);
                        double sum = 0;
                        for (int i = 0; i < spaceDim; ++ i)
                        {
                                double size = spaceOffset [i] +
                                        (c [i] + 0.5) * boxSize [i];
                                sum += size;
                                popSizes [i] += size;
                        }
                        popSize += sum;
                }
                popSize /= morseCount;
                sbug << "Avg pop sizes: ";
                for (int i = 0; i < spaceDim; ++ i)
                {
                        popSizes [i] /= morseCount;
                        sbug << (i ? " + " : "") << popSizes [i];
                }
                sbug << " = " << popSize << ".\n";
 
                // determine if this Morse set is at the bottom of the CMG
                bool bottom = true;
                for (int m = 0; m < nSets; ++ m)
                {
                        if (morseDec. connected (n, m))
                        {
                                bottom = false;
                                break;
                        }
                }
 
                // store the index of the attractor and the population size
                if ((attractor < 0) && bottom)
                {
                        attrPopSize = popSize;
                        for (int i = 0; i < spaceDim; ++ i)
                                attrPopSizes [i] = popSizes [i];
                        attractor = exitSet. empty () ? n : (-1 - n);
                }
        }
 
        // output the index of the attractor and the attr. population size
        int attrIndex = (attractor < 0) ? (-1 - attractor) : attractor;
        sbug << ((attractor < 0) ? "Non-isol a" : "A") << "ttractor " <<
                attrIndex << " (" << morseDec [attrIndex]. size () <<
                " cubes), avg pop sizes: ";
        for (int i = 0; i < spaceDim; ++ i)
                sbug << (i ? " + " : "") << attrPopSizes [i];
        sbug << " = " << attrPopSize << ". Coord ranges: ";
        CoordMinMax coordMinMax;
        coordMinMax (morseDec [attrIndex]);
        sbug << coordMinMax << ".\n";
 
        sbug << "===== Post-processing completed in " << timePostProc <<
                ". =====\n";
 
        return;
} /* processMorseDec */
 
 
} // namespace procMorseDecHarv
} // namespace custom
 
 
#endif // _CMGRAPHS_PROCHARV_H_