Aquila-Consortium/libsharp2
0
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

fix MPI mode

Browse source 
rename symbolic transform constants
This commit is contained in:
Martin Reinecke 2012-07-13 14:46:46 +02:00
parent 7928e13156
commit 4b0f2f1f54
9 changed files with 54 additions and 48 deletions
Download patch file Download diff file Expand all files Collapse all files

22
libsharp/sharp.c
View file

@ -360,7 +360,7 @@ static void fill_alm (const sharp_alm_info *ainfo, void *alm, dcmplx value,
static void init_output (sharp_job *job)
{
if (job->add_output) return;
if (job->type == MAP2ALM)
if (job->type == SHARP_MAP2ALM)
for (int i=0; i<job->ntrans*job->nalm; ++i)
fill_alm (job->ainfo,job->alm[i],0.,job->fde);
else
@ -374,10 +374,10 @@ static void alloc_phase (sharp_job *job, int nm, int ntheta)
static void dealloc_phase (sharp_job *job)
{ DEALLOC(job->phase); }
//FIXME: set phase to zero if not MAP2ALM?
//FIXME: set phase to zero if not SHARP_MAP2ALM?
static void map2phase (sharp_job *job, int mmax, int llim, int ulim)
{
if (job->type != MAP2ALM) return;
if (job->type != SHARP_MAP2ALM) return;
int pstride = 2*job->ntrans*job->nmaps;
#pragma omp parallel
{
@ -403,7 +403,7 @@ static void dealloc_almtmp (sharp_job *job)
static void alm2almtmp (sharp_job *job, int lmax, int mi)
{
if (job->type!=MAP2ALM)
if (job->type!=SHARP_MAP2ALM)
for (int l=job->ainfo->mval[mi]; l<=lmax; ++l)
{
ptrdiff_t aidx = sharp_alm_index(job->ainfo,l,mi);
@ -423,7 +423,7 @@ static void alm2almtmp (sharp_job *job, int lmax, int mi)
static void almtmp2alm (sharp_job *job, int lmax, int mi)
{
if (job->type != MAP2ALM) return;
if (job->type != SHARP_MAP2ALM) return;
for (int l=job->ainfo->mval[mi]; l<=lmax; ++l)
{
ptrdiff_t aidx = sharp_alm_index(job->ainfo,l,mi);
@ -439,7 +439,7 @@ static void almtmp2alm (sharp_job *job, int lmax, int mi)
static void phase2map (sharp_job *job, int mmax, int llim, int ulim)
{
if (job->type == MAP2ALM) return;
if (job->type == SHARP_MAP2ALM) return;
int pstride = 2*job->ntrans*job->nmaps;
#pragma omp parallel
{
@ -465,7 +465,7 @@ void sharp_execute_job (sharp_job *job)
int lmax = job->ainfo->lmax,
mmax=sharp_get_mmax(job->ainfo->mval, job->ainfo->nm);
job->norm_l = (job->type==ALM2MAP_DERIV1) ?
job->norm_l = (job->type==SHARP_ALM2MAP_DERIV1) ?
sharp_Ylmgen_get_d1norm (lmax) :
sharp_Ylmgen_get_norm (lmax, job->spin);
@ -546,14 +546,14 @@ static void sharp_build_job_common (sharp_job *job, sharp_jobtype type, int spin
const sharp_alm_info *alm_info, int ntrans)
{
UTIL_ASSERT((ntrans>0),"bad number of simultaneous transforms");
if (type==ALM2MAP_DERIV1) spin=1;
if (type==SHARP_ALM2MAP_DERIV1) spin=1;
UTIL_ASSERT((spin>=0)&&(spin<=30), "bad spin");
job->type = type;
job->spin = spin;
job->norm_l = NULL;
job->add_output = add_output;
job->nmaps = (type==ALM2MAP_DERIV1) ? 2 : ((spin>0) ? 2 : 1);
job->nalm = (type==ALM2MAP_DERIV1) ? 1 : ((spin>0) ? 2 : 1);
job->nmaps = (type==SHARP_ALM2MAP_DERIV1) ? 2 : ((spin>0) ? 2 : 1);
job->nalm = (type==SHARP_ALM2MAP_DERIV1) ? 1 : ((spin>0) ? 2 : 1);
job->ginfo = geom_info;
job->ainfo = alm_info;
job->nv = sharp_nv_oracle (type, spin, ntrans);
@ -653,7 +653,7 @@ int sharp_nv_oracle (sharp_jobtype type, int spin, int ntrans)
static int in_oracle=0;
if (in_oracle) return -20;
if (type==ALM2MAP_DERIV1) spin=1;
if (type==SHARP_ALM2MAP_DERIV1) spin=1;
UTIL_ASSERT((ntrans>0),"bad number of simultaneous transforms");
UTIL_ASSERT((spin>=0)&&(spin<=30), "bad spin");
ntrans=IMIN(ntrans,maxtr);

8
libsharp/sharp.h
View file

@ -140,9 +140,9 @@ void sharp_destroy_geom_info (sharp_geom_info *info);
/*! \{ */
/*! Enumeration of SHARP job types. */
typedef enum { MAP2ALM, /*!< analysis */
ALM2MAP, /*!< synthesis */
ALM2MAP_DERIV1 /*!< currently unused */
typedef enum { SHARP_MAP2ALM, /*!< analysis */
SHARP_ALM2MAP, /*!< synthesis */
SHARP_ALM2MAP_DERIV1 /*!< synthesis of first derivatives */
} sharp_jobtype;
typedef enum { FLOAT, DOUBLE } sharp_fde;
@ -171,7 +171,7 @@ typedef struct
/*! Initializes \a job with the appropriate parameters to perform the required
SHT.
\param type the type of SHT (currently ALM2MAP and MAP2ALM)
\param type the type of SHT
\param spin the spin of the quantities to be transformed
\param add_output if 0, the output arrays will be overwritten,
else the result will be added to the output arrays.

13
libsharp/sharp_acctest.c
View file

@ -112,7 +112,7 @@ static void check_sign_scale(void)
alm[i][j]=1.+_Complex_I;
sharp_job job;
sharpd_build_job(&job,ALM2MAP,0,0,&alm[0],&map[0],tinfo,alms,ntrans);
sharpd_build_job(&job,SHARP_ALM2MAP,0,0,&alm[0],&map[0],tinfo,alms,ntrans);
sharp_execute_job(&job);
for (int it=0; it<ntrans; ++it)
{
@ -123,7 +123,7 @@ static void check_sign_scale(void)
UTIL_ASSERT(FAPPROX(map[it][npix-1],-1.234675107554816442e+01,1e-12),
"error");
}
sharpd_build_job(&job,ALM2MAP,1,0,&alm[0],&map[0],tinfo,alms,ntrans);
sharpd_build_job(&job,SHARP_ALM2MAP,1,0,&alm[0],&map[0],tinfo,alms,ntrans);
sharp_execute_job(&job);
for (int it=0; it<ntrans; ++it)
{
@ -141,7 +141,7 @@ static void check_sign_scale(void)
"error");
}
sharpd_build_job(&job,ALM2MAP,2,0,&alm[0],&map[0],tinfo,alms,ntrans);
sharpd_build_job(&job,SHARP_ALM2MAP,2,0,&alm[0],&map[0],tinfo,alms,ntrans);
sharp_execute_job(&job);
for (int it=0; it<ntrans; ++it)
{
@ -159,7 +159,8 @@ static void check_sign_scale(void)
"error");
}
sharpd_build_job(&job,ALM2MAP_DERIV1,1,0,&alm[0],&map[0],tinfo,alms,ntrans);
sharpd_build_job(&job,SHARP_ALM2MAP_DERIV1,1,0,&alm[0],&map[0],tinfo,alms,
ntrans);
sharp_execute_job(&job);
for (int it=0; it<ntrans; ++it)
{
@ -207,11 +208,11 @@ static void check_accuracy (sharp_geom_info *tinfo, ptrdiff_t lmax,
ALLOC2D(alm2,dcmplx,ncomp,nalms);
sharp_job job;
sharpd_build_job(&job,ALM2MAP,spin,0,&alm[0],&map[0],tinfo,alms,ntrans);
sharpd_build_job(&job,SHARP_ALM2MAP,spin,0,&alm[0],&map[0],tinfo,alms,ntrans);
job.nv=nv;
sharp_execute_job(&job);
sharpd_build_job(&job,MAP2ALM,spin,0,&alm2[0],&map[0],tinfo,alms,ntrans);
sharpd_build_job(&job,SHARP_MAP2ALM,spin,0,&alm2[0],&map[0],tinfo,alms,ntrans);
job.nv=nv;
sharp_execute_job(&job);
measure_errors(alm,alm2,nalms,ncomp);

6
libsharp/sharp_bench.c
View file

@ -111,9 +111,9 @@ int main(void)
for (int spin=0; spin<=2; spin+=2)
{
fprintf(fp,"{");
for (sharp_jobtype type=MAP2ALM; type<=ALM2MAP_DERIV1; ++type)
for (sharp_jobtype type=SHARP_MAP2ALM; type<=SHARP_ALM2MAP_DERIV1; ++type)
{
if ((type==ALM2MAP_DERIV1) && (spin==0))
if ((type==SHARP_ALM2MAP_DERIV1) && (spin==0))
fprintf(fp,"-1");
else
{
@ -135,7 +135,7 @@ int main(void)
DEALLOC(time);
fprintf(fp,"%d",nvbest);
}
if (type!=ALM2MAP_DERIV1) fprintf(fp,",");
if (type!=SHARP_ALM2MAP_DERIV1) fprintf(fp,",");
}
fprintf(fp,(spin==0)?"},":"}");
printf("\n");

8
libsharp/sharp_core_inc2.c
View file

@ -641,8 +641,8 @@ static void Z(inner_loop) (sharp_job *job, const int *ispair,
switch (job->type)
{
case ALM2MAP:
case ALM2MAP_DERIV1:
case SHARP_ALM2MAP:
case SHARP_ALM2MAP_DERIV1:
{
if (job->spin==0)
{
@ -700,7 +700,7 @@ static void Z(inner_loop) (sharp_job *job, const int *ispair,
itot=idx[itot];
cth.s[i]=cth_[itot];
}
(job->type==ALM2MAP) ?
(job->type==SHARP_ALM2MAP) ?
Z(calc_alm2map_spin ) (cth.b,gen,job,&p1.b,&p2.b,&done) :
Z(calc_alm2map_deriv1) (cth.b,gen,job,&p1.b,&p2.b,&done);
}
@ -736,7 +736,7 @@ static void Z(inner_loop) (sharp_job *job, const int *ispair,
}
break;
}
case MAP2ALM:
case SHARP_MAP2ALM:
{
if (job->spin==0)
{

8
libsharp/sharp_core_inc3.c
View file

@ -629,8 +629,8 @@ static void Y(inner_loop) (sharp_job *job, const int *ispair,
switch (job->type)
{
case ALM2MAP:
case ALM2MAP_DERIV1:
case SHARP_ALM2MAP:
case SHARP_ALM2MAP_DERIV1:
{
if (job->spin==0)
{
@ -688,7 +688,7 @@ static void Y(inner_loop) (sharp_job *job, const int *ispair,
itot=idx[itot];
cth.s[i]=cth_[itot];
}
(job->type==ALM2MAP) ?
(job->type==SHARP_ALM2MAP) ?
Y(calc_alm2map_spin )
(cth.b,gen,job,&p1[0].b,&p2[0].b,njobs,&done) :
Y(calc_alm2map_deriv1)
@ -726,7 +726,7 @@ static void Y(inner_loop) (sharp_job *job, const int *ispair,
}
break;
}
case MAP2ALM:
case SHARP_MAP2ALM:
{
if (job->spin==0)
{

14
libsharp/sharp_mpi.c
View file

@ -182,20 +182,20 @@ static void sharp_communicate_map2alm (const sharp_mpi_info *minfo, dcmplx **ph)
static void alloc_phase_mpi (sharp_job *job, int nm, int ntheta,
int nmfull, int nthetafull)
{
ptrdiff_t phase_size = (job->type==MAP2ALM) ?
ptrdiff_t phase_size = (job->type==SHARP_MAP2ALM) ?
(ptrdiff_t)(nmfull)*ntheta : (ptrdiff_t)(nm)*nthetafull;
job->phase=RALLOC(dcmplx,2*job->ntrans*job->nmaps*phase_size);
}
static void alm2map_comm (sharp_job *job, const sharp_mpi_info *minfo)
{
if (job->type != MAP2ALM)
if (job->type != SHARP_MAP2ALM)
sharp_communicate_alm2map (minfo,&job->phase);
}
static void map2alm_comm (sharp_job *job, const sharp_mpi_info *minfo)
{
if (job->type == MAP2ALM)
if (job->type == SHARP_MAP2ALM)
sharp_communicate_map2alm (minfo,&job->phase);
}
@ -209,7 +209,7 @@ void sharp_execute_job_mpi (sharp_job *job, MPI_Comm comm)
int lmax = job->ainfo->lmax;
job->norm_l = Ylmgen_get_norm (lmax, job->spin);
job->norm_l = sharp_Ylmgen_get_norm (lmax, job->spin);
sharp_mpi_info minfo;
sharp_make_mpi_info(comm, job, &minfo);
@ -244,8 +244,8 @@ void sharp_execute_job_mpi (sharp_job *job, MPI_Comm comm)
#pragma omp parallel
{
sharp_job ljob = *job;
Ylmgen_C generator;
Ylmgen_init (&generator,lmax,minfo.mmax,ljob.spin);
sharp_Ylmgen_C generator;
sharp_Ylmgen_init (&generator,lmax,minfo.mmax,ljob.spin);
alloc_almtmp(&ljob,lmax);
#pragma omp for schedule(dynamic,1)
@ -262,7 +262,7 @@ void sharp_execute_job_mpi (sharp_job *job, MPI_Comm comm)
almtmp2alm (&ljob, lmax, mi);
}
Ylmgen_destroy(&generator);
sharp_Ylmgen_destroy(&generator);
dealloc_almtmp(&ljob);
#pragma omp critical

10
libsharp/sharp_test.c
View file

@ -117,7 +117,7 @@ static void map2alm_iter (sharp_geom_info *tinfo, double **map,
sharp_make_triangular_alm_info(lmax,mmax,1,&alms);
sharp_job job;
sharpd_build_job(&job,MAP2ALM,spin,0,&alm[0],&map[0],tinfo,alms,ntrans);
sharpd_build_job(&job,SHARP_MAP2ALM,spin,0,&alm[0],&map[0],tinfo,alms,ntrans);
sharp_execute_job(&job);
printf("wall time for map2alm: %fs\n",job.time);
printf("Performance: %fGFLOPs/s\n",1e-9*job.opcnt/job.time);
@ -128,7 +128,8 @@ static void map2alm_iter (sharp_geom_info *tinfo, double **map,
double **map2;
ALLOC2D(map2,double,ncomp,npix);
printf ("\niteration %i:\n", iter+1);
sharpd_build_job(&job,ALM2MAP,spin,0,&alm[0],&map2[0],tinfo,alms,ntrans);
sharpd_build_job(&job,SHARP_ALM2MAP,spin,0,&alm[0],&map2[0],tinfo,alms,
ntrans);
sharp_execute_job(&job);
printf("wall time for alm2map: %fs\n",job.time);
printf("Performance: %fGFLOPs/s\n",1e-9*job.opcnt/job.time);
@ -136,7 +137,8 @@ static void map2alm_iter (sharp_geom_info *tinfo, double **map,
for (ptrdiff_t m=0; m<npix; ++m)
map2[i][m] = map[i][m]-map2[i][m];
sharpd_build_job(&job,MAP2ALM,spin,1,&alm[0],&map2[0],tinfo,alms,ntrans);
sharpd_build_job(&job,SHARP_MAP2ALM,spin,1,&alm[0],&map2[0],tinfo,alms,
ntrans);
sharp_execute_job(&job);
printf("wall time for map2alm: %fs\n",job.time);
printf("Performance: %fGFLOPs/s\n",1e-9*job.opcnt/job.time);
@ -170,7 +172,7 @@ static void check_accuracy (sharp_geom_info *tinfo, ptrdiff_t lmax,
sharp_job job;
printf ("\niteration 0:\n");
sharpd_build_job(&job,ALM2MAP,spin,0,&alm[0],&map[0],tinfo,alms,ntrans);
sharpd_build_job(&job,SHARP_ALM2MAP,spin,0,&alm[0],&map[0],tinfo,alms,ntrans);
sharp_execute_job(&job);
printf("wall time for alm2map: %fs\n",job.time);
printf("Performance: %fGFLOPs/s\n",1e-9*job.opcnt/job.time);

13
libsharp/sharp_test_mpi.c
View file

@ -59,6 +59,7 @@
#include "sharp_almhelpers.h"
#include "c_utils.h"
#include "walltime_c.h"
#include "sharp_announce.h"
#include "sharp_core.h"
typedef complex double dcmplx;
@ -197,7 +198,7 @@ static void map2alm_iter (sharp_geom_info *tinfo, double **map,
reduce_alm_info(alms);
sharp_job job;
sharpd_build_job(&job,MAP2ALM,spin,0,&alm[0],&map[0],tinfo,alms,ntrans);
sharpd_build_job(&job,SHARP_MAP2ALM,spin,0,&alm[0],&map[0],tinfo,alms,ntrans);
sharp_execute_job_mpi(&job,MPI_COMM_WORLD);
unsigned long long opcnt=totalops(job.opcnt);
double timer=maxTime(job.time);
@ -210,7 +211,8 @@ static void map2alm_iter (sharp_geom_info *tinfo, double **map,
double **map2;
ALLOC2D(map2,double,ncomp,npix);
if (mytask==0) printf ("\niteration %i:\n", iter+1);
sharpd_build_job(&job,ALM2MAP,spin,0,&alm[0],&map2[0],tinfo,alms,ntrans);
sharpd_build_job(&job,SHARP_ALM2MAP,spin,0,&alm[0],&map2[0],tinfo,alms,
ntrans);
sharp_execute_job_mpi(&job,MPI_COMM_WORLD);
opcnt=totalops(job.opcnt);
timer=maxTime(job.time);
@ -220,7 +222,8 @@ static void map2alm_iter (sharp_geom_info *tinfo, double **map,
for (ptrdiff_t m=0; m<npix; ++m)
map2[i][m] = map[i][m]-map2[i][m];
sharpd_build_job(&job,MAP2ALM,spin,1,&alm[0],&map2[0],tinfo,alms,ntrans);
sharpd_build_job(&job,SHARP_MAP2ALM,spin,1,&alm[0],&map2[0],tinfo,alms,
ntrans);
sharp_execute_job_mpi(&job,MPI_COMM_WORLD);
opcnt=totalops(job.opcnt);
timer=maxTime(job.time);
@ -258,7 +261,7 @@ static void check_accuracy (sharp_geom_info *tinfo, ptrdiff_t lmax,
if (mytask==0) printf ("\niteration 0:\n");
sharp_job job;
sharpd_build_job(&job,ALM2MAP,spin,0,&alm[0],&map[0],tinfo,alms,ntrans);
sharpd_build_job(&job,SHARP_ALM2MAP,spin,0,&alm[0],&map[0],tinfo,alms,ntrans);
sharp_execute_job_mpi(&job,MPI_COMM_WORLD);
unsigned long long opcnt=totalops(job.opcnt);
double timer=maxTime(job.time);
@ -280,7 +283,7 @@ int main(int argc, char **argv)
MPI_Comm_size(MPI_COMM_WORLD,&ntasks);
MPI_Comm_rank(MPI_COMM_WORLD,&mytask);
module_startup_c("sharp_test_mpi",argc,7,
sharp_module_startup("sharp_test_mpi",argc,7,
"<healpix|ecp|gauss> <lmax> <nside|nphi> <niter> <spin> <ntrans>",
mytask==0);
int lmax=atoi(argv[2]);