diff --git a/src/AnalysisList.h b/src/AnalysisList.h index fd9470058c..2decf17c52 100644 --- a/src/AnalysisList.h +++ b/src/AnalysisList.h @@ -4,6 +4,9 @@ /// Hold all Analyses to be performed. class AnalysisList { public: + /// Analysis setup status + enum AnalysisStatusType { NO_SETUP = 0, SETUP, INACTIVE }; + AnalysisList(); ~AnalysisList(); void Clear(); @@ -13,9 +16,13 @@ class AnalysisList { int DoAnalyses(); void List() const; bool Empty() const { return analysisList_.empty(); } + unsigned int size() const { return analysisList_.size(); } +# ifdef MPI + Analysis& Ana(int i) { return *(analysisList_[i].ptr_); } + AnalysisStatusType Status(int i) const { return analysisList_[i].status_; } + ArgList const& Args(int i) const { return analysisList_[i].args_; } +# endif private: - /// Analysis setup status - enum AnalysisStatusType { NO_SETUP = 0, SETUP, INACTIVE }; struct AnaHolder { Analysis* ptr_; ///< Pointer to Analysis ArgList args_; ///< Arguments associated with Analysis diff --git a/src/Analysis_HausdorffDistance.cpp b/src/Analysis_HausdorffDistance.cpp new file mode 100644 index 0000000000..905ea50511 --- /dev/null +++ b/src/Analysis_HausdorffDistance.cpp @@ -0,0 +1,229 @@ +#include // std::min, max +#include "Analysis_HausdorffDistance.h" +#include "CpptrajStdio.h" +#include "DataSet_MatrixFlt.h" + +Analysis_HausdorffDistance::Analysis_HausdorffDistance() : + outType_(BASIC), + out_(0), + ab_out_(0), + ba_out_(0) +{} + +/** Assume input matrix contains all distances beween sets A and B. + * The directed Hausdorff distance from A to B will be the maximum of all + * minimums in each row. + * The directed Hausdorff distance from B to A will be the maximum of all + * minimums in each column. + * The symmetric Hausdorff distance is the max of the two directed distances. + * \param m1 The matrix containing distances from A to B. + * \param hd_ab The directed Hausdorff distance from A to B. + * \param hd_ba The directed Hausdorff distance from B to A. + * \return the symmetric Hausdorff distance. + */ +double Analysis_HausdorffDistance::CalcHausdorffFromMatrix(DataSet_2D const& m1, + double& hd_ab, double& hd_ba) +{ +// if (m1 == 0) { +// mprinterr("Internal Error: Analysis_HausdorffDistance::(): Matrix set is null.\n"); +// return -1.0; +// } + if (m1.Size() < 1) { + mprinterr("Error: '%s' is empty.\n", m1.legend()); + return -1.0; + } + // Hausdorff distance from A to B. + hd_ab = 0.0; + for (unsigned int row = 0; row != m1.Nrows(); row++) + { + double minRow = m1.GetElement(0, row); + for (unsigned int col = 1; col != m1.Ncols(); col++) + minRow = std::min( minRow, m1.GetElement(col, row) ); + //mprintf("DEBUG: Min row %6u is %12.4f\n", row, minRow); + hd_ab = std::max( hd_ab, minRow ); + } + //mprintf("DEBUG: Hausdorff A to B= %12.4f\n", hd_ab); + // Hausdorff distance from B to A. + hd_ba = 0.0; + for (unsigned int col = 0; col != m1.Ncols(); col++) + { + double minCol = m1.GetElement(col, 0); + for (unsigned int row = 1; row != m1.Nrows(); row++) + minCol = std::min( minCol, m1.GetElement(col, row) ); + //mprintf("DEBUG: Min col %6u is %12.4f\n", col, minCol); + hd_ba = std::max( hd_ba, minCol); + } + //mprintf("DEBUG: Hausdorff B to A= %12.4f\n", hd_ba); + // Symmetric Hausdorff distance + double hd = std::max( hd_ab, hd_ba ); + + return hd; +} + +// Analysis_HausdorffDistance::Help() +void Analysis_HausdorffDistance::Help() const { + mprintf("\t [ ...]\n" + "\t[outtype {basic|trimatrix nrows <#>|fullmatrix nrows <#> [ncols <#>]}]\n" + "\t[name ] [out ] [outab ] [outba ]\n" + " Given 1 or more 2D matrices containing distances between two sets\n" + " A and B, calculate the symmetric Hausdorff distance for each matrix.\n" + " The results can be saved as an array or as an upper-triangular\n" + " matrix with the specified number of rows.\n"); +} + + +// Analysis_HausdorffDistance::Setup() +Analysis::RetType Analysis_HausdorffDistance::Setup(ArgList& analyzeArgs, AnalysisSetup& setup, int debugIn) +{ + // Keywords + int nrows = -1; + int ncols = -1; + std::string outtypearg = analyzeArgs.GetStringKey("outtype"); + if (!outtypearg.empty()) { + if (outtypearg == "basic") + outType_ = BASIC; + else if (outtypearg == "trimatrix") { + outType_ = UPPER_TRI_MATRIX; + nrows = analyzeArgs.getKeyInt("nrows", -1); + if (nrows < 1) { + mprinterr("Error: 'nrows' must be specified and > 0 for 'trimatrix'\n"); + return Analysis::ERR; + } + } else if (outtypearg == "fullmatrix") { + outType_ = FULL_MATRIX; + nrows = analyzeArgs.getKeyInt("nrows", -1); + if (nrows < 1) { + mprinterr("Error: 'nrows' must be specified and > 0 for 'fullmatrix'\n"); + return Analysis::ERR; + } + ncols = analyzeArgs.getKeyInt("ncols", nrows); + if (ncols < 1) { + mprinterr("Error: 'ncols' must be > 0 for 'fullmatrix'\n"); + return Analysis::ERR; + } + } else { + mprinterr("Error: Unrecognized keyword for 'outtype': %s\n", outtypearg.c_str()); + return Analysis::ERR; + } + } else + outType_ = BASIC; + std::string dsname = analyzeArgs.GetStringKey("name"); + DataFile* df = setup.DFL().AddDataFile( analyzeArgs.GetStringKey("out"), analyzeArgs ); + DataFile* dfab = setup.DFL().AddDataFile( analyzeArgs.GetStringKey("outab"), analyzeArgs ); + DataFile* dfba = setup.DFL().AddDataFile( analyzeArgs.GetStringKey("outba"), analyzeArgs ); + // Get input data sets + std::string dsarg = analyzeArgs.GetStringNext(); + while (!dsarg.empty()) { + DataSetList selected = setup.DSL().GetMultipleSets( dsarg ); + for (DataSetList::const_iterator set = selected.begin(); set != selected.end(); ++set) + { + if ((*set)->Group() == DataSet::MATRIX_2D) + inputSets_.AddCopyOfSet( *set ); + else + mprintf("Warning: Currently only 2D matrices supported; skipping set '%s'\n", + (*set)->legend()); + } + //inputSets_ += setup.DSL().GetMultipleSets( dsarg ); + dsarg = analyzeArgs.GetStringNext(); + } + if (inputSets_.empty()) { + mprinterr("Error: No data sets specified.\n"); + return Analysis::ERR; + } + // Output data set + out_ = 0; + if (outType_ == BASIC) { + out_ = setup.DSL().AddSet(DataSet::FLOAT, dsname, "HAUSDORFF"); + if (out_ == 0) return Analysis::ERR; + // Directed sets + ab_out_ = setup.DSL().AddSet(DataSet::FLOAT, MetaData(out_->Meta().Name(),"AB")); + if (ab_out_ == 0) return Analysis::ERR; + ba_out_ = setup.DSL().AddSet(DataSet::FLOAT, MetaData(out_->Meta().Name(),"BA")); + if (ba_out_ == 0) return Analysis::ERR; + } else if (outType_ == UPPER_TRI_MATRIX || outType_ == FULL_MATRIX) { + out_ = setup.DSL().AddSet(DataSet::MATRIX_FLT, dsname, "HAUSDORFF"); + ab_out_ = setup.DSL().AddSet(DataSet::MATRIX_FLT, MetaData(out_->Meta().Name(),"AB")); + ba_out_ = setup.DSL().AddSet(DataSet::MATRIX_FLT, MetaData(out_->Meta().Name(),"BA")); + if (out_ == 0 || ab_out_ == 0 || ba_out_ == 0) return Analysis::ERR; + if (outType_ == UPPER_TRI_MATRIX) { + if (((DataSet_2D*)out_)->AllocateTriangle( nrows )) return Analysis::ERR; + if (((DataSet_2D*)ab_out_)->AllocateTriangle( nrows )) return Analysis::ERR; + if (((DataSet_2D*)ba_out_)->AllocateTriangle( nrows )) return Analysis::ERR; + } else if (outType_ == FULL_MATRIX) { + if (((DataSet_2D*)out_)->Allocate2D( nrows,ncols )) return Analysis::ERR; + if (((DataSet_2D*)ab_out_)->Allocate2D( nrows,ncols )) return Analysis::ERR; + if (((DataSet_2D*)ba_out_)->Allocate2D( nrows,ncols )) return Analysis::ERR; + } + if (out_->Size() != inputSets_.size()) { + mprinterr("Warning: Number of input data sets (%zu) != number of expected" + " sets in matrix (%zu)\n", inputSets_.size(), out_->Size()); + return Analysis::ERR; + } + // Directed sets + + } + if (df != 0) + df->AddDataSet( out_ ); + if (dfab != 0) + df->AddDataSet( ab_out_ ); + if (dfba != 0) + df->AddDataSet( ba_out_ ); + + mprintf(" HAUSDORFF:\n"); + mprintf("\tCalculating Hausdorff distances from the following 2D distance matrices:\n\t "); + for (DataSetList::const_iterator it = inputSets_.begin(); it != inputSets_.end(); ++it) + mprintf(" %s", (*it)->legend()); + mprintf("\n"); + if (outType_ == BASIC) + mprintf("\tOutput will be stored in 1D array set '%s'\n", out_->legend()); + else if (outType_ == UPPER_TRI_MATRIX) + mprintf("\tOutput will be stored in upper-triangular matrix set '%s' with %i rows.\n", + out_->legend(), nrows); + else if (outType_ == FULL_MATRIX) + mprintf("\tOutput will be stored in matrix set '%s' with %i rows, %i columns.\n", + out_->legend(), nrows, ncols); + mprintf("\tDirected A->B distance output set: %s\n", ab_out_->legend()); + mprintf("\tDirected B->A distance output set: %s\n", ba_out_->legend()); + if (df != 0) mprintf("\tOutput set written to '%s'\n", df->DataFilename().full()); + if (dfab != 0) mprintf("\tA->B output set written to '%s'\n", dfab->DataFilename().full()); + if (dfba != 0) mprintf("\tB->A output set written to '%s'\n", dfba->DataFilename().full()); + + return Analysis::OK; +} + +// Analysis_HausdorffDistance::Analyze() +Analysis::RetType Analysis_HausdorffDistance::Analyze() { + // Get the Hausdorff distance for each set. + double hd = 0.0; + double hd_ab = 0.0; + double hd_ba = 0.0; + int idx = 0; + for (DataSetList::const_iterator it = inputSets_.begin(); it != inputSets_.end(); ++it) + { + hd = -1.0; + hd_ab = -1.0; + hd_ba = -1.0; + if ( (*it)->Group() == DataSet::MATRIX_2D ) + hd = CalcHausdorffFromMatrix( static_cast( *(*it) ), hd_ab, hd_ba ); + else + mprintf("Warning: '%s' type not yet supported for Hausdorff\n", (*it)->legend()); + mprintf("%12.4f %s\n", hd, (*it)->legend()); + float fhd = (float)hd; + float fab = (float)hd_ab; + float fba = (float)hd_ba; + switch (outType_) { + case BASIC : + out_->Add(idx, &fhd); + ab_out_->Add(idx, &fab); + ba_out_->Add(idx++, &fba); + break; + case UPPER_TRI_MATRIX : + case FULL_MATRIX : + ((DataSet_MatrixFlt*)out_)->AddElement( fhd ); + ((DataSet_MatrixFlt*)ab_out_)->AddElement( fhd ); + ((DataSet_MatrixFlt*)ba_out_)->AddElement( fhd ); + break; + } + } + return Analysis::OK; +} diff --git a/src/Analysis_HausdorffDistance.h b/src/Analysis_HausdorffDistance.h new file mode 100644 index 0000000000..15770643d6 --- /dev/null +++ b/src/Analysis_HausdorffDistance.h @@ -0,0 +1,28 @@ +#ifndef INC_ANALYSIS_HAUSDORFFDISTANCE_H +#define INC_ANALYSIS_HAUSDORFFDISTANCE_H +#include "Analysis.h" +#include "DataSet_2D.h" +/// Compute the Hausdorff distance between 2 or more data sets. +/** The Hausdorff distance between two sets of points A and B is defined as: + * h(A,B) = max[a in A]{ min[b in B]( d(a,b) } } + */ +class Analysis_HausdorffDistance : public Analysis { + public: + Analysis_HausdorffDistance(); + DispatchObject* Alloc() const { return (DispatchObject*)new Analysis_HausdorffDistance(); } + void Help() const; + + Analysis::RetType Setup(ArgList&, AnalysisSetup&, int); + Analysis::RetType Analyze(); + private: + static double CalcHausdorffFromMatrix(DataSet_2D const&, double&, double&); + + enum OutType { BASIC = 0, UPPER_TRI_MATRIX, FULL_MATRIX }; + + DataSetList inputSets_; + OutType outType_; + DataSet* out_; + DataSet* ab_out_; + DataSet* ba_out_; +}; +#endif diff --git a/src/Command.cpp b/src/Command.cpp index f71b5693aa..7ecf86bfe4 100644 --- a/src/Command.cpp +++ b/src/Command.cpp @@ -19,6 +19,7 @@ #include "Exec_DataSetCmd.h" #include "Exec_GenerateAmberRst.h" #include "Exec_Help.h" +#include "Exec_ParallelAnalysis.h" #include "Exec_Precision.h" #include "Exec_PrintData.h" #include "Exec_ReadData.h" @@ -174,6 +175,7 @@ #include "Analysis_Multicurve.h" #include "Analysis_TI.h" #include "Analysis_ConstantPHStats.h" +#include "Analysis_HausdorffDistance.h" CmdList Command::commands_ = CmdList(); @@ -206,6 +208,7 @@ void Command::Init() { Command::AddCmd( new Exec_ListAll(), Cmd::EXE, 1, "list" ); Command::AddCmd( new Exec_NoExitOnError(), Cmd::EXE, 1, "noexitonerror" ); Command::AddCmd( new Exec_NoProgress(), Cmd::EXE, 1, "noprogress" ); + Command::AddCmd( new Exec_ParallelAnalysis(),Cmd::EXE, 1, "parallelanalysis" ); Command::AddCmd( new Exec_Precision(), Cmd::EXE, 1, "precision" ); Command::AddCmd( new Exec_PrintData(), Cmd::EXE, 1, "printdata" ); Command::AddCmd( new Exec_QuietBlocks(), Cmd::EXE, 1, "quietblocks" ); @@ -364,6 +367,7 @@ void Command::Init() { Command::AddCmd( new Analysis_Matrix(), Cmd::ANA, 2, "diagmatrix", "matrix" ); Command::AddCmd( new Analysis_Divergence(), Cmd::ANA, 1, "divergence" ); Command::AddCmd( new Analysis_FFT(), Cmd::ANA, 1, "fft" ); + Command::AddCmd( new Analysis_HausdorffDistance,Cmd::ANA,1,"hausdorff" ); Command::AddCmd( new Analysis_Hist(), Cmd::ANA, 2, "hist", "histogram" ); Command::AddCmd( new Analysis_Integrate(), Cmd::ANA, 1, "integrate" ); Command::AddCmd( new Analysis_IRED(), Cmd::ANA, 1, "ired" ); diff --git a/src/CpptrajState.h b/src/CpptrajState.h index 7d274d2540..13371c1744 100644 --- a/src/CpptrajState.h +++ b/src/CpptrajState.h @@ -28,6 +28,8 @@ class CpptrajState { DataSetList& DSL() { return DSL_; } DataFileList const& DFL() const { return DFL_; } DataFileList& DFL() { return DFL_; } + AnalysisList const& Analyses() const { return analysisList_; } + AnalysisList& Analyses() { return analysisList_; } TrajModeType Mode() const { return mode_; } int Debug() const { return debug_; } bool ShowProgress() const { return showProgress_;} diff --git a/src/DataFileList.cpp b/src/DataFileList.cpp index 44c62e9764..325235bcf5 100644 --- a/src/DataFileList.cpp +++ b/src/DataFileList.cpp @@ -339,6 +339,34 @@ void DataFileList::WriteAllDF() { # endif } +#ifdef MPI +// DataFileList::WriteAllDF() +/** Call write for all DataFiles in list for which writeFile is true. Once + * a file has been written set writeFile to false; it can be reset to + * true if new DataSets are added to it. All threads are allowed to try + * writing files. FIXME this is a hack until a better way of determining + * write threads is found. + */ +void DataFileList::AllThreads_WriteAllDF() { + if (fileList_.empty()) return; +# ifdef TIMER + Timer datafile_time; + datafile_time.Start(); +# endif + for (DFarray::iterator df = fileList_.begin(); df != fileList_.end(); ++df) { + if ( (*df)->DFLwrite() ) { + (*df)->SetThreadCanWrite( true ); + (*df)->WriteDataOut(); + (*df)->SetDFLwrite( false ); + } + } +# ifdef TIMER + datafile_time.Stop(); + mprintf("TIME: Write of all data files took %.4f seconds.\n", datafile_time.Total()); +# endif +} +#endif + // DataFileList::UnwrittenData() bool DataFileList::UnwrittenData() const { for (DFarray::const_iterator df = fileList_.begin(); df != fileList_.end(); ++df) diff --git a/src/DataFileList.h b/src/DataFileList.h index c49ff60127..84c2ce1e01 100644 --- a/src/DataFileList.h +++ b/src/DataFileList.h @@ -44,6 +44,7 @@ class DataFileList { CpptrajFile* AddCpptrajFile(FileName const&,std::string const&,CFtype,bool); # ifdef MPI CpptrajFile* AddCpptrajFile(FileName const&,std::string const&,CFtype,bool,Parallel::Comm const&); + void AllThreads_WriteAllDF(); # endif /// List DataFiles and CpptrajFiles. void List() const; diff --git a/src/DataIO_Gnuplot.cpp b/src/DataIO_Gnuplot.cpp index f8054cb917..6826ed5182 100644 --- a/src/DataIO_Gnuplot.cpp +++ b/src/DataIO_Gnuplot.cpp @@ -220,6 +220,7 @@ void DataIO_Gnuplot::WriteHelp() { "\tpm3d: Normal pm3d map output.\n" "\tnopm3d: Turn off pm3d\n" "\tjpeg: Plot will write to a JPEG file when used with gnuplot.\n" + "\ttitle: Plot title. Default is file name.\n" // "\tbinary: Use binary output\n" "\tnoheader: Do not format plot; data output only.\n" "\tpalette : Change gnuplot pm3d palette to :\n" @@ -242,6 +243,7 @@ int DataIO_Gnuplot::processWriteArgs(ArgList &argIn) { if (argIn.hasKey("jpeg")) jpegout_ = true; if (argIn.hasKey("binary")) binary_ = true; if (argIn.hasKey("noheader")) writeHeader_ = false; + title_ = argIn.GetStringKey("title"); if (!writeHeader_ && jpegout_) { mprintf("Warning: jpeg output not supported with 'noheader' option.\n"); jpegout_ = false; @@ -309,7 +311,13 @@ void DataIO_Gnuplot::WriteRangeAndHeader(Dimension const& Xdim, size_t Xmax, file_.Printf("set yrange [%8.3f:%8.3f]\nset xrange [%8.3f:%8.3f]\n", Ydim.Coord(0) - Ydim.Step(), Ydim.Coord(Ymax + 1), Xdim.Coord(0) - Xdim.Step(), Xdim.Coord(Xmax + 1)); - file_.Printf("splot \"%s\"%s%s title \"%s\"\n", data_fname_.full(), binaryFlag, pm3dstr.c_str(), file_.Filename().base()); + const char* tout; + if (!title_.empty()) + tout = title_.c_str(); + else + tout = file_.Filename().base(); + file_.Printf("splot \"%s\"%s%s title \"%s\"\n", data_fname_.full(), binaryFlag, pm3dstr.c_str(), + tout); } // DataIO_Gnuplot::Finish() diff --git a/src/DataIO_Gnuplot.h b/src/DataIO_Gnuplot.h index 2b21df8745..dc966b7152 100644 --- a/src/DataIO_Gnuplot.h +++ b/src/DataIO_Gnuplot.h @@ -16,6 +16,7 @@ class DataIO_Gnuplot : public DataIO { private: CpptrajFile file_; FileName data_fname_; ///< Data file name + std::string title_; ///< Plot title (output) typedef std::vector LabelArray; LabelArray Xlabels_; LabelArray Ylabels_; diff --git a/src/DataSet.h b/src/DataSet.h index fc3c971e26..4e80b7d914 100644 --- a/src/DataSet.h +++ b/src/DataSet.h @@ -67,6 +67,9 @@ class DataSet { # ifdef MPI /// Piece this DataSet together from multiple threads. virtual int Sync(size_t, std::vector const&, Parallel::Comm const&) = 0; + // TODO pure virtual + virtual int SendSet(int, Parallel::Comm const&) { return 1; } + virtual int RecvSet(int, Parallel::Comm const&) { return 1; } # endif // ----------------------------------------------------- /// Associate additional data with this set. diff --git a/src/DataSet_MatrixFlt.cpp b/src/DataSet_MatrixFlt.cpp index 5a6dcfd8a2..9eabeeac10 100644 --- a/src/DataSet_MatrixFlt.cpp +++ b/src/DataSet_MatrixFlt.cpp @@ -26,4 +26,32 @@ int DataSet_MatrixFlt::Sync(size_t total, std::vector const& rank_frames, commIn.ReduceMaster( 0, &(mat_[0]), mat_.size(), MPI_FLOAT, MPI_SUM ); return 0; } + +int DataSet_MatrixFlt::SendSet(int dest, Parallel::Comm const& commIn) { + // First send the size + type + int size[3]; + size[0] = mat_.Ncols(); + size[1] = mat_.Nrows(); + size[2] = (int)mat_.Type(); + commIn.Send( &size, 3, MPI_INT, dest, 1700 ); + // Then send the matrix + commIn.Send( &(mat_[0]), mat_.size(), MPI_FLOAT, dest, 1701 ); + return 0; +} + +int DataSet_MatrixFlt::RecvSet(int src, Parallel::Comm const& commIn) { + // First receive the size + type + int size[3]; + commIn.Recv( &size, 3, MPI_INT, src, 1700 ); + if (size[2] == 0) // FULL + mat_.resize(size[0], size[1]); + else if (size[2] == 1) // HALF + mat_.resize(size[0], 0); + else // TRI + mat_.resize(0, size[1]); + // Then receive the matrix + commIn.Recv( &(mat_[0]), mat_.size(), MPI_FLOAT, src, 1701); + return 0; +} + #endif diff --git a/src/DataSet_MatrixFlt.h b/src/DataSet_MatrixFlt.h index 9a9362eec0..90ab289fce 100644 --- a/src/DataSet_MatrixFlt.h +++ b/src/DataSet_MatrixFlt.h @@ -13,6 +13,8 @@ class DataSet_MatrixFlt : public DataSet_2D { # ifdef MPI // FIXME: Currently just sums up. Should this be a separate Sync function? int Sync(size_t, std::vector const&, Parallel::Comm const&); + int SendSet(int, Parallel::Comm const&); + int RecvSet(int, Parallel::Comm const&); # endif void Info() const { return; } void WriteBuffer(CpptrajFile&, SizeArray const&) const; diff --git a/src/Exec_ParallelAnalysis.cpp b/src/Exec_ParallelAnalysis.cpp new file mode 100644 index 0000000000..436e068024 --- /dev/null +++ b/src/Exec_ParallelAnalysis.cpp @@ -0,0 +1,127 @@ +#include +#include "Exec_ParallelAnalysis.h" +#include "CpptrajStdio.h" +#ifdef MPI +#include "Timer.h" + +// Exec_ParallelAnalysis::Help() +void Exec_ParallelAnalysis::Help() const +{ + mprintf("\t[sync]\n" + " Divide all currently set up Analyses among available MPI processes and run.\n" + " If 'sync' is specifed, sync all results back to master process; this may be\n" + " required if performing subsequent analyses.\n"); +} + +// Exec_ParallelAnalysis::Execute() +Exec::RetType Exec_ParallelAnalysis::Execute(CpptrajState& State, ArgList& argIn) +{ + Timer t_total; + t_total.Start(); + Timer t_sync; + bool syncToMaster = argIn.hasKey("sync"); + std::vector setSizesBefore; + if (syncToMaster) { + t_sync.Start(); + setSizesBefore.reserve( State.DSL().size() ); + for (DataSetList::const_iterator it = State.DSL().begin(); it != State.DSL().end(); ++it) + setSizesBefore.push_back( (*it)->Size() ); + t_sync.Stop(); + } + // DEBUG - Have each thread report what analyses it knows about and what + // data sets it has. +/* + for (int rank = 0; rank < Parallel::World().Size(); rank++) { + if (rank == Parallel::World().Rank()) { + printf("Rank %i, %u analyses, %zu data sets:\n", rank, State.Analyses().size(), + State.DSL().size()); + for (DataSetList::const_iterator it = State.DSL().begin(); it != State.DSL().end(); ++it) + printf("\t'%s' (%zu)\n", (*it)->Meta().PrintName().c_str(), (*it)->Size()); + } + Parallel::World().Barrier(); + } + Parallel::World().Barrier(); +*/ + mprintf(" PARALLELANALYSIS: Will attempt to run current Analyses in parallel.\n\n" + "*** THIS COMMAND IS STILL EXPERIMENTAL! ***\n\n"); + if (syncToMaster) + mprintf("\tResulting data sets will be synced back to master.\n"); + // Naively divide up all analyses among threads. + int my_start, my_stop; + int nelts = Parallel::World().DivideAmongProcesses( my_start, my_stop, State.Analyses().size() ); + rprintf("Dividing %zu analyses among %i threads: %i to %i (%i)\n", + State.Analyses().size(), Parallel::World().Size(), my_start, my_stop, nelts); + + // Each thread runs the analyses they are responsible for. + int nerr = 0; + for (int na = my_start; na != my_stop; na++) { + // TODO check setup status + if (State.Analyses().Ana(na).Analyze() != Analysis::OK) { + rprinterr("Error: Analysis failed: '%s'\n", State.Analyses().Args(na).ArgLine()); + nerr++; + } + } + // This error check serves as a barrier + if (Parallel::World().CheckError( nerr )) return CpptrajState::ERR; + State.DFL().AllThreads_WriteAllDF(); + State.Analyses().Clear(); + if (syncToMaster) { + t_sync.Start(); + // Check which sizes have changed. + if (setSizesBefore.size() != State.DSL().size()) { + mprintf("Warning: Number of sets have changed. Not attempting to sync sets to master.\n"); + } else { + for (unsigned int idx = 0; idx != State.DSL().size(); idx++) { + int setHasChanged = 0; + if (!Parallel::World().Master()) { + if (setSizesBefore[idx] != State.DSL()[idx]->Size()) { + if (State.Debug() > 0) + rprintf("Set '%s' size has changed from %u to %zu\n", + State.DSL()[idx]->legend(), setSizesBefore[idx], State.DSL()[idx]->Size()); + setHasChanged = 1; + } + } + int totalChanged; + Parallel::World().AllReduce(&totalChanged, &setHasChanged, 1, MPI_INT, MPI_SUM); + if (totalChanged > 0) { + if (totalChanged == 1) { + int sourceRank = 0; + if (setHasChanged == 1) + setHasChanged = Parallel::World().Rank(); + Parallel::World().ReduceMaster(&sourceRank, &setHasChanged, 1, MPI_INT, MPI_SUM); + if (State.Debug() > 0) + mprintf("DEBUG: Need to sync '%s' from %i\n", State.DSL()[idx]->legend(), sourceRank); + if (Parallel::World().Master()) + State.DSL()[idx]->RecvSet( sourceRank, Parallel::World() ); + else if (setHasChanged == Parallel::World().Rank()) + State.DSL()[idx]->SendSet( 0, Parallel::World() ); + } else + mprintf("Warning: '%s' exists on multiple threads. Not syncing.\n", + State.DSL()[idx]->legend()); + } + } + } + t_sync.Stop(); + } + t_total.Stop(); + if (syncToMaster) + t_sync.WriteTiming(2, "Sync:", t_total.Total()); + t_total.WriteTiming(1, "Total:"); + return CpptrajState::OK; +} +#else /* MPI */ +// ============================================================================= +// NON-MPI CODE BELOW +// Exec_ParallelAnalysis::Help() +void Exec_ParallelAnalysis::Help() const +{ + mprintf(" This command is only available in MPI-enabled builds.\n"); +} + +// Exec_ParallelAnalysis::Execute() +Exec::RetType Exec_ParallelAnalysis::Execute(CpptrajState& State, ArgList& argIn) +{ + mprinterr("Error: This command is only available in MPI-enabled builds.\n"); + return CpptrajState::ERR; +} +#endif /* MPI */ diff --git a/src/Exec_ParallelAnalysis.h b/src/Exec_ParallelAnalysis.h new file mode 100644 index 0000000000..94f05af319 --- /dev/null +++ b/src/Exec_ParallelAnalysis.h @@ -0,0 +1,12 @@ +#ifndef INC_EXEC_PARALLELANALYSIS_H +#define INC_EXEC_PARALLELANALYSIS_H +#include "Exec.h" +/// +class Exec_ParallelAnalysis : public Exec { + public: + Exec_ParallelAnalysis() : Exec(GENERAL) {} + void Help() const; + DispatchObject* Alloc() const { return (DispatchObject*)new Exec_ParallelAnalysis(); } + RetType Execute(CpptrajState&, ArgList&); +}; +#endif diff --git a/src/Matrix.h b/src/Matrix.h index 2443e88402..3245381253 100644 --- a/src/Matrix.h +++ b/src/Matrix.h @@ -24,7 +24,7 @@ template class Matrix { /// \return total number of elements in the matrix. size_t size() const { return nelements_; } /// \return current matrix type. - //MType Type() const { return type_; } + MType Type() const { return type_; } /// \return estimated size in bytes. static size_t sizeInBytes(size_t,size_t); /// Set up matrix for given number of cols and rows. diff --git a/src/Parallel.h b/src/Parallel.h index 39157673af..81cf0624f1 100644 --- a/src/Parallel.h +++ b/src/Parallel.h @@ -34,6 +34,9 @@ * 1400+X: Action_AtomicCorr: Atomic movement vectors * 1500 : Action_NAstruct: Array containing BP step info on rank. * 1501+X: Array of step series data from rank. + * 1600+X: Ensemble sort + * 1700 : DataSet_MatrixFlt size + * 1701 : DataSet_MatrixFlt buffer */ class Parallel { public: diff --git a/src/Version.h b/src/Version.h index 8e03c23d5f..c711599bcd 100644 --- a/src/Version.h +++ b/src/Version.h @@ -20,5 +20,5 @@ * Whenever a number that precedes is incremented, all subsequent * numbers should be reset to 0. */ -#define CPPTRAJ_INTERNAL_VERSION "V4.11.0" +#define CPPTRAJ_INTERNAL_VERSION "V4.11.1" #endif diff --git a/src/cpptrajdepend b/src/cpptrajdepend index da9e8ee14b..178553e614 100644 --- a/src/cpptrajdepend +++ b/src/cpptrajdepend @@ -97,6 +97,7 @@ Analysis_CrossCorr.o : Analysis_CrossCorr.cpp ActionState.h Analysis.h AnalysisS Analysis_CurveFit.o : Analysis_CurveFit.cpp ActionState.h Analysis.h AnalysisState.h Analysis_CurveFit.h ArgList.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h CoordinateInfo.h CpptrajFile.h CpptrajStdio.h CurveFit.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_1D.h DataSet_Coords.h DataSet_Coords_REF.h DataSet_Mesh.h Dimension.h DispatchObject.h FileIO.h FileName.h FileTypes.h Frame.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h Parallel.h ParameterTypes.h RPNcalc.h Range.h ReferenceFrame.h ReplicaDimArray.h Residue.h Spline.h StringRoutines.h SymbolExporting.h TextFormat.h Timer.h Topology.h Vec3.h Analysis_Divergence.o : Analysis_Divergence.cpp ActionState.h Analysis.h AnalysisState.h Analysis_Divergence.h ArgList.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h CoordinateInfo.h CpptrajFile.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_1D.h DataSet_Coords.h DataSet_Coords_REF.h Dimension.h DispatchObject.h FileIO.h FileName.h FileTypes.h Frame.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h Residue.h SymbolExporting.h TextFormat.h Timer.h Topology.h Vec3.h Analysis_FFT.o : Analysis_FFT.cpp ActionState.h Analysis.h AnalysisState.h Analysis_FFT.h ArgList.h Array1D.h ArrayIterator.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h ComplexArray.h CoordinateInfo.h CpptrajFile.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_1D.h DataSet_Coords.h DataSet_Coords_REF.h Dimension.h DispatchObject.h FileIO.h FileName.h FileTypes.h Frame.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h Parallel.h ParameterTypes.h PubFFT.h Range.h ReferenceFrame.h ReplicaDimArray.h Residue.h SymbolExporting.h TextFormat.h Timer.h Topology.h Vec3.h +Analysis_HausdorffDistance.o : Analysis_HausdorffDistance.cpp ActionState.h Analysis.h AnalysisState.h Analysis_HausdorffDistance.h ArgList.h ArrayIterator.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h CoordinateInfo.h CpptrajFile.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_2D.h DataSet_Coords.h DataSet_Coords_REF.h DataSet_MatrixFlt.h Dimension.h DispatchObject.h FileIO.h FileName.h FileTypes.h Frame.h MaskToken.h Matrix.h Matrix_3x3.h MetaData.h Molecule.h NameType.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h Residue.h SymbolExporting.h TextFormat.h Timer.h Topology.h Vec3.h Analysis_Hist.o : Analysis_Hist.cpp ActionFrameCounter.h ActionState.h Analysis.h AnalysisState.h Analysis_Hist.h ArgList.h ArrayIterator.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h BondSearch.h Box.h CharMask.h Constants.h CoordinateInfo.h CpptrajFile.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_1D.h DataSet_2D.h DataSet_3D.h DataSet_Coords.h DataSet_Coords_REF.h DataSet_GridFlt.h DataSet_MatrixDbl.h DataSet_double.h Dimension.h DispatchObject.h FileIO.h FileName.h FileTypes.h Frame.h FrameArray.h FramePtrArray.h Grid.h GridBin.h HistBin.h MaskToken.h Matrix.h Matrix_3x3.h MetaData.h Molecule.h NameType.h OutputTrajCommon.h Parallel.h ParameterTypes.h ParmFile.h ParmIO.h Range.h ReferenceFrame.h ReplicaDimArray.h Residue.h StringRoutines.h SymbolExporting.h TextFormat.h Timer.h Topology.h TrajectoryFile.h TrajectoryIO.h Trajout_Single.h Vec3.h Analysis_IRED.o : Analysis_IRED.cpp ActionState.h Analysis.h AnalysisState.h Analysis_IRED.h ArgList.h ArrayIterator.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h ComplexArray.h Constants.h CoordinateInfo.h Corr.h CpptrajFile.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_1D.h DataSet_2D.h DataSet_Coords.h DataSet_Coords_REF.h DataSet_MatrixDbl.h DataSet_Modes.h DataSet_Vector.h DataSet_double.h Dimension.h DispatchObject.h FileIO.h FileName.h FileTypes.h Frame.h MaskToken.h Matrix.h Matrix_3x3.h MetaData.h Molecule.h NameType.h Parallel.h ParameterTypes.h PubFFT.h Range.h ReferenceFrame.h ReplicaDimArray.h Residue.h SymbolExporting.h TextFormat.h Timer.h Topology.h Vec3.h Analysis_Integrate.o : Analysis_Integrate.cpp ActionState.h Analysis.h AnalysisState.h Analysis_Integrate.h ArgList.h Array1D.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h CoordinateInfo.h CpptrajFile.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_1D.h DataSet_Coords.h DataSet_Coords_REF.h DataSet_Mesh.h Dimension.h DispatchObject.h FileIO.h FileName.h FileTypes.h Frame.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h Residue.h Spline.h SymbolExporting.h TextFormat.h Timer.h Topology.h Vec3.h @@ -153,7 +154,7 @@ Cluster_ReadInfo.o : Cluster_ReadInfo.cpp ArgList.h ArrayIterator.h AssociatedDa Cmd.o : Cmd.cpp Cmd.h DispatchObject.h CmdInput.o : CmdInput.cpp CmdInput.h StringRoutines.h CmdList.o : CmdList.cpp Cmd.h CmdList.h DispatchObject.h -Command.o : Command.cpp Action.h ActionFrameCounter.h ActionList.h ActionState.h Action_Align.h Action_Angle.h Action_AreaPerMol.h Action_AtomMap.h Action_AtomicCorr.h Action_AtomicFluct.h Action_AutoImage.h Action_Average.h Action_Bounds.h Action_Box.h Action_Center.h Action_Channel.h Action_CheckChirality.h Action_CheckStructure.h Action_Closest.h Action_ClusterDihedral.h Action_Contacts.h Action_CreateCrd.h Action_CreateReservoir.h Action_DNAionTracker.h Action_DSSP.h Action_Density.h Action_Diffusion.h Action_Dihedral.h Action_Dipole.h Action_DistRmsd.h Action_Distance.h Action_Energy.h Action_Esander.h Action_FilterByData.h Action_FixAtomOrder.h Action_FixImagedBonds.h Action_GIST.h Action_Grid.h Action_GridFreeEnergy.h Action_HydrogenBond.h Action_Image.h Action_InfraredSpectrum.h Action_Jcoupling.h Action_LESsplit.h Action_LIE.h Action_LipidOrder.h Action_MakeStructure.h Action_Mask.h Action_Matrix.h Action_MinImage.h Action_Molsurf.h Action_MultiDihedral.h Action_MultiVector.h Action_NAstruct.h Action_NMRrst.h Action_NativeContacts.h Action_OrderParameter.h Action_Outtraj.h Action_PairDist.h Action_Pairwise.h Action_Principal.h Action_Projection.h Action_Pucker.h Action_Radgyr.h Action_Radial.h Action_RandomizeIons.h Action_Remap.h Action_ReplicateCell.h Action_Rmsd.h Action_Rotate.h Action_RunningAvg.h Action_STFC_Diffusion.h Action_Scale.h Action_SetVelocity.h Action_Spam.h Action_Strip.h Action_Surf.h Action_SymmetricRmsd.h Action_Temperature.h Action_Translate.h Action_Unstrip.h Action_Unwrap.h Action_Vector.h Action_VelocityAutoCorr.h Action_Volmap.h Action_Volume.h Action_Watershell.h Analysis.h AnalysisList.h AnalysisState.h Analysis_AmdBias.h Analysis_AutoCorr.h Analysis_Average.h Analysis_Clustering.h Analysis_ConstantPHStats.h Analysis_Corr.h Analysis_CrankShaft.h Analysis_CrdFluct.h Analysis_CrossCorr.h Analysis_CurveFit.h Analysis_Divergence.h Analysis_FFT.h Analysis_Hist.h Analysis_IRED.h Analysis_Integrate.h Analysis_KDE.h Analysis_Lifetime.h Analysis_LowestCurve.h Analysis_Matrix.h Analysis_MeltCurve.h Analysis_Modes.h Analysis_MultiHist.h Analysis_Multicurve.h Analysis_Overlap.h Analysis_PhiPsi.h Analysis_Regression.h Analysis_RemLog.h Analysis_Rms2d.h Analysis_RmsAvgCorr.h Analysis_Rotdif.h Analysis_RunningAvg.h Analysis_Spline.h Analysis_State.h Analysis_Statistics.h Analysis_TI.h Analysis_Timecorr.h Analysis_VectorMath.h Analysis_Wavelet.h ArgList.h Array1D.h ArrayIterator.h AssociatedData.h Atom.h AtomExtra.h AtomMap.h AtomMask.h AxisType.h BaseIOtype.h Box.h BufferedLine.h CharMask.h ClusterDist.h ClusterList.h ClusterMap.h ClusterNode.h ClusterSieve.h Cmd.h CmdInput.h CmdList.h Command.h ComplexArray.h Constraints.h Control.h CoordinateInfo.h Corr.h Cph.h CpptrajFile.h CpptrajState.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_1D.h DataSet_2D.h DataSet_3D.h DataSet_Cmatrix.h DataSet_Coords.h DataSet_Coords_CRD.h DataSet_Coords_REF.h DataSet_GridFlt.h DataSet_Mat3x3.h DataSet_MatrixDbl.h DataSet_MatrixFlt.h DataSet_Mesh.h DataSet_Modes.h DataSet_RemLog.h DataSet_Vector.h DataSet_double.h DataSet_float.h DataSet_integer.h DataSet_pH.h DataSet_string.h Deprecated.h DihedralSearch.h Dimension.h DispatchObject.h DistRoutines.h Energy.h Energy_Sander.h EnsembleIn.h EnsembleOut.h EnsembleOutList.h Ewald.h Exec.h Exec_Analyze.h Exec_Calc.h Exec_Change.h Exec_ClusterMap.h Exec_CombineCoords.h Exec_Commands.h Exec_CompareTop.h Exec_CrdAction.h Exec_CrdOut.h Exec_CreateSet.h Exec_DataFile.h Exec_DataFilter.h Exec_DataSetCmd.h Exec_GenerateAmberRst.h Exec_Help.h Exec_LoadCrd.h Exec_LoadTraj.h Exec_ParmBox.h Exec_ParmSolvent.h Exec_ParmStrip.h Exec_ParmWrite.h Exec_PermuteDihedrals.h Exec_Precision.h Exec_PrintData.h Exec_ReadData.h Exec_ReadEnsembleData.h Exec_ReadInput.h Exec_RotateDihedral.h Exec_RunAnalysis.h Exec_ScaleDihedralK.h Exec_SequenceAlign.h Exec_SortEnsembleData.h Exec_SplitCoords.h Exec_System.h Exec_Top.h Exec_Traj.h Exec_UpdateParameters.h Exec_ViewRst.h FileIO.h FileName.h FileTypes.h Frame.h FrameArray.h FramePtrArray.h Grid.h GridAction.h GridBin.h HistBin.h Hungarian.h ImageTypes.h ImagedAction.h InputTrajCommon.h MapAtom.h MaskToken.h Matrix.h Matrix_3x3.h MetaData.h Molecule.h NameType.h NetcdfFile.h OnlineVarT.h OutputTrajCommon.h PDBfile.h PairList.h Parallel.h ParameterHolders.h ParameterTypes.h PubFFT.h RPNcalc.h Random.h Range.h ReferenceAction.h ReferenceFrame.h RemdReservoirNC.h ReplicaDimArray.h ReplicaInfo.h Residue.h Spline.h StructureCheck.h SymbolExporting.h SymmetricRmsdCalc.h TextFormat.h Timer.h Topology.h TrajFrameCounter.h TrajectoryFile.h TrajectoryIO.h Trajin.h TrajinList.h TrajoutList.h Trajout_Single.h VariableArray.h Vec3.h molsurf.h +Command.o : Command.cpp Action.h ActionFrameCounter.h ActionList.h ActionState.h Action_Align.h Action_Angle.h Action_AreaPerMol.h Action_AtomMap.h Action_AtomicCorr.h Action_AtomicFluct.h Action_AutoImage.h Action_Average.h Action_Bounds.h Action_Box.h Action_Center.h Action_Channel.h Action_CheckChirality.h Action_CheckStructure.h Action_Closest.h Action_ClusterDihedral.h Action_Contacts.h Action_CreateCrd.h Action_CreateReservoir.h Action_DNAionTracker.h Action_DSSP.h Action_Density.h Action_Diffusion.h Action_Dihedral.h Action_Dipole.h Action_DistRmsd.h Action_Distance.h Action_Energy.h Action_Esander.h Action_FilterByData.h Action_FixAtomOrder.h Action_FixImagedBonds.h Action_GIST.h Action_Grid.h Action_GridFreeEnergy.h Action_HydrogenBond.h Action_Image.h Action_InfraredSpectrum.h Action_Jcoupling.h Action_LESsplit.h Action_LIE.h Action_LipidOrder.h Action_MakeStructure.h Action_Mask.h Action_Matrix.h Action_MinImage.h Action_Molsurf.h Action_MultiDihedral.h Action_MultiVector.h Action_NAstruct.h Action_NMRrst.h Action_NativeContacts.h Action_OrderParameter.h Action_Outtraj.h Action_PairDist.h Action_Pairwise.h Action_Principal.h Action_Projection.h Action_Pucker.h Action_Radgyr.h Action_Radial.h Action_RandomizeIons.h Action_Remap.h Action_ReplicateCell.h Action_Rmsd.h Action_Rotate.h Action_RunningAvg.h Action_STFC_Diffusion.h Action_Scale.h Action_SetVelocity.h Action_Spam.h Action_Strip.h Action_Surf.h Action_SymmetricRmsd.h Action_Temperature.h Action_Translate.h Action_Unstrip.h Action_Unwrap.h Action_Vector.h Action_VelocityAutoCorr.h Action_Volmap.h Action_Volume.h Action_Watershell.h Analysis.h AnalysisList.h AnalysisState.h Analysis_AmdBias.h Analysis_AutoCorr.h Analysis_Average.h Analysis_Clustering.h Analysis_ConstantPHStats.h Analysis_Corr.h Analysis_CrankShaft.h Analysis_CrdFluct.h Analysis_CrossCorr.h Analysis_CurveFit.h Analysis_Divergence.h Analysis_FFT.h Analysis_HausdorffDistance.h Analysis_Hist.h Analysis_IRED.h Analysis_Integrate.h Analysis_KDE.h Analysis_Lifetime.h Analysis_LowestCurve.h Analysis_Matrix.h Analysis_MeltCurve.h Analysis_Modes.h Analysis_MultiHist.h Analysis_Multicurve.h Analysis_Overlap.h Analysis_PhiPsi.h Analysis_Regression.h Analysis_RemLog.h Analysis_Rms2d.h Analysis_RmsAvgCorr.h Analysis_Rotdif.h Analysis_RunningAvg.h Analysis_Spline.h Analysis_State.h Analysis_Statistics.h Analysis_TI.h Analysis_Timecorr.h Analysis_VectorMath.h Analysis_Wavelet.h ArgList.h Array1D.h ArrayIterator.h AssociatedData.h Atom.h AtomExtra.h AtomMap.h AtomMask.h AxisType.h BaseIOtype.h Box.h BufferedLine.h CharMask.h ClusterDist.h ClusterList.h ClusterMap.h ClusterNode.h ClusterSieve.h Cmd.h CmdInput.h CmdList.h Command.h ComplexArray.h Constraints.h Control.h CoordinateInfo.h Corr.h Cph.h CpptrajFile.h CpptrajState.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_1D.h DataSet_2D.h DataSet_3D.h DataSet_Cmatrix.h DataSet_Coords.h DataSet_Coords_CRD.h DataSet_Coords_REF.h DataSet_GridFlt.h DataSet_Mat3x3.h DataSet_MatrixDbl.h DataSet_MatrixFlt.h DataSet_Mesh.h DataSet_Modes.h DataSet_RemLog.h DataSet_Vector.h DataSet_double.h DataSet_float.h DataSet_integer.h DataSet_pH.h DataSet_string.h Deprecated.h DihedralSearch.h Dimension.h DispatchObject.h DistRoutines.h Energy.h Energy_Sander.h EnsembleIn.h EnsembleOut.h EnsembleOutList.h Ewald.h Exec.h Exec_Analyze.h Exec_Calc.h Exec_Change.h Exec_ClusterMap.h Exec_CombineCoords.h Exec_Commands.h Exec_CompareTop.h Exec_CrdAction.h Exec_CrdOut.h Exec_CreateSet.h Exec_DataFile.h Exec_DataFilter.h Exec_DataSetCmd.h Exec_GenerateAmberRst.h Exec_Help.h Exec_LoadCrd.h Exec_LoadTraj.h Exec_ParallelAnalysis.h Exec_ParmBox.h Exec_ParmSolvent.h Exec_ParmStrip.h Exec_ParmWrite.h Exec_PermuteDihedrals.h Exec_Precision.h Exec_PrintData.h Exec_ReadData.h Exec_ReadEnsembleData.h Exec_ReadInput.h Exec_RotateDihedral.h Exec_RunAnalysis.h Exec_ScaleDihedralK.h Exec_SequenceAlign.h Exec_SortEnsembleData.h Exec_SplitCoords.h Exec_System.h Exec_Top.h Exec_Traj.h Exec_UpdateParameters.h Exec_ViewRst.h FileIO.h FileName.h FileTypes.h Frame.h FrameArray.h FramePtrArray.h Grid.h GridAction.h GridBin.h HistBin.h Hungarian.h ImageTypes.h ImagedAction.h InputTrajCommon.h MapAtom.h MaskToken.h Matrix.h Matrix_3x3.h MetaData.h Molecule.h NameType.h NetcdfFile.h OnlineVarT.h OutputTrajCommon.h PDBfile.h PairList.h Parallel.h ParameterHolders.h ParameterTypes.h PubFFT.h RPNcalc.h Random.h Range.h ReferenceAction.h ReferenceFrame.h RemdReservoirNC.h ReplicaDimArray.h ReplicaInfo.h Residue.h Spline.h StructureCheck.h SymbolExporting.h SymmetricRmsdCalc.h TextFormat.h Timer.h Topology.h TrajFrameCounter.h TrajectoryFile.h TrajectoryIO.h Trajin.h TrajinList.h TrajoutList.h Trajout_Single.h VariableArray.h Vec3.h molsurf.h ComplexArray.o : ComplexArray.cpp ArrayIterator.h ComplexArray.h Constraints.o : Constraints.cpp ArgList.h Atom.h AtomExtra.h AtomMask.h Box.h CharMask.h Constants.h Constraints.h CoordinateInfo.h CpptrajStdio.h FileName.h Frame.h MaskToken.h Matrix_3x3.h Molecule.h NameType.h Parallel.h ParameterTypes.h Range.h ReplicaDimArray.h Residue.h SymbolExporting.h Topology.h Vec3.h Control.o : Control.cpp Action.h ActionFrameCounter.h ActionList.h ActionState.h Analysis.h AnalysisList.h AnalysisState.h ArgList.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h Control.h CoordinateInfo.h CpptrajFile.h CpptrajState.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_Coords.h DataSet_Coords_REF.h Dimension.h DispatchObject.h EnsembleIn.h EnsembleOut.h EnsembleOutList.h FileIO.h FileName.h FileTypes.h Frame.h FrameArray.h FramePtrArray.h InputTrajCommon.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h OutputTrajCommon.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h ReplicaInfo.h Residue.h StringRoutines.h SymbolExporting.h TextFormat.h Timer.h Topology.h TrajFrameCounter.h TrajectoryFile.h TrajectoryIO.h Trajin.h TrajinList.h TrajoutList.h Trajout_Single.h VariableArray.h Vec3.h @@ -249,6 +250,7 @@ Exec_GenerateAmberRst.o : Exec_GenerateAmberRst.cpp Action.h ActionFrameCounter. Exec_Help.o : Exec_Help.cpp Action.h ActionFrameCounter.h ActionList.h ActionState.h Analysis.h AnalysisList.h AnalysisState.h ArgList.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h BondSearch.h Box.h CharMask.h Cmd.h CmdList.h Command.h Control.h CoordinateInfo.h CpptrajFile.h CpptrajState.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_Coords.h DataSet_Coords_REF.h Dimension.h DispatchObject.h EnsembleIn.h EnsembleOut.h EnsembleOutList.h Exec.h Exec_Help.h FileIO.h FileName.h FileTypes.h Frame.h FrameArray.h FramePtrArray.h InputTrajCommon.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h OutputTrajCommon.h Parallel.h ParameterTypes.h ParmFile.h ParmIO.h Range.h ReferenceFrame.h ReplicaDimArray.h ReplicaInfo.h Residue.h SymbolExporting.h TextFormat.h Timer.h Topology.h TrajFrameCounter.h TrajectoryFile.h TrajectoryIO.h Trajin.h TrajinList.h TrajoutList.h Trajout_Single.h VariableArray.h Vec3.h Exec_LoadCrd.o : Exec_LoadCrd.cpp Action.h ActionFrameCounter.h ActionList.h ActionState.h Analysis.h AnalysisList.h AnalysisState.h ArgList.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h CoordinateInfo.h CpptrajFile.h CpptrajState.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_Coords.h DataSet_Coords_REF.h Dimension.h DispatchObject.h EnsembleIn.h EnsembleOut.h EnsembleOutList.h Exec.h Exec_LoadCrd.h FileIO.h FileName.h FileTypes.h Frame.h FrameArray.h FramePtrArray.h InputTrajCommon.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h OutputTrajCommon.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h ReplicaInfo.h Residue.h SymbolExporting.h TextFormat.h Timer.h Topology.h TrajFrameCounter.h TrajectoryFile.h TrajectoryIO.h Trajin.h TrajinList.h Trajin_Single.h TrajoutList.h Trajout_Single.h Vec3.h Exec_LoadTraj.o : Exec_LoadTraj.cpp Action.h ActionFrameCounter.h ActionList.h ActionState.h Analysis.h AnalysisList.h AnalysisState.h ArgList.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h CoordinateInfo.h CpptrajFile.h CpptrajState.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_Coords.h DataSet_Coords_REF.h DataSet_Coords_TRJ.h Dimension.h DispatchObject.h EnsembleIn.h EnsembleOut.h EnsembleOutList.h Exec.h Exec_LoadTraj.h FileIO.h FileName.h FileTypes.h Frame.h FrameArray.h FramePtrArray.h InputTrajCommon.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h OutputTrajCommon.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h ReplicaInfo.h Residue.h SymbolExporting.h TextFormat.h Timer.h Topology.h TrajFrameCounter.h TrajFrameIndex.h TrajectoryFile.h TrajectoryIO.h Trajin.h TrajinList.h TrajoutList.h Trajout_Single.h Vec3.h +Exec_ParallelAnalysis.o : Exec_ParallelAnalysis.cpp Action.h ActionFrameCounter.h ActionList.h ActionState.h Analysis.h AnalysisList.h AnalysisState.h ArgList.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h CoordinateInfo.h CpptrajFile.h CpptrajState.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_Coords.h DataSet_Coords_REF.h Dimension.h DispatchObject.h EnsembleIn.h EnsembleOut.h EnsembleOutList.h Exec.h Exec_ParallelAnalysis.h FileIO.h FileName.h FileTypes.h Frame.h FrameArray.h FramePtrArray.h InputTrajCommon.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h OutputTrajCommon.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h ReplicaInfo.h Residue.h SymbolExporting.h TextFormat.h Timer.h Topology.h TrajFrameCounter.h TrajectoryFile.h TrajectoryIO.h Trajin.h TrajinList.h TrajoutList.h Trajout_Single.h Vec3.h Exec_ParmBox.o : Exec_ParmBox.cpp Action.h ActionFrameCounter.h ActionList.h ActionState.h Analysis.h AnalysisList.h AnalysisState.h ArgList.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h CoordinateInfo.h CpptrajFile.h CpptrajState.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_Coords.h DataSet_Coords_REF.h Dimension.h DispatchObject.h EnsembleIn.h EnsembleOut.h EnsembleOutList.h Exec.h Exec_ParmBox.h FileIO.h FileName.h FileTypes.h Frame.h FrameArray.h FramePtrArray.h InputTrajCommon.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h OutputTrajCommon.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h ReplicaInfo.h Residue.h SymbolExporting.h TextFormat.h Timer.h Topology.h TrajFrameCounter.h TrajectoryFile.h TrajectoryIO.h Trajin.h TrajinList.h TrajoutList.h Trajout_Single.h Vec3.h Exec_ParmSolvent.o : Exec_ParmSolvent.cpp Action.h ActionFrameCounter.h ActionList.h ActionState.h Analysis.h AnalysisList.h AnalysisState.h ArgList.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h CoordinateInfo.h CpptrajFile.h CpptrajState.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_Coords.h DataSet_Coords_REF.h Dimension.h DispatchObject.h EnsembleIn.h EnsembleOut.h EnsembleOutList.h Exec.h Exec_ParmSolvent.h FileIO.h FileName.h FileTypes.h Frame.h FrameArray.h FramePtrArray.h InputTrajCommon.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h OutputTrajCommon.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h ReplicaInfo.h Residue.h SymbolExporting.h TextFormat.h Timer.h Topology.h TrajFrameCounter.h TrajectoryFile.h TrajectoryIO.h Trajin.h TrajinList.h TrajoutList.h Trajout_Single.h Vec3.h Exec_ParmStrip.o : Exec_ParmStrip.cpp Action.h ActionFrameCounter.h ActionList.h ActionState.h Analysis.h AnalysisList.h AnalysisState.h ArgList.h AssociatedData.h Atom.h AtomExtra.h AtomMask.h BaseIOtype.h Box.h CharMask.h CoordinateInfo.h CpptrajFile.h CpptrajState.h CpptrajStdio.h DataFile.h DataFileList.h DataIO.h DataSet.h DataSetList.h DataSet_Coords.h DataSet_Coords_REF.h Dimension.h DispatchObject.h EnsembleIn.h EnsembleOut.h EnsembleOutList.h Exec.h Exec_ParmStrip.h FileIO.h FileName.h FileTypes.h Frame.h FrameArray.h FramePtrArray.h InputTrajCommon.h MaskToken.h Matrix_3x3.h MetaData.h Molecule.h NameType.h OutputTrajCommon.h Parallel.h ParameterTypes.h Range.h ReferenceFrame.h ReplicaDimArray.h ReplicaInfo.h Residue.h SymbolExporting.h TextFormat.h Timer.h Topology.h TrajFrameCounter.h TrajectoryFile.h TrajectoryIO.h Trajin.h TrajinList.h TrajoutList.h Trajout_Single.h Vec3.h diff --git a/src/cpptrajfiles b/src/cpptrajfiles index c595430641..f9395fd816 100644 --- a/src/cpptrajfiles +++ b/src/cpptrajfiles @@ -97,6 +97,7 @@ COMMON_SOURCES=ActionFrameCounter.cpp \ Analysis_CurveFit.cpp \ Analysis_Divergence.cpp \ Analysis_FFT.cpp \ + Analysis_HausdorffDistance.cpp \ Analysis_Hist.cpp \ Analysis_Integrate.cpp \ Analysis_IRED.cpp \ @@ -246,6 +247,7 @@ COMMON_SOURCES=ActionFrameCounter.cpp \ Exec_Help.cpp \ Exec_LoadCrd.cpp \ Exec_LoadTraj.cpp \ + Exec_ParallelAnalysis.cpp \ Exec_ParmBox.cpp \ Exec_ParmSolvent.cpp \ Exec_ParmStrip.cpp \ diff --git a/test/Makefile b/test/Makefile index 86043702da..f2e7ef3939 100644 --- a/test/Makefile +++ b/test/Makefile @@ -440,6 +440,9 @@ test.cmdline: test.xyzfmt: @-cd Test_XYZformat && ./RunTest.sh $(OPT) +test.hausdorff: + @-cd Test_Hausdorff && ./RunTest.sh $(OPT) + # Every test target should go here. COMPLETETESTS=test.general \ test.strip \ @@ -581,7 +584,8 @@ COMPLETETESTS=test.general \ test.infraredspec \ test.cphstats \ test.cmdline \ - test.xyzfmt + test.xyzfmt \ + test.hausdorff test.all: $(MAKE) test.complete summary diff --git a/test/Test_Hausdorff/RunTest.sh b/test/Test_Hausdorff/RunTest.sh new file mode 100755 index 0000000000..072765b699 --- /dev/null +++ b/test/Test_Hausdorff/RunTest.sh @@ -0,0 +1,93 @@ +#!/bin/bash + +. ../MasterTest.sh + +CleanFiles matrix.dat hd.in hd.dat rms2d.dat rms2d.gnu hausdorff.matrix.dat \ + hausdorff.matrix.gnu hausdorff.fullmatrix.gnu + +TESTNAME='Hausdorff distance tests.' +Requires netcdf + +INPUT='-i hd.in' + +# Simple Hausdorff distance from matrix. +cat > matrix.dat < hd.in < hd.in < hd.in <> hd.in + done +done +cat >> hd.in < hd.in <> hd.in <> hd.in <