diff --git a/external/cosmotool/build_tools/gather_sources.py b/external/cosmotool/build_tools/gather_sources.py new file mode 100644 index 0000000..ca8c3d1 --- /dev/null +++ b/external/cosmotool/build_tools/gather_sources.py @@ -0,0 +1,161 @@ +#+ +# This is CosmoTool (./build_tools/gather_sources.py) -- Copyright (C) Guilhem Lavaux (2007-2013) +# +# guilhem.lavaux@gmail.com +# +# This software is a computer program whose purpose is to provide a toolbox for cosmological +# data analysis (e.g. filters, generalized Fourier transforms, power spectra, ...) +# +# This software is governed by the CeCILL license under French law and +# abiding by the rules of distribution of free software. You can use, +# modify and/ or redistribute the software under the terms of the CeCILL +# license as circulated by CEA, CNRS and INRIA at the following URL +# "http://www.cecill.info". +# +# As a counterpart to the access to the source code and rights to copy, +# modify and redistribute granted by the license, users are provided only +# with a limited warranty and the software's author, the holder of the +# economic rights, and the successive licensors have only limited +# liability. +# +# In this respect, the user's attention is drawn to the risks associated +# with loading, using, modifying and/or developing or reproducing the +# software by the user in light of its specific status of free software, +# that may mean that it is complicated to manipulate, and that also +# therefore means that it is reserved for developers and experienced +# professionals having in-depth computer knowledge. Users are therefore +# encouraged to load and test the software's suitability as regards their +# requirements in conditions enabling the security of their systems and/or +# data to be ensured and, more generally, to use and operate it in the +# same conditions as regards security. +# +# The fact that you are presently reading this means that you have had +# knowledge of the CeCILL license and that you accept its terms. +#+ +import shutil +import tempfile +import re +from git import Repo,Tree,Blob + +def apply_license(license, relimit, filename): + header = re.sub(r'@FILENAME@', filename, license) + + f = file(filename) + lines = f.read() + f.close() + + lines = re.sub(relimit, '', lines) + lines = header + lines + + with tempfile.NamedTemporaryFile(delete=False) as tmp_sources: + tmp_sources.write(lines) + + shutil.move(tmp_sources.name, filename) + + +def apply_python_license(filename): + license="""#+ +# This is CosmoTool (@FILENAME@) -- Copyright (C) Guilhem Lavaux (2007-2013) +# +# guilhem.lavaux@gmail.com +# +# This software is a computer program whose purpose is to provide a toolbox for cosmological +# data analysis (e.g. filters, generalized Fourier transforms, power spectra, ...) +# +# This software is governed by the CeCILL license under French law and +# abiding by the rules of distribution of free software. You can use, +# modify and/ or redistribute the software under the terms of the CeCILL +# license as circulated by CEA, CNRS and INRIA at the following URL +# "http://www.cecill.info". +# +# As a counterpart to the access to the source code and rights to copy, +# modify and redistribute granted by the license, users are provided only +# with a limited warranty and the software's author, the holder of the +# economic rights, and the successive licensors have only limited +# liability. +# +# In this respect, the user's attention is drawn to the risks associated +# with loading, using, modifying and/or developing or reproducing the +# software by the user in light of its specific status of free software, +# that may mean that it is complicated to manipulate, and that also +# therefore means that it is reserved for developers and experienced +# professionals having in-depth computer knowledge. Users are therefore +# encouraged to load and test the software's suitability as regards their +# requirements in conditions enabling the security of their systems and/or +# data to be ensured and, more generally, to use and operate it in the +# same conditions as regards security. +# +# The fact that you are presently reading this means that you have had +# knowledge of the CeCILL license and that you accept its terms. +#+ +""" + + print("Shell/Python file: %s" % filename) + relimit=r'^#\+\n(#.*\n)*#\+\n' + apply_license(license, relimit, filename) + + +def apply_cpp_license(filename): + license="""/*+ +This is CosmoTool (@FILENAME@) -- Copyright (C) Guilhem Lavaux (2007-2013) + +guilhem.lavaux@gmail.com + +This software is a computer program whose purpose is to provide a toolbox for cosmological +data analysis (e.g. filters, generalized Fourier transforms, power spectra, ...) + +This software is governed by the CeCILL license under French law and +abiding by the rules of distribution of free software. You can use, +modify and/ or redistribute the software under the terms of the CeCILL +license as circulated by CEA, CNRS and INRIA at the following URL +"http://www.cecill.info". + +As a counterpart to the access to the source code and rights to copy, +modify and redistribute granted by the license, users are provided only +with a limited warranty and the software's author, the holder of the +economic rights, and the successive licensors have only limited +liability. + +In this respect, the user's attention is drawn to the risks associated +with loading, using, modifying and/or developing or reproducing the +software by the user in light of its specific status of free software, +that may mean that it is complicated to manipulate, and that also +therefore means that it is reserved for developers and experienced +professionals having in-depth computer knowledge. Users are therefore +encouraged to load and test the software's suitability as regards their +requirements in conditions enabling the security of their systems and/or +data to be ensured and, more generally, to use and operate it in the +same conditions as regards security. + +The fact that you are presently reading this means that you have had +knowledge of the CeCILL license and that you accept its terms. ++*/ +""" + relimit = r'^(?s)/\*\+.*\+\*/' + print("C++ file: %s" % filename) + apply_license(license, relimit, filename) + + +def analyze_tree(prefix, t): + for ename,entry in t.items(): + if ename == 'external' or ename == 'zobov': + continue + if type(entry) == Tree: + analyze_tree(prefix + "/" + ename, entry) + elif type(entry) == Blob: + if ename == './src/hdf5_flash.h' or ename == './src/h5_readFlash.cpp' or ename == './src/h5_readFlash.hpp': + continue + + if re.match(".*\.(sh|py|pyx)$",ename) != None: + fname=prefix+"/"+ename + apply_python_license(fname) + if re.match('.*\.(cpp|hpp|h)$', ename) != None: + fname=prefix+"/"+ename + apply_cpp_license(fname) + + +if __name__=="__main__": + repo = Repo(".") + assert repo.bare == False + t = repo.tree() + analyze_tree(".", t) diff --git a/external/cosmotool/sample/testkd3.cpp b/external/cosmotool/sample/testkd3.cpp new file mode 100644 index 0000000..50eba07 --- /dev/null +++ b/external/cosmotool/sample/testkd3.cpp @@ -0,0 +1,187 @@ +/*+ +This is CosmoTool (./sample/testkd3.cpp) -- Copyright (C) Guilhem Lavaux (2007-2013) + +guilhem.lavaux@gmail.com + +This software is a computer program whose purpose is to provide a toolbox for cosmological +data analysis (e.g. filters, generalized Fourier transforms, power spectra, ...) + +This software is governed by the CeCILL license under French law and +abiding by the rules of distribution of free software. You can use, +modify and/ or redistribute the software under the terms of the CeCILL +license as circulated by CEA, CNRS and INRIA at the following URL +"http://www.cecill.info". + +As a counterpart to the access to the source code and rights to copy, +modify and redistribute granted by the license, users are provided only +with a limited warranty and the software's author, the holder of the +economic rights, and the successive licensors have only limited +liability. + +In this respect, the user's attention is drawn to the risks associated +with loading, using, modifying and/or developing or reproducing the +software by the user in light of its specific status of free software, +that may mean that it is complicated to manipulate, and that also +therefore means that it is reserved for developers and experienced +professionals having in-depth computer knowledge. Users are therefore +encouraged to load and test the software's suitability as regards their +requirements in conditions enabling the security of their systems and/or +data to be ensured and, more generally, to use and operate it in the +same conditions as regards security. + +The fact that you are presently reading this means that you have had +knowledge of the CeCILL license and that you accept its terms. ++*/ + +#include +#include +#include +#include +#include +#define __KD_TREE_NUMNODES +#include "mykdtree.hpp" +#include "kdtree_splitters.hpp" +#include + +#define NTRY 30 +#define NGB 24 +#define ND 3 + +using namespace std; +using namespace CosmoTool; + +typedef KDTree > MyTree; +//typedef KDTree MyTree; +typedef KDCell MyCell; + +static double periodic(double a) +{ + while (a < -0.5) + a+=1; + while (a > 0.5) + a -= 1; + return a; +} + +MyCell *findNearest(MyTree::coords& xc, MyCell *cells, uint32_t Ncells) +{ + MyCell *near2 = 0; + double R2 = INFINITY; + for (int i = 0; i < Ncells; i++) + { + double d2 = 0; + for (int j = 0; j < ND; j++) + { + double delta = periodic(xc[j]-cells[i].coord[j]); + d2 += delta*delta; + } + if (d2 < R2) + { + near2 = &cells[i]; + R2 = d2; + } + } + return near2; +} + +int main() +{ + uint32_t Ncells = 3000000; + MyCell *cells = new MyCell[Ncells]; + + for (int i = 0; i < Ncells; i++) + { + cells[i].active = true; + for (int l = 0; l < ND; l++) + cells[i].coord[l] = drand48(); + } + + // Check timing + clock_t startTimer = clock(); + MyTree tree(cells, Ncells); + clock_t endTimer = clock(); + + tree.setPeriodic(true, 1.0); + + clock_t delta = endTimer-startTimer; + double myTime = delta*1.0/CLOCKS_PER_SEC * 1.0; + + cout << "KDTree build = " << myTime << " s" << endl; + + MyTree::coords *xc = new MyTree::coords[NTRY]; + + cout << "Generating seeds..." << endl; + for (int k = 0; k < NTRY; k++) + { + for (int l = 0; l < ND; l++) + xc[k][l] = drand48(); + } + + // Check consistency + cout << "Check consistency..." << endl; + MyCell **ngb = new MyCell *[NGB]; + double *distances = new double[NGB]; + + ofstream fngb("nearest.txt"); + for (int k = 0; k < NTRY; k++) { + cout << "Seed = " << xc[k][0] << " " << xc[k][1] << " " << xc[k][2] << endl; + tree.getNearestNeighbours(xc[k], NGB, ngb, distances); + int last = -1; + + for (uint32_t i = 0; i < NGB; i++) + { + if (ngb[i] == 0) + continue; + + last = i; + + double d2 = 0; + for (int l = 0; l < 3; l++) + d2 += ({double delta = periodic(xc[k][l] - ngb[i]->coord[l]); delta*delta;}); + fngb << ngb[i]->coord[0] << " " << ngb[i]->coord[1] << " " << ngb[i]->coord[2] << " " << sqrt(d2) << " " << sqrt(distances[i]) << endl; + } + fngb << endl << endl; + double farther_dist = distances[last]; + for (uint32_t i = 0; i < Ncells; i++) + { + bool found = false; + // If the points is not in the list, it means it is farther than the farthest point + for (int j =0; j < NGB; j++) + { + if (&cells[i] == ngb[j]) { + found = true; + break; + } + } + double dist_to_seed = 0; + for (int l = 0; l < 3; l++) + { + double delta = periodic(xc[k][l]-cells[i].coord[l]); + dist_to_seed += delta*delta; + } + double delta_epsilon = fabs(farther_dist/dist_to_seed-1); + if (!found) + { + if (dist_to_seed <= farther_dist && + delta_epsilon > 1000*sqrt(numeric_limits::epsilon())) + { + cout << "SHOULD BE IN (d=" << dist_to_seed << "): Failed iteration " << k << " particle " << i << " (farthest=" << sqrt(farther_dist) << ")" << endl; + cout << "delta_epsilon = " << delta_epsilon << endl; + abort(); + } + } + else + { + if (dist_to_seed > farther_dist && + delta_epsilon > 1000*sqrt(numeric_limits::epsilon())) + { + cout << "SHOULD BE OUT (d=" << sqrt(dist_to_seed) << "): Failed iteration " << k << " particle " << i << " (farthest=" << sqrt(farther_dist) << ")" << endl; + cout << "delta_epsilon = " << delta_epsilon << " epsilon = " << 100*sqrt(numeric_limits::epsilon()) << endl; + abort(); + } + } + } + } + + return 0; +} diff --git a/external/cosmotool/src/replicateGenerator.hpp b/external/cosmotool/src/replicateGenerator.hpp new file mode 100644 index 0000000..fd7035f --- /dev/null +++ b/external/cosmotool/src/replicateGenerator.hpp @@ -0,0 +1,100 @@ +/*+ +This is CosmoTool (./src/replicateGenerator.hpp) -- Copyright (C) Guilhem Lavaux (2007-2013) + +guilhem.lavaux@gmail.com + +This software is a computer program whose purpose is to provide a toolbox for cosmological +data analysis (e.g. filters, generalized Fourier transforms, power spectra, ...) + +This software is governed by the CeCILL license under French law and +abiding by the rules of distribution of free software. You can use, +modify and/ or redistribute the software under the terms of the CeCILL +license as circulated by CEA, CNRS and INRIA at the following URL +"http://www.cecill.info". + +As a counterpart to the access to the source code and rights to copy, +modify and redistribute granted by the license, users are provided only +with a limited warranty and the software's author, the holder of the +economic rights, and the successive licensors have only limited +liability. + +In this respect, the user's attention is drawn to the risks associated +with loading, using, modifying and/or developing or reproducing the +software by the user in light of its specific status of free software, +that may mean that it is complicated to manipulate, and that also +therefore means that it is reserved for developers and experienced +professionals having in-depth computer knowledge. Users are therefore +encouraged to load and test the software's suitability as regards their +requirements in conditions enabling the security of their systems and/or +data to be ensured and, more generally, to use and operate it in the +same conditions as regards security. + +The fact that you are presently reading this means that you have had +knowledge of the CeCILL license and that you accept its terms. ++*/ +#ifndef __REPLICATE_GENERATOR_HPP +#define __REPLICATE_GENERATOR_HPP + +#include +#include "algo.hpp" + +namespace CosmoTool +{ + + template + class ReplicateGenerator + { + public: + typedef Coord Coords[N]; + Coord replicate; + + ReplicateGenerator(const Coords& x, Coord shift) + { + face = 0; + replicate = shift; + numFaces = spower(3); + std::copy(x, x+N, x_base); + if (!next()) + abort(); + } + + bool next() + { + if (face == numFaces) + return false; + + face++; + + bool no_move = true; + int q_face = face; + for (int i = 0; i < N; i++) + { + int c_face; + c_face = q_face % 3; + q_face /= 3; + x_shifted[i] = x_base[i] + (c_face-1)*replicate; + no_move = no_move && (c_face == 1); + } + if (no_move) + return next(); + return true; + } + + const Coord *getPosition() + { + return x_shifted; + } + + void getPosition(Coords& x_out) + { + std::copy(x_shifted, x_shifted+N, x_out); + } + + private: + Coord x_shifted[N], x_base[N]; + long face, numFaces; + }; + +}; + +#endif