This commit is contained in:
Martin Reinecke 2019-01-10 13:59:36 +01:00
parent c9684732b8
commit ef2907f050
13 changed files with 91 additions and 183 deletions

View file

@ -9,12 +9,10 @@ src_sharp = \
pocketfft/pocketfft.h \
libsharp/sharp.c \
libsharp/sharp_almhelpers.c \
libsharp/sharp_announce.c \
libsharp/sharp_core.c \
libsharp/sharp_geomhelpers.c \
libsharp/sharp_legendre_roots.c \
libsharp/sharp_ylmgen_c.c \
libsharp/sharp_announce.h \
libsharp/sharp_internal.h \
libsharp/sharp_legendre_roots.h \
libsharp/sharp_vecsupport.h \

View file

@ -853,7 +853,8 @@ NOINLINE static void sharp_execute_job (sharp_job *job)
init_output (job);
int nchunks, chunksize;
get_chunk_info(job->ginfo->npairs,sharp_veclen()*sharp_max_nvec(),&nchunks,&chunksize);
get_chunk_info(job->ginfo->npairs,sharp_veclen()*sharp_max_nvec(job->spin),
&nchunks,&chunksize);
//FIXME: needs to be changed to "nm"
alloc_phase (job,mmax+1,chunksize);

View file

@ -29,8 +29,8 @@
* \author Martin Reinecke \author Dag Sverre Seljebotn
*/
#ifndef PLANCK_SHARP_LOWLEVEL_H
#define PLANCK_SHARP_LOWLEVEL_H
#ifndef PLANCK_SHARP_H
#define PLANCK_SHARP_H
#include <stddef.h>
@ -207,16 +207,13 @@ typedef enum { SHARP_DP = 1<<4,
\param type the type of SHT
\param spin the spin of the quantities to be transformed
\param alm contains pointers to the a_lm coefficients. If \a spin==0,
alm[0] points to the a_lm of the first SHT, alm[1] to those of the second
etc. If \a spin>0, alm[0] and alm[1] point to the a_lm of the first SHT,
alm[2] and alm[3] to those of the second, etc. The exact data type of \a alm
alm[0] points to the a_lm of the SHT. If \a spin>0, alm[0] and alm[1]
point to the two a_lm sets of the SHT. The exact data type of \a alm
depends on whether the SHARP_DP flag is set.
\param map contains pointers to the maps. If \a spin==0,
map[0] points to the map of the first SHT, map[1] to that of the second
etc. If \a spin>0, or \a type is SHARP_ALM2MAP_DERIV1, map[0] and map[1]
point to the maps of the first SHT, map[2] and map[3] to those of the
second, etc. The exact data type of \a map depends on whether the SHARP_DP
flag is set.
map[0] points to the map of the SHT. If \a spin>0, or \a type is
SHARP_ALM2MAP_DERIV1, map[0] and map[1] point to the two maps of the SHT.
The exact data type of \a map depends on whether the SHARP_DP flag is set.
\param geom_info A \c sharp_geom_info object compatible with the provided
\a map arrays.
\param alm_info A \c sharp_alm_info object compatible with the provided

View file

@ -1,98 +0,0 @@
/*
* This file is part of libc_utils.
*
* libc_utils is free software; you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation; either version 2 of the License, or
* (at your option) any later version.
*
* libc_utils 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 libc_utils; if not, write to the Free Software
* Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
*/
/*
* libc_utils is being developed at the Max-Planck-Institut fuer Astrophysik
* and financially supported by the Deutsches Zentrum fuer Luft- und Raumfahrt
* (DLR).
*/
/*! \file sharp_announce.c
* Banner for module startup
*
* Copyright (C) 2012 Max-Planck-Society
* \author Martin Reinecke
*/
#include <stdio.h>
#include <string.h>
#include <stdlib.h>
#ifdef _OPENMP
#include <omp.h>
#endif
#ifdef USE_MPI
#include <mpi.h>
#endif
#include "sharp_announce.h"
#include "sharp_internal.h"
static void OpenMP_status(void)
{
#ifndef _OPENMP
printf("OpenMP: not supported by this binary\n");
#else
int threads = omp_get_max_threads();
if (threads>1)
printf("OpenMP active: max. %d threads.\n",threads);
else
printf("OpenMP active, but running with 1 thread only.\n");
#endif
}
static void MPI_status(void)
{
#ifndef USE_MPI
printf("MPI: not supported by this binary\n");
#else
int tasks;
MPI_Comm_size(MPI_COMM_WORLD,&tasks);
if (tasks>1)
printf("MPI active with %d tasks.\n",tasks);
else
printf("MPI active, but running with 1 task only.\n");
#endif
}
static void vecmath_status(void)
{ printf("Supported vector length: %d\n",sharp_veclen()); }
void sharp_announce (const char *name)
{
size_t m, nlen=strlen(name);
printf("\n+-");
for (m=0; m<nlen; ++m) printf("-");
printf("-+\n");
printf("| %s |\n", name);
printf("+-");
for (m=0; m<nlen; ++m) printf("-");
printf("-+\n\n");
vecmath_status();
OpenMP_status();
MPI_status();
printf("\n");
}
void sharp_module_startup (const char *name, int argc, int argc_expected,
const char *argv_expected, int verbose)
{
if (verbose) sharp_announce (name);
if (argc==argc_expected) return;
if (verbose) fprintf(stderr, "Usage: %s %s\n", name, argv_expected);
exit(1);
}

View file

@ -1,39 +0,0 @@
/*
* This file is part of libc_utils.
*
* libc_utils is free software; you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation; either version 2 of the License, or
* (at your option) any later version.
*
* libc_utils 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 libc_utils; if not, write to the Free Software
* Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
*/
/*
* libc_utils is being developed at the Max-Planck-Institut fuer Astrophysik
* and financially supported by the Deutsches Zentrum fuer Luft- und Raumfahrt
* (DLR).
*/
/*! \file sharp_announce.h
* Banner for module startup
*
* Copyright (C) 2012 Max-Planck-Society
* \author Martin Reinecke
*/
#ifndef SHARP_ANNOUNCE_H
#define SHARP_ANNOUNCE_H
void sharp_announce (const char *name);
void sharp_module_startup (const char *name, int argc, int argc_expected,
const char *argv_expected, int verbose);
#endif

View file

@ -1096,7 +1096,7 @@ int sharp_veclen(void)
return VLEN;
}
int sharp_max_nvec(void)
int sharp_max_nvec(int spin)
{
return nv0;
return (spin==0) ? nv0 : nvx;
}

View file

@ -25,7 +25,7 @@
/*! \file sharp_cxx.h
* Spherical transform library
*
* Copyright (C) 2012-2017 Max-Planck-Society
* Copyright (C) 2012-2019 Max-Planck-Society
* \author Martin Reinecke
*/
@ -120,15 +120,13 @@ template<typename T> class sharp_cxxjob: public sharp_base
{
void *aptr=conv(alm), *mptr=conv(map);
int flags=cxxjobhelper__<T>::val | (add ? SHARP_ADD : 0);
sharp_execute (SHARP_ALM2MAP, 0, &aptr, &mptr, ginfo, ainfo, 1,
flags,0,0);
sharp_execute (SHARP_ALM2MAP, 0, &aptr, &mptr, ginfo, ainfo, flags, 0, 0);
}
void alm2map (const std::complex<T> *alm, T *map, bool add) const
{
void *aptr=conv(alm), *mptr=conv(map);
int flags=cxxjobhelper__<T>::val | (add ? SHARP_ADD : 0);
sharp_execute (SHARP_ALM2MAP, 0, &aptr, &mptr, ginfo, ainfo, 1,
flags,0,0);
sharp_execute (SHARP_ALM2MAP, 0, &aptr, &mptr, ginfo, ainfo, flags, 0, 0);
}
void alm2map_spin (const T *alm1, const T *alm2,
T *map1, T *map2, int spin, bool add) const
@ -137,7 +135,7 @@ template<typename T> class sharp_cxxjob: public sharp_base
aptr[0]=conv(alm1); aptr[1]=conv(alm2);
mptr[0]=conv(map1); mptr[1]=conv(map2);
int flags=cxxjobhelper__<T>::val | (add ? SHARP_ADD : 0);
sharp_execute (SHARP_ALM2MAP,spin,aptr,mptr,ginfo,ainfo,1,flags,0,0);
sharp_execute (SHARP_ALM2MAP,spin,aptr,mptr,ginfo,ainfo,flags, 0, 0);
}
void alm2map_spin (const std::complex<T> *alm1, const std::complex<T> *alm2,
T *map1, T *map2, int spin, bool add) const
@ -146,14 +144,14 @@ template<typename T> class sharp_cxxjob: public sharp_base
aptr[0]=conv(alm1); aptr[1]=conv(alm2);
mptr[0]=conv(map1); mptr[1]=conv(map2);
int flags=cxxjobhelper__<T>::val | (add ? SHARP_ADD : 0);
sharp_execute (SHARP_ALM2MAP,spin,aptr,mptr,ginfo,ainfo,1,flags,0,0);
sharp_execute (SHARP_ALM2MAP, spin, aptr, mptr, ginfo, ainfo, flags,0,0);
}
void alm2map_der1 (const T *alm, T *map1, T *map2, bool add) const
{
void *aptr=conv(alm), *mptr[2];
mptr[0]=conv(map1); mptr[1]=conv(map2);
int flags=cxxjobhelper__<T>::val | (add ? SHARP_ADD : 0);
sharp_execute (SHARP_ALM2MAP_DERIV1,1,&aptr,mptr,ginfo,ainfo,1,flags,0,0);
sharp_execute (SHARP_ALM2MAP_DERIV1,1,&aptr,mptr,ginfo,ainfo,flags,0,0);
}
void alm2map_der1 (const std::complex<T> *alm, T *map1, T *map2, bool add)
const
@ -161,7 +159,7 @@ template<typename T> class sharp_cxxjob: public sharp_base
void *aptr=conv(alm), *mptr[2];
mptr[0]=conv(map1); mptr[1]=conv(map2);
int flags=cxxjobhelper__<T>::val | (add ? SHARP_ADD : 0);
sharp_execute (SHARP_ALM2MAP_DERIV1,1,&aptr,mptr,ginfo,ainfo,1,flags,0,0);
sharp_execute (SHARP_ALM2MAP_DERIV1,1,&aptr,mptr,ginfo,ainfo,flags,0,0);
}
void alm2map_adjoint (const T *map, T *alm, bool add) const
{
@ -173,19 +171,19 @@ template<typename T> class sharp_cxxjob: public sharp_base
{
void *aptr=conv(alm), *mptr=conv(map);
int flags=cxxjobhelper__<T>::val | (add ? SHARP_ADD : 0);
sharp_execute (SHARP_Yt,0,&aptr,&mptr,ginfo,ainfo,1,flags,0,0);
sharp_execute (SHARP_Yt,0,&aptr,&mptr,ginfo,ainfo,flags,0,0);
}
void map2alm (const T *map, T *alm, bool add) const
{
void *aptr=conv(alm), *mptr=conv(map);
int flags=cxxjobhelper__<T>::val | (add ? SHARP_ADD : 0);
sharp_execute (SHARP_MAP2ALM,0,&aptr,&mptr,ginfo,ainfo,1,flags,0,0);
sharp_execute (SHARP_MAP2ALM,0,&aptr,&mptr,ginfo,ainfo,flags,0,0);
}
void map2alm (const T *map, std::complex<T> *alm, bool add) const
{
void *aptr=conv(alm), *mptr=conv(map);
int flags=cxxjobhelper__<T>::val | (add ? SHARP_ADD : 0);
sharp_execute (SHARP_MAP2ALM,0,&aptr,&mptr,ginfo,ainfo,1,flags,0,0);
sharp_execute (SHARP_MAP2ALM,0,&aptr,&mptr,ginfo,ainfo,flags,0,0);
}
void map2alm_spin (const T *map1, const T *map2, T *alm1, T *alm2,
int spin, bool add) const
@ -194,7 +192,7 @@ template<typename T> class sharp_cxxjob: public sharp_base
aptr[0]=conv(alm1); aptr[1]=conv(alm2);
mptr[0]=conv(map1); mptr[1]=conv(map2);
int flags=cxxjobhelper__<T>::val | (add ? SHARP_ADD : 0);
sharp_execute (SHARP_MAP2ALM,spin,aptr,mptr,ginfo,ainfo,1,flags,0,0);
sharp_execute (SHARP_MAP2ALM,spin,aptr,mptr,ginfo,ainfo,flags,0,0);
}
void map2alm_spin (const T *map1, const T *map2, std::complex<T> *alm1,
std::complex<T> *alm2, int spin, bool add) const
@ -203,7 +201,7 @@ template<typename T> class sharp_cxxjob: public sharp_base
aptr[0]=conv(alm1); aptr[1]=conv(alm2);
mptr[0]=conv(map1); mptr[1]=conv(map2);
int flags=cxxjobhelper__<T>::val | (add ? SHARP_ADD : 0);
sharp_execute (SHARP_MAP2ALM,spin,aptr,mptr,ginfo,ainfo,1,flags,0,0);
sharp_execute (SHARP_MAP2ALM,spin,aptr,mptr,ginfo,ainfo,flags,0,0);
}
};

View file

@ -25,7 +25,7 @@
/*! \file sharp_geomhelpers.h
* SHARP helper function for the creation of grid geometries
*
* Copyright (C) 2006-2013 Max-Planck-Society
* Copyright (C) 2006-2019 Max-Planck-Society
* \author Martin Reinecke
*/

View file

@ -25,7 +25,7 @@
/*! \file sharp_internal.h
* Internally used functionality for the spherical transform library.
*
* Copyright (C) 2006-2018 Max-Planck-Society
* Copyright (C) 2006-2019 Max-Planck-Society
* \author Martin Reinecke \author Dag Sverre Seljebotn
*/
@ -67,6 +67,6 @@ void inner_loop (sharp_job *job, const int *ispair,const double *cth,
const int *mlim);
int sharp_veclen(void);
int sharp_max_nvec(void);
int sharp_max_nvec(int spin);
#endif

View file

@ -24,7 +24,7 @@
/*! \file sharp_legendre_roots.h
*
* Copyright (C) 2006-2012 Max-Planck-Society
* Copyright (C) 2006-2019 Max-Planck-Society
* \author Martin Reinecke
*/

View file

@ -25,7 +25,7 @@
/*! \file sharp_mpi.c
* Functionality only needed for MPI-parallel transforms
*
* Copyright (C) 2012-2013 Max-Planck-Society
* Copyright (C) 2012-2019 Max-Planck-Society
* \author Martin Reinecke \author Dag Sverre Seljebotn
*/

View file

@ -25,7 +25,7 @@
/*! \file sharp_mpi.h
* Interface for the spherical transform library with MPI support.
*
* Copyright (C) 2011,2012 Max-Planck-Society
* Copyright (C) 2011-2019 Max-Planck-Society
* \author Martin Reinecke \author Dag Sverre Seljebotn
*/
@ -40,21 +40,18 @@ extern "C" {
#endif
/*! Performs an MPI parallel libsharp SHT job. The interface deliberately does
not use the C99 "complex" data type, in order to be callable from C.
not use the C99 "complex" data type, in order to be callable from C89 and C++.
\param comm the MPI communicator to be used for this SHT
\param type the type of SHT
\param spin the spin of the quantities to be transformed
\param alm contains pointers to the a_lm coefficients. If \a spin==0,
alm[0] points to the a_lm of the first SHT, alm[1] to those of the second
etc. If \a spin>0, alm[0] and alm[1] point to the a_lm of the first SHT,
alm[2] and alm[3] to those of the second, etc. The exact data type of \a alm
alm[0] points to the a_lm of the SHT. If \a spin>0, alm[0] and alm[1]
point to the two a_lm sets of the SHT. The exact data type of \a alm
depends on whether the SHARP_DP flag is set.
\param map contains pointers to the maps. If \a spin==0,
map[0] points to the map of the first SHT, map[1] to that of the second
etc. If \a spin>0, or \a type is SHARP_ALM2MAP_DERIV1, map[0] and map[1]
point to the maps of the first SHT, map[2] and map[3] to those of the
second, etc. The exact data type of \a map depends on whether the SHARP_DP
flag is set.
map[0] points to the map of the SHT. If \a spin>0, or \a type is
SHARP_ALM2MAP_DERIV1, map[0] and map[1] point to the two maps of the SHT.
The exact data type of \a map depends on whether the SHARP_DP flag is set.
\param geom_info A \c sharp_geom_info object compatible with the provided
\a map arrays. The total map geometry is the union of all \a geom_info
objects over the participating MPI tasks.

View file

@ -24,7 +24,7 @@
/* \file sharp_testsuite.c
*
* Copyright (C) 2012-2018 Max-Planck-Society
* Copyright (C) 2012-2019 Max-Planck-Society
* \author Martin Reinecke
*/
@ -42,9 +42,63 @@
#include "sharp_geomhelpers.h"
#include "sharp_almhelpers.h"
#include "c_utils.h"
#include "sharp_announce.h"
#include "memusage.h"
static void OpenMP_status(void)
{
#ifndef _OPENMP
printf("OpenMP: not supported by this binary\n");
#else
int threads = omp_get_max_threads();
if (threads>1)
printf("OpenMP active: max. %d threads.\n",threads);
else
printf("OpenMP active, but running with 1 thread only.\n");
#endif
}
static void MPI_status(void)
{
#ifndef USE_MPI
printf("MPI: not supported by this binary\n");
#else
int tasks;
MPI_Comm_size(MPI_COMM_WORLD,&tasks);
if (tasks>1)
printf("MPI active with %d tasks.\n",tasks);
else
printf("MPI active, but running with 1 task only.\n");
#endif
}
static void vecmath_status(void)
{ printf("Supported vector length: %d\n",sharp_veclen()); }
static void sharp_announce (const char *name)
{
size_t m, nlen=strlen(name);
printf("\n+-");
for (m=0; m<nlen; ++m) printf("-");
printf("-+\n");
printf("| %s |\n", name);
printf("+-");
for (m=0; m<nlen; ++m) printf("-");
printf("-+\n\n");
vecmath_status();
OpenMP_status();
MPI_status();
printf("\n");
}
static void sharp_module_startup (const char *name, int argc, int argc_expected,
const char *argv_expected, int verbose)
{
if (verbose) sharp_announce (name);
if (argc==argc_expected) return;
if (verbose) fprintf(stderr, "Usage: %s %s\n", name, argv_expected);
exit(1);
}
typedef complex double dcmplx;
int ntasks, mytask;