toeplitz.h

Go to the documentation of this file.

#ifndef         TOEPLITZ_H_
#define         TOEPLITZ_H_

#ifdef W_MPI
#include <mpi.h>
#endif

#ifdef W_OPENMP
#include <omp.h>
#endif

#include <fftw3.h>
#include <stdlib.h>
#include <stdio.h>
#include <math.h>
#include <string.h>

//=========================================================================
//Basic functions definition
#define max(a,b) (((a) > (b)) ? (a) : (b))
#define min(a,b) (((a) < (b)) ? (a) : (b))
 
//=========================================================================
//Fixed parameters



#ifndef MPI_USER_TAG
#define MPI_USER_TAG 123
#endif

//Define this parameter to use fftw multithreading
//This is not fully tested
#ifndef fftw_MULTITHREADING
#define fftw_MULTITHREADING 
#endif


#ifndef NFFT_DEFAULT
#define NFFT_DEFAULT 1 /*1*/
#endif



#ifndef FFTW_FLAG_AUTO
#define FFTW_FLAG_AUTO FFTW_ESTIMATE
#endif


//Parameters to define the computational strategy
#ifndef FLAG_STGY
#define FLAG_STGY

#define FLAG_BS 0  //0:auto  1:fixed  2:zero  3:3lambda  4:4lambda  5:formula2 
#define FLAG_NFFT 0  //0:auto  1:fixed  2:numthreads  3:fftwthreads 
#define FLAG_FFTW FFTW_FLAG_AUTO  //ESTIMATE, MEASURE, PATIENT, EXHAUSTIVE. Default is MEASURE
#define FLAG_NO_RSHP 0  //0:auto  1:yes  1:no
#define FLAG_NOFFT 0 //0:auto  1:yes  1:no
#define FLAG_BLOCKINGCOMM 0  //0:auto  1:noblocking  2:blocking 
#define FIXED_NFFT 0  //fixed init value for nfft
#define FIXED_BS 0   //fixed init value for blockside
#define FLAG_VERBOSE 0
#define FLAG_SKIP_BUILD_GAPPY_BLOCKS 0
#define FLAG_PARAM_DISTMIN_FIXED 0
#define FLAG_PRECOMPUTE_LVL 0   //0: no precompute  1: precompute plans  2: precomputes Toeplitz and plans

#endif

//=========================================================================
//Global parameters

extern int VERBOSE;
extern int PRINT_RANK;

//=========================================================================
//Strutures definition


typedef struct Block {
     int64_t idv;
     double *T_block;  //pointer of the Toeplitz data
     int lambda;
     int n;
/* For precomputed fftw
     int bs;
     int nfft;
     fftw_complex *T_fft;
     fftw_complex *V_fft;
     double *V_rfft;
     fftw_plan plan_f;
     fftw_plan plan_b;
     fftw_plan plan_f_T;
*/
} Block;


typedef struct Flag {
       int flag_bs;  //bs used formula
       int flag_nfft;   
       int flag_fftw;
       int flag_no_rshp;  //with or without
       int flag_nofft;
       int flag_blockingcomm;
       int fixed_nfft;  //init value for nfft
       int fixed_bs;   //long long int
       int flag_verbose;
       int flag_skip_build_gappy_blocks;
       int flag_param_distmin_fixed;
       int flag_precompute_lvl;
} Flag;

typedef struct Gap {
  int64_t *id0gap;
  int *lgap;
  int ngap;
} Gap;


typedef struct Tpltz {
     int64_t nrow;  //n total
     int m_cw;  //V column number in the linear row-wise order (vect row-wise order)
     int m_rw;  //V column number in the uniform row-wise order (matrix row-wise order)
     Block *tpltzblocks;
     int nb_blocks_loc;
     int nb_blocks_tot;
     int64_t idp;
     int local_V_size;
     Flag flag_stgy;
     MPI_Comm comm;
} Tpltz;


//=========================================================================
//Groups definition for documentation

//=========================================================================
// User routines definition (API)

//Wizard routines
int stbmmProd( Tpltz Nm1, double *V);


//Sequential routines (group 11)
int tpltz_init(int n, int lambda, int *nfft, int *blocksize, fftw_complex **T_fft, double *T, fftw_complex **V_fft, double **V_rfft, fftw_plan *plan_f, fftw_plan *plan_b, Flag flag_stgy);

int tpltz_cleanup(fftw_complex **T_fft, fftw_complex **V_fft, double **V_rfft,fftw_plan *plan_f, fftw_plan *plan_b);

int stmm_core(double **V, int n, int m, double *T, fftw_complex *T_fft, int blocksize, int lambda, fftw_complex *V_fft, double *V_rfft, int nfft, fftw_plan plan_f, fftw_plan plan_b, int flag_offset, int flag_nofft);

int stmm_main(double **V, int n, int m, int id0, int l, double *T, fftw_complex *T_fft, int lambda, fftw_complex *V_fft, double *V_rfft, fftw_plan plan_f, fftw_plan plan_b, int blocksize, int nfft, Flag flag_stgy);

int stmm(double **V, int n, int m, double *T, int lambda, Flag flag_stgy);

//int stbmm(double **V, int nrow, int m, int m_rowwise, Block *tpltzblocks, int nb_blocks_local, int nb_blocks_all, int idp, int local_V_size, Flag flag_stgy);

//int gstbmm(double **V, int nrow, int m, int m_rowwise, Block *tpltzblocks, int nb_blocks_local, int nb_blocks_all, int id0p, int local_V_size, int *id0gap, int *lgap, int ngap,Flag flag_stgy);

int stbmm(double **V, int nrow, int m_cw, int m_rw, Block *tpltzblocks, int nb_blocks, int64_t idp, int local_V_size, Flag flag_stgy);

int gstbmm(double **V, int nrow, int m_cw, int m_rw, Block *tpltzblocks, int nb_blocks, int64_t idp, int local_V_size, int64_t *id0gap, int *lgap, int ngap, Flag flag_stgy);


int reset_gaps(double **V, int id0,int local_V_size, int m, int nrow, int m_rowwise, int64_t *id0gap, int *lgap, int ngap);


//Mpi routines (group 12)
#ifdef W_MPI
int mpi_stmm(double **V, int n, int m, int id0, int l, double *T, int lambda, Flag flag_stgy, MPI_Comm comm);

int mpi_stbmm(double **V, int64_t nrow, int m, int m_rowwise, Block *tpltzblocks, int nb_blocks_local, int nb_blocks_all, int64_t idp, int local_V_size, Flag flag_stgy, MPI_Comm comm);

int mpi_gstbmm(double **V, int nrow, int m, int m_rowwise, Block *tpltzblocks, int nb_blocks_local, int nb_blocks_all, int id0p, int local_V_size, int64_t *id0gap, int *lgap, int ngap, Flag flag_stgy, MPI_Comm comm);

#endif


//=========================================================================
// User routines definition

//Low level routines (group 21)
int flag_stgy_init_auto(Flag *flag_stgy);

int flag_stgy_init_zeros(Flag *flag_stgy);

int flag_stgy_init_defined(Flag *flag_stgy);

int print_flag_stgy_init(Flag flag_stgy);

int define_blocksize(int n, int lambda, int bs_flag, int fixed_bs);

int define_nfft(int n_thread, int flag_nfft, int fixed_nfft);

int fftw_init_omp_threads();

int rhs_init_fftw(int *nfft, int fft_size, fftw_complex **V_fft, double **V_rfft, fftw_plan *plan_f, fftw_plan *plan_b, int fftw_flag);

int circ_init_fftw(double *T, int fft_size, int lambda, fftw_complex **T_fft);

int scmm_direct(int fft_size, int nfft, fftw_complex *C_fft, int ncol, double *V_rfft, double **CV, fftw_complex *V_fft, fftw_plan plan_f_V, fftw_plan plan_b_CV);

int scmm_basic(double **V, int blocksize, int m, fftw_complex *C_fft, double **CV, fftw_complex *V_fft, double *V_rfft, int nfft, fftw_plan plan_f_V, fftw_plan plan_b_CV);

int stmm_simple_basic(double **V, int n, int m, double *T, int lambda, double **TV);

int build_gappy_blocks(int nrow, int m, Block *tpltzblocks, int nb_blocks_local, int nb_blocks_all, int64_t *id0gap, int *lgap, int ngap, Block *tpltzblocks_gappy, int *nb_blocks_gappy_final, int flag_param_distmin_fixed);


//Internal routines (group 22)
int print_error_message(int error_number, char const *file, int line);

int copy_block(int ninrow, int nincol, double *Vin, int noutrow, int noutcol, double *Vout, int inrow, int incol, int nblockrow, int nblockcol, int outrow, int outcol, double norm, int set_zero_flag);

int vect2nfftblock(double *V1, int v1_size, double *V2, int fft_size, int nfft, int lambda);

int nfftblock2vect(double *V2, int fft_size, int nfft, int lambda, double *V1, int v1_size);

int get_overlapping_blocks_params(int nbloc, Block *tpltzblocks, int local_V_size, int64_t nrow, int64_t idp, int64_t *idpnew, int *local_V_size_new, int *nnew, int *ifirstBlock, int *ilastBlock);



//Wizard routines
int stbmmProd( Tpltz Nm1, double *V);

//=========================================================================
#endif      /* !TOEPLITZ_H_ */