QMCPACK
TraceSamples< T > Struct Template Reference
+ Inheritance diagram for TraceSamples< T >:
+ Collaboration diagram for TraceSamples< T >:

Public Member Functions

 TraceSamples ()
 
 ~TraceSamples ()
 
void set_verbose (bool v)
 
int size ()
 
void assign_sample_index (const std::string &domain, const std::string &name, int index, std::string label="")
 
template<int D>
Array< T, D > * checkout_array (const std::string &domain, const std::string &name, TinyVector< int, DMAX > shape)
 
template<int D>
Array< T, D > * checkout_array (const ParticleSet &P, const std::string &name, TinyVector< int, DMAX > shape)
 
TraceSample< T > * get_trace (const std::string &domain, const std::string &name)
 
CombinedTraceSample< T > * get_combined_trace (const std::string &domain, const std::string &name)
 
bool make_combined_trace (const std::string &name, std::vector< std::string > &names, std::vector< TraceReal > &weights)
 
void set_unit_size (int usize)
 
void screen_writes (TraceRequest &request)
 
void order_by_size ()
 
void set_buffer_ranges (int &starting_index)
 
int total_size ()
 
int min_buffer_index ()
 
int max_buffer_index ()
 
void combine_samples ()
 
void reset_combined_samples ()
 
void finalize ()
 
void register_hdf_data (hdf_archive &f)
 
void write_summary (std::string type, std::string pad=" ")
 
void user_report (const std::string &type, const std::string &pad=" ")
 

Public Attributes

std::vector< TraceSample< T > * > samples
 
std::map< std::string, std::map< std::string, int > > sample_indices
 
std::vector< TraceSample< T > * > ordered_samples
 
std::vector< CombinedTraceSample< T > * > combined_samples
 
std::vector< Vector< T > * > combined_sample_vectors
 
bool verbose
 

Detailed Description

template<typename T>
struct qmcplusplus::TraceSamples< T >

Definition at line 751 of file TraceManager.h.

Constructor & Destructor Documentation

◆ TraceSamples()

TraceSamples ( )
inline

Definition at line 760 of file TraceManager.h.

760 : verbose(false) {}

◆ ~TraceSamples()

~TraceSamples ( )
inline

Definition at line 762 of file TraceManager.h.

762 { finalize(); }

Member Function Documentation

◆ assign_sample_index()

void assign_sample_index ( const std::string &  domain,
const std::string &  name,
int  index,
std::string  label = "" 
)
inline

Definition at line 771 of file TraceManager.h.

Referenced by TraceSamples< std::complex< TraceReal > >::checkout_array(), and TraceSamples< std::complex< TraceReal > >::make_combined_trace().

772  {
773  if (sample_indices.count(domain) > 0 && sample_indices[domain].count(name) > 0)
774  {
775  APP_ABORT("TraceSamples::checkout " + label + " variable " + name + " already exists in domain " + domain);
776  }
777  else
778  {
779  sample_indices[domain][name] = index;
780  }
781  }
#define APP_ABORT(msg)
Widely used but deprecated fatal error macros from legacy code.
Definition: AppAbort.h:27
std::map< std::string, std::map< std::string, int > > sample_indices
Definition: TraceManager.h:754

◆ checkout_array() [1/2]

Array<T, D>* checkout_array ( const std::string &  domain,
const std::string &  name,
TinyVector< int, DMAX shape 
)
inline

Definition at line 784 of file TraceManager.h.

Referenced by TraceManager::checkout_complex(), TraceManager::checkout_int(), TraceManager::checkout_real(), and qmcplusplus::TEST_CASE().

785  {
786  int index = samples.size();
787  assign_sample_index(domain, name, index, "array");
788  std::array<size_t, D> subshape;
789  for (int idim = 0; idim < D; idim++)
790  subshape[idim] = shape[idim];
791  Array<T, D>* a = new Array<T, D>(subshape);
792  TraceSample<T>* s = new TraceSample<T>(domain, name, index, D, shape, a->storage());
793  samples.push_back(s);
794  if (verbose)
795  app_log() << "TraceSamples::checkout_array " << domain << " " << name << " " << index << std::endl;
796  return a;
797  }
std::ostream & app_log()
Definition: OutputManager.h:65
Container_t & storage()
Definition: OhmmsArray.h:60
std::vector< TraceSample< T > * > samples
Definition: TraceManager.h:753
void assign_sample_index(const std::string &domain, const std::string &name, int index, std::string label="")
Definition: TraceManager.h:771

◆ checkout_array() [2/2]

Array<T, D>* checkout_array ( const ParticleSet P,
const std::string &  name,
TinyVector< int, DMAX shape 
)
inline

Definition at line 801 of file TraceManager.h.

802  {
803  const std::string& domain = P.parentName();
804  int index = samples.size();
805  assign_sample_index(domain, name, index, "array");
806  std::array<size_t, D> subshape;
807  for (int idim = 0; idim < D; idim++)
808  subshape[idim] = shape[idim];
809  Array<T, D>* a = new Array<T, D>(subshape);
810  TraceSample<T>* s = new TraceSample<T>(domain, name, index, D, shape, a->storage());
811  samples.push_back(s);
812  s->array_trace = true;
813  if (verbose)
814  app_log() << "TraceSamples::checkout_array " << domain << " " << name << " " << index << std::endl;
815  return a;
816  }
std::ostream & app_log()
Definition: OutputManager.h:65
Container_t & storage()
Definition: OhmmsArray.h:60
std::vector< TraceSample< T > * > samples
Definition: TraceManager.h:753
void assign_sample_index(const std::string &domain, const std::string &name, int index, std::string label="")
Definition: TraceManager.h:771

◆ combine_samples()

void combine_samples ( )
inline

Definition at line 963 of file TraceManager.h.

Referenced by TraceBuffer< TraceInt >::collect_sample().

964  {
965  for (int i = 0; i < combined_samples.size(); ++i)
966  combined_samples[i]->combine();
967  }
std::vector< CombinedTraceSample< T > * > combined_samples
Definition: TraceManager.h:756

◆ finalize()

void finalize ( )
inline

Definition at line 977 of file TraceManager.h.

Referenced by TraceManager::finalize_traces(), and TraceSamples< std::complex< TraceReal > >::~TraceSamples().

978  {
979  //note combined_samples pointers are a subset of those in samples
980  // and so only samples needs to be deleted
981  delete_iter(samples.begin(), samples.end());
983  samples.resize(0);
984  ordered_samples.resize(0);
985  combined_samples.resize(0);
986  combined_sample_vectors.resize(0);
987  sample_indices.clear();
988  }
void delete_iter(IT first, IT last)
delete the pointers in [first,last)
std::vector< Vector< T > * > combined_sample_vectors
Definition: TraceManager.h:757
std::vector< TraceSample< T > * > ordered_samples
Definition: TraceManager.h:755
std::vector< CombinedTraceSample< T > * > combined_samples
Definition: TraceManager.h:756
std::vector< TraceSample< T > * > samples
Definition: TraceManager.h:753
std::map< std::string, std::map< std::string, int > > sample_indices
Definition: TraceManager.h:754

◆ get_combined_trace()

CombinedTraceSample<T>* get_combined_trace ( const std::string &  domain,
const std::string &  name 
)
inline

Definition at line 837 of file TraceManager.h.

Referenced by TraceManager::get_complex_combined_trace(), TraceManager::get_int_combined_trace(), and TraceManager::get_real_combined_trace().

838  {
839  CombinedTraceSample<T>* ts = NULL;
840  for (int i = 0; i < combined_samples.size(); ++i)
841  {
842  CombinedTraceSample<T>& tsc = *combined_samples[i];
843  if (tsc.domain == domain && tsc.name == name)
844  {
845  ts = &tsc;
846  break;
847  }
848  }
849  return ts;
850  }
std::vector< CombinedTraceSample< T > * > combined_samples
Definition: TraceManager.h:756

◆ get_trace()

TraceSample<T>* get_trace ( const std::string &  domain,
const std::string &  name 
)
inline

Definition at line 819 of file TraceManager.h.

Referenced by TraceManager::get_complex_trace(), TraceManager::get_int_trace(), and TraceManager::get_real_trace().

820  {
821  TraceSample<T>* ts = NULL;
822  for (int i = 0; i < samples.size(); ++i)
823  {
824  TraceSample<T>& tsc = *samples[i];
825  if (tsc.domain == domain && tsc.name == name)
826  {
827  ts = &tsc;
828  break;
829  }
830  }
831  if (ts == NULL)
832  APP_ABORT("TraceSamples::get_trace failed to get trace for quantity " + name + " in domain " + domain);
833  return ts;
834  }
#define APP_ABORT(msg)
Widely used but deprecated fatal error macros from legacy code.
Definition: AppAbort.h:27
std::vector< TraceSample< T > * > samples
Definition: TraceManager.h:753

◆ make_combined_trace()

bool make_combined_trace ( const std::string &  name,
std::vector< std::string > &  names,
std::vector< TraceReal > &  weights 
)
inline

Definition at line 853 of file TraceManager.h.

Referenced by TraceBuffer< TraceInt >::make_combined_trace().

856  {
857  bool created = false;
858  if (names.size() != weights.size())
859  APP_ABORT("TraceSamples::make_combined_trace names and weights must be the same size");
860  std::map<std::string, std::map<std::string, int>>::iterator it;
861  for (it = sample_indices.begin(); it != sample_indices.end(); it++)
862  {
863  std::string domain = it->first;
864  std::map<std::string, int>& indices = it->second;
865  bool any_present = false;
866  for (int i = 0; i < names.size(); ++i)
867  any_present = any_present || indices.count(names[i]) > 0;
868  if (any_present)
869  {
870  int index = samples.size();
871  auto* sample = new Vector<T>;
872  CombinedTraceSample<T>* combined = new CombinedTraceSample<T>(domain, name, index, 0, *sample);
873  for (int i = 0; i < names.size(); ++i)
874  {
875  if (indices.count(names[i]) > 0)
876  {
877  TraceSample<T>* component = samples[indices[names[i]]];
878  combined->add_component(component, weights[i]);
879  }
880  }
881  assign_sample_index(domain, name, index);
882  samples.push_back(combined);
883  combined_samples.push_back(combined);
884  combined_sample_vectors.push_back(sample);
885  created = true;
886  }
887  }
888  return created;
889  }
std::vector< Vector< T > * > combined_sample_vectors
Definition: TraceManager.h:757
#define APP_ABORT(msg)
Widely used but deprecated fatal error macros from legacy code.
Definition: AppAbort.h:27
std::vector< CombinedTraceSample< T > * > combined_samples
Definition: TraceManager.h:756
std::vector< TraceSample< T > * > samples
Definition: TraceManager.h:753
void assign_sample_index(const std::string &domain, const std::string &name, int index, std::string label="")
Definition: TraceManager.h:771
std::map< std::string, std::map< std::string, int > > sample_indices
Definition: TraceManager.h:754

◆ max_buffer_index()

int max_buffer_index ( )
inline

Definition at line 954 of file TraceManager.h.

Referenced by TraceBuffer< TraceInt >::test_buffer_write().

955  {
956  int max_index = -1;
957  for (int i = 0; i < samples.size(); i++)
958  max_index = std::max(max_index, samples[i]->buffer_end);
959  return max_index;
960  }
std::vector< TraceSample< T > * > samples
Definition: TraceManager.h:753

◆ min_buffer_index()

int min_buffer_index ( )
inline

Definition at line 945 of file TraceManager.h.

Referenced by TraceBuffer< TraceInt >::test_buffer_write().

946  {
947  int min_index = 2000000000;
948  for (int i = 0; i < samples.size(); i++)
949  min_index = std::min(min_index, samples[i]->buffer_start);
950  return min_index;
951  }
T min(T a, T b)
std::vector< TraceSample< T > * > samples
Definition: TraceManager.h:753

◆ order_by_size()

void order_by_size ( )
inline

Definition at line 916 of file TraceManager.h.

Referenced by TraceBuffer< TraceInt >::order_and_resize().

917  {
918  for (int i = 0; i < samples.size(); i++)
919  samples[i]->set_data_size();
920  ordered_samples.resize(samples.size());
921  copy(samples.begin(), samples.end(), ordered_samples.begin());
922  sort(ordered_samples.begin(), ordered_samples.end(), TraceSample_comp<T>);
923  }
std::vector< TraceSample< T > * > ordered_samples
Definition: TraceManager.h:755
void copy(const Array< T1, 3 > &src, Array< T2, 3 > &dest)
Definition: Blitz.h:639
std::vector< TraceSample< T > * > samples
Definition: TraceManager.h:753

◆ register_hdf_data()

void register_hdf_data ( hdf_archive f)
inline

Definition at line 991 of file TraceManager.h.

Referenced by TraceBuffer< TraceInt >::register_hdf_data().

992  {
993  std::map<std::string, std::map<std::string, int>>::iterator it;
994  std::map<std::string, int>::iterator it2;
995  for (it = sample_indices.begin(); it != sample_indices.end(); it++)
996  {
997  const std::string& domain = it->first;
998  std::map<std::string, int>& indices = it->second;
999  f.push(domain);
1000  for (it2 = indices.begin(); it2 != indices.end(); ++it2)
1001  {
1002  const std::string& quantity = it2->first;
1003  TraceSample<T>& sample = *samples[it2->second];
1004  if (sample.write)
1005  {
1006  f.push(quantity);
1007  f.write(sample.dimension, "dimension");
1008  f.write(sample.shape, "shape");
1009  f.write(sample.size, "size");
1010  f.write(sample.unit_size, "unit_size");
1011  f.write(sample.buffer_start, "row_start");
1012  f.write(sample.buffer_end, "row_end");
1013  f.pop();
1014  }
1015  }
1016  f.pop();
1017  }
1018  }
std::vector< TraceSample< T > * > samples
Definition: TraceManager.h:753
std::map< std::string, std::map< std::string, int > > sample_indices
Definition: TraceManager.h:754

◆ reset_combined_samples()

void reset_combined_samples ( )
inline

Definition at line 970 of file TraceManager.h.

Referenced by TraceBuffer< TraceInt >::collect_sample().

971  {
972  for (int i = 0; i < combined_samples.size(); ++i)
973  combined_samples[i]->reset();
974  }
std::vector< CombinedTraceSample< T > * > combined_samples
Definition: TraceManager.h:756

◆ screen_writes()

void screen_writes ( TraceRequest request)
inline

Definition at line 899 of file TraceManager.h.

Referenced by TraceManager::screen_writes().

900  {
901  for (int i = 0; i < samples.size(); i++)
902  {
903  TraceSample<T>& s = *samples[i];
904  bool stream = request.screen_sample(s.domain, s.name, s.write);
905  if (verbose)
906  app_log() << "TraceRequest screening " << s.name << " in domain " << s.domain << ". stream: " << stream
907  << " write: " << s.write << std::endl;
908  if (!stream && !s.is_combined())
909  app_log() << "warning: quantity " + s.name + " in domain " + s.domain +
910  " was not requested but is streaming anyway"
911  << std::endl;
912  }
913  }
std::ostream & app_log()
Definition: OutputManager.h:65
std::vector< TraceSample< T > * > samples
Definition: TraceManager.h:753

◆ set_buffer_ranges()

void set_buffer_ranges ( int &  starting_index)
inline

Definition at line 926 of file TraceManager.h.

Referenced by TraceBuffer< TraceInt >::order_and_resize().

927  {
928  for (int i = 0; i < ordered_samples.size(); i++)
929  {
930  TraceSample<T>& sample = *ordered_samples[i];
931  sample.set_buffer_range(starting_index);
932  }
933  }
std::vector< TraceSample< T > * > ordered_samples
Definition: TraceManager.h:755

◆ set_unit_size()

void set_unit_size ( int  usize)
inline

Definition at line 892 of file TraceManager.h.

Referenced by TraceBuffer< TraceInt >::order_and_resize(), and TraceSamples< std::complex< TraceReal > >::set_unit_size().

893  {
894  for (int i = 0; i < samples.size(); i++)
895  samples[i]->set_unit_size(usize);
896  }
void set_unit_size(int usize)
Definition: TraceManager.h:892
std::vector< TraceSample< T > * > samples
Definition: TraceManager.h:753

◆ set_verbose()

void set_verbose ( bool  v)
inline

Definition at line 765 of file TraceManager.h.

Referenced by TraceManager::distribute().

765 { verbose = v; }

◆ size()

int size ( void  )
inline

Definition at line 768 of file TraceManager.h.

Referenced by TraceBuffer< TraceInt >::test_buffer_write().

768 { return samples.size(); }
std::vector< TraceSample< T > * > samples
Definition: TraceManager.h:753

◆ total_size()

int total_size ( )
inline

Definition at line 936 of file TraceManager.h.

Referenced by TraceBuffer< TraceInt >::test_buffer_write().

937  {
938  int s = 0;
939  for (int i = 0; i < samples.size(); i++)
940  s += samples[i]->sample.size() * samples[i]->unit_size;
941  return s;
942  }
std::vector< TraceSample< T > * > samples
Definition: TraceManager.h:753

◆ user_report()

void user_report ( const std::string &  type,
const std::string &  pad = "  " 
)
inline

Definition at line 1059 of file TraceManager.h.

Referenced by TraceBuffer< TraceInt >::user_report().

1060  {
1061  std::string pad2 = pad + " ";
1062  app_log() << pad << type << " traces provided by estimators" << std::endl;
1063  std::map<std::string, std::map<std::string, int>>::iterator it;
1064  std::map<std::string, int>::iterator it2;
1065  for (it = sample_indices.begin(); it != sample_indices.end(); it++)
1066  {
1067  const std::string& domain = it->first;
1068  std::map<std::string, int>& indices = it->second;
1069  app_log() << pad2 << "domain " << domain << ": ";
1070  int n = 0;
1071  for (it2 = indices.begin(); it2 != indices.end(); ++it2)
1072  {
1073  if (n % 5 == 0)
1074  app_log() << std::endl << pad2 << " ";
1075  n++;
1076  const std::string& quantity = it2->first;
1077  app_log() << quantity << " ";
1078  }
1079  app_log() << std::endl;
1080  }
1081  }
std::ostream & app_log()
Definition: OutputManager.h:65
std::map< std::string, std::map< std::string, int > > sample_indices
Definition: TraceManager.h:754

◆ write_summary()

void write_summary ( std::string  type,
std::string  pad = "  " 
)
inline

Definition at line 1021 of file TraceManager.h.

Referenced by TraceSamples< std::complex< TraceReal > >::write_summary(), and TraceBuffer< TraceInt >::write_summary().

1022  {
1023  std::string pad2 = pad + " ";
1024  std::string pad3 = pad2 + " ";
1025  std::string pad4 = pad3 + " ";
1026  app_log() << pad << "TraceSamples<" << type << ">" << std::endl;
1027  app_log() << pad2 << "nsamples = " << samples.size() << std::endl;
1028  app_log() << pad2 << "ncombined_samples = " << combined_samples.size() << std::endl;
1029  app_log() << pad2 << "sample_indices" << std::endl;
1030  std::map<std::string, std::map<std::string, int>>::iterator it;
1031  std::map<std::string, int>::iterator it2;
1032  for (it = sample_indices.begin(); it != sample_indices.end(); it++)
1033  {
1034  const std::string& domain = it->first;
1035  std::map<std::string, int>& indices = it->second;
1036  app_log() << pad3 << domain << std::endl;
1037  for (it2 = indices.begin(); it2 != indices.end(); ++it2)
1038  app_log() << pad4 << it2->first << " = " << it2->second << std::endl;
1039  }
1040  app_log() << pad2 << "end sample_indices" << std::endl;
1041  app_log() << pad2 << "combined_sample_vectors = ";
1042  for (int i = 0; i < combined_sample_vectors.size(); ++i)
1043  app_log() << (size_t)combined_sample_vectors[i] << " ";
1044  app_log() << std::endl;
1045  app_log() << pad2 << "combined_samples" << std::endl;
1046  for (int i = 0; i < combined_samples.size(); ++i)
1047  combined_samples[i]->write_summary_combined(i, pad3);
1048  app_log() << pad2 << "end combined_samples" << std::endl;
1049  app_log() << pad2 << "samples" << std::endl;
1050  for (int i = 0; i < ordered_samples.size(); ++i)
1051  ordered_samples[i]->write_summary(i, pad3);
1052  //for(int i=0; i<samples.size(); ++i)
1053  // samples[i]->write_summary(i,pad3);
1054  app_log() << pad2 << "end samples" << std::endl;
1055  app_log() << pad << "end TraceSamples<" << type << ">" << std::endl;
1056  }
std::ostream & app_log()
Definition: OutputManager.h:65
std::vector< Vector< T > * > combined_sample_vectors
Definition: TraceManager.h:757
std::vector< TraceSample< T > * > ordered_samples
Definition: TraceManager.h:755
void write_summary(std::string type, std::string pad=" ")
std::vector< CombinedTraceSample< T > * > combined_samples
Definition: TraceManager.h:756
std::vector< TraceSample< T > * > samples
Definition: TraceManager.h:753
std::map< std::string, std::map< std::string, int > > sample_indices
Definition: TraceManager.h:754

Member Data Documentation

◆ combined_sample_vectors

◆ combined_samples

◆ ordered_samples

◆ sample_indices

◆ samples

◆ verbose


The documentation for this struct was generated from the following file: