Setup.

CMakeLists.txt

+# CONFIGURE CMAKE
+cmake_minimum_required(VERSION 2.8.3)
+set(CMAKE_INSTALL_PREFIX ${CMAKE_SOURCE_DIR})
+set(LOCAL_INSTALL_DIR ${CMAKE_INSTALL_PREFIX}/soap)
+set(CMAKE_MODULE_PATH ${CMAKE_CURRENT_SOURCE_DIR}/cmake)
+
+# BUILD OPTIONS
+option(BUILD_SHARED_LIBS "Build shared libs" ON)
+
+# LOCAL PACKAGES
+#set(BOOST_ROOT "/home/cp605/packages/install_boost/boost_1_60_0/libboost")
+#set(ESPP_ROOT "/home/cp605/nomad/espressopp")
+#set(ESPP_LIBRARIES ${ESPP_ROOT}/_espressopp.so)
+#include_directories(${ESPP_INCLUDE_DIRS})
+
+# FIND PACKAGES
+find_package(PythonLibs)
+include_directories(${PYTHON_INCLUDE_DIRS})
+
+find_package(Boost COMPONENTS python mpi filesystem serialization)
+include_directories(${Boost_INCLUDE_DIRS})
+
+find_package(MPI REQUIRED)
+include_directories(${MPI_INCLUDE_PATH})
+
+find_package(FFTW3 REQUIRED)
+include_directories(${FFTW3_INCLUDES})
+
+find_package(GSL REQUIRED)
+include_directories(${GSL_INCLUDE_DIRS})
+
+#link_libraries(${ESPP_LIBRARIES} ${Boost_LIBRARIES} ${PYTHON_LIBRARIES} ${MPI_LIBRARIES} ${FFTW3_LIBRARIES})
+
+# SOME FRIENDLY MESSAGING
+get_property(local_dirs DIRECTORY ${CMAKE_CURRENT_SOURCE_DIR} PROPERTY INCLUDE_DIRECTORIES)
+get_property(local_libs DIRECTORY ${CMAKE_CURRENT_SOURCE_DIR} PROPERTY LINK_LIBRARIES)
+
+message(STATUS "Include directories: ")
+foreach(dir ${local_dirs})
+  message(STATUS " o ${dir}")
+endforeach()
+message(STATUS "Linked libraries:    ")
+foreach(lib ${local_dirs})
+  message(STATUS " o ${lib}")
+endforeach()
+
+#message(STATUS "Espresso++ lib:     ")
+#message(STATUS " o ${ESPP_LIBRARIES}")
+
+# SUBDIRECTORIES
+add_subdirectory(src)
+

README.md

+# soapxx

cmake/FindFFTW3.cmake

+# - Find FFTW3
+# Find the native FFTW3 includes and library
+#
+# FFTW3_INCLUDES - where to find fftw3.h
+# FFTW3_LIBRARIES - List of libraries when using FFTW3.
+# FFTW3_FOUND - True if FFTW3 found.
+
+if (FFTW3_INCLUDES)
+  # Already in cache, be silent
+  set (FFTW3_FIND_QUIETLY TRUE)
+endif (FFTW3_INCLUDES)
+
+find_path (FFTW3_INCLUDES fftw3.h)
+
+find_library (FFTW3_LIBRARIES NAMES fftw3)
+
+# handle the QUIETLY and REQUIRED arguments and set FFTW3_FOUND to TRUE if
+# all listed variables are TRUE
+include (FindPackageHandleStandardArgs)
+find_package_handle_standard_args (FFTW3 DEFAULT_MSG FFTW3_LIBRARIES FFTW3_INCLUDES)
+
+mark_as_advanced (FFTW3_LIBRARIES FFTW3_INCLUDES)

cmake/FindGSL.cmake

