mirror/MUSIC
1
0
Fork 0
mirror of https://github.com/cosmo-sims/MUSIC.git synced 2024-09-16 13:33:46 +02:00
Code Releases Activity

WIP commit, transition to monofonic logger and config_file

Browse source
This commit is contained in:
Oliver Hahn 2023-02-14 10:34:58 -08:00
parent 5304b8141a
commit 0fa33b9ea9
52 changed files with 2802 additions and 2864 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

546
src/config_file.hh
View file

@ -1,67 +1,73 @@
/*
config_file.hh - This file is part of MUSIC -
a code to generate multi-scale initial conditions
for cosmological simulations
Copyright (C) 2010 Oliver Hahn
*/
// This file is part of monofonIC (MUSIC2)
// A software package to generate ICs for cosmological simulations
// Copyright (C) 2020 by Oliver Hahn
//
// monofonIC 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.
//
// monofonIC 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 program. If not, see <http://www.gnu.org/licenses/>.
#ifndef __CONFIG_FILE_HH
#define __CONFIG_FILE_HH
#pragma once
#include <string>
#include <sstream>
#include <map>
#include <fstream>
#include <iostream>
#include <iomanip>
#include <iostream>
#include <map>
#include <sstream>
#include <stdexcept>
#include <string>
#include <typeinfo>
#include "log.hh"
#include <logger.hh>
/*!
* @class config_file
* @brief provides read/write access to configuration options
*
* This class provides access to the configuration file. The
* This class provides access to the configuration file. The
* configuration is stored in hash-pairs and can be queried and
* validated by the responsible class/routine
*/
class config_file {
//! current line number
unsigned m_iLine;
unsigned iline_;
//! hash table for key/value pairs, stored as strings
std::map<std::string, std::string> m_Items;
std::map<std::string, std::string> items_;
public:
//! removes all white space from string source
/*!
* @param source the string to be trimmed
* @param delims a string of delimiting characters
* @return trimmed string
*/
std::string trim(std::string const& source, char const* delims = " \t\r\n") const{
std::string trim(std::string const &source,
char const *delims = " \t\r\n") const {
std::string result(source);
//... skip initial whitespace ...
std::string::size_type index = result.find_last_not_of(delims);
if(index != std::string::npos)
if (index != std::string::npos)
result.erase(++index);
//... find beginning of trailing whitespace ...
index = result.find_first_not_of(delims);
//... remove trailing whitespace ...
if(index != std::string::npos)
if (index != std::string::npos)
result.erase(0, index);
else
result.erase();
return result;
}
//! converts between different variable types
/*!
* The main purpose of this function is to parse and convert
@ -70,286 +76,306 @@ public:
* @param oval the interpreted/converted value
*/
template <class in_value, class out_value>
void convert( const in_value & ival, out_value & oval) const
{
std::stringstream ss;
ss << ival; //.. insert value into stream
ss >> oval; //.. retrieve value from stream
void convert(const in_value &ival, out_value &oval) const {
std::stringstream ss;
ss << ival; //.. insert value into stream
ss >> oval; //.. retrieve value from stream
if (! ss.eof()) {
//.. conversion error
std::cerr << "Error: conversion of \'" << ival << "\' failed." << std::endl;
throw ErrInvalidConversion(std::string("invalid conversion to ")+typeid(out_value).name()+'.');
}
if (!ss.eof()) {
//.. conversion error
music::elog << "Error: conversion of \'" << ival << "\' failed."
<< std::endl;
throw except_invalid_conversion(std::string("invalid conversion to ") +
typeid(out_value).name() + '.');
}
}
//! constructor of class config_file
/*! @param FileName the path/name of the configuration file to be parsed
/*! @param filename the path/name of the configuration file to be parsed
*/
config_file( std::string const& FileName )
: m_iLine(0), m_Items()
{
std::ifstream file(FileName.c_str());
if( !file.is_open() )
throw std::runtime_error(std::string("Error: Could not open config file \'")+FileName+std::string("\'"));
explicit config_file(std::string const &filename) : iline_(0), items_() {
std::ifstream file(filename.c_str());
if (!file.is_open()){
music::elog << "Could not open config file \'" << filename << "\'." << std::endl;
throw std::runtime_error(
std::string("Error: Could not open config file \'") + filename +
std::string("\'"));
}
std::string line;
std::string name;
std::string value;
std::string inSection;
int posEqual;
m_iLine=0;
std::string in_section;
int pos_equal;
iline_ = 0;
//.. walk through all lines ..
while (std::getline(file,line)) {
++m_iLine;
while (std::getline(file, line)) {
++iline_;
//.. encounterd EOL ?
if (! line.length()) continue;
if (!line.length())
continue;
//.. encountered comment ?
unsigned long idx;
if( (idx=line.find_first_of("#;%")) != std::string::npos )
if ((idx = line.find_first_of("#;%")) != std::string::npos)
line.erase(idx);
//.. encountered section tag ?
if (line[0] == '[') {
inSection=trim(line.substr(1,line.find(']')-1));
in_section = trim(line.substr(1, line.find(']') - 1));
continue;
}
//.. seek end of entry name ..
posEqual=line.find('=');
name = trim(line.substr(0,posEqual));
value = trim(line.substr(posEqual+1));
if( (size_t)posEqual==std::string::npos && (name.size()!=0||value.size()!=0) )
{
LOGWARN("Ignoring non-assignment in %s:%d",FileName.c_str(),m_iLine);
continue;
}
if(name.length()==0&&value.size()!=0)
{
LOGWARN("Ignoring assignment missing entry name in %s:%d",FileName.c_str(),m_iLine);
continue;
}
if(value.length()==0&&name.size()!=0)
{
LOGWARN("Empty entry will be ignored in %s:%d",FileName.c_str(),m_iLine);
continue;
}
if( value.length()==0&&name.size()==0)
continue;
pos_equal = line.find('=');
name = trim(line.substr(0, pos_equal));
value = trim(line.substr(pos_equal + 1));
if ((size_t)pos_equal == std::string::npos &&
(name.size() != 0 || value.size() != 0)) {
music::wlog << "Ignoring non-assignment in " << filename << ":"
<< iline_ << std::endl;
continue;
}
if (name.length() == 0 && value.size() != 0) {
music::wlog << "Ignoring assignment missing entry name in "
<< filename << ":" << iline_ << std::endl;
continue;
}
if (value.length() == 0 && name.size() != 0) {
music::wlog << "Empty entry will be ignored in " << filename << ":"
<< iline_ << std::endl;
continue;
}
if (value.length() == 0 && name.size() == 0)
continue;
//.. add key/value pair to hash table ..
if( m_Items.find(inSection+'/'+name) != m_Items.end() )
LOGWARN("Redeclaration overwrites previous value in %s:%d",FileName.c_str(),m_iLine);
m_Items[inSection+'/'+name] = value;
if (items_.find(in_section + '/' + name) != items_.end()) {
music::wlog << "Redeclaration overwrites previous value in "
<< filename << ":" << iline_ << std::endl;
}
items_[in_section + '/' + name] = value;
}
}
//! inserts a key/value pair in the hash map
/*! @param key the key value, usually "section/key"
* @param value the value of the key, also a string
*/
void insertValue( std::string const& key, std::string const& value )
{
m_Items[key] = value;
}
//! inserts a key/value pair in the hash map
/*! @param section section name. values are stored under "section/key"
* @param key the key value usually "section/key"
* @param value the value of the key, also a string
*/
void insertValue( std::string const& section, std::string const& key, std::string const& value )
{
m_Items[section+'/'+key] = value;
}
//! checks if a key is part of the hash map
/*! @param section the section name of the key
* @param key the key name to be checked
* @return true if the key is present, false otherwise
*/
bool containsKey( std::string const& section, std::string const& key )
{
std::map<std::string,std::string>::const_iterator i = m_Items.find(section+'/'+key);
if ( i == m_Items.end() )
return false;
return true;
}
//! checks if a key is part of the hash map
/*! @param key the key name to be checked
* @return true if the key is present, false otherwise
*/
bool containsKey( std::string const& key )
{
std::map<std::string,std::string>::const_iterator i = m_Items.find(key);
if ( i == m_Items.end() )
return false;
return true;
}
//! return value of a key
/*! returns the value of a given key, throws a ErrItemNotFound
* exception if the key is not available in the hash map.
* @param key the key name
* @return the value of the key
* @sa ErrItemNotFound
*/
template<class T> T getValue( std::string const& key ) const{
return getValue<T>( "", key );
}
//! return value of a key
/*! returns the value of a given key, throws a ErrItemNotFound
* exception if the key is not available in the hash map.
* @param section the section name for the key
* @param key the key name
* @return the value of the key
* @sa ErrItemNotFound
*/
template<class T> T getValue( std::string const& section, std::string const& key ) const
{
T r;
std::map<std::string,std::string>::const_iterator i = m_Items.find(section + '/' + key);
if ( i == m_Items.end() )
throw ErrItemNotFound('\'' + section + '/' + key + std::string("\' not found."));
convert(i->second,r);
return r;
}
//! inserts a key/value pair in the hash map
/*! @param key the key value, usually "section/key"
* @param value the value of the key, also a string
*/
void insert_value(std::string const &key, std::string const &value) {
items_[key] = value;
}
//! exception safe version of getValue
/*! returns the value of a given key, returns a default value rather
* than a ErrItemNotFound exception if the key is not found.
* @param section the section name for the key
* @param key the key name
* @param default_value the value that is returned if the key is not found
* @return the key value (if key found) otherwise default_value
*/
template<class T> T getValueSafe( std::string const& section, std::string const& key, T default_value ) const
{
T r;
try{
r = getValue<T>( section, key );
} catch( ErrItemNotFound& ) {
r = default_value;
}
return r;
}
//! inserts a key/value pair in the hash map
/*! @param section section name. values are stored under "section/key"
* @param key the key value usually "section/key"
* @param value the value of the key, also a string
*/
void insert_value(std::string const &section, std::string const &key,
std::string const &value) {
items_[section + '/' + key] = value;
}
//! exception safe version of getValue
/*! returns the value of a given key, returns a default value rather
* than a ErrItemNotFound exception if the key is not found.
* @param key the key name
* @param default_value the value that is returned if the key is not found
* @return the key value (if key found) otherwise default_value
*/
template<class T> T getValueSafe( std::string const& key, T default_value ) const
{
return getValueSafe( "", key, default_value );
}
//! checks if a key is part of the hash map
/*! @param section the section name of the key
* @param key the key name to be checked
* @return true if the key is present, false otherwise
*/
bool contains_key(std::string const &section, std::string const &key) {
std::map<std::string, std::string>::const_iterator i =
items_.find(section + '/' + key);
if (i == items_.end())
return false;
return true;
}
//! dumps all key-value pairs to a std::ostream
void dump( std::ostream& out )
{
std::map<std::string,std::string>::const_iterator i = m_Items.begin();
while( i!=m_Items.end() )
{
if( i->second.length() > 0 )
out << std::setw(24) << std::left << i->first << " = " << i->second << std::endl;
++i;
}
}
void log_dump( void )
{
LOGUSER("List of all configuration options:");
std::map<std::string,std::string>::const_iterator i = m_Items.begin();
while( i!=m_Items.end() )
{
if( i->second.length() > 0 )
LOGUSER(" %24s = %s",(i->first).c_str(),(i->second).c_str());//out << std::setw(24) << std::left << i->first << " = " << i->second << std::endl;
++i;
}
}
//! checks if a key is part of the hash map
/*! @param key the key name to be checked
* @return true if the key is present, false otherwise
*/
bool contains_key(std::string const &key) {
std::map<std::string, std::string>::const_iterator i = items_.find(key);
if (i == items_.end())
return false;
return true;
}
//! return value of a key
/*! returns the value of a given key, throws a except_item_not_found
* exception if the key is not available in the hash map.
* @param key the key name
* @return the value of the key
* @sa except_item_not_found
*/
template <class T> T get_value(std::string const &key) const {
return get_value<T>("", key);
}
//! return value of a key
/*! returns the value of a given key, throws a except_item_not_found
* exception if the key is not available in the hash map.
* @param section the section name for the key
* @param key the key name
* @return the value of the key
* @sa except_item_not_found
*/
template <class T>
T get_value_basic(std::string const &section, std::string const &key) const {
T r;
std::map<std::string, std::string>::const_iterator i =
items_.find(section + '/' + key);
if (i == items_.end()){
throw except_item_not_found('\'' + section + '/' + key +
std::string("\' not found."));
}
convert(i->second, r);
return r;
}
template <class T>
T get_value(std::string const &section, std::string const &key) const
{
T r;
try
{
r = get_value_basic<T>(section, key);
}
catch (except_item_not_found& e)
{
music::elog << e.what() << std::endl;
throw;
}
return r;
}
//! exception safe version of get_value
/*! returns the value of a given key, returns a default value rather
* than a except_item_not_found exception if the key is not found.
* @param section the section name for the key
* @param key the key name
* @param default_value the value that is returned if the key is not found
* @return the key value (if key found) otherwise default_value
*/
template <class T>
T get_value_safe(std::string const &section, std::string const &key,
T default_value) const {
T r;
try {
r = get_value_basic<T>(section, key);
} catch (except_item_not_found&) {
r = default_value;
music::dlog << "Item \'" << section << "/" << key << " not found in config. Default = \'" << default_value << "\'" << std::endl;
}
return r;
}
//! exception safe version of get_value
/*! returns the value of a given key, returns a default value rather
* than a except_item_not_found exception if the key is not found.
* @param key the key name
* @param default_value the value that is returned if the key is not found
* @return the key value (if key found) otherwise default_value
*/
template <class T>
T get_value_safe(std::string const &key, T default_value) const {
return get_value_safe("", key, default_value);
}
//! dumps all key-value pairs to a std::ostream
void dump(std::ostream &out) {
std::map<std::string, std::string>::const_iterator i = items_.begin();
while (i != items_.end()) {
if (i->second.length() > 0)
out << std::setw(24) << std::left << i->first << " = " << i->second
<< std::endl;
++i;
}
}
void dump_to_log(void) {
music::ulog << "List of all configuration options:" << std::endl;
std::map<std::string, std::string>::const_iterator i = items_.begin();
while (i != items_.end()) {
if (i->second.length() > 0)
music::ulog << std::setw(28) << i->first << " = " << i->second
<< std::endl;
++i;
}
}
//--- EXCEPTIONS ---
//! runtime error that is thrown if key is not found in getValue
class ErrItemNotFound : public std::runtime_error{
public:
ErrItemNotFound( std::string itemname )
: std::runtime_error( itemname.c_str() )
{}
};
//! runtime error that is thrown if key is not found in get_value
class except_item_not_found : public std::runtime_error {
public:
except_item_not_found(std::string itemname)
: std::runtime_error(itemname.c_str()) {}
};
//! runtime error that is thrown if type conversion fails
class ErrInvalidConversion : public std::runtime_error{
public:
ErrInvalidConversion( std::string errmsg )
: std::runtime_error( errmsg )
{}
};
//! runtime error that is thrown if type conversion fails
class except_invalid_conversion : public std::runtime_error {
public:
except_invalid_conversion(std::string errmsg) : std::runtime_error(errmsg) {}
};
//! runtime error that is thrown if identifier is not found in keys
class ErrIllegalIdentifier : public std::runtime_error{
public:
ErrIllegalIdentifier( std::string errmsg )
: std::runtime_error( errmsg )
{}
};
//! runtime error that is thrown if identifier is not found in keys
class ErrIllegalIdentifier : public std::runtime_error {
public:
ErrIllegalIdentifier(std::string errmsg) : std::runtime_error(errmsg) {}
};
};
//==== below are template specialisations =======================================//
//==== below are template specialisations
//=======================================//
//... Function: getValue( strSection, strEntry ) ...
//... Descript: specialization of getValue for type boolean to interpret strings ...
//... Function: get_value( strSection, strEntry ) ...
//... Descript: specialization of get_value for type boolean to interpret strings
//...
//... like "true" and "false" etc.
//... converts the string to type bool, returns type bool ...
template<>
inline bool config_file::getValue<bool>( std::string const& strSection, std::string const& strEntry ) const{
std::string r1 = getValue<std::string>( strSection, strEntry );
if( r1=="true" || r1=="yes" || r1=="on" || r1=="1" )
return true;
if( r1=="false" || r1=="no" || r1=="off" || r1=="0" )
return false;
throw ErrIllegalIdentifier(std::string("Illegal identifier \'")+r1+std::string("\' in \'")+strEntry+std::string("\'."));
//return false;
}
template <>
inline bool config_file::get_value<bool>(std::string const &strSection,
std::string const &strEntry) const {
std::string r1 = get_value<std::string>(strSection, strEntry);
if (r1 == "true" || r1 == "yes" || r1 == "on" || r1 == "1")
return true;
if (r1 == "false" || r1 == "no" || r1 == "off" || r1 == "0")
return false;
music::elog << "Illegal identifier \'" << r1 << "\' in \'" << strEntry << "\'." << std::endl;
throw ErrIllegalIdentifier(std::string("Illegal identifier \'") + r1 +
std::string("\' in \'") + strEntry +
std::string("\'."));
// return false;
}
template<>
inline bool config_file::getValueSafe<bool>( std::string const& strSection, std::string const& strEntry, bool defaultValue ) const{
template <>
inline bool config_file::get_value_safe<bool>(std::string const &strSection,
std::string const &strEntry,
bool defaultValue) const {
std::string r1;
try{
r1 = getValue<std::string>( strSection, strEntry );
if( r1=="true" || r1=="yes" || r1=="on" || r1=="1" )
try {
r1 = get_value_basic<std::string>(strSection, strEntry);
std::transform(r1.begin(), r1.end(),r1.begin(), ::toupper);
if (r1 == "YAY" || r1 == "TRUE" || r1 == "YES" || r1 == "ON" || r1 == "1")
return true;
if( r1=="false" || r1=="no" || r1=="off" || r1=="0" )
if (r1 == "NAY" || r1 == "FALSE" || r1 == "NO" || r1 == "OFF" || r1 == "0")
return false;
} catch( ErrItemNotFound& ) {
} catch (except_item_not_found&) {
return defaultValue;
}
return defaultValue;
}
template<>
inline void config_file::convert<std::string,std::string>( const std::string & ival, std::string & oval) const
{
oval = ival;
template <>
inline void
config_file::convert<std::string, std::string>(const std::string &ival,
std::string &oval) const {
oval = ival;
}
#endif //__CONFIG_FILE_HH

66
src/constraints.cc
View file

@ -86,8 +86,8 @@ void compute_sigma_tophat( config_file& cf, transfer_function *ptf, double R, st
double D0 = ccalc.CalcGrowthFactor(1.0);
double sigma8 = cf.getValue<double>("cosmology","sigma_8");
double nspec = cf.getValue<double>("cosmology","nspec");
double sigma8 = cf.get_value<double>("cosmology","sigma_8");
double nspec = cf.get_value<double>("cosmology","nspec");
double sigma0 = 0.0;
{
@ -129,8 +129,8 @@ void compute_sigma_gauss( config_file& cf, transfer_function *ptf, double R, std
double D0 = ccalc.CalcGrowthFactor(1.0);
double sigma8 = cf.getValue<double>("cosmology","sigma_8");
double nspec = cf.getValue<double>("cosmology","nspec");
double sigma8 = cf.get_value<double>("cosmology","sigma_8");
double nspec = cf.get_value<double>("cosmology","nspec");
double sigma0 = 0.0;
{
@ -170,14 +170,14 @@ constraint_set::constraint_set( config_file& cf, transfer_function *ptf )
unsigned i=0;
double astart = 1.0/(1.0+pcf_->getValue<double>("setup","zstart"));
unsigned levelmin = pcf_->getValue<unsigned>("setup","levelmin");
unsigned levelmin_TF = pcf_->getValueSafe<unsigned>("setup","levelmin_TF",levelmin);
constr_level_ = pcf_->getValueSafe<unsigned>("constraints","level",levelmin_TF);
double astart = 1.0/(1.0+pcf_->get_value<double>("setup","zstart"));
unsigned levelmin = pcf_->get_value<unsigned>("setup","levelmin");
unsigned levelmin_TF = pcf_->get_value_safe<unsigned>("setup","levelmin_TF",levelmin);
constr_level_ = pcf_->get_value_safe<unsigned>("constraints","level",levelmin_TF);
constr_level_ = std::max(constr_level_,levelmin_TF);
double omegam = pcf_->getValue<double>("cosmology","Omega_m");
double omegam = pcf_->get_value<double>("cosmology","Omega_m");
double rhom = omegam*2.77519737e11; //... mean matter density in Msun/Mpc^3
//... use EdS density for estimation
@ -192,9 +192,9 @@ constraint_set::constraint_set( config_file& cf, transfer_function *ptf )
char temp1[128];
std::string temp2;
sprintf(temp1,"constraint[%u].type",i);
if( cf.containsKey( "constraints", temp1 ) )
if( cf.contains_key( "constraints", temp1 ) )
{
std::string str_type = cf.getValue<std::string>( "constraints", temp1 );
std::string str_type = cf.get_value<std::string>( "constraints", temp1 );
if( constr_type_map.find(str_type) == constr_type_map.end() )
throw std::runtime_error("Unknown constraint type!\n");
@ -205,51 +205,51 @@ constraint_set::constraint_set( config_file& cf, transfer_function *ptf )
//... read position of constraint
sprintf(temp1,"constraint[%u].pos",i);
temp2 = cf.getValue<std::string>( "constraints", temp1 );
temp2 = cf.get_value<std::string>( "constraints", temp1 );
sscanf(temp2.c_str(), "%lf,%lf,%lf", &new_c.x, &new_c.y, &new_c.z);
if( new_c.type == halo)
{
//.. halo type constraints take mass and collapse redshift
sprintf(temp1,"constraint[%u].mass",i);
double mass = cf.getValue<double>( "constraints", temp1 );
double mass = cf.get_value<double>( "constraints", temp1 );
sprintf(temp1,"constraint[%u].zform",i);
double zcoll = cf.getValue<double>( "constraints", temp1 );
double zcoll = cf.get_value<double>( "constraints", temp1 );
new_c.Rg = pow((mass/pow(2.*M_PI,1.5)/rhom),1./3.);
new_c.sigma = 1.686/(pccalc_->CalcGrowthFactor(1./(1.+zcoll))/pccalc_->CalcGrowthFactor(1.0));
LOGINFO("sigma of constraint : %g", new_c.sigma );
music::ilog.Print("sigma of constraint : %g", new_c.sigma );
new_c.sigma *=pccalc_->CalcGrowthFactor(astart)/pccalc_->CalcGrowthFactor(1.0);
LOGINFO("Constraint %d : halo with %g h-1 M_o",i,pow(2.*M_PI,1.5)*rhom*pow(new_c.Rg,3));
music::ilog.Print("Constraint %d : halo with %g h-1 M_o",i,pow(2.*M_PI,1.5)*rhom*pow(new_c.Rg,3));
}
else if( new_c.type == peak )
{
//... peak type constraints take a scale and a peak height
//sprintf(temp1,"constraint[%u].Rg",i);
//new_c.Rg = cf.getValue<double>( "constraints", temp1 );
//new_c.Rg = cf.get_value<double>( "constraints", temp1 );
//double mass = pow(new_c.Rg,3.0)*rhom*pow(2.*M_PI,1.5);
sprintf(temp1,"constraint[%u].mass",i);
double mass = cf.getValue<double>( "constraints", temp1 );
double mass = cf.get_value<double>( "constraints", temp1 );
new_c.Rg = pow((mass/pow(2.*M_PI,1.5)/rhom),1./3.);
double Rtophat = pow(mass/4.0*3.0/M_PI/rhom,1./3.);
sprintf(temp1,"constraint[%u].nu",i);
double nu = cf.getValue<double>( "constraints", temp1 );
double nu = cf.get_value<double>( "constraints", temp1 );
std::vector<double> z,sigma;
compute_sigma_tophat( cf, ptf, Rtophat, z, sigma );
double zcoll = find_coll_z( z, sigma, nu );
//LOGINFO("Probable collapse redshift for constraint %d : z = %f @ M = %g", i, zcoll,mass );
//music::ilog.Print("Probable collapse redshift for constraint %d : z = %f @ M = %g", i, zcoll,mass );
compute_sigma_gauss( cf, ptf, new_c.Rg, z, sigma );
new_c.sigma = nu*sigma.back();
//LOGINFO("Constraint %d : peak with Rg=%g h-1 Mpc and nu = %g",i,new_c.Rg,new_c.sigma);
LOGINFO("Constraint %3d : peak",i);
LOGINFO(" M = %g h-1 M_o, nu = %.2f sigma", mass, nu );
LOGINFO(" estimated z_coll = %f, sigma = %f", zcoll, new_c.sigma );
//music::ilog.Print("Constraint %d : peak with Rg=%g h-1 Mpc and nu = %g",i,new_c.Rg,new_c.sigma);
music::ilog.Print("Constraint %3d : peak",i);
music::ilog.Print(" M = %g h-1 M_o, nu = %.2f sigma", mass, nu );
music::ilog.Print(" estimated z_coll = %f, sigma = %f", zcoll, new_c.sigma );
}
@ -263,7 +263,7 @@ constraint_set::constraint_set( config_file& cf, transfer_function *ptf )
++i;
}
LOGINFO("Found %d density constraint(s) to be obeyed.",cset_.size());
music::ilog.Print("Found %d density constraint(s) to be obeyed.",cset_.size());
}
@ -272,8 +272,8 @@ void constraint_set::wnoise_constr_corr( double dx, size_t nx, size_t ny, size_t
double lsub = nx*dx;
double dk = 2.0*M_PI/lsub, d3k=dk*dk*dk;
double pnorm = pcf_->getValue<double>("cosmology","pnorm");
double nspec = pcf_->getValue<double>("cosmology","nspec");
double pnorm = pcf_->get_value<double>("cosmology","pnorm");
double nspec = pcf_->get_value<double>("cosmology","nspec");
pnorm *= dplus0_*dplus0_;
size_t nconstr = cset_.size();
@ -299,7 +299,7 @@ void constraint_set::wnoise_constr_corr( double dx, size_t nx, size_t ny, size_t
chisq += cset_[i].sigma*cinv(i,j)*cset_[j].sigma;
chisq0 += g0[i]*cinv(i,j)*g0[j];
}
LOGINFO("Chi squared for the constraints:\n sampled = %f, desired = %f", chisq0, chisq );
music::ilog.Print("Chi squared for the constraints:\n sampled = %f, desired = %f", chisq0, chisq );
std::vector<double> sigma(nconstr,0.0);
@ -369,7 +369,7 @@ void constraint_set::wnoise_constr_corr( double dx, size_t nx, size_t ny, size_t
}
for(int i=0; i<(int)nconstr; ++i )
LOGINFO("Constraint %3d : sigma = %+6f (%+6f)",i,sigma[i],cset_[i].sigma);
music::ilog.Print("Constraint %3d : sigma = %+6f (%+6f)",i,sigma[i],cset_[i].sigma);
}
@ -381,8 +381,8 @@ void constraint_set::wnoise_constr_corr( double dx, fftw_complex* cw, size_t nx,
g0.assign(nconstr,0.0);
double pnorm = pcf_->getValue<double>("cosmology","pnorm");
double nspec = pcf_->getValue<double>("cosmology","nspec");
double pnorm = pcf_->get_value<double>("cosmology","pnorm");
double nspec = pcf_->get_value<double>("cosmology","nspec");
pnorm *= dplus0_*dplus0_;
double lsub = nx*dx;
double dk = 2.0*M_PI/lsub, d3k=dk*dk*dk;
@ -439,8 +439,8 @@ void constraint_set::icov_constr( double dx, size_t nx, size_t ny, size_t nz, ma
size_t nconstr = cset_.size();
size_t nzp=nz/2+1;
double pnorm = pcf_->getValue<double>("cosmology","pnorm");
double nspec = pcf_->getValue<double>("cosmology","nspec");
double pnorm = pcf_->get_value<double>("cosmology","pnorm");
double nspec = pcf_->get_value<double>("cosmology","nspec");
pnorm *= dplus0_*dplus0_;
cij = matrix(nconstr,nconstr);

16
src/constraints.hh
View file

@ -179,7 +179,7 @@ public:
return;
unsigned nlvl = 1<<ilevel;
double boxlength = pcf_->getValue<double>("setup","boxlength");
double boxlength = pcf_->get_value<double>("setup","boxlength");
//... compute constraint coordinates for grid
for( size_t i=0; i<cset_.size(); ++i )
@ -191,20 +191,20 @@ public:
cset_[i].gRg2 = cset_[i].gRg*cset_[i].gRg;
if(cset_[i].gRg > 0.5*lx[0])
LOGWARN("Constraint %d appears to be too large scale",i);
music::wlog.Print("Constraint %d appears to be too large scale",i);
}
std::vector<double> g0;
// unsigned levelmax = pcf_->getValue<unsigned>("setup","levelmax");
unsigned levelmin = pcf_->getValue<unsigned>("setup","levelmin_TF");
// unsigned levelmax = pcf_->get_value<unsigned>("setup","levelmax");
unsigned levelmin = pcf_->get_value<unsigned>("setup","levelmin_TF");
bool bperiodic = ilevel==levelmin;
double dx = pcf_->getValue<double>("setup","boxlength")/(1<<ilevel);
double dx = pcf_->get_value<double>("setup","boxlength")/(1<<ilevel);
LOGINFO("Computing constrained realization...");
music::ilog.Print("Computing constrained realization...");
if( bperiodic )
{
@ -284,7 +284,7 @@ public:
(*wnoise)((x0[0]+i),(x0[1]+j),(x0[2]+k)) = w[q]*fftnorm;
}
LOGINFO("Applied constraints to level %d.",ilevel);
music::ilog.Print("Applied constraints to level %d.",ilevel);
delete[] w;
@ -397,7 +397,7 @@ public:
}
LOGINFO("Applied constraints to level %d.",ilevel);
music::ilog.Print("Applied constraints to level %d.",ilevel);
delete[] w;

20
src/convolution_kernel.cc
View file

@ -45,8 +45,8 @@ void perform(kernel *pk, void *pd, bool shift, bool fix, bool flip)
std::cout << " - Performing density convolution... ("
<< cparam_.nx << ", " << cparam_.ny << ", " << cparam_.nz << ")\n";
LOGUSER("Performing kernel convolution on (%5d,%5d,%5d) grid", cparam_.nx, cparam_.ny, cparam_.nz);
LOGUSER("Performing forward FFT...");
music::ulog.Print("Performing kernel convolution on (%5d,%5d,%5d) grid", cparam_.nx, cparam_.ny, cparam_.nz);
music::ulog.Print("Performing forward FFT...");
#ifdef FFTW3
#ifdef SINGLE_PRECISION
fftwf_plan plan, iplan;
@ -82,13 +82,13 @@ void perform(kernel *pk, void *pd, bool shift, bool fix, bool flip)
if (shift)
{
double boxlength = pk->pcf_->getValue<double>("setup", "boxlength");
double stagfact = pk->pcf_->getValueSafe<double>("setup", "baryon_staggering", 0.5);
int lmax = pk->pcf_->getValue<int>("setup", "levelmax");
double boxlength = pk->pcf_->get_value<double>("setup", "boxlength");
double stagfact = pk->pcf_->get_value_safe<double>("setup", "baryon_staggering", 0.5);
int lmax = pk->pcf_->get_value<int>("setup", "levelmax");
double dxmax = boxlength / (1 << lmax);
double dxcur = cparam_.lx / cparam_.nx;
//std::cerr << "Performing staggering shift for SPH\n";
LOGUSER("Performing staggering shift for SPH");
music::ulog.Print("Performing staggering shift for SPH");
dstag = stagfact * 2.0 * M_PI / cparam_.nx * dxmax / dxcur;
}
@ -161,7 +161,7 @@ void perform(kernel *pk, void *pd, bool shift, bool fix, bool flip)
IM(cdata[0]) = 0.0;
LOGUSER("Performing backward FFT...");
music::ulog.Print("Performing backward FFT...");
#ifdef FFTW3
#ifdef SINGLE_PRECISION
@ -215,9 +215,9 @@ public:
kernel_k(config_file &cf, transfer_function *ptf, refinement_hierarchy &refh, tf_type type)
: kernel(cf, ptf, refh, type)
{
boxlength_ = pcf_->getValue<double>("setup", "boxlength");
nspec_ = pcf_->getValue<double>("cosmology", "nspec");
pnorm_ = pcf_->getValue<double>("cosmology", "pnorm");
boxlength_ = pcf_->get_value<double>("setup", "boxlength");
nspec_ = pcf_->get_value<double>("cosmology", "nspec");
pnorm_ = pcf_->get_value<double>("cosmology", "pnorm");
volfac_ = 1.0; //pow(boxlength,3)/pow(2.0*M_PI,3);
kfac_ = 2.0 * M_PI / boxlength_;
kmax_ = kfac_ / 2;

58
src/densities.cc
View file

@ -184,7 +184,7 @@ void fft_interpolate(m1 &V, m2 &v, bool from_basegrid = false)
ozf -= mzf/2; //nzf / 8;
}
LOGUSER("FFT interpolate: offset=%d,%d,%d size=%d,%d,%d", oxf, oyf, ozf, nxf, nyf, nzf);
music::ulog.Print("FFT interpolate: offset=%d,%d,%d size=%d,%d,%d", oxf, oyf, ozf, nxf, nyf, nzf);
// cut out piece of coarse grid that overlaps the fine:
assert(nxf % 2 == 0 && nyf % 2 == 0 && nzf % 2 == 0);
@ -345,19 +345,19 @@ void GenerateDensityUnigrid(config_file &cf, transfer_function *ptf, tf_type typ
{
unsigned levelmin, levelmax, levelminPoisson;
levelminPoisson = cf.getValue<unsigned>("setup", "levelmin");
levelmin = cf.getValueSafe<unsigned>("setup", "levelmin_TF", levelminPoisson);
levelmax = cf.getValue<unsigned>("setup", "levelmax");
levelminPoisson = cf.get_value<unsigned>("setup", "levelmin");
levelmin = cf.get_value_safe<unsigned>("setup", "levelmin_TF", levelminPoisson);
levelmax = cf.get_value<unsigned>("setup", "levelmax");
bool kspace = cf.getValue<bool>("setup", "kspace_TF");
bool kspace = cf.get_value<bool>("setup", "kspace_TF");
bool fix = cf.getValueSafe<bool>("setup","fix_mode_amplitude",false);
bool flip = cf.getValueSafe<bool>("setup","flip_mode_amplitude",false);
bool fix = cf.get_value_safe<bool>("setup","fix_mode_amplitude",false);
bool flip = cf.get_value_safe<bool>("setup","flip_mode_amplitude",false);
unsigned nbase = 1 << levelmin;
std::cerr << " - Running unigrid version\n";
LOGUSER("Running unigrid density convolution...");
music::ulog.Print("Running unigrid density convolution...");
//... select the transfer function to be used
convolution::kernel_creator *the_kernel_creator;
@ -365,7 +365,7 @@ void GenerateDensityUnigrid(config_file &cf, transfer_function *ptf, tf_type typ
if (kspace)
{
std::cout << " - Using k-space transfer function kernel.\n";
LOGUSER("Using k-space transfer function kernel.");
music::ulog.Print("Using k-space transfer function kernel.");
#ifdef SINGLE_PRECISION
the_kernel_creator = convolution::get_kernel_map()["tf_kernel_k_float"];
@ -376,7 +376,7 @@ void GenerateDensityUnigrid(config_file &cf, transfer_function *ptf, tf_type typ
else
{
std::cout << " - Using real-space transfer function kernel.\n";
LOGUSER("Using real-space transfer function kernel.");
music::ulog.Print("Using real-space transfer function kernel.");
#ifdef SINGLE_PRECISION
the_kernel_creator = convolution::get_kernel_map()["tf_kernel_real_float"];
@ -390,7 +390,7 @@ void GenerateDensityUnigrid(config_file &cf, transfer_function *ptf, tf_type typ
//...
std::cout << " - Performing noise convolution on level " << std::setw(2) << levelmax << " ..." << std::endl;
LOGUSER("Performing noise convolution on level %3d", levelmax);
music::ulog.Print("Performing noise convolution on level %3d", levelmax);
//... create convolution mesh
DensityGrid<real_t> *top = new DensityGrid<real_t>(nbase, nbase, nbase);
@ -428,7 +428,6 @@ void GenerateDensityHierarchy(config_file &cf, transfer_function *ptf, tf_type t
unsigned levelmin, levelmax, levelminPoisson;
std::vector<long> rngseeds;
std::vector<std::string> rngfnames;
bool kspaceTF;
double tstart, tend;
@ -438,17 +437,16 @@ void GenerateDensityHierarchy(config_file &cf, transfer_function *ptf, tf_type t
tstart = (double)clock() / CLOCKS_PER_SEC;
#endif
levelminPoisson = cf.getValue<unsigned>("setup", "levelmin");
levelmin = cf.getValueSafe<unsigned>("setup", "levelmin_TF", levelminPoisson);
levelmax = cf.getValue<unsigned>("setup", "levelmax");
kspaceTF = cf.getValue<bool>("setup", "kspace_TF");
levelminPoisson = cf.get_value<unsigned>("setup", "levelmin");
levelmin = cf.get_value_safe<unsigned>("setup", "levelmin_TF", levelminPoisson);
levelmax = cf.get_value<unsigned>("setup", "levelmax");
bool fix = cf.getValueSafe<bool>("setup","fix_mode_amplitude",false);
bool flip = cf.getValueSafe<bool>("setup","flip_mode_amplitude",false);
bool fourier_splicing = cf.getValueSafe<bool>("setup","fourier_splicing",true);
bool fix = cf.get_value_safe<bool>("setup","fix_mode_amplitude",false);
bool flip = cf.get_value_safe<bool>("setup","flip_mode_amplitude",false);
bool fourier_splicing = cf.get_value_safe<bool>("setup","fourier_splicing",true);
if( fix && levelmin != levelmax ){
LOGWARN("You have chosen mode fixing for a zoom. This is not well tested,\n please proceed at your own risk...");
music::wlog.Print("You have chosen mode fixing for a zoom. This is not well tested,\n please proceed at your own risk...");
}
unsigned nbase = 1 << levelmin;
@ -456,7 +454,7 @@ void GenerateDensityHierarchy(config_file &cf, transfer_function *ptf, tf_type t
convolution::kernel_creator *the_kernel_creator;
std::cout << " - Using k-space transfer function kernel.\n";
LOGUSER("Using k-space transfer function kernel.");
music::ulog.Print("Using k-space transfer function kernel.");
#ifdef SINGLE_PRECISION
the_kernel_creator = convolution::get_kernel_map()["tf_kernel_k_float"];
@ -476,7 +474,7 @@ void GenerateDensityHierarchy(config_file &cf, transfer_function *ptf, tf_type t
// do coarse level
top = new DensityGrid<real_t>(nbase, nbase, nbase);
LOGINFO("Performing noise convolution on level %3d", levelmin);
music::ilog.Print("Performing noise convolution on level %3d", levelmin);
rand.load(*top, levelmin);
convolution::perform<real_t>(the_tf_kernel->fetch_kernel(levelmin, false), reinterpret_cast<void *>(top->get_data_ptr()), shift, fix, flip);
@ -485,24 +483,24 @@ void GenerateDensityHierarchy(config_file &cf, transfer_function *ptf, tf_type t
for (int i = 1; i < nlevels; ++i)
{
LOGINFO("Performing noise convolution on level %3d...", levelmin + i);
music::ilog.Print("Performing noise convolution on level %3d...", levelmin + i);
/////////////////////////////////////////////////////////////////////////
//... add new refinement patch
LOGUSER("Allocating refinement patch");
LOGUSER(" offset=(%5d,%5d,%5d)", refh.offset(levelmin + i, 0),
music::ulog.Print("Allocating refinement patch");
music::ulog.Print(" offset=(%5d,%5d,%5d)", refh.offset(levelmin + i, 0),
refh.offset(levelmin + i, 1), refh.offset(levelmin + i, 2));
LOGUSER(" size =(%5d,%5d,%5d)", refh.size(levelmin + i, 0),
music::ulog.Print(" size =(%5d,%5d,%5d)", refh.size(levelmin + i, 0),
refh.size(levelmin + i, 1), refh.size(levelmin + i, 2));
if( refh.get_margin() > 0 ){
fine = new PaddedDensitySubGrid<real_t>( refh.offset(levelmin + i, 0), refh.offset(levelmin + i, 1), refh.offset(levelmin + i, 2),
refh.size(levelmin + i, 0), refh.size(levelmin + i, 1), refh.size(levelmin + i, 2),
refh.get_margin(), refh.get_margin(), refh.get_margin() );
LOGUSER(" margin = %d",refh.get_margin());
music::ulog.Print(" margin = %d",refh.get_margin());
}else{
fine = new PaddedDensitySubGrid<real_t>( refh.offset(levelmin + i, 0), refh.offset(levelmin + i, 1), refh.offset(levelmin + i, 2),
refh.size(levelmin + i, 0), refh.size(levelmin + i, 1), refh.size(levelmin + i, 2));
LOGUSER(" margin = %d",refh.size(levelmin + i, 0)/2);
music::ulog.Print(" margin = %d",refh.size(levelmin + i, 0)/2);
}
/////////////////////////////////////////////////////////////////////////
@ -554,7 +552,7 @@ void GenerateDensityHierarchy(config_file &cf, transfer_function *ptf, tf_type t
if( !fourier_splicing ){
coarsen_density(refh,delta,false);
}
LOGUSER("Finished computing the density field in %fs", tend - tstart);
music::ulog.Print("Finished computing the density field in %fs", tend - tstart);
}
/*******************************************************************************************/
@ -586,7 +584,7 @@ void normalize_density(grid_hierarchy &delta)
}
std::cout << " - Top grid mean density is off by " << sum << ", correcting..." << std::endl;
LOGUSER("Grid mean density is %g. Correcting...", sum);
music::ulog.Print("Grid mean density is %g. Correcting...", sum);
for (unsigned i = levelmin; i <= levelmax; ++i)
{

26
src/general.hh
View file

@ -11,7 +11,7 @@
#ifndef __GENERAL_HH
#define __GENERAL_HH
#include "log.hh"
#include "logger.hh"
#include <cassert>
#include "omp.h"
@ -125,23 +125,23 @@ typedef struct cosmology{
cosmology( config_file cf )
{
double zstart = cf.getValue<double>( "setup", "zstart" );
double zstart = cf.get_value<double>( "setup", "zstart" );
astart = 1.0/(1.0+zstart);
Omega_b = cf.getValue<double>( "cosmology", "Omega_b" );
Omega_m = cf.getValue<double>( "cosmology", "Omega_m" );
Omega_DE = cf.getValue<double>( "cosmology", "Omega_L" );
w_0 = cf.getValueSafe<double>( "cosmology", "w0", -1.0 );
w_a = cf.getValueSafe<double>( "cosmology", "wa", 0.0 );
Omega_b = cf.get_value<double>( "cosmology", "Omega_b" );
Omega_m = cf.get_value<double>( "cosmology", "Omega_m" );
Omega_DE = cf.get_value<double>( "cosmology", "Omega_L" );
w_0 = cf.get_value_safe<double>( "cosmology", "w0", -1.0 );
w_a = cf.get_value_safe<double>( "cosmology", "wa", 0.0 );
Omega_r = cf.getValueSafe<double>( "cosmology", "Omega_r", 0.0 ); // no longer default to nonzero (8.3e-5)
Omega_r = cf.get_value_safe<double>( "cosmology", "Omega_r", 0.0 ); // no longer default to nonzero (8.3e-5)
Omega_k = 1.0 - Omega_m - Omega_DE - Omega_r;
H0 = cf.getValue<double>( "cosmology", "H0" );
sigma8 = cf.getValue<double>( "cosmology", "sigma_8" );
nspect = cf.getValue<double>( "cosmology", "nspec" );
WDMg_x = cf.getValueSafe<double>( "cosmology", "WDMg_x", 1.5 );
WDMmass = cf.getValueSafe<double>( "cosmology", "WDMmass", 0.0 );
H0 = cf.get_value<double>( "cosmology", "H0" );
sigma8 = cf.get_value<double>( "cosmology", "sigma_8" );
nspect = cf.get_value<double>( "cosmology", "nspec" );
WDMg_x = cf.get_value_safe<double>( "cosmology", "WDMg_x", 1.5 );
WDMmass = cf.get_value_safe<double>( "cosmology", "WDMmass", 0.0 );
dplus = 0.0;
pnorm = 0.0;

124
src/log.cc
View file

@ -1,124 +0,0 @@
/*
log.cc - This file is part of MUSIC -
a code to generate multi-scale initial conditions
for cosmological simulations
Copyright (C) 2010 Oliver Hahn
*/
#include "log.hh"
#include <iostream>
#include <algorithm>
std::string RemoveMultipleWhiteSpaces( std::string s );
std::string MUSIC::log::outputFile_;
std::ofstream MUSIC::log::outputStream_;
std::list<MUSIC::log::message> MUSIC::log::messages_;
void (*MUSIC::log::receiver)(const message&) = NULL;
MUSIC::log::messageType MUSIC::log::logLevel_;
std::string RemoveMultipleWhiteSpaces( std::string s )
{
std::string search = " "; // this is 2 spaces
size_t index;
while( (index = s.find(search)) != std::string::npos )
{ // remove 1 character from the string at index
s.erase(index,1);
}
return s;
}
void MUSIC::log::send(messageType type, const std::string& text_)
//void MUSIC::log::send(messageType type, std::stringstream& textstr)
{
std::string text(text_);// = textstr.str();
// Skip logging if minimum level is higher
if (logLevel_)
if (type < logLevel_) return;
// log message
MUSIC::log::message m;
m.type = type;
m.text = text;
time_t t = time(NULL);
m.when = localtime(&t);
messages_.push_back(m);
if( type==Info||type==Warning||type==Error||type==FatalError )
{
std::cout << " - ";
if(type==Warning)
std::cout << "WARNING: ";
if(type==Error)
std::cout << "ERROR: ";
if(type==FatalError)
std::cout << "FATAL: ";
std::cout << text << std::endl;
}
std::replace(text.begin(),text.end(),'\n',' ');
RemoveMultipleWhiteSpaces(text);
// if enabled logging to file
if(outputStream_.is_open())
{
// print time
char buffer[9];
strftime(buffer, 9, "%X", m.when);
outputStream_ << buffer;
// print type
switch(type)
{
case Info: outputStream_ << " | info | "; break;
case DebugInfo: outputStream_ << " | debug | "; break;
case Warning: outputStream_ << " | warning | "; break;
case Error: outputStream_ << " | ERROR | "; break;
case FatalError:outputStream_ << " | FATAL | "; break;
case User: outputStream_ << " | info | "; break;
default: outputStream_ << " | ";
}
// print description
outputStream_ << text << std::endl;
}
// if user wants to catch messages, send it to him
if(receiver)
receiver(m);
}
void MUSIC::log::setOutput(const std::string& filename)
{
//logDebug("Setting output log file: " + filename);
outputFile_ = filename;
// close old one
if(outputStream_.is_open())
outputStream_.close();
// create file
outputStream_.open(filename.c_str());
if(!outputStream_.is_open())
LOGERR("Cannot create/open logfile \'%s\'.",filename.c_str());
}
void MUSIC::log::setLevel(const MUSIC::log::messageType level)
{
logLevel_ = level;
}
MUSIC::log::~log()
{
if(outputStream_.is_open())
outputStream_.close();
}

173
src/log.hh
View file

@ -1,173 +0,0 @@
/*
log.hh - This file is part of MUSIC -
a code to generate multi-scale initial conditions
for cosmological simulations
Copyright (C) 2010 Oliver Hahn
*/
#ifndef __LOG_HH
#define __LOG_HH
#include <string>
#include <list>
#include <fstream>
#include <ctime>
#include <cstdarg>
#include <sstream>
/*!
* \brief System for logging runtime library errors, warnings, etc.
*
* This is the class that catches every (debug) info, warning, error, or user message and
* processes it. Messages can be written to files and/or forwarded to user function for
* processing messages.
*/
namespace MUSIC
{
class log
{
public:
log(){}
~log();
/*!
* \brief Types of logged messages.
*/
enum messageType
{
Info,
DebugInfo,
Warning,
Error,
FatalError,
User
};
/*!
* \brief Logged message of type MessageType with some info.
*/
struct message
{
messageType type;
std::string text;
tm* when;
};
/*!
* \brief Open file where to log the messages.
*/
static void setOutput(const std::string& filename);
/*!
* \brief Get the filename of log.
*/
static const std::string& output() { return outputFile_; }
/*!
* \brief Add a new message to log.
* \param type Type of the new message.
* \param text Message.
* \remarks Message is directly passes to user reciever if one is set.
*/
static void send(messageType type, const std::string& text);
//static void send(messageType type, std::string& text);
/*!
* \brief Get the list of all of the logged messages.
*/
static const std::list<message>& messages() { return messages_; }
/*!
* \brief Get the last logged message.
*/
static const message& lastMessage() { return messages_.back(); }
/*!
* \brief Set user function to receive newly sent messages to logger.
*/
static void setUserReceiver(void (*userFunc)(const message&)) { receiver = userFunc; }
/*!
* \brief Set minimum level of message to be logged.
*/
static void setLevel(const log::messageType level);
private:
static std::string outputFile_;
static std::ofstream outputStream_;
static std::list<message> messages_;
static messageType logLevel_;
static void (*receiver)(const message&);
};
}
inline void LOGERR( const char* str, ... )
{
char out[1024];
va_list argptr;
va_start(argptr,str);
va_end(argptr);
vsprintf(out,str,argptr);
MUSIC::log::send(MUSIC::log::Error, std::string(out));
}
inline void LOGWARN( const char* str, ... )
{
char out[1024];
va_list argptr;
va_start(argptr,str);
va_end(argptr);
vsprintf(out,str,argptr);
MUSIC::log::send(MUSIC::log::Warning, std::string(out));
}
inline void LOGFATAL( const char* str, ... )
{
char out[1024];
va_list argptr;
va_start(argptr,str);
va_end(argptr);
vsprintf(out,str,argptr);
MUSIC::log::send(MUSIC::log::FatalError, std::string(out));
}
inline void LOGDEBUG( const char* str, ... )
{
char out[1024];
va_list argptr;
va_start(argptr,str);
va_end(argptr);
vsprintf(out,str,argptr);
MUSIC::log::send(MUSIC::log::DebugInfo, std::string(out));
}
inline void LOGUSER( const char* str, ... )
{
char out[1024];
va_list argptr;
va_start(argptr,str);
va_end(argptr);
vsprintf(out,str,argptr);
MUSIC::log::send(MUSIC::log::User, std::string(out));
}
inline void LOGINFO( const char* str, ... )
{
char out[1024];
va_list argptr;
va_start(argptr,str);
va_end(argptr);
vsprintf(out,str,argptr);
MUSIC::log::send(MUSIC::log::Info, std::string(out));
}
#endif //__LOG_HH

60
src/logger.cc Normal file
View file

@ -0,0 +1,60 @@
// This file is part of monofonIC (MUSIC2)
// A software package to generate ICs for cosmological simulations
// Copyright (C) 2020 by Oliver Hahn & Michael Michaux (this file)
//
// monofonIC 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.
//
// monofonIC 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 program. If not, see <http://www.gnu.org/licenses/>.
#include <logger.hh>
namespace music {
std::ofstream logger::output_file_;
log_level logger::log_level_ = log_level::off;
void logger::set_level(const log_level &level) {
log_level_ = level;
}
log_level logger::get_level() {
return log_level_;
}
void logger::set_output(const std::string filename) {
if (output_file_.is_open()) {
output_file_.close();
}
output_file_.open(filename, std::ofstream::out);
assert(output_file_.is_open());
}
void logger::unset_output() {
if (output_file_.is_open()) {
output_file_.close();
}
}
std::ofstream &logger::get_output() {
return output_file_;
}
// global instantiations for different levels
logger the_logger;
log_stream flog(the_logger, log_level::fatal);
log_stream elog(the_logger, log_level::error);
log_stream wlog(the_logger, log_level::warning);
log_stream ilog(the_logger, log_level::info);
log_stream ulog(the_logger, log_level::user);
log_stream dlog(the_logger, log_level::debug);
} // namespace music

154
src/logger.hh Normal file
View file

@ -0,0 +1,154 @@
// This file is part of monofonIC (MUSIC2)
// A software package to generate ICs for cosmological simulations
// Copyright (C) 2020 by Oliver Hahn & Michael Michaux (this file)
//
// monofonIC 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.
//
// monofonIC 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 program. If not, see <http://www.gnu.org/licenses/>.
#pragma once
#include <algorithm>
#include <cassert>
#include <cstdarg>
#include <fstream>
#include <iostream>
namespace music {
enum log_level : int {
off = 0,
fatal = 1,
error = 2,
warning = 3,
info = 4,
user = 5,
debug = 6,
};
class logger {
private:
static log_level log_level_;
static std::ofstream output_file_;
public:
logger() = default;
~logger() = default;
static void set_level(const log_level &level);
static log_level get_level();
static void set_output(const std::string filename);
static void unset_output();
static std::ofstream &get_output();
template <typename T> logger &operator<<(const T &item) {
std::cout << item;
if (output_file_.is_open()) {
output_file_ << item;
}
return *this;
}
logger &operator<<(std::ostream &(*fp)(std::ostream &)) {
std::cout << fp;
if (output_file_.is_open()) {
output_file_ << fp;
}
return *this;
}
};
class log_stream {
private:
logger &logger_;
log_level stream_level_;
std::string line_prefix_, line_postfix_;
bool newline;
public:
log_stream(logger &logger, const log_level &level)
: logger_(logger), stream_level_(level), newline(true) {
switch (stream_level_) {
case log_level::fatal:
line_prefix_ = "\033[31mFatal : ";
break;
case log_level::error:
line_prefix_ = "\033[31mError : ";
break;
case log_level::warning:
line_prefix_ = "\033[33mWarning : ";
break;
case log_level::info:
case log_level::user:
//line_prefix_ = " | Info | ";
line_prefix_ = " \033[0m";
break;
case log_level::debug:
line_prefix_ = "Debug : \033[0m";
break;
default:
line_prefix_ = "\033[0m";
break;
}
line_postfix_ = "\033[0m";
}
~log_stream() = default;
inline std::string GetPrefix() const {
return line_prefix_;
}
template <typename T> log_stream &operator<<(const T &item) {
if (logger::get_level() >= stream_level_) {
if (newline) {
logger_ << line_prefix_;
newline = false;
}
logger_ << item;
}
return *this;
}
log_stream &operator<<(std::ostream &(*fp)(std::ostream &)) {
if (logger::get_level() >= stream_level_) {
logger_ << fp;
logger_ << line_postfix_;
newline = true;
}
return *this;
}
inline void Print(const char *str, ...) {
char out[1024];
va_list argptr;
va_start(argptr, str);
vsprintf(out, str, argptr);
va_end(argptr);
std::string out_string = std::string(out);
out_string.erase(std::remove(out_string.begin(), out_string.end(), '\n'),
out_string.end());
(*this) << out_string << std::endl;
}
};
// global instantiations for different levels
extern logger glogger;
extern log_stream flog;
extern log_stream elog;
extern log_stream wlog;
extern log_stream ilog;
extern log_stream ulog;
extern log_stream dlog;
} // namespace music

242
src/main.cc
View file

@ -86,12 +86,12 @@ void splash(void)
<< " this is " << THE_CODE_NAME << " version " << THE_CODE_VERSION << "\n\n";
#if defined(CMAKE_BUILD)
LOGINFO("Version built from git rev.: %s, tag: %s, branch: %s", GIT_REV, GIT_TAG, GIT_BRANCH);
music::ilog.Print("Version built from git rev.: %s, tag: %s, branch: %s", GIT_REV, GIT_TAG, GIT_BRANCH);
#endif
#if defined(SINGLE_PRECISION)
LOGINFO("Version was compiled for single precision.");
music::ilog.Print("Version was compiled for single precision.");
#else
LOGINFO("Version was compiled for double precision.");
music::ilog.Print("Version was compiled for double precision.");
#endif
std::cout << "\n\n";
}
@ -100,10 +100,10 @@ void modify_grid_for_TF(const refinement_hierarchy &rh_full, refinement_hierarch
{
unsigned lbase, lbaseTF, lmax, overlap;
lbase = cf.getValue<unsigned>("setup", "levelmin");
lmax = cf.getValue<unsigned>("setup", "levelmax");
lbaseTF = cf.getValueSafe<unsigned>("setup", "levelmin_TF", lbase);
overlap = cf.getValueSafe<unsigned>("setup", "overlap", 4);
lbase = cf.get_value<unsigned>("setup", "levelmin");
lmax = cf.get_value<unsigned>("setup", "levelmax");
lbaseTF = cf.get_value_safe<unsigned>("setup", "levelmin_TF", lbase);
overlap = cf.get_value_safe<unsigned>("setup", "overlap", 4);
rh_TF = rh_full;
unsigned pad = overlap;
@ -151,10 +151,10 @@ void modify_grid_for_TF(const refinement_hierarchy &rh_full, refinement_hierarch
void print_hierarchy_stats(config_file &cf, const refinement_hierarchy &rh)
{
double omegam = cf.getValue<double>("cosmology", "Omega_m");
double omegab = cf.getValue<double>("cosmology", "Omega_b");
bool bbaryons = cf.getValue<bool>("setup", "baryons");
double boxlength = cf.getValue<double>("setup", "boxlength");
double omegam = cf.get_value<double>("cosmology", "Omega_m");
double omegab = cf.get_value<double>("cosmology", "Omega_b");
bool bbaryons = cf.get_value<bool>("setup", "baryons");
double boxlength = cf.get_value<double>("setup", "boxlength");
unsigned levelmin = rh.levelmin();
double dx = boxlength / (double)(1 << levelmin), dx3 = dx * dx * dx;
@ -221,11 +221,11 @@ void store_grid_structure(config_file &cf, const refinement_hierarchy &rh)
{
sprintf(str1, "offset(%d,%d)", i, j);
sprintf(str2, "%ld", rh.offset(i, j));
cf.insertValue("setup", str1, str2);
cf.insert_value("setup", str1, str2);
sprintf(str1, "size(%d,%d)", i, j);
sprintf(str2, "%ld", rh.size(i, j));
cf.insertValue("setup", str1, str2);
cf.insert_value("setup", str1, str2);
}
}
}
@ -306,6 +306,12 @@ int main(int argc, const char *argv[])
unsigned lbase, lmax, lbaseTF;
#if defined(NDEBUG)
music::logger::set_level(music::log_level::info);
#else
music::logger::set_level(music::log_level::debug);
#endif
//------------------------------------------------------------------------------
//... parse command line options
//------------------------------------------------------------------------------
@ -330,29 +336,29 @@ int main(int argc, const char *argv[])
char logfname[128];
sprintf(logfname, "%s_log.txt", argv[1]);
MUSIC::log::setOutput(logfname);
music::logger::set_output(logfname);
time_t ltime = time(NULL);
LOGINFO("Opening log file \'%s\'.", logfname);
LOGUSER("Running %s, version %s", THE_CODE_NAME, THE_CODE_VERSION);
LOGUSER("Log is for run started %s", asctime(localtime(&ltime)));
music::ilog.Print("Opening log file \'%s\'.", logfname);
music::ulog.Print("Running %s, version %s", THE_CODE_NAME, THE_CODE_VERSION);
music::ulog.Print("Log is for run started %s", asctime(localtime(&ltime)));
#ifdef FFTW3
LOGUSER("Code was compiled using FFTW version 3.x");
music::ulog.Print("Code was compiled using FFTW version 3.x");
#else
LOGUSER("Code was compiled using FFTW version 2.x");
music::ulog.Print("Code was compiled using FFTW version 2.x");
#endif
#ifdef SINGLETHREAD_FFTW
LOGUSER("Code was compiled for single-threaded FFTW");
music::ulog.Print("Code was compiled for single-threaded FFTW");
#else
LOGUSER("Code was compiled for multi-threaded FFTW");
LOGUSER("Running with a maximum of %d OpenMP threads", omp_get_max_threads());
music::ulog.Print("Code was compiled for multi-threaded FFTW");
music::ulog.Print("Running with a maximum of %d OpenMP threads", omp_get_max_threads());
#endif
#ifdef SINGLE_PRECISION
LOGUSER("Code was compiled for single precision.");
music::ulog.Print("Code was compiled for single precision.");
#else
LOGUSER("Code was compiled for double precision.");
music::ulog.Print("Code was compiled for double precision.");
#endif
//------------------------------------------------------------------------------
@ -367,35 +373,35 @@ int main(int argc, const char *argv[])
//... initialize some parameters about grid set-up
//------------------------------------------------------------------------------
boxlength = cf.getValue<double>("setup", "boxlength");
lbase = cf.getValue<unsigned>("setup", "levelmin");
lmax = cf.getValue<unsigned>("setup", "levelmax");
lbaseTF = cf.getValueSafe<unsigned>("setup", "levelmin_TF", lbase);
boxlength = cf.get_value<double>("setup", "boxlength");
lbase = cf.get_value<unsigned>("setup", "levelmin");
lmax = cf.get_value<unsigned>("setup", "levelmax");
lbaseTF = cf.get_value_safe<unsigned>("setup", "levelmin_TF", lbase);
if (lbase == lmax && !force_shift)
cf.insertValue("setup", "no_shift", "yes");
cf.insert_value("setup", "no_shift", "yes");
if (lbaseTF < lbase)
{
std::cout << " - WARNING: levelminTF < levelmin. This is not good!\n"
<< " I will set levelminTF = levelmin.\n";
LOGUSER("levelminTF < levelmin. set levelminTF = levelmin.");
music::ulog.Print("levelminTF < levelmin. set levelminTF = levelmin.");
lbaseTF = lbase;
cf.insertValue("setup", "levelmin_TF", cf.getValue<std::string>("setup", "levelmin"));
cf.insert_value("setup", "levelmin_TF", cf.get_value<std::string>("setup", "levelmin"));
}
// .. determine if spectral sampling should be used
if (!cf.containsKey("setup", "kspace_TF"))
cf.insertValue("setup", "kspace_TF", "yes");
if (!cf.contains_key("setup", "kspace_TF"))
cf.insert_value("setup", "kspace_TF", "yes");
bool bspectral_sampling = cf.getValue<bool>("setup", "kspace_TF");
bool bspectral_sampling = cf.get_value<bool>("setup", "kspace_TF");
if (bspectral_sampling)
LOGINFO("Using k-space sampled transfer functions...");
music::ilog.Print("Using k-space sampled transfer functions...");
else
LOGINFO("Using real space sampled transfer functions...");
music::ilog.Print("Using real space sampled transfer functions...");
//------------------------------------------------------------------------------
//... initialize multithread FFTW
@ -419,10 +425,10 @@ int main(int argc, const char *argv[])
//... initialize cosmology
//------------------------------------------------------------------------------
bool
do_baryons = cf.getValue<bool>("setup", "baryons"),
do_2LPT = cf.getValueSafe<bool>("setup", "use_2LPT", false),
do_LLA = cf.getValueSafe<bool>("setup", "use_LLA", false),
do_counter_mode = cf.getValueSafe<bool>("setup", "zero_zoom_velocity", false);
do_baryons = cf.get_value<bool>("setup", "baryons"),
do_2LPT = cf.get_value_safe<bool>("setup", "use_2LPT", false),
do_LLA = cf.get_value_safe<bool>("setup", "use_LLA", false),
do_counter_mode = cf.get_value_safe<bool>("setup", "zero_zoom_velocity", false);
transfer_function_plugin *the_transfer_function_plugin = select_transfer_function_plugin(cf);
@ -439,7 +445,7 @@ int main(int argc, const char *argv[])
cosmo.pnorm *= cosmo.dplus * cosmo.dplus;
//... directly use the normalisation via a parameter rather than the calculated one
cosmo.pnorm = cf.getValueSafe<double>("setup", "force_pnorm", cosmo.pnorm);
cosmo.pnorm = cf.get_value_safe<double>("setup", "force_pnorm", cosmo.pnorm);
double vfac2lpt = 1.0;
@ -453,11 +459,11 @@ int main(int argc, const char *argv[])
{
char tmpstr[128];
sprintf(tmpstr, "%.12g", cosmo.pnorm);
cf.insertValue("cosmology", "pnorm", tmpstr);
cf.insert_value("cosmology", "pnorm", tmpstr);
sprintf(tmpstr, "%.12g", cosmo.dplus);
cf.insertValue("cosmology", "dplus", tmpstr);
cf.insert_value("cosmology", "dplus", tmpstr);
sprintf(tmpstr, "%.12g", cosmo.vfact);
cf.insertValue("cosmology", "vfact", tmpstr);
cf.insert_value("cosmology", "vfact", tmpstr);
}
the_region_generator = select_region_generator_plugin(cf);
@ -490,17 +496,17 @@ int main(int argc, const char *argv[])
modify_grid_for_TF(rh_Poisson, rh_TF, cf);
// rh_TF.output();
LOGUSER("Grid structure for Poisson solver:");
music::ulog.Print("Grid structure for Poisson solver:");
rh_Poisson.output_log();
LOGUSER("Grid structure for density convolution:");
music::ulog.Print("Grid structure for density convolution:");
rh_TF.output_log();
//------------------------------------------------------------------------------
//... initialize the output plug-in
//------------------------------------------------------------------------------
std::string outformat, outfname;
outformat = cf.getValue<std::string>("output", "format");
outfname = cf.getValue<std::string>("output", "filename");
outformat = cf.get_value<std::string>("output", "format");
outfname = cf.get_value<std::string>("output", "filename");
output_plugin *the_output_plugin = select_output_plugin(cf);
//------------------------------------------------------------------------------
@ -509,24 +515,24 @@ int main(int argc, const char *argv[])
std::cout << "=============================================================\n";
std::cout << " GENERATING WHITE NOISE\n";
std::cout << "-------------------------------------------------------------\n";
LOGUSER("Computing white noise...");
music::ulog.Print("Computing white noise...");
// rand_gen rand(cf, rh_TF, the_transfer_function_plugin);
rand.initialize_for_grid_structure( rh_TF );
//------------------------------------------------------------------------------
//... initialize the Poisson solver
//------------------------------------------------------------------------------
// bool bdefd = cf.getValueSafe<bool> ( "poisson" , "fft_fine", true );
// bool bdefd = cf.get_value_safe<bool> ( "poisson" , "fft_fine", true );
bool bdefd = true; // we set this by default and don't allow it to be changed outside any more
bool bglass = cf.getValueSafe<bool>("output", "glass", false);
bool bsph = cf.getValueSafe<bool>("setup", "do_SPH", false) && do_baryons;
bool bglass = cf.get_value_safe<bool>("output", "glass", false);
bool bsph = cf.get_value_safe<bool>("setup", "do_SPH", false) && do_baryons;
bool bbshift = bsph && !bglass;
bool kspace = cf.getValueSafe<bool>("poisson", "kspace", false);
bool kspace = cf.get_value_safe<bool>("poisson", "kspace", false);
bool kspace2LPT = kspace;
bool decic_DM = cf.getValueSafe<bool>("output", "glass_cicdeconvolve", false);
bool decic_baryons = cf.getValueSafe<bool>("output", "glass_cicdeconvolve", false) & bsph;
bool decic_DM = cf.get_value_safe<bool>("output", "glass_cicdeconvolve", false);
bool decic_baryons = cf.get_value_safe<bool>("output", "glass_cicdeconvolve", false) & bsph;
std::array<double,3> counter_mode_amp;
@ -544,7 +550,7 @@ int main(int argc, const char *argv[])
else
poisson_solver_name = std::string("mg_poisson");
unsigned grad_order = cf.getValueSafe<unsigned>("poisson", "grad_order", 4);
unsigned grad_order = cf.get_value_safe<unsigned>("poisson", "grad_order", 4);
//... switch off if using kspace anyway
// bdefd &= !kspace;
@ -560,7 +566,7 @@ int main(int argc, const char *argv[])
{
if (!do_2LPT)
{
LOGUSER("Entering 1LPT branch");
music::ulog.Print("Entering 1LPT branch");
//------------------------------------------------------------------------------
//... cdm density and displacements
@ -568,7 +574,7 @@ int main(int argc, const char *argv[])
std::cout << "=============================================================\n";
std::cout << " COMPUTING DARK MATTER DISPLACEMENTS\n";
std::cout << "-------------------------------------------------------------\n";
LOGUSER("Computing dark matter displacements...");
music::ulog.Print("Computing dark matter displacements...");
grid_hierarchy f(nbnd); //, u(nbnd);
tf_type my_tf_type = cdm;
@ -581,7 +587,7 @@ int main(int argc, const char *argv[])
normalize_density(f);
LOGUSER("Writing CDM data");
music::ulog.Print("Writing CDM data");
the_output_plugin->write_dm_mass(f);
the_output_plugin->write_dm_density(f);
@ -592,7 +598,7 @@ int main(int argc, const char *argv[])
if (!bdefd)
f.deallocate();
LOGUSER("Writing CDM potential");
music::ulog.Print("Writing CDM potential");
the_output_plugin->write_dm_potential(u);
//------------------------------------------------------------------------------
@ -615,14 +621,14 @@ int main(int argc, const char *argv[])
//... displacement
the_poisson_solver->gradient(icoord, u, data_forIO);
double dispmax = compute_finest_absmax(data_forIO);
LOGINFO("max. %c-displacement of HR particles is %f [mean dx]", 'x' + icoord, dispmax * (double)(1ll << data_forIO.levelmax()));
music::ilog.Print("max. %c-displacement of HR particles is %f [mean dx]", 'x' + icoord, dispmax * (double)(1ll << data_forIO.levelmax()));
coarsen_density(rh_Poisson, data_forIO, false);
//... compute counter-mode to minimize advection errors
counter_mode_amp[icoord] = compute_finest_mean(data_forIO);
if( do_counter_mode ) add_constant_value( data_forIO, -counter_mode_amp[icoord] );
LOGUSER("Writing CDM displacements");
music::ulog.Print("Writing CDM displacements");
the_output_plugin->write_dm_position(icoord, data_forIO);
}
if (do_baryons)
@ -638,7 +644,7 @@ int main(int argc, const char *argv[])
std::cout << "=============================================================\n";
std::cout << " COMPUTING BARYON DENSITY\n";
std::cout << "-------------------------------------------------------------\n";
LOGUSER("Computing baryon density...");
music::ulog.Print("Computing baryon density...");
GenerateDensityHierarchy(cf, the_transfer_function_plugin, baryon, rh_TF, rand, f, false, bbshift);
coarsen_density(rh_Poisson, f, bspectral_sampling);
f.add_refinement_mask(rh_Poisson.get_coord_shift());
@ -646,7 +652,7 @@ int main(int argc, const char *argv[])
if (!do_LLA)
{
LOGUSER("Writing baryon density");
music::ulog.Print("Writing baryon density");
the_output_plugin->write_gas_density(f);
}
@ -676,7 +682,7 @@ int main(int argc, const char *argv[])
the_poisson_solver->gradient(icoord, u, data_forIO);
coarsen_density(rh_Poisson, data_forIO, false);
LOGUSER("Writing baryon displacements");
music::ulog.Print("Writing baryon displacements");
the_output_plugin->write_gas_position(icoord, data_forIO);
}
u.deallocate();
@ -692,7 +698,7 @@ int main(int argc, const char *argv[])
compute_LLA_density(u, f, grad_order);
u.deallocate();
normalize_density(f);
LOGUSER("Writing baryon density");
music::ulog.Print("Writing baryon density");
the_output_plugin->write_gas_density(f);
}
@ -707,11 +713,11 @@ int main(int argc, const char *argv[])
std::cout << "=============================================================\n";
std::cout << " COMPUTING VELOCITIES\n";
std::cout << "-------------------------------------------------------------\n";
LOGUSER("Computing velocitites...");
music::ulog.Print("Computing velocitites...");
if (do_baryons || the_transfer_function_plugin->tf_has_velocities())
{
LOGUSER("Generating velocity perturbations...");
music::ulog.Print("Generating velocity perturbations...");
GenerateDensityHierarchy(cf, the_transfer_function_plugin, vtotal, rh_TF, rand, f, false, false);
coarsen_density(rh_Poisson, f, bspectral_sampling);
f.add_refinement_mask(rh_Poisson.get_coord_shift());
@ -745,26 +751,26 @@ int main(int argc, const char *argv[])
//... velocity kick to keep refined region centered?
double sigv = compute_finest_sigma(data_forIO);
LOGINFO("sigma of %c-velocity of high-res particles is %f", 'x' + icoord, sigv);
music::ilog.Print("sigma of %c-velocity of high-res particles is %f", 'x' + icoord, sigv);
double meanv = compute_finest_mean(data_forIO);
LOGINFO("mean of %c-velocity of high-res particles is %f", 'x' + icoord, meanv);
LOGUSER("mean of %c-velocity of high-res particles is %f", 'x' + icoord, meanv);
music::ilog.Print("mean of %c-velocity of high-res particles is %f", 'x' + icoord, meanv);
music::ulog.Print("mean of %c-velocity of high-res particles is %f", 'x' + icoord, meanv);
double maxv = compute_finest_absmax(data_forIO);
LOGINFO("max of abs of %c-velocity of high-res particles is %f", 'x' + icoord, maxv);
music::ilog.Print("max of abs of %c-velocity of high-res particles is %f", 'x' + icoord, maxv);
coarsen_density(rh_Poisson, data_forIO, false);
// add counter velocity-mode
if( do_counter_mode ) add_constant_value( data_forIO, -counter_mode_amp[icoord]*cosmo.vfact );
LOGUSER("Writing CDM velocities");
music::ulog.Print("Writing CDM velocities");
the_output_plugin->write_dm_velocity(icoord, data_forIO);
if (do_baryons)
{
LOGUSER("Writing baryon velocities");
music::ulog.Print("Writing baryon velocities");
the_output_plugin->write_gas_velocity(icoord, data_forIO);
}
}
@ -774,11 +780,11 @@ int main(int argc, const char *argv[])
}
else
{
LOGINFO("Computing separate velocities for CDM and baryons:");
music::ilog.Print("Computing separate velocities for CDM and baryons:");
std::cout << "=============================================================\n";
std::cout << " COMPUTING DARK MATTER VELOCITIES\n";
std::cout << "-------------------------------------------------------------\n";
LOGUSER("Computing dark matter velocitites...");
music::ulog.Print("Computing dark matter velocitites...");
//... we do baryons and have velocity transfer functions, or we do SPH and not to shift
//... do DM first
@ -815,21 +821,21 @@ int main(int argc, const char *argv[])
data_forIO *= cosmo.vfact;
double sigv = compute_finest_sigma(data_forIO);
LOGINFO("sigma of %c-velocity of high-res DM is %f", 'x' + icoord, sigv);
music::ilog.Print("sigma of %c-velocity of high-res DM is %f", 'x' + icoord, sigv);
double meanv = compute_finest_mean(data_forIO);
LOGINFO("mean of %c-velocity of high-res particles is %f", 'x' + icoord, meanv);
LOGUSER("mean of %c-velocity of high-res particles is %f", 'x' + icoord, meanv);
music::ilog.Print("mean of %c-velocity of high-res particles is %f", 'x' + icoord, meanv);
music::ulog.Print("mean of %c-velocity of high-res particles is %f", 'x' + icoord, meanv);
double maxv = compute_finest_absmax(data_forIO);
LOGINFO("max of abs of %c-velocity of high-res particles is %f", 'x' + icoord, maxv);
music::ilog.Print("max of abs of %c-velocity of high-res particles is %f", 'x' + icoord, maxv);
coarsen_density(rh_Poisson, data_forIO, false);
// add counter velocity mode
if( do_counter_mode ) add_constant_value( data_forIO, -counter_mode_amp[icoord]*cosmo.vfact );
LOGUSER("Writing CDM velocities");
music::ulog.Print("Writing CDM velocities");
the_output_plugin->write_dm_velocity(icoord, data_forIO);
}
u.deallocate();
@ -839,7 +845,7 @@ int main(int argc, const char *argv[])
std::cout << "=============================================================\n";
std::cout << " COMPUTING BARYON VELOCITIES\n";
std::cout << "-------------------------------------------------------------\n";
LOGUSER("Computing baryon velocitites...");
music::ulog.Print("Computing baryon velocitites...");
//... do baryons
GenerateDensityHierarchy(cf, the_transfer_function_plugin, vbaryon, rh_TF, rand, f, false, bbshift);
coarsen_density(rh_Poisson, f, bspectral_sampling);
@ -874,21 +880,21 @@ int main(int argc, const char *argv[])
data_forIO *= cosmo.vfact;
double sigv = compute_finest_sigma(data_forIO);
LOGINFO("sigma of %c-velocity of high-res baryons is %f", 'x' + icoord, sigv);
music::ilog.Print("sigma of %c-velocity of high-res baryons is %f", 'x' + icoord, sigv);
double meanv = compute_finest_mean(data_forIO);
LOGINFO("mean of %c-velocity of high-res baryons is %f", 'x' + icoord, meanv);
LOGUSER("mean of %c-velocity of high-res baryons is %f", 'x' + icoord, meanv);
music::ilog.Print("mean of %c-velocity of high-res baryons is %f", 'x' + icoord, meanv);
music::ulog.Print("mean of %c-velocity of high-res baryons is %f", 'x' + icoord, meanv);
double maxv = compute_finest_absmax(data_forIO);
LOGINFO("max of abs of %c-velocity of high-res baryons is %f", 'x' + icoord, maxv);
music::ilog.Print("max of abs of %c-velocity of high-res baryons is %f", 'x' + icoord, maxv);
coarsen_density(rh_Poisson, data_forIO, false);
// add counter velocity mode
if( do_counter_mode ) add_constant_value( data_forIO, -counter_mode_amp[icoord]*cosmo.vfact );
LOGUSER("Writing baryon velocities");
music::ulog.Print("Writing baryon velocities");
the_output_plugin->write_gas_velocity(icoord, data_forIO);
}
u.deallocate();
@ -903,7 +909,7 @@ int main(int argc, const char *argv[])
else
{
//.. use 2LPT ...
LOGUSER("Entering 2LPT branch");
music::ulog.Print("Entering 2LPT branch");
grid_hierarchy f(nbnd), u1(nbnd), u2LPT(nbnd), f2LPT(nbnd);
@ -916,12 +922,12 @@ int main(int argc, const char *argv[])
if (my_tf_type == total)
{
std::cout << " COMPUTING VELOCITIES\n";
LOGUSER("Computing velocities...");
music::ulog.Print("Computing velocities...");
}
else
{
std::cout << " COMPUTING DARK MATTER VELOCITIES\n";
LOGUSER("Computing dark matter velocities...");
music::ulog.Print("Computing dark matter velocities...");
}
std::cout << "-------------------------------------------------------------\n";
@ -948,16 +954,16 @@ int main(int argc, const char *argv[])
else
f.deallocate();
LOGINFO("Computing 2LPT term....");
music::ilog.Print("Computing 2LPT term....");
if (!kspace2LPT)
compute_2LPT_source(u1, f2LPT, grad_order);
else
{
LOGUSER("computing term using FFT");
music::ulog.Print("computing term using FFT");
compute_2LPT_source_FFT(cf, u1, f2LPT);
}
LOGINFO("Solving 2LPT Poisson equation");
music::ilog.Print("Solving 2LPT Poisson equation");
u2LPT = u1;
u2LPT.zero();
the_poisson_solver->solve(f2LPT, u2LPT);
@ -996,11 +1002,11 @@ int main(int argc, const char *argv[])
double sigv = compute_finest_sigma(data_forIO);
double meanv = compute_finest_mean(data_forIO);
LOGINFO("mean of %c-velocity of high-res particles is %f", 'x' + icoord, meanv);
LOGUSER("mean of %c-velocity of high-res particles is %f", 'x' + icoord, meanv);
music::ilog.Print("mean of %c-velocity of high-res particles is %f", 'x' + icoord, meanv);
music::ulog.Print("mean of %c-velocity of high-res particles is %f", 'x' + icoord, meanv);
double maxv = compute_finest_absmax(data_forIO);
LOGINFO("max of abs of %c-velocity of high-res particles is %f", 'x' + icoord, maxv);
music::ilog.Print("max of abs of %c-velocity of high-res particles is %f", 'x' + icoord, maxv);
std::cerr << " - velocity component " << icoord << " : sigma = " << sigv << std::endl;
std::cerr << " - velocity component " << icoord << " : mean = " << meanv << std::endl;
@ -1011,12 +1017,12 @@ int main(int argc, const char *argv[])
counter_mode_amp[icoord] = compute_finest_mean(data_forIO);
if( do_counter_mode ) add_constant_value( data_forIO, -counter_mode_amp[icoord] );
LOGUSER("Writing CDM velocities");
music::ulog.Print("Writing CDM velocities");
the_output_plugin->write_dm_velocity(icoord, data_forIO);
if (do_baryons && !the_transfer_function_plugin->tf_has_velocities() && !bsph)
{
LOGUSER("Writing baryon velocities");
music::ulog.Print("Writing baryon velocities");
the_output_plugin->write_gas_velocity(icoord, data_forIO);
}
}
@ -1029,7 +1035,7 @@ int main(int argc, const char *argv[])
std::cout << "=============================================================\n";
std::cout << " COMPUTING BARYON VELOCITIES\n";
std::cout << "-------------------------------------------------------------\n";
LOGUSER("Computing baryon displacements...");
music::ulog.Print("Computing baryon displacements...");
GenerateDensityHierarchy(cf, the_transfer_function_plugin, vbaryon, rh_TF, rand, f, false, bbshift);
coarsen_density(rh_Poisson, f, bspectral_sampling);
@ -1045,7 +1051,7 @@ int main(int argc, const char *argv[])
//... compute 1LPT term
the_poisson_solver->solve(f, u1);
LOGINFO("Writing baryon potential");
music::ilog.Print("Writing baryon potential");
the_output_plugin->write_gas_potential(u1);
//... compute 2LPT term
@ -1094,11 +1100,11 @@ int main(int argc, const char *argv[])
double sigv = compute_finest_sigma(data_forIO);
double meanv = compute_finest_mean(data_forIO);
LOGINFO("mean of %c-velocity of high-res baryons is %f", 'x' + icoord, meanv);
LOGUSER("mean of %c-velocity of high-res baryons is %f", 'x' + icoord, meanv);
music::ilog.Print("mean of %c-velocity of high-res baryons is %f", 'x' + icoord, meanv);
music::ulog.Print("mean of %c-velocity of high-res baryons is %f", 'x' + icoord, meanv);
double maxv = compute_finest_absmax(data_forIO);
LOGINFO("max of abs of %c-velocity of high-res baryons is %f", 'x' + icoord, maxv);
music::ilog.Print("max of abs of %c-velocity of high-res baryons is %f", 'x' + icoord, maxv);
std::cerr << " - velocity component " << icoord << " : sigma = " << sigv << std::endl;
std::cerr << " - velocity component " << icoord << " : mean = " << meanv << std::endl;
@ -1108,7 +1114,7 @@ int main(int argc, const char *argv[])
// add counter velocity mode
if( do_counter_mode ) add_constant_value( data_forIO, -counter_mode_amp[icoord] );
LOGUSER("Writing baryon velocities");
music::ulog.Print("Writing baryon velocities");
the_output_plugin->write_gas_velocity(icoord, data_forIO);
}
data_forIO.deallocate();
@ -1118,7 +1124,7 @@ int main(int argc, const char *argv[])
std::cout << "=============================================================\n";
std::cout << " COMPUTING DARK MATTER DISPLACEMENTS\n";
std::cout << "-------------------------------------------------------------\n";
LOGUSER("Computing dark matter displacements...");
music::ulog.Print("Computing dark matter displacements...");
//... if baryons are enabled, the displacements have to be recomputed
//... otherwise we can compute them directly from the velocities
@ -1134,7 +1140,7 @@ int main(int argc, const char *argv[])
f.add_refinement_mask(rh_Poisson.get_coord_shift());
normalize_density(f);
LOGUSER("Writing CDM data");
music::ulog.Print("Writing CDM data");
the_output_plugin->write_dm_density(f);
the_output_plugin->write_dm_mass(f);
u1 = f;
@ -1206,14 +1212,14 @@ int main(int argc, const char *argv[])
the_poisson_solver->gradient(icoord, u1, data_forIO);
double dispmax = compute_finest_absmax(data_forIO);
LOGINFO("max. %c-displacement of HR particles is %f [mean dx]", 'x' + icoord, dispmax * (double)(1ll << data_forIO.levelmax()));
music::ilog.Print("max. %c-displacement of HR particles is %f [mean dx]", 'x' + icoord, dispmax * (double)(1ll << data_forIO.levelmax()));
coarsen_density(rh_Poisson, data_forIO, false);
// add counter mode
if( do_counter_mode ) add_constant_value( data_forIO, -counter_mode_amp[icoord]/cosmo.vfact );
LOGUSER("Writing CDM displacements");
music::ulog.Print("Writing CDM displacements");
the_output_plugin->write_dm_position(icoord, data_forIO);
}
@ -1225,7 +1231,7 @@ int main(int argc, const char *argv[])
std::cout << "=============================================================\n";
std::cout << " COMPUTING BARYON DENSITY\n";
std::cout << "-------------------------------------------------------------\n";
LOGUSER("Computing baryon density...");
music::ulog.Print("Computing baryon density...");
GenerateDensityHierarchy(cf, the_transfer_function_plugin, baryon, rh_TF, rand, f, true, false);
coarsen_density(rh_Poisson, f, bspectral_sampling);
@ -1259,7 +1265,7 @@ int main(int argc, const char *argv[])
compute_LLA_density(u1, f, grad_order);
normalize_density(f);
LOGUSER("Writing baryon density");
music::ulog.Print("Writing baryon density");
the_output_plugin->write_gas_density(f);
}
}
@ -1268,14 +1274,14 @@ int main(int argc, const char *argv[])
std::cout << "=============================================================\n";
std::cout << " COMPUTING BARYON DISPLACEMENTS\n";
std::cout << "-------------------------------------------------------------\n";
LOGUSER("Computing baryon displacements...");
music::ulog.Print("Computing baryon displacements...");
GenerateDensityHierarchy(cf, the_transfer_function_plugin, baryon, rh_TF, rand, f, false, bbshift);
coarsen_density(rh_Poisson, f, bspectral_sampling);
f.add_refinement_mask(rh_Poisson.get_coord_shift());
normalize_density(f);
LOGUSER("Writing baryon density");
music::ulog.Print("Writing baryon density");
the_output_plugin->write_gas_density(f);
u1 = f;
u1.zero();
@ -1331,7 +1337,7 @@ int main(int argc, const char *argv[])
if( do_counter_mode ) add_constant_value( data_forIO, -counter_mode_amp[icoord]/cosmo.vfact );
LOGUSER("Writing baryon displacements");
music::ulog.Print("Writing baryon displacements");
the_output_plugin->write_gas_position(icoord, data_forIO);
}
}
@ -1346,8 +1352,8 @@ int main(int argc, const char *argv[])
}
catch (std::runtime_error &excp)
{
LOGERR("Fatal error occured. Code will exit:");
LOGERR("Exception: %s", excp.what());
music::elog.Print("Fatal error occured. Code will exit:");
music::elog.Print("Exception: %s", excp.what());
std::cerr << " - " << excp.what() << std::endl;
std::cerr << " - A fatal error occured. We need to exit...\n";
bfatal = true;
@ -1358,7 +1364,7 @@ int main(int argc, const char *argv[])
if (!bfatal)
{
std::cout << " - Wrote output file \'" << outfname << "\'\n using plugin \'" << outformat << "\'...\n";
LOGUSER("Wrote output file \'%s\'.", outfname.c_str());
music::ulog.Print("Wrote output file \'%s\'.", outfname.c_str());
}
//------------------------------------------------------------------------------
@ -1383,9 +1389,9 @@ int main(int argc, const char *argv[])
ltime = time(NULL);
LOGUSER("Run finished succesfully on %s", asctime(localtime(&ltime)));
music::ulog.Print("Run finished succesfully on %s", asctime(localtime(&ltime)));
cf.log_dump();
cf.dump_to_log();
return 0;
}

87
src/mesh.hh
View file

@ -19,8 +19,6 @@
#include <math.h>
#include "config_file.hh"
#include "log.hh"
#include "region_generator.hh"
#include <array>
@ -239,7 +237,7 @@ public:
{
#ifdef DEBUG
if (ix < 0 || ix >= (int)m_nx || iy < 0 || iy >= (int)m_ny || iz < 0 || iz >= (int)m_nz)
LOGERR("Array index (%d,%d,%d) out of bounds", ix, iy, iz);
music::elog.Print("Array index (%d,%d,%d) out of bounds", ix, iy, iz);
#endif
return m_pdata[((size_t)ix * m_ny + (size_t)iy) * m_nz + (size_t)iz];
@ -250,7 +248,7 @@ public:
{
#ifdef DEBUG
if (ix < 0 || ix >= (int)m_nx || iy < 0 || iy >= (int)m_ny || iz < 0 || iz >= (int)m_nz)
LOGERR("Array index (%d,%d,%d) out of bounds", ix, iy, iz);
music::elog.Print("Array index (%d,%d,%d) out of bounds", ix, iy, iz);
#endif
return m_pdata[((size_t)ix * m_ny + (size_t)iy) * m_nz + (size_t)iz];
@ -293,7 +291,7 @@ public:
{
if (v.m_nx * v.m_ny * v.m_nz != m_nx * m_ny * m_nz)
{
LOGERR("Meshvar::operator*= : attempt to operate on incompatible data");
music::elog.Print("Meshvar::operator*= : attempt to operate on incompatible data");
throw std::runtime_error("Meshvar::operator*= : attempt to operate on incompatible data");
}
for (size_t i = 0; i < m_nx * m_ny * m_nz; ++i)
@ -307,7 +305,7 @@ public:
{
if (v.m_nx * v.m_ny * v.m_nz != m_nx * m_ny * m_nz)
{
LOGERR("Meshvar::operator/= : attempt to operate on incompatible data");
music::elog.Print("Meshvar::operator/= : attempt to operate on incompatible data");
throw std::runtime_error("Meshvar::operator/= : attempt to operate on incompatible data");
}
@ -322,7 +320,7 @@ public:
{
if (v.m_nx * v.m_ny * v.m_nz != m_nx * m_ny * m_nz)
{
LOGERR("Meshvar::operator+= : attempt to operate on incompatible data");
music::elog.Print("Meshvar::operator+= : attempt to operate on incompatible data");
throw std::runtime_error("Meshvar::operator+= : attempt to operate on incompatible data");
}
for (size_t i = 0; i < m_nx * m_ny * m_nz; ++i)
@ -336,7 +334,7 @@ public:
{
if (v.m_nx * v.m_ny * v.m_nz != m_nx * m_ny * m_nz)
{
LOGERR("Meshvar::operator-= : attempt to operate on incompatible data");
music::elog.Print("Meshvar::operator-= : attempt to operate on incompatible data");
throw std::runtime_error("Meshvar::operator-= : attempt to operate on incompatible data");
}
for (size_t i = 0; i < m_nx * m_ny * m_nz; ++i)
@ -590,7 +588,7 @@ public:
if (ilevel >= m_pgrids.size())
{
LOGERR("Attempt to access level %d but maxlevel = %d", ilevel, m_pgrids.size() - 1);
music::elog.Print("Attempt to access level %d but maxlevel = %d", ilevel, m_pgrids.size() - 1);
throw std::runtime_error("Fatal: attempt to access non-existent grid");
}
return m_pgrids[ilevel];
@ -601,7 +599,7 @@ public:
{
if (ilevel >= m_pgrids.size())
{
LOGERR("Attempt to access level %d but maxlevel = %d", ilevel, m_pgrids.size() - 1);
music::elog.Print("Attempt to access level %d but maxlevel = %d", ilevel, m_pgrids.size() - 1);
throw std::runtime_error("Fatal: attempt to access non-existent grid");
}
@ -992,7 +990,7 @@ public:
{
if (!is_consistent(gh))
{
LOGERR("GridHierarchy::operator*= : attempt to operate on incompatible data");
music::elog.Print("GridHierarchy::operator*= : attempt to operate on incompatible data");
throw std::runtime_error("GridHierarchy::operator*= : attempt to operate on incompatible data");
}
for (unsigned i = 0; i < m_pgrids.size(); ++i)
@ -1005,7 +1003,7 @@ public:
{
if (!is_consistent(gh))
{
LOGERR("GridHierarchy::operator/= : attempt to operate on incompatible data");
music::elog.Print("GridHierarchy::operator/= : attempt to operate on incompatible data");
throw std::runtime_error("GridHierarchy::operator/= : attempt to operate on incompatible data");
}
for (unsigned i = 0; i < m_pgrids.size(); ++i)
@ -1029,7 +1027,7 @@ public:
{
if (!is_consistent(gh))
{
LOGERR("GridHierarchy::operator-= : attempt to operate on incompatible data");
music::elog.Print("GridHierarchy::operator-= : attempt to operate on incompatible data");
throw std::runtime_error("GridHierarchy::operator-= : attempt to operate on incompatible data");
}
for (unsigned i = 0; i < m_pgrids.size(); ++i)
@ -1269,24 +1267,24 @@ public:
: cf_(cf)
{
//... query the parameter data we need
levelmin_ = cf_.getValue<unsigned>("setup", "levelmin");
levelmax_ = cf_.getValue<unsigned>("setup", "levelmax");
levelmin_tf_ = cf_.getValueSafe<unsigned>("setup", "levelmin_TF", levelmin_);
align_top_ = cf_.getValueSafe<bool>("setup", "align_top", false);
preserve_dims_ = cf_.getValueSafe<bool>("setup", "preserve_dims", false);
equal_extent_ = cf_.getValueSafe<bool>("setup", "force_equal_extent", false);
blocking_factor_ = cf.getValueSafe<unsigned>("setup", "blocking_factor", 0);
margin_ = cf.getValueSafe<int>("setup","convolution_margin",32);
levelmin_ = cf_.get_value<unsigned>("setup", "levelmin");
levelmax_ = cf_.get_value<unsigned>("setup", "levelmax");
levelmin_tf_ = cf_.get_value_safe<unsigned>("setup", "levelmin_TF", levelmin_);
align_top_ = cf_.get_value_safe<bool>("setup", "align_top", false);
preserve_dims_ = cf_.get_value_safe<bool>("setup", "preserve_dims", false);
equal_extent_ = cf_.get_value_safe<bool>("setup", "force_equal_extent", false);
blocking_factor_ = cf.get_value_safe<unsigned>("setup", "blocking_factor", 0);
margin_ = cf.get_value_safe<int>("setup","convolution_margin",32);
bool bnoshift = cf_.getValueSafe<bool>("setup", "no_shift", false);
bool force_shift = cf_.getValueSafe<bool>("setup", "force_shift", false);
bool bnoshift = cf_.get_value_safe<bool>("setup", "no_shift", false);
bool force_shift = cf_.get_value_safe<bool>("setup", "force_shift", false);
gridding_unit_ = cf.getValueSafe<unsigned>("setup", "gridding_unit", 2);
gridding_unit_ = cf.get_value_safe<unsigned>("setup", "gridding_unit", 2);
if (gridding_unit_ != 2 && blocking_factor_==0) {
blocking_factor_ = gridding_unit_; // THIS WILL LIKELY CAUSE PROBLEMS WITH NYX
}else if (gridding_unit_ != 2 && blocking_factor_!=0 && gridding_unit_!=blocking_factor_ ) {
LOGERR("incompatible gridding unit %d and blocking factor specified", gridding_unit_, blocking_factor_ );
music::elog.Print("incompatible gridding unit %d and blocking factor specified", gridding_unit_, blocking_factor_ );
throw std::runtime_error("Incompatible gridding unit and blocking factor!");
}
@ -1301,10 +1299,11 @@ public:
lxref_[i] = x1ref[i] - x0ref_[i];
bhave_nref = false;
std::string region_type = cf.getValueSafe<std::string>("setup", "region", "box");
std::string region_type = cf.get_value_safe<std::string>("setup", "region", "box");
LOGINFO("refinement region is \'%s\', w/ bounding box\n left = [%f,%f,%f]\n right = [%f,%f,%f]",
region_type.c_str(), x0ref_[0], x0ref_[1], x0ref_[2], x1ref[0], x1ref[1], x1ref[2]);
music::ilog << " refinement region is \'" << region_type.c_str() << "\', w/ bounding box" << std::endl;
music::ilog << " left = [" << x0ref_[0] << "," << x0ref_[1] << "," << x0ref_[2] << "]" << std::endl;
music::ilog << " right = [" << x1ref[0] << "," << x1ref[1] << "," << x1ref[2] << "]" << std::endl;
bhave_nref = the_region_generator->is_grid_dim_forced(lnref_);
}
@ -1327,10 +1326,10 @@ public:
if ((levelmin_ != levelmax_) && (!bnoshift || force_shift))
{
int random_base_grid_unit = cf.getValueSafe<int>("random","base_unit",1);
int random_base_grid_unit = cf.get_value_safe<int>("random","base_unit",1);
int shift_unit = get_shift_unit( random_base_grid_unit, levelmin_ );
if( shift_unit != 1 ){
LOGINFO("volume can only be shifted by multiples of %d coarse cells.",shift_unit);
music::ilog.Print("volume can only be shifted by multiples of %d coarse cells.",shift_unit);
}
xshift_[0] = (int)((0.5-xc[0]) * (double)ncoarse / shift_unit + 0.5) * shift_unit;//ARJ(int)((0.5 - xc[0]) * ncoarse);
xshift_[1] = (int)((0.5-xc[1]) * (double)ncoarse / shift_unit + 0.5) * shift_unit;//ARJ(int)((0.5 - xc[1]) * ncoarse);
@ -1349,11 +1348,11 @@ public:
char strtmp[32];
sprintf(strtmp, "%ld", xshift_[0]);
cf_.insertValue("setup", "shift_x", strtmp);
cf_.insert_value("setup", "shift_x", strtmp);
sprintf(strtmp, "%ld", xshift_[1]);
cf_.insertValue("setup", "shift_y", strtmp);
cf_.insert_value("setup", "shift_y", strtmp);
sprintf(strtmp, "%ld", xshift_[2]);
cf_.insertValue("setup", "shift_z", strtmp);
cf_.insert_value("setup", "shift_z", strtmp);
rshift_[0] = -(double)xshift_[0] / ncoarse;
rshift_[1] = -(double)xshift_[1] / ncoarse;
@ -1414,7 +1413,7 @@ public:
lnref_[1] % (1ul << (levelmax_ - levelmin_)) != 0 ||
lnref_[2] % (1ul << (levelmax_ - levelmin_)) != 0)
{
LOGERR("specified ref_dims and align_top=yes but cannot be aligned with coarse grid!");
music::elog.Print("specified ref_dims and align_top=yes but cannot be aligned with coarse grid!");
throw std::runtime_error("specified ref_dims and align_top=yes but cannot be aligned with coarse grid!");
}
}
@ -1458,7 +1457,7 @@ public:
else
{
//... require alignment with coarser grid
LOGINFO("Internal refinement bounding box error: [%d,%d]x[%d,%d]x[%d,%d]", il, ir, jl, jr, kl, kr);
music::ilog.Print("- Internal refinement bounding box: [%d,%d]x[%d,%d]x[%d,%d]", il, ir, jl, jr, kl, kr);
il -= il % gridding_unit_;
jl -= jl % gridding_unit_;
@ -1492,7 +1491,7 @@ public:
if (il >= ir || jl >= jr || kl >= kr)
{
LOGERR("Internal refinement bounding box error: [%d,%d]x[%d,%d]x[%d,%d]", il, ir, jl, jr, kl, kr);
music::elog.Print("Internal refinement bounding box error: [%d,%d]x[%d,%d]x[%d,%d]", il, ir, jl, jr, kl, kr);
throw std::runtime_error("refinement_hierarchy: Internal refinement bounding box error 1");
}
//... determine offsets
@ -1508,7 +1507,7 @@ public:
if (bhave_nref && (lnref_[0] != lnref_[1] || lnref_[0] != lnref_[2]))
{
LOGERR("Specified equal_extent=yes conflicting with ref_dims which are not equal.");
music::elog.Print("Specified equal_extent=yes conflicting with ref_dims which are not equal.");
throw std::runtime_error("Specified equal_extent=yes conflicting with ref_dims which are not equal.");
}
size_t ilevel = levelmax_;
@ -1531,7 +1530,7 @@ public:
}
}
padding_ = cf_.getValueSafe<unsigned>("setup", "padding", 8);
padding_ = cf_.get_value_safe<unsigned>("setup", "padding", 8);
//... determine position of coarser grids
for (unsigned ilevel = levelmax_ - 1; ilevel > levelmin_; --ilevel)
@ -1585,7 +1584,7 @@ public:
if (il >= ir || jl >= jr || kl >= kr || il < 0 || jl < 0 || kl < 0)
{
LOGERR("Internal refinement bounding box error: [%d,%d]x[%d,%d]x[%d,%d], level=%d", il, ir, jl, jr, kl, kr, ilevel);
music::elog.Print("Internal refinement bounding box error: [%d,%d]x[%d,%d]x[%d,%d], level=%d", il, ir, jl, jr, kl, kr, ilevel);
throw std::runtime_error("refinement_hierarchy: Internal refinement bounding box error 2");
}
absoffsets_[ilevel] = { il, jl, kl };
@ -1650,7 +1649,7 @@ public:
len_[ilevel][1] > index_t(1ul << (ilevel - 1)) ||
len_[ilevel][2] > index_t(1ul << (ilevel - 1)))
{
LOGERR("On level %d, subgrid is larger than half the box. This is not allowed!", ilevel);
music::elog.Print("On level %d, subgrid is larger than half the box. This is not allowed!", ilevel);
throw std::runtime_error("Fatal: Subgrid larger than half boxin zoom.");
}
}
@ -1758,7 +1757,7 @@ public:
}
if ((old_levelmin != levelmin_) && print)
LOGINFO("refinement_hierarchy: set new levelmin to %d", levelmin_);
music::ilog.Print("refinement_hierarchy: set new levelmin to %d", levelmin_);
}
//! get absolute grid offset for a specified level along a specified dimension (in fine grid units)
@ -1832,11 +1831,11 @@ public:
void output_log(void) const
{
LOGUSER(" Domain shifted by (%5d,%5d,%5d)", xshift_[0], xshift_[1], xshift_[2]);
music::ulog.Print(" Domain shifted by (%5d,%5d,%5d)", xshift_[0], xshift_[1], xshift_[2]);
for (unsigned ilevel = levelmin_; ilevel <= levelmax_; ++ilevel)
{
LOGUSER(" Level %3d : offset = (%5d,%5d,%5d)", ilevel, offsets_[ilevel][0], offsets_[ilevel][1], offsets_[ilevel][2]);
LOGUSER(" size = (%5d,%5d,%5d)", len_[ilevel][0], len_[ilevel][1], len_[ilevel][2]);
music::ulog.Print(" Level %3d : offset = (%5d,%5d,%5d)", ilevel, offsets_[ilevel][0], offsets_[ilevel][1], offsets_[ilevel][2]);
music::ulog.Print(" size = (%5d,%5d,%5d)", len_[ilevel][0], len_[ilevel][1], len_[ilevel][2]);
}
}
};

12
src/mg_solver.hh
View file

@ -454,7 +454,7 @@ double solver<S,I,O,T>::compute_error( const GridHierarchy<T>& uh, const GridHie
if( verbose )
std::cout << " Level " << std::setw(6) << ilevel << ", Error = " << err << std::endl;
LOGDEBUG("[mg] level %3d, residual %g, rel. error %g",ilevel, mean_res, err);
music::dlog.Print("[mg] level %3d, residual %g, rel. error %g",ilevel, mean_res, err);
maxerr = std::max(maxerr,err);
@ -504,7 +504,7 @@ double solver<S,I,O,T>::compute_RMS_resid( const GridHierarchy<T>& uh, const Gri
if( verbose && !m_is_ini )
std::cout << " Level " << std::setw(6) << ilevel << ", Error = " << err_rel << std::endl;
LOGDEBUG("[mg] level %3d, rms residual %g, rel. error %g",ilevel, err_abs, err_rel);
music::dlog.Print("[mg] level %3d, rms residual %g, rel. error %g",ilevel, err_abs, err_rel);
if( err_rel > maxerr )
maxerr = err_rel;
@ -535,7 +535,7 @@ double solver<S,I,O,T>::solve( GridHierarchy<T>& uh, double acc, double h, bool
while (true)
{
LOGUSER("Performing multi-grid V-cycle...");
music::ulog.Print("Performing multi-grid V-cycle...");
twoGrid( uh.levelmax() );
//err = compute_RMS_resid( *m_pu, *m_pf, fullverbose );
@ -543,7 +543,7 @@ double solver<S,I,O,T>::solve( GridHierarchy<T>& uh, double acc, double h, bool
++niter;
if( fullverbose ){
LOGUSER(" multigrid iteration %3d, maximum RMS residual = %g", niter, err );
music::ulog.Print(" multigrid iteration %3d, maximum RMS residual = %g", niter, err );
std::cout << " - Step No. " << std::setw(3) << niter << ", Max Err = " << err << std::endl;
std::cout << " ---------------------------------------------------\n";
}
@ -558,12 +558,12 @@ double solver<S,I,O,T>::solve( GridHierarchy<T>& uh, double acc, double h, bool
if( err > acc )
{
std::cout << "Error : no convergence in Poisson solver" << std::endl;
LOGERR("No convergence in Poisson solver, final error: %g.",err);
music::elog.Print("No convergence in Poisson solver, final error: %g.",err);
}
else if( verbose )
{
std::cout << " - Converged in " << niter << " steps to " << maxerr << std::endl;
LOGUSER("Poisson solver converged to max. error of %g in %d steps.",err,niter);
music::ulog.Print("Poisson solver converged to max. error of %g in %d steps.",err,niter);
}

2
src/output.cc
View file

@ -36,7 +36,7 @@ void print_output_plugins()
output_plugin *select_output_plugin( config_file& cf )
{
std::string formatname = cf.getValue<std::string>( "output", "format" );
std::string formatname = cf.get_value<std::string>( "output", "format" );
output_plugin_creator *the_output_plugin_creator
= get_output_plugin_map()[ formatname ];

8
src/output.hh
View file

@ -62,7 +62,7 @@ protected:
for( unsigned i=levelmin_; i<=levelmax_; ++i )
{
sprintf( str, "%s(%u,%d)", name.c_str(), i, icomp );
*oit = cf_.getValue<unsigned>( "setup", str );
*oit = cf_.get_value<unsigned>( "setup", str );
++oit;
}
}
@ -73,9 +73,9 @@ public:
explicit output_plugin( config_file& cf )
: cf_(cf)
{
fname_ = cf.getValue<std::string>("output","filename");
levelmin_ = cf.getValue<unsigned>( "setup", "levelmin" );
levelmax_ = cf.getValue<unsigned>( "setup", "levelmax" );
fname_ = cf.get_value<std::string>("output","filename");
levelmin_ = cf.get_value<unsigned>( "setup", "levelmin" );
levelmax_ = cf.get_value<unsigned>( "setup", "levelmax" );
query_grid_prop( "offset", 0, std::back_inserter(offx_) );
query_grid_prop( "offset", 1, std::back_inserter(offy_) );

2
src/plugins/convex_hull.hh
View file

@ -9,7 +9,7 @@
#include <omp.h>
#endif
#include "log.hh"
#include "logger.hh"
#include <array>

32
src/plugins/nyx_plugin/output_nyx.cpp
View file

@ -143,7 +143,7 @@ public:
char **argv;
BoxLib::Initialize(argc,argv);
bool bhave_hydro = cf_.getValue<bool>("setup","baryons");
bool bhave_hydro = cf_.get_value<bool>("setup","baryons");
if (bhave_hydro)
n_data_items = 10;
@ -249,11 +249,11 @@ public:
// throw std::runtime_error("Error in nyx_output_plugin!");
// }
bool haveblockingfactor = cf.containsKey( "setup", "blocking_factor");
bool haveblockingfactor = cf.contains_key( "setup", "blocking_factor");
if( !haveblockingfactor )
{
LOGERR("nyx output plug-in requires that \'blocking_factor\' is set!");
music::elog.Print("nyx output plug-in requires that \'blocking_factor\' is set!");
throw std::runtime_error("nyx output plug-in requires that \'blocking_factor\' is set!");
}
@ -265,19 +265,19 @@ public:
the_sim_header.offset.push_back( 0 );
the_sim_header.offset.push_back( 0 );
the_sim_header.a_start = 1.0/(1.0+cf.getValue<double>("setup","zstart"));
the_sim_header.dx = cf.getValue<double>("setup","boxlength")/the_sim_header.dimensions[0]/(cf.getValue<double>("cosmology","H0")*0.01); // not sure?!?
the_sim_header.boxlength=cf.getValue<double>("setup","boxlength");
the_sim_header.h0 = cf.getValue<double>("cosmology","H0")*0.01;
the_sim_header.a_start = 1.0/(1.0+cf.get_value<double>("setup","zstart"));
the_sim_header.dx = cf.get_value<double>("setup","boxlength")/the_sim_header.dimensions[0]/(cf.get_value<double>("cosmology","H0")*0.01); // not sure?!?
the_sim_header.boxlength=cf.get_value<double>("setup","boxlength");
the_sim_header.h0 = cf.get_value<double>("cosmology","H0")*0.01;
if( bhave_hydro )
the_sim_header.omega_b = cf.getValue<double>("cosmology","Omega_b");
the_sim_header.omega_b = cf.get_value<double>("cosmology","Omega_b");
else
the_sim_header.omega_b = 0.0;
the_sim_header.omega_m = cf.getValue<double>("cosmology","Omega_m");
the_sim_header.omega_v = cf.getValue<double>("cosmology","Omega_L");
the_sim_header.vfact = cf.getValue<double>("cosmology","vfact")*the_sim_header.h0; //.. need to multiply by h, nyx wants this factor for non h-1 units
the_sim_header.omega_m = cf.get_value<double>("cosmology","Omega_m");
the_sim_header.omega_v = cf.get_value<double>("cosmology","Omega_L");
the_sim_header.vfact = cf.get_value<double>("cosmology","vfact")*the_sim_header.h0; //.. need to multiply by h, nyx wants this factor for non h-1 units
std::cout << "creating output object" << std::endl;
}
@ -432,7 +432,7 @@ public:
inputs << "cosmo.ic-source = MUSIC " << std::endl;
inputs << "amr.blocking_factor = " << cf_.getValue<double>("setup","blocking_factor") << std::endl;
inputs << "amr.blocking_factor = " << cf_.get_value<double>("setup","blocking_factor") << std::endl;
inputs << "nyx.do_hydro = "<< (the_sim_header.omega_b>0?1:0) << std::endl;
inputs << "amr.max_level = " << levelmax_-levelmin_ << std::endl;
@ -467,7 +467,7 @@ public:
const Real cur_time = 0.0;
std::cout << "in writeLevelPlotFile" << std::endl;
double h0 = cf_.getValue<double>("cosmology", "H0")*0.01;
double h0 = cf_.get_value<double>("cosmology", "H0")*0.01;
// for (MFIter mfi(mf); mfi.isValid(); ++mfi)
// {
@ -495,7 +495,7 @@ public:
for (i = 0; i < BL_SPACEDIM; i++)
os << 0 << ' '; //ProbLo
os << '\n';
double boxlength = cf_.getValue<double>("setup","boxlength");
double boxlength = cf_.get_value<double>("setup","boxlength");
for (i = 0; i < BL_SPACEDIM; i++)
os << boxlength/h0 << ' '; //ProbHi
os << '\n';
@ -522,7 +522,7 @@ public:
os << 0 << ' ';
os << '\n';
double dx = cf_.getValue<double>("setup","boxlength")/gridp/h0;
double dx = cf_.get_value<double>("setup","boxlength")/gridp/h0;
for (i = 0; i <= f_lev; i++)
{
for (int k = 0; k < BL_SPACEDIM; k++)
@ -558,7 +558,7 @@ public:
os << 0 << '\n';
double cellsize[3];
double dx = cf_.getValue<double>("setup","boxlength")/gridp/h0;
double dx = cf_.get_value<double>("setup","boxlength")/gridp/h0;
for (n = 0; n < BL_SPACEDIM; n++)
{
cellsize[n] = dx;

32
src/plugins/output_arepo.cc
View file

@ -523,26 +523,26 @@ public:
{
// ensure that everyone knows we want to do SPH, implies: bsph=1, bbshift=1, decic_baryons=1
// -> instead of just writing gas densities (which are here ignored), the gas displacements are also written
cf.insertValue("setup", "do_SPH", "yes");
cf.insert_value("setup", "do_SPH", "yes");
// init header and config parameters
nPartTotal = std::vector<long long>(NTYPES, 0);
massTable = std::vector<double>(NTYPES, 0.0);
coarsePartType = cf.getValueSafe<unsigned>("output", "arepo_coarsetype", COARSE_DM_DEFAULT_PARTTYPE);
UnitLength_in_cm = cf.getValueSafe<double>("output", "arepo_unitlength", 3.085678e21); // 1.0 kpc
UnitMass_in_g = cf.getValueSafe<double>("output", "arepo_unitmass", 1.989e43); // 1.0e10 solar masses
UnitVelocity_in_cm_per_s = cf.getValueSafe<double>("output", "arepo_unitvel", 1e5); // 1 km/sec
coarsePartType = cf.get_value_safe<unsigned>("output", "arepo_coarsetype", COARSE_DM_DEFAULT_PARTTYPE);
UnitLength_in_cm = cf.get_value_safe<double>("output", "arepo_unitlength", 3.085678e21); // 1.0 kpc
UnitMass_in_g = cf.get_value_safe<double>("output", "arepo_unitmass", 1.989e43); // 1.0e10 solar masses
UnitVelocity_in_cm_per_s = cf.get_value_safe<double>("output", "arepo_unitvel", 1e5); // 1 km/sec
omega0 = cf.getValue<double>("cosmology", "Omega_m");
omega_b = cf.getValue<double>("cosmology", "Omega_b");
omega_L = cf.getValue<double>("cosmology", "Omega_L");
redshift = cf.getValue<double>("setup", "zstart");
boxSize = cf.getValue<double>("setup", "boxlength");
doBaryons = cf.getValueSafe<bool>("setup", "baryons", false);
useLongIDs = cf.getValueSafe<bool>("output", "arepo_longids", false);
numFiles = cf.getValueSafe<unsigned>("output", "arepo_num_files", 1);
doublePrec = cf.getValueSafe<bool>("output", "arepo_doubleprec", 0);
omega0 = cf.get_value<double>("cosmology", "Omega_m");
omega_b = cf.get_value<double>("cosmology", "Omega_b");
omega_L = cf.get_value<double>("cosmology", "Omega_L");
redshift = cf.get_value<double>("setup", "zstart");
boxSize = cf.get_value<double>("setup", "boxlength");
doBaryons = cf.get_value_safe<bool>("setup", "baryons", false);
useLongIDs = cf.get_value_safe<bool>("output", "arepo_longids", false);
numFiles = cf.get_value_safe<unsigned>("output", "arepo_num_files", 1);
doublePrec = cf.get_value_safe<bool>("output", "arepo_doubleprec", 0);
for (unsigned i = 0; i < numFiles; i++)
nPart.push_back(std::vector<unsigned int>(NTYPES, 0));
@ -583,7 +583,7 @@ public:
}
// calculate Tini for gas
hubbleParam = cf.getValue<double>("cosmology", "H0") / 100.0;
hubbleParam = cf.get_value<double>("cosmology", "H0") / 100.0;
double astart = 1.0 / (1.0 + redshift);
double h2 = hubbleParam * hubbleParam;
@ -604,7 +604,7 @@ public:
if (coarsePartType == GAS_PARTTYPE || coarsePartType == HIGHRES_DM_PARTTYPE)
throw std::runtime_error("Error: Specified illegal Arepo particle type for coarse particles.");
if (coarsePartType == STAR_PARTTYPE)
LOGWARN("WARNING: Specified coarse particle type will collide with stars if USE_SFR enabled.");
music::wlog.Print("WARNING: Specified coarse particle type will collide with stars if USE_SFR enabled.");
// create file(s)
for (unsigned i = 0; i < numFiles; i++)

70
src/plugins/output_art.cc
View file

@ -122,7 +122,7 @@ protected:
if (!this->good())
{
LOGERR("Could not open buffer file in ART output plug-in");
music::elog.Print("Could not open buffer file in ART output plug-in");
throw std::runtime_error("Could not open buffer file in ART output plug-in");
}
@ -130,8 +130,8 @@ protected:
if (blk != (size_t)(npart * sizeof(T_store)))
{
LOGERR("Internal consistency error in ART output plug-in");
LOGERR("Expected %d bytes in temp file but found %d", npart * (unsigned)sizeof(T_store), blk);
music::elog.Print("Internal consistency error in ART output plug-in");
music::elog.Print("Expected %d bytes in temp file but found %d", npart * (unsigned)sizeof(T_store), blk);
throw std::runtime_error("Internal consistency error in ART output plug-in");
}
}
@ -147,7 +147,7 @@ protected:
if (!this->good())
{
LOGERR("Could not open buffer file \'%s\' in ART output plug-in", fname.c_str());
music::elog.Print("Could not open buffer file \'%s\' in ART output plug-in", fname.c_str());
throw std::runtime_error("Could not open buffer file in ART output plug-in");
}
@ -155,8 +155,8 @@ protected:
if (blk != (size_t)(npart * sizeof(T_store)))
{
LOGERR("Internal consistency error in ART output plug-in");
LOGERR("Expected %d bytes in temp file but found %d", npart * (unsigned)sizeof(T_store), blk);
music::elog.Print("Internal consistency error in ART output plug-in");
music::elog.Print("Expected %d bytes in temp file but found %d", npart * (unsigned)sizeof(T_store), blk);
throw std::runtime_error("Internal consistency error in ART output plug-in");
}
}
@ -180,7 +180,7 @@ protected:
int blksize = hsize_;
if (swap_endianness_)
{
LOGINFO("ART : swap_endianness option enabled");
music::ilog.Print("ART : swap_endianness option enabled");
blksize = bytereorder(blksize);
this_header.aexpN = bytereorder(this_header.aexpN);
this_header.aexp0 = bytereorder(this_header.aexp0);
@ -243,7 +243,7 @@ protected:
ofs.write((char *)&this_header.extras, sizeof(this_header.extras));
ofs.write((char *)&blksize, sizeof(int));
ofs.close();
LOGINFO("ART : done writing header file.");
music::ilog.Print("ART : done writing header file.");
}
void write_pt_file(void) //pt.dat
@ -262,7 +262,7 @@ protected:
ofs.write((char *)&this_ptf, sizeof(ptf));
ofs.write((char *)&blksize, sizeof(int));
ofs.close();
LOGINFO("ART : done writing pt file.");
music::ilog.Print("ART : done writing pt file.");
}
void adjust_buf_endianness(T_store *buf)
@ -321,7 +321,7 @@ protected:
size_t npleft, n2read;
size_t npcdm = npcdm_;
LOGINFO("writing DM data to ART format file");
music::ilog.Print("writing DM data to ART format file");
//ofs.open(fname_.c_str(), std::ios::binary|std::ios::trunc );
pistream ifs_x, ifs_y, ifs_z, ifs_vx, ifs_vy, ifs_vz;
@ -401,7 +401,7 @@ protected:
delete[] tmp5;
delete[] tmp6;
LOGINFO("ART : done writing DM file.");
music::ilog.Print("ART : done writing DM file.");
}
/*
@ -441,7 +441,7 @@ protected:
size_t npleft, n2read;
size_t npcgas = npcdm_; // # of gas elemets should be equal to # of dm elements
LOGINFO("writing gas data to ART format file");
music::ilog.Print("writing gas data to ART format file");
//ofs.open(fname_.c_str(), std::ios::binary|std::ios::trunc );
pistream ifs_x, ifs_y, ifs_z, ifs_vx, ifs_vy, ifs_vz;
@ -521,14 +521,14 @@ protected:
delete[] tmp5;
delete[] tmp6;
LOGINFO("ART : done writing gas file.");
music::ilog.Print("ART : done writing gas file.");
// Temperature
const double Tcmb0 = 2.726;
const double h2 = header_.hubble * header_.hubble;
const double adec = 1.0 / (160. * pow(omegab_ * h2 / 0.022, 2.0 / 5.0));
const double Tini = astart_ < adec ? Tcmb0 / astart_ : Tcmb0 / astart_ / astart_ * adec;
const double mu = (Tini > 1.e4) ? 4.0 / (8. - 5. * YHe_) : 4.0 / (1. + 3. * (1. - YHe_));
LOGINFO("ART : set initial gas temperature to %.3f K (%.3f K/mu)", Tini, Tini / mu);
music::ilog.Print("ART : set initial gas temperature to %.3f K (%.3f K/mu)", Tini, Tini / mu);
}
public:
@ -538,29 +538,29 @@ public:
if (mkdir(fname_.c_str(), 0777))
;
do_baryons_ = cf.getValueSafe<bool>("setup", "baryons", false);
do_baryons_ = cf.get_value_safe<bool>("setup", "baryons", false);
// We need to say that we want to do SPH for baryons
// because if not MUSIC does not calculate/write gas positions
cf.insertValue("setup", "do_SPH", "yes");
cf.insert_value("setup", "do_SPH", "yes");
// header size (alignment problem)
hsize_ = 529; // dm & hydro run
omegab_ = cf.getValueSafe<double>("cosmology", "Omega_b", 0.0);
omegam_ = cf.getValue<double>("cosmology", "Omega_m");
zstart_ = cf.getValue<double>("setup", "zstart");
omegab_ = cf.get_value_safe<double>("cosmology", "Omega_b", 0.0);
omegam_ = cf.get_value<double>("cosmology", "Omega_m");
zstart_ = cf.get_value<double>("setup", "zstart");
astart_ = 1.0 / (1.0 + zstart_);
swap_endianness_ = cf.getValueSafe<bool>("output", "art_swap_endian", true);
swap_endianness_ = cf.get_value_safe<bool>("output", "art_swap_endian", true);
int levelmin = cf.getValue<unsigned>("setup", "levelmin");
int levelmax = cf.getValue<unsigned>("setup", "levelmax");
int levelmin = cf.get_value<unsigned>("setup", "levelmin");
int levelmax = cf.get_value<unsigned>("setup", "levelmax");
block_buf_size_ = (size_t)(pow(pow(2, levelmax), 2)); //Npage=nrow^2; Number of particles in each page
YHe_ = cf.getValueSafe<double>("cosmology", "YHe", 0.248);
gamma_ = cf.getValueSafe<double>("cosmology", "gamma", 5.0 / 3.0);
YHe_ = cf.get_value_safe<double>("cosmology", "YHe", 0.248);
gamma_ = cf.get_value_safe<double>("cosmology", "gamma", 5.0 / 3.0);
// Set header
std::string thead;
thead = cf.getValueSafe<std::string>("output", "header", "ICs generated using MUSIC");
thead = cf.get_value_safe<std::string>("output", "header", "ICs generated using MUSIC");
strcpy(header_.head, thead.c_str()); // text for the header; any easy way to add also the version?
std::string ws = " "; // Filling with blanks. Any better way?
for (int i = thead.size(); i < 45; i++)
@ -570,7 +570,7 @@ public:
header_.aexpN = astart_;
header_.aexp0 = header_.aexpN;
header_.amplt = 0.0; // Amplitude of density fluctuations
header_.astep = cf.getValue<double>("output", "astep"); // Seems that this must also be in the config file
header_.astep = cf.get_value<double>("output", "astep"); // Seems that this must also be in the config file
ptf_.astep = header_.astep; // to write pt file
header_.istep = 0; // step (=0 in IC)
header_.partw = 0.0; // mass of highest res particle. SEE BELOW
@ -590,12 +590,12 @@ public:
//header_.partw SEE BELOW
header_.Nseed = 0; // random number used ( 0 for MUSIC? or set the random number used in the lowest level?)
header_.Om0 = cf.getValue<double>("cosmology", "Omega_m"); //Omega_m
header_.Oml0 = cf.getValue<double>("cosmology", "Omega_L"); //Omega_L
header_.hubble = cf.getValue<double>("cosmology", "H0") / 100; //hubble constant h=H/100
header_.Om0 = cf.get_value<double>("cosmology", "Omega_m"); //Omega_m
header_.Oml0 = cf.get_value<double>("cosmology", "Omega_L"); //Omega_L
header_.hubble = cf.get_value<double>("cosmology", "H0") / 100; //hubble constant h=H/100
header_.Wp5 = 0.0; // 0.0
header_.Ocurv = 1.0 - header_.Oml0 - header_.Om0; //
header_.Omb0 = cf.getValue<double>("cosmology", "Omega_b");
header_.Omb0 = cf.get_value<double>("cosmology", "Omega_b");
; // this parameter only appears in header in hydro runs
for (int i = 0; i < 10; i++)
{
@ -611,12 +611,12 @@ public:
{
header_.extras[i] = 0.0; //extras[20-99]
}
header_.extras[13] = cf.getValueSafe<double>("cosmology", "Omega_b", 0.0);
header_.extras[14] = cf.getValue<double>("cosmology", "sigma_8");
header_.extras[15] = cf.getValue<double>("cosmology", "nspec"); //Slope of the Power spectrum
header_.extras[79] = cf.getValue<double>("setup", "boxlength");
header_.extras[13] = cf.get_value_safe<double>("cosmology", "Omega_b", 0.0);
header_.extras[14] = cf.get_value<double>("cosmology", "sigma_8");
header_.extras[15] = cf.get_value<double>("cosmology", "nspec"); //Slope of the Power spectrum
header_.extras[79] = cf.get_value<double>("setup", "boxlength");
LOGINFO("ART : done header info.");
music::ilog.Print("ART : done header info.");
}
void write_dm_mass(const grid_hierarchy &gh)

66
src/plugins/output_cart.cc
View file

@ -114,7 +114,7 @@ class cart_output_plugin : public output_plugin
if( !this->good() )
{
LOGERR("Could not open buffer file in CART output plug-in");
music::elog.Print("Could not open buffer file in CART output plug-in");
throw std::runtime_error("Could not open buffer file in CART output plug-in");
}
@ -122,8 +122,8 @@ class cart_output_plugin : public output_plugin
if( blk != (size_t)(npart*sizeof(T_store)) )
{
LOGERR("Internal consistency error in CART output plug-in");
LOGERR("Expected %d bytes in temp file but found %d",npart*(unsigned)sizeof(T_store),blk);
music::elog.Print("Internal consistency error in CART output plug-in");
music::elog.Print("Expected %d bytes in temp file but found %d",npart*(unsigned)sizeof(T_store),blk);
throw std::runtime_error("Internal consistency error in CART output plug-in");
}
}
@ -137,7 +137,7 @@ class cart_output_plugin : public output_plugin
if( !this->good() )
{
LOGERR("Could not open buffer file \'%s\' in CART output plug-in",fname.c_str());
music::elog.Print("Could not open buffer file \'%s\' in CART output plug-in",fname.c_str());
throw std::runtime_error("Could not open buffer file in CART output plug-in");
}
@ -145,8 +145,8 @@ class cart_output_plugin : public output_plugin
if( blk != (size_t)(npart*sizeof(T_store)) )
{
LOGERR("Internal consistency error in CART output plug-in");
LOGERR("Expected %d bytes in temp file but found %d",npart*(unsigned)sizeof(T_store),blk);
music::elog.Print("Internal consistency error in CART output plug-in");
music::elog.Print("Expected %d bytes in temp file but found %d",npart*(unsigned)sizeof(T_store),blk);
throw std::runtime_error("Internal consistency error in CART output plug-in");
}
}
@ -171,7 +171,7 @@ class cart_output_plugin : public output_plugin
int blksize = hsize_;
if( swap_endianness_ )
{
LOGINFO("CART : swap_endianness option enabled");
music::ilog.Print("CART : swap_endianness option enabled");
blksize = bytereorder( blksize );
this_header.aexpN = bytereorder( this_header.aexpN );
this_header.aexp0 = bytereorder( this_header.aexp0 );
@ -240,7 +240,7 @@ class cart_output_plugin : public output_plugin
ofs.write( (char *)&this_header.extras,sizeof(this_header.extras));
ofs.write( (char *)&blksize, sizeof(int) );
ofs.close();
LOGINFO("CART : done writing header file.");
music::ilog.Print("CART : done writing header file.");
}
void write_pt_file( void ) //pt.dat
@ -259,7 +259,7 @@ class cart_output_plugin : public output_plugin
// ofs.write( (char *)&this_ptf,sizeof(ptf));
// ofs.write( (char *)&blksize, sizeof(int) );
// ofs.close();
// LOGINFO("CART : done writing pt file.");
// music::ilog.Print("CART : done writing pt file.");
}
@ -319,7 +319,7 @@ class cart_output_plugin : public output_plugin
size_t npleft, n2read;
size_t npcdm = npcdm_;
LOGINFO("writing DM data to CART format file");
music::ilog.Print("writing DM data to CART format file");
//ofs.open(fname_.c_str(), std::ios::binary|std::ios::trunc );
pistream ifs_x, ifs_y, ifs_z, ifs_vx, ifs_vy, ifs_vz;
@ -395,7 +395,7 @@ class cart_output_plugin : public output_plugin
delete[] tmp5;
delete[] tmp6;
LOGINFO("CART : done writing DM file.");
music::ilog.Print("CART : done writing DM file.");
}
@ -446,7 +446,7 @@ class cart_output_plugin : public output_plugin
size_t npleft, n2read;
size_t npcgas = npcdm_; // # of gas elemets should be equal to # of dm elements
LOGINFO("writing gas data to CART format file");
music::ilog.Print("writing gas data to CART format file");
//ofs.open(fname_.c_str(), std::ios::binary|std::ios::trunc );
pistream ifs_x, ifs_y, ifs_z, ifs_vx, ifs_vy, ifs_vz, ifs_pma;
@ -530,7 +530,7 @@ class cart_output_plugin : public output_plugin
delete[] tmp6;
delete[] tmp7;
LOGINFO("CART : done writing gas file.");
music::ilog.Print("CART : done writing gas file.");
}
@ -542,26 +542,26 @@ class cart_output_plugin : public output_plugin
{
mkdir( fname_.c_str(), 0777 );
do_baryons_ = cf.getValue<bool>("setup","baryons");
do_baryons_ = cf.get_value<bool>("setup","baryons");
hsize_ = 533; // dm & hydro run (omega_b is included in header -- 529 for oldstyle)
omegab_ = cf.getValue<double>("cosmology","Omega_b");
omegam_ = cf.getValue<double>("cosmology","Omega_m");
zstart_ = cf.getValue<double>("setup","zstart");
omegab_ = cf.get_value<double>("cosmology","Omega_b");
omegam_ = cf.get_value<double>("cosmology","Omega_m");
zstart_ = cf.get_value<double>("setup","zstart");
astart_ = 1.0/(1.0+zstart_);
//snl this doesn't corrently swap particle endianness and you check on the CART end anyway
swap_endianness_ = cf.getValueSafe<bool>("output","art_swap_endian",false);
swap_endianness_ = cf.get_value_safe<bool>("output","art_swap_endian",false);
int levelmin = cf.getValue<unsigned>("setup","levelmin");
int levelmax = cf.getValue<unsigned>("setup","levelmax");
int levelmin = cf.get_value<unsigned>("setup","levelmin");
int levelmax = cf.get_value<unsigned>("setup","levelmax");
block_buf_size_ = (size_t) (pow(pow(2,levelmax),2)); //Npage=nrow^2; Number of particles in each page
YHe_ = cf.getValueSafe<double>("cosmology","YHe",0.248);
gamma_ = cf.getValueSafe<double>("cosmology","gamma",5.0/3.0);
YHe_ = cf.get_value_safe<double>("cosmology","YHe",0.248);
gamma_ = cf.get_value_safe<double>("cosmology","gamma",5.0/3.0);
// Set header
std::string thead;
thead=cf.getValueSafe<std::string>("output","header","ICs generated using MUSIC");
thead=cf.get_value_safe<std::string>("output","header","ICs generated using MUSIC");
strcpy(header_.head,thead.c_str()); // text for the header; any easy way to add also the version?
std::string ws = " "; // Filling with blanks. Any better way?
for (int i=thead.size(); i<45;i++)
@ -571,7 +571,7 @@ class cart_output_plugin : public output_plugin
header_.aexpN = astart_;
header_.aexp0 = header_.aexpN;
header_.amplt = 0.0; // Amplitude of density fluctuations
header_.astep = 0.0; //cf.getValue<double>("output","astep");
header_.astep = 0.0; //cf.get_value<double>("output","astep");
ptf_.astep=header_.astep; // to write pt file
header_.istep = 0; // step (=0 in IC)
header_.partw = 0.0; // mass of highest res particle. SEE BELOW
@ -603,10 +603,10 @@ class cart_output_plugin : public output_plugin
}
header_.Nseed = 0; // random number used ( 0 for MUSIC? or set the random number used in the lowest level?)
header_.Omb0 = cf.getValue<double>("cosmology","Omega_b");; // this parameter only appears in header in hydro runs
header_.Om0 = cf.getValue<double>("cosmology","Omega_m"); //Omega_m
header_.Oml0 = cf.getValue<double>("cosmology","Omega_L"); //Omega_L
header_.hubble = cf.getValue<double>("cosmology","H0")/100.; //hubble constant h=H/100
header_.Omb0 = cf.get_value<double>("cosmology","Omega_b");; // this parameter only appears in header in hydro runs
header_.Om0 = cf.get_value<double>("cosmology","Omega_m"); //Omega_m
header_.Oml0 = cf.get_value<double>("cosmology","Omega_L"); //Omega_L
header_.hubble = cf.get_value<double>("cosmology","H0")/100.; //hubble constant h=H/100
header_.Wp5 = 0.0;
header_.Ocurv = 1.0 - header_.Oml0 - header_.Om0;
@ -627,9 +627,9 @@ class cart_output_plugin : public output_plugin
{
header_.extras[i] = 0.0; //extras[20-99]
}
header_.extras[13] = cf.getValueSafe<double>("cosmology","Omega_b",0.0);
header_.extras[14] = cf.getValue<double>("cosmology","sigma_8");
header_.extras[15] = cf.getValue<double>("cosmology","nspec"); //Slope of the Power spectrum
header_.extras[13] = cf.get_value_safe<double>("cosmology","Omega_b",0.0);
header_.extras[14] = cf.get_value<double>("cosmology","sigma_8");
header_.extras[15] = cf.get_value<double>("cosmology","nspec"); //Slope of the Power spectrum
header_.magic1 = 0.1234f ;
header_.DelDC = 0.0;
header_.abox = astart_;
@ -638,9 +638,9 @@ class cart_output_plugin : public output_plugin
//header_.extras[NFILL-2] = Tini;
header_.extras[NFILL-1] = cf.getValue<double>("setup","boxlength");
header_.extras[NFILL-1] = cf.get_value<double>("setup","boxlength");
LOGINFO("CART : done header info.");
music::ilog.Print("CART : done header info.");
}

64
src/plugins/output_enzo.cc
View file

@ -314,11 +314,11 @@ public:
throw std::runtime_error("Error in enzo_output_plugin!");
}
bool bhave_hydro = cf_.getValue<bool>("setup", "baryons");
bool align_top = cf.getValueSafe<bool>("setup", "align_top", false);
bool bhave_hydro = cf_.get_value<bool>("setup", "baryons");
bool align_top = cf.get_value_safe<bool>("setup", "align_top", false);
if (!align_top)
LOGWARN("Old ENZO versions may require \'align_top=true\'!");
music::wlog.Print("Old ENZO versions may require \'align_top=true\'!");
the_sim_header.dimensions.push_back(1 << levelmin_);
the_sim_header.dimensions.push_back(1 << levelmin_);
@ -328,18 +328,18 @@ public:
the_sim_header.offset.push_back(0);
the_sim_header.offset.push_back(0);
the_sim_header.a_start = 1.0 / (1.0 + cf.getValue<double>("setup", "zstart"));
the_sim_header.dx = cf.getValue<double>("setup", "boxlength") / the_sim_header.dimensions[0] / (cf.getValue<double>("cosmology", "H0") * 0.01); // not sure?!?
the_sim_header.h0 = cf.getValue<double>("cosmology", "H0") * 0.01;
the_sim_header.a_start = 1.0 / (1.0 + cf.get_value<double>("setup", "zstart"));
the_sim_header.dx = cf.get_value<double>("setup", "boxlength") / the_sim_header.dimensions[0] / (cf.get_value<double>("cosmology", "H0") * 0.01); // not sure?!?
the_sim_header.h0 = cf.get_value<double>("cosmology", "H0") * 0.01;
if (bhave_hydro)
the_sim_header.omega_b = cf.getValue<double>("cosmology", "Omega_b");
the_sim_header.omega_b = cf.get_value<double>("cosmology", "Omega_b");
else
the_sim_header.omega_b = 0.0;
the_sim_header.omega_m = cf.getValue<double>("cosmology", "Omega_m");
the_sim_header.omega_v = cf.getValue<double>("cosmology", "Omega_L");
the_sim_header.vfact = cf.getValue<double>("cosmology", "vfact") * the_sim_header.h0; //.. need to multiply by h, ENZO wants this factor for non h-1 units
the_sim_header.omega_m = cf.get_value<double>("cosmology", "Omega_m");
the_sim_header.omega_v = cf.get_value<double>("cosmology", "Omega_L");
the_sim_header.vfact = cf.get_value<double>("cosmology", "vfact") * the_sim_header.h0; //.. need to multiply by h, ENZO wants this factor for non h-1 units
}
~enzo_output_plugin()
@ -353,8 +353,8 @@ public:
void write_dm_density(const grid_hierarchy &gh)
{ /* write the parameter file data */
bool bhave_hydro = cf_.getValue<bool>("setup", "baryons");
double refine_region_fraction = cf_.getValueSafe<double>("output", "enzo_refine_region_fraction", 0.8);
bool bhave_hydro = cf_.get_value<bool>("setup", "baryons");
double refine_region_fraction = cf_.get_value_safe<double>("output", "enzo_refine_region_fraction", 0.8);
char filename[256];
unsigned nbase = (unsigned)pow(2, levelmin_);
@ -412,9 +412,9 @@ public:
<< "CosmologyOmegaMatterNow = " << the_sim_header.omega_m << "\n"
<< "CosmologyOmegaLambdaNow = " << the_sim_header.omega_v << "\n"
<< "CosmologyHubbleConstantNow = " << the_sim_header.h0 << " // in 100 km/s/Mpc\n"
<< "CosmologyComovingBoxSize = " << cf_.getValue<double>("setup", "boxlength") << " // in Mpc/h\n"
<< "CosmologyComovingBoxSize = " << cf_.get_value<double>("setup", "boxlength") << " // in Mpc/h\n"
<< "CosmologyMaxExpansionRate = 0.015 // maximum allowed delta(a)/a\n"
<< "CosmologyInitialRedshift = " << cf_.getValue<double>("setup", "zstart") << " //\n"
<< "CosmologyInitialRedshift = " << cf_.get_value<double>("setup", "zstart") << " //\n"
<< "CosmologyFinalRedshift = 0 //\n"
<< "GravitationalConstant = 1 // this must be true for cosmology\n"
<< "#\n"
@ -528,31 +528,31 @@ public:
double h = 1.0 / (1 << levelmin_);
double shift[3];
shift[0] = -(double)cf_.getValue<int>("setup", "shift_x") * h;
shift[1] = -(double)cf_.getValue<int>("setup", "shift_y") * h;
shift[2] = -(double)cf_.getValue<int>("setup", "shift_z") * h;
shift[0] = -(double)cf_.get_value<int>("setup", "shift_x") * h;
shift[1] = -(double)cf_.get_value<int>("setup", "shift_y") * h;
shift[2] = -(double)cf_.get_value<int>("setup", "shift_z") * h;
if (gh.levelmin() != gh.levelmax())
{
LOGINFO("Global density extrema: ");
LOGINFO(" minimum: delta=%f at (%f,%f,%f) (level=%d)", rhomin, loc_rhomin[0], loc_rhomin[1], loc_rhomin[2], lvl_rhomin);
LOGINFO(" shifted back at (%f,%f,%f)", loc_rhomin[0] + shift[0], loc_rhomin[1] + shift[1], loc_rhomin[2] + shift[2]);
LOGINFO(" maximum: delta=%f at (%f,%f,%f) (level=%d)", rhomax, loc_rhomax[0], loc_rhomax[1], loc_rhomax[2], lvl_rhomax);
LOGINFO(" shifted back at (%f,%f,%f)", loc_rhomax[0] + shift[0], loc_rhomax[1] + shift[1], loc_rhomax[2] + shift[2]);
music::ilog.Print("Global density extrema: ");
music::ilog.Print(" minimum: delta=%f at (%f,%f,%f) (level=%d)", rhomin, loc_rhomin[0], loc_rhomin[1], loc_rhomin[2], lvl_rhomin);
music::ilog.Print(" shifted back at (%f,%f,%f)", loc_rhomin[0] + shift[0], loc_rhomin[1] + shift[1], loc_rhomin[2] + shift[2]);
music::ilog.Print(" maximum: delta=%f at (%f,%f,%f) (level=%d)", rhomax, loc_rhomax[0], loc_rhomax[1], loc_rhomax[2], lvl_rhomax);
music::ilog.Print(" shifted back at (%f,%f,%f)", loc_rhomax[0] + shift[0], loc_rhomax[1] + shift[1], loc_rhomax[2] + shift[2]);
LOGINFO("Density extrema on finest level: ");
LOGINFO(" minimum: delta=%f at (%f,%f,%f)", rhomin_lm, loc_rhomin_lm[0], loc_rhomin_lm[1], loc_rhomin_lm[2]);
LOGINFO(" shifted back at (%f,%f,%f)", loc_rhomin_lm[0] + shift[0], loc_rhomin_lm[1] + shift[1], loc_rhomin_lm[2] + shift[2]);
LOGINFO(" maximum: delta=%f at (%f,%f,%f)", rhomax_lm, loc_rhomax_lm[0], loc_rhomax_lm[1], loc_rhomax_lm[2]);
LOGINFO(" shifted back at (%f,%f,%f)", loc_rhomax_lm[0] + shift[0], loc_rhomax_lm[1] + shift[1], loc_rhomax_lm[2] + shift[2]);
music::ilog.Print("Density extrema on finest level: ");
music::ilog.Print(" minimum: delta=%f at (%f,%f,%f)", rhomin_lm, loc_rhomin_lm[0], loc_rhomin_lm[1], loc_rhomin_lm[2]);
music::ilog.Print(" shifted back at (%f,%f,%f)", loc_rhomin_lm[0] + shift[0], loc_rhomin_lm[1] + shift[1], loc_rhomin_lm[2] + shift[2]);
music::ilog.Print(" maximum: delta=%f at (%f,%f,%f)", rhomax_lm, loc_rhomax_lm[0], loc_rhomax_lm[1], loc_rhomax_lm[2]);
music::ilog.Print(" shifted back at (%f,%f,%f)", loc_rhomax_lm[0] + shift[0], loc_rhomax_lm[1] + shift[1], loc_rhomax_lm[2] + shift[2]);
}
else
{
LOGINFO("Global density extrema: ");
LOGINFO(" minimum: delta=%f at (%f,%f,%f)", rhomin, loc_rhomin[0], loc_rhomin[1], loc_rhomin[2]);
LOGINFO(" shifted back at (%f,%f,%f)", loc_rhomin[0] + shift[0], loc_rhomin[1] + shift[1], loc_rhomin[2] + shift[2]);
LOGINFO(" maximum: delta=%f at (%f,%f,%f)", rhomax, loc_rhomax[0], loc_rhomax[1], loc_rhomax[2]);
LOGINFO(" shifted back at (%f,%f,%f)", loc_rhomax[0] + shift[0], loc_rhomax[1] + shift[1], loc_rhomax[2] + shift[2]);
music::ilog.Print("Global density extrema: ");
music::ilog.Print(" minimum: delta=%f at (%f,%f,%f)", rhomin, loc_rhomin[0], loc_rhomin[1], loc_rhomin[2]);
music::ilog.Print(" shifted back at (%f,%f,%f)", loc_rhomin[0] + shift[0], loc_rhomin[1] + shift[1], loc_rhomin[2] + shift[2]);
music::ilog.Print(" maximum: delta=%f at (%f,%f,%f)", rhomax, loc_rhomax[0], loc_rhomax[1], loc_rhomax[2]);
music::ilog.Print(" shifted back at (%f,%f,%f)", loc_rhomax[0] + shift[0], loc_rhomax[1] + shift[1], loc_rhomax[2] + shift[2]);
}
}

2727
src/plugins/output_gadget2.cc
View file

File diff suppressed because it is too large Load diff

104
src/plugins/output_gadget2_2comp.cc
View file

@ -9,7 +9,7 @@
*/
#include <fstream>
#include "log.hh"
#include "logger.hh"
#include "output.hh"
#include "mg_interp.hh"
#include "mesh.hh"
@ -120,8 +120,8 @@ protected:
expected = ((unsigned long long) npart*(unsigned long long)sizeof(T_store));
if( blk != expected )
{
LOGERR("Internal consistency error in gadget2 output plug-in, open_and_check");
LOGERR("Expected %d particles (%lld bytes) in temp file %s but found %lld",npart, expected ,ffname.c_str(), blk);
music::elog.Print("Internal consistency error in gadget2 output plug-in, open_and_check");
music::elog.Print("Expected %d particles (%lld bytes) in temp file %s but found %lld",npart, expected ,ffname.c_str(), blk);
//throw std::runtime_error("Internal consistency error in gadget2 output plug-in");
}
@ -138,7 +138,7 @@ protected:
if( !this->good() )
{
LOGERR("Could not open buffer file in gadget2 output plug-in");
music::elog.Print("Could not open buffer file in gadget2 output plug-in");
throw std::runtime_error("Could not open buffer file in gadget2 output plug-in");
}
@ -146,8 +146,8 @@ protected:
if( blk != npart*sizeof(T_store) )
{
LOGERR("Internal consistency error in gadget2 output plug-in");
LOGERR("Expected %ld bytes in temp file but found %ld",npart*sizeof(T_store),blk);
music::elog.Print("Internal consistency error in gadget2 output plug-in");
music::elog.Print("Expected %ld bytes in temp file but found %ld",npart*sizeof(T_store),blk);
throw std::runtime_error("Internal consistency error in gadget2 output plug-in");
}
@ -166,7 +166,7 @@ protected:
if( !this->good() )
{
LOGERR("Could not open buffer file \'%s\' in gadget2 output plug-in",fname.c_str());
music::elog.Print("Could not open buffer file \'%s\' in gadget2 output plug-in",fname.c_str());
throw std::runtime_error("Could not open buffer file in gadget2 output plug-in");
}
@ -174,8 +174,8 @@ protected:
if( blk != npart*sizeof(T_store) )
{
LOGERR("Internal consistency error in gadget2 output plug-in");
LOGERR("Expected %ld bytes in temp file but found %ld",npart*sizeof(T_store),blk);
music::elog.Print("Internal consistency error in gadget2 output plug-in");
music::elog.Print("Expected %ld bytes in temp file but found %ld",npart*sizeof(T_store),blk);
throw std::runtime_error("Internal consistency error in gadget2 output plug-in");
}
@ -193,7 +193,7 @@ protected:
if( !this->good() )
{
LOGERR("Could not open buffer file in gadget2 output plug-in");
music::elog.Print("Could not open buffer file in gadget2 output plug-in");
throw std::runtime_error("Could not open buffer file in gadget2 output plug-in");
}
@ -201,8 +201,8 @@ protected:
if( blk != npart*sizeof(T_store) )
{
LOGERR("Internal consistency error in gadget2 output plug-in");
LOGERR("Expected %ld bytes in temp file but found %ld",npart*sizeof(T_store),blk);
music::elog.Print("Internal consistency error in gadget2 output plug-in");
music::elog.Print("Expected %ld bytes in temp file but found %ld",npart*sizeof(T_store),blk);
throw std::runtime_error("Internal consistency error in gadget2 output plug-in");
}
@ -225,7 +225,7 @@ protected:
if( !this->good() )
{
LOGERR("Could not open buffer file \'%s\' in gadget2 output plug-in",fname.c_str());
music::elog.Print("Could not open buffer file \'%s\' in gadget2 output plug-in",fname.c_str());
throw std::runtime_error("Could not open buffer file in gadget2 output plug-in");
}
@ -233,8 +233,8 @@ protected:
if( blk != npart*sizeof(T_store) )
{
LOGERR("Internal consistency error in gadget2 output plug-in");
LOGERR("Expected %ld bytes in temp file but found %ld",npart*sizeof(T_store),blk);
music::elog.Print("Internal consistency error in gadget2 output plug-in");
music::elog.Print("Expected %ld bytes in temp file but found %ld",npart*sizeof(T_store),blk);
throw std::runtime_error("Internal consistency error in gadget2 output plug-in");
}
@ -290,7 +290,7 @@ protected:
n2read = std::min(block_buf_size_,npleft);
if( header_.npart[5] > 0 )
LOGERR("Multi-resolution setup not supported for 2comp hack");
music::elog.Print("Multi-resolution setup not supported for 2comp hack");
std::cout << " - Writing " << nptot << " particles to Gadget file...\n"
<< " type 1 : " << header_.npart[1] << "\n"
@ -333,7 +333,7 @@ protected:
if( nptot >= 1ul<<32 )
{
bneed_long_ids = true;
LOGWARN("Need long particle IDs, make sure to enable in Gadget!");
music::wlog.Print("Need long particle IDs, make sure to enable in Gadget!");
}
@ -722,15 +722,15 @@ public:
gadget2_2comp_output_plugin( config_file& cf )
: output_plugin( cf )//, ofs_( fname_.c_str(), std::ios::binary|std::ios::trunc )
{
block_buf_size_ = cf_.getValueSafe<unsigned>("output","gadget_blksize",2*1048576);
block_buf_size_ = cf_.get_value_safe<unsigned>("output","gadget_blksize",2*1048576);
//... ensure that everyone knows we want to do SPH
cf.insertValue("setup","do_SPH","yes");
cf.insert_value("setup","do_SPH","yes");
//bbndparticles_ = !cf_.getValueSafe<bool>("output","gadget_nobndpart",false);
//bbndparticles_ = !cf_.get_value_safe<bool>("output","gadget_nobndpart",false);
npartmax_ = 1<<30;
nfiles_ = cf.getValueSafe<unsigned>("output","gadget_num_files",1);
nfiles_ = cf.get_value_safe<unsigned>("output","gadget_num_files",1);
@ -744,7 +744,7 @@ public:
ofs_.open(ffname, std::ios::binary|std::ios::trunc );
if(!ofs_.good())
{
LOGERR("gadget-2 output plug-in could not open output file \'%s\' for writing!",ffname);
music::elog.Print("gadget-2 output plug-in could not open output file \'%s\' for writing!",ffname);
throw std::runtime_error(std::string("gadget-2 output plug-in could not open output file \'")+std::string(ffname)+"\' for writing!\n");
}
ofs_.close();
@ -753,7 +753,7 @@ public:
ofs_.open(fname_.c_str(), std::ios::binary|std::ios::trunc );
if(!ofs_.good())
{
LOGERR("gadget-2 output plug-in could not open output file \'%s\' for writing!",fname_.c_str());
music::elog.Print("gadget-2 output plug-in could not open output file \'%s\' for writing!",fname_.c_str());
throw std::runtime_error(std::string("gadget-2 output plug-in could not open output file \'")+fname_+"\' for writing!\n");
}
ofs_.close();
@ -776,29 +776,29 @@ public:
header_.mass[i] = 0.0;
}
YHe_ = cf.getValueSafe<double>("cosmology","YHe",0.248);
gamma_ = cf.getValueSafe<double>("cosmology","gamma",5.0/3.0);
YHe_ = cf.get_value_safe<double>("cosmology","YHe",0.248);
gamma_ = cf.get_value_safe<double>("cosmology","gamma",5.0/3.0);
do_baryons_ = cf.getValueSafe<bool>("setup","baryons",false);
omegab_ = cf.getValueSafe<double>("cosmology","Omega_b",0.045);
do_baryons_ = cf.get_value_safe<bool>("setup","baryons",false);
omegab_ = cf.get_value_safe<double>("cosmology","Omega_b",0.045);
//... write displacements in kpc/h rather than Mpc/h?
kpcunits_ = cf.getValueSafe<bool>("output","gadget_usekpc",false);
kpcunits_ = cf.get_value_safe<bool>("output","gadget_usekpc",false);
do_glass_ = cf.getValueSafe<bool>("output","glass", false);
do_glass_ = cf.get_value_safe<bool>("output","glass", false);
if( do_glass_ )
{
LOGINFO("Will use provided glass rather than Cartesian mesh for particle placement.");
music::ilog.Print("Will use provided glass rather than Cartesian mesh for particle placement.");
fname_glass_cdm_ = cf.getValue<std::string>("output","glass_file_cdm");
fname_glass_cdm_ = cf.get_value<std::string>("output","glass_file_cdm");
if( do_baryons_ )
fname_glass_baryon_ = fname_glass_cdm_;//cf.getValue<std::string>("output","glass_file_baryon");
fname_glass_baryon_ = fname_glass_cdm_;//cf.get_value<std::string>("output","glass_file_baryon");
}
//... set time ......................................................
header_.redshift = cf.getValue<double>("setup","zstart");
header_.redshift = cf.get_value<double>("setup","zstart");
header_.time = 1.0/(1.0+header_.redshift);
//... SF flags
@ -808,13 +808,13 @@ public:
//...
header_.num_files = nfiles_;
header_.BoxSize = cf.getValue<double>("setup","boxlength");
header_.Omega0 = cf.getValue<double>("cosmology","Omega_m");
header_.BoxSize = cf.get_value<double>("setup","boxlength");
header_.Omega0 = cf.get_value<double>("cosmology","Omega_m");
omegam_ = header_.Omega0;
omegac_ = omegam_ - omegab_;
header_.OmegaLambda = cf.getValue<double>("cosmology","Omega_L");
header_.HubbleParam = cf.getValue<double>("cosmology","H0");
header_.OmegaLambda = cf.get_value<double>("cosmology","Omega_L");
header_.HubbleParam = cf.get_value<double>("cosmology","H0");
header_.flag_stellarage = 0;
header_.flag_metals = 0;
@ -930,16 +930,16 @@ public:
//... determine if we need to shift the coordinates back
double *shift = NULL;
if( cf_.getValueSafe<bool>("output","shift_back",false ) )
if( cf_.get_value_safe<bool>("output","shift_back",false ) )
{
if( coord == 0 )
std::cout << " - gadget2 output plug-in will shift particle positions back...\n";
double h = 1.0/(1<<levelmin_);
shift = new double[3];
shift[0] = -(double)cf_.getValue<int>( "setup", "shift_x" )*h;
shift[1] = -(double)cf_.getValue<int>( "setup", "shift_y" )*h;
shift[2] = -(double)cf_.getValue<int>( "setup", "shift_z" )*h;
shift[0] = -(double)cf_.get_value<int>( "setup", "shift_x" )*h;
shift[1] = -(double)cf_.get_value<int>( "setup", "shift_y" )*h;
shift[2] = -(double)cf_.get_value<int>( "setup", "shift_z" )*h;
}
size_t npart = npfine+npcoarse;
@ -1023,7 +1023,7 @@ public:
std::ifstream ofg( fname_glass_cdm_.c_str(), std::ios::binary );
if( !ofg.good() )
LOGERR("could not open glass input file \'%s\'",fname_glass_cdm_.c_str());
music::elog.Print("could not open glass input file \'%s\'",fname_glass_cdm_.c_str());
io_header glasshead;
unsigned blksz;
@ -1036,7 +1036,7 @@ public:
//size_t nreq = gh.size(gh.levelmax(), 0)*gh.size(gh.levelmax(), 1)*gh.size(gh.levelmax(), 2);
/*if( nreq != (size_t)glasshead.npart[1] )
{
LOGERR("glass file contains %d particles, but should contain %ld",glasshead.npart[1],nreq);
music::elog.Print("glass file contains %d particles, but should contain %ld",glasshead.npart[1],nreq);
throw std::runtime_error("glass file does not contain the right amount of particles");
}*/
@ -1212,7 +1212,7 @@ public:
std::ifstream ofg( fname_glass_cdm_.c_str(), std::ios::binary );
if( !ofg.good() )
LOGERR("could not open glass input file \'%s\'",fname_glass_cdm_.c_str());
music::elog.Print("could not open glass input file \'%s\'",fname_glass_cdm_.c_str());
io_header glasshead;
unsigned blksz;
@ -1375,7 +1375,7 @@ public:
std::ifstream ofg( fname_glass_baryon_.c_str(), std::ios::binary );
if( !ofg.good() )
LOGERR("could not open glass input file \'%s\'",fname_glass_cdm_.c_str());
music::elog.Print("could not open glass input file \'%s\'",fname_glass_cdm_.c_str());
io_header glasshead;
unsigned blksz;
@ -1388,7 +1388,7 @@ public:
//size_t nreq = gh.size(gh.levelmax(), 0)*gh.size(gh.levelmax(), 1)*gh.size(gh.levelmax(), 2);
/*if( nreq != (size_t)glasshead.npart[1] )
{
LOGERR("glass file contains %d particles, but should contain %ld",glasshead.npart[1],nreq);
music::elog.Print("glass file contains %d particles, but should contain %ld",glasshead.npart[1],nreq);
throw std::runtime_error("glass file does not contain the right amount of particles");
}*/
@ -1478,16 +1478,16 @@ public:
//... determine if we need to shift the coordinates back
double *shift = NULL;
if( cf_.getValueSafe<bool>("output","shift_back",false ) )
if( cf_.get_value_safe<bool>("output","shift_back",false ) )
{
if( coord == 0 )
std::cout << " - gadget2 output plug-in will shift particle positions back...\n";
double h = 1.0/(1<<levelmin_);
shift = new double[3];
shift[0] = -(double)cf_.getValue<int>( "setup", "shift_x" )*h;
shift[1] = -(double)cf_.getValue<int>( "setup", "shift_y" )*h;
shift[2] = -(double)cf_.getValue<int>( "setup", "shift_z" )*h;
shift[0] = -(double)cf_.get_value<int>( "setup", "shift_x" )*h;
shift[1] = -(double)cf_.get_value<int>( "setup", "shift_y" )*h;
shift[2] = -(double)cf_.get_value<int>( "setup", "shift_z" )*h;
}
unsigned long long npart = npfine;
@ -1575,7 +1575,7 @@ public:
std::ifstream ofg( fname_glass_baryon_.c_str(), std::ios::binary );
if( !ofg.good() )
LOGERR("could not open glass input file \'%s\'",fname_glass_cdm_.c_str());
music::elog.Print("could not open glass input file \'%s\'",fname_glass_cdm_.c_str());
io_header glasshead;
unsigned blksz;
@ -1588,7 +1588,7 @@ public:
//size_t nreq = gh.size(gh.levelmax(), 0)*gh.size(gh.levelmax(), 1)*gh.size(gh.levelmax(), 2);
/*if( nreq != (size_t)glasshead.npart[1] )
{
LOGERR("glass file contains %d particles, but should contain %ld",glasshead.npart[1],nreq);
music::elog.Print("glass file contains %d particles, but should contain %ld",glasshead.npart[1],nreq);
throw std::runtime_error("glass file does not contain the right amount of particles");
}*/

102
src/plugins/output_gadget_tetmesh.cc
View file

@ -16,7 +16,7 @@
#include <cstring>
#include <algorithm>
#include "log.hh"
#include "logger.hh"
#include "region_generator.hh"
#include "output.hh"
#include "mg_interp.hh"
@ -520,8 +520,8 @@ protected:
ifs.read( (char*)&blk, sizeof(size_t) );
if( blk != blksize )
{
LOGERR("Internal consistency error in gadget2 output plug-in");
LOGERR("Expected %ld bytes in temp file but found %ld",blksize,blk);
music::elog.Print("Internal consistency error in gadget2 output plug-in");
music::elog.Print("Expected %ld bytes in temp file but found %ld",blksize,blk);
throw std::runtime_error("Internal consistency error in gadget2 output plug-in");
}
ifs.seekg( offset, std::ios::cur );
@ -539,7 +539,7 @@ protected:
if( !this->good() )
{
LOGERR("Could not open buffer file in gadget2 output plug-in");
music::elog.Print("Could not open buffer file in gadget2 output plug-in");
throw std::runtime_error("Could not open buffer file in gadget2 output plug-in");
}
@ -547,8 +547,8 @@ protected:
if( blk != blksize )
{
LOGERR("Internal consistency error in gadget2 output plug-in");
LOGERR("Expected %ld bytes in temp file but found %ld",blksize,blk);
music::elog.Print("Internal consistency error in gadget2 output plug-in");
music::elog.Print("Expected %ld bytes in temp file but found %ld",blksize,blk);
throw std::runtime_error("Internal consistency error in gadget2 output plug-in");
}
@ -567,7 +567,7 @@ protected:
if( !this->good() )
{
LOGERR("Could not open buffer file \'%s\' in gadget2 output plug-in",fname.c_str());
music::elog.Print("Could not open buffer file \'%s\' in gadget2 output plug-in",fname.c_str());
throw std::runtime_error("Could not open buffer file in gadget2 output plug-in");
}
@ -575,8 +575,8 @@ protected:
if( blk != blksize )
{
LOGERR("Internal consistency error in gadget2 output plug-in");
LOGERR("Expected %ld bytes in temp file but found %ld",blksize,blk);
music::elog.Print("Internal consistency error in gadget2 output plug-in");
music::elog.Print("Expected %ld bytes in temp file but found %ld",blksize,blk);
throw std::runtime_error("Internal consistency error in gadget2 output plug-in");
}
@ -655,7 +655,7 @@ protected:
size_t curr_block_buf_size = block_buf_size_;
size_t idcount = 0;
LOGWARN("Need long particle IDs, will write 64bit, make sure to enable in Gadget!");
music::wlog.Print("Need long particle IDs, will write 64bit, make sure to enable in Gadget!");
for( unsigned ifile=0; ifile<nfiles_; ++ifile )
{
@ -917,22 +917,22 @@ public:
gadget_tetmesh_output_plugin( config_file& cf )
: output_plugin( cf )
{
block_buf_size_ = cf_.getValueSafe<unsigned>("output","gadget_blksize",1048576);
block_buf_size_ = cf_.get_value_safe<unsigned>("output","gadget_blksize",1048576);
//... ensure that everyone knows we want to do SPH
cf.insertValue("setup","do_SPH","yes");
cf.insert_value("setup","do_SPH","yes");
//bbndparticles_ = !cf_.getValueSafe<bool>("output","gadget_nobndpart",false);
//bbndparticles_ = !cf_.get_value_safe<bool>("output","gadget_nobndpart",false);
npartmax_ = 1<<30;
nfiles_ = cf.getValueSafe<unsigned>("output","gadget_num_files",1);
nfiles_ = cf.get_value_safe<unsigned>("output","gadget_num_files",1);
blongids_ = cf.getValueSafe<bool>("output","gadget_longids",false);
blongids_ = cf.get_value_safe<bool>("output","gadget_longids",false);
shift_halfcell_ = cf.getValueSafe<bool>("output","gadget_cell_centered",false);
shift_halfcell_ = cf.get_value_safe<bool>("output","gadget_cell_centered",false);
//if( nfiles_ < (int)ceil((double)npart/(double)npartmax_) )
// LOGWARN("Should use more files.");
// music::wlog.Print("Should use more files.");
if (nfiles_ > 1 )
{
@ -943,7 +943,7 @@ public:
ofs_.open(ffname, std::ios::binary|std::ios::trunc );
if(!ofs_.good())
{
LOGERR("gadget-2 output plug-in could not open output file \'%s\' for writing!",ffname);
music::elog.Print("gadget-2 output plug-in could not open output file \'%s\' for writing!",ffname);
throw std::runtime_error(std::string("gadget-2 output plug-in could not open output file \'")+std::string(ffname)+"\' for writing!\n");
}
ofs_.close();
@ -952,7 +952,7 @@ public:
ofs_.open(fname_.c_str(), std::ios::binary|std::ios::trunc );
if(!ofs_.good())
{
LOGERR("gadget-2 output plug-in could not open output file \'%s\' for writing!",fname_.c_str());
music::elog.Print("gadget-2 output plug-in could not open output file \'%s\' for writing!",fname_.c_str());
throw std::runtime_error(std::string("gadget-2 output plug-in could not open output file \'")+fname_+"\' for writing!\n");
}
ofs_.close();
@ -978,30 +978,30 @@ public:
header_.mass[i] = 0.0;
}
YHe_ = cf.getValueSafe<double>("cosmology","YHe",0.248);
gamma_ = cf.getValueSafe<double>("cosmology","gamma",5.0/3.0);
YHe_ = cf.get_value_safe<double>("cosmology","YHe",0.248);
gamma_ = cf.get_value_safe<double>("cosmology","gamma",5.0/3.0);
do_baryons_ = cf.getValueSafe<bool>("setup","baryons",false);
omegab_ = cf.getValueSafe<double>("cosmology","Omega_b",0.045);
do_baryons_ = cf.get_value_safe<bool>("setup","baryons",false);
omegab_ = cf.get_value_safe<double>("cosmology","Omega_b",0.045);
//... write displacements in kpc/h rather than Mpc/h?
kpcunits_ = cf.getValueSafe<bool>("output","gadget_usekpc",false);
msolunits_ = cf.getValueSafe<bool>("output","gadget_usemsol",false);
kpcunits_ = cf.get_value_safe<bool>("output","gadget_usekpc",false);
msolunits_ = cf.get_value_safe<bool>("output","gadget_usemsol",false);
blagrangeids_as_vertids_ = cf.getValueSafe<bool>("output","gadget_lagrangevertid",true);
blagrangeids_as_vertids_ = cf.get_value_safe<bool>("output","gadget_lagrangevertid",true);
/*bndparticletype_ = cf.getValueSafe<unsigned>("output","gadget_coarsetype",5);
/*bndparticletype_ = cf.get_value_safe<unsigned>("output","gadget_coarsetype",5);
if( bndparticletype_ == 0 || bndparticletype_ == 1 || bndparticletype_ == 4 ||
bndparticletype_ > 5 )
{
LOGERR("Coarse particles cannot be of Gadget particle type %d in output plugin.", bndparticletype_);
music::elog.Print("Coarse particles cannot be of Gadget particle type %d in output plugin.", bndparticletype_);
throw std::runtime_error("Specified illegal Gadget particle type for coarse particles");
}*/
//... set time ......................................................
header_.redshift = cf.getValue<double>("setup","zstart");
header_.redshift = cf.get_value<double>("setup","zstart");
header_.time = 1.0/(1.0+header_.redshift);
//... SF flags
@ -1011,10 +1011,10 @@ public:
//...
header_.num_files = nfiles_;//1;
header_.BoxSize = cf.getValue<double>("setup","boxlength");
header_.Omega0 = cf.getValue<double>("cosmology","Omega_m");
header_.OmegaLambda = cf.getValue<double>("cosmology","Omega_L");
header_.HubbleParam = cf.getValue<double>("cosmology","H0")/100.0;
header_.BoxSize = cf.get_value<double>("setup","boxlength");
header_.Omega0 = cf.get_value<double>("cosmology","Omega_m");
header_.OmegaLambda = cf.get_value<double>("cosmology","Omega_L");
header_.HubbleParam = cf.get_value<double>("cosmology","H0")/100.0;
header_.flag_stellarage = 0;
header_.flag_metals = 0;
@ -1075,7 +1075,7 @@ public:
if( it==idmap.end() ){
foundall = false;
LOGERR("This should not happen : Lagrange ID %llu not found!", REMOVE_DECORATION_BITS( P[ip].get_vertex(i) ));
music::elog.Print("This should not happen : Lagrange ID %llu not found!", REMOVE_DECORATION_BITS( P[ip].get_vertex(i) ));
throw std::runtime_error("FATAL");
break;
}
@ -1125,13 +1125,13 @@ public:
if( num_p + newnum_p_per_split > num_p_alloc )
{
P = reinterpret_cast<particle*>( realloc( reinterpret_cast<void*>(P), (num_p_alloc+=num_p_realloc_blocksize)*sizeof(particle) ) );
LOGINFO("reallocated particle buffer. new size = %llu MBytes.", num_p_alloc * sizeof(particle)/1024/1024 );
music::ilog.Print("reallocated particle buffer. new size = %llu MBytes.", num_p_alloc * sizeof(particle)/1024/1024 );
}
split_lagrange_cube( ip );
}
delete_duplicates();
LOGINFO("refined tet mesh to level %d : now have %lld particles", ilevel+1, num_p );
music::ilog.Print("refined tet mesh to level %d : now have %lld particles", ilevel+1, num_p );
}
@ -1194,9 +1194,9 @@ public:
header_.npart[5] = num_p_t5;
header_.npartTotal[5] = num_p_t5;
LOGINFO(" active HR particles (type 1) : %llu", num_p_t1 );
LOGINFO(" active LR particles (type 5) : %llu", num_p_t5 );
LOGINFO(" passive particles (type 2) : %llu", num_p_t2 );
music::ilog.Print(" active HR particles (type 1) : %llu", num_p_t1 );
music::ilog.Print(" active LR particles (type 5) : %llu", num_p_t5 );
music::ilog.Print(" passive particles (type 2) : %llu", num_p_t2 );
// write all particle masses
@ -1240,14 +1240,14 @@ public:
if( nwritten != num_p )
{
LOGERR("Internal consistency error while writing temporary file for masses");
music::elog.Print("Internal consistency error while writing temporary file for masses");
throw std::runtime_error("Internal consistency error while writing temporary file for masses");
}
ofs_temp.write( (char *)&blksize, sizeof(size_t) );
if( ofs_temp.bad() ){
LOGERR("I/O error while writing temporary file for masses");
music::elog.Print("I/O error while writing temporary file for masses");
throw std::runtime_error("I/O error while writing temporary file for masses");
}
}
@ -1304,14 +1304,14 @@ public:
if( nwritten != 8*num_p )
{
LOGERR("Internal consistency error while writing temporary file for connectivities");
music::elog.Print("Internal consistency error while writing temporary file for connectivities");
throw std::runtime_error("Internal consistency error while writing temporary file for connectivities");
}
ofs_temp.write( (char *)&blksize, sizeof(size_t) );
if( ofs_temp.bad() ){
LOGERR("I/O error while writing temporary file for connectivities");
music::elog.Print("I/O error while writing temporary file for connectivities");
throw std::runtime_error("I/O error while writing temporary file for connectivities");
}
}
@ -1350,14 +1350,14 @@ public:
if( nwritten != num_p )
{
LOGERR("Internal consistency error while writing temporary file for masses");
music::elog.Print("Internal consistency error while writing temporary file for masses");
throw std::runtime_error("Internal consistency error while writing temporary file for masses");
}
ofs_temp.write( (char *)&blksize, sizeof(size_t) );
if( ofs_temp.bad() ){
LOGERR("I/O error while writing temporary file for masses");
music::elog.Print("I/O error while writing temporary file for masses");
throw std::runtime_error("I/O error while writing temporary file for masses");
}
}
@ -1396,14 +1396,14 @@ public:
if( nwritten != num_p )
{
LOGERR("Internal consistency error while writing temporary file for masses");
music::elog.Print("Internal consistency error while writing temporary file for masses");
throw std::runtime_error("Internal consistency error while writing temporary file for masses");
}
ofs_temp.write( (char *)&blksize, sizeof(size_t) );
if( ofs_temp.bad() ){
LOGERR("I/O error while writing temporary file for masses");
music::elog.Print("I/O error while writing temporary file for masses");
throw std::runtime_error("I/O error while writing temporary file for masses");
}
}
@ -1460,14 +1460,14 @@ public:
if( nwritten != num_p )
{
LOGERR("Internal consistency error while writing temporary file for masses");
music::elog.Print("Internal consistency error while writing temporary file for masses");
throw std::runtime_error("Internal consistency error while writing temporary file for masses");
}
ofs_temp.write( (char *)&blksize, sizeof(size_t) );
if( ofs_temp.bad() ){
LOGERR("I/O error while writing temporary file for masses");
music::elog.Print("I/O error while writing temporary file for masses");
throw std::runtime_error("I/O error while writing temporary file for masses");
}
}
@ -1526,14 +1526,14 @@ public:
if( nwritten != num_p )
{
LOGERR("Internal consistency error while writing temporary file for velocities");
music::elog.Print("Internal consistency error while writing temporary file for velocities");
throw std::runtime_error("Internal consistency error while writing temporary file for vleocities");
}
ofs_temp.write( (char *)&blksize, sizeof(size_t) );
if( ofs_temp.bad() ){
LOGERR("I/O error while writing temporary file for velocities");
music::elog.Print("I/O error while writing temporary file for velocities");
throw std::runtime_error("I/O error while writing temporary file for velocities");
}
}

36
src/plugins/output_generic.cc
View file

@ -157,30 +157,30 @@ public:
double h = 1.0/(1<<levelmin_);
double shift[3];
shift[0] = -(double)cf_.getValue<int>( "setup", "shift_x" )*h;
shift[1] = -(double)cf_.getValue<int>( "setup", "shift_y" )*h;
shift[2] = -(double)cf_.getValue<int>( "setup", "shift_z" )*h;
shift[0] = -(double)cf_.get_value<int>( "setup", "shift_x" )*h;
shift[1] = -(double)cf_.get_value<int>( "setup", "shift_y" )*h;
shift[2] = -(double)cf_.get_value<int>( "setup", "shift_z" )*h;
if( gh.levelmin() != gh.levelmax() )
{
LOGINFO("Global density extrema: ");
LOGINFO(" minimum: delta=%f at (%f,%f,%f) (level=%d)",rhomin,loc_rhomin[0],loc_rhomin[1],loc_rhomin[2],lvl_rhomin);
LOGINFO(" shifted back at (%f,%f,%f)",loc_rhomin[0]+shift[0],loc_rhomin[1]+shift[1],loc_rhomin[2]+shift[2]);
LOGINFO(" maximum: delta=%f at (%f,%f,%f) (level=%d)",rhomax,loc_rhomax[0],loc_rhomax[1],loc_rhomax[2],lvl_rhomax);
LOGINFO(" shifted back at (%f,%f,%f)",loc_rhomax[0]+shift[0],loc_rhomax[1]+shift[1],loc_rhomax[2]+shift[2]);
music::ilog.Print("Global density extrema: ");
music::ilog.Print(" minimum: delta=%f at (%f,%f,%f) (level=%d)",rhomin,loc_rhomin[0],loc_rhomin[1],loc_rhomin[2],lvl_rhomin);
music::ilog.Print(" shifted back at (%f,%f,%f)",loc_rhomin[0]+shift[0],loc_rhomin[1]+shift[1],loc_rhomin[2]+shift[2]);
music::ilog.Print(" maximum: delta=%f at (%f,%f,%f) (level=%d)",rhomax,loc_rhomax[0],loc_rhomax[1],loc_rhomax[2],lvl_rhomax);
music::ilog.Print(" shifted back at (%f,%f,%f)",loc_rhomax[0]+shift[0],loc_rhomax[1]+shift[1],loc_rhomax[2]+shift[2]);
LOGINFO("Density extrema on finest level: ");
LOGINFO(" minimum: delta=%f at (%f,%f,%f)",rhomin_lm,loc_rhomin_lm[0],loc_rhomin_lm[1],loc_rhomin_lm[2]);
LOGINFO(" shifted back at (%f,%f,%f)",loc_rhomin_lm[0]+shift[0],loc_rhomin_lm[1]+shift[1],loc_rhomin_lm[2]+shift[2]);
LOGINFO(" maximum: delta=%f at (%f,%f,%f)",rhomax_lm,loc_rhomax_lm[0],loc_rhomax_lm[1],loc_rhomax_lm[2]);
LOGINFO(" shifted back at (%f,%f,%f)",loc_rhomax_lm[0]+shift[0],loc_rhomax_lm[1]+shift[1],loc_rhomax_lm[2]+shift[2]);
music::ilog.Print("Density extrema on finest level: ");
music::ilog.Print(" minimum: delta=%f at (%f,%f,%f)",rhomin_lm,loc_rhomin_lm[0],loc_rhomin_lm[1],loc_rhomin_lm[2]);
music::ilog.Print(" shifted back at (%f,%f,%f)",loc_rhomin_lm[0]+shift[0],loc_rhomin_lm[1]+shift[1],loc_rhomin_lm[2]+shift[2]);
music::ilog.Print(" maximum: delta=%f at (%f,%f,%f)",rhomax_lm,loc_rhomax_lm[0],loc_rhomax_lm[1],loc_rhomax_lm[2]);
music::ilog.Print(" shifted back at (%f,%f,%f)",loc_rhomax_lm[0]+shift[0],loc_rhomax_lm[1]+shift[1],loc_rhomax_lm[2]+shift[2]);
}else{
LOGINFO("Global density extrema: ");
LOGINFO(" minimum: delta=%f at (%f,%f,%f)",rhomin,loc_rhomin[0],loc_rhomin[1],loc_rhomin[2]);
LOGINFO(" shifted back at (%f,%f,%f)",loc_rhomin[0]+shift[0],loc_rhomin[1]+shift[1],loc_rhomin[2]+shift[2]);
LOGINFO(" maximum: delta=%f at (%f,%f,%f)",rhomax,loc_rhomax[0],loc_rhomax[1],loc_rhomax[2]);
LOGINFO(" shifted back at (%f,%f,%f)",loc_rhomax[0]+shift[0],loc_rhomax[1]+shift[1],loc_rhomax[2]+shift[2]);
music::ilog.Print("Global density extrema: ");
music::ilog.Print(" minimum: delta=%f at (%f,%f,%f)",rhomin,loc_rhomin[0],loc_rhomin[1],loc_rhomin[2]);
music::ilog.Print(" shifted back at (%f,%f,%f)",loc_rhomin[0]+shift[0],loc_rhomin[1]+shift[1],loc_rhomin[2]+shift[2]);
music::ilog.Print(" maximum: delta=%f at (%f,%f,%f)",rhomax,loc_rhomax[0],loc_rhomax[1],loc_rhomax[2]);
music::ilog.Print(" shifted back at (%f,%f,%f)",loc_rhomax[0]+shift[0],loc_rhomax[1]+shift[1],loc_rhomax[2]+shift[2]);
}
}

36
src/plugins/output_grafic2.cc
View file

@ -42,12 +42,12 @@ protected:
header loc_head;
double
boxlength = cf_.getValue<double>("setup","boxlength"),
H0 = cf_.getValue<double>("cosmology","H0"),
zstart = cf_.getValue<double>("setup","zstart"),
boxlength = cf_.get_value<double>("setup","boxlength"),
H0 = cf_.get_value<double>("cosmology","H0"),
zstart = cf_.get_value<double>("setup","zstart"),
astart = 1.0/(1.0+zstart),
omegam = cf_.getValue<double>("cosmology","Omega_m"),
omegaL = cf_.getValue<double>("cosmology","Omega_L");
omegam = cf_.get_value<double>("cosmology","Omega_m"),
omegaL = cf_.get_value<double>("cosmology","Omega_L");
loc_head.n1 = gh.get_grid(ilevel)->size(0);
loc_head.n2 = gh.get_grid(ilevel)->size(1);
@ -235,7 +235,7 @@ protected:
size_t nref;
nref = restrict_mask( n1, n2, n3, o1, o2, o3, n1c, n2c, n3c, &data[0], &data_coarse[0] );
LOGINFO("%f of cells on level %d are refined",(double)nref/(n1c*n2c*n3c),ilevel);
music::ilog.Print("%f of cells on level %d are refined",(double)nref/(n1c*n2c*n3c),ilevel);
sprintf(ff,"%s/level_%03d/ic_refmap",fname_.c_str(), ilevel );
std::ofstream ofs(ff,std::ios::binary|std::ios::trunc);
@ -328,7 +328,7 @@ protected:
unsigned naddref = 8; // initialize with settings for 10 additional levels of refinement
unsigned nexpand = (cf_.getValue<unsigned>("setup","padding")-1)/2;
unsigned nexpand = (cf_.get_value<unsigned>("setup","padding")-1)/2;
// -- AMR_PARAMS -- //
ofst << "&AMR_PARAMS\n"
@ -366,8 +366,8 @@ protected:
<< "/\n\n";
LOGINFO("The grafic2 output plug-in wrote the grid data to a partial");
LOGINFO(" RAMSES namelist file \'%s\'",fname_.c_str() );
music::ilog.Print("The grafic2 output plug-in wrote the grid data to a partial");
music::ilog.Print(" RAMSES namelist file \'%s\'",fname_.c_str() );
}
void write_ramses_namelist_old( const grid_hierarchy& gh )
@ -474,10 +474,10 @@ public:
}
bhavehydro_ = cf.getValue<bool>("setup","baryons");
//metal_floor_ = cf.getValueSafe<float>("output","ramses_metal_floor",1e-5);
passive_variable_index_ = cf.getValueSafe<int>("output","ramses_pvar_idx",1);
passive_variable_value_ = cf.getValueSafe<float>("output","ramses_pvar_val",1.0f);
bhavehydro_ = cf.get_value<bool>("setup","baryons");
//metal_floor_ = cf.get_value_safe<float>("output","ramses_metal_floor",1e-5);
passive_variable_index_ = cf.get_value_safe<int>("output","ramses_pvar_idx",1);
passive_variable_value_ = cf.get_value_safe<float>("output","ramses_pvar_val",1.0f);
}
/*~grafic2_output_plugin()
@ -487,7 +487,7 @@ public:
void write_dm_position( int coord, const grid_hierarchy& gh )
{
double
boxlength = cf_.getValue<double>("setup","boxlength");
boxlength = cf_.get_value<double>("setup","boxlength");
for(unsigned ilevel=levelmin_; ilevel<=levelmax_; ++ilevel )
{
@ -505,7 +505,7 @@ public:
void write_dm_velocity( int coord, const grid_hierarchy& gh )
{
double
boxlength = cf_.getValue<double>("setup","boxlength");
boxlength = cf_.get_value<double>("setup","boxlength");
for(unsigned ilevel=levelmin_; ilevel<=levelmax_; ++ilevel )
{
@ -523,7 +523,7 @@ public:
void write_gas_velocity( int coord, const grid_hierarchy& gh )
{
double
boxlength = cf_.getValue<double>("setup","boxlength");
boxlength = cf_.get_value<double>("setup","boxlength");
for(unsigned ilevel=levelmin_; ilevel<=levelmax_; ++ilevel )
{
@ -560,9 +560,9 @@ public:
if(! bhavehydro_ )
write_gas_density(gh);
if( cf_.getValueSafe<bool>("output","ramses_nml",true) )
if( cf_.get_value_safe<bool>("output","ramses_nml",true) )
write_ramses_namelist(gh);
else if( cf_.getValueSafe<bool>("output","ramses_old_nml",false) )
else if( cf_.get_value_safe<bool>("output","ramses_old_nml",false) )
write_ramses_namelist_old(gh);
if( gh.levelmin() != gh.levelmax() )

32
src/plugins/output_swift.cc
View file

@ -523,26 +523,26 @@ public:
{
// ensure that everyone knows we want to do SPH, implies: bsph=1, bbshift=1, decic_baryons=1
// -> instead of just writing gas densities (which are here ignored), the gas displacements are also written
cf.insertValue("setup", "do_SPH", "yes");
cf.insert_value("setup", "do_SPH", "yes");
// init header and config parameters
nPartTotal = std::vector<long long>(NTYPES, 0);
massTable = std::vector<double>(NTYPES, 0.0);
coarsePartType = cf.getValueSafe<unsigned>("output", "arepo_coarsetype", COARSE_DM_DEFAULT_PARTTYPE);
UnitLength_in_cm = cf.getValueSafe<double>("output", "arepo_unitlength", 3.085678e21); // 1.0 kpc
UnitMass_in_g = cf.getValueSafe<double>("output", "arepo_unitmass", 1.989e43); // 1.0e10 solar masses
UnitVelocity_in_cm_per_s = cf.getValueSafe<double>("output", "arepo_unitvel", 1e5); // 1 km/sec
coarsePartType = cf.get_value_safe<unsigned>("output", "arepo_coarsetype", COARSE_DM_DEFAULT_PARTTYPE);
UnitLength_in_cm = cf.get_value_safe<double>("output", "arepo_unitlength", 3.085678e21); // 1.0 kpc
UnitMass_in_g = cf.get_value_safe<double>("output", "arepo_unitmass", 1.989e43); // 1.0e10 solar masses
UnitVelocity_in_cm_per_s = cf.get_value_safe<double>("output", "arepo_unitvel", 1e5); // 1 km/sec
omega0 = cf.getValue<double>("cosmology", "Omega_m");
omega_b = cf.getValue<double>("cosmology", "Omega_b");
omega_L = cf.getValue<double>("cosmology", "Omega_L");
redshift = cf.getValue<double>("setup", "zstart");
boxSize = cf.getValue<double>("setup", "boxlength");
doBaryons = cf.getValueSafe<bool>("setup", "baryons", false);
useLongIDs = cf.getValueSafe<bool>("output", "arepo_longids", false);
numFiles = cf.getValueSafe<unsigned>("output", "arepo_num_files", 1);
doublePrec = cf.getValueSafe<bool>("output", "arepo_doubleprec", 0);
omega0 = cf.get_value<double>("cosmology", "Omega_m");
omega_b = cf.get_value<double>("cosmology", "Omega_b");
omega_L = cf.get_value<double>("cosmology", "Omega_L");
redshift = cf.get_value<double>("setup", "zstart");
boxSize = cf.get_value<double>("setup", "boxlength");
doBaryons = cf.get_value_safe<bool>("setup", "baryons", false);
useLongIDs = cf.get_value_safe<bool>("output", "arepo_longids", false);
numFiles = cf.get_value_safe<unsigned>("output", "arepo_num_files", 1);
doublePrec = cf.get_value_safe<bool>("output", "arepo_doubleprec", 0);
for (unsigned i = 0; i < numFiles; i++)
nPart.push_back(std::vector<unsigned int>(NTYPES, 0));
@ -583,7 +583,7 @@ public:
}
// calculate Tini for gas
hubbleParam = cf.getValue<double>("cosmology", "H0") / 100.0;
hubbleParam = cf.get_value<double>("cosmology", "H0") / 100.0;
double astart = 1.0 / (1.0 + redshift);
double h2 = hubbleParam * hubbleParam;
@ -604,7 +604,7 @@ public:
if (coarsePartType == GAS_PARTTYPE || coarsePartType == HIGHRES_DM_PARTTYPE)
throw std::runtime_error("Error: Specified illegal Arepo particle type for coarse particles.");
if (coarsePartType == STAR_PARTTYPE)
LOGWARN("WARNING: Specified coarse particle type will collide with stars if USE_SFR enabled.");
music::wlog.Print("WARNING: Specified coarse particle type will collide with stars if USE_SFR enabled.");
// create file(s)
for (unsigned i = 0; i < numFiles; i++)

80
src/plugins/output_tipsy.cc
View file

@ -107,7 +107,7 @@ protected:
if( !this->good() )
{
LOGERR("Could not open buffer file in TIPSY output plug-in");
music::elog.Print("Could not open buffer file in TIPSY output plug-in");
throw std::runtime_error("Could not open buffer file in TIPSY output plug-in");
}
@ -115,8 +115,8 @@ protected:
if( blk != (size_t)(npart*sizeof(T_store)) )
{
LOGERR("Internal consistency error in TIPSY output plug-in");
LOGERR("Expected %d bytes in temp file but found %d",npart*(unsigned)sizeof(T_store),blk);
music::elog.Print("Internal consistency error in TIPSY output plug-in");
music::elog.Print("Expected %d bytes in temp file but found %d",npart*(unsigned)sizeof(T_store),blk);
throw std::runtime_error("Internal consistency error in TIPSY output plug-in");
}
}
@ -133,7 +133,7 @@ protected:
if( !this->good() )
{
LOGERR("Could not open buffer file \'%s\' in TIPSY output plug-in",fname.c_str());
music::elog.Print("Could not open buffer file \'%s\' in TIPSY output plug-in",fname.c_str());
throw std::runtime_error("Could not open buffer file in TIPSY output plug-in");
}
@ -141,8 +141,8 @@ protected:
if( blk != (size_t)(npart*sizeof(T_store)) )
{
LOGERR("Internal consistency error in TIPSY output plug-in");
LOGERR("Expected %d bytes in temp file but found %d",npart*(unsigned)sizeof(T_store),blk);
music::elog.Print("Internal consistency error in TIPSY output plug-in");
music::elog.Print("Expected %d bytes in temp file but found %d",npart*(unsigned)sizeof(T_store),blk);
throw std::runtime_error("Internal consistency error in TIPSY output plug-in");
}
}
@ -158,7 +158,7 @@ protected:
if( !this->good() )
{
LOGERR("Could not open buffer file in TIPSY output plug-in");
music::elog.Print("Could not open buffer file in TIPSY output plug-in");
throw std::runtime_error("Could not open buffer file in TIPSY output plug-in");
}
@ -166,8 +166,8 @@ protected:
if( blk != npart*sizeof(T_store) )
{
LOGERR("Internal consistency error in TIPSY output plug-in");
LOGERR("Expected %ld bytes in temp file but found %ld",npart*sizeof(T_store),blk);
music::elog.Print("Internal consistency error in TIPSY output plug-in");
music::elog.Print("Expected %ld bytes in temp file but found %ld",npart*sizeof(T_store),blk);
throw std::runtime_error("Internal consistency error in TIPSY output plug-in");
}
@ -190,7 +190,7 @@ protected:
if( !this->good() )
{
LOGERR("Could not open buffer file \'%s\' in TIPSY output plug-in",fname.c_str());
music::elog.Print("Could not open buffer file \'%s\' in TIPSY output plug-in",fname.c_str());
throw std::runtime_error("Could not open buffer file in TIPSY output plug-in");
}
@ -198,8 +198,8 @@ protected:
if( blk != npart*sizeof(T_store) )
{
LOGERR("Internal consistency error in TIPSY output plug-in");
LOGERR("Expected %ld bytes in temp file but found %ld",npart*sizeof(T_store),blk);
music::elog.Print("Internal consistency error in TIPSY output plug-in");
music::elog.Print("Expected %ld bytes in temp file but found %ld",npart*sizeof(T_store),blk);
throw std::runtime_error("Internal consistency error in TIPSY output plug-in");
}
@ -373,7 +373,7 @@ protected:
//std::cout << " - Writing " << nptot << " particles to tipsy file...\n";
LOGINFO("TIPSY : output plugin will write:\n DM particles : %d\n SPH particles : %d",header_.ndark, header_.nsph);
music::ilog.Print("TIPSY : output plugin will write:\n DM particles : %d\n SPH particles : %d",header_.ndark, header_.nsph);
@ -413,7 +413,7 @@ protected:
if( with_baryons_ )
{
LOGINFO("TIPSY : writing baryon data");
music::ilog.Print("TIPSY : writing baryon data");
// compute gas temperature
@ -428,7 +428,7 @@ protected:
//const double ceint = 1.3806e-16/1.6726e-24 * Tini * npol / mu / unitv / unitv;
T_store temperature = (T_store) Tini;
LOGINFO("TIPSY : set initial gas temperature to %.2f K (mu = %.2f)",Tini,mu);
music::ilog.Print("TIPSY : set initial gas temperature to %.2f K (mu = %.2f)",Tini,mu);
// write
@ -514,7 +514,7 @@ protected:
//... dark matter particles ..................................................
LOGINFO("TIPSY : writing DM data");
music::ilog.Print("TIPSY : writing DM data");
ifs_x.open( fnx, npcdm );
ifs_y.open( fny, npcdm );
@ -620,7 +620,7 @@ protected:
LOGINFO("TIPSY : done writing.");
music::ilog.Print("TIPSY : done writing.");
}
@ -630,37 +630,37 @@ public:
tipsy_output_plugin( config_file& cf )
: output_plugin( cf ), ofs_( fname_.c_str(), std::ios::binary|std::ios::trunc )
{
block_buf_size_ = cf_.getValueSafe<unsigned>("output","tipsy_blksize",10485760); // default buffer size is 10 MB
block_buf_size_ = cf_.get_value_safe<unsigned>("output","tipsy_blksize",10485760); // default buffer size is 10 MB
//... ensure that everyone knows we want to do SPH
cf.insertValue("setup","do_SPH","yes");
with_baryons_ = cf_.getValue<bool>("setup","baryons");
cf.insert_value("setup","do_SPH","yes");
with_baryons_ = cf_.get_value<bool>("setup","baryons");
//bbndparticles_ = !cf_.getValueSafe<bool>("output","gadget_nobndpart",false);
//bbndparticles_ = !cf_.get_value_safe<bool>("output","gadget_nobndpart",false);
npartmax_ = 1<<30;
if(!ofs_.good())
{
LOGERR("tipsy output plug-in could not open output file \'%s\' for writing!",fname_.c_str());
music::elog.Print("tipsy output plug-in could not open output file \'%s\' for writing!",fname_.c_str());
throw std::runtime_error(std::string("tipsy output plug-in could not open output file \'")+fname_+"\' for writing!\n");
}
ofs_.close();
double zstart = cf.getValue<double>("setup","zstart");
double zstart = cf.get_value<double>("setup","zstart");
astart_ = 1.0/(1.0+zstart);
omegam_ = cf.getValue<double>("cosmology","Omega_m");
omegab_ = cf.getValue<double>("cosmology","Omega_b");
boxsize_ = cf.getValue<double>("setup","boxlength");
omegam_ = cf.get_value<double>("cosmology","Omega_m");
omegab_ = cf.get_value<double>("cosmology","Omega_b");
boxsize_ = cf.get_value<double>("setup","boxlength");
epsfac_ = cf.getValueSafe<double>("output","tipsy_eps",0.05);
epsfac_coarse_ = cf.getValueSafe<double>("output","tipsy_eps_coarse",epsfac_);
epsfac_gas_ = cf.getValueSafe<double>("output","tipsy_eps_gas",epsfac_);
epsfac_ = cf.get_value_safe<double>("output","tipsy_eps",0.05);
epsfac_coarse_ = cf.get_value_safe<double>("output","tipsy_eps_coarse",epsfac_);
epsfac_gas_ = cf.get_value_safe<double>("output","tipsy_eps_gas",epsfac_);
H0_ = cf.getValue<double>("cosmology","H0");
YHe_ = cf.getValueSafe<double>("cosmology","YHe",0.248);
gamma_ = cf.getValueSafe<double>("cosmology","gamma",5.0/3.0);
H0_ = cf.get_value<double>("cosmology","H0");
YHe_ = cf.get_value_safe<double>("cosmology","YHe",0.248);
gamma_ = cf.get_value_safe<double>("cosmology","gamma",5.0/3.0);
native_ = cf.getValueSafe<bool>("output","tipsy_native",false);
native_ = cf.get_value_safe<bool>("output","tipsy_native",false);
}
@ -746,7 +746,7 @@ public:
if( nwritten != nptot )
{
LOGERR("TIPSY output plugin wrote %ld, should have %ld", nwritten, nptot);
music::elog.Print("TIPSY output plugin wrote %ld, should have %ld", nwritten, nptot);
throw std::runtime_error("Internal consistency error while writing temporary file for DM masses");
}
ofs_temp.write( (char *)&blksize, sizeof(size_t) );
@ -804,14 +804,14 @@ public:
}
if( nwritten != nptot ){
LOGERR("TIPSY output plugin wrote %ld gas particles, should have %ld", nwritten, nptot);
music::elog.Print("TIPSY output plugin wrote %ld gas particles, should have %ld", nwritten, nptot);
throw std::runtime_error("Internal consistency error while writing temporary file for baryon masses");
}
ofs_temp.write( (char *)&blksize, sizeof(size_t) );
if( ofs_temp.bad() ){
LOGERR("I/O error while writing temporary file for baryon masse");
music::elog.Print("I/O error while writing temporary file for baryon masse");
throw std::runtime_error("I/O error while writing temporary file for baryon masses");
}
}
@ -866,7 +866,7 @@ public:
if( nwritten != nptot )
{
LOGERR("TIPSY output plugin wrote %ld, should have %ld", nwritten, nptot);
music::elog.Print("TIPSY output plugin wrote %ld, should have %ld", nwritten, nptot);
throw std::runtime_error("Internal consistency error while writing temporary file for positions");
}
@ -922,7 +922,7 @@ public:
if( nwritten != nptot )
{
LOGERR("TIPSY output plugin wrote %ld, should have %ld", nwritten, nptot);
music::elog.Print("TIPSY output plugin wrote %ld, should have %ld", nwritten, nptot);
throw std::runtime_error("Internal consistency error while writing temporary file for DM velocities");
}
@ -996,7 +996,7 @@ public:
if( nwritten != npart )
{
LOGERR("TIPSY output plugin wrote %ld, should have %ld", nwritten, npart);
music::elog.Print("TIPSY output plugin wrote %ld, should have %ld", nwritten, npart);
throw std::runtime_error("Internal consistency error while writing temporary file for baryon velocities");
}
@ -1068,7 +1068,7 @@ public:
if( nwritten != npart )
{
LOGERR("TIPSY output plugin wrote %ld, should have %ld", nwritten, npart);
music::elog.Print("TIPSY output plugin wrote %ld, should have %ld", nwritten, npart);
throw std::runtime_error("Internal consistency error while writing temporary file for baryon positions");
}

100
src/plugins/output_tipsy_resample.cc
View file

@ -114,7 +114,7 @@ protected:
if (!this->good ())
{
LOGERR ("Could not open buffer file in TIPSY output plug-in");
music::elog.Print ("Could not open buffer file in TIPSY output plug-in");
throw std::runtime_error("Could not open buffer file in TIPSY output plug-in");
}
@ -122,8 +122,8 @@ protected:
if (blk != (size_t) (npart * sizeof (T_store)))
{
LOGERR ("Internal consistency error in TIPSY output plug-in");
LOGERR ("Expected %d bytes in temp file but found %d", npart * (unsigned) sizeof (T_store), blk);
music::elog.Print ("Internal consistency error in TIPSY output plug-in");
music::elog.Print ("Expected %d bytes in temp file but found %d", npart * (unsigned) sizeof (T_store), blk);
throw std::runtime_error("Internal consistency error in TIPSY output plug-in");
}
}
@ -140,7 +140,7 @@ protected:
if (!this->good ())
{
LOGERR ("Could not open buffer file \'%s\' in TIPSY output plug-in", fname.c_str ());
music::elog.Print ("Could not open buffer file \'%s\' in TIPSY output plug-in", fname.c_str ());
throw std::runtime_error("Could not open buffer file in TIPSY output plug-in");
}
@ -148,8 +148,8 @@ protected:
if (blk != (size_t) (npart * sizeof (T_store)))
{
LOGERR ("Internal consistency error in TIPSY output plug-in");
LOGERR ("Expected %d bytes in temp file but found %d", npart * (unsigned) sizeof (T_store), blk);
music::elog.Print ("Internal consistency error in TIPSY output plug-in");
music::elog.Print ("Expected %d bytes in temp file but found %d", npart * (unsigned) sizeof (T_store), blk);
throw std::runtime_error("Internal consistency error in TIPSY output plug-in");
}
}
@ -167,7 +167,7 @@ protected:
if (!this->good ())
{
LOGERR ("Could not open buffer file in TIPSY output plug-in");
music::elog.Print ("Could not open buffer file in TIPSY output plug-in");
throw std::runtime_error("Could not open buffer file in TIPSY output plug-in");
}
@ -175,8 +175,8 @@ protected:
if (blk != npart * sizeof (T_store))
{
LOGERR ("Internal consistency error in TIPSY output plug-in");
LOGERR ("Expected %ld bytes in temp file but found %ld", npart * sizeof (T_store), blk);
music::elog.Print ("Internal consistency error in TIPSY output plug-in");
music::elog.Print ("Expected %ld bytes in temp file but found %ld", npart * sizeof (T_store), blk);
throw std::runtime_error("Internal consistency error in TIPSY output plug-in");
}
@ -199,7 +199,7 @@ protected:
if (!this->good ())
{
LOGERR ("Could not open buffer file \'%s\' in TIPSY output plug-in", fname.c_str ());
music::elog.Print ("Could not open buffer file \'%s\' in TIPSY output plug-in", fname.c_str ());
throw std::runtime_error("Could not open buffer file in TIPSY output plug-in");
}
@ -207,8 +207,8 @@ protected:
if (blk != npart * sizeof (T_store))
{
LOGERR ("Internal consistency error in TIPSY output plug-in");
LOGERR ("Expected %ld bytes in temp file but found %ld", npart * sizeof (T_store), blk);
music::elog.Print ("Internal consistency error in TIPSY output plug-in");
music::elog.Print ("Expected %ld bytes in temp file but found %ld", npart * sizeof (T_store), blk);
throw std::runtime_error("Internal consistency error in TIPSY output plug-in");
}
@ -380,7 +380,7 @@ protected:
//std::cout << " - Writing " << nptot << " particles to tipsy file...\n";
LOGINFO
music::ilog.Print
("TIPSY : output plugin will write:\n DM particles : %d\n SPH particles : %d",
header_.ndark, header_.nsph);
@ -415,7 +415,7 @@ protected:
if (with_baryons_)
{
LOGINFO ("TIPSY : writing baryon data");
music::ilog.Print ("TIPSY : writing baryon data");
// compute gas temperature
@ -434,7 +434,7 @@ protected:
//const double ceint = 1.3806e-16/1.6726e-24 * Tini * npol / mu / unitv / unitv;
T_store temperature = (T_store) Tini;
LOGINFO("TIPSY : set initial gas temperature to %.2f K (mu = %.2f)", Tini, mu);
music::ilog.Print("TIPSY : set initial gas temperature to %.2f K (mu = %.2f)", Tini, mu);
// write
@ -503,7 +503,7 @@ protected:
//... dark matter particles ..................................................
LOGINFO ("TIPSY : writing DM data");
music::ilog.Print ("TIPSY : writing DM data");
ifs_x.open (fnx, npcdm);
ifs_y.open (fny, npcdm);
@ -592,7 +592,7 @@ protected:
LOGINFO ("TIPSY : done writing.");
music::ilog.Print ("TIPSY : done writing.");
}
@ -602,66 +602,66 @@ public:
tipsy_output_plugin_res (config_file & cf)
: output_plugin (cf), ofs_ (fname_.c_str (), std::ios::binary | std::ios::trunc)
{
block_buf_size_ = cf_.getValueSafe < unsigned >("output", "tipsy_blksize", 10485760); // default buffer size is 10 MB
block_buf_size_ = cf_.get_value_safe < unsigned >("output", "tipsy_blksize", 10485760); // default buffer size is 10 MB
//... ensure that everyone knows we want to do SPH
cf.insertValue ("setup", "do_SPH", "yes");
with_baryons_ = cf_.getValue < bool > ("setup", "baryons");
cf.insert_value ("setup", "do_SPH", "yes");
with_baryons_ = cf_.get_value < bool > ("setup", "baryons");
//bbndparticles_ = !cf_.getValueSafe<bool>("output","gadget_nobndpart",false);
//bbndparticles_ = !cf_.get_value_safe<bool>("output","gadget_nobndpart",false);
npartmax_ = 1 << 30;
if (!ofs_.good ())
{
LOGERR("tipsy output plug-in could not open output file \'%s\' for writing!",fname_.c_str ());
music::elog.Print("tipsy output plug-in could not open output file \'%s\' for writing!",fname_.c_str ());
throw std::runtime_error (std::string("tipsy output plug-in could not open output file \'")+ fname_ + "\' for writing!\n");
}
ofs_.close ();
double zstart = cf.getValue < double >("setup", "zstart");
double zstart = cf.get_value < double >("setup", "zstart");
astart_ = 1.0 / (1.0 + zstart);
omegam_ = cf.getValue < double >("cosmology", "Omega_m");
omegab_ = cf.getValue < double >("cosmology", "Omega_b");
boxsize_ = cf.getValue < double >("setup", "boxlength");
epsfac_ = cf.getValueSafe < double >("output", "tipsy_eps", 0.05);
H0_ = cf.getValue < double >("cosmology", "H0");
YHe_ = cf.getValueSafe < double >("cosmology", "YHe", 0.248);
gamma_ = cf.getValueSafe < double >("cosmology", "gamma", 5.0 / 3.0);
omegam_ = cf.get_value < double >("cosmology", "Omega_m");
omegab_ = cf.get_value < double >("cosmology", "Omega_b");
boxsize_ = cf.get_value < double >("setup", "boxlength");
epsfac_ = cf.get_value_safe < double >("output", "tipsy_eps", 0.05);
H0_ = cf.get_value < double >("cosmology", "H0");
YHe_ = cf.get_value_safe < double >("cosmology", "YHe", 0.248);
gamma_ = cf.get_value_safe < double >("cosmology", "gamma", 5.0 / 3.0);
bresample_ = cf_.containsKey("output","tipsy_resfine");
bresample_ = cf_.contains_key("output","tipsy_resfine");
if( bresample_ )
{
LOGINFO("Resampling in high-res region enabled for TIPSY output,\n" \
music::ilog.Print("Resampling in high-res region enabled for TIPSY output,\n" \
" Setting option \'[output]/glass_cicdeconvolve=yes\'.");
np_resample_ = cf_.getValue<size_t> ("output","tipsy_resfine");
np_resample_ = cf_.get_value<size_t> ("output","tipsy_resfine");
unsigned nfine[3];
char tempstr[256];
unsigned levelmax = cf_.getValue<unsigned>("setup","levelmax");
unsigned levelmax = cf_.get_value<unsigned>("setup","levelmax");
sprintf(tempstr,"size(%d,0)",levelmax);
nfine[0] = cf_.getValue<unsigned>("setup",tempstr);
nfine[0] = cf_.get_value<unsigned>("setup",tempstr);
sprintf(tempstr,"size(%d,1)",levelmax);
nfine[1] = cf_.getValue<unsigned>("setup",tempstr);
nfine[1] = cf_.get_value<unsigned>("setup",tempstr);
sprintf(tempstr,"size(%d,2)",levelmax);
nfine[2] = cf_.getValue<unsigned>("setup",tempstr);
nfine[2] = cf_.get_value<unsigned>("setup",tempstr);
if( nfine[0]!=nfine[1] || nfine[0]!=nfine[2] )
{
LOGERR("Need to set \'[setup]/force_equal_extent=yes\' when using \'tipsy_refine=yes\'!");
music::elog.Print("Need to set \'[setup]/force_equal_extent=yes\' when using \'tipsy_refine=yes\'!");
throw std::runtime_error("Need to set \'[setup]/force_equal_extent=yes\' when using \'tipsy_refine=yes\'!");
}
double resfac = (double)nfine[0]/(double)np_resample_;
sprintf(tempstr,"%g",resfac*0.5);
cf_.insertValue("setup","baryon_staggering",std::string(tempstr));
cf_.insert_value("setup","baryon_staggering",std::string(tempstr));
cf_.insertValue("output","glass_cicdeconvolve","yes");
cf_.insert_value("output","glass_cicdeconvolve","yes");
}
}
@ -781,7 +781,7 @@ public:
if (nwritten != nptot)
{
LOGERR ("TIPSY output plugin wrote %ld, should have %ld", nwritten, nptot);
music::elog.Print ("TIPSY output plugin wrote %ld, should have %ld", nwritten, nptot);
throw std::runtime_error("Internal consistency error while writing temporary file for DM masses");
}
ofs_temp.write ((char *) &blksize, sizeof (size_t));
@ -865,19 +865,19 @@ public:
double cmass, bmass(0.0);
std::cout << "-------------------------------------------------------------" << std::endl;
LOGINFO("TIPSY: particle resampling is enabled");
LOGINFO("TIPSY: new high-res particles have the masses:");
music::ilog.Print("TIPSY: particle resampling is enabled");
music::ilog.Print("TIPSY: new high-res particles have the masses:");
if( with_baryons_ )
{
cmass = (omegam_-omegab_)*rhom*dx3;
bmass = omegab_*rhom*dx3;
LOGINFO("TIPSY: DM particle mass = %g h-1 M_o",cmass);
LOGINFO("TIPSY: baryon particle mass = %g h-1 M_o",bmass);
music::ilog.Print("TIPSY: DM particle mass = %g h-1 M_o",cmass);
music::ilog.Print("TIPSY: baryon particle mass = %g h-1 M_o",bmass);
}
else
{
cmass = omegam_*rhom*dx3;
LOGINFO("TIPSY: particle mass = %g h-1 M_o",cmass);
music::ilog.Print("TIPSY: particle mass = %g h-1 M_o",cmass);
}
}
}
@ -1030,7 +1030,7 @@ public:
if (nwritten != nptot)
{
LOGERR ("TIPSY output plugin wrote %ld, should have %ld", nwritten,
music::elog.Print ("TIPSY output plugin wrote %ld, should have %ld", nwritten,
nptot);
throw std::
runtime_error
@ -1131,7 +1131,7 @@ public:
if (nwritten != nptot)
{
LOGERR ("TIPSY output plugin wrote %ld, should have %ld", nwritten, nptot);
music::elog.Print ("TIPSY output plugin wrote %ld, should have %ld", nwritten, nptot);
throw std::runtime_error("Internal consistency error while writing temporary file for DM velocities");
}
@ -1240,7 +1240,7 @@ public:
if (nwritten != npart)
{
LOGERR ("TIPSY output plugin wrote %ld, should have %ld", nwritten, npart);
music::elog.Print ("TIPSY output plugin wrote %ld, should have %ld", nwritten, npart);
throw std::runtime_error("Internal consistency error while writing temporary file for baryon velocities");
}
@ -1356,7 +1356,7 @@ public:
if (nwritten != npart)
{
LOGERR ("TIPSY output plugin wrote %ld, should have %ld", nwritten, npart);
music::elog.Print ("TIPSY output plugin wrote %ld, should have %ld", nwritten, npart);
throw std::runtime_error("Internal consistency error while writing temporary file for baryon positions");
}

12
src/plugins/point_file_reader.hh
View file

@ -5,7 +5,7 @@
#include <sstream>
#include <string>
#include <vector>
#include "log.hh"
#include <logger.hh>
struct point_reader{
@ -31,7 +31,7 @@ struct point_reader{
std::ifstream ifs(fname.c_str());
if( !ifs )
{
LOGERR("region_ellipsoid_plugin::read_points_from_file : Could not open file \'%s\'",fname.c_str());
music::elog.Print("region_ellipsoid_plugin::read_points_from_file : Could not open file \'%s\'",fname.c_str());
throw std::runtime_error("region_ellipsoid_plugin::read_points_from_file : cannot open point file.");
}
@ -58,16 +58,16 @@ struct point_reader{
++row;
if( row>1 && colcount != colcount1 )
LOGERR("error on line %d of input file",row);
music::elog.Print("error on line %d of input file",row);
//std::cout << std::endl;
}
LOGINFO("region point file appears to contain %d columns",colcount);
music::ilog.Print("region point file appears to contain %d columns",colcount);
if( p.size()%3 != 0 && p.size()%6 != 0 )
{
LOGERR("Region point file \'%s\' does not contain triplets (%d elems)",fname.c_str(),p.size());
music::elog.Print("Region point file \'%s\' does not contain triplets (%d elems)",fname.c_str(),p.size());
throw std::runtime_error("region_ellipsoid_plugin::read_points_from_file : file does not contain triplets.");
}
@ -123,7 +123,7 @@ struct point_reader{
}
}
else
LOGERR("Problem interpreting the region point file \'%s\'", fname.c_str() );
music::elog.Print("Problem interpreting the region point file \'%s\'", fname.c_str() );
num_columns = colcount;
}

68
src/plugins/random_music.cc
View file

@ -46,17 +46,17 @@ public:
levelmin_ = prefh_->levelmin();
levelmax_ = prefh_->levelmax();
ran_cube_size_ = pcf_->getValueSafe<unsigned>("random", "cubesize", DEF_RAN_CUBE_SIZE);
disk_cached_ = pcf_->getValueSafe<bool>("random", "disk_cached", true);
restart_ = pcf_->getValueSafe<bool>("random", "restart", false);
ran_cube_size_ = pcf_->get_value_safe<unsigned>("random", "cubesize", DEF_RAN_CUBE_SIZE);
disk_cached_ = pcf_->get_value_safe<bool>("random", "disk_cached", true);
restart_ = pcf_->get_value_safe<bool>("random", "restart", false);
pcf_->insertValue("setup","fourier_splicing","true");
pcf_->insert_value("setup","fourier_splicing","true");
mem_cache_.assign(levelmax_ - levelmin_ + 1, (std::vector<real_t> *)NULL);
if (restart_ && !disk_cached_)
{
LOGERR("Cannot restart from mem cached random numbers.");
music::elog.Print("Cannot restart from mem cached random numbers.");
throw std::runtime_error("Cannot restart from mem cached random numbers.");
}
@ -93,8 +93,8 @@ void RNG_music::parse_random_parameters(void)
std::string tempstr;
bool noseed = false;
sprintf(seedstr, "seed[%d]", i);
if (pcf_->containsKey("random", seedstr))
tempstr = pcf_->getValue<std::string>("random", seedstr);
if (pcf_->contains_key("random", seedstr))
tempstr = pcf_->get_value<std::string>("random", seedstr);
else
{
// "-2" means that no seed entry was found for that level
@ -116,7 +116,7 @@ void RNG_music::parse_random_parameters(void)
{
if (ltemp <= 0)
{
LOGERR("Specified seed [random]/%s needs to be a number >0!", seedstr);
music::elog.Print("Specified seed [random]/%s needs to be a number >0!", seedstr);
throw std::runtime_error("Seed values need to be >0");
}
rngseeds_.push_back(ltemp);
@ -126,7 +126,7 @@ void RNG_music::parse_random_parameters(void)
{
rngfnames_.push_back(tempstr);
rngseeds_.push_back(-1);
LOGINFO("Random numbers for level %3d will be read from file.", i);
music::ilog.Print("Random numbers for level %3d will be read from file.", i);
}
}
@ -141,7 +141,7 @@ void RNG_music::parse_random_parameters(void)
void RNG_music::compute_random_numbers(void)
{
bool rndsign = pcf_->getValueSafe<bool>("random", "grafic_sign", false);
bool rndsign = pcf_->get_value_safe<bool>("random", "grafic_sign", false);
std::vector<rng *> randc(std::max(levelmax_, levelmin_seed_) + 1, (rng *)NULL);
@ -160,7 +160,7 @@ void RNG_music::compute_random_numbers(void)
// #warning add possibility to read noise from file also here!
if (rngfnames_[i].size() > 0)
LOGINFO("Warning: Cannot use filenames for higher levels currently! Ignoring!");
music::ilog.Print("Warning: Cannot use filenames for higher levels currently! Ignoring!");
randc[i] = new rng(*randc[i - 1], ran_cube_size_, rngseeds_[i]);
delete randc[i - 1];
@ -180,7 +180,7 @@ void RNG_music::compute_random_numbers(void)
for (int ilevel = levelmin_seed_ - 1; ilevel >= (int)levelmin_; --ilevel)
{
if (rngseeds_[ilevel - levelmin_] > 0)
LOGINFO("Warning: random seed for level %d will be ignored.\n"
music::ilog.Print("Warning: random seed for level %d will be ignored.\n"
" consistency requires that it is obtained by restriction from level %d",
ilevel, levelmin_seed_);
@ -212,11 +212,11 @@ void RNG_music::compute_random_numbers(void)
{
int lx[3], x0[3];
int shift[3], levelmin_poisson;
shift[0] = pcf_->getValue<int>("setup", "shift_x");
shift[1] = pcf_->getValue<int>("setup", "shift_y");
shift[2] = pcf_->getValue<int>("setup", "shift_z");
shift[0] = pcf_->get_value<int>("setup", "shift_x");
shift[1] = pcf_->get_value<int>("setup", "shift_y");
shift[2] = pcf_->get_value<int>("setup", "shift_z");
levelmin_poisson = pcf_->getValue<unsigned>("setup", "levelmin");
levelmin_poisson = pcf_->get_value<unsigned>("setup", "levelmin");
int lfac = 1 << (ilevel - levelmin_poisson);
@ -255,11 +255,11 @@ void RNG_music::compute_random_numbers(void)
void RNG_music::store_rnd(int ilevel, rng *prng)
{
int shift[3], levelmin_poisson;
shift[0] = pcf_->getValue<int>("setup", "shift_x");
shift[1] = pcf_->getValue<int>("setup", "shift_y");
shift[2] = pcf_->getValue<int>("setup", "shift_z");
shift[0] = pcf_->get_value<int>("setup", "shift_x");
shift[1] = pcf_->get_value<int>("setup", "shift_y");
shift[2] = pcf_->get_value<int>("setup", "shift_z");
levelmin_poisson = pcf_->getValue<unsigned>("setup", "levelmin");
levelmin_poisson = pcf_->get_value<unsigned>("setup", "levelmin");
int lfac = 1 << (ilevel - levelmin_poisson);
@ -279,7 +279,7 @@ void RNG_music::store_rnd(int ilevel, rng *prng)
char fname[128];
sprintf(fname, "grafic_wnoise_%04d.bin", ilevel);
LOGUSER("Storing white noise field for grafic in file \'%s\'...", fname);
music::ulog.Print("Storing white noise field for grafic in file \'%s\'...", fname);
std::ofstream ofs(fname, std::ios::binary | std::ios::trunc);
data.assign(N * N, 0.0);
@ -325,8 +325,8 @@ void RNG_music::store_rnd(int ilevel, rng *prng)
char fname[128];
sprintf(fname, "grafic_wnoise_%04d.bin", ilevel);
LOGUSER("Storing white noise field for grafic in file \'%s\'...", fname);
LOGDEBUG("(%d,%d,%d) -- (%d,%d,%d) -- lfac = %d", nx, ny, nz, i0, j0, k0, lfac);
music::ulog.Print("Storing white noise field for grafic in file \'%s\'...", fname);
music::dlog.Print("(%d,%d,%d) -- (%d,%d,%d) -- lfac = %d", nx, ny, nz, i0, j0, k0, lfac);
std::ofstream ofs(fname, std::ios::binary | std::ios::trunc);
data.assign(nx * ny, 0.0);
@ -372,7 +372,7 @@ void RNG_music::store_rnd(int ilevel, rng *prng)
char fname[128];
sprintf(fname, "wnoise_%04d.bin", ilevel);
LOGUSER("Storing white noise field in file \'%s\'...", fname);
music::ulog.Print("Storing white noise field in file \'%s\'...", fname);
std::ofstream ofs(fname, std::ios::binary | std::ios::trunc);
@ -418,7 +418,7 @@ void RNG_music::store_rnd(int ilevel, rng *prng)
char fname[128];
sprintf(fname, "wnoise_%04d.bin", ilevel);
LOGUSER("Storing white noise field in file \'%s\'...", fname);
music::ulog.Print("Storing white noise field in file \'%s\'...", fname);
std::ofstream ofs(fname, std::ios::binary | std::ios::trunc);
@ -474,7 +474,7 @@ void RNG_music::store_rnd(int ilevel, rng *prng)
mem_cache_[ilevel - levelmin_] = new std::vector<real_t>(nx * ny * nz, 0.0);
LOGUSER("Copying white noise to mem cache...");
music::ulog.Print("Copying white noise to mem cache...");
#pragma omp parallel for
for (int i = 0; i < nx; ++i)
@ -489,12 +489,12 @@ void RNG_music::fill_grid(int ilevel, DensityGrid<real_t> &A)
{
if (!initialized_)
{
LOGERR("Call to RNG_music::fill_grid before call to RNG_music::initialize_for_grid_structure");
music::elog.Print("Call to RNG_music::fill_grid before call to RNG_music::initialize_for_grid_structure");
throw std::runtime_error("invalid call order for random number generator");
}
if (restart_)
LOGINFO("Attempting to restart using random numbers for level %d\n from file \'wnoise_%04d.bin\'.", ilevel,
music::ilog.Print("Attempting to restart using random numbers for level %d\n from file \'wnoise_%04d.bin\'.", ilevel,
ilevel);
if (disk_cached_)
@ -502,12 +502,12 @@ void RNG_music::fill_grid(int ilevel, DensityGrid<real_t> &A)
char fname[128];
sprintf(fname, "wnoise_%04d.bin", ilevel);
LOGUSER("Loading white noise from file \'%s\'...", fname);
music::ulog.Print("Loading white noise from file \'%s\'...", fname);
std::ifstream ifs(fname, std::ios::binary);
if (!ifs.good())
{
LOGERR("White noise file \'%s\'was not found.", fname);
music::elog.Print("White noise file \'%s\'was not found.", fname);
throw std::runtime_error("A white noise file was not found. This is an internal inconsistency and bad.");
}
@ -553,7 +553,7 @@ void RNG_music::fill_grid(int ilevel, DensityGrid<real_t> &A)
}
else
{
LOGERR("White noise file is not aligned with array. File: [%d,%d,%d]. Mem: [%d,%d,%d].", nx, ny, nz, A.size(0),
music::elog.Print("White noise file is not aligned with array. File: [%d,%d,%d]. Mem: [%d,%d,%d].", nx, ny, nz, A.size(0),
A.size(1), A.size(2));
throw std::runtime_error(
"White noise file is not aligned with array. This is an internal inconsistency and bad.");
@ -577,16 +577,16 @@ void RNG_music::fill_grid(int ilevel, DensityGrid<real_t> &A)
}
else
{
LOGUSER("Copying white noise from memory cache...");
music::ulog.Print("Copying white noise from memory cache...");
if (mem_cache_[ilevel - levelmin_] == NULL)
LOGERR("Tried to access mem-cached random numbers for level %d. But these are not available!\n", ilevel);
music::elog.Print("Tried to access mem-cached random numbers for level %d. But these are not available!\n", ilevel);
int nx(A.size(0)), ny(A.size(1)), nz(A.size(2));
if ((size_t)nx * (size_t)ny * (size_t)nz != mem_cache_[ilevel - levelmin_]->size())
{
LOGERR("White noise file is not aligned with array. File: [%d,%d,%d]. Mem: [%d,%d,%d].", nx, ny, nz, A.size(0),
music::elog.Print("White noise file is not aligned with array. File: [%d,%d,%d]. Mem: [%d,%d,%d].", nx, ny, nz, A.size(0),
A.size(1), A.size(2));
throw std::runtime_error("White noise file is not aligned with array. This is an internal inconsistency and bad");
}

40
src/plugins/random_music_wnoise_generator.cc
View file

@ -11,7 +11,7 @@
template <typename T>
void rapid_proto_ngenic_rng(size_t res, long baseseed, music_wnoise_generator<T> &R)
{
LOGUSER("Invoking the N-GenIC random number generator");
music::ulog.Print("Invoking the N-GenIC random number generator");
unsigned *seedtable = new unsigned[res * res];
@ -162,7 +162,7 @@ template <typename T>
music_wnoise_generator<T>::music_wnoise_generator(unsigned res, unsigned cubesize, long baseseed, int *x0, int *lx)
: res_(res), cubesize_(cubesize), ncubes_(1), baseseed_(baseseed)
{
LOGINFO("Generating random numbers (1) with seed %ld", baseseed);
music::ilog.Print("Generating random numbers (1) with seed %ld", baseseed);
initialize();
fill_subvolume(x0, lx);
@ -172,7 +172,7 @@ template <typename T>
music_wnoise_generator<T>::music_wnoise_generator(unsigned res, unsigned cubesize, long baseseed, bool zeromean)
: res_(res), cubesize_(cubesize), ncubes_(1), baseseed_(baseseed)
{
LOGINFO("Generating random numbers (2) with seed %ld", baseseed);
music::ilog.Print("Generating random numbers (2) with seed %ld", baseseed);
double mean = 0.0;
size_t res_l = res;
@ -182,7 +182,7 @@ music_wnoise_generator<T>::music_wnoise_generator(unsigned res, unsigned cubesiz
cubesize_ = res_;
if (!musicnoise)
LOGERR("This currently breaks compatibility. Need to disable by hand! Make sure to not check into repo");
music::elog.Print("This currently breaks compatibility. Need to disable by hand! Make sure to not check into repo");
initialize();
@ -227,7 +227,7 @@ music_wnoise_generator<T>::music_wnoise_generator(unsigned res, std::string rand
std::ifstream ifs(randfname.c_str(), std::ios::binary);
if (!ifs)
{
LOGERR("Could not open random number file \'%s\'!", randfname.c_str());
music::elog.Print("Could not open random number file \'%s\'!", randfname.c_str());
throw std::runtime_error(std::string("Could not open random number file \'") + randfname + std::string("\'!"));
}
@ -323,7 +323,7 @@ music_wnoise_generator<T>::music_wnoise_generator(unsigned res, std::string rand
std::vector<float> in_float;
std::vector<double> in_double;
LOGINFO("Random number file \'%s\'\n contains %ld numbers. Reading...", randfname.c_str(), nx * ny * nz);
music::ilog.Print("Random number file \'%s\'\n contains %ld numbers. Reading...", randfname.c_str(), nx * ny * nz);
long double sum = 0.0, sum2 = 0.0;
size_t count = 0;
@ -422,7 +422,7 @@ music_wnoise_generator<T>::music_wnoise_generator(unsigned res, std::string rand
mean = sum / count;
var = sum2 / count - mean * mean;
LOGINFO("Random numbers in file have \n mean = %f and var = %f", mean, var);
music::ilog.Print("Random numbers in file have \n mean = %f and var = %f", mean, var);
}
//... copy construct by averaging down
@ -436,7 +436,7 @@ music_wnoise_generator<T>::music_wnoise_generator(/*const*/ music_wnoise_generat
size_t count = 0;
LOGINFO("Generating a coarse white noise field by k-space degrading");
music::ilog.Print("Generating a coarse white noise field by k-space degrading");
//... initialize properties of container
res_ = rc.res_ / 2;
cubesize_ = res_;
@ -445,7 +445,7 @@ music_wnoise_generator<T>::music_wnoise_generator(/*const*/ music_wnoise_generat
if (sizeof(fftw_real) != sizeof(T))
{
LOGERR("type mismatch with fftw_real in k-space averaging");
music::elog.Print("type mismatch with fftw_real in k-space averaging");
throw std::runtime_error("type mismatch with fftw_real in k-space averaging");
}
@ -580,7 +580,7 @@ music_wnoise_generator<T>::music_wnoise_generator(/*const*/ music_wnoise_generat
rmean = sum / count;
rvar = sum2 / count - rmean * rmean;
LOGINFO("Restricted random numbers have\n mean = %f, var = %f", rmean, rvar);
music::ilog.Print("Restricted random numbers have\n mean = %f, var = %f", rmean, rvar);
}
template <typename T>
@ -610,7 +610,9 @@ music_wnoise_generator<T>::music_wnoise_generator(music_wnoise_generator<T> &rc,
}
LOGINFO("Generating a constrained random number set with seed %ld\n using coarse mode replacement...", baseseed);
music::ilog << "Generating a constrained random number set with seed " << baseseed << std::endl;
music::ilog << " using coarse mode replacement..." << std::endl;
assert(lx[0] % 2 == 0 && lx[1] % 2 == 0 && lx[2] % 2 == 0);
size_t nx = lx[0], ny = lx[1], nz = lx[2],
nxc = lx[0] / 2, nyc = lx[1] / 2, nzc = lx[2] / 2;
@ -803,7 +805,7 @@ void music_wnoise_generator<T>::register_cube(int i, int j, int k)
rnums_.push_back(NULL);
cubemap_[icube] = rnums_.size() - 1;
#ifdef DEBUG
LOGDEBUG("registering new cube %d,%d,%d . ID = %ld, memloc = %ld", i, j, k, icube, cubemap_[icube]);
music::dlog.Print("registering new cube %d,%d,%d . ID = %ld, memloc = %ld", i, j, k, icube, cubemap_[icube]);
#endif
}
}
@ -827,7 +829,7 @@ double music_wnoise_generator<T>::fill_cube(int i, int j, int k)
if (it == cubemap_.end())
{
LOGERR("Attempt to access non-registered random number cube!");
music::elog.Print("Attempt to access non-registered random number cube!");
throw std::runtime_error("Attempt to access non-registered random number cube!");
}
@ -864,7 +866,7 @@ void music_wnoise_generator<T>::subtract_from_cube(int i, int j, int k, double v
if (it == cubemap_.end())
{
LOGERR("Attempt to access unallocated RND cube %d,%d,%d in music_wnoise_generator::subtract_from_cube", i, j, k);
music::elog.Print("Attempt to access unallocated RND cube %d,%d,%d in music_wnoise_generator::subtract_from_cube", i, j, k);
throw std::runtime_error("Attempt to access unallocated RND cube in music_wnoise_generator::subtract_from_cube");
}
@ -890,7 +892,7 @@ void music_wnoise_generator<T>::free_cube(int i, int j, int k)
if (it == cubemap_.end())
{
LOGERR("Attempt to access unallocated RND cube %d,%d,%d in music_wnoise_generator::free_cube", i, j, k);
music::elog.Print("Attempt to access unallocated RND cube %d,%d,%d in music_wnoise_generator::free_cube", i, j, k);
throw std::runtime_error("Attempt to access unallocated RND cube in music_wnoise_generator::free_cube");
}
@ -914,7 +916,7 @@ void music_wnoise_generator<T>::initialize(void)
cubesize_ = res_;
}
LOGINFO("Generating random numbers w/ sample cube size of %d", cubesize_);
music::ilog.Print("Generating random numbers w/ sample cube size of %d", cubesize_);
}
template <typename T>
@ -931,8 +933,8 @@ double music_wnoise_generator<T>::fill_subvolume(int *i0, int *n)
ncube[2] = (int)(n[2] / cubesize_) + 2;
#ifdef DEBUG
LOGDEBUG("random numbers needed for region %d,%d,%d ..+ %d,%d,%d", i0[0], i0[1], i0[2], n[0], n[1], n[2]);
LOGDEBUG("filling cubes %d,%d,%d ..+ %d,%d,%d", i0cube[0], i0cube[1], i0cube[2], ncube[0], ncube[1], ncube[2]);
music::dlog.Print("random numbers needed for region %d,%d,%d ..+ %d,%d,%d", i0[0], i0[1], i0[2], n[0], n[1], n[2]);
music::dlog.Print("filling cubes %d,%d,%d ..+ %d,%d,%d", i0cube[0], i0cube[1], i0cube[2], ncube[0], ncube[1], ncube[2]);
#endif
double mean = 0.0;
@ -1023,7 +1025,7 @@ void music_wnoise_generator<T>::print_allocated(void)
if (rnums_[i] != NULL)
ncount++;
LOGINFO(" -> %d of %d random number cubes currently allocated", ncount, ntot);
music::ilog.Print(" -> %d of %d random number cubes currently allocated", ncount, ntot);
}
template class music_wnoise_generator<float>;

6
src/plugins/random_music_wnoise_generator.hh
View file

@ -58,7 +58,7 @@ protected:
if (it == cubemap_.end())
{
LOGERR("attempting to copy data from non-existing RND cube %d,%d,%d", i, j, k);
music::elog.Print("attempting to copy data from non-existing RND cube %d,%d,%d", i, j, k);
throw std::runtime_error("attempting to copy data from non-existing RND cube");
}
@ -163,7 +163,7 @@ public:
if (it == cubemap_.end())
{
LOGERR("Attempting to copy data from non-existing RND cube %d,%d,%d @ %d,%d,%d", ic, jc, kc, i, j, k);
music::elog.Print("Attempting to copy data from non-existing RND cube %d,%d,%d @ %d,%d,%d", ic, jc, kc, i, j, k);
throw std::runtime_error("attempting to copy data from non-existing RND cube");
}
@ -171,7 +171,7 @@ public:
if (rnums_[cubeidx] == NULL)
{
LOGERR("Attempting to access data from non-allocated RND cube %d,%d,%d", ic, jc, kc);
music::elog.Print("Attempting to access data from non-allocated RND cube %d,%d,%d", ic, jc, kc);
throw std::runtime_error("attempting to access data from non-allocated RND cube");
}

83
src/plugins/random_panphasia.cc
View file

@ -156,7 +156,7 @@ protected:
public:
explicit RNG_panphasia(config_file &cf) : RNG_plugin(cf)
{
descriptor_string_ = pcf_->getValue<std::string>("random", "descriptor");
descriptor_string_ = pcf_->get_value<std::string>("random", "descriptor");
#ifdef _OPENMP
num_threads_ = omp_get_max_threads();
@ -169,7 +169,7 @@ public:
// parse the descriptor for its properties
pdescriptor_ = new panphasia_descriptor(descriptor_string_);
LOGINFO("PANPHASIA: descriptor \'%s\' is base %d,", pdescriptor_->name.c_str(), pdescriptor_->i_base);
music::ilog.Print("PANPHASIA: descriptor \'%s\' is base %d,", pdescriptor_->name.c_str(), pdescriptor_->i_base);
// write panphasia base size into config file for the grid construction
// as the gridding unit we use the least common multiple of 2 and i_base
@ -177,19 +177,19 @@ public:
// ARJ ss << lcm(2, pdescriptor_->i_base);
// ss << two_or_largest_power_two_less_than(pdescriptor_->i_base);//ARJ
ss << 2; // ARJ - set gridding unit to two
pcf_->insertValue("setup", "gridding_unit", ss.str());
pcf_->insert_value("setup", "gridding_unit", ss.str());
ss.str(std::string());
ss << pdescriptor_->i_base;
pcf_->insertValue("random", "base_unit", ss.str());
pcf_->insert_value("random", "base_unit", ss.str());
pcf_->insertValue("setup","fourier_splicing","false");
pcf_->insert_value("setup","fourier_splicing","false");
}
void initialize_for_grid_structure(const refinement_hierarchy &refh)
{
prefh_ = &refh;
levelmin_ = prefh_->levelmin();
levelmin_final_ = pcf_->getValue<unsigned>("setup", "levelmin");
levelmin_final_ = pcf_->get_value<unsigned>("setup", "levelmin");
levelmax_ = prefh_->levelmax();
if( refh.get_margin() < 0 ){
@ -199,7 +199,7 @@ public:
}
clear_panphasia_thread_states();
LOGINFO("PANPHASIA: running with %d threads", num_threads_);
music::ilog.Print("PANPHASIA: running with %d threads", num_threads_);
// if ngrid is not a multiple of i_base, then we need to enlarge and then sample down
ngrid_ = 1 << levelmin_;
@ -211,9 +211,9 @@ public:
int ratio = 1 << lextra_;
grid_rescale_fac_ = 1.0;
coordinate_system_shift_[0] = -pcf_->getValue<int>("setup", "shift_x");
coordinate_system_shift_[1] = -pcf_->getValue<int>("setup", "shift_y");
coordinate_system_shift_[2] = -pcf_->getValue<int>("setup", "shift_z");
coordinate_system_shift_[0] = -pcf_->get_value<int>("setup", "shift_x");
coordinate_system_shift_[1] = -pcf_->get_value<int>("setup", "shift_y");
coordinate_system_shift_[2] = -pcf_->get_value<int>("setup", "shift_z");
incongruent_fields_ = false;
if (ngrid_ != ratio * pdescriptor_->i_base)
@ -221,7 +221,7 @@ public:
incongruent_fields_ = true;
ngrid_ = 2 * ratio * pdescriptor_->i_base;
grid_rescale_fac_ = (double)ngrid_ / (1 << levelmin_);
LOGINFO("PANPHASIA: will use a higher resolution:\n"
music::ilog.Print("PANPHASIA: will use a higher resolution:\n"
" (%d -> %d) * 2**ref compatible with PANPHASIA\n"
" will Fourier interpolate after.",
grid_m_, grid_p_);
@ -361,21 +361,21 @@ void RNG_panphasia::fill_grid(int level, DensityGrid<real_t> &R)
}
if ((nx_m[0] != nx_m[1]) || (nx_m[0] != nx_m[2]))
LOGERR("Fatal error: non-cubic refinement being requested");
music::elog.Print("Fatal error: non-cubic refinement being requested");
inter_grid_phase_adjustment_ = M_PI * (1.0 / (double)nx_m[0] - 1.0 / (double)nxremap[0]);
// LOGINFO("The value of the phase adjustement is %f\n", inter_grid_phase_adjustment_);
// music::ilog.Print("The value of the phase adjustement is %f\n", inter_grid_phase_adjustment_);
// LOGINFO("ileft[0],ileft[1],ileft[2] %d %d %d", ileft[0], ileft[1], ileft[2]);
// LOGINFO("ileft_corner[0,1,2] %d %d %d", ileft_corner[0], ileft_corner[1], ileft_corner[2]);
// music::ilog.Print("ileft[0],ileft[1],ileft[2] %d %d %d", ileft[0], ileft[1], ileft[2]);
// music::ilog.Print("ileft_corner[0,1,2] %d %d %d", ileft_corner[0], ileft_corner[1], ileft_corner[2]);
// LOGINFO("iexpand_left[1,2,3] = (%d, %d, %d) Max %d ",iexpand_left[0],iexpand_left[1],iexpand_left[2], ileft_max_expand);
// music::ilog.Print("iexpand_left[1,2,3] = (%d, %d, %d) Max %d ",iexpand_left[0],iexpand_left[1],iexpand_left[2], ileft_max_expand);
// LOGINFO("ileft_corner_m[0,1,2] = (%d,%d,%d)",ileft_corner_m[0],ileft_corner_m[1],ileft_corner_m[2]);
// LOGINFO("grid_m_ %d grid_p_ %d",grid_m_,grid_p_);
// LOGINFO("nx_m[0,1,2] = (%d,%d,%d)",nx_m[0],nx_m[1],nx_m[2]);
// LOGINFO("ileft_corner_p[0,1,2] = (%d,%d,%d)",ileft_corner_p[0],ileft_corner_p[1],ileft_corner_p[2]);
// LOGINFO("nxremap[0,1,2] = (%d,%d,%d)",nxremap[0],nxremap[1],nxremap[2]);
// music::ilog.Print("ileft_corner_m[0,1,2] = (%d,%d,%d)",ileft_corner_m[0],ileft_corner_m[1],ileft_corner_m[2]);
// music::ilog.Print("grid_m_ %d grid_p_ %d",grid_m_,grid_p_);
// music::ilog.Print("nx_m[0,1,2] = (%d,%d,%d)",nx_m[0],nx_m[1],nx_m[2]);
// music::ilog.Print("ileft_corner_p[0,1,2] = (%d,%d,%d)",ileft_corner_p[0],ileft_corner_p[1],ileft_corner_p[2]);
// music::ilog.Print("nxremap[0,1,2] = (%d,%d,%d)",nxremap[0],nxremap[1],nxremap[2]);
size_t ngp = nxremap[0] * nxremap[1] * (nxremap[2] + 2);
@ -391,7 +391,7 @@ void RNG_panphasia::fill_grid(int level, DensityGrid<real_t> &R)
pc3 = reinterpret_cast<fftw_complex *>(pr3);
pc4 = reinterpret_cast<fftw_complex *>(pr4);
LOGINFO("calculating PANPHASIA random numbers for level %d...", level);
music::ilog.Print("calculating PANPHASIA random numbers for level %d...", level);
clear_panphasia_thread_states();
double t1 = get_wtime();
@ -424,8 +424,7 @@ void RNG_panphasia::fill_grid(int level, DensityGrid<real_t> &R)
lextra = (log10((double)ng_level / (double)d.i_base) + 0.001) / log10(2.0);
level_p = d.wn_level_base + lextra;
int ratio = 1 << lextra;
assert(ng_level == ratio * d.i_base);
assert(ng_level == (1 << lextra) * d.i_base);
ix_rel[0] = ileft_corner_p[0];
ix_rel[1] = ileft_corner_p[1];
@ -440,7 +439,7 @@ void RNG_panphasia::fill_grid(int level, DensityGrid<real_t> &R)
set_phases_and_rel_origin_(&lstate[mythread], descriptor, &level_p, &ix_rel[0], &ix_rel[1], &ix_rel[2],
&verbosity);
// LOGUSER(" called set_phases_and_rel_origin level %d ix_rel iy_rel iz_rel %d %d %d\n", level_p, ix_rel[0], ix_rel[1], ix_rel[2]);
// music::ulog.Print(" called set_phases_and_rel_origin level %d ix_rel iy_rel iz_rel %d %d %d\n", level_p, ix_rel[0], ix_rel[1], ix_rel[2]);
odd_x = ix_rel[0] % 2;
odd_y = ix_rel[1] % 2;
@ -493,14 +492,14 @@ void RNG_panphasia::fill_grid(int level, DensityGrid<real_t> &R)
}
}
}
LOGUSER("time for calculating PANPHASIA for level %d : %f s, %f µs/cell", level, get_wtime() - t1,
music::ulog.Print("time for calculating PANPHASIA for level %d : %f s, %f µs/cell", level, get_wtime() - t1,
1e6 * (get_wtime() - t1) / ((double)nxremap[2] * (double)nxremap[1] * (double)nxremap[0]));
LOGUSER("time for calculating PANPHASIA for level %d : %f s, %f µs/cell", level, get_wtime() - t1,
music::ulog.Print("time for calculating PANPHASIA for level %d : %f s, %f µs/cell", level, get_wtime() - t1,
1e6 * (get_wtime() - t1) / ((double)nxremap[2] * (double)nxremap[1] * (double)nxremap[0]));
//////////////////////////////////////////////////////////////////////////////////////////////
LOGUSER("\033[31mtiming level %d [adv_panphasia_cell_properties]: %f s\033[0m", level, get_wtime() - tp);
music::ulog.Print("\033[31mtiming level %d [adv_panphasia_cell_properties]: %f s\033[0m", level, get_wtime() - tp);
tp = get_wtime();
/////////////////////////////////////////////////////////////////////////
@ -585,7 +584,7 @@ void RNG_panphasia::fill_grid(int level, DensityGrid<real_t> &R)
// END
LOGUSER("\033[31mtiming level %d [build panphasia field]: %f s\033[0m", level, get_wtime() - tp);
music::ulog.Print("\033[31mtiming level %d [build panphasia field]: %f s\033[0m", level, get_wtime() - tp);
tp = get_wtime();
//***************************************************************
@ -618,8 +617,8 @@ void RNG_panphasia::fill_grid(int level, DensityGrid<real_t> &R)
lextra = (log10((double)ng_level / (double)d.i_base) + 0.001) / log10(2.0);
level_p = d.wn_level_base + lextra;
int ratio = 1 << lextra;
assert(ng_level == ratio * d.i_base);
assert(ng_level == (1 << lextra) * d.i_base);
ix_rel[0] = ileft_corner_p[0];
ix_rel[1] = ileft_corner_p[1];
@ -634,7 +633,7 @@ void RNG_panphasia::fill_grid(int level, DensityGrid<real_t> &R)
set_phases_and_rel_origin_(&lstate[mythread], descriptor, &level_p, &ix_rel[0], &ix_rel[1], &ix_rel[2],
&verbosity);
// LOGINFO(" called set_phases_and_rel_origin level %d ix_rel iy_rel iz_rel %d %d %d\n", level_p, ix_rel[0], ix_rel[1], ix_rel[2]);
// music::ilog.Print(" called set_phases_and_rel_origin level %d ix_rel iy_rel iz_rel %d %d %d\n", level_p, ix_rel[0], ix_rel[1], ix_rel[2]);
odd_x = ix_rel[0] % 2;
odd_y = ix_rel[1] % 2;
@ -685,10 +684,10 @@ void RNG_panphasia::fill_grid(int level, DensityGrid<real_t> &R)
}
}
}
LOGINFO("time for calculating PANPHASIA for level %d : %f s, %f µs/cell", level, get_wtime() - t1,
music::ilog.Print("time for calculating PANPHASIA for level %d : %f s, %f µs/cell", level, get_wtime() - t1,
1e6 * (get_wtime() - t1) / ((double)nxremap[2] * (double)nxremap[1] * (double)nxremap[0]));
LOGUSER("\033[31mtiming level %d [adv_panphasia_cell_properties2]: %f s \033[0m", level, get_wtime() - tp);
music::ulog.Print("\033[31mtiming level %d [adv_panphasia_cell_properties2]: %f s \033[0m", level, get_wtime() - tp);
tp = get_wtime();
/////////////////////////////////////////////////////////////////////////
@ -766,7 +765,7 @@ void RNG_panphasia::fill_grid(int level, DensityGrid<real_t> &R)
}
}
LOGUSER("\033[31mtiming level %d [build panphasia field2]: %f s\033[0m", level, get_wtime() - tp);
music::ulog.Print("\033[31mtiming level %d [build panphasia field2]: %f s\033[0m", level, get_wtime() - tp);
tp = get_wtime();
// END
@ -782,7 +781,7 @@ void RNG_panphasia::fill_grid(int level, DensityGrid<real_t> &R)
if (incongruent_fields_)
{
LOGUSER("Remapping fields from dimension %d -> %d", nxremap[0], nx_m[0]);
music::ulog.Print("Remapping fields from dimension %d -> %d", nxremap[0], nx_m[0]);
memset(pr1, 0, ngp * sizeof(fftw_real));
#pragma omp parallel for
@ -818,12 +817,12 @@ void RNG_panphasia::fill_grid(int level, DensityGrid<real_t> &R)
// if (level == 9)
// {
// LOGUSER("DC mode of level is %g", RE(pc0[0]));
// music::ulog.Print("DC mode of level is %g", RE(pc0[0]));
// // RE(pc0[0]) = 1e8;
// // IM(pc0[0]) = 0.0;
// }
LOGUSER("\033[31mtiming level %d [remap noncongruent]: %f s\033[0m", level, get_wtime() - tp);
music::ulog.Print("\033[31mtiming level %d [remap noncongruent]: %f s\033[0m", level, get_wtime() - tp);
tp = get_wtime();
/////////////////////////////////////////////////////////////////////////
// transform back
@ -837,7 +836,7 @@ void RNG_panphasia::fill_grid(int level, DensityGrid<real_t> &R)
delete[] pr3;
delete[] pr4;
LOGUSER("Copying random field data %d,%d,%d -> %d,%d,%d", nxremap[0], nxremap[1], nxremap[2], nx[0], nx[1], nx[2]);
music::ulog.Print("Copying random field data %d,%d,%d -> %d,%d,%d", nxremap[0], nxremap[1], nxremap[2], nx[0], nx[1], nx[2]);
// n = 1<<level;
// ng = n;
@ -873,10 +872,10 @@ void RNG_panphasia::fill_grid(int level, DensityGrid<real_t> &R)
sum2 = (sum2 - sum * sum);
LOGUSER("done with PANPHASIA for level %d:\n mean=%g, var=%g", level, sum, sum2);
LOGUSER("Copying into R array: nx[0],nx[1],nx[2] %d %d %d \n", nx[0], nx[1], nx[2]);
music::ulog.Print("done with PANPHASIA for level %d:\n mean=%g, var=%g", level, sum, sum2);
music::ulog.Print("Copying into R array: nx[0],nx[1],nx[2] %d %d %d \n", nx[0], nx[1], nx[2]);
LOGINFO("PANPHASIA level %d mean and variance are\n <p> = %g | var(p) = %g", level, sum, sum2);
music::ilog.Print("PANPHASIA level %d mean and variance are\n <p> = %g | var(p) = %g", level, sum, sum2);
}
namespace

22
src/plugins/region_convex_hull.cc
View file

@ -38,7 +38,7 @@ private:
void apply_shift( size_t Np, double *p, int *shift, int levelmin )
{
double dx = 1.0/(1<<levelmin);
LOGINFO("unapplying shift of previous zoom region to region particles :\n" \
music::ilog.Print("unapplying shift of previous zoom region to region particles :\n" \
"\t [%d,%d,%d] = (%f,%f,%f)",shift[0],shift[1],shift[2],shift[0]*dx,shift[1]*dx,shift[2]*dx);
for( size_t i=0,i3=0; i<Np; i++,i3+=3 )
@ -52,23 +52,23 @@ public:
{
std::vector<double> pp;
vfac_ = cf.getValue<double>("cosmology","vfact");
padding_ = cf.getValue<int>("setup","padding");
vfac_ = cf.get_value<double>("cosmology","vfact");
padding_ = cf.get_value<int>("setup","padding");
std::string point_file = cf.getValue<std::string>("setup","region_point_file");
std::string point_file = cf.get_value<std::string>("setup","region_point_file");
point_reader pfr;
pfr.read_points_from_file( point_file, vfac_, pp );
if( cf.containsKey("setup","region_point_shift") )
if( cf.contains_key("setup","region_point_shift") )
{
std::string point_shift = cf.getValue<std::string>("setup","region_point_shift");
std::string point_shift = cf.get_value<std::string>("setup","region_point_shift");
if(sscanf( point_shift.c_str(), "%d,%d,%d", &shift[0],&shift[1],&shift[2] )!=3){
LOGERR("Error parsing triple for region_point_shift");
music::elog.Print("Error parsing triple for region_point_shift");
throw std::runtime_error("Error parsing triple for region_point_shift");
}
unsigned point_levelmin = cf.getValue<unsigned>("setup","region_point_levelmin");
unsigned point_levelmin = cf.get_value<unsigned>("setup","region_point_levelmin");
apply_shift( pp.size()/3, &pp[0], shift, point_levelmin );
shift_level = point_levelmin;
@ -90,13 +90,13 @@ public:
phull_ = new convex_hull<double>( &pp[0], pp.size()/3, anchor_pt_ );
//expand the ellipsoid by one grid cell
unsigned levelmax = cf.getValue<unsigned>("setup","levelmax");
unsigned levelmax = cf.get_value<unsigned>("setup","levelmax");
double dx = 1.0/(1ul<<levelmax);
phull_->expand( sqrt(3.)*dx );
// output the center
double c[3] = { phull_->centroid_[0], phull_->centroid_[1], phull_->centroid_[2] };
LOGINFO("Region center from convex hull centroid determined at\n\t (%f,%f,%f)",c[0],c[1],c[2]);
music::ilog.Print("Region center from convex hull centroid determined at\n\t (%f,%f,%f)",c[0],c[1],c[2]);
//-----------------------------------------------------------------
// when querying the bounding box, do we need extra padding?
@ -104,7 +104,7 @@ public:
// conditions should be added here
{
std::string output_plugin = cf.getValue<std::string>("output","format");
std::string output_plugin = cf.get_value<std::string>("output","format");
if( output_plugin == std::string("grafic2") )
do_extra_padding_ = true;
}

66
src/plugins/region_ellipsoid.cc
View file

@ -295,7 +295,7 @@ protected:
}
if (count >= maxit)
LOGERR("No convergence in min_ellipsoid::compute: maximum number of iterations reached!");
music::elog.Print("No convergence in min_ellipsoid::compute: maximum number of iterations reached!");
delete[] unew;
}
@ -305,7 +305,7 @@ public:
: N(N_), axes_computed(false), hold_point_data(true)
{
// --- initialize ---
LOGINFO("computing minimum bounding ellipsoid from %lld points", N);
music::ilog.Print("computing minimum bounding ellipsoid from %lld points", N);
Q = new float[4 * N];
u = new float[N];
@ -536,7 +536,7 @@ public:
{
if (!axes_computed)
{
LOGUSER("computing ellipsoid axes.....");
music::ulog.Print("computing ellipsoid axes.....");
compute_axes();
}
// float muold[3] = {mu[0],mu[1],mu[2]};
@ -558,7 +558,7 @@ public:
Inverse_3x3(A, Ainv);
// LOGUSER("computing ellipsoid axes.....");
// music::ulog.Print("computing ellipsoid axes.....");
compute_axes();
// print();
@ -583,7 +583,7 @@ private:
void apply_shift(size_t Np, double *p, int *shift, int levelmin)
{
double dx = 1.0 / (1 << levelmin);
LOGINFO("unapplying shift of previous zoom region to region particles :\n"
music::ilog.Print("unapplying shift of previous zoom region to region particles :\n"
"\t [%d,%d,%d] = (%f,%f,%f)",
shift[0], shift[1], shift[2], shift[0] * dx, shift[1] * dx, shift[2] * dx);
@ -602,25 +602,25 @@ public:
pellip_.push_back(NULL);
// sanity check
if (!cf.containsKey("setup", "region_point_file") &&
!(cf.containsKey("setup", "region_ellipsoid_matrix[0]") &&
cf.containsKey("setup", "region_ellipsoid_matrix[1]") &&
cf.containsKey("setup", "region_ellipsoid_matrix[2]") &&
cf.containsKey("setup", "region_ellipsoid_center")))
if (!cf.contains_key("setup", "region_point_file") &&
!(cf.contains_key("setup", "region_ellipsoid_matrix[0]") &&
cf.contains_key("setup", "region_ellipsoid_matrix[1]") &&
cf.contains_key("setup", "region_ellipsoid_matrix[2]") &&
cf.contains_key("setup", "region_ellipsoid_center")))
{
LOGERR("Insufficient data for region=ellipsoid\n Need to specify either \'region_point_file\' or the ellipsoid equation.");
music::elog.Print("Insufficient data for region=ellipsoid\n Need to specify either \'region_point_file\' or the ellipsoid equation.");
throw std::runtime_error("Insufficient data for region=ellipsoid");
}
//
vfac_ = cf.getValue<double>("cosmology", "vfact");
padding_ = cf.getValue<int>("setup", "padding");
vfac_ = cf.get_value<double>("cosmology", "vfact");
padding_ = cf.get_value<int>("setup", "padding");
std::string point_file;
if (cf.containsKey("setup", "region_point_file"))
if (cf.contains_key("setup", "region_point_file"))
{
point_file = cf.getValue<std::string>("setup", "region_point_file");
point_file = cf.get_value<std::string>("setup", "region_point_file");
point_reader pfr;
pfr.read_points_from_file(point_file, vfac_, pp);
@ -639,15 +639,15 @@ public:
pp.swap(xx);
}
if (cf.containsKey("setup", "region_point_shift"))
if (cf.contains_key("setup", "region_point_shift"))
{
std::string point_shift = cf.getValue<std::string>("setup", "region_point_shift");
std::string point_shift = cf.get_value<std::string>("setup", "region_point_shift");
if (sscanf(point_shift.c_str(), "%d,%d,%d", &shift[0], &shift[1], &shift[2]) != 3)
{
LOGERR("Error parsing triple for region_point_shift");
music::elog.Print("Error parsing triple for region_point_shift");
throw std::runtime_error("Error parsing triple for region_point_shift");
}
unsigned point_levelmin = cf.getValue<unsigned>("setup", "region_point_levelmin");
unsigned point_levelmin = cf.get_value<unsigned>("setup", "region_point_levelmin");
apply_shift(pp.size() / 3, &pp[0], shift, point_levelmin);
shift_level = point_levelmin;
@ -660,28 +660,28 @@ public:
double A[9] = {0}, c[3] = {0};
std::string strtmp;
strtmp = cf.getValue<std::string>("setup", "region_ellipsoid_matrix[0]");
strtmp = cf.get_value<std::string>("setup", "region_ellipsoid_matrix[0]");
if (sscanf(strtmp.c_str(), "%lf,%lf,%lf", &A[0], &A[1], &A[2]) != 3)
{
LOGERR("Error parsing triple for region_ellipsoid_matrix[0]");
music::elog.Print("Error parsing triple for region_ellipsoid_matrix[0]");
throw std::runtime_error("Error parsing triple for region_ellipsoid_matrix[0]");
}
strtmp = cf.getValue<std::string>("setup", "region_ellipsoid_matrix[1]");
strtmp = cf.get_value<std::string>("setup", "region_ellipsoid_matrix[1]");
if (sscanf(strtmp.c_str(), "%lf,%lf,%lf", &A[3], &A[4], &A[5]) != 3)
{
LOGERR("Error parsing triple for region_ellipsoid_matrix[1]");
music::elog.Print("Error parsing triple for region_ellipsoid_matrix[1]");
throw std::runtime_error("Error parsing triple for region_ellipsoid_matrix[1]");
}
strtmp = cf.getValue<std::string>("setup", "region_ellipsoid_matrix[2]");
strtmp = cf.get_value<std::string>("setup", "region_ellipsoid_matrix[2]");
if (sscanf(strtmp.c_str(), "%lf,%lf,%lf", &A[6], &A[7], &A[8]) != 3)
{
LOGERR("Error parsing triple for region_ellipsoid_matrix[2]");
music::elog.Print("Error parsing triple for region_ellipsoid_matrix[2]");
throw std::runtime_error("Error parsing triple for region_ellipsoid_matrix[2]");
}
strtmp = cf.getValue<std::string>("setup", "region_ellipsoid_center");
strtmp = cf.get_value<std::string>("setup", "region_ellipsoid_center");
if (sscanf(strtmp.c_str(), "%lf,%lf,%lf", &c[0], &c[1], &c[2]) != 3)
{
LOGERR("Error parsing triple for region_ellipsoid_center");
music::elog.Print("Error parsing triple for region_ellipsoid_center");
throw std::runtime_error("Error parsing triple for region_ellipsoid_center");
}
@ -693,7 +693,7 @@ public:
{
// compute convex hull and use only hull points to speed things up
// BUT THIS NEEDS MORE TESTING BEFORE IT GOES IN THE REPO
LOGINFO("Computing convex hull for %llu points", pp.size()/3 );
music::ilog.Print("Computing convex hull for %llu points", pp.size()/3 );
convex_hull<double> ch( &pp[0], pp.size()/3 );
std::set<int> unique;
ch.get_defining_indices( unique );
@ -723,14 +723,14 @@ public:
pellip_[levelmax_]->get_center(c);
pellip_[levelmax_]->get_matrix(A);
LOGINFO("Region center for ellipsoid determined at\n\t xc = ( %f %f %f )", c[0], c[1], c[2]);
LOGINFO("Ellipsoid matrix determined as\n\t ( %f %f %f )\n\t A = ( %f %f %f )\n\t ( %f %f %f )",
music::ilog.Print("Region center for ellipsoid determined at\n\t xc = ( %f %f %f )", c[0], c[1], c[2]);
music::ilog.Print("Ellipsoid matrix determined as\n\t ( %f %f %f )\n\t A = ( %f %f %f )\n\t ( %f %f %f )",
A[0], A[1], A[2], A[3], A[4], A[5], A[6], A[7], A[8]);
// expand the ellipsoid by one grid cell
unsigned levelmax = cf.getValue<unsigned>("setup", "levelmax");
unsigned npad = cf.getValue<unsigned>("setup", "padding");
unsigned levelmax = cf.get_value<unsigned>("setup", "levelmax");
unsigned npad = cf.get_value<unsigned>("setup", "padding");
double dx = 1.0 / (1ul << levelmax);
pellip_[levelmax_]->expand_ellipsoid(dx);
@ -748,7 +748,7 @@ public:
// conditions should be added here
{
std::string output_plugin = cf.getValue<std::string>("output", "format");
std::string output_plugin = cf.get_value<std::string>("output", "format");
if (output_plugin == std::string("grafic2"))
do_extra_padding_ = true;
}

20
src/plugins/region_multibox.cc
View file

@ -136,34 +136,34 @@ public:
{
res = 1<<(levelmin_-1);
//check parameters
if ( !cf.containsKey("setup", "region_point_file"))
if ( !cf.contains_key("setup", "region_point_file"))
{
LOGERR("Missing parameter \'region_point_file\' needed for region=multibox");
music::elog.Print("Missing parameter \'region_point_file\' needed for region=multibox");
throw std::runtime_error("Missing parameter for region=multibox");
}
if ( cf.containsKey("setup", "region_point_level"))
if ( cf.contains_key("setup", "region_point_level"))
{
res = 1<<(levelmin_-1);
res = 1<<(cf.getValue<int>("setup","region_point_level")-1);
res = 1<<(cf.get_value<int>("setup","region_point_level")-1);
}
vfac_ = cf.getValue<double>("cosmology","vfact");
vfac_ = cf.get_value<double>("cosmology","vfact");
if (levelmin_ >= levelmax_)
{
LOGERR("Why are you specifying a region if you aren't refining?");
music::elog.Print("Why are you specifying a region if you aren't refining?");
throw std::runtime_error("region=multibox needs levelmax>levelmin");
}
std::vector<double> pp;
point_reader pfr;
pfr.read_points_from_file(cf.getValue<std::string>("setup", "region_point_file"), vfac_, pp);
LOGINFO("Multibox Resolution: %d", res);
pfr.read_points_from_file(cf.get_value<std::string>("setup", "region_point_file"), vfac_, pp);
music::ilog.Print("Multibox Resolution: %d", res);
//Initialize the grid with zeros, the base level
refgrid = grid(res,slice(res,col(res,levelmin_)));
//Build the highest refinement region
deposit_particles(pp, res);
LOGINFO("Deposited Multigrid Particles");
music::ilog.Print("Deposited Multigrid Particles");
//Build the intermediate refinement regions
build_refgrid();
LOGINFO("Built Multigrid");
music::ilog.Print("Built Multigrid");
get_center(cen);
dump_grid();
}

4
src/plugins/transfer_bbks.cc
View file

@ -35,12 +35,12 @@ public:
bool bSugiyama(true);
try{
bSugiyama= pcf_->getValue<bool>( "cosmology", "sugiyama_corr" );
bSugiyama= pcf_->get_value<bool>( "cosmology", "sugiyama_corr" );
}catch(...){
throw std::runtime_error("Error in \'tranfer_bbks_plugin\': need to specify \'[cosmology]/sugiyama_corr = [true/false]");
}
FreeGamma = pcf_->getValueSafe<double>( "cosmology", "gamma", FreeGamma );
FreeGamma = pcf_->get_value_safe<double>( "cosmology", "gamma", FreeGamma );
if( FreeGamma <= 0.0 ){
m_Gamma = Omega0*0.01*cosmo_.H0;

18
src/plugins/transfer_camb.cc
View file

@ -35,10 +35,10 @@ private:
#ifdef WITH_MPI
if (MPI::COMM_WORLD.Get_rank() == 0) {
#endif
LOGINFO("Reading tabulated transfer function data from file \n \'%s\'",
music::ilog.Print("Reading tabulated transfer function data from file \n \'%s\'",
m_filename_Tk.c_str());
LOGINFO("CAUTION: make sure that this transfer function \n\t has been output for z=%f!",m_zstart);
music::ilog.Print("CAUTION: make sure that this transfer function \n\t has been output for z=%f!",m_zstart);
std::string line;
std::ifstream ifs(m_filename_Tk.c_str());
@ -85,7 +85,7 @@ private:
ss >> dummy; // v_b-v_cdm
if( ss.bad() || ss.fail() ){
LOGERR("Error reading the transfer function file (corrupt or not in expected format)!");
music::elog.Print("Error reading the transfer function file (corrupt or not in expected format)!");
throw std::runtime_error("Error reading transfer function file \'" +
m_filename_Tk + "\'");
}
@ -126,10 +126,10 @@ private:
ifs.close();
LOGINFO("Read CAMB transfer function table with %d rows", m_nlines);
music::ilog.Print("Read CAMB transfer function table with %d rows", m_nlines);
if (m_linbaryoninterp)
LOGINFO("Using log-lin interpolation for baryons\n (TF is not "
music::ilog.Print("Using log-lin interpolation for baryons\n (TF is not "
"positive definite)");
#ifdef WITH_MPI
@ -165,10 +165,10 @@ public:
transfer_CAMB_plugin(config_file &cf)
: transfer_function_plugin(cf)
{
m_filename_Tk = pcf_->getValue<std::string>("cosmology", "transfer_file");
m_Omega_m=cf.getValue<double>("cosmology","Omega_m"); //MvD
m_Omega_b=cf.getValue<double>("cosmology","Omega_b"); //MvD
m_zstart =cf.getValue<double>("setup","zstart"); //MvD
m_filename_Tk = pcf_->get_value<std::string>("cosmology", "transfer_file");
m_Omega_m=cf.get_value<double>("cosmology","Omega_m"); //MvD
m_Omega_b=cf.get_value<double>("cosmology","Omega_b"); //MvD
m_zstart =cf.get_value<double>("setup","zstart"); //MvD
read_table();

42
src/plugins/transfer_eisenstein.cc
View file

@ -198,7 +198,7 @@ public:
transfer_eisenstein_plugin( config_file &cf )//Cosmology aCosm, double Tcmb = 2.726 )
: transfer_function_plugin(cf)
{
double Tcmb = pcf_->getValueSafe("cosmology","Tcmb",2.726);
double Tcmb = pcf_->get_value_safe("cosmology","Tcmb",2.726);
etf_.set_parameters( cosmo_, Tcmb );
@ -239,20 +239,20 @@ public:
transfer_eisenstein_wdm_plugin( config_file &cf )
: transfer_function_plugin(cf), m_h0( cosmo_.H0*0.01 )
{
double Tcmb = pcf_->getValueSafe("cosmology","Tcmb",2.726);
double Tcmb = pcf_->get_value_safe("cosmology","Tcmb",2.726);
etf_.set_parameters( cosmo_, Tcmb );
typemap_.insert( std::pair<std::string,int>( "BODE", wdm_bode ) );
typemap_.insert( std::pair<std::string,int>( "VIEL", wdm_viel ) ); // add the other types
typemap_.insert( std::pair<std::string,int>( "BODE_WRONG", wdm_bode_wrong ) ); // add the other types
omegam_ = cf.getValue<double>("cosmology","Omega_m");
omegab_ = cf.getValue<double>("cosmology","Omega_b");
wdmm_ = cf.getValue<double>("cosmology","WDMmass");
omegam_ = cf.get_value<double>("cosmology","Omega_m");
omegab_ = cf.get_value<double>("cosmology","Omega_b");
wdmm_ = cf.get_value<double>("cosmology","WDMmass");
H0_ = cf.getValue<double>("cosmology","H0");
type_ = cf.getValueSafe<std::string>("cosmology","WDMtftype","BODE");
H0_ = cf.get_value<double>("cosmology","H0");
type_ = cf.get_value_safe<std::string>("cosmology","WDMtftype","BODE");
//type_ = std::string( toupper( type_.c_str() ) );
@ -265,8 +265,8 @@ public:
{
//... parameterisation from Bode et al. (2001), ApJ, 556, 93
case wdm_bode:
wdmnu_ = cf.getValueSafe<double>("cosmology","WDMnu",1.0);
wdmgx_ = cf.getValueSafe<double>("cosmology","WDMg_x",1.5);
wdmnu_ = cf.get_value_safe<double>("cosmology","WDMnu",1.0);
wdmgx_ = cf.get_value_safe<double>("cosmology","WDMg_x",1.5);
m_WDMalpha = 0.05 * pow( omegam_/0.4,0.15)
*pow(H0_*0.01/0.65,1.3)*pow(wdmm_,-1.15)
*pow(1.5/wdmgx_,0.29);
@ -275,7 +275,7 @@ public:
//... parameterisation from Viel et al. (2005), Phys Rev D, 71
case wdm_viel:
wdmnu_ = cf.getValueSafe<double>("cosmology","WDMnu",1.12);
wdmnu_ = cf.get_value_safe<double>("cosmology","WDMnu",1.12);
m_WDMalpha = 0.049 * pow( omegam_/0.25,0.11)
*pow(H0_*0.01/0.7,1.22)*pow(wdmm_,-1.11);
break;
@ -284,16 +284,16 @@ public:
//.... below is for historical reasons due to the buggy parameterisation
//.... in early versions of MUSIC, but apart from H instead of h, Bode et al.
case wdm_bode_wrong:
wdmnu_ = cf.getValueSafe<double>("cosmology","WDMnu",1.0);
wdmgx_ = cf.getValueSafe<double>("cosmology","WDMg_x",1.5);
wdmnu_ = cf.get_value_safe<double>("cosmology","WDMnu",1.0);
wdmgx_ = cf.get_value_safe<double>("cosmology","WDMg_x",1.5);
m_WDMalpha = 0.05 * pow( omegam_/0.4,0.15)
*pow(H0_/0.65,1.3)*pow(wdmm_,-1.15)
*pow(1.5/wdmgx_,0.29);
break;
default:
wdmnu_ = cf.getValueSafe<double>("cosmology","WDMnu",1.0);
wdmgx_ = cf.getValueSafe<double>("cosmology","WDMg_x",1.5);
wdmnu_ = cf.get_value_safe<double>("cosmology","WDMnu",1.0);
wdmgx_ = cf.get_value_safe<double>("cosmology","WDMg_x",1.5);
m_WDMalpha = 0.05 * pow( omegam_/0.4,0.15)
*pow(H0_*0.01/0.65,1.3)*pow(wdmm_,-1.15)
*pow(1.5/wdmgx_,0.29);
@ -329,20 +329,20 @@ public:
transfer_eisenstein_cdmbino_plugin( config_file &cf )
: transfer_function_plugin(cf), m_h0( cosmo_.H0*0.01 )
{
double Tcmb = pcf_->getValueSafe("cosmology","Tcmb",2.726);
double Tcmb = pcf_->get_value_safe("cosmology","Tcmb",2.726);
etf_.set_parameters( cosmo_, Tcmb );
omegam_ = cf.getValue<double>("cosmology","Omega_m");
omegab_ = cf.getValue<double>("cosmology","Omega_b");
H0_ = cf.getValue<double>("cosmology","H0");
omegam_ = cf.get_value<double>("cosmology","Omega_m");
omegab_ = cf.get_value<double>("cosmology","Omega_b");
H0_ = cf.get_value<double>("cosmology","H0");
mcdm_ = cf.getValueSafe<double>("cosmology","CDM_mass", 100.0); // bino particle mass in GeV
Tkd_ = cf.getValueSafe<double>("cosmology","CDM_Tkd", 33.0); // temperature at which CDM particle kinetically decouples (in MeV)
mcdm_ = cf.get_value_safe<double>("cosmology","CDM_mass", 100.0); // bino particle mass in GeV
Tkd_ = cf.get_value_safe<double>("cosmology","CDM_Tkd", 33.0); // temperature at which CDM particle kinetically decouples (in MeV)
kfs_ = 1.7e6 / m_h0 * sqrt( mcdm_ / 100. * Tkd_ / 30. ) / (1.0 + log( Tkd_ / 30. ) / 19.2 );
kd_ = 3.8e7 / m_h0 * sqrt( mcdm_ / 100. * Tkd_ / 30. );
LOGINFO(" bino CDM: k_fs = %g, k_d = %g", kfs_, kd_ );
music::ilog.Print(" bino CDM: k_fs = %g, k_d = %g", kfs_, kd_ );
}

6
src/plugins/transfer_eisenstein_suppressLSS.cc
View file

@ -201,9 +201,9 @@ public:
transfer_eisensteinS_plugin( config_file &cf )//Cosmology aCosm, double Tcmb = 2.726 )
: transfer_function_plugin(cf)
{
double Tcmb = pcf_->getValueSafe<double>("cosmology","Tcmb",2.726);
//double boxlength = pcf_->getValue<double>("setup","boxlength");
ktrunc_ = pcf_->getValue<double>("cosmology","ktrunc");
double Tcmb = pcf_->get_value_safe<double>("cosmology","Tcmb",2.726);
//double boxlength = pcf_->get_value<double>("setup","boxlength");
ktrunc_ = pcf_->get_value<double>("cosmology","ktrunc");
normfac_ = 2.0/M_PI;
etf_.set_parameters( cosmo_, Tcmb );

2
src/plugins/transfer_inflation.cc
View file

@ -34,7 +34,7 @@ public:
transfer_inflation_plugin( config_file& cf )
: transfer_function_plugin( cf )
{
ns2_ = 0.5*cf.getValue<double>("cosmology","nspec");
ns2_ = 0.5*cf.get_value<double>("cosmology","nspec");
tf_distinct_ = true;
}

6
src/plugins/transfer_linger++.cc
View file

@ -124,15 +124,15 @@ public:
transfer_LINGERpp_plugin( config_file& cf )
: transfer_function_plugin( cf )
{
m_filename_Tk = pcf_->getValue<std::string>("cosmology","transfer_file");
m_filename_Tk = pcf_->get_value<std::string>("cosmology","transfer_file");
//.. disable the baryon-CDM relative velocity (both follow the total matter potential)
m_bnovrel = pcf_->getValueSafe<bool>("cosmology","no_vrel",false);
m_bnovrel = pcf_->get_value_safe<bool>("cosmology","no_vrel",false);
//.. normalize at z=0 rather than using the linearly scaled zini spectrum
//.. this can be different due to radiation still being non-negligible at
//.. high redshifts
m_bz0norm = pcf_->getValueSafe<bool>("cosmology","z0norm",true);
m_bz0norm = pcf_->get_value_safe<bool>("cosmology","z0norm",true);
tf_distinct_ = true;
tf_withvel_ = true;

2
src/plugins/transfer_music.cc
View file

@ -120,7 +120,7 @@ public:
transfer_MUSIC_plugin( config_file& cf )
: transfer_function_plugin( cf )
{
m_filename_Tk = pcf_->getValue<std::string>("cosmology","transfer_file");
m_filename_Tk = pcf_->get_value<std::string>("cosmology","transfer_file");
read_table( );

96
src/poisson.cc
View file

@ -32,7 +32,7 @@ void print_poisson_plugins()
if ((*it).second)
std::cout << "\t\'" << (*it).first << "\'\n";
LOGINFO("Poisson plug-in :: %s", std::string((*it).first).c_str());
music::ilog.Print("Poisson plug-in :: %s", std::string((*it).first).c_str());
++it;
}
@ -58,9 +58,9 @@ typedef multigrid::solver<stencil_19P<double>, interp_O7_fluxcorr, mg_straight,
double multigrid_poisson_plugin::solve(grid_hierarchy &f, grid_hierarchy &u)
{
LOGUSER("Initializing multi-grid Poisson solver...");
music::ulog.Print("Initializing multi-grid Poisson solver...");
unsigned verbosity = cf_.getValueSafe<unsigned>("setup", "verbosity", 2);
unsigned verbosity = cf_.get_value_safe<unsigned>("setup", "verbosity", 2);
if (verbosity > 0)
{
@ -72,32 +72,32 @@ double multigrid_poisson_plugin::solve(grid_hierarchy &f, grid_hierarchy &u)
std::string ps_smoother_name;
unsigned ps_presmooth, ps_postsmooth, order;
acc = cf_.getValueSafe<double>("poisson", "accuracy", acc);
ps_presmooth = cf_.getValueSafe<unsigned>("poisson", "pre_smooth", 3);
ps_postsmooth = cf_.getValueSafe<unsigned>("poisson", "post_smooth", 3);
ps_smoother_name = cf_.getValueSafe<std::string>("poisson", "smoother", "gs");
order = cf_.getValueSafe<unsigned>("poisson", "laplace_order", 4);
acc = cf_.get_value_safe<double>("poisson", "accuracy", acc);
ps_presmooth = cf_.get_value_safe<unsigned>("poisson", "pre_smooth", 3);
ps_postsmooth = cf_.get_value_safe<unsigned>("poisson", "post_smooth", 3);
ps_smoother_name = cf_.get_value_safe<std::string>("poisson", "smoother", "gs");
order = cf_.get_value_safe<unsigned>("poisson", "laplace_order", 4);
multigrid::opt::smtype ps_smtype = multigrid::opt::sm_gauss_seidel;
if (ps_smoother_name == std::string("gs"))
{
ps_smtype = multigrid::opt::sm_gauss_seidel;
LOGUSER("Selected Gauss-Seidel multigrid smoother");
music::ulog.Print("Selected Gauss-Seidel multigrid smoother");
}
else if (ps_smoother_name == std::string("jacobi"))
{
ps_smtype = multigrid::opt::sm_jacobi;
LOGUSER("Selected Jacobi multigrid smoother");
music::ulog.Print("Selected Jacobi multigrid smoother");
}
else if (ps_smoother_name == std::string("sor"))
{
ps_smtype = multigrid::opt::sm_sor;
LOGUSER("Selected SOR multigrid smoother");
music::ulog.Print("Selected SOR multigrid smoother");
}
else
{
LOGWARN("Unknown multigrid smoother \'%s\' specified. Reverting to Gauss-Seidel.", ps_smoother_name.c_str());
music::wlog.Print("Unknown multigrid smoother \'%s\' specified. Reverting to Gauss-Seidel.", ps_smoother_name.c_str());
std::cerr << " - Warning: unknown smoother \'" << ps_smoother_name << "\' for multigrid solver!\n"
<< " reverting to \'gs\' (Gauss-Seidel)" << std::endl;
}
@ -113,25 +113,25 @@ double multigrid_poisson_plugin::solve(grid_hierarchy &f, grid_hierarchy &u)
//----- run Poisson solver -----//
if (order == 2)
{
LOGUSER("Running multigrid solver with 2nd order Laplacian...");
music::ulog.Print("Running multigrid solver with 2nd order Laplacian...");
poisson_solver_O2 ps(f, ps_smtype, ps_presmooth, ps_postsmooth);
err = ps.solve(u, acc, true);
}
else if (order == 4)
{
LOGUSER("Running multigrid solver with 4th order Laplacian...");
music::ulog.Print("Running multigrid solver with 4th order Laplacian...");
poisson_solver_O4 ps(f, ps_smtype, ps_presmooth, ps_postsmooth);
err = ps.solve(u, acc, true);
}
else if (order == 6)
{
LOGUSER("Running multigrid solver with 6th order Laplacian..");
music::ulog.Print("Running multigrid solver with 6th order Laplacian..");
poisson_solver_O6 ps(f, ps_smtype, ps_presmooth, ps_postsmooth);
err = ps.solve(u, acc, true);
}
else
{
LOGERR("Invalid order specified for Laplace operator");
music::elog.Print("Invalid order specified for Laplace operator");
throw std::runtime_error("Invalid order specified for Laplace operator");
}
@ -155,7 +155,7 @@ double multigrid_poisson_plugin::gradient(int dir, grid_hierarchy &u, grid_hiera
{
Du = u;
unsigned order = cf_.getValueSafe<unsigned>("poisson", "grad_order", 4);
unsigned order = cf_.get_value_safe<unsigned>("poisson", "grad_order", 4);
switch (order)
{
@ -169,7 +169,7 @@ double multigrid_poisson_plugin::gradient(int dir, grid_hierarchy &u, grid_hiera
implementation().gradient_O6(dir, u, Du);
break;
default:
LOGERR("Invalid order %d specified for gradient operator", order);
music::elog.Print("Invalid order %d specified for gradient operator", order);
throw std::runtime_error("Invalid order specified for gradient operator!");
}
@ -180,7 +180,7 @@ double multigrid_poisson_plugin::gradient_add(int dir, grid_hierarchy &u, grid_h
{
// Du = u;
unsigned order = cf_.getValueSafe<unsigned>("poisson", "grad_order", 4);
unsigned order = cf_.get_value_safe<unsigned>("poisson", "grad_order", 4);
switch (order)
{
@ -194,7 +194,7 @@ double multigrid_poisson_plugin::gradient_add(int dir, grid_hierarchy &u, grid_h
implementation().gradient_add_O6(dir, u, Du);
break;
default:
LOGERR("Invalid order %d specified for gradient operator!", order);
music::elog.Print("Invalid order %d specified for gradient operator!", order);
throw std::runtime_error("Invalid order specified for gradient operator!");
}
@ -203,7 +203,7 @@ double multigrid_poisson_plugin::gradient_add(int dir, grid_hierarchy &u, grid_h
void multigrid_poisson_plugin::implementation::gradient_O2(int dir, grid_hierarchy &u, grid_hierarchy &Du)
{
LOGUSER("Computing a 2nd order finite difference gradient...");
music::ulog.Print("Computing a 2nd order finite difference gradient...");
for (unsigned ilevel = u.levelmin(); ilevel <= u.levelmax(); ++ilevel)
{
@ -232,12 +232,12 @@ void multigrid_poisson_plugin::implementation::gradient_O2(int dir, grid_hierarc
(*pvar)(ix, iy, iz) = 0.5 * ((*u.get_grid(ilevel))(ix, iy, iz + 1) - (*u.get_grid(ilevel))(ix, iy, iz - 1)) * h;
}
LOGUSER("Done computing a 2nd order finite difference gradient.");
music::ulog.Print("Done computing a 2nd order finite difference gradient.");
}
void multigrid_poisson_plugin::implementation::gradient_add_O2(int dir, grid_hierarchy &u, grid_hierarchy &Du)
{
LOGUSER("Computing a 2nd order finite difference gradient...");
music::ulog.Print("Computing a 2nd order finite difference gradient...");
for (unsigned ilevel = u.levelmin(); ilevel <= u.levelmax(); ++ilevel)
{
@ -266,12 +266,12 @@ void multigrid_poisson_plugin::implementation::gradient_add_O2(int dir, grid_hie
(*pvar)(ix, iy, iz) += 0.5 * ((*u.get_grid(ilevel))(ix, iy, iz + 1) - (*u.get_grid(ilevel))(ix, iy, iz - 1)) * h;
}
LOGUSER("Done computing a 4th order finite difference gradient.");
music::ulog.Print("Done computing a 4th order finite difference gradient.");
}
void multigrid_poisson_plugin::implementation::gradient_O4(int dir, grid_hierarchy &u, grid_hierarchy &Du)
{
LOGUSER("Computing a 4th order finite difference gradient...");
music::ulog.Print("Computing a 4th order finite difference gradient...");
for (unsigned ilevel = u.levelmin(); ilevel <= u.levelmax(); ++ilevel)
{
@ -302,12 +302,12 @@ void multigrid_poisson_plugin::implementation::gradient_O4(int dir, grid_hierarc
(*pvar)(ix, iy, iz) = ((*u.get_grid(ilevel))(ix, iy, iz - 2) - 8.0 * (*u.get_grid(ilevel))(ix, iy, iz - 1) + 8.0 * (*u.get_grid(ilevel))(ix, iy, iz + 1) - (*u.get_grid(ilevel))(ix, iy, iz + 2)) * h;
}
LOGUSER("Done computing a 4th order finite difference gradient.");
music::ulog.Print("Done computing a 4th order finite difference gradient.");
}
void multigrid_poisson_plugin::implementation::gradient_add_O4(int dir, grid_hierarchy &u, grid_hierarchy &Du)
{
LOGUSER("Computing a 4th order finite difference gradient...");
music::ulog.Print("Computing a 4th order finite difference gradient...");
for (unsigned ilevel = u.levelmin(); ilevel <= u.levelmax(); ++ilevel)
{
@ -338,12 +338,12 @@ void multigrid_poisson_plugin::implementation::gradient_add_O4(int dir, grid_hie
(*pvar)(ix, iy, iz) += ((*u.get_grid(ilevel))(ix, iy, iz - 2) - 8.0 * (*u.get_grid(ilevel))(ix, iy, iz - 1) + 8.0 * (*u.get_grid(ilevel))(ix, iy, iz + 1) - (*u.get_grid(ilevel))(ix, iy, iz + 2)) * h;
}
LOGUSER("Done computing a 4th order finite difference gradient.");
music::ulog.Print("Done computing a 4th order finite difference gradient.");
}
void multigrid_poisson_plugin::implementation::gradient_O6(int dir, grid_hierarchy &u, grid_hierarchy &Du)
{
LOGUSER("Computing a 6th order finite difference gradient...");
music::ulog.Print("Computing a 6th order finite difference gradient...");
for (unsigned ilevel = u.levelmin(); ilevel <= u.levelmax(); ++ilevel)
{
@ -376,12 +376,12 @@ void multigrid_poisson_plugin::implementation::gradient_O6(int dir, grid_hierarc
(-(*u.get_grid(ilevel))(ix, iy, iz - 3) + 9.0 * (*u.get_grid(ilevel))(ix, iy, iz - 2) - 45.0 * (*u.get_grid(ilevel))(ix, iy, iz - 1) + 45.0 * (*u.get_grid(ilevel))(ix, iy, iz + 1) - 9.0 * (*u.get_grid(ilevel))(ix, iy, iz + 2) + (*u.get_grid(ilevel))(ix, iy, iz + 3)) * h;
}
LOGUSER("Done computing a 6th order finite difference gradient.");
music::ulog.Print("Done computing a 6th order finite difference gradient.");
}
void multigrid_poisson_plugin::implementation::gradient_add_O6(int dir, grid_hierarchy &u, grid_hierarchy &Du)
{
LOGUSER("Computing a 6th order finite difference gradient...");
music::ulog.Print("Computing a 6th order finite difference gradient...");
for (unsigned ilevel = u.levelmin(); ilevel <= u.levelmax(); ++ilevel)
{
@ -414,7 +414,7 @@ void multigrid_poisson_plugin::implementation::gradient_add_O6(int dir, grid_hie
(-(*u.get_grid(ilevel))(ix, iy, iz - 3) + 9.0 * (*u.get_grid(ilevel))(ix, iy, iz - 2) - 45.0 * (*u.get_grid(ilevel))(ix, iy, iz - 1) + 45.0 * (*u.get_grid(ilevel))(ix, iy, iz + 1) - 9.0 * (*u.get_grid(ilevel))(ix, iy, iz + 2) + (*u.get_grid(ilevel))(ix, iy, iz + 3)) * h;
}
LOGUSER("Done computing a 6th order finite difference gradient.");
music::ulog.Print("Done computing a 6th order finite difference gradient.");
}
/**************************************************************************************/
@ -423,13 +423,13 @@ void multigrid_poisson_plugin::implementation::gradient_add_O6(int dir, grid_hie
double fft_poisson_plugin::solve(grid_hierarchy &f, grid_hierarchy &u)
{
LOGUSER("Entering k-space Poisson solver...");
music::ulog.Print("Entering k-space Poisson solver...");
unsigned verbosity = cf_.getValueSafe<unsigned>("setup", "verbosity", 2);
unsigned verbosity = cf_.get_value_safe<unsigned>("setup", "verbosity", 2);
if (f.levelmin() != f.levelmax())
{
LOGERR("Attempt to run k-space Poisson solver on non unigrid mesh.");
music::elog.Print("Attempt to run k-space Poisson solver on non unigrid mesh.");
throw std::runtime_error("fft_poisson_plugin::solve : k-space method can only be used in unigrid mode (levelmin=levelmax)");
}
@ -459,7 +459,7 @@ double fft_poisson_plugin::solve(grid_hierarchy &f, grid_hierarchy &u)
}
//... perform FFT and Poisson solve................................
LOGUSER("Performing forward transform.");
music::ulog.Print("Performing forward transform.");
#ifdef FFTW3
#ifdef SINGLE_PRECISION
@ -519,7 +519,7 @@ double fft_poisson_plugin::solve(grid_hierarchy &f, grid_hierarchy &u)
RE(cdata[0]) = 0.0;
IM(cdata[0]) = 0.0;
LOGUSER("Performing backward transform.");
music::ulog.Print("Performing backward transform.");
#ifdef FFTW3
#ifdef SINGLE_PRECISION
@ -592,14 +592,14 @@ double fft_poisson_plugin::solve(grid_hierarchy &f, grid_hierarchy &u)
}
}
LOGUSER("Done with k-space Poisson solver.");
music::ulog.Print("Done with k-space Poisson solver.");
return 0.0;
}
double fft_poisson_plugin::gradient(int dir, grid_hierarchy &u, grid_hierarchy &Du)
{
LOGUSER("Computing a gradient in k-space...\n");
music::ulog.Print("Computing a gradient in k-space...\n");
if (u.levelmin() != u.levelmax())
throw std::runtime_error("fft_poisson_plugin::gradient : k-space method can only be used in unigrid mode (levelmin=levelmax)");
@ -658,11 +658,11 @@ double fft_poisson_plugin::gradient(int dir, grid_hierarchy &u, grid_hierarchy &
double fac = -1.0 / (double)((size_t)nx * (size_t)ny * (size_t)nz);
double kfac = 2.0 * M_PI;
bool do_glass = cf_.getValueSafe<bool>("output", "glass", false);
bool deconvolve_cic = do_glass | cf_.getValueSafe<bool>("output", "glass_cicdeconvolve", false);
bool do_glass = cf_.get_value_safe<bool>("output", "glass", false);
bool deconvolve_cic = do_glass | cf_.get_value_safe<bool>("output", "glass_cicdeconvolve", false);
if (deconvolve_cic)
LOGINFO("CIC deconvolution is enabled for kernel!");
music::ilog.Print("CIC deconvolution is enabled for kernel!");
#pragma omp parallel for
for (int i = 0; i < nx; ++i)
@ -775,7 +775,7 @@ double fft_poisson_plugin::gradient(int dir, grid_hierarchy &u, grid_hierarchy &
delete[] data;
LOGUSER("Done with k-space gradient.\n");
music::ulog.Print("Done with k-space gradient.\n");
return 0.0;
}
@ -903,7 +903,7 @@ void do_poisson_hybrid(fftw_real *data, int idir, int nxp, int nyp, int nzp, boo
fftw_complex *cdata = reinterpret_cast<fftw_complex *>(data);
if (deconvolve_cic)
LOGINFO("CIC deconvolution step is enabled.");
music::ilog.Print("CIC deconvolution step is enabled.");
#ifdef FFTW3
#ifdef SINGLE_PRECISION
@ -1012,7 +1012,7 @@ void poisson_hybrid(T &f, int idir, int order, bool periodic, bool deconvolve_ci
int xo = 0, yo = 0, zo = 0;
int nmax = std::max(nx, std::max(ny, nz));
LOGUSER("Entering hybrid Poisson solver...");
music::ulog.Print("Entering hybrid Poisson solver...");
const int boundary = 32;
@ -1074,11 +1074,11 @@ void poisson_hybrid(T &f, int idir, int order, bool periodic, bool deconvolve_ci
break;
default:
std::cerr << " - ERROR: invalid operator order specified in deconvolution.";
LOGERR("Invalid operator order specified in deconvolution.");
music::elog.Print("Invalid operator order specified in deconvolution.");
break;
}
LOGUSER("Copying hybrid correction factor...");
music::ulog.Print("Copying hybrid correction factor...");
#pragma omp parallel for
for (int i = 0; i < nx; ++i)
@ -1091,7 +1091,7 @@ void poisson_hybrid(T &f, int idir, int order, bool periodic, bool deconvolve_ci
delete[] data;
LOGUSER("Done with hybrid Poisson solve.");
music::ulog.Print("Done with hybrid Poisson solve.");
}
/**************************************************************************************/

6
src/random.cc
View file

@ -32,21 +32,21 @@ void print_RNG_plugins()
RNG_plugin *select_RNG_plugin(config_file &cf)
{
std::string rngname = cf.getValueSafe<std::string>("random", "generator", "MUSIC");
std::string rngname = cf.get_value_safe<std::string>("random", "generator", "MUSIC");
RNG_plugin_creator *the_RNG_plugin_creator = get_RNG_plugin_map()[rngname];
if (!the_RNG_plugin_creator)
{
std::cerr << " - Error: random number generator plug-in \'" << rngname << "\' not found." << std::endl;
LOGERR("Invalid/Unregistered random number generator plug-in encountered : %s", rngname.c_str());
music::elog.Print("Invalid/Unregistered random number generator plug-in encountered : %s", rngname.c_str());
print_RNG_plugins();
throw std::runtime_error("Unknown random number generator plug-in");
}
else
{
std::cout << " - Selecting random number generator plug-in \'" << rngname << "\'..." << std::endl;
LOGUSER("Selecting random number generator plug-in : %s", rngname.c_str());
music::ulog.Print("Selecting random number generator plug-in : %s", rngname.c_str());
}
RNG_plugin *the_RNG_plugin = the_RNG_plugin_creator->create(cf);

4
src/random.hh
View file

@ -98,8 +98,8 @@ public:
random_number_generator(config_file &cf, transfer_function *ptf = NULL)
: pcf_(&cf) //, prefh_( &refh )
{
levelmin_ = pcf_->getValue<int>("setup", "levelmin");
levelmax_ = pcf_->getValue<int>("setup", "levelmax");
levelmin_ = pcf_->get_value<int>("setup", "levelmin");
levelmax_ = pcf_->get_value<int>("setup", "levelmax");
generator_ = select_RNG_plugin(cf);
}

46
src/region_generator.cc
View file

@ -24,21 +24,21 @@ void print_region_generator_plugins()
region_generator_plugin *select_region_generator_plugin(config_file &cf)
{
std::string rgname = cf.getValueSafe<std::string>("setup", "region", "box");
std::string rgname = cf.get_value_safe<std::string>("setup", "region", "box");
region_generator_plugin_creator *the_region_generator_plugin_creator = get_region_generator_plugin_map()[rgname];
if (!the_region_generator_plugin_creator)
{
std::cerr << " - Error: region generator plug-in \'" << rgname << "\' not found." << std::endl;
LOGERR("Invalid/Unregistered region generator plug-in encountered : %s", rgname.c_str());
music::elog.Print("Invalid/Unregistered region generator plug-in encountered : %s", rgname.c_str());
print_region_generator_plugins();
throw std::runtime_error("Unknown region generator plug-in");
}
else
{
std::cout << " - Selecting region generator plug-in \'" << rgname << "\'..." << std::endl;
LOGUSER("Selecting region generator plug-in : %s", rgname.c_str());
music::ulog.Print("Selecting region generator plug-in : %s", rgname.c_str());
}
region_generator_plugin *the_region_generator_plugin = the_region_generator_plugin_creator->create(cf);
@ -70,43 +70,43 @@ public:
region_box_plugin(config_file &cf)
: region_generator_plugin(cf)
{
levelmin_ = pcf_->getValue<unsigned>("setup", "levelmin");
levelmax_ = pcf_->getValue<unsigned>("setup", "levelmax");
levelmin_ = pcf_->get_value<unsigned>("setup", "levelmin");
levelmax_ = pcf_->get_value<unsigned>("setup", "levelmax");
if (levelmin_ != levelmax_)
{
padding_ = cf.getValue<int>("setup", "padding");
padding_ = cf.get_value<int>("setup", "padding");
std::string temp;
if (!pcf_->containsKey("setup", "ref_offset") && !pcf_->containsKey("setup", "ref_center"))
if (!pcf_->contains_key("setup", "ref_offset") && !pcf_->contains_key("setup", "ref_center"))
{
LOGERR("Found levelmin!=levelmax but neither ref_offset nor ref_center was specified.");
music::elog.Print("Found levelmin!=levelmax but neither ref_offset nor ref_center was specified.");
throw std::runtime_error("Found levelmin!=levelmax but neither ref_offset nor ref_center was specified.");
}
if (!pcf_->containsKey("setup", "ref_extent") && !pcf_->containsKey("setup", "ref_dims"))
if (!pcf_->contains_key("setup", "ref_extent") && !pcf_->contains_key("setup", "ref_dims"))
{
LOGERR("Found levelmin!=levelmax but neither ref_extent nor ref_dims was specified.");
music::elog.Print("Found levelmin!=levelmax but neither ref_extent nor ref_dims was specified.");
throw std::runtime_error("Found levelmin!=levelmax but neither ref_extent nor ref_dims was specified.");
}
if (pcf_->containsKey("setup", "ref_extent"))
if (pcf_->contains_key("setup", "ref_extent"))
{
temp = pcf_->getValue<std::string>("setup", "ref_extent");
temp = pcf_->get_value<std::string>("setup", "ref_extent");
std::remove_if(temp.begin(), temp.end(), isspace);
if (sscanf(temp.c_str(), "%lf,%lf,%lf", &lxref_[0], &lxref_[1], &lxref_[2]) != 3)
{
LOGERR("Error parsing triple for ref_extent");
music::elog.Print("Error parsing triple for ref_extent");
throw std::runtime_error("Error parsing triple for ref_extent");
}
bhave_nref_ = false;
}
else if (pcf_->containsKey("setup", "ref_dims"))
else if (pcf_->contains_key("setup", "ref_dims"))
{
temp = pcf_->getValue<std::string>("setup", "ref_dims");
temp = pcf_->get_value<std::string>("setup", "ref_dims");
std::remove_if(temp.begin(), temp.end(), isspace);
if (sscanf(temp.c_str(), "%lu,%lu,%lu", &lnref_[0], &lnref_[1], &lnref_[2]) != 3)
{
LOGERR("Error parsing triple for ref_dims");
music::elog.Print("Error parsing triple for ref_dims");
throw std::runtime_error("Error parsing triple for ref_dims");
}
bhave_nref_ = true;
@ -116,26 +116,26 @@ public:
lxref_[2] = lnref_[2] * 1.0 / (double)(1 << levelmax_);
}
if (pcf_->containsKey("setup", "ref_center"))
if (pcf_->contains_key("setup", "ref_center"))
{
temp = pcf_->getValue<std::string>("setup", "ref_center");
temp = pcf_->get_value<std::string>("setup", "ref_center");
std::remove_if(temp.begin(), temp.end(), isspace);
if (sscanf(temp.c_str(), "%lf,%lf,%lf", &xcref_[0], &xcref_[1], &xcref_[2]) != 3)
{
LOGERR("Error parsing triple for ref_center");
music::elog.Print("Error parsing triple for ref_center");
throw std::runtime_error("Error parsing triple for ref_center");
}
x0ref_[0] = fmod(xcref_[0] - 0.5 * lxref_[0] + 1.0, 1.0);
x0ref_[1] = fmod(xcref_[1] - 0.5 * lxref_[1] + 1.0, 1.0);
x0ref_[2] = fmod(xcref_[2] - 0.5 * lxref_[2] + 1.0, 1.0);
}
else if (pcf_->containsKey("setup", "ref_offset"))
else if (pcf_->contains_key("setup", "ref_offset"))
{
temp = pcf_->getValue<std::string>("setup", "ref_offset");
temp = pcf_->get_value<std::string>("setup", "ref_offset");
std::remove_if(temp.begin(), temp.end(), isspace);
if (sscanf(temp.c_str(), "%lf,%lf,%lf", &x0ref_[0], &x0ref_[1], &x0ref_[2]) != 3)
{
LOGERR("Error parsing triple for ref_offset");
music::elog.Print("Error parsing triple for ref_offset");
throw std::runtime_error("Error parsing triple for ref_offset");
}
xcref_[0] = fmod(x0ref_[0] + 0.5 * lxref_[0], 1.0);
@ -146,7 +146,7 @@ public:
// conditions should be added here
{
do_extra_padding_ = false;
std::string output_plugin = cf.getValue<std::string>("output", "format");
std::string output_plugin = cf.get_value<std::string>("output", "format");
if (output_plugin == std::string("grafic2"))
do_extra_padding_ = true;
padding_fine_ = 0.0;

6
src/region_generator.hh
View file

@ -2,7 +2,7 @@
#define __REGION_GENERATOR_HH
#include <vector>
#include "config_file.hh"
#include <config_file.hh>
#include <array>
using vec3_t = std::array<double,3>;
@ -22,8 +22,8 @@ public:
region_generator_plugin( config_file& cf )
: pcf_( &cf )
{
levelmin_ = cf.getValue<int>("setup","levelmin");
levelmax_ = cf.getValue<int>("setup","levelmax");
levelmin_ = cf.get_value<int>("setup","levelmin");
levelmax_ = cf.get_value<int>("setup","levelmax");
}
//! destructor

6
src/transfer_function.cc
View file

@ -37,7 +37,7 @@ void print_transfer_function_plugins()
transfer_function_plugin *select_transfer_function_plugin( config_file& cf )
{
std::string tfname = cf.getValue<std::string>( "cosmology", "transfer" );
std::string tfname = cf.get_value<std::string>( "cosmology", "transfer" );
transfer_function_plugin_creator *the_transfer_function_plugin_creator
= get_transfer_function_plugin_map()[ tfname ];
@ -45,14 +45,14 @@ transfer_function_plugin *select_transfer_function_plugin( config_file& cf )
if( !the_transfer_function_plugin_creator )
{
std::cerr << " - Error: transfer function plug-in \'" << tfname << "\' not found." << std::endl;
LOGERR("Invalid/Unregistered transfer function plug-in encountered : %s",tfname.c_str() );
music::elog.Print("Invalid/Unregistered transfer function plug-in encountered : %s",tfname.c_str() );
print_transfer_function_plugins();
throw std::runtime_error("Unknown transfer function plug-in");
}else
{
std::cout << " - Selecting transfer function plug-in \'" << tfname << "\'..." << std::endl;
LOGUSER("Selecting transfer function plug-in : %s",tfname.c_str() );
music::ulog.Print("Selecting transfer function plug-in : %s",tfname.c_str() );
}
transfer_function_plugin *the_transfer_function_plugin

14
src/transfer_function.hh
View file

@ -53,14 +53,14 @@ public:
: pcf_( &cf ), tf_distinct_(false), tf_withvel_(false), tf_withtotal0_(false), tf_velunits_(false)
{
real_t zstart;
zstart = pcf_->getValue<real_t>( "setup", "zstart" );
zstart = pcf_->get_value<real_t>( "setup", "zstart" );
cosmo_.astart = 1.0/(1.0+zstart);
cosmo_.Omega_b = pcf_->getValue<real_t>( "cosmology", "Omega_b" );
cosmo_.Omega_m = pcf_->getValue<real_t>( "cosmology", "Omega_m" );
cosmo_.Omega_DE = pcf_->getValue<real_t>( "cosmology", "Omega_L" );
cosmo_.H0 = pcf_->getValue<real_t>( "cosmology", "H0" );
cosmo_.sigma8 = pcf_->getValue<real_t>( "cosmology", "sigma_8" );
cosmo_.nspect = pcf_->getValue<real_t>( "cosmology", "nspec" );
cosmo_.Omega_b = pcf_->get_value<real_t>( "cosmology", "Omega_b" );
cosmo_.Omega_m = pcf_->get_value<real_t>( "cosmology", "Omega_m" );
cosmo_.Omega_DE = pcf_->get_value<real_t>( "cosmology", "Omega_L" );
cosmo_.H0 = pcf_->get_value<real_t>( "cosmology", "H0" );
cosmo_.sigma8 = pcf_->get_value<real_t>( "cosmology", "sigma_8" );
cosmo_.nspect = pcf_->get_value<real_t>( "cosmology", "nspec" );
}
//! destructor

14
tools/compute_ellipsoid.cc
View file

@ -22,9 +22,9 @@
#include <gsl/gsl_eigen.h>
#define LOGERR printf
#define LOGINFO printf
#define LOGUSER printf
#define music::elog.Print printf
#define music::ilog.Print printf
#define music::ulog.Print printf
/***** Math helper functions ******/
@ -295,7 +295,7 @@ protected:
}
if( count >= maxit )
LOGERR("No convergence in min_ellipsoid::compute: maximum number of iterations reached!");
music::elog.Print("No convergence in min_ellipsoid::compute: maximum number of iterations reached!");
delete[] unew;
}
@ -305,7 +305,7 @@ public:
: N( N_ ), axes_computed( false ), hold_point_data( true )
{
// --- initialize ---
LOGINFO("computing minimum bounding ellipsoid from %lld points",N);
music::ilog.Print("computing minimum bounding ellipsoid from %lld points",N);
Q = new float[4*N];
u = new float[N];
@ -470,7 +470,7 @@ public:
{
if( !axes_computed )
{
LOGUSER("computing ellipsoid axes.....");
music::ulog.Print("computing ellipsoid axes.....");
compute_axes();
}
float muold[3] = {mu[0],mu[1],mu[2]};
@ -493,7 +493,7 @@ public:
Inverse_3x3( A, Ainv );
//LOGUSER("computing ellipsoid axes.....");
//music::ulog.Print("computing ellipsoid axes.....");
compute_axes();
//print();