Cyclops Tensor Framework
parallel arithmetic on multidimensional arrays
CTF::Tensor< dtype > Class Template Reference

an instance of a tensor within a CTF world More...

#include <tensor.h>

Inheritance diagram for CTF::Tensor< dtype >:
Collaboration diagram for CTF::Tensor< dtype >:

Public Member Functions

 Tensor ()
 default constructor More...
 
 Tensor (int order, int const *len, int const *sym, World &wrld=get_universe(), char const *name=NULL, bool profile=0, CTF_int::algstrct const &sr=Ring< dtype >())
 defines tensor filled with zeros on the default algstrct More...
 
 Tensor (int order, int const *len, int const *sym, World &wrld, CTF_int::algstrct const &sr, char const *name=NULL, bool profile=0)
 defines a tensor filled with zeros on a specified algstrct More...
 
 Tensor (int order, int const *len, World &wrld=get_universe(), CTF_int::algstrct const &sr=Ring< dtype >(), char const *name=NULL, bool profile=0)
 defines a nonsymmetric tensor filled with zeros on a specified algstrct More...
 
 Tensor (int order, bool is_sparse, int const *len, int const *sym, World &wrld=get_universe(), CTF_int::algstrct const &sr=Ring< dtype >(), char const *name=NULL, bool profile=0)
 defines a (sparse) tensor on a specified algstrct More...
 
 Tensor (int order, bool is_sparse, int const *len, World &wrld=get_universe(), CTF_int::algstrct const &sr=Ring< dtype >(), char const *name=NULL, bool profile=0)
 defines a nonsymmetric tensor filled with zeros on a specified algstrct More...
 
 Tensor (Tensor< dtype > const &A)
 copies a tensor, copying the data of A More...
 
 Tensor (tensor const &A)
 copies a tensor, copying the data of A (same as above) More...
 
 Tensor (bool copy, tensor const &A)
 copies a tensor (setting data to zero or copying A) More...
 
 Tensor (tensor &A, int const *new_sym)
 repacks the tensor other to a different symmetry (assumes existing data contains the symmetry and keeps only values with indices in increasing order) WARN: LIMITATION: new_sym must cause unidirectional structural changes, i.e. {NS,NS}->{SY,NS} OK, {SY,NS}->{NS,NS} OK, {NS,NS,SY,NS}->{SY,NS,NS,NS} NOT OK! More...
 
 Tensor (tensor const &A, World &wrld)
 creates a zeroed out copy (data not copied) of a tensor in a different world More...
 
 Tensor (int order, int const *len, int const *sym, World &wrld, char const *idx, Idx_Partition const &prl, Idx_Partition const &blk=Idx_Partition(), char const *name=NULL, bool profile=0, CTF_int::algstrct const &sr=Ring< dtype >())
 defines tensor filled with zeros on the default algstrct on a user-specified distributed layout More...
 
 Tensor (int order, bool is_sparse, int const *len, int const *sym, World &wrld, char const *idx, Idx_Partition const &prl, Idx_Partition const &blk=Idx_Partition(), char const *name=NULL, bool profile=0, CTF_int::algstrct const &sr=Ring< dtype >())
 defines tensor filled with zeros on the default algstrct on a user-specified distributed layout More...
 
Typ_Idx_Tensor< dtype > operator[] (char const *idx_map)
 associated an index map with the tensor for future operation More...
 
Typ_Idx_Tensor< dtype > i (char const *idx_map)
 
void read (int64_t npair, Pair< dtype > *pairs)
 Gives the values associated with any set of indices. More...
 
void read (int64_t npair, int64_t const *global_idx, dtype *data)
 Gives the values associated with any set of indices The sparse data is defined in coordinate format. The tensor index (i,j,k,l) of a tensor with edge lengths {m,n,p,q} is associated with the global index g via the formula g=i+j*m+k*m*n+l*m*n*p. The row index is first and the column index is second for matrices, which means they are column major. More...
 
void read (int64_t npair, dtype alpha, dtype beta, Pair< dtype > *pairs)
 sparse read with accumulation: pairs[i].d = alpha*A[pairs[i].k]+beta*pairs[i].d More...
 
void read (int64_t npair, dtype alpha, dtype beta, int64_t const *global_idx, dtype *data)
 sparse read with accumulation: data[i] = alpha*A[global_idx[i]]+beta*data[i] More...
 
void get_local_data (int64_t *npair, int64_t **global_idx, dtype **data, bool nonzeros_only=false, bool unpack_sym=false) const
 Gives the global indices and values associated with the local data. More...
 
void read_local (int64_t *npair, int64_t **global_idx, dtype **data, bool unpack_sym=false) const
 Using get_local_data(), which returns an array that must be freed with delete [], is more efficient, this version exists for backwards compatibility. gives the global indices and values associated with the local data WARNING-1: for sparse tensors this includes the zeros to maintain consistency with the behavior for dense tensors, use get_local_pairs to get only nonzeros. More...
 
void get_local_pairs (int64_t *npair, Pair< dtype > **pairs, bool nonzeros_only=false, bool unpack_sym=false) const
 gives the global indices and values associated with the local data More...
 
void read_local (int64_t *npair, Pair< dtype > **pairs, bool unpack_sym=false) const
 gives the global indices and values associated with the local data WARNING: for sparse tensors this includes the zeros to maintain consistency with the behavior for dense tensors, use get_lcoal_pairs to get only nonzeros More...
 
void get_all_data (int64_t *npair, dtype **data, bool unpack=false)
 collects the entire tensor data on each process (not memory scalable) More...
 
void read_all (int64_t *npair, dtype **data, bool unpack=false)
 collects the entire tensor data on each process (not memory scalable) More...
 
int64_t read_all (dtype *data, bool unpack=false)
 collects the entire tensor data on each process (not memory scalable) More...
 
void write (int64_t npair, int64_t const *global_idx, dtype const *data)
 writes in values associated with any set of indices The sparse data is defined in coordinate format. The tensor index (i,j,k,l) of a tensor with edge lengths {m,n,p,q} is associated with the global index g via the formula g=i+j*m+k*m*n+l*m*n*p. The row index is first and the column index is second for matrices, which means they are column major. More...
 
void write (int64_t npair, Pair< dtype > const *pairs)
 writes in values associated with any set of indices More...
 
void write (int64_t npair, dtype alpha, dtype beta, int64_t const *global_idx, dtype const *data)
 sparse add: A[global_idx[i]] = beta*A[global_idx[i]]+alpha*data[i] More...
 
void write (int64_t npair, dtype alpha, dtype beta, Pair< dtype > const *pairs)
 sparse add: A[pairs[i].k] = alpha*A[pairs[i].k]+beta*pairs[i].d More...
 
void contract (dtype alpha, CTF_int::tensor &A, char const *idx_A, CTF_int::tensor &B, char const *idx_B, dtype beta, char const *idx_C)
 contracts C[idx_C] = beta*C[idx_C] + alpha*A[idx_A]*B[idx_B] More...
 
void contract (dtype alpha, CTF_int::tensor &A, char const *idx_A, CTF_int::tensor &B, char const *idx_B, dtype beta, char const *idx_C, Bivar_Function< dtype > fseq)
 contracts computes C[idx_C] = beta*C[idx_C] fseq(alpha*A[idx_A],B[idx_B]) More...
 
void sum (dtype alpha, CTF_int::tensor &A, char const *idx_A, dtype beta, char const *idx_B)
 sums B[idx_B] = beta*B[idx_B] + alpha*A[idx_A] More...
 
void sum (dtype alpha, CTF_int::tensor &A, char const *idx_A, dtype beta, char const *idx_B, Univar_Function< dtype > fseq)
 sums B[idx_B] = beta*B[idx_B] + fseq(alpha*A[idx_A]) More...
 
void scale (dtype alpha, char const *idx_A)
 scales A[idx_A] = alpha*A[idx_A] More...
 
void scale (dtype alpha, char const *idx_A, Endomorphism< dtype > fseq)
 scales A[idx_A] = fseq(alpha*A[idx_A]) More...
 
dtype * get_mapped_data (char const *idx, Idx_Partition const &prl, Idx_Partition const &blk=Idx_Partition(), bool unpack=true)
 returns local data of tensor with parallel distribution prl and local blocking blk More...
 
void write (char const *idx, dtype const *data, Idx_Partition const &prl, Idx_Partition const &blk=Idx_Partition(), bool unpack=true)
 writes data to tensor from parallel distribution prl and local blocking blk More...
 
double estimate_time (CTF_int::tensor &A, char const *idx_A, CTF_int::tensor &B, char const *idx_B, char const *idx_C)
 estimate the time of a contraction C[idx_C] = A[idx_A]*B[idx_B] More...
 
double estimate_time (CTF_int::tensor &A, char const *idx_A, char const *idx_B)
 estimate the time of a sum B[idx_B] = A[idx_A] More...
 
Tensor< dtype > slice (int const *offsets, int const *ends) const
 cuts out a slice (block) of this tensor A[offsets,ends) result will always be fully nonsymmetric More...
 
Tensor< dtype > slice (int64_t corner_off, int64_t corner_end) const
 cuts out a slice (block) of this tensor with corners specified by global index result will always be fully nonsymmetric More...
 
Tensor< dtype > slice (int const *offsets, int const *ends, World *oworld) const
 cuts out a slice (block) of this tensor A[offsets,ends) result will always be fully nonsymmetric More...
 
Tensor< dtype > slice (int64_t corner_off, int64_t corner_end, World *oworld) const
 cuts out a slice (block) of this tensor with corners specified by global index result will always be fully nonsymmetric More...
 
void slice (int const *offsets, int const *ends, dtype beta, CTF_int::tensor const &A, int const *offsets_A, int const *ends_A, dtype alpha)
 adds to a slice (block) of this tensor = B B[offsets,ends)=beta*B[offsets,ends) + alpha*A[offsets_A,ends_A) More...
 
void slice (int64_t corner_off, int64_t corner_end, dtype beta, CTF_int::tensor const &A, int64_t corner_off_A, int64_t corner_end_A, dtype alpha)
 adds to a slice (block) of this tensor = B B[offsets,ends)=beta*B[offsets,ends) + alpha*A[offsets_A,ends_A) More...
 
void permute (dtype beta, CTF_int::tensor &A, int *const *perms_A, dtype alpha)
 Apply permutation to matrix, potentially extracting a slice B[i,j,...] = beta*B[...] + alpha*A[perms_A[0][i],perms_A[1][j],...]. More...
 
void permute (int *const *perms_B, dtype beta, CTF_int::tensor &A, dtype alpha)
 Apply permutation to matrix, potentially extracting a slice B[perms_B[0][i],perms_B[0][j],...] = beta*B[...] + alpha*A[i,j,...]. More...
 
void sparsify ()
 reduce tensor to sparse format, storing only nonzero data, or data above a specified threshold. makes dense tensors sparse. cleans sparse tensors of any 'computed' zeros. More...
 
void sparsify (dtype threshold, bool take_abs=true)
 reduce tensor to sparse format, storing only nonzero data, or data above a specified threshold. makes dense tensors sparse. cleans sparse tensors of any small values More...
 
void sparsify (std::function< bool(dtype)> filter)
 sparsifies tensor keeping only values v such that filter(v) = true More...
 
void read_sparse_from_file (const char *fpath, bool with_vals=true)
 read sparse tensor from file, entries of tensor must be stored one per line, as i_1 ... i_order v, to create entry T[i_1, ..., i_order] = v or as i_1 ... i_order, to create entry T[i_1, ..., i_order] = mulid More...
 
void write_sparse_to_file (const char *fpath, bool with_vals=true)
 write sparse tensor to file, entries of tensor will be stored one per line, as i_1 ... i_order v, corresponding to entry T[i_1, ..., i_order] = v or as i_1 ... i_order if with_vals =false More...
 
void add_to_subworld (Tensor< dtype > *tsr, dtype alpha, dtype beta)
 accumulates this tensor to a tensor object defined on a different world More...
 
void add_to_subworld (Tensor< dtype > *tsr)
 accumulates this tensor to a tensor object defined on a different world More...
 
void add_from_subworld (Tensor< dtype > *tsr, dtype alpha, dtype beta)
 accumulates this tensor from a tensor object defined on a different world More...
 
void add_from_subworld (Tensor< dtype > *tsr)
 accumulates this tensor from a tensor object defined on a different world More...
 
void align (CTF_int::tensor const &A)
 aligns data mapping with tensor A More...
 
dtype reduce (OP op)
 performs a reduction on the tensor More...
 
dtype norm1 ()
 computes the entrywise 1-norm of the tensor More...
 
dtype norm2 ()
 computes the frobenius norm of the tensor (needs sqrt()!) More...
 
dtype norm_infty ()
 finds the max absolute value element of the tensor More...
 
void norm1 (double &nrm)
 computes the entrywise 1-norm of the tensor More...
 
void norm2 (double &nrm)
 computes the frobenius norm of the tensor More...
 
void norm_infty (double &nrm)
 finds the max absolute value element of the tensor More...
 
dtype * get_raw_data (int64_t *size) const
 gives the raw current local data with padding included More...
 
const dtype * raw_data (int64_t *size) const
 gives a read-only copy of the raw current local data with padding included More...
 
void get_max_abs (int n, dtype *data) const
 obtains a small number of the biggest elements of the tensor in sorted order (e.g. eigenvalues) More...
 
void fill_random (dtype rmin, dtype rmax)
 fills local unique tensor elements to random values in the range [min,max] works only for dtype in {float,double,int,int64_t}, for others you can use Transform() uses mersenne twister seeded every time a global world is created (reseeded if all worlds destroyed) More...
 
void fill_sp_random (dtype rmin, dtype rmax, double frac_sp)
 generate roughly frac_sp*dense_tensor_size nonzeros between rmin and rmax, works only for dtype in {float,double,int,int64_t}, for others you can use Transform() uses mersenne twister seeded every time a global world is created (reseeded if all worlds destroyed) More...
 
void profile_on ()
 turns on profiling for tensor More...
 
void profile_off ()
 turns off profiling for tensor More...
 
void set_name (char const *name)
 sets tensor name More...
 
Tensor< dtype > & operator= (dtype val)
 sets all values in the tensor to val More...
 
Tensor< dtype > & operator= (const Tensor< dtype > A)
 sets the tensor More...
 
Sparse_Tensor< dtype > operator[] (std::vector< int64_t > indices)
 gives handle to sparse index subset of tensors More...
 
void print (FILE *fp, dtype cutoff) const
 prints tensor data to file using process 0 (modify print(...) overload in set.h if you would like a different print format) More...
 
void print (FILE *fp=stdout) const
 
void prnt () const
 
void compare (const Tensor< dtype > &A, FILE *fp=stdout, double cutoff=-1.0)
 prints two sets of tensor data side-by-side to file using process 0 More...
 
 ~Tensor ()
 frees CTF tensor More...
 
template<>
void read_sparse_from_file (const char *fpath, bool with_vals)
 
template<>
void read_sparse_from_file (const char *fpath, bool with_vals)
 
template<>
void read_sparse_from_file (const char *fpath, bool with_vals)
 
template<>
void read_sparse_from_file (const char *fpath, bool with_vals)
 
template<>
void write_sparse_to_file (const char *fpath, bool with_vals)
 
template<>
void write_sparse_to_file (const char *fpath, bool with_vals)
 
template<>
void write_sparse_to_file (const char *fpath, bool with_vals)
 
template<>
void write_sparse_to_file (const char *fpath, bool with_vals)
 
template<>
void fill_random (double rmin, double rmax)
 
template<>
void fill_random (float rmin, float rmax)
 
template<>
void fill_random (int64_t rmin, int64_t rmax)
 
template<>
void fill_random (int rmin, int rmax)
 
template<>
void fill_sp_random (double rmin, double rmax, double frac_sp)
 
template<>
void fill_sp_random (float rmin, float rmax, double frac_sp)
 
template<>
void fill_sp_random (int rmin, int rmax, double frac_sp)
 
template<>
void fill_sp_random (int64_t rmin, int64_t rmax, double frac_sp)
 
- Public Member Functions inherited from CTF_int::tensor
CTF::Idx_Tensor operator[] (char const *idx_map)
 associated an index map with the tensor for future operation More...
 
 tensor ()
 
 ~tensor ()
 class free self More...
 
void free_self ()
 destructor More...
 
 tensor (algstrct const *sr, int order, int const *edge_len, int const *sym, CTF::World *wrld, bool alloc_data=true, char const *name=NULL, bool profile=1, bool is_sparse=0)
 defines a tensor object with some mapping (if alloc_data) More...
 
 tensor (algstrct const *sr, int order, bool is_sparse, int const *edge_len, int const *sym, CTF::World *wrld, char const *idx, CTF::Idx_Partition const &prl, CTF::Idx_Partition const &blk, char const *name=NULL, bool profile=1)
 defines a tensor object with some mapping (if alloc_data) More...
 
 tensor (tensor const *other, bool copy=1, bool alloc_data=1)
 creates tensor copy, unfolds other if other is folded More...
 
 tensor (tensor *other, int const *new_sym)
 repacks the tensor other to a different symmetry (assumes existing data contains the symmetry and keeps only values with indices in increasing order) WARN: LIMITATION: new_sym must cause unidirectional structural changes, i.e. {NS,NS}->{SY,NS} OK, {SY,NS}->{NS,NS} OK, {NS,NS,SY,NS}->{SY,NS,NS,NS} NOT OK! More...
 
int * calc_phase () const
 compute the cyclic phase of each tensor dimension More...
 
int calc_tot_phase () const
 calculate the total number of blocks of the tensor More...
 
int64_t calc_nvirt () const
 calculate virtualization factor of tensor return virtualization factor More...
 
int64_t calc_npe () const
 calculate the number of processes this tensor is distributed over return number of processes owning a block of the tensor More...
 
void set_padding ()
 sets padding and local size of a tensor given a mapping More...
 
int set_zero ()
 elects a mapping and sets tensor data to zero More...
 
int set (char const *val)
 sets tensor data to val More...
 
int zero_out_padding ()
 sets padded portion of tensor to zero (this should be maintained internally) More...
 
void scale_diagonals (int const *sym_mask)
 scales each element by 1/(number of entries equivalent to it after permutation of indices for which sym_mask is 1) More...
 
void addinv ()
 
void print_map (FILE *stream=stdout, bool allcall=1) const
 displays mapping information More...
 
void set_name (char const *name)
 set the tensor name More...
 
char const * get_name () const
 get the tensor name More...
 
void profile_on ()
 turn on profiling More...
 
void profile_off ()
 turn off profiling More...
 
void get_raw_data (char **data, int64_t *size) const
 get raw data pointer without copy WARNING: includes padding More...
 
int write (int64_t num_pair, char const *alpha, char const *beta, char *mapped_data, char const rw='w')
 Add tensor data new=alpha*new+beta*old with <key, value> pairs where key is the global index for the value. More...
 
void write (int64_t num_pair, char const *alpha, char const *beta, int64_t const *inds, char const *data)
 Add tensor data new=alpha*new+beta*old with <key, value> pairs where key is the global index for the value. More...
 
void read (int64_t num_pair, char const *alpha, char const *beta, int64_t const *inds, char *data)
 read tensor data with <key, value> pairs where key is the global index for the value, which gets filled in with beta times the old values plus alpha times the values read from the tensor. More...
 
int read (int64_t num_pair, char const *alpha, char const *beta, char *mapped_data)
 read tensor data with <key, value> pairs where key is the global index for the value, which gets filled in with beta times the old values plus alpha times the values read from the tensor. More...
 
char * read (char const *idx, CTF::Idx_Partition const &prl, CTF::Idx_Partition const &blk, bool unpack)
 returns local data of tensor with parallel distribution prl and local blocking blk More...
 
int read (int64_t num_pair, char *mapped_data)
 read tensor data with <key, value> pairs where key is the global index for the value, which gets filled in. More...
 
int64_t get_tot_size (bool packed)
 get number of elements in whole tensor More...
 
int allread (int64_t *num_pair, char **all_data, bool unpack)
 read entire tensor with each processor (in packed layout). WARNING: will use an 'unscalable' amount of memory. More...
 
int allread (int64_t *num_pair, char *all_data, bool unpack=true)
 read entire tensor with each processor (in packed layout). WARNING: will use an 'unscalable' amount of memory. More...
 
void slice (int const *offsets_B, int const *ends_B, char const *beta, tensor *A, int const *offsets_A, int const *ends_A, char const *alpha)
 accumulates out a slice (block) of this tensor = B B[offsets,ends)=beta*B[offsets,ends) + alpha*A[offsets_A,ends_A) More...
 
int permute (tensor *A, int *const *permutation_A, char const *alpha, int *const *permutation_B, char const *beta)
 
int sparsify (char const *threshold=NULL, bool take_abs=true)
 reduce tensor to sparse format, storing only nonzero data, or data above a specified threshold. makes dense tensors sparse. cleans sparse tensors of any 'computed' zeros. More...
 
int sparsify (std::function< bool(char const *)> f)
 sparsifies tensor keeping only values v such that filter(v) = true More...
 
int read_local (int64_t *num_pair, char **mapped_data, bool unpack_sym=false) const
 read tensor data pairs local to processor including those with zero values WARNING: for sparse tensors this includes the zeros to maintain consistency with the behavior for dense tensors, use read_local_nnz to get only nonzeros More...
 
int read_local_nnz (int64_t *num_pair, char **mapped_data, bool unpack_sym=false) const
 read tensor data pairs local to processor that have nonzero values More...
 
void read_local (int64_t *num_pair, int64_t **inds, char **data, bool unpack_sym=false) const
 read tensor data pairs local to processor including those with zero values WARNING: for sparse tensors this includes the zeros to maintain consistency with the behavior for dense tensors, use read_local_nnz to get only nonzeros More...
 
void read_local_nnz (int64_t *num_pair, int64_t **inds, char **data, bool unpack_sym=false) const
 read tensor data pairs local to processor that have nonzero values More...
 
int align (tensor const *B)
 align mapping of thisa tensor to that of B More...
 
int reduce_sum (char *result)
 Performs an elementwise summation reduction on a tensor. More...
 
int reduce_sum (char *result, algstrct const *sr_other)
 Performs an elementwise summation reduction on a tensor with summation defined by sr_other. More...
 
int reduce_sumabs (char *result)
 Performs an elementwise absolute value summation reduction on a tensor. More...
 
int reduce_sumabs (char *result, algstrct const *sr_other)
 Performs an elementwise absolute value summation reduction on a tensor. More...
 
int reduce_sumsq (char *result)
 computes the sum of squares of the elements More...
 
int get_max_abs (int n, char *data) const
 obtains the largest n elements (in absolute value) of the tensor More...
 
void print (FILE *fp=stdout, char const *cutoff=NULL) const
 prints tensor data to file using process 0 More...
 
void prnt () const
 
void compare (const tensor *A, FILE *fp, char const *cutoff)
 prints two sets of tensor data side-by-side to file using process 0 More...
 
void orient_subworld (CTF::World *greater_world, int &bw_mirror_rank, int &fw_mirror_rank, distribution *&odst, char **sub_buffer_)
 maps data from this world (subcomm) to the correct order of processors with respect to a parent (greater_world) comm More...
 
void add_to_subworld (tensor *tsr_sub, char const *alpha, char const *beta)
 accumulates this tensor to a tensor object defined on a different world More...
 
void add_from_subworld (tensor *tsr_sub, char const *alpha, char const *beta)
 accumulates this tensor from a tensor object defined on a different world More...
 
void unfold (bool was_mod=0)
 undo the folding of a local tensor block unsets is_folded and deletes rec_tsr More...
 
void remove_fold ()
 removes folding without doing transpose unsets is_folded and deletes rec_tsr More...
 
double est_time_unfold ()
 estimate cost of potential transpose involved in undoing the folding of a local tensor block More...
 
void fold (int nfold, int const *fold_idx, int const *idx_map, int *all_fdim, int **all_flen)
 fold a tensor by putting the symmetry-preserved portion in the leading dimensions of the tensor sets is_folded and creates rec_tsr with aliased data More...
 
void pull_alias (tensor const *other)
 pulls data from an tensor with an aliased buffer More...
 
void clear_mapping ()
 zeros out mapping More...
 
int redistribute (distribution const &old_dist, int const *old_offsets=NULL, int *const *old_permutation=NULL, int const *new_offsets=NULL, int *const *new_permutation=NULL)
 permutes the data of a tensor to its new layout More...
 
double est_redist_time (distribution const &old_dist, double nnz_frac)
 
int64_t get_redist_mem (distribution const &old_dist, double nnz_frac)
 
int map_tensor_rem (int num_phys_dims, CommData *phys_comm, int fill=0)
 map the remainder of a tensor More...
 
int extract_diag (int const *idx_map, int rw, tensor *&new_tsr, int **idx_map_new)
 extracts the diagonal of a tensor if the index map specifies to do so More...
 
void set_sym (int const *sym)
 sets symmetry, WARNING: for internal use only !!!! More...
 
void set_new_nnz_glb (int64_t const *nnz_blk)
 sets the number of nonzeros both locally (nnz_loc) and overall globally (nnz_tot) More...
 
void spmatricize (int m, int n, int nrow_idx, bool csr)
 transposes local data in preparation for summation or contraction, transforms to COO or CSR format for sparse More...
 
void despmatricize (int nrow_idx, bool csr)
 transposes back local data from sparse matrix format to key-value pair format More...
 
void leave_home_with_buffer ()
 degister home buffer More...
 
void register_size (int64_t size)
 register buffer allocation for this tensor More...
 
void deregister_size ()
 deregister buffer allocation for this tensor More...
 
void write_dense_to_file (MPI_File &file, int64_t offset=0)
 write all tensor data to binary file in element order, unpacking from sparse or symmetric formats More...
 
void write_dense_to_file (char const *filename)
 write all tensor data to binary file in element order, unpacking from sparse or symmetric formats More...
 
void read_dense_from_file (MPI_File &file, int64_t offset=0)
 read all tensor data from binary file in element order, which should be stored as nonsymmetric and dense as done in write_dense_to_file() More...
 
void read_dense_from_file (char const *filename)
 read all tensor data from binary file in element order, which should be stored as nonsymmetric and dense as done in write_dense_to_file() More...
 
template<typename dtype_A , typename dtype_B >
void conv_type (tensor *B)
 convert this tensor from dtype_A to dtype_B and store the result in B (primarily needed for python interface) More...
 
template<typename dtype_A , typename dtype_B >
void exp_helper (tensor *A)
 
template<typename dtype >
void compare_elementwise (tensor *A, tensor *B)
 do an elementwise comparison(==) of two tensors with elements of type dtype (primarily needed for python interface), store result in this tensor (has to be boolean tensor) More...
 
template<typename dtype >
void not_equals (tensor *A, tensor *B)
 do an elementwise comparison(!=) of two tensors with elements of type dtype (primarily needed for python interface), store result in this tensor (has to be boolean tensor) More...
 
template<typename dtype >
void smaller_than (tensor *A, tensor *B)
 do an elementwise comparison(<) of two tensors with elements of type dtype (primarily needed for python interface), store result in this tensor (has to be boolean tensor) More...
 
template<typename dtype >
void smaller_equal_than (tensor *A, tensor *B)
 do an elementwise comparison(<=) of two tensors with elements of type dtype (primarily needed for python interface), store result in this tensor (has to be boolean tensor) More...
 
template<typename dtype >
void larger_than (tensor *A, tensor *B)
 do an elementwise comparison(>) of two tensors with elements of type dtype (primarily needed for python interface), store result in this tensor (has to be boolean tensor) More...
 
template<typename dtype >
void larger_equal_than (tensor *A, tensor *B)
 do an elementwise comparison(>=) of two tensors with elements of type dtype (primarily needed for python interface), store result in this tensor (has to be boolean tensor) More...
 
template<typename dtype >
void true_divide (tensor *A)
 
tensorself_reduce (int const *idx_A, int **new_idx_A, int order_B, int const *idx_B, int **new_idx_B, int order_C=0, int const *idx_C=NULL, int **new_idx_C=NULL)
 performs a partial reduction on the tensor (used in summation and contraction) More...
 

Additional Inherited Members

- Data Fields inherited from CTF_int::tensor
CTF::Worldwrld
 distributed processor context on which tensor is defined More...
 
algstrctsr
 algstrct on which tensor elements and operations are defined More...
 
int * sym
 symmetries among tensor dimensions More...
 
int order
 number of tensor dimensions More...
 
int * lens
 unpadded tensor edge lengths More...
 
int * pad_edge_len
 padded tensor edge lengths More...
 
int * padding
 padding along each edge length (less than distribution phase) More...
 
char * name
 name given to tensor More...
 
int is_scp_padded
 whether tensor data has additional padding More...
 
int * scp_padding
 additional padding, may be greater than ScaLAPACK phase More...
 
int * sym_table
 order-by-order table of dimensional symmetry relations More...
 
bool is_mapped
 whether a mapping has been selected More...
 
topologytopo
 topology to which the tensor is mapped More...
 
mappingedge_map
 mappings of each tensor dimension onto topology dimensions More...
 
int64_t size
 current size of local tensor data chunk (mapping-dependent) More...
 
int64_t registered_alloc_size
 size CTF keeps track of for memory usage More...
 
bool is_folded
 whether the data is folded/transposed into a (lower-order) tensor More...
 
int * inner_ordering
 ordering of the dimensions according to which the tensori s folded More...
 
tensorrec_tsr
 representation of folded tensor (shares data pointer) More...
 
bool is_cyclic
 whether the tensor data is cyclically distributed (blocked if false) More...
 
bool is_data_aliased
 whether the tensor data is an alias of another tensor object's data More...
 
tensorslay
 tensor object associated with tensor object whose data pointer needs to be preserved, needed for ScaLAPACK wrapper FIXME: home buffer should presumably take care of this... More...
 
bool has_zero_edge_len
 if true tensor has a zero edge length, so is zero, which short-cuts stuff More...
 
char * data
 tensor data, either the data or the key-value pairs should exist at any given time More...
 
bool has_home
 whether the tensor has a home mapping/buffer More...
 
char * home_buffer
 buffer associated with home mapping of tensor, to which it is returned More...
 
int64_t home_size
 size of home buffer More...
 
bool is_home
 whether the latest tensor data is in the home buffer More...
 
bool left_home_transp
 whether the tensor left home to transpose More...
 
bool profile
 whether profiling should be done for contractions/sums involving this tensor More...
 
bool is_sparse
 whether only the non-zero elements of the tensor are stored More...
 
bool is_csr
 whether CSR or COO if folded More...
 
int nrow_idx
 how many modes are folded into matricized row More...
 
int64_t nnz_loc
 number of local nonzero elements More...
 
int64_t nnz_tot
 maximum number of local nonzero elements over all procs More...
 
int64_t * nnz_blk
 nonzero elements in each block owned locally More...
 
- Protected Member Functions inherited from CTF_int::tensor
void init (algstrct const *sr, int order, int const *edge_len, int const *sym, CTF::World *wrld, bool alloc_data, char const *name, bool profile, bool is_sparse)
 initializes tensor data More...
 
PairIterator read_all_pairs (int64_t *num_pair, bool unpack)
 read all pairs with each processor (packed) More...
 
void copy_tensor_data (tensor const *other)
 copies all tensor data from other More...
 
void set_distribution (char const *idx, CTF::Idx_Partition const &prl, CTF::Idx_Partition const &blk)
 set edge mappings as specified More...
 

Detailed Description

template<typename dtype = double>
class CTF::Tensor< dtype >

an instance of a tensor within a CTF world

Definition at line 74 of file tensor.h.

Constructor & Destructor Documentation

template<typename dtype >
CTF::Tensor< dtype >::Tensor ( )

default constructor

Definition at line 12 of file tensor.cxx.

template<typename dtype >
CTF::Tensor< dtype >::Tensor ( int  order,
int const *  len,
int const *  sym,
World wrld = get_universe(),
char const *  name = NULL,
bool  profile = 0,
CTF_int::algstrct const &  sr = Ring<dtype>() 
)

defines tensor filled with zeros on the default algstrct

Parameters
[in]ordernumber of dimensions of tensor
[in]lenedge lengths of tensor
[in]symsymmetries of tensor (e.g. symmetric matrix -> sym={SY, NS})
[in]wrlda world for the tensor to live in
[in]namean optionary name for the tensor
[in]profileset to 1 to profile contractions involving this tensor
[in]srdefines the tensor arithmetic for this tensor

Definition at line 16 of file tensor.cxx.

References ctf.core::dtype, CTF_int::accumulatable::el_size, and IASSERT.

template<typename dtype >
CTF::Tensor< dtype >::Tensor ( int  order,
int const *  len,
int const *  sym,
World wrld,
CTF_int::algstrct const &  sr,
char const *  name = NULL,
bool  profile = 0 
)

defines a tensor filled with zeros on a specified algstrct

Parameters
[in]ordernumber of dimensions of tensor
[in]lenedge lengths of tensor
[in]symsymmetries of tensor (e.g. symmetric matrix -> sym={SY, NS})
[in]wrlda world for the tensor to live in
[in]srdefines the tensor arithmetic for this tensor
[in]namean optionary name for the tensor
[in]profileset to 1 to profile contractions involving this tensor

Definition at line 28 of file tensor.cxx.

References ctf.core::dtype, CTF_int::accumulatable::el_size, and IASSERT.

template<typename dtype >
CTF::Tensor< dtype >::Tensor ( int  order,
int const *  len,
World wrld = get_universe(),
CTF_int::algstrct const &  sr = Ring<dtype>(),
char const *  name = NULL,
bool  profile = 0 
)

defines a nonsymmetric tensor filled with zeros on a specified algstrct

Parameters
[in]ordernumber of dimensions of tensor
[in]lenedge lengths of tensor
[in]wrlda world for the tensor to live in
[in]srdefines the tensor arithmetic for this tensor
[in]namean optionary name for the tensor
[in]profileset to 1 to profile contractions involving this tensor

Definition at line 67 of file tensor.cxx.

References ctf.core::dtype, CTF_int::accumulatable::el_size, and IASSERT.

template<typename dtype >
CTF::Tensor< dtype >::Tensor ( int  order,
bool  is_sparse,
int const *  len,
int const *  sym,
World wrld = get_universe(),
CTF_int::algstrct const &  sr = Ring<dtype>(),
char const *  name = NULL,
bool  profile = 0 
)

defines a (sparse) tensor on a specified algstrct

Parameters
[in]ordernumber of dimensions of tensor
[in]is_sparseif 1 then tensor will be sparse and non-trivial elements won't be stored
[in]lenedge lengths of tensor
[in]symsymmetries of tensor (e.g. symmetric matrix -> sym={SY, NS})
[in]wrlda world for the tensor to live in
[in]srdefines the tensor arithmetic for this tensor
[in]namean optionary name for the tensor
[in]profileset to 1 to profile contractions involving this tensor

Definition at line 41 of file tensor.cxx.

References ctf.core::dtype, CTF_int::accumulatable::el_size, and IASSERT.

template<typename dtype >
CTF::Tensor< dtype >::Tensor ( int  order,
bool  is_sparse,
int const *  len,
World wrld = get_universe(),
CTF_int::algstrct const &  sr = Ring<dtype>(),
char const *  name = NULL,
bool  profile = 0 
)

defines a nonsymmetric tensor filled with zeros on a specified algstrct

Parameters
[in]ordernumber of dimensions of tensor
[in]is_sparseif 1 then tensor will be sparse and non-trivial elements won't be stored
[in]lenedge lengths of tensor
[in]wrlda world for the tensor to live in
[in]srdefines the tensor arithmetic for this tensor
[in]namean optionary name for the tensor
[in]profileset to 1 to profile contractions involving this tensor

Definition at line 54 of file tensor.cxx.

References ctf.core::dtype, CTF_int::accumulatable::el_size, and IASSERT.

template<typename dtype >
CTF::Tensor< dtype >::Tensor ( Tensor< dtype > const &  A)

copies a tensor, copying the data of A

Parameters
[in]Atensor to copy

Definition at line 116 of file tensor.cxx.

template<typename dtype >
CTF::Tensor< dtype >::Tensor ( tensor const &  A)

copies a tensor, copying the data of A (same as above)

Parameters
[in]Atensor to copy

Definition at line 120 of file tensor.cxx.

template<typename dtype >
CTF::Tensor< dtype >::Tensor ( bool  copy,
tensor const &  A 
)

copies a tensor (setting data to zero or copying A)

Parameters
[in]copywhether to copy the data of A into the new tensor
[in]Atensor to copy

Definition at line 111 of file tensor.cxx.

template<typename dtype >
CTF::Tensor< dtype >::Tensor ( tensor A,
int const *  new_sym 
)

repacks the tensor other to a different symmetry (assumes existing data contains the symmetry and keeps only values with indices in increasing order) WARN: LIMITATION: new_sym must cause unidirectional structural changes, i.e. {NS,NS}->{SY,NS} OK, {SY,NS}->{NS,NS} OK, {NS,NS,SY,NS}->{SY,NS,NS,NS} NOT OK!

Parameters
[in]Atensor to copy
[in]new_symnew symmetry array (replaces this->sym)

Definition at line 129 of file tensor.cxx.

References CTF::Tensor< dtype >::operator[]().

template<typename dtype >
CTF::Tensor< dtype >::Tensor ( tensor const &  A,
World wrld 
)

creates a zeroed out copy (data not copied) of a tensor in a different world

Parameters
[in]Atensor whose characteristics to copy
[in]wrlda world for the tensor we are creating to live in, can be different from A

Definition at line 124 of file tensor.cxx.

template<typename dtype >
CTF::Tensor< dtype >::Tensor ( int  order,
int const *  len,
int const *  sym,
World wrld,
char const *  idx,
Idx_Partition const &  prl,
Idx_Partition const &  blk = Idx_Partition(),
char const *  name = NULL,
bool  profile = 0,
CTF_int::algstrct const &  sr = Ring<dtype>() 
)

defines tensor filled with zeros on the default algstrct on a user-specified distributed layout

Parameters
[in]ordernumber of dimensions of tensor
[in]lenedge lengths of tensor
[in]symsymmetries of tensor (e.g. symmetric matrix -> sym={SY, NS})
[in]wrlda world for the tensor to live in
[in]idxassignment of characters to each dim
[in]prlmesh processor topology with character labels
[in]blklocal blocking with processor labels
[in]namean optionary name for the tensor
[in]profileset to 1 to profile contractions involving this tensor
[in]srdefines the tensor arithmetic for this tensor

Definition at line 79 of file tensor.cxx.

References ctf.core::dtype, CTF_int::accumulatable::el_size, IASSERT, and CTF_int::tensor::sr.

template<typename dtype >
CTF::Tensor< dtype >::Tensor ( int  order,
bool  is_sparse,
int const *  len,
int const *  sym,
World wrld,
char const *  idx,
Idx_Partition const &  prl,
Idx_Partition const &  blk = Idx_Partition(),
char const *  name = NULL,
bool  profile = 0,
CTF_int::algstrct const &  sr = Ring<dtype>() 
)

defines tensor filled with zeros on the default algstrct on a user-specified distributed layout

Parameters
[in]ordernumber of dimensions of tensor
[in]is_sparsewhether tensor is sparse
[in]lenedge lengths of tensor
[in]symsymmetries of tensor (e.g. symmetric matrix -> sym={SY, NS})
[in]wrlda world for the tensor to live in
[in]idxassignment of characters to each dim
[in]prlmesh processor topology with character labels
[in]blklocal blocking with processor labels
[in]namean optionary name for the tensor
[in]profileset to 1 to profile contractions involving this tensor
[in]srdefines the tensor arithmetic for this tensor

Definition at line 94 of file tensor.cxx.

References ctf.core::dtype, CTF_int::accumulatable::el_size, IASSERT, and CTF_int::tensor::sr.

template<typename dtype >
CTF::Tensor< dtype >::~Tensor ( )

frees CTF tensor

Definition at line 149 of file tensor.cxx.

Member Function Documentation

template<typename dtype >
void CTF::Tensor< dtype >::add_from_subworld ( Tensor< dtype > *  tsr,
dtype  alpha,
dtype  beta 
)

accumulates this tensor from a tensor object defined on a different world

Parameters
[in]tsra tensor object of the same characteristic as this tensor, but on a different world/MPI_comm
[in]alphascaling factor for tensor tsr (default 1.0)
[in]betascaling factor for this tensor (default 1.0)

Definition at line 560 of file tensor.cxx.

References CTF_int::tensor::add_from_subworld(), CTF_int::algstrct::clone(), CTF_int::tensor::sr, and CTF_int::tensor::tensor().

Referenced by readall_test(), and test_subworld_gemm().

template<typename dtype >
void CTF::Tensor< dtype >::add_from_subworld ( Tensor< dtype > *  tsr)

accumulates this tensor from a tensor object defined on a different world

Parameters
[in]tsra tensor object of the same characteristic as this tensor, but on a different world/MPI_comm

Definition at line 574 of file tensor.cxx.

References CTF_int::tensor::add_from_subworld(), CTF_int::algstrct::clone(), CTF_int::algstrct::mulid(), CTF_int::tensor::sr, and CTF_int::tensor::tensor().

template<typename dtype >
void CTF::Tensor< dtype >::add_to_subworld ( Tensor< dtype > *  tsr,
dtype  alpha,
dtype  beta 
)

accumulates this tensor to a tensor object defined on a different world

Parameters
[in]tsra tensor object of the same characteristic as this tensor, but on a different world/MPI_comm
[in]alphascaling factor for this tensor (default 1.0)
[in]betascaling factor for tensor tsr (default 1.0)

Definition at line 540 of file tensor.cxx.

References CTF_int::tensor::add_to_subworld(), CTF_int::algstrct::clone(), CTF_int::tensor::sr, and CTF_int::tensor::tensor().

Referenced by CTF::Tensor< dtype >::add_to_subworld(), and test_subworld_gemm().

template<typename dtype >
void CTF::Tensor< dtype >::add_to_subworld ( Tensor< dtype > *  tsr)

accumulates this tensor to a tensor object defined on a different world

Parameters
[in]tsra tensor object of the same characteristic as this tensor, but on a different world/MPI_comm

Definition at line 554 of file tensor.cxx.

References CTF::Tensor< dtype >::add_to_subworld(), CTF_int::algstrct::mulid(), and CTF_int::tensor::sr.

template<typename dtype = double>
void CTF::Tensor< dtype >::align ( CTF_int::tensor const &  A)

aligns data mapping with tensor A

Parameters
[in]Aalign with this tensor

Definition at line 720 of file tensor.cxx.

References CTF_int::tensor::align(), CTF::World::cdt, CTF_int::CommData::cm, IASSERT, CTF_int::SUCCESS, and CTF_int::tensor::wrld.

Referenced by gemm_4D().

template<typename dtype >
void CTF::Tensor< dtype >::compare ( const Tensor< dtype > &  A,
FILE *  fp = stdout,
double  cutoff = -1.0 
)

prints two sets of tensor data side-by-side to file using process 0

Parameters
[in]fpfile to print to e.g. stdout
[in]Atensor to compare against
[in]cutoffdo not print values of absolute value smaller than this

Definition at line 424 of file tensor.cxx.

References CTF_int::tensor::compare().

template<typename dtype = double>
void CTF::Tensor< dtype >::contract ( dtype  alpha,
CTF_int::tensor A,
char const *  idx_A,
CTF_int::tensor B,
char const *  idx_B,
dtype  beta,
char const *  idx_C 
)

contracts C[idx_C] = beta*C[idx_C] + alpha*A[idx_A]*B[idx_B]

Parameters
[in]alphaA*B scaling factor
[in]Afirst operand tensor
[in]idx_Aindices of A in contraction, e.g. "ik" -> A_{ik}
[in]Bsecond operand tensor
[in]idx_Bindices of B in contraction, e.g. "kj" -> B_{kj}
[in]betaC scaling factor
[in]idx_Cindices of C (this tensor), e.g. "ij" -> C_{ij}

Referenced by CTF::Tensor< dtype >::fill_sp_random(), and weigh_4D().

template<typename dtype = double>
void CTF::Tensor< dtype >::contract ( dtype  alpha,
CTF_int::tensor A,
char const *  idx_A,
CTF_int::tensor B,
char const *  idx_B,
dtype  beta,
char const *  idx_C,
Bivar_Function< dtype >  fseq 
)

contracts computes C[idx_C] = beta*C[idx_C] fseq(alpha*A[idx_A],B[idx_B])

Parameters
[in]alphaA*B scaling factor
[in]Afirst operand tensor
[in]idx_Aindices of A in contraction, e.g. "ik" -> A_{ik}
[in]Bsecond operand tensor
[in]idx_Bindices of B in contraction, e.g. "kj" -> B_{kj}
[in]betaC scaling factor
[in]idx_Cindices of C (this tensor), e.g. "ij" -> C_{ij}
[in]fseqsequential operation to execute, default is multiply-add
template<typename dtype = double>
double CTF::Tensor< dtype >::estimate_time ( CTF_int::tensor A,
char const *  idx_A,
CTF_int::tensor B,
char const *  idx_B,
char const *  idx_C 
)

estimate the time of a contraction C[idx_C] = A[idx_A]*B[idx_B]

Parameters
[in]Afirst operand tensor
[in]idx_Aindices of A in contraction, e.g. "ik" -> A_{ik}
[in]Bsecond operand tensor
[in]idx_Bindices of B in contraction, e.g. "kj" -> B_{kj}
[in]idx_Cindices of C (this tensor), e.g. "ij" -> C_{ij}
Returns
time in seconds, at the moment not at all precise

Referenced by CTF::Tensor< dtype >::operator=().

template<typename dtype = double>
double CTF::Tensor< dtype >::estimate_time ( CTF_int::tensor A,
char const *  idx_A,
char const *  idx_B 
)

estimate the time of a sum B[idx_B] = A[idx_A]

Parameters
[in]Afirst operand tensor
[in]idx_Aindices of A in contraction, e.g. "ik" -> A_{ik}
[in]idx_Bindices of B in contraction, e.g. "kj" -> B_{kj}
Returns
time in seconds, at the moment not at all precise
template<typename dtype >
void CTF::Tensor< dtype >::fill_random ( dtype  rmin,
dtype  rmax 
)

fills local unique tensor elements to random values in the range [min,max] works only for dtype in {float,double,int,int64_t}, for others you can use Transform() uses mersenne twister seeded every time a global world is created (reseeded if all worlds destroyed)

Parameters
[in]rminminimum random value
[in]rmaxmaximum random value

Definition at line 928 of file tensor.cxx.

References IASSERT, CTF_int::tensor::name, CTF::World::rank, and CTF_int::tensor::wrld.

Referenced by bench_ao_mo_transf(), bench_redistribution(), bitonic(), bivar_function(), bivar_transform(), block_sparse(), ccsdt_t3_to_t2(), checkpoint(), endomorphism(), fast_sym_4D(), jacobi(), main(), matmul(), neural(), reduce_bcast(), scan_test(), sparse_mp3(), spectral(), spmv(), test_ao_mo_transf(), test_qr(), test_svd(), train_dns_vec_mat(), train_off_vec_mat(), train_ttm(), and univar_function().

template<>
void CTF::Tensor< double >::fill_random ( double  rmin,
double  rmax 
)
inline

Definition at line 948 of file tensor.cxx.

template<>
void CTF::Tensor< float >::fill_random ( float  rmin,
float  rmax 
)
inline

Definition at line 953 of file tensor.cxx.

template<>
void CTF::Tensor< int64_t >::fill_random ( int64_t  rmin,
int64_t  rmax 
)
inline

Definition at line 958 of file tensor.cxx.

template<>
void CTF::Tensor< int >::fill_random ( int  rmin,
int  rmax 
)
inline

Definition at line 963 of file tensor.cxx.

template<typename dtype >
void CTF::Tensor< dtype >::fill_sp_random ( dtype  rmin,
dtype  rmax,
double  frac_sp 
)

generate roughly frac_sp*dense_tensor_size nonzeros between rmin and rmax, works only for dtype in {float,double,int,int64_t}, for others you can use Transform() uses mersenne twister seeded every time a global world is created (reseeded if all worlds destroyed)

Parameters
[in]rminminimum random value
[in]rmaxmaximum random value
[in]frac_spdesired expected nonzero fraction

Definition at line 969 of file tensor.cxx.

References IASSERT, CTF_int::tensor::name, CTF::World::rank, and CTF_int::tensor::wrld.

Referenced by matmul(), neural(), setup_unstructured(), sparse_checkpoint(), test_mis(), test_mis2(), and train_sps_vec_mat().

template<>
void CTF::Tensor< double >::fill_sp_random ( double  rmin,
double  rmax,
double  frac_sp 
)
inline

Definition at line 1018 of file tensor.cxx.

template<>
void CTF::Tensor< float >::fill_sp_random ( float  rmin,
float  rmax,
double  frac_sp 
)
inline

Definition at line 1023 of file tensor.cxx.

template<>
void CTF::Tensor< int >::fill_sp_random ( int  rmin,
int  rmax,
double  frac_sp 
)
inline

Definition at line 1028 of file tensor.cxx.

template<typename dtype = double>
void CTF::Tensor< dtype >::get_all_data ( int64_t *  npair,
dtype **  data,
bool  unpack = false 
)

collects the entire tensor data on each process (not memory scalable)

Parameters
[out]npairnumber of values in the tensor
[out]datapointer to the data of the entire tensor, should be released with delete []
[in]unpackif true any symmetric tensor is unpacked, otherwise only unique elements are read
template<typename dtype >
void CTF::Tensor< dtype >::get_local_data ( int64_t *  npair,
int64_t **  global_idx,
dtype **  data,
bool  nonzeros_only = false,
bool  unpack_sym = false 
) const

Gives the global indices and values associated with the local data.

Parameters
[out]npairnumber of local values
[out]global_idxindex within global tensor of each data value
[out]datapointer to local values in the order of the indices, should be released with delete []
[in]nonzeros_onlyif true, outputs all tensor elements, if false ignores those equivalent to the additive identity (zero)
[in]unpack_symif true, outputs all tensor elements, if false only those unique with respect to symmetry

Definition at line 159 of file tensor.cxx.

References CTF_int::alloc(), CTF_int::algstrct::alloc(), ctf.core::dtype, CTF::Tensor< dtype >::i(), IASSERT, CTF_int::PairIterator::k(), CTF_int::algstrct::pair_dealloc(), CTF_int::tensor::read_local(), CTF_int::tensor::read_local_nnz(), CTF_int::PairIterator::read_val(), CTF_int::tensor::sr, and CTF_int::SUCCESS.