+# - Find gsl
+# Find the native GSL headers and libraries.
+#
+#  GSL_INCLUDE_DIRS - where to find gsl/gsl_linalg.h, etc.
+#  GSL_LIBRARIES    - List of libraries when using gsl.
+#  GSL_FOUND        - True if gsl found.
+#
+# Copyright 2009-2011 The VOTCA Development Team (http://www.votca.org)
+#
+# Licensed under the Apache License, Version 2.0 (the "License");
+# you may not use this file except in compliance with the License.
+# You may obtain a copy of the License at
+#
+#     http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+#
+
+find_package(PkgConfig)
+
+pkg_check_modules(PC_GSL gsl)
+find_path(GSL_INCLUDE_DIR gsl/gsl_linalg.h HINTS ${PC_GSL_INCLUDE_DIRS})
+find_library(GSL_LIBRARY NAMES gsl HINTS ${PC_GSL_LIBRARY_DIRS} )
+
+include(FindPackageHandleStandardArgs)
+# handle the QUIETLY and REQUIRED arguments and set GSL_FOUND to TRUE
+# if all listed variables are TRUE
+find_package_handle_standard_args(GSL DEFAULT_MSG GSL_LIBRARY GSL_INCLUDE_DIR)
+
+set(GSL_LIBRARIES "${GSL_LIBRARY}")
+set(GSL_INCLUDE_DIRS ${GSL_INCLUDE_DIR} )
+
+if (GSL_FOUND)
+  include(CheckLibraryExists)
+  #adding MATH_LIBRARIES here to allow static libs, this does not harm us as we are anyway using it
+  check_library_exists("${GSL_LIBRARIES};${MATH_LIBRARIES}" gsl_linalg_QR_decomp "" FOUND_QR_DECOMP)
+  if(NOT FOUND_QR_DECOMP)
+    #Some distributions like OpenSUSE need cblas in the linker flags, too
+    find_library(GSLCBLAS_LIBRARY NAMES gslcblas cblas HINTS ${PC_GSL_LIBRARY_DIRS} )
+    find_package_handle_standard_args(GSLCBLAS DEFAULT_MSG GSLCBLAS_LIBRARY)
+    if (GSLCBLAS_FOUND)
+      check_library_exists("${GSLCBLAS_LIBRARY};${MATH_LIBRARIES}" cblas_dsyrk "" FOUND_DSYRK)
+      if(NOT FOUND_DSYRK)
+        message(FATAL_ERROR "Could not find cblas_dsyrk in ${GSLCBLAS_LIBRARY}, take a look at the error message in ${CMAKE_BINARY_DIR}${CMAKE_FILES_DIRECTORY}/CMakeError.log to find out what was going wrong. If you are using a static lib (.a) make sure you have specified all dependencies of libcblas in GSLCBLAS_LIBRARY by hand (i.e. -DGSLCBLAS_LIBRARY='/path/to/libcblas.so;/path/to/libm.so')! If your gsl was build against an different version of cblas, specify it in GSLCBLAS_LIBRARY")
+      endif(NOT FOUND_DSYRK)
+      set(GSL_LIBRARIES "${GSL_LIBRARY};${GSLCBLAS_LIBRARY}")
+      check_library_exists("${GSL_LIBRARIES};${MATH_LIBRARIES}" gsl_linalg_QR_lssolve "" FOUND_QR_DECOMP_AGAIN)
+    endif(GSLCBLAS_FOUND)
+  endif(NOT FOUND_QR_DECOMP)
+  if(NOT FOUND_QR_DECOMP AND NOT FOUND_QR_DECOMP_AGAIN)
+      message(FATAL_ERROR "Could not find gsl_linalg_QR_decompx in ${GSL_LIBRARY}, take a look at the error message in ${CMAKE_BINARY_DIR}${CMAKE_FILES_DIRECTORY}/CMakeError.log to find out what was going wrong.  If you are using a static lib (.a) make sure you have specified all dependencies of libgsl in GSL_LIBRARY by hand (i.e. -DGSL_LIBRARY='/path/to/libgsl.so;/path/to/libm.so') !")
+  endif(NOT FOUND_QR_DECOMP AND NOT FOUND_QR_DECOMP_AGAIN)
+endif (GSL_FOUND)
+
+mark_as_advanced(GSL_INCLUDE_DIR GSL_LIBRARY)

src/CMakeLists.txt

+file(GLOB_RECURSE LOCAL_SOURCES *.cpp)
+
+add_library(_soapxx ${LOCAL_SOURCES})
+target_link_libraries(_soapxx ${Boost_LIBRARIES} ${PYTHON_LIBRARIES} ${MPI_LIBRARIES} ${GSL_LIBRARIES})
+set_target_properties(_soapxx PROPERTIES PREFIX "" SUFFIX ".so" LIBRARY_OUTPUT_DIRECTORY .)
+
+install(TARGETS _soapxx LIBRARY DESTINATION ${LOCAL_INSTALL_DIR})
+install(FILES __init__.py DESTINATION ${LOCAL_INSTALL_DIR})
+#install(DIRECTORY linalg DESTINATION ${LOCAL_INSTALL_DIR})
+
+add_subdirectory(linalg)
+add_subdirectory(tools)
+

src/__init__.py

+from _soapxx import *
+from linalg import *
+

