/* * Project: MoleCuilder * Description: creates and alters molecular systems * Copyright (C) 2013 University of Bonn. All rights reserved. * * * This file is part of MoleCuilder. * * MoleCuilder 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 2 of the License, or * (at your option) any later version. * * MoleCuilder 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 MoleCuilder. If not, see . */ /* * AnalyseFragmentationResultsAction.cpp * * Created on: Mar 8, 2013 * Author: heber */ // include config.h #ifdef HAVE_CONFIG_H #include #endif // include headers that implement a archive in simple text format // and before MemDebug due to placement new #include #include #include "CodePatterns/MemDebug.hpp" #include #include #include #include #include #include #include "CodePatterns/Assert.hpp" #include "CodePatterns/Log.hpp" #ifdef HAVE_JOBMARKET #include "JobMarket/types.hpp" #else typedef size_t JobId_t; #endif #include "Fragmentation/Summation/Containers/FragmentationChargeDensity.hpp" #include "Fragmentation/Summation/Containers/FragmentationLongRangeResults.hpp" #include "Fragmentation/Summation/Containers/FragmentationResultContainer.hpp" #include "Fragmentation/Summation/Containers/FragmentationShortRangeResults.hpp" #include "Fragmentation/Summation/Containers/MPQCData.hpp" #include "Fragmentation/Summation/Containers/MPQCData_printKeyNames.hpp" #include "Fragmentation/Homology/HomologyContainer.hpp" #include "Fragmentation/Homology/HomologyGraph.hpp" #include "Fragmentation/KeySetsContainer.hpp" #include "Fragmentation/Summation/SetValues/Fragment.hpp" #include "Fragmentation/Summation/SetValues/Histogram.hpp" #include "Fragmentation/Summation/SetValues/IndexedVectors.hpp" #include "Fragmentation/Summation/IndexSetContainer.hpp" #include "Fragmentation/Summation/writeTable.hpp" #ifdef HAVE_VMG #include "Fragmentation/Summation/Containers/VMGData.hpp" #include "Fragmentation/Summation/Containers/VMGDataFused.hpp" #include "Fragmentation/Summation/Containers/VMGDataMap.hpp" #include "Fragmentation/Summation/Containers/VMGData_printKeyNames.hpp" #endif #include "World.hpp" #include "Actions/FragmentationAction/AnalyseFragmentationResultsAction.hpp" using namespace MoleCuilder; // and construct the stuff #include "AnalyseFragmentationResultsAction.def" #include "Action_impl_pre.hpp" /** =========== define the function ====================== */ void writeToFile(const std::string &filename, const std::string contents) { std::ofstream tablefile(filename.c_str()); tablefile << contents; tablefile.close(); } /** Print (short range) energy, forces, and timings from received results. * * @param results summed up results container */ void printReceivedShortResults( const FragmentationShortRangeResults &results) { // print tables (without eigenvalues, they go extra) { typedef boost::mpl::remove::type MPQCDataEnergyVector_noeigenvalues_t; const std::string energyresult = writeTable()( results.Result_Energy_fused, results.getMaxLevel()); LOG(0, "Energy table is \n" << energyresult); std::string filename; filename += FRAGMENTPREFIX + std::string("_Energy.dat"); writeToFile(filename, energyresult); } { typedef boost::mpl::list< MPQCDataFused::energy_eigenvalues > MPQCDataEigenvalues_t; const std::string eigenvalueresult = writeTable()( results.Result_Energy_fused, results.getMaxLevel()); LOG(0, "Eigenvalue table is \n" << eigenvalueresult); std::string filename; filename += FRAGMENTPREFIX + std::string("_Eigenvalues.dat"); writeToFile(filename, eigenvalueresult); } { const std::string forceresult = writeTable()( results.Result_Force_fused, results.getMaxLevel()); LOG(0, "Force table is \n" << forceresult); std::string filename; filename += FRAGMENTPREFIX + std::string("_Forces.dat"); writeToFile(filename, forceresult); } // we don't want to print grid to a table { // print times (without flops for now) typedef boost::mpl::remove< boost::mpl::remove::type, MPQCDataFused::times_gather_flops>::type MPQCDataTimeVector_noflops_t; const std::string timesresult = writeTable()( results.Result_Time_fused, results.getMaxLevel()); LOG(0, "Times table is \n" << timesresult); std::string filename; filename += FRAGMENTPREFIX + std::string("_Times.dat"); writeToFile(filename, timesresult); } } /** Print long range energy from received results. * * @param results summed up results container */ void printReceivedFullResults( const FragmentationLongRangeResults &results) { // print tables (without eigenvalues, they go extra) if (results.Result_LongRange_fused.size() >= (results.getMaxLevel()-2)) { const std::string gridresult = writeTable()( results.Result_LongRange_fused, results.getMaxLevel(), 2); LOG(0, "VMG table is \n" << gridresult); std::string filename; filename += FRAGMENTPREFIX + std::string("_VMGEnergy.dat"); writeToFile(filename, gridresult); } if (results.Result_LongRangeIntegrated_fused.size() >= (results.getMaxLevel()-2)) { const std::string gridresult = writeTable()( results.Result_LongRangeIntegrated_fused, results.getMaxLevel(), 2); LOG(0, "LongRange table is \n" << gridresult); std::string filename; filename += FRAGMENTPREFIX + std::string("_LongRangeEnergy.dat"); writeToFile(filename, gridresult); } } bool appendToHomologyFile( const boost::filesystem::path &homology_file, const FragmentationShortRangeResults &shortrangeresults, const FragmentationLongRangeResults &longrangeresults ) { /// read homology container (if present) HomologyContainer homology_container; if (boost::filesystem::exists(homology_file)) { std::ifstream returnstream(homology_file.string().c_str()); if (returnstream.good()) { boost::archive::text_iarchive ia(returnstream); ia >> homology_container; } else { ELOG(2, "Failed to parse from " << homology_file.string() << "."); } returnstream.close(); } else { LOG(2, "Could not open " << homology_file.string() << ", creating empty container."); } /// append all fragments to a HomologyContainer HomologyContainer::container_t values; // convert KeySetContainer to IndexSetContainer IndexSetContainer::ptr ForceContainer(new IndexSetContainer(shortrangeresults.getForceKeySet())); const IndexSetContainer::Container_t &Indices = shortrangeresults.getContainer(); const IndexSetContainer::Container_t &ForceIndices = ForceContainer->getContainer(); IndexSetContainer::Container_t::const_iterator iter = Indices.begin(); IndexSetContainer::Container_t::const_iterator forceiter = ForceIndices.begin(); /// go through all fragments for (;iter != Indices.end(); ++iter, ++forceiter) // go through each IndexSet { /// create new graph entry in HomologyContainer which is (key,value) type LOG(1, "INFO: Creating new graph with " << **forceiter << "."); HomologyGraph graph(**forceiter); LOG(2, "DEBUG: Created graph " << graph << "."); const IndexSet::ptr &index = *iter; HomologyContainer::value_t value; // obtain fragment as key std::map >::const_iterator fragmentiter = longrangeresults.Result_perIndexSet_Fragment.find(index); ASSERT( fragmentiter != longrangeresults.Result_perIndexSet_Fragment.end(), "appendToHomologyFile() - cannot find index "+toString(*index) +" in FragmentResults."); value.first = boost::fusion::at_key(fragmentiter->second.first); // obtain energy as value std::map >::const_iterator energyiter = shortrangeresults.Result_perIndexSet_Energy.find(index); ASSERT( energyiter != shortrangeresults.Result_perIndexSet_Energy.end(), "appendToHomologyFile() - cannot find index "+toString(*index) +" in FragmentResults."); // value.second = boost::fusion::at_key(energyiter->second.first); // values value.second = boost::fusion::at_key(energyiter->second.second); // contributions values.insert( std::make_pair( graph, value) ); } homology_container.insert(values); LOG(1, "INFO: Listing all present atomic ids ..."); std::stringstream output; for (World::AtomIterator iter = World::getInstance().getAtomIter(); iter != World::getInstance().atomEnd(); ++iter) output << (*iter)->getId() << " "; LOG(1, "INFO: { " << output.str() << "} ."); // for debugging: print container LOG(1, "INFO: Listing all present homologies ..."); for (HomologyContainer::container_t::const_iterator iter = homology_container.begin(); iter != homology_container.end(); ++iter) { LOG(1, "INFO: graph " << iter->first << " has Fragment " << iter->second.first << " and associated energy " << iter->second.second << "."); } /// store homology container again std::ofstream outputstream(homology_file.string().c_str()); if (outputstream.good()) { // check if opened boost::archive::text_oarchive oa(outputstream); oa << homology_container; if (outputstream.fail()) { // check if correctly written LOG(1, "Failed to write to file " << homology_file.string() << "."); return false; } else outputstream.close(); } else { LOG(1, "Failed to open file " << homology_file.string() << " for writing."); return false; } return true; } Action::state_ptr FragmentationAnalyseFragmentationResultsAction::performCall() { // if file is given, parse from file into resultscontainer FragmentationResultContainer& container = FragmentationResultContainer::getInstance(); if (!params.resultsfile.get().empty()) { boost::filesystem::path resultsfile = params.resultsfile.get(); if (boost::filesystem::exists(resultsfile)) { LOG(1, "INFO: Parsing results from " << resultsfile.string() << "."); std::ifstream returnstream(resultsfile.string().c_str()); if (returnstream.good()) { boost::archive::text_iarchive ia(returnstream); ia >> container; } } else { ELOG(1, "Given file" << resultsfile.string() << " does not exist."); } } // get data and keysets from ResultsContainer const std::map &shortrangedata = container.getShortRangeResults(); const KeySetsContainer &keysets = container.getKeySets(); const KeySetsContainer &forcekeysets = container.getForceKeySets(); const bool DoLongrange = container.areFullRangeResultsPresent(); if (keysets.KeySets.empty()) { ELOG(2, "There are no results in the container."); return Action::failure; } FragmentationShortRangeResults shortrangeresults(shortrangedata, keysets, forcekeysets); shortrangeresults(shortrangedata); printReceivedShortResults(shortrangeresults); #ifdef HAVE_VMG if (DoLongrange) { if ( World::getInstance().getAllAtoms().size() == 0) { ELOG(1, "Please load the full molecule into the world before starting this action."); return Action::failure; } FragmentationChargeDensity summedChargeDensity(shortrangedata,keysets); const std::vector full_sample = summedChargeDensity.getFullSampledGrid(); std::map longrangeData = container.getLongRangeResults(); // remove full solution corresponding to full_sample from map (must be highest ids), has to be treated extra std::map::iterator iter = longrangeData.end(); for (size_t i=0;i::iterator remove_iter = iter; std::vector fullsolutionData; for (; iter != longrangeData.end(); ++iter) fullsolutionData.push_back(iter->second); longrangeData.erase(remove_iter, longrangeData.end()); // Final phase: sum up and print result FragmentationLongRangeResults longrangeresults( shortrangedata,longrangeData,keysets, forcekeysets); longrangeresults( shortrangedata, longrangeData, fullsolutionData, full_sample); printReceivedFullResults(longrangeresults); // append all keysets to homology file if (!params.homology_file.get().empty()) { const boost::filesystem::path &homology_file = params.homology_file.get(); if (homology_file.string() != "") { LOG(1, "INFO: Appending HomologyGraphs to file " << homology_file.string() << "."); if (!appendToHomologyFile(homology_file, shortrangeresults, longrangeresults)) return Action::failure; } } } #endif container.clear(); return Action::success; } Action::state_ptr FragmentationAnalyseFragmentationResultsAction::performUndo(Action::state_ptr _state) { return Action::success; } Action::state_ptr FragmentationAnalyseFragmentationResultsAction::performRedo(Action::state_ptr _state){ return Action::success; } bool FragmentationAnalyseFragmentationResultsAction::canUndo() { return false; } bool FragmentationAnalyseFragmentationResultsAction::shouldUndo() { return false; } /** =========== end of function ====================== */