Skip to content

Commit e4af5bf

Browse files
author
Martin D. Weinberg
committed
Added particle index to accumulate() members. Fixed Cylindrical calls with pcavar
1 parent ae7941f commit e4af5bf

9 files changed

Lines changed: 130 additions & 95 deletions

File tree

‎expui/BasisFactory.H‎

Lines changed: 4 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -184,13 +184,13 @@ namespace BasisClasses
184184
{ return coefret; }
185185

186186
//! Accumulate new coefficients from full phase space
187-
virtual void accumulate(double mass,
188-
double x, double y, double z,
189-
double u, double v, double w) = 0;
187+
virtual void accumulate(double x, double y, double z,
188+
double u, double v, double w,
189+
double mass, unsigned long int indx) = 0;
190190

191191
//! Accumulate new coefficients from coordinates only
192192
virtual void accumulate(double x, double y, double z,
193-
double mass) = 0;
193+
double mass, unsigned long int indx) = 0;
194194

195195
//! Get mass on grid
196196
double getMass(void) { return totalMass; }

‎expui/BiorthBasis.H‎

Lines changed: 31 additions & 12 deletions
Original file line numberDiff line numberDiff line change
@@ -110,6 +110,14 @@ namespace BasisClasses
110110
Eigen::VectorXd sampleMasses;
111111
//@}
112112

113+
//! Round time key to emulated fixed-point arithmetic
114+
double roundTime(double time)
115+
{
116+
// Eight decimal places should be enough here...
117+
const double multiplier = 1.0e+08; // std::pow(10.0, 8);
118+
return std::floor(time * multiplier + 0.5) / multiplier;
119+
}
120+
113121
//! Write H5 covariance; returns updated group count
114122
virtual unsigned writeCovarH5(HighFive::Group& group, unsigned count, double time);
115123

@@ -174,13 +182,14 @@ namespace BasisClasses
174182
//@}
175183

176184
//! Accumulate new coefficients
177-
virtual void accumulate(double x, double y, double z, double mass) = 0;
185+
virtual void accumulate(double x, double y, double z, double mass,
186+
unsigned long indx = 0) = 0;
178187

179188
//! Accumulate new coefficients
180-
virtual void accumulate(double mass,
181-
double x, double y, double z,
182-
double u, double v, double w)
183-
{ accumulate(x, y, z, mass); }
189+
virtual void accumulate(double x, double y, double z,
190+
double u, double v, double w,
191+
double mass, unsigned long indx = 0)
192+
{ accumulate(x, y, z, mass, indx); }
184193

185194
//! Get mass on grid
186195
double getMass(void) { return totalMass; }
@@ -264,6 +273,9 @@ namespace BasisClasses
264273
const std::string& runtag,
265274
double time=0.0);
266275

276+
//! Make covariance after accumulation
277+
virtual void makeCoefCovariance(void) {}
278+
267279
};
268280

269281
/**
@@ -398,7 +410,8 @@ namespace BasisClasses
398410
void make_coefs(void);
399411

400412
//! Accumulate new coefficients
401-
virtual void accumulate(double x, double y, double z, double mass);
413+
virtual void accumulate(double x, double y, double z, double mass,
414+
unsigned long indx = 0);
402415

403416
//! Return current maximum harmonic order in expansion
404417
int getLmax() { return lmax; }
@@ -643,7 +656,8 @@ namespace BasisClasses
643656
void make_coefs(void);
644657

645658
//! Accumulate new coefficients
646-
virtual void accumulate(double x, double y, double z, double mass);
659+
virtual void accumulate(double x, double y, double z, double mass,
660+
unsigned long indx = 0);
647661

648662
//! Return current maximum harmonic order in expansion
649663
int getMmax() { return mmax; }
@@ -800,7 +814,8 @@ namespace BasisClasses
800814
void make_coefs(void);
801815

802816
//! Accumulate new coefficients
803-
virtual void accumulate(double x, double y, double z, double mass);
817+
virtual void accumulate(double x, double y, double z, double mass,
818+
unsigned long indx = 0);
804819

805820
//! Return current maximum harmonic order in expansion
806821
int getMmax() { return mmax; }
@@ -958,7 +973,8 @@ namespace BasisClasses
958973
void make_coefs(void);
959974

960975
//! Accumulate new coefficients
961-
virtual void accumulate(double x, double y, double z, double mass);
976+
virtual void accumulate(double x, double y, double z, double mass,
977+
unsigned long indx = 0);
962978

963979
//! Return current maximum harmonic order in expansion
964980
int getMmax() { return mmax; }
@@ -1004,6 +1020,7 @@ namespace BasisClasses
10041020
std::tie(sampleCounts, sampleMasses) = sl->getCovarSamples();
10051021
return {sampleCounts, sampleMasses};
10061022
}
1023+
10071024
};
10081025

10091026
/**
@@ -1120,7 +1137,8 @@ namespace BasisClasses
11201137
void make_coefs(void);
11211138

11221139
//! Accumulate new coefficients
1123-
virtual void accumulate(double x, double y, double z, double mass);
1140+
virtual void accumulate(double x, double y, double z, double mass,
1141+
unsigned long indx = 0);
11241142

11251143
//! Return current maximum harmonic order in expansion
11261144
Eigen::Vector3i getNmax() { return {nmaxx, nmaxy, nmaxz}; }
@@ -1258,7 +1276,8 @@ namespace BasisClasses
12581276
void make_coefs(void);
12591277

12601278
//! Accumulate new coefficients
1261-
virtual void accumulate(double x, double y, double z, double mass);
1279+
virtual void accumulate(double x, double y, double z, double mass,
1280+
unsigned long indx = 0);
12621281

12631282
//! Return current maximum harmonic order in expansion
12641283
Eigen::Vector3i getNmax() { return {nmaxx, nmaxy, nmaxz}; }
@@ -1302,7 +1321,7 @@ namespace BasisClasses
13021321
std::vector<Eigen::VectorXd> sampleMasses;
13031322
std::vector<std::vector<std::vector<CoefCovarType>>> covarData;
13041323

1305-
const double fixedPointPrecision = 1.0e04;
1324+
const double fixedPointPrecision = 1.0e08;
13061325
std::map<int, unsigned> timeMap;
13071326
unsigned timeIndex(double time)
13081327
{

‎expui/BiorthBasis.cc‎

Lines changed: 48 additions & 35 deletions
Original file line numberDiff line numberDiff line change
@@ -265,7 +265,7 @@ namespace BasisClasses
265265
if (conf["EVEN_M"]) EVEN_M = conf["EVEN_M"].as<bool>();
266266
if (conf["M0_ONLY"]) M0_only = conf["M0_ONLY"].as<bool>();
267267
if (conf["pcavar"]) pcavar = conf["pcavar"].as<bool>();
268-
if (conf["sampT"]) sampT = conf["sampT"].as<int>();
268+
if (conf["subsamp"]) sampT = conf["subsamp"].as<int>();
269269
}
270270
catch (YAML::Exception & error) {
271271
if (myid==0) std::cout << "Error parsing parameter stanza for <"
@@ -568,7 +568,8 @@ namespace BasisClasses
568568
coefctr = {0.0, 0.0, 0.0};
569569
}
570570

571-
void Spherical::accumulate(double x, double y, double z, double mass)
571+
void Spherical::accumulate(double x, double y, double z, double mass,
572+
unsigned long indx)
572573
{
573574
double fac, fac1, fac2, fac4;
574575
double norm = -4.0*M_PI;
@@ -1444,6 +1445,8 @@ namespace BasisClasses
14441445
if (unmatched.size())
14451446
throw YamlConfigError("Basis::Basis::Cylindrical", "parameter", unmatched, __FILE__, __LINE__);
14461447

1448+
int sampT = 0;
1449+
14471450
// Assign values from YAML
14481451
//
14491452
try {
@@ -1495,6 +1498,7 @@ namespace BasisClasses
14951498
if (conf["vflag" ]) vflag = conf["vflag" ].as<int>();
14961499
if (conf["pyname" ]) pyname = conf["pyname" ].as<std::string>();
14971500
if (conf["pcavar"] ) pcavar = conf["pcavar" ].as<bool>();
1501+
if (conf["subsamp"] ) sampT = conf["subsamp" ].as<int>();
14981502

14991503
// Deprecation warning
15001504
if (conf["density" ]) {
@@ -1544,6 +1548,7 @@ namespace BasisClasses
15441548
EmpCylSL::CMAPZ = cmapZ;
15451549
EmpCylSL::logarithmic = logarithmic;
15461550
EmpCylSL::VFLAG = vflag;
1551+
EmpCylSL::PCAVAR = pcavar;
15471552

15481553
// Check for non-null cache file name. This must be specified
15491554
// to prevent recomputation and unexpected behavior.
@@ -1566,7 +1571,11 @@ namespace BasisClasses
15661571
//
15671572
if (mlim>=0) sl->set_mlim(mlim);
15681573
if (EVEN_M) sl->setEven(EVEN_M);
1569-
1574+
if (sampT>0) {
1575+
sl->setSampT(sampT);
1576+
sl->init_pca();
1577+
}
1578+
15701579
// Cache override for old Eigen cache
15711580
//
15721581
if (oldcache) sl->AllowOldCache();
@@ -1790,11 +1799,12 @@ namespace BasisClasses
17901799
}
17911800
}
17921801

1793-
void Cylindrical::accumulate(double x, double y, double z, double mass)
1802+
void Cylindrical::accumulate(double x, double y, double z, double mass,
1803+
unsigned long indx)
17941804
{
17951805
double R = sqrt(x*x + y*y);
17961806
double phi = atan2(y, x);
1797-
sl->accumulate(R, z, phi, mass, 0, 0);
1807+
sl->accumulate(R, z, phi, mass, indx, 0, 0, pcavar);
17981808
}
17991809

18001810
void Cylindrical::reset_coefs(void)
@@ -1857,7 +1867,7 @@ namespace BasisClasses
18571867

18581868
void Cylindrical::make_coefs(void)
18591869
{
1860-
sl->make_coefficients();
1870+
sl->make_coefficients(pcavar);
18611871
}
18621872

18631873

@@ -2170,7 +2180,8 @@ namespace BasisClasses
21702180
coefctr = {0.0, 0.0, 0.0};
21712181
}
21722182

2173-
void FlatDisk::accumulate(double x, double y, double z, double mass)
2183+
void FlatDisk::accumulate(double x, double y, double z, double mass,
2184+
unsigned long indx)
21742185
{
21752186
// Normalization factors
21762187
//
@@ -2974,7 +2985,8 @@ namespace BasisClasses
29742985
coefctr = {0.0, 0.0, 0.0};
29752986
}
29762987

2977-
void CBDisk::accumulate(double x, double y, double z, double mass)
2988+
void CBDisk::accumulate(double x, double y, double z, double mass,
2989+
unsigned long int idx)
29782990
{
29792991
// Normalization factors
29802992
//
@@ -3432,7 +3444,8 @@ namespace BasisClasses
34323444
coefctr = {0.0, 0.0, 0.0};
34333445
}
34343446

3435-
void Slab::accumulate(double x, double y, double z, double mass)
3447+
void Slab::accumulate(double x, double y, double z, double mass,
3448+
unsigned long int idx)
34363449
{
34373450
// Truncate to slab with sides in [0,1]
34383451
if (x<0.0)
@@ -3949,7 +3962,8 @@ namespace BasisClasses
39493962
coefctr = {0.0, 0.0, 0.0};
39503963
}
39513964

3952-
void Cube::accumulate(double x, double y, double z, double mass)
3965+
void Cube::accumulate(double x, double y, double z, double mass,
3966+
unsigned long int indx)
39533967
{
39543968
// Truncate to cube with sides in [0,1]
39553969
if (x<0.0)
@@ -4179,7 +4193,7 @@ namespace BasisClasses
41794193
}
41804194

41814195
if (use) {
4182-
accumulate(pp(0), pp(1), pp(2), p->mass);
4196+
accumulate(pp(0), pp(1), pp(2), p->mass, p->indx);
41834197
}
41844198
}
41854199
make_coefs();
@@ -4298,7 +4312,7 @@ namespace BasisClasses
42984312
coefindx++;
42994313

43004314
if (use) {
4301-
accumulate(pp(0), pp(1), pp(2), m(n));
4315+
accumulate(pp(0), pp(1), pp(2), m(n), n);
43024316
}
43034317
}
43044318
}
@@ -4336,7 +4350,7 @@ namespace BasisClasses
43364350
coefindx++;
43374351

43384352
if (use) {
4339-
accumulate(pp(0), pp(1), pp(2), m(n));
4353+
accumulate(pp(0), pp(1), pp(2), m(n), n);
43404354
}
43414355
}
43424356
}
@@ -4826,6 +4840,7 @@ namespace BasisClasses
48264840

48274841
// Add a time attribute
48284842
//
4843+
time = roundTime(time);
48294844
stanza.createAttribute<double>("Time", HighFive::DataSpace::From(time)).write(time);
48304845

48314846
// Add the sample statistics
@@ -4903,14 +4918,22 @@ namespace BasisClasses
49034918
return;
49044919
}
49054920

4921+
// Round time
4922+
//
4923+
time = roundTime(time);
4924+
49064925
// The H5 filename
4926+
//
49074927
std::string fname = "coefcovar." + compname + "." + runtag + ".h5";
49084928

49094929
// Check if file exists?
49104930
//
49114931
try {
4912-
// Open the HDF5 file in read-only mode
4913-
HighFive::File file(fname, HighFive::File::ReadOnly);
4932+
// Open the HDF5 file in read-write mode, creating if it doesn't
4933+
// exist
4934+
HighFive::File file(fname,
4935+
HighFive::File::ReadWrite |
4936+
HighFive::File::Create);
49144937

49154938
// Check for version string
49164939
std::string path = "CovarianceFileVersion";
@@ -4921,23 +4944,6 @@ namespace BasisClasses
49214944
return;
49224945
}
49234946

4924-
} catch (const HighFive::Exception& err) {
4925-
// Handle HighFive specific errors (e.g., file not found)
4926-
std::cerr << "HighFive Error: " << err.what() << std::endl;
4927-
} catch (const std::exception& err) {
4928-
// Handle other general exceptions
4929-
std::cerr << "Error: " << err.what() << std::endl;
4930-
}
4931-
4932-
// Write coefficients
4933-
//
4934-
try {
4935-
// Create a new hdf5 file
4936-
//
4937-
HighFive::File file(fname,
4938-
HighFive::File::ReadWrite |
4939-
HighFive::File::Create);
4940-
49414947
// Write the Version string
49424948
//
49434949
file.createAttribute<std::string>("CovarianceFileVersion", HighFive::DataSpace::From(CovarianceFileVersion)).write(CovarianceFileVersion);
@@ -4967,8 +4973,14 @@ namespace BasisClasses
49674973
//
49684974
dataset.write(count);
49694975

4970-
} catch (HighFive::Exception& err) {
4971-
std::cerr << err.what() << std::endl;
4976+
} catch (const HighFive::Exception& err) {
4977+
// Handle HighFive specific errors (e.g., file not found)
4978+
throw std::runtime_error
4979+
(std::string("BiorthBasis::writeCoefCovariance HighFive Error: ") + err.what());
4980+
} catch (const std::exception& err) {
4981+
// Handle other general exceptions
4982+
throw std::runtime_error
4983+
(std::string("BiorthBasis::writeCoefCovariance Error: ") + err.what());
49724984
}
49734985
}
49744986

@@ -4996,7 +5008,8 @@ namespace BasisClasses
49965008
dataset.write(count);
49975009

49985010
} catch (HighFive::Exception& err) {
4999-
std::cerr << err.what() << std::endl;
5011+
throw std::runtime_error
5012+
(std::string("BiorthBasis::extendCoefCovariance: HighFive error: ") + err.what());
50005013
}
50015014
}
50025015

‎expui/FieldBasis.H‎

Lines changed: 5 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -156,12 +156,13 @@ namespace BasisClasses
156156
(Eigen::VectorXd& m, RowMatrixXd& p, bool roundrobin=true, bool posvelrows=false);
157157

158158
//! Accumulate new coefficients
159-
virtual void accumulate(double mass,
160-
double x, double y, double z,
161-
double u, double v, double w);
159+
virtual void accumulate(double x, double y, double z,
160+
double u, double v, double w,
161+
double mass, unsigned long int indx = 0);
162162

163163
[[deprecated("not relevant for this class")]]
164-
virtual void accumulate(double x, double y, double z, double mass)
164+
virtual void accumulate(double x, double y, double z, double mass,
165+
unsigned long indx = 0)
165166
{
166167
}
167168

0 commit comments

Comments
 (0)