Referenced by bitonic_sort(), DFT_matrix(), diag_ctr(), diag_sym(), fast_sym(), fft(), Integrals::fill_rand(), Amplitudes::fill_rand(), flatten_block_sparse_matrix(), fold_unfold(), force_integration(), gemm_4D(), main(), multi_tsr_sym(), permute_multiworld(), repack(), scalar(), scan(), setup_unstructured(), sparse_permuted_slice(), sptensor_sum(), strassen(), sy_times_ns(), test_dft(), test_dft_3D(), test_recursive_matmul(), test_subworld_gemm(), trace(), twiddle_matrix(), and weigh_4D().

template<typename dtype >
void CTF::Tensor< dtype >::get_local_pairs ( int64_t *  npair,
Pair< dtype > **  pairs,
bool  nonzeros_only = false,
bool  unpack_sym = false 
) const

gives the global indices and values associated with the local data

Parameters
[out]npairnumber of local values
[out]pairspointer to local key-value pairs, should be released with delete []
[in]nonzeros_onlyif true, outputs all tensor elements, if false ignores those equivalent to the additive identity (zero)
[in]unpack_symif true, outputs all tensor elements, if false only those unique with respect to symmetry

Definition at line 201 of file tensor.cxx.

References IASSERT, CTF_int::tensor::read_local(), CTF_int::tensor::read_local_nnz(), and CTF_int::SUCCESS.

Referenced by apsp(), mis2(), and CTF::write_sparse_to_file_base().

template<typename dtype >
dtype * CTF::Tensor< dtype >::get_mapped_data ( char const *  idx,
Idx_Partition const &  prl,
Idx_Partition const &  blk = Idx_Partition(),
bool  unpack = true 
)

returns local data of tensor with parallel distribution prl and local blocking blk

Parameters
[in]idxassignment of characters to each dim
[in]prlmesh processor topology with character labels
[in]blklocal blocking with processor labels
[in]unpackwhether to unpack from symmetric layout
Returns
local piece of data of tensor in this distribution, release via delete []

Definition at line 1129 of file tensor.cxx.

References ctf.core::dtype, and CTF_int::tensor::read().

template<typename dtype >
void CTF::Tensor< dtype >::get_max_abs ( int  n,
dtype *  data 
) const

obtains a small number of the biggest elements of the tensor in sorted order (e.g. eigenvalues)

Parameters
[in]nnumber of elements to collect
[in]dataoutput data (should be preallocated to size at least n)

WARNING: currently functional only for dtype=double

Definition at line 920 of file tensor.cxx.

References CTF_int::tensor::get_max_abs(), IASSERT, and CTF_int::SUCCESS.

template<typename dtype >
dtype * CTF::Tensor< dtype >::get_raw_data ( int64_t *  size) const

gives the raw current local data with padding included

Parameters
[out]sizeof local data chunk
Returns
pointer to local data

Definition at line 152 of file tensor.cxx.

References CTF_int::tensor::data, and ctf.core::dtype.

Referenced by CTF::Scalar< dtype >::get_val(), CTF::Scalar< dtype >::Scalar(), and CTF::Matrix< dtype >::svd().

template<typename dtype = double>
dtype CTF::Tensor< dtype >::norm1 ( )
inline

computes the entrywise 1-norm of the tensor

Definition at line 806 of file tensor.h.

References CTF::OP_SUMABS.

template<typename dtype >
void CTF::Tensor< dtype >::norm1 ( double &  nrm)

computes the entrywise 1-norm of the tensor

Definition at line 823 of file tensor.cxx.

References IASSERT, CTF_int::tensor::name, CTF::World::rank, and CTF_int::tensor::wrld.

template<typename dtype >
void CTF::Tensor< dtype >::norm2 ( double &  nrm)

computes the frobenius norm of the tensor

Definition at line 864 of file tensor.cxx.

References IASSERT, CTF_int::tensor::name, CTF::World::rank, and CTF_int::tensor::wrld.

template<typename dtype = double>
dtype CTF::Tensor< dtype >::norm_infty ( )
inline

finds the max absolute value element of the tensor

Definition at line 816 of file tensor.h.

References ctf.core::dtype, and CTF::OP_MAXABS.

template<typename dtype >
void CTF::Tensor< dtype >::norm_infty ( double &  nrm)

finds the max absolute value element of the tensor

Definition at line 894 of file tensor.cxx.

References IASSERT, CTF_int::tensor::name, CTF::World::rank, and CTF_int::tensor::wrld.

template<typename dtype >
Tensor< dtype > & CTF::Tensor< dtype >::operator= ( dtype  val)
template<typename dtype = double>
Typ_Idx_Tensor<dtype> CTF::Tensor< dtype >::operator[] ( char const *  idx_map)

associated an index map with the tensor for future operation

Parameters
[in]idx_mapindex assignment for this tensor

Referenced by CTF::Tensor< dtype >::Tensor().

template<typename dtype >
Sparse_Tensor< dtype > CTF::Tensor< dtype >::operator[] ( std::vector< int64_t >  indices)

gives handle to sparse index subset of tensors

Parameters
[in]indicesvector of indices to sparse tensor

Definition at line 1199 of file tensor.cxx.

template<typename dtype >
void CTF::Tensor< dtype >::permute ( dtype  beta,
CTF_int::tensor A,
int *const *  perms_A,
dtype  alpha 
)

Apply permutation to matrix, potentially extracting a slice B[i,j,...] = beta*B[...] + alpha*A[perms_A[0][i],perms_A[1][j],...].

Parameters
[in]betascaling factor for values of tensor B (this)
[in]Aspecification of operand tensor A must live on the same World or a subset of the World on which B lives
[in]perms_Aspecifies permutations for tensor A, e.g. A[perms_A[0][i],perms_A[1][j]] if a subarray NULL, no permutation applied to this index, if an entry is -1, the corresponding entries of the tensor are skipped (A must then be smaller than B)
[in]alphascaling factor for A tensor

Definition at line 429 of file tensor.cxx.

References IASSERT, CTF_int::tensor::permute(), and CTF_int::SUCCESS.

Referenced by ctf.core.tensor::__setitem__(), permute_multiworld(), and sparse_permuted_slice().

template<typename dtype >
void CTF::Tensor< dtype >::permute ( int *const *  perms_B,
dtype  beta,
CTF_int::tensor A,
dtype  alpha 
)

Apply permutation to matrix, potentially extracting a slice B[perms_B[0][i],perms_B[0][j],...] = beta*B[...] + alpha*A[i,j,...].

Parameters
[in]perms_Bspecifies permutations for tensor B, e.g. B[perms_B[0][i],perms_B[1][j]] if a subarray NULL, no permutation applied to this index, if an entry is -1, the corresponding entries of the tensor are skipped (A must then be smaller than B)
[in]betascaling factor for values of tensor B (this)
[in]Aspecification of operand tensor A must live on the same World or a superset of the World on which B lives
[in]alphascaling factor for A tensor

Definition at line 439 of file tensor.cxx.

References IASSERT, CTF_int::tensor::permute(), and CTF_int::SUCCESS.

Referenced by ctf.core.tensor::__setitem__().

template<typename dtype >
void CTF::Tensor< dtype >::print ( FILE *  fp,
dtype  cutoff 
) const

prints tensor data to file using process 0 (modify print(...) overload in set.h if you would like a different print format)

Parameters
[in]fpfile to print to e.g. stdout
[in]cutoffdo not print values of absolute value smaller than this

Definition at line 407 of file tensor.cxx.

References CTF_int::tensor::print().

Referenced by checkpoint(), fast_sym(), fast_sym_4D(), main(), scalar(), and smooth_jacobi().

template<typename dtype >
void CTF::Tensor< dtype >::print ( FILE *  fp = stdout) const

Definition at line 412 of file tensor.cxx.

References CTF_int::tensor::print().

template<typename dtype >
void CTF::Tensor< dtype >::prnt ( ) const

Definition at line 417 of file tensor.cxx.

References CTF_int::tensor::print().

template<typename dtype >
void CTF::Tensor< dtype >::profile_off ( )

turns off profiling for tensor

Definition at line 402 of file tensor.cxx.

References CTF_int::tensor::profile_off().

template<typename dtype >
void CTF::Tensor< dtype >::profile_on ( )

turns on profiling for tensor

Definition at line 397 of file tensor.cxx.

References CTF_int::tensor::profile_on().

template<typename dtype = double>
const dtype* CTF::Tensor< dtype >::raw_data ( int64_t *  size) const

gives a read-only copy of the raw current local data with padding included

Parameters
[out]sizeof local data chunk
Returns
pointer to read-only copy of local data
template<typename dtype >
void CTF::Tensor< dtype >::read ( int64_t  npair,
Pair< dtype > *  pairs 
)

Gives the values associated with any set of indices.

Parameters
[in]npairnumber of values to fetch
[in,out]pairsa prealloced pointer to key-value pairs

Definition at line 246 of file tensor.cxx.

References IASSERT, CTF_int::tensor::read(), and CTF_int::SUCCESS.

Referenced by ctf.core.tensor::__getitem__(), bench_redistribution(), CTF::Matrix< dtype >::read_mat(), and weigh_4D().

template<typename dtype >
void CTF::Tensor< dtype >::read ( int64_t  npair,
int64_t const *  global_idx,
dtype *  data 
)

Gives the values associated with any set of indices The sparse data is defined in coordinate format. The tensor index (i,j,k,l) of a tensor with edge lengths {m,n,p,q} is associated with the global index g via the formula g=i+j*m+k*m*n+l*m*n*p. The row index is first and the column index is second for matrices, which means they are column major.

Parameters
[in]npairnumber of values to fetch
[in]global_idxindex within global tensor of each value to fetch
[in,out]dataa prealloced pointer to the data with the specified indices

Definition at line 226 of file tensor.cxx.

References CTF::Pair< dtype >::d, CTF::Tensor< dtype >::i(), IASSERT, CTF::Pair< dtype >::k, CTF_int::algstrct::pair_alloc(), CTF_int::algstrct::pair_dealloc(), CTF_int::tensor::read(), CTF_int::tensor::sr, and CTF_int::SUCCESS.

Referenced by ctf.core.tensor::__getitem__().

template<typename dtype >
void CTF::Tensor< dtype >::read ( int64_t  npair,
dtype  alpha,
dtype  beta,
Pair< dtype > *  pairs 
)

sparse read with accumulation: pairs[i].d = alpha*A[pairs[i].k]+beta*pairs[i].d

Parameters
[in]npairnumber of values to read into tensor
[in]alphascaling factor on read data
[in]betascaling factor on value in initial pairs vector
[in]pairskey-value pairs to which tensor values should be accumulated

Definition at line 357 of file tensor.cxx.

References IASSERT, CTF_int::tensor::read(), and CTF_int::SUCCESS.

Referenced by ctf.core.tensor::__getitem__().

template<typename dtype >
void CTF::Tensor< dtype >::read ( int64_t  npair,
dtype  alpha,
dtype  beta,
int64_t const *  global_idx,
dtype *  data 
)

sparse read with accumulation: data[i] = alpha*A[global_idx[i]]+beta*data[i]

Parameters
[in]npairnumber of values to read into tensor
[in]alphascaling factor on read data
[in]betascaling factor on value in initial values vector
[in]global_idxglobal index within tensor of value to add
[in]datavalues to which tensor values should be accumulated

Definition at line 336 of file tensor.cxx.

References CTF::Pair< dtype >::d, CTF::Tensor< dtype >::i(), IASSERT, CTF::Pair< dtype >::k, CTF_int::algstrct::pair_alloc(), CTF_int::algstrct::pair_dealloc(), CTF_int::tensor::read(), CTF_int::tensor::sr, and CTF_int::SUCCESS.

Referenced by ctf.core.tensor::__getitem__().

template<typename dtype >
void CTF::Tensor< dtype >::read_all ( int64_t *  npair,
dtype **  data,
bool  unpack = false 
)

collects the entire tensor data on each process (not memory scalable)

Parameters
[out]npairnumber of values in the tensor
[out]datapointer to the data of the entire tensor, should be released with free
[in]unpackif true any symmetric tensor is unpacked, otherwise only unique elements are read

Definition at line 377 of file tensor.cxx.

References CTF_int::tensor::allread(), IASSERT, and CTF_int::SUCCESS.

Referenced by bitonic(), bivar_function(), bivar_transform(), endomorphism(), CTF::Matrix< dtype >::print_matrix(), readall_test(), scan_test(), ctf.core.tensor::to_nparray(), and univar_function().

template<typename dtype >
int64_t CTF::Tensor< dtype >::read_all ( dtype *  data,
bool  unpack = false 
)

collects the entire tensor data on each process (not memory scalable)

Parameters
[in,out]datapreallocated pointer to the data of the entire tensor
[in]unpackif true any symmetric tensor is unpacked, otherwise only unique elements are read

Definition at line 384 of file tensor.cxx.

References CTF_int::tensor::allread(), IASSERT, and CTF_int::SUCCESS.

Referenced by ctf.core.tensor::to_nparray().

template<typename dtype >
void CTF::Tensor< dtype >::read_local ( int64_t *  npair,
int64_t **  global_idx,
dtype **  data,
bool  unpack_sym = false 
) const

Using get_local_data(), which returns an array that must be freed with delete [], is more efficient, this version exists for backwards compatibility. gives the global indices and values associated with the local data WARNING-1: for sparse tensors this includes the zeros to maintain consistency with the behavior for dense tensors, use get_local_pairs to get only nonzeros.

Parameters
[out]npairnumber of local values
[out]global_idxindex within global tensor of each data value
[out]datapointer to local values in the order of the indices, should be released with free
[in]unpack_symif true, outputs all tensor elements, if false only those unique with respect to symmetry

Definition at line 182 of file tensor.cxx.

References CTF_int::alloc(), ctf.core::dtype, CTF::Tensor< dtype >::i(), IASSERT, CTF_int::PairIterator::k(), CTF_int::algstrct::pair_dealloc(), CTF_int::tensor::read_local(), CTF_int::PairIterator::read_val(), CTF_int::tensor::sr, and CTF_int::SUCCESS.

Referenced by fast_diagram(), fast_tensor_ctr(), CTF::fold_unfold(), fold_unfold(), and get_rand_as_tsr().

template<typename dtype >
void CTF::Tensor< dtype >::read_local ( int64_t *  npair,
Pair< dtype > **  pairs,
bool  unpack_sym = false 
) const

gives the global indices and values associated with the local data WARNING: for sparse tensors this includes the zeros to maintain consistency with the behavior for dense tensors, use get_lcoal_pairs to get only nonzeros

Parameters
[out]npairnumber of local values
[out]pairspointer to local key-value pairs, should be released with free
[in]unpack_symif true, outputs all tensor elements, if false only those unique with respect to symmetry

Definition at line 216 of file tensor.cxx.

References IASSERT, CTF_int::tensor::read_local(), and CTF_int::SUCCESS.

template<>
void CTF::Tensor< int >::read_sparse_from_file ( const char *  fpath,
bool  with_vals 
)
inline

Definition at line 485 of file tensor.cxx.

template<>
void CTF::Tensor< double >::read_sparse_from_file ( const char *  fpath,
bool  with_vals 
)
inline

Definition at line 490 of file tensor.cxx.

template<>
void CTF::Tensor< float >::read_sparse_from_file ( const char *  fpath,
bool  with_vals 
)
inline

Definition at line 495 of file tensor.cxx.

template<>
void CTF::Tensor< int64_t >::read_sparse_from_file ( const char *  fpath,
bool  with_vals 
)
inline

Definition at line 500 of file tensor.cxx.

template<typename dtype = double>
void CTF::Tensor< dtype >::read_sparse_from_file ( const char *  fpath,
bool  with_vals = true 
)

read sparse tensor from file, entries of tensor must be stored one per line, as i_1 ... i_order v, to create entry T[i_1, ..., i_order] = v or as i_1 ... i_order, to create entry T[i_1, ..., i_order] = mulid

Parameters
[in]fpathstring of file name to read from
[in]with_valswhether vs are provided in file

Referenced by sparse_checkpoint().

template<typename dtype >
dtype CTF::Tensor< dtype >::reduce ( OP  op)
template<typename dtype = double>
void CTF::Tensor< dtype >::scale ( dtype  alpha,
char const *  idx_A 
)

scales A[idx_A] = alpha*A[idx_A]

Parameters
[in]alphaA scaling factor
[in]idx_Aindices of A (this tensor), e.g. "ij" -> A_{ij}

Referenced by CTF::Tensor< dtype >::fill_sp_random().

template<typename dtype = double>
void CTF::Tensor< dtype >::scale ( dtype  alpha,
char const *  idx_A,
Endomorphism< dtype >  fseq 
)

scales A[idx_A] = fseq(alpha*A[idx_A])

Parameters
[in]alphaA scaling factor
[in]idx_Aindices of A (this tensor), e.g. "ij" -> A_{ij}
[in]fsequser defined function
template<typename dtype >
void CTF::Tensor< dtype >::set_name ( char const *  name)

sets tensor name

Parameters
[in]namenew for tensor

Definition at line 392 of file tensor.cxx.

References CTF_int::tensor::set_name().

template<typename dtype >
Tensor< dtype > CTF::Tensor< dtype >::slice ( int const *  offsets,
int const *  ends 
) const

cuts out a slice (block) of this tensor A[offsets,ends) result will always be fully nonsymmetric

Parameters
[in]offsetsbottom left corner of block
[in]endstop right corner of block
Returns
new tensor corresponding to requested slice

Definition at line 643 of file tensor.cxx.

References CTF_int::tensor::wrld.

Referenced by btwn_cnt_fast(), flatten_block_sparse_matrix(), CTF::Matrix< dtype >::qr(), recursive_matmul(), CTF::Tensor< dtype >::slice(), strassen(), and CTF::Matrix< dtype >::svd().

template<typename dtype >
Tensor< dtype > CTF::Tensor< dtype >::slice ( int64_t  corner_off,
int64_t  corner_end 
) const

cuts out a slice (block) of this tensor with corners specified by global index result will always be fully nonsymmetric

Parameters
[in]corner_offtop left corner of block
[in]corner_endbottom right corner of block
Returns
new tensor corresponding to requested slice

Definition at line 650 of file tensor.cxx.

References CTF::Tensor< dtype >::slice(), and CTF_int::tensor::wrld.

template<typename dtype >
Tensor< dtype > CTF::Tensor< dtype >::slice ( int const *  offsets,
int const *  ends,
World oworld 
) const

cuts out a slice (block) of this tensor A[offsets,ends) result will always be fully nonsymmetric

Parameters
[in]offsetsbottom left corner of block
[in]endstop right corner of block
[in]oworldthe world in which the new tensor should be defined
Returns
new tensor corresponding to requested slice which lives on oworld

Definition at line 657 of file tensor.cxx.

References CTF_int::algstrct::addid(), CTF_int::alloc(), CTF_int::cdealloc(), ctf.core::dtype, CTF::Tensor< dtype >::i(), IASSERT, CTF_int::tensor::lens, CTF_int::algstrct::mulid(), NS, CTF_int::tensor::order, CTF::Tensor< dtype >::slice(), CTF_int::tensor::sr, and CTF_int::tensor::sym.

template<typename dtype >
Tensor< dtype > CTF::Tensor< dtype >::slice ( int64_t  corner_off,
int64_t  corner_end,
World oworld 
) const

cuts out a slice (block) of this tensor with corners specified by global index result will always be fully nonsymmetric

Parameters
[in]corner_offtop left corner of block
[in]corner_endbottom right corner of block
[in]oworldthe world in which the new tensor should be defined
Returns
new tensor corresponding to requested slice which lives on oworld

Definition at line 699 of file tensor.cxx.

References CTF_int::cdealloc(), CTF_int::cvrt_idx(), CTF::Tensor< dtype >::i(), CTF_int::tensor::lens, CTF_int::tensor::order, CTF::Tensor< dtype >::slice(), and ctf.core::tsr.

template<typename dtype >
void CTF::Tensor< dtype >::slice ( int const *  offsets,
int const *  ends,
dtype  beta,
CTF_int::tensor const &  A,
int const *  offsets_A,
int const *  ends_A,
dtype  alpha 
)

adds to a slice (block) of this tensor = B B[offsets,ends)=beta*B[offsets,ends) + alpha*A[offsets_A,ends_A)

Parameters
[in]offsetsbottom left corner of block
[in]endstop right corner of block
[in]betascaling factor of this tensor
[in]Atensor who owns pure-operand slice
[in]offsets_Abottom left corner of block of A
[in]ends_Atop right corner of block of A
[in]alphascaling factor of tensor A

Definition at line 586 of file tensor.cxx.

References CTF::World::comm, IASSERT, CTF_int::tensor::slice(), and CTF_int::tensor::wrld.

template<typename dtype >
void CTF::Tensor< dtype >::slice ( int64_t  corner_off,
int64_t  corner_end,
dtype  beta,
CTF_int::tensor const &  A,
int64_t  corner_off_A,
int64_t  corner_end_A,
dtype  alpha 
)

adds to a slice (block) of this tensor = B B[offsets,ends)=beta*B[offsets,ends) + alpha*A[offsets_A,ends_A)

Parameters
[in]corner_offtop left corner of block
[in]corner_endbottom right corner of block
[in]betascaling factor of this tensor
[in]Atensor who owns pure-operand slice
[in]corner_off_Atop left corner of block of A
[in]corner_end_Abottom right corner of block of A
[in]alphascaling factor of tensor A

Definition at line 614 of file tensor.cxx.

References CTF_int::cdealloc(), CTF_int::cvrt_idx(), CTF::Tensor< dtype >::i(), CTF_int::tensor::lens, CTF_int::tensor::order, and CTF_int::tensor::slice().

template<typename dtype >
void CTF::Tensor< dtype >::sparsify ( )

reduce tensor to sparse format, storing only nonzero data, or data above a specified threshold. makes dense tensors sparse. cleans sparse tensors of any 'computed' zeros.

Definition at line 449 of file tensor.cxx.

References IASSERT, CTF_int::tensor::sparsify(), and CTF_int::SUCCESS.

Referenced by apsp(), btwn_cnt_fast(), jacobi(), mis(), mis2(), sparse_mp3(), spmv(), and train_off_vec_mat().

template<typename dtype >
void CTF::Tensor< dtype >::sparsify ( dtype  threshold,
bool  take_abs = true 
)

reduce tensor to sparse format, storing only nonzero data, or data above a specified threshold. makes dense tensors sparse. cleans sparse tensors of any small values

Parameters
[in]thresholdall values smaller or equal to than this one will be removed/not stored (by default is NULL, meaning only zeros are removed, so same as threshold=additive identity)
[in]take_abswhether to take absolute value when comparing to threshold

Definition at line 455 of file tensor.cxx.

References IASSERT, CTF_int::tensor::sparsify(), and CTF_int::SUCCESS.

template<typename dtype >
void CTF::Tensor< dtype >::sparsify ( std::function< bool(dtype)>  filter)

sparsifies tensor keeping only values v such that filter(v) = true

Parameters
[in]filterboolean function to apply to values to determine whether to keep them

Definition at line 461 of file tensor.cxx.

References ctf.core::dtype, IASSERT, CTF_int::tensor::sparsify(), and CTF_int::SUCCESS.

template<typename dtype = double>
void CTF::Tensor< dtype >::sum ( dtype  alpha,
CTF_int::tensor A,
char const *  idx_A,
dtype  beta,
char const *  idx_B 
)

sums B[idx_B] = beta*B[idx_B] + alpha*A[idx_A]

Parameters
[in]alphaA scaling factor
[in]Afirst operand tensor
[in]idx_Aindices of A in sum, e.g. "ij" -> A_{ij}
[in]betaB scaling factor
[in]idx_Bindices of B (this tensor), e.g. "ij" -> B_{ij}

Referenced by CTF::Tensor< dtype >::fill_sp_random(), and CTF::Tensor< dtype >::operator=().

template<typename dtype = double>
void CTF::Tensor< dtype >::sum ( dtype  alpha,
CTF_int::tensor A,
char const *  idx_A,
dtype  beta,
char const *  idx_B,
Univar_Function< dtype >  fseq 
)

sums B[idx_B] = beta*B[idx_B] + fseq(alpha*A[idx_A])

Parameters
[in]alphaA scaling factor
[in]Afirst operand tensor
[in]idx_Aindices of A in sum, e.g. "ij" -> A_{ij}
[in]betaB scaling factor
[in]idx_Bindices of B (this tensor), e.g. "ij" -> B_{ij}
[in]fseqsequential operation to execute, default is multiply-add
template<typename dtype >
void CTF::Tensor< dtype >::write ( int64_t  npair,
int64_t const *  global_idx,
dtype const *  data 
)

writes in values associated with any set of indices The sparse data is defined in coordinate format. The tensor index (i,j,k,l) of a tensor with edge lengths {m,n,p,q} is associated with the global index g via the formula g=i+j*m+k*m*n+l*m*n*p. The row index is first and the column index is second for matrices, which means they are column major.

Parameters
[in]npairnumber of values to write into tensor
[in]global_idxglobal index within tensor of value to write
[in]datavalues to write to the indices

Definition at line 264 of file tensor.cxx.

References CTF_int::algstrct::addid(), CTF::Pair< dtype >::d, CTF::Tensor< dtype >::i(), IASSERT, CTF::Pair< dtype >::k, CTF_int::algstrct::mulid(), CTF_int::algstrct::pair_alloc(), CTF_int::algstrct::pair_dealloc(), CTF_int::tensor::sr, CTF_int::SUCCESS, and CTF_int::tensor::write().

Referenced by ctf.core.tensor::__setitem__(), bitonic_sort(), block_sparse(), DFT_matrix(), diag_ctr(), diag_sym(), fast_diagram(), fast_sym(), fast_tensor_ctr(), fft(), Integrals::fill_rand(), Amplitudes::fill_rand(), CTF::fill_sp_random_base(), CTF::fold_unfold(), fold_unfold(), force_integration(), force_integration_sparse(), ctf.core.tensor::from_nparray(), gemm_4D(), get_rand_as_tsr(), main(), mis2(), multi_tsr_sym(), neural(), permute_multiworld(), CTF::read_sparse_from_file_base(), repack(), scalar(), scan(), setup_3d_Poisson(), setup_unstructured(), sparse_permuted_slice(), sptensor_sum(), sssp(), strassen(), sy_times_ns(), test_dft(), test_dft_3D(), test_recursive_matmul(), test_subworld_gemm(), trace(), twiddle_matrix(), weigh_4D(), and CTF::Matrix< dtype >::write_mat().

template<typename dtype >
void CTF::Tensor< dtype >::write ( int64_t  npair,
Pair< dtype > const *  pairs 
)

writes in values associated with any set of indices

Parameters
[in]npairnumber of values to write into tensor
[in]pairskey-value pairs to write to the tensor

Definition at line 286 of file tensor.cxx.

References CTF_int::algstrct::addid(), IASSERT, CTF_int::algstrct::mulid(), CTF_int::tensor::sr, CTF_int::SUCCESS, and CTF_int::tensor::write().

Referenced by ctf.core.tensor::__setitem__(), and ctf.core.tensor::from_nparray().

template<typename dtype >
void CTF::Tensor< dtype >::write ( int64_t  npair,
dtype  alpha,
dtype  beta,
int64_t const *  global_idx,
dtype const *  data 
)

sparse add: A[global_idx[i]] = beta*A[global_idx[i]]+alpha*data[i]

Parameters
[in]npairnumber of values to write into tensor
[in]alphascaling factor on value to add
[in]betascaling factor on original data
[in]global_idxglobal index within tensor of value to add
[in]datavalues to add to the tensor

Definition at line 298 of file tensor.cxx.

References CTF::Pair< dtype >::d, CTF::Tensor< dtype >::i(), IASSERT, CTF::Pair< dtype >::k, CTF_int::algstrct::pair_alloc(), CTF_int::algstrct::pair_dealloc(), CTF_int::tensor::sr, CTF_int::SUCCESS, and CTF_int::tensor::write().

Referenced by ctf.core.tensor::__setitem__(), and ctf.core.tensor::from_nparray().

template<typename dtype >
void CTF::Tensor< dtype >::write ( int64_t  npair,
dtype  alpha,
dtype  beta,
Pair< dtype > const *  pairs 
)

sparse add: A[pairs[i].k] = alpha*A[pairs[i].k]+beta*pairs[i].d

Parameters
[in]npairnumber of values to write into tensor
[in]alphascaling factor on value to add
[in]betascaling factor on original data
[in]pairskey-value pairs to add to the tensor

Definition at line 324 of file tensor.cxx.

References IASSERT, CTF_int::SUCCESS, and CTF_int::tensor::write().

Referenced by ctf.core.tensor::__setitem__(), and ctf.core.tensor::from_nparray().

template<typename dtype = double>
void CTF::Tensor< dtype >::write ( char const *  idx,
dtype const *  data,
Idx_Partition const &  prl,
Idx_Partition const &  blk = Idx_Partition(),
bool  unpack = true 
)

writes data to tensor from parallel distribution prl and local blocking blk

Parameters
[in]idxassignment of characters to each dim
[in]datato write from this distribution
[in]prlmesh processor topology with character labels
[in]blklocal blocking with processor labels
[in]unpackwhether written data is unpacked from symmetric layout

Referenced by ctf.core.tensor::__setitem__(), and ctf.core.tensor::from_nparray().

template<>
void CTF::Tensor< int >::write_sparse_to_file ( const char *  fpath,
bool  with_vals 
)
inline

Definition at line 519 of file tensor.cxx.

template<>
void CTF::Tensor< double >::write_sparse_to_file ( const char *  fpath,
bool  with_vals 
)
inline

Definition at line 524 of file tensor.cxx.

template<>
void CTF::Tensor< float >::write_sparse_to_file ( const char *  fpath,
bool  with_vals 
)
inline

Definition at line 529 of file tensor.cxx.

template<>
void CTF::Tensor< int64_t >::write_sparse_to_file ( const char *  fpath,
bool  with_vals 
)
inline

Definition at line 534 of file tensor.cxx.

template<typename dtype = double>
void CTF::Tensor< dtype >::write_sparse_to_file ( const char *  fpath,
bool  with_vals = true 
)

write sparse tensor to file, entries of tensor will be stored one per line, as i_1 ... i_order v, corresponding to entry T[i_1, ..., i_order] = v or as i_1 ... i_order if with_vals =false

Parameters
[in]fpathstring of file name to read from
[in]with_valswhether vs should be written to file

Referenced by sparse_checkpoint().


The documentation for this class was generated from the following files: