Cyclops Tensor Framework
parallel arithmetic on multidimensional arrays
CTF_int::tensor Class Reference

internal distributed tensor class More...

#include <untyped_tensor.h>

Inheritance diagram for CTF_int::tensor:
Collaboration diagram for CTF_int::tensor:

Public Member Functions

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...
 

Data Fields

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

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

internal distributed tensor class

Definition at line 20 of file untyped_tensor.h.

Constructor & Destructor Documentation

CTF_int::tensor::tensor ( )
CTF_int::tensor::~tensor ( )

class free self

Definition at line 78 of file untyped_tensor.cxx.

CTF_int::tensor::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)

Parameters
[in]srdefines the tensor arithmetic for this tensor
[in]ordernumber of dimensions of tensor
[in]edge_lenedge lengths of tensor
[in]symsymmetries of tensor (e.g. symmetric matrix -> sym={SY, NS})
[in]wrlda distributed context for the tensor to live in
[in]alloc_datawhether to allocate and set data to zero immediately
[in]namean optionary name for the tensor
[in]profileset to 1 to profile contractions involving this tensor
[in]is_sparseset to 1 to store only nontrivial tensor elements

Definition at line 82 of file untyped_tensor.cxx.

CTF_int::tensor::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)

Parameters
[in]srdefines the tensor arithmetic for this tensor
[in]ordernumber of dimensions of tensor
[in]is_sparsewhether to make tensor sparse
[in]edge_lenedge lengths of tensor
[in]symsymmetries of tensor (e.g. symmetric matrix -> sym={SY, NS})
[in]wrlda distributed context 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

Definition at line 94 of file untyped_tensor.cxx.

References CTF_int::algstrct::addid(), CTF_int::alloc(), CTF_int::algstrct::alloc(), CTF_int::accumulatable::el_size, and CTF_int::algstrct::set().

CTF_int::tensor::tensor ( tensor const *  other,
bool  copy = 1,
bool  alloc_data = 1 
)

creates tensor copy, unfolds other if other is folded

Parameters
[in]othertensor to copy
[in]copywhether to copy mapping and data
[in]alloc_datawhether th allocate data

Definition at line 136 of file untyped_tensor.cxx.

References CTF_int::alloc(), CTF_int::cdealloc(), DPRINTF, has_zero_edge_len, is_sparse, lens, name, order, profile, CTF::World::rank, sr, sym, and wrld.

CTF_int::tensor::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!

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

Definition at line 157 of file untyped_tensor.cxx.

References CTF_int::alloc(), CTF_int::cdealloc(), DPRINTF, has_zero_edge_len, is_sparse, lens, name, NS, order, profile, CTF::World::rank, sr, CTF_int::summation::sum_tensors(), sym, and wrld.

Member Function Documentation

void CTF_int::tensor::add_from_subworld ( tensor tsr_sub,
char const *  alpha,
char const *  beta 
)

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

Parameters
[in]tsr_subid of tensor on a subcomm of this CTF inst
[in]alphascaling factor for this tensor
[in]betascaling factor for tensor tsr

Definition at line 1037 of file untyped_tensor.cxx.

References CTF_int::cyclic_reshuffle(), order, and orient_subworld().

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

void CTF_int::tensor::add_to_subworld ( tensor tsr_sub,
char const *  alpha,
char const *  beta 
)

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

Parameters
[in]tsr_subtensor on a subcomm of this world
[in]alphascaling factor for this tensor
[in]betascaling factor for tensor tsr

Definition at line 991 of file untyped_tensor.cxx.

References ASSERT, CTF_int::cyclic_reshuffle(), data, CTF_int::algstrct::mdtype(), order, orient_subworld(), CTF_int::distribution::size, slice(), and sr.

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

void CTF_int::tensor::addinv ( )
int CTF_int::tensor::align ( tensor const *  B)

align mapping of thisa tensor to that of B

brief copy A into this (B). Realloc if necessary param[in] A tensor to copy

Parameters
[in]Btensor handle of B

Definition at line 1749 of file untyped_tensor.cxx.

References ASSERT, CTF_int::comp_dim_map(), CTF_int::copy_mapping(), edge_map, is_folded, order, CTF_int::SUCCESS, topo, and wrld.

Referenced by CTF::Tensor< dtype >::align(), and compare().

int CTF_int::tensor::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.

Parameters
[out]num_pairnumber of values read
[in,out]all_datavalues read (allocated by library)
[in]unpackif true any symmetric tensor is unpacked, otherwise only unique elements are read

Definition at line 1724 of file untyped_tensor.cxx.

References CTF_int::PairIterator::ptr, CTF_int::PairIterator::read_val(), and CTF_int::SUCCESS.

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

int CTF_int::tensor::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.

Parameters
[out]num_pairnumber of values read
[in,out]all_datapreallocated mapped_data values read
[in]unpackif true any symmetric tensor is unpacked, otherwise only unique elements are read

Definition at line 1738 of file untyped_tensor.cxx.

References CTF_int::PairIterator::read_val(), and CTF_int::SUCCESS.

int64_t CTF_int::tensor::calc_npe ( ) const

calculate the number of processes this tensor is distributed over return number of processes owning a block of the tensor

Definition at line 435 of file untyped_tensor.cxx.

References CTF_int::mapping::child, CTF_int::mapping::has_child, CTF_int::mapping::np, CTF_int::PHYSICAL_MAP, and CTF_int::mapping::type.

Referenced by CTF_int::get_len_ordering().

int * CTF_int::tensor::calc_phase ( ) const

compute the cyclic phase of each tensor dimension

Returns
int * of cyclic phases

Definition at line 394 of file untyped_tensor.cxx.

References CTF_int::alloc_ptr(), and CTF_int::mapping::calc_phase().

int CTF_int::tensor::calc_tot_phase ( ) const

calculate the total number of blocks of the tensor

Returns
int total phase factor

Definition at line 406 of file untyped_tensor.cxx.

References CTF_int::cdealloc().

void CTF_int::tensor::clear_mapping ( )

zeros out mapping

Definition at line 2189 of file untyped_tensor.cxx.

References CTF_int::mapping::clear().

Referenced by CTF_int::desymmetrize(), CTF_int::scaling::execute(), CTF_int::get_len_ordering(), and read().

void CTF_int::tensor::compare ( const tensor A,
FILE *  fp,
char const *  cutoff 
)

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

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

Definition at line 1930 of file untyped_tensor.cxx.

References CTF_int::algstrct::abs, align(), CTF_int::alloc_ptr(), CTF_int::cdealloc(), CTF::World::cdt, CTF_int::CommData::cm, CTF_int::accumulatable::el_size, CTF_int::algstrct::isequal(), lens, CTF_int::algstrct::min(), CTF_int::CommData::np, order, CTF_int::algstrct::pair_alloc(), CTF_int::algstrct::pair_size(), CTF_int::algstrct::print(), print_map(), CTF_int::CommData::rank, read_local(), CTF_int::PairIterator::sort(), sr, and wrld.

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

template<typename dtype >
template void CTF_int::tensor::compare_elementwise< bool > ( 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)

Parameters
[in]Afirst operand
[in]Bsecond operand

Definition at line 27 of file untyped_tensor_tmpl.h.

References ctf.core::a, ctf.core::b, ctf.core::dtype, CTF_int::accumulatable::el_size, operator[](), order, and sr.

template<typename dtype_A , typename dtype_B >
void CTF_int::tensor::conv_type ( tensor B)

convert this tensor from dtype_A to dtype_B and store the result in B (primarily needed for python interface)

Parameters
[in]Boutput tensor

Definition at line 7 of file untyped_tensor_tmpl.h.

References ctf.core::a, and order.

void CTF_int::tensor::copy_tensor_data ( tensor const *  other)
protected
void CTF_int::tensor::deregister_size ( )

deregister buffer allocation for this tensor

Definition at line 2889 of file untyped_tensor.cxx.

References CTF_int::inc_tot_mem_used().

void CTF_int::tensor::despmatricize ( int  nrow_idx,
bool  csr 
)

transposes back local data from sparse matrix format to key-value pair format

Parameters
[in]nrow_idxnumber of indices to fold into column
[in]csrwhether to go from csr (1) or coo (0) layout

Definition at line 2790 of file untyped_tensor.cxx.

References CTF_int::COO_Matrix::all_data, ASSERT, CTF_int::cdealloc(), CTF_int::COO_Matrix::get_data(), CTF_int::COO_Matrix::nnz(), CTF_int::CSR_Matrix::nnz(), TAU_FSTART, and TAU_FSTOP.

double CTF_int::tensor::est_time_unfold ( )

estimate cost of potential transpose involved in undoing the folding of a local tensor block

Returns
estimated time for transpose

Definition at line 2075 of file untyped_tensor.cxx.

References CTF_int::alloc_ptr(), CTF_int::calc_dim(), CTF_int::cdealloc(), CTF_int::est_time_transp(), NS, and CTF_int::sy_packed_size().

template<typename dtype_A , typename dtype_B >
template void CTF_int::tensor::exp_helper< int32_t, float > ( tensor A)

function store the e**value in tensor A into this (primarily needed for python interface)

Definition at line 18 of file untyped_tensor_tmpl.h.

References ctf.core::a, ctf.core::exp(), operator[](), and order.

int CTF_int::tensor::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

Parameters
[in]idx_mapindex map of tensor for this operation
[in]rwif 1 this writes to the diagonal, if 0 it reads the diagonal
[in,out]new_tsrif rw=1 this will be output as new tensor if rw=0 this should be input as the tensor of the extracted diagonal
[out]idx_map_newif rw=1 this will be the new index map

Definition at line 2456 of file untyped_tensor.cxx.

References CTF_int::alloc_ptr(), CTF_int::cdealloc(), CTF_int::summation::execute(), CTF_int::PairIterator::k(), CTF_int::NEGATIVE, NS, read_local_nnz(), CTF_int::SUCCESS, ctf.core::sum(), write(), and CTF_int::PairIterator::write_val().

Referenced by CTF_int::summation::estimate_time(), CTF_int::get_len_ordering(), and CTF_int::summation::sum_tensors().

void CTF_int::tensor::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

Parameters
[in]nfoldnumber of global indices we are folding
[in]fold_idxwhich global indices we are folding
[in]idx_maphow this tensor indices map to the global indices
[out]all_fdimnumber of dimensions including unfolded dimensions
[out]all_flenedge lengths including unfolded dimensions

Definition at line 2101 of file untyped_tensor.cxx.

References CTF_int::alloc_ptr(), CTF_int::calc_dim(), CTF_int::cdealloc(), NS, and CTF_int::sy_packed_size().

void CTF_int::tensor::free_self ( )
int CTF_int::tensor::get_max_abs ( int  n,
char *  data 
) const

obtains the largest n elements (in absolute value) of the tensor

Parameters
[in]nnumber of elements to fill
[in,out]datapreallocated array of size at least n, in which to put the elements

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

char const * CTF_int::tensor::get_name ( ) const

get the tensor name

Returns
tensor name

Definition at line 703 of file untyped_tensor.cxx.

void CTF_int::tensor::get_raw_data ( char **  data,
int64_t *  size 
) const

get raw data pointer without copy WARNING: includes padding

Parameters
[out]dataraw local data in char * format
[out]sizenumber of elements in data

Definition at line 715 of file untyped_tensor.cxx.

Referenced by CTF::Idx_Tensor::execute().

int64_t CTF_int::tensor::get_redist_mem ( distribution const &  old_dist,
double  nnz_frac 
)
int64_t CTF_int::tensor::get_tot_size ( bool  packed = false)

get number of elements in whole tensor