src/base/logger.hpp

+#ifndef _SOAP_LOGGER_H
+#define	_SOAP_LOGGER_H
+
+#include <sstream>
+#include <iostream>
+
+namespace soap { 
+
+
+
+enum TLogLevel {logERROR, logWARNING, logINFO, logDEBUG};
+ 
+/*
+ * Macros to use the Logger: LOG(level,logger) << message
+ */
+#define LOGIF(level, log) \
+if ( &log != NULL && level > (log).getReportLevel() ) ; \
+else (log)(level)
+
+#define LOG(plog) \
+if (plog == NULL) ; \
+else (*plog)()
+
+/*
+ * Custom buffer to store messages
+ */
+class LogBuffer : public std::stringbuf {
+
+public:
+	LogBuffer() : std::stringbuf(),
+                _errorPreface("ERROR "), _warnPreface("WARNING "),
+                _infoPreface(""), _dbgPreface("DEBUG "),
+                _writePreface(true) {}
+        
+        // sets the log level (needed for output)
+	void setLogLevel(TLogLevel LogLevel) { _LogLevel = LogLevel; }
+        
+        // sets Multithreading (buffering required)
+        void setMultithreading( bool maverick ) { _maverick = maverick; }
+        
+        // sets preface strings for logERROR, logWARNING, ...
+        void setPreface(TLogLevel level, std::string preface) {            
+            switch ( level )
+            {
+
+                case logERROR: 
+                    _errorPreface = preface;
+                    break;
+                case logWARNING:
+                    _warnPreface = preface;
+                    break;
+                case logINFO:
+                    _infoPreface = preface;
+                    break;
+                case logDEBUG:
+                    _dbgPreface = preface;
+                    break;
+            }
+        }
+        
+        void EnablePreface() { _writePreface = true; }
+        void DisablePreface() { _writePreface = false; }
+        
+        // flushes all collected messages
+        void FlushBuffer(){ std::cout << _stringStream.str(); _stringStream.str(""); }
+        
+        // returns the pointer to the collected messages
+        std::string Messages() { 
+            std::string _messages = _stringStream.str(); 
+            _stringStream.str("");
+            return _messages; 
+        }
+
+private:
+  
+  // Log Level (WARNING, INFO, etc) 
+  TLogLevel _LogLevel;
+  
+  // temporary buffer to store messages
+  std::ostringstream _stringStream;
+  
+  // Multithreading
+  bool _maverick;
+  bool _writePreface;
+  
+  std::string _timePreface;
+  std::string _errorPreface;
+  std::string _warnPreface;
+  std::string _infoPreface;
+  std::string _dbgPreface;
+  
+
+protected:
+	virtual int sync() {
+            
+            std::ostringstream _message;
+
+            if (_writePreface) {
+                switch ( _LogLevel )
+                {
+                    case logERROR: 
+                        _message << _errorPreface;
+                        break;
+                    case logWARNING:
+                        _message << _warnPreface;
+                        break;      
+                    case logINFO:
+                        _message << _infoPreface;
+                        break;      
+                    case logDEBUG:
+                        _message << _dbgPreface;
+                        break;      
+                }
+            }
+
+            if ( !_maverick ) {
+                // collect all messages of one thread
+                _stringStream << _message.str()  << "" << str();
+            } else {
+                // if only one thread outputs, flush immediately
+                std::cout << _message.str() << "" << str() << std::flush;
+            }
+            _message.str("");
+	    str("");
+	    return 0;
+	}
+      
+};
+
+
+/** \class Logger
+*   \brief Logger is used for thread-safe output of messages
+*
+*  Logger writes messages into LogBuffer. 
+*  Inheritance from ostream allows to use overloaded << and >> for writing.
+*  Example:  
+*
+*  \code
+*  #include <votca/ctp/logger.h>
+*  Logger* log = new Logger(); // create a logger object
+*  log->setReportLevel(logDEBUG); // output only log messages starting from a DEBUG level
+*  LOG(logERROR,*log) << "Error detected" << flush; // write to the logger at an ERROR level
+*  cout << log; // output logger content to standard output
+*  \endcode
+*
+*  Logger has four predefined log levels: logERROR, logWARNING, logINFO, logDEBUG.
+*/
+class Logger : public std::ostream {
+
+       friend std::ostream& operator<<( std::ostream& out, Logger&  logger ) {
+           out << logger.Messages();
+           return out;
+       }
+    
+public:
+	Logger( TLogLevel ReportLevel) : std::ostream(new LogBuffer()) {
+            _ReportLevel = ReportLevel; 
+            _maverick = true;
+     }
+
+	 Logger() : std::ostream(new LogBuffer()) {
+		 _ReportLevel = logDEBUG;
+		 _maverick = true;
+		 dynamic_cast<LogBuffer *>( rdbuf() )->setMultithreading(_maverick);
+	 }
+        
+	~Logger() {
+            //dynamic_cast<LogBuffer *>( rdbuf())->FlushBuffer();
+            delete rdbuf();
+            rdbuf(NULL);
+	}
+        
+	Logger &operator()( TLogLevel LogLevel = logINFO) {
+		//rdbuf()->pubsync();
+		dynamic_cast<LogBuffer *>( rdbuf() )->setLogLevel(LogLevel);
+		return *this;
+	}
+        
+        void setReportLevel( TLogLevel ReportLevel ) { _ReportLevel = ReportLevel; }
+        void setMultithreading( bool maverick ) { 
+            _maverick = maverick;
+            dynamic_cast<LogBuffer *>( rdbuf() )->setMultithreading( _maverick );
+        }
+        bool isMaverick() { return _maverick; }
+        
+        TLogLevel getReportLevel( ) { return _ReportLevel; }
+        
+        void setPreface(TLogLevel level, std::string preface) {            
+            dynamic_cast<LogBuffer *>( rdbuf() )->setPreface(level, preface);
+        }
+        
+        void EnablePreface() {
+            dynamic_cast<LogBuffer *>( rdbuf() )->EnablePreface();
+        }
+        
+        void DisablePreface() {
+            dynamic_cast<LogBuffer *>( rdbuf() )->DisablePreface();
+        }
+        
+private:
+    // at what level of detail output messages
+    TLogLevel _ReportLevel;
+    
+    // if true, only a single processor job is executed
+    bool      _maverick;
+    
+    std::string Messages() {
+        return dynamic_cast<LogBuffer *>( rdbuf() )->Messages();
+    }
+        
+};
+
+/**
+*   \brief Timestamp returns the current time as a string
+*  Example: cout << TimeStamp()
+*/
+class TimeStamp 
+{
+  public:
+    friend std::ostream & operator<<(std::ostream &os, const TimeStamp& ts)
+    {
+        time_t rawtime;
+        tm * timeinfo;
+        time(&rawtime);
+        timeinfo = localtime( &rawtime );
+        os  << (timeinfo->tm_year)+1900
+            << "-" << timeinfo->tm_mon + 1
+            << "-" << timeinfo->tm_mday 
+            << " " << timeinfo->tm_hour
+            << ":" << timeinfo->tm_min 
+            << ":"  << timeinfo->tm_sec;
+         return os;    
+    }
+    
+    explicit TimeStamp() {};
+    
+};
+
+}
+
+#endif

src/base/objectfactory.hpp

+#ifndef _objectfactory_H
+#define	_objectfactory_H
+
+#include <map>
+#include <list>
+#include <iostream>
+#include <stdexcept>
+#include <boost/lexical_cast.hpp>
+
+namespace soap { namespace base {
+
+/**
+    The REGISTER_OBJECT macro allows to easily register an object in an object factory.
+ */
+#define REGISTER_OBJECT(factory, object, key) \
+    namespace { \
+        ObjectFactoryRegister<object> \
+            _register_##object(factory, key); \
+    }
+
+/**
+    \brief template class for object factory
+
+    This class is a template for an object factory. The factory creates an instance of an derived class,
+    be giving a key (e.g. a string) for the object which identifies it uniquely. This allows the implementation
+    of new features (e.g. new file formats, new mapping algorithms) without touching or recompiling existing bits of code.
+
+    If you don't understand this, read the book by Alexandresku (Modern C++ design)
+    everything explained there in detail!
+*/
+template<typename key_t, typename T>
+class ObjectFactory
+{
+private:
+    typedef T* (*creator_t)();
+public:
+
+    typedef T abstract_type;
+    typedef std::map<key_t, creator_t> assoc_map;
+    
+    ObjectFactory() {}
+    ~ObjectFactory() {};
+    
+    /**
+     * \brief register an object
+     * \param key identifier
+     * \param creator create policy
+     *
+     * This function is called to register an object in the factory. After an object is registered,
+     * an instance of it can be created by calling Create specifying the corresponding key.
+     */
+    void Register(const key_t &key, creator_t creator );
+
+    template< typename obj_t >
+    void Register(const key_t &key);
+
+    /**
+       Create an instance of the object identified by key.
+    */
+    T *Create(const key_t &key);
+    bool IsRegistered(const key_t & _id) const;
+
+    static ObjectFactory<key_t, T>& Instance()
+    {
+        static ObjectFactory<key_t, T> _this;
+        return _this;
+    }
+
+    const assoc_map &getObjects() { return _objects; }
+private:
+    assoc_map _objects;
+};
+
+template<class parent, class T> parent* create_policy_new()
+{
+    return new T();
+}
+
+template<typename key_t, typename T>
+inline void ObjectFactory<key_t, T>::Register(const key_t &key, creator_t creator)
+{
+    (void)_objects.insert(typename assoc_map::value_type(key, creator)).second;
+}
+
+template<typename key_t, typename T>
+template< typename obj_t >
+inline void ObjectFactory<key_t, T>::Register(const key_t &key)
+{
+    Register(key, create_policy_new<abstract_type, obj_t>);
+}
+
+
+template<typename key_t, typename T>
+inline T* ObjectFactory<key_t, T>::Create(const key_t &key)
+{
+    typename assoc_map::const_iterator it(_objects.find(key));
+    if (it != _objects.end())
+        return (it->second)();
+    else
+        throw std::runtime_error("factory key " + boost::lexical_cast<std::string>(key) + " not found.");
+}
+
+/*template<typename key_t, typename T>
+inline static ObjectFactory<key_t, T>& ObjectFactory<key_t, T>::Instance()
+{
+    static ObjectFactory<key_t, T> _this;
+    return _this;
+}*/
+
+template<typename key_t, typename T>
+inline bool ObjectFactory<key_t, T>::IsRegistered(const key_t & _id) const
+{
+	return ( _objects.find(_id)!= _objects.end() );
+}
+
+/*std::string list_keys() const {
+    std::stringstream _str;
+    for (typename assoc_map::const_iterator it(map.begin()); it != map.end(); it++) {
+        _str << (*it).first << "\n";
+    }
+    return _str.str();
+}
+*/
+
+template<typename object_type >
+class ObjectFactoryRegister {
+public:
+    template<typename factory_type, typename key_type>
+    ObjectFactoryRegister(factory_type &factory, key_type &key) {
+        factory.Register(key, &create_policy_new<typename factory_type::abstract_type, object_type>);
+    }
+};
+
+}}
+
+#endif	/* _objectfactory_H */

src/bindings.cpp

+#include "bindings.hpp"
+
+namespace soap {
+
+}
+
+BOOST_PYTHON_MODULE(_soapxx)
+{
+    using namespace boost::python;
+    soap::Particle::registerPython();
+    soap::Segment::registerPython();
+    soap::Structure::registerPython();
+    soap::Options::registerPython();
+    soap::Spectrum::registerPython();
+
+    soap::RadialBasisFactory::registerAll();
+}
+

src/bindings.hpp

+#ifndef _SOAP_LINALG_BINDINGS_H
+#define	_SOAP_LINALG_BINDINGS_H
+#include <boost/python.hpp>
+#include "structure.hpp"
+#include "options.hpp"
+#include "spectrum.hpp"
+
+namespace soap {
+
+}
+
+
+#endif
+

src/boundary.hpp

+#ifndef _SOAP_BOUNDARY_HPP
+#define	_SOAP_BOUNDARY_HPP
+
+#include "linalg/matrix.hpp"
+
+namespace soap {
+
+class BoundaryOpen;
+class BoundaryOrthorhombic;
+class BoundaryTriclinic;
+
+class Boundary
+{
+public:
+
+	enum eBoxType {
+		typeOpen = 0, typeOrthorhombic, typeTriclinic
+	};
+
+    virtual ~Boundary() {;}
+
+    void setBox(const matrix &box) {
+        _box = box;
+    }
+
+    const matrix &getBox() { return _box; };
+    virtual double BoxVolume() {
+		vec a = _box.getCol(0);
+		vec b = _box.getCol(1);
+		vec c = _box.getCol(2);
+		return (a^b)*c;
+    }
+    virtual vec connect(const vec &r_i, const vec &r_j) const = 0;
+
+    virtual eBoxType getBoxType() { return _type; }
+
+protected:
+    matrix _box;
+    eBoxType _type;
+};
+
+
+class BoundaryOpen : public Boundary
+{
+public:
+	BoundaryOpen(const matrix &box) {
+		 _type = Boundary::typeOpen;
+		 _box = box;
+	}
+    vec connect(const vec &r_i, const vec &r_j) const {
+    	return r_j - r_i;
+    }
+};