Parameters
[in]packedif false (default) ignore symmetry
Returns
number of elements (including zeros)

Definition at line 1712 of file untyped_tensor.cxx.

References CTF_int::packed_size().

void CTF_int::tensor::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 
)
protected

initializes tensor data

Parameters
[in]srdefines the tensor arithmetic for this tensor
[in]ordernumber of dimensions of tensor
[in]edge_lenedge lengths of tensor
[in]symsymmetries of tensor (e.g. symmetric matrix -> sym={SY, NS})
[in]wrlda distributed context for the tensor to live in
[in]alloc_dataset to 1 if tensor should be mapped and data buffer allocated
[in]namean optionary name for the tensor
[in]profileset to 1 to profile contractions involving this tensor
[in]is_sparseset to 1 to store only nontrivial tensor elements

Definition at line 306 of file untyped_tensor.cxx.

References CTF_int::alloc(), CTF_int::alloc_ptr(), ASSERT, CTF_int::algstrct::clone(), DPRINTF, CTF_int::NOT_MAPPED, CTF_int::SUCCESS, TAU_FSTART, TAU_FSTOP, and CTF_int::mapping::type.

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

template<typename dtype >
template void CTF_int::tensor::larger_equal_than< bool > ( 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)

Parameters
[in]Afirst operand
[in]Bsecond operand

Definition at line 87 of file untyped_tensor_tmpl.h.

References ctf.core::a, ctf.core::b, ctf.core::dtype, CTF_int::accumulatable::el_size, operator[](), order, and sr.

template<typename dtype >
template void CTF_int::tensor::larger_than< bool > ( 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)

Parameters
[in]Afirst operand
[in]Bsecond operand

Definition at line 75 of file untyped_tensor_tmpl.h.

References ctf.core::a, ctf.core::b, ctf.core::dtype, CTF_int::accumulatable::el_size, operator[](), order, and sr.

void CTF_int::tensor::leave_home_with_buffer ( )

degister home buffer

Definition at line 2863 of file untyped_tensor.cxx.

References DPRINTF.

int CTF_int::tensor::map_tensor_rem ( int  num_phys_dims,
CommData phys_comm,
int  fill = 0 
)

map the remainder of a tensor

Parameters
[in]num_phys_dimsnumber of physical processor grid dimensions
[in]phys_commdimensional communicators
[in]fillwhether to map everything

Definition at line 2405 of file untyped_tensor.cxx.

References CTF_int::alloc_ptr(), CTF_int::cdealloc(), CTF_int::mapping::cdt, CTF_int::mapping::child, CTF_int::mapping::has_child, CTF_int::map_tensor(), CTF_int::NOT_MAPPED, CTF_int::PHYSICAL_MAP, and CTF_int::mapping::type.

Referenced by CTF_int::scaling::execute(), and CTF_int::get_len_ordering().

template<typename dtype >
template void CTF_int::tensor::not_equals< bool > ( 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)

Parameters
[in]Afirst operand
[in]Bsecond operand

Definition at line 39 of file untyped_tensor_tmpl.h.

References ctf.core::a, ctf.core::b, ctf.core::dtype, CTF_int::accumulatable::el_size, operator[](), order, and sr.

Idx_Tensor CTF_int::tensor::operator[] ( char const *  idx_map)

associated an index map with the tensor for future operation

Parameters
[in]idx_mapindex assignment for this tensor

Definition at line 31 of file untyped_tensor.cxx.

Referenced by compare_elementwise(), exp_helper(), larger_equal_than(), larger_than(), not_equals(), smaller_equal_than(), smaller_than(), and true_divide().

void CTF_int::tensor::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

Parameters
[in]greater_worldcomm with respect to which the data needs to be ordered
[out]bw_mirror_rankprocessor rank in greater_world from which data is received
[out]fw_mirror_rankprocessor rank in greater_world to which data is sent
[out]odstdistribution mapping of data on output defined on oriented subworld
[out]sub_buffer_allocated buffer of received data on oriented subworld

Definition at line 777 of file untyped_tensor.cxx.

References CTF_int::alloc(), CTF_int::CommData::allred(), ASSERT, CTF_int::CommData::bcast(), CTF_int::cdealloc(), CTF::World::cdt, CTF::World::comm, CTF_int::get_distribution_size(), CTF::World::rank, CTF_int::distribution::serialize(), and CTF_int::distribution::size.

Referenced by add_from_subworld(), and add_to_subworld().

int CTF_int::tensor::permute ( tensor A,
int *const *  permutation_A,
char const *  alpha,
int *const *  permutation_B,
char const *  beta 
)

Permutes a tensor along each dimension skips if perm set to -1, generalizes slice. one of permutation_A or permutation_B has to be set to NULL, if multiworld read, then the parent world tensor should not be being permuted

Parameters
[in]Apure-operand tensor
[in]permutation_Amappings for each dimension of A indices
[in]alphascaling factor for A
[in]permutation_Bmappings for each dimension of B (this) indices
[in]betascaling factor for current values of B

Definition at line 720 of file untyped_tensor.cxx.

References ASSERT, CTF_int::depermute_keys(), has_zero_edge_len, is_sparse, lens, CTF::World::np, order, CTF_int::algstrct::pair_dealloc(), CTF_int::permute_keys(), read_local(), sr, write(), and wrld.

Referenced by ctf.core.tensor::__setitem__(), and CTF::Tensor< dtype >::permute().

void CTF_int::tensor::print ( FILE *  fp = stdout,
char const *  cutoff = NULL 
) const

prints tensor data to file using process 0

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

Definition at line 1825 of file untyped_tensor.cxx.

References CTF_int::alloc_ptr(), CTF_int::cdealloc(), read_local_nnz(), and sparsify().

Referenced by CTF::Tensor< dtype >::print(), and CTF::Tensor< dtype >::prnt().

void CTF_int::tensor::print_map ( FILE *  stream = stdout,
bool  allcall = 1 
) const

displays mapping information

Parameters
[in]streamoutput log (e.g. stdout)
[in]allcall(if 1 print only with proc 0)

Definition at line 652 of file untyped_tensor.cxx.

References ctf.core::dim, ctf.core::np(), and CTF_int::PHYSICAL_MAP.

Referenced by compare(), CTF_int::scaling::execute(), CTF_int::get_len_ordering(), read(), and CTF_int::summation::sum_tensors().

void CTF_int::tensor::prnt ( ) const

Definition at line 1822 of file untyped_tensor.cxx.

void CTF_int::tensor::profile_off ( )

turn off profiling

Definition at line 711 of file untyped_tensor.cxx.

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

void CTF_int::tensor::profile_on ( )

turn on profiling

Definition at line 707 of file untyped_tensor.cxx.

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

void CTF_int::tensor::pull_alias ( tensor const *  other)

pulls data from an tensor with an aliased buffer

Parameters
[in]othertensor with aliased data to pull from

Definition at line 2176 of file untyped_tensor.cxx.

References ASSERT, CTF_int::copy_mapping(), data, edge_map, has_home, home_buffer, is_data_aliased, is_home, order, and topo.

Referenced by CTF_int::summation::estimate_time(), and CTF_int::get_len_ordering().

void CTF_int::tensor::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.

Parameters
[in]num_pairnumber of pairs to read
[in]alphascaling factor of read value
[in]betascaling factor of old value
[in]indices64-bit global indices
[in]datavalues (num_pair of them to read)

Definition at line 1203 of file untyped_tensor.cxx.

References CTF_int::PairIterator::read_val(), CTF_int::PairIterator::write_key(), and CTF_int::PairIterator::write_val().

Referenced by ctf.core.tensor::__getitem__(), CTF::Tensor< dtype >::get_mapped_data(), CTF::Tensor< dtype >::read(), and read().

int CTF_int::tensor::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.

Parameters
[in]num_pairnumber of pairs to read
[in]alphascaling factor of read value
[in]betascaling factor of old value
[in]mapped_datapairs to write

Definition at line 1221 of file untyped_tensor.cxx.

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

char * CTF_int::tensor::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

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

Definition at line 1299 of file untyped_tensor.cxx.

References clear_mapping(), data, has_home, is_data_aliased, is_home, NS, print_map(), read(), redistribute(), set_distribution(), and ctf.core::tsr.

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

int CTF_int::tensor::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.

Parameters
[in]num_pairnumber of pairs to read
[in,out]mapped_datapairs to read

Definition at line 1228 of file untyped_tensor.cxx.

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

PairIterator CTF_int::tensor::read_all_pairs ( int64_t *  num_pair,
bool  unpack 
)
protected

read all pairs with each processor (packed)

Parameters
[out]num_pairnumber of values read
[in]unpackwhether to read all or unique pairs return pair iterator with allocated all pairs read

Definition at line 1648 of file untyped_tensor.cxx.

References CTF_int::alloc_ptr(), CTF_int::cdealloc(), and CTF_int::PairIterator::sort().

void CTF_int::tensor::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()

Parameters
[in]filestream to read from, the user should open, (optionally) set view, and close after function
[in]offsetdisplacement in bytes at which to start in the file (ought ot be the same on all processors)

Definition at line 2941 of file untyped_tensor.cxx.

References CTF_int::alloc(), CTF_int::cdealloc(), NS, CTF_int::packed_size(), read_dense_from_file(), and CTF_int::summation::sum_tensors().

Referenced by checkpoint(), and read_dense_from_file().

void CTF_int::tensor::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()

Parameters
[in]filenamestream to read from

Definition at line 2982 of file untyped_tensor.cxx.

int CTF_int::tensor::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

Parameters
[out]num_pairnumber of values read
[out]mapped_datavalues read

Definition at line 1553 of file untyped_tensor.cxx.

References CTF_int::alloc_ptr(), ASSERT, CTF_int::mapping::calc_phase(), CTF_int::mapping::calc_phys_phase(), CTF_int::mapping::calc_phys_rank(), CTF_int::cdealloc(), CTF_int::mapping::cdt, data, edge_map, CTF_int::summation::execute(), has_zero_edge_len, is_folded, is_mapped, is_sparse, CTF_int::topology::lda, ctf.core::np(), NS, order, pad_edge_len, padding, CTF_int::PHYSICAL_MAP, CTF_int::read_loc_pairs(), read_local(), size, CTF_int::SUCCESS, sym, TAU_FSTART, TAU_FSTOP, topo, ctf.core::tsr, CTF_int::mapping::type, and write().

Referenced by compare(), CTF::Tensor< dtype >::get_local_data(), CTF::Tensor< dtype >::get_local_pairs(), permute(), CTF::Tensor< dtype >::read_local(), read_local(), and slice().

void CTF_int::tensor::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

Parameters
[out]num_pairnumber of values read
[out]indices64-bit global indices
[out]datavalues (num_pair of them to read)

Definition at line 1505 of file untyped_tensor.cxx.

References CTF_int::alloc(), CTF_int::ConstPairIterator::k(), and CTF_int::ConstPairIterator::read_val().

int CTF_int::tensor::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

Parameters
[out]num_pairnumber of values read
[out]mapped_datavalues read

Definition at line 1538 of file untyped_tensor.cxx.

References data, is_data_aliased, nnz_loc, sparsify(), and CTF_int::SUCCESS.

Referenced by extract_diag(), CTF::Tensor< dtype >::get_local_data(), CTF::Tensor< dtype >::get_local_pairs(), print(), and ctf.core.tensor::reshape().

void CTF_int::tensor::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

Parameters
[out]num_pairnumber of values read
[out]indices64-bit global indices
[out]datavalues (num_pair of them to read)

Definition at line 1522 of file untyped_tensor.cxx.

References CTF_int::alloc(), CTF_int::ConstPairIterator::k(), and CTF_int::ConstPairIterator::read_val().

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

int CTF_int::tensor::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

Parameters
[in]old_distprevious distribution to remap data from
[in]old_offsetsoffsets from corner of tensor
[in]old_permutationpermutation of rows/cols/...
[in]new_offsetsoffsets from corner of tensor
[in]new_permutationpermutation of rows/cols/...

Definition at line 2201 of file untyped_tensor.cxx.

References ABORT, CTF_int::alloc(), CTF_int::block_reshuffle(), CTF_int::can_block_reshuffle(), CTF_int::cdealloc(), dgtog_reshuffle(), DPRINTF, CTF_int::padded_reshuffle(), CTF_int::distribution::phase, CTF_int::distribution::size, CTF_int::spredist_mdl, CTF::Timer::start(), CTF::Timer::stop(), CTF_int::SUCCESS, and VPRINTF.

Referenced by CTF_int::scaling::execute(), CTF_int::get_len_ordering(), and read().

int CTF_int::tensor::reduce_sum ( char *  result)

Performs an elementwise summation reduction on a tensor.

Parameters
[out]resultresult of reduction operation

Definition at line 1771 of file untyped_tensor.cxx.

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

int CTF_int::tensor::reduce_sum ( char *  result,
algstrct const *  sr_other 
)

Performs an elementwise summation reduction on a tensor with summation defined by sr_other.

Parameters
[out]resultresult of reduction operation
[in]sr_otheran algebraic structure (at least a monoid) defining the summation operation

Definition at line 1775 of file untyped_tensor.cxx.

References ASSERT, data, CTF_int::accumulatable::el_size, CTF_int::summation::execute(), CTF_int::algstrct::mulid(), and CTF_int::SUCCESS.

int CTF_int::tensor::reduce_sumabs ( char *  result)

Performs an elementwise absolute value summation reduction on a tensor.

Parameters
[out]resultresult of reduction operation

Definition at line 1789 of file untyped_tensor.cxx.

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

int CTF_int::tensor::reduce_sumabs ( char *  result,
algstrct const *  sr_other 
)

Performs an elementwise absolute value summation reduction on a tensor.

Parameters
[out]resultresult of reduction operation
[in]sr_otheran algebraic structure (at least a monoid) defining the summation operation

Definition at line 1793 of file untyped_tensor.cxx.

References CTF_int::algstrct::abs, ASSERT, data, CTF_int::summation::execute(), CTF_int::algstrct::mulid(), CTF_int::SUCCESS, and univar_function().

int CTF_int::tensor::reduce_sumsq ( char *  result)

computes the sum of squares of the elements

Parameters
[out]resultresult of reduction operation

Definition at line 1808 of file untyped_tensor.cxx.

References ASSERT, data, CTF_int::contraction::execute(), and CTF_int::SUCCESS.

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

void CTF_int::tensor::register_size ( int64_t  size)

register buffer allocation for this tensor

Definition at line 2883 of file untyped_tensor.cxx.

References CTF_int::inc_tot_mem_used().

void CTF_int::tensor::remove_fold ( )

removes folding without doing transpose unsets is_folded and deletes rec_tsr

Definition at line 2067 of file untyped_tensor.cxx.

References CTF_int::cdealloc().

Referenced by CTF_int::get_len_ordering().

void CTF_int::tensor::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)

Parameters
[in]sym_maskidentifies which tensor indices are part of the symmetric group which diagonals we want to scale (i.e. sym_mask [1,1] does A["ii"]= (1./2.)*A["ii"])

Definition at line 2629 of file untyped_tensor.cxx.

References CTF_int::alloc_ptr(), ASSERT, CTF_int::mapping::calc_phase(), CTF_int::mapping::calc_phys_phase(), CTF_int::mapping::calc_phys_rank(), CTF_int::cdealloc(), CTF_int::mapping::cdt, ctf.core::np(), CTF_int::PHYSICAL_MAP, CTF_int::scal_diag(), TAU_FSTART, TAU_FSTOP, and CTF_int::mapping::type.

Referenced by CTF_int::get_len_ordering().

tensor * CTF_int::tensor::self_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)

Parameters
[in]idx_Aindex map of this tensor as defined by summation/contraction
[out]new_idx_Ahow idx_A needs to be transformed
[in]order_Bnumber of modes in tensor B
[in]idx_Bindex map containing all indices in another tensor involved in the operaiton
[out]new_idx_Bhow idx_B needs to be transformed
[in]order_Cnumber of modes in tensor C
[in]idx_Cindex map containing all indices in another tensor involved in the operaiton (should be NULL for summation)
[out]new_idx_Chow idx_C needs to be transformed, untouched if idx_C is NULL
[in]idx_Cindex map containing all indices in another tensor involved in the operaiton (should be NULL for summation)
[out]new_idx_Chow idx_C needs to be transformed, untouched if idx_C is NULL

Definition at line 2990 of file untyped_tensor.cxx.

References CTF_int::alloc(), CTF_int::summation::execute(), and NS.

Referenced by CTF_int::summation::estimate_time(), and CTF_int::get_len_ordering().

int CTF_int::tensor::set ( char const *  val)

sets tensor data to val

Definition at line 493 of file untyped_tensor.cxx.

void CTF_int::tensor::set_distribution ( char const *  idx,
CTF::Idx_Partition const &  prl,
CTF::Idx_Partition const &  blk 
)
protected

set edge mappings as specified

Parameters
[in]idxassignment of characters to each dim
[in]prlmesh processor topology with character labels
[in]blklocal blocking with processor labels

Definition at line 1233 of file untyped_tensor.cxx.

References ASSERT, CTF_int::cdealloc(), CTF_int::mapping::cdt, CTF_int::check_self_mapping(), CTF_int::mapping::child, CTF_int::conv_idx(), CTF_int::find_topology(), CTF_int::mapping::has_child, CTF::Idx_Partition::idx, CTF::Partition::lens, CTF_int::NOT_MAPPED, CTF_int::mapping::np, CTF::Partition::order, CTF::Idx_Partition::part, CTF_int::PHYSICAL_MAP, CTF::Idx_Partition::reduce_order(), CTF_int::mapping::type, and CTF_int::VIRTUAL_MAP.

Referenced by read().

void CTF_int::tensor::set_name ( char const *  name)

set the tensor name

Parameters
[in]nameto set

Definition at line 697 of file untyped_tensor.cxx.

References CTF_int::alloc(), and CTF_int::cdealloc().

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

void CTF_int::tensor::set_new_nnz_glb ( int64_t const *  nnz_blk)

sets the number of nonzeros both locally (nnz_loc) and overall globally (nnz_tot)

Parameters
[in]nnz_blknumber of nonzeros in each block

Definition at line 2716 of file untyped_tensor.cxx.

Referenced by CTF_int::desymmetrize(), CTF_int::summation::estimate_time(), CTF_int::scaling::execute(), CTF_int::get_len_ordering(), CTF_int::summation::sum_tensors(), and CTF_int::symmetrize().

void CTF_int::tensor::set_padding ( )
void CTF_int::tensor::set_sym ( int const *  sym)

sets symmetry, WARNING: for internal use only !!!!

Parameters
[in]sym

Definition at line 2701 of file untyped_tensor.cxx.

References NS.

Referenced by CTF_int::get_len_ordering(), and CTF_int::summation::sum_tensors().

void CTF_int::tensor::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)

Parameters
[in]offsets_Bbottom left corner of block
[in]ends_Btop 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 861 of file untyped_tensor.cxx.

References CTF_int::alloc(), ASSERT, CTF_int::cdealloc(), CTF_int::depad_tsr(), has_zero_edge_len, lens, CTF::World::np, NS, order, CTF_int::pad_key(), CTF_int::algstrct::pair_alloc(), CTF_int::algstrct::pair_dealloc(), read_local(), sr, sym, write(), and wrld.

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

template<typename dtype >
template void CTF_int::tensor::smaller_equal_than< bool > ( 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)

Parameters
[in]Afirst operand
[in]Bsecond operand

Definition at line 63 of file untyped_tensor_tmpl.h.

References ctf.core::a, ctf.core::b, ctf.core::dtype, CTF_int::accumulatable::el_size, operator[](), order, and sr.

template<typename dtype >
template void CTF_int::tensor::smaller_than< bool > ( 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)

Parameters
[in]Afirst operand
[in]Bsecond operand

Definition at line 51 of file untyped_tensor_tmpl.h.

References ctf.core::a, ctf.core::b, ctf.core::dtype, CTF_int::accumulatable::el_size, operator[](), order, and sr.

int CTF_int::tensor::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.

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 1332 of file untyped_tensor.cxx.

References CTF_int::SUCCESS.

Referenced by CTF_int::get_len_ordering(), print(), read_local_nnz(), CTF::Tensor< dtype >::sparsify(), CTF_int::subsample(), and CTF_int::summation::sum_tensors().

int CTF_int::tensor::sparsify ( std::function< bool(char const *)>  f)

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

Parameters
[in]fboolean function to apply to values to determine whether to keep them, must be deterministic

Definition at line 1355 of file untyped_tensor.cxx.

References CTF_int::alloc(), CTF_int::alloc_ptr(), ASSERT, CTF_int::mapping::calc_phase(), CTF_int::mapping::calc_phys_phase(), CTF_int::mapping::calc_phys_rank(), CTF_int::cdealloc(), CTF_int::mapping::cdt, CTF_int::depad_tsr(), CTF_int::pad_key(), CTF_int::PHYSICAL_MAP, CTF_int::spsfy_tsr(), CTF_int::SUCCESS, TAU_FSTART, TAU_FSTOP, CTF_int::mapping::type, and CTF_int::PairIterator::write().

void CTF_int::tensor::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

Parameters
[in]mnumber of rows in matrix
[in]nnumber of columns in matrix
[in]nrow_idxnumber of indices to fold into column
[in]csrwhether to do csr (1) or coo (0) layout

Definition at line 2728 of file untyped_tensor.cxx.

References CTF_int::COO_Matrix::all_data, CTF_int::alloc(), CTF_int::alloc_ptr(), ASSERT, CTF_int::cdealloc(), CTF_int::get_coo_size(), CTF_int::get_csr_size(), CTF_int::COO_Matrix::set_data(), TAU_FSTART, and TAU_FSTOP.

template<typename dtype >
template void CTF_int::tensor::true_divide< bool > ( tensor A)

Definition at line 99 of file untyped_tensor_tmpl.h.

References ctf.core::a, ctf.core::dtype, operator[](), and order.

void CTF_int::tensor::unfold ( bool  was_mod = 0)

undo the folding of a local tensor block unsets is_folded and deletes rec_tsr

Parameters
[in]was_modtrue if data was modified, controls whether to discard sparse data

Definition at line 2029 of file untyped_tensor.cxx.

References CTF_int::alloc_ptr(), ASSERT, CTF_int::calc_dim(), CTF_int::cdealloc(), CTF_int::nosym_transpose(), NS, and CTF_int::sy_packed_size().

Referenced by CTF_int::summation::estimate_time(), CTF_int::scaling::execute(), CTF_int::get_len_ordering(), and CTF_int::summation::sum_tensors().

int CTF_int::tensor::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.

Parameters
[in]num_pairnumber of pairs to write
[in]alphascaling factor of written (read) value
[in]betascaling factor of old (existing) value
[in]mapped_datapairs to write
[in]rwweather to read (r) or write (w)

Definition at line 1082 of file untyped_tensor.cxx.

References CTF_int::alloc_ptr(), ASSERT, CTF_int::mapping::calc_phase(), CTF_int::mapping::calc_phys_phase(), CTF_int::mapping::calc_phys_rank(), CTF_int::cdealloc(), CTF_int::mapping::cdt, data, DEBUG_PRINTF, edge_map, CTF_int::ERROR, has_zero_edge_len, is_mapped, CTF_int::topology::lda, order, pad_edge_len, padding, CTF_int::PHYSICAL_MAP, CTF_int::proc_bytes_available(), set_padding(), CTF_int::SUCCESS, sym, TAU_FSTART, TAU_FSTOP, topo, ctf.core::tsr, CTF_int::mapping::type, and CTF_int::wr_pairs_layout().

Referenced by ctf.core.tensor::__setitem__(), extract_diag(), ctf.core.tensor::from_nparray(), permute(), read_local(), slice(), and CTF::Tensor< dtype >::write().

void CTF_int::tensor::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.

Parameters
[in]num_pairnumber of pairs to write
[in]alphascaling factor of written (read) value
[in]betascaling factor of old (existing) value
[in]indices64-bit global indices
[in]datavalues (num_pair of them)

Definition at line 1067 of file untyped_tensor.cxx.

References CTF_int::PairIterator::write_key(), and CTF_int::PairIterator::write_val().

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

void CTF_int::tensor::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

Parameters
[in,out]filestream to write to, the user should open, (optionally) set view, and close after function
[in]offsetdisplacement in bytes at which to start in the file (ought ot be the same on all processors)

Definition at line 2894 of file untyped_tensor.cxx.

References CTF_int::alloc(), CTF_int::cdealloc(), NS, CTF_int::packed_size(), CTF_int::PairIterator::read_val(), write_dense_to_file(), and CTF_int::PairIterator::write_key().

Referenced by checkpoint(), and write_dense_to_file().

void CTF_int::tensor::write_dense_to_file ( char const *  filename)

write all tensor data to binary file in element order, unpacking from sparse or symmetric formats

Parameters
[in]filenamestream to write to

Definition at line 2933 of file untyped_tensor.cxx.

Field Documentation

bool CTF_int::tensor::has_home
bool CTF_int::tensor::has_zero_edge_len

if true tensor has a zero edge length, so is zero, which short-cuts stuff

Definition at line 115 of file untyped_tensor.h.

Referenced by CTF_int::scaling::execute(), CTF_int::get_len_ordering(), permute(), read_local(), slice(), tensor(), and write().

char* CTF_int::tensor::home_buffer

buffer associated with home mapping of tensor, to which it is returned

Definition at line 121 of file untyped_tensor.h.

Referenced by copy_tensor_data(), CTF_int::desymmetrize(), CTF_int::summation::estimate_time(), CTF_int::scaling::execute(), CTF_int::get_len_ordering(), CTF_int::nosym_transpose(), and pull_alias().

int64_t CTF_int::tensor::home_size
int* CTF_int::tensor::inner_ordering

ordering of the dimensions according to which the tensori s folded

Definition at line 104 of file untyped_tensor.h.

Referenced by copy_tensor_data(), CTF_int::get_len_ordering(), CTF_int::seq_tsr_ctr::seq_tsr_ctr(), and CTF_int::seq_tsr_spctr::seq_tsr_spctr().

bool CTF_int::tensor::is_csr

whether CSR or COO if folded

Definition at line 133 of file untyped_tensor.h.

bool CTF_int::tensor::is_cyclic

whether the tensor data is cyclically distributed (blocked if false)

Definition at line 108 of file untyped_tensor.h.

Referenced by copy_tensor_data(), CTF_int::distribution::distribution(), CTF_int::scaling::execute(), and CTF_int::get_len_ordering().

bool CTF_int::tensor::is_data_aliased

whether the tensor data is an alias of another tensor object's data

Definition at line 110 of file untyped_tensor.h.

Referenced by CTF_int::desymmetrize(), CTF_int::summation::estimate_time(), CTF_int::scaling::execute(), CTF_int::get_len_ordering(), pull_alias(), read(), and read_local_nnz().

bool CTF_int::tensor::is_folded

whether the data is folded/transposed into a (lower-order) tensor

Definition at line 102 of file untyped_tensor.h.

Referenced by align(), copy_tensor_data(), and read_local().

bool CTF_int::tensor::is_home
bool CTF_int::tensor::is_mapped
int CTF_int::tensor::is_scp_padded

whether tensor data has additional padding

Definition at line 86 of file untyped_tensor.h.

bool CTF_int::tensor::left_home_transp

whether the tensor left home to transpose

Definition at line 127 of file untyped_tensor.h.

Referenced by CTF_int::nosym_transpose().

int64_t CTF_int::tensor::nnz_loc
int64_t CTF_int::tensor::nnz_tot

maximum number of local nonzero elements over all procs

Definition at line 139 of file untyped_tensor.h.

Referenced by btwn_cnt_fast(), copy_tensor_data(), CTF_int::summation::estimate_time(), CTF_int::get_len_ordering(), mis(), mis2(), setup_3d_Poisson(), setup_unstructured(), test_mis(), and test_mis2().

int CTF_int::tensor::nrow_idx

how many modes are folded into matricized row

Definition at line 135 of file untyped_tensor.h.

int CTF_int::tensor::order

number of tensor dimensions

Definition at line 76 of file untyped_tensor.h.

Referenced by ctf.core.tensor::__get__(), ctf.core.tensor::__getitem__(), ctf.core.tensor::__setitem__(), CTF_int::abs_helper(), add_from_subworld(), CTF_int::add_sym_perm(), add_to_subworld(), align(), block_sparse(), CTF_int::calc_fold_lnmk(), check_asym(), CTF_int::check_self_mapping(), check_sym(), compare(), compare_elementwise(), CTF_int::conj_helper(), CTF_int::contraction::contraction(), conv_type(), copy_tensor_data(), CTF_int::ctr_2d_gen_build(), CTF_int::ctr_virt::ctr_virt(), CTF_int::desymmetrize(), CTF_int::distribution::distribution(), CTF_int::summation::estimate_time(), CTF_int::scaling::execute(), exp_helper(), fast_tensor_ctr(), CTF::fill_sp_random_base(), CTF::get_core_tensor(), get_core_tensor(), get_core_tensor_hooi(), get_factor_matrices(), CTF::get_factor_matrices(), CTF_int::get_full_intm(), CTF_int::get_imag(), CTF_int::get_len_ordering(), get_rand_as_tsr(), CTF_int::get_real(), CTF_int::get_sym_perms(), CTF::Idx_Tensor::get_uniq_inds(), hooi(), ctf.core.tensor::i(), CTF::Idx_Tensor::Idx_Tensor(), larger_equal_than(), larger_than(), CTF_int::map_self_indices(), CTF::Matrix< dtype >::Matrix(), not_equals(), CTF::Scalar< dtype >::operator=(), CTF::Tensor< dtype >::operator=(), CTF_int::order_perm(), permute(), ctf.core.tensor::permute(), pull_alias(), read_local(), CTF::read_sparse_from_file_base(), CTF::real_norm1(), rec_scan(), CTF_int::seq_tsr_ctr::seq_tsr_ctr(), CTF_int::seq_tsr_spctr::seq_tsr_spctr(), CTF_int::seq_tsr_spsum::seq_tsr_spsum(), CTF_int::seq_tsr_sum::seq_tsr_sum(), CTF_int::set_imag(), CTF_int::set_real(), slice(), CTF::Tensor< dtype >::slice(), smaller_equal_than(), smaller_than(), CTF_int::spctr_pin_keys::spctr_pin_keys(), CTF_int::spctr_virt::spctr_virt(), CTF_int::sum_bool_tsr(), CTF_int::summation::sum_tensors(), CTF_int::summation::summation(), CTF_int::symmetrize(), tensor(), ctf.core.tensor::transpose(), true_divide(), CTF_int::tspsum_map::tspsum_map(), CTF_int::tspsum_permute::tspsum_permute(), CTF_int::tspsum_pin_keys::tspsum_pin_keys(), CTF_int::tspsum_virt::tspsum_virt(), CTF_int::tsum_virt::tsum_virt(), CTF::Vector< dtype >::Vector(), write(), and CTF::write_sparse_to_file_base().

int* CTF_int::tensor::padding

padding along each edge length (less than distribution phase)

Definition at line 82 of file untyped_tensor.h.

Referenced by copy_tensor_data(), CTF_int::distribution::distribution(), read_local(), and write().

bool CTF_int::tensor::profile

whether profiling should be done for contractions/sums involving this tensor

Definition at line 129 of file untyped_tensor.h.

Referenced by CTF_int::desymmetrize(), CTF_int::get_len_ordering(), CTF::Scalar< dtype >::operator=(), CTF::Tensor< dtype >::operator=(), CTF_int::symmetrize(), and tensor().

tensor* CTF_int::tensor::rec_tsr
int64_t CTF_int::tensor::registered_alloc_size

size CTF keeps track of for memory usage

Definition at line 100 of file untyped_tensor.h.

int* CTF_int::tensor::scp_padding

additional padding, may be greater than ScaLAPACK phase

Definition at line 88 of file untyped_tensor.h.

tensor* CTF_int::tensor::slay

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...

Definition at line 113 of file untyped_tensor.h.

int* CTF_int::tensor::sym_table

order-by-order table of dimensional symmetry relations

Definition at line 90 of file untyped_tensor.h.

Referenced by CTF_int::get_len_ordering(), and CTF_int::map_self_indices().

CTF::World* CTF_int::tensor::wrld

distributed processor context on which tensor is defined

Definition at line 70 of file untyped_tensor.h.

Referenced by align(), CTF::Tensor< dtype >::align(), ao_mo_transf_naive(), ao_mo_transf_slice(), btwn_cnt_fast(), btwn_cnt_naive(), check_asym(), check_sym(), compare(), copy_tensor_data(), CTF_int::desymmetrize(), divide_EaEi(), CTF_int::Bifun_Term::execute(), CTF_int::scaling::execute(), CTF::Tensor< dtype >::fill_random(), CTF::Tensor< dtype >::fill_sp_random(), CTF::fill_sp_random_base(), CTF::get_core_tensor(), CTF::Matrix< dtype >::get_desc(), CTF::get_factor_matrices(), CTF_int::get_full_intm(), CTF_int::get_len_ordering(), get_rand_as_tsr(), CTF::Scalar< dtype >::get_val(), mis2(), mp3(), CTF::Tensor< dtype >::norm1(), CTF::Tensor< dtype >::norm2(), CTF::Tensor< dtype >::norm_infty(), CTF_int::Term::operator double(), CTF_int::Term::operator float(), CTF_int::Term::operator int(), CTF_int::Term::operator int64_t(), CTF::Scalar< dtype >::operator=(), CTF::Tensor< dtype >::operator=(), permute(), CTF::Matrix< dtype >::print_matrix(), CTF::Matrix< dtype >::qr(), CTF::Matrix< dtype >::read_mat(), CTF::read_sparse_from_file_base(), rec_scan(), recursive_matmul(), scan(), CTF_int::seq_tsr_ctr::seq_tsr_ctr(), CTF_int::seq_tsr_spctr::seq_tsr_spctr(), setup(), slice(), CTF::Tensor< dtype >::slice(), smooth_jacobi(), sparse_mp3(), CTF_int::sum_bool_tsr(), CTF_int::summation::sum_tensors(), CTF::Matrix< dtype >::svd(), CTF_int::symmetrize(), tensor(), test_alg_multigrid(), test_mis(), test_mis2(), vcycle(), CTF::Idx_Tensor::where_am_i(), CTF::Matrix< dtype >::write_mat(), and CTF::write_sparse_to_file_base().


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