You signed in with another tab or window. Reload to refresh your session.You signed out in another tab or window. Reload to refresh your session.You switched accounts on another tab or window. Reload to refresh your session.Dismiss alert
{{ message }}
Repository navigation
Commit 1ce87f3
Browse filesBrowse the repository at this point in the historyBrowse files
TPC SC: review fixes and further fixes of the space-charge calculation
Fixes from the review and of issues found afterwards. The distortions and
corrections of the default calculation are unchanged (bit-identical), apart
from the sub-um fix of the cartesian getCorrections() (see Other).
File IO
- no std::string_view::data() for functions expecting a C string: file
names, tree names and options are copied to a std::string before they are
passed to TFile, TString, gSystem and the column lists of RDataFrame
(writeToFile, initFromFile, dumpSlice, dumpInterpolation, dumpToFile,
dumpMetaData, readMetaData, DataContainer3D::getVertices,
setAnalyticalCorrectionsDistortionsFromFile, calculateElectronDriftPath).
A file name given as a view into a longer string created a second file
with the longer name in writeToFile() with option UPDATE.
- the tree names of the set*FromFile() functions are std::string: they were
string_views to the temporary returned by fmt::format(), used after the
end of its lifetime in checkGridFromFile() and initFromFile().
Poisson solver
- poissonSolver(side) and poissonSolver() use the stopping criterion set for
the PoissonSolver (setConvergenceError()) instead of overwriting it with
the default; the overloads with a stopping criterion set it as before.
Default unchanged: 1e-8 for double, 1e-6 for float.
- the factorisation of the periodic tridiagonal systems of the line
relaxation in phi (Thomas algorithm + Sherman-Morrison) is calculated once
per multi grid level (factorisePhiLines(), PhiLineFactorisation) instead
of for each sweep; the coefficients of a level are set in
setLevelCoefficients(). Same potential (bit-identical). Poisson solver per
side (realistic density, 32 threads): 129x129x1800 1.74 s -> 1.57 s,
257x257x180 and 129x129x180 unchanged.
- the parameter of relaxLinePhi() is not shadowed by a local variable.
Density from histograms
- the fill buffer of the histograms is flushed before the bin array is read:
a histogram with entries still in the buffer gave a zero density.
- fillChargeDensityFromHisto(hA, hC): a different number of bins of the A-
and C-side histograms is an error (the C-side histogram was read outside
of its bin array); the density is not filled and the function returns
false, also if a histogram is not found in the file. Different ranges of
the C-side histogram and a z axis of the A-side histogram not starting at
0 give a warning, the density is unchanged.
- one accessor type for the bin content (bin array of TH3F/TH3D or
GetBinContent() for other classes): 2 instead of 12 instantiations of the
rebinning loop per DataT, same density and time.
Tricubic interpolation
- the combined interpolation of three containers (interpolate(z, r, phi,
other1, other2, ...)) is only used for the same data layout, grid and
extrapolation settings, otherwise the three values are interpolated
separately: a container with the same number of points but a different
shape was read outside of its data at the boundary of the grid. The check
is inline (the combined interpolation is ~1% slower than without the
check). Same distortions and corrections.
Other
- cartesian getCorrections(): the phi in [0, 2 pi) is only used for the look
up, the corrected position is calculated with the phi of the point.
bringTo02PiGen() adds the single precision 2 pi (larger than 2 pi by
1.7e-7), which rotated the corrected position of points with y < 0 by
1.7e-7 rad (up to 0.45 um at the outer radius), also for
SpaceCharge<double>.
Tests: stopping criterion of poissonSolver(), density with a pending fill
buffer, different binning and ranges of the A- and C-side histograms, return
value of fillChargeDensityFromHisto(), combined = separate tricubic
interpolation for different extrapolation settings, grids and shapes.
Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
/// relaxation of the given level: line relaxation in phi (MGLevel::lineRelaxPhi) with the factorisation phiLineFact of the level or point relaxation (relax3D())
/// red-black Gauss-Seidel line relaxation in phi direction: the potential of all phi slices on one line (r, z) is solved exactly (periodic tridiagonal system)
123
141
/// for given potential on the neighbouring lines. The lines with (ir + iz) even and odd are relaxed one after the other.
124
142
/// Robust if the phi direction is strongly coupled compared to r and z (fine grid in phi), where the point relaxation smooths the error badly.
143
+
/// \param phiLineFact factorisation of the periodic tridiagonal systems of the level (see factorisePhiLines())
/// \returns the grids of the multi grid levels. From one level to the next the r and z direction are coarsened and the phi direction is coarsened (if the number of slices is even)
129
148
/// only if it is not weakly coupled compared to r and z (see MGParameters::phiCoarseningThreshold). Line relaxation in phi is used on the levels where phi is strongly coupled (see MGParameters::phiLineRelaxThreshold)
@@ -273,6 +292,7 @@ class PoissonSolver
273
292
/// \param nPost number of smoothing after coarsening
274
293
/// \param ratioZ ratio between square of grid r and grid z
275
294
/// \param levels grids of the multi grid levels
295
+
/// \param phiLineFact factorisation of the line relaxation in phi of the levels (empty for the levels without line relaxation)
276
296
/// \param tvArrayV vector of V potential in different grids
277
297
/// \param tvCharge vector of charge distribution in different grids
278
298
/// \param tvResidue vector of residue calculation in different grids
@@ -281,9 +301,9 @@ class PoissonSolver
281
301
/// \param coefficient3 coefficient for relaxation (ratio r/z)
282
302
/// \param coefficient4 coefficient for relaxation (ratio for grid_r)
283
303
/// \param inverseCoefficient4 coefficient for relaxation (inverse coefficient4)
/// step 0: set the space-charge density from two TH3 histograms containing the space-charge density for A and C side seperately
139
+
/// step 0: set the space-charge density from two TH3 histograms containing the space-charge density for A and C side seperately.
140
+
/// Both histograms need the same number of bins (phi, r, z), the C-side histogram is expected in phi and r with the same ranges as the A-side histogram
141
+
/// and in z with the mirrored range (e.g. A side: z in [0, zMax], C side: z in [-zMax, 0]). The z axis of the A-side histogram is expected to start at 0.
142
+
/// Different ranges give a warning: the bins are used with the expected ranges
140
143
/// \param hisSCDensity3D_A histogram for the space charge density for A-side
141
144
/// \param hisSCDensity3D_C histogram for the space charge density for C-side
/// step 0: set the space-charge density from two TH3 histograms containing the space charge density for A and C side separately which are stored in a ROOT file
145
149
/// \param file path to root file containing the space-charge density
146
150
/// \param nameA name of the space-charge density histogram for the A-side
147
151
/// \param nameC name of the space-charge density histogram for the C-side
/// \return true if the density was filled, false if a histogram was not found or the number of bins of the histograms differs (the density is not changed)
/// step 0: set the space charge density from std::vector<CalDet> containing the space charge density. Each entry in the object corresponds to one z slice
151
156
/// \param calSCDensity3D pad-by-pad CalDet for the space charge density
@@ -282,14 +287,25 @@ class SpaceCharge
282
287
/// \param Side of the TPC
283
288
voidsetPotentialBoundaryOuterRadius(const std::function<DataT(DataT)>& potentialFunc, const Side side);
284
289
285
-
/// step 1: use the O2TPCPoissonSolver class to numerically calculate the potential with set space charge density and boundary conditions from potential
290
+
/// step 1: use the O2TPCPoissonSolver class to numerically calculate the potential with set space charge density and boundary conditions from potential.
291
+
/// The stopping criterion of the poisson solver is the one set with PoissonSolver<DataT>::setConvergenceError() (default: PoissonSolver<DataT>::getDefaultConvergenceError())
286
292
/// \param side side of the TPC
287
-
/// \param stoppingConvergence stopping criterion used in the poisson solver (relative change of the potential between two multi grid V cycles, see MGParameters::relativeConvergence; default 1e-8 for double, 1e-6 for float)
/// step 1: use the O2TPCPoissonSolver class to numerically calculate the potential with set space charge density and boundary conditions from potential for A and C side in parallel
291
-
/// \param stoppingConvergence stopping criterion used in the poisson solver (relative change of the potential between two multi grid V cycles, see MGParameters::relativeConvergence; default 1e-8 for double, 1e-6 for float)
/// step 1: use the O2TPCPoissonSolver class to numerically calculate the potential with set space charge density and boundary conditions from potential
296
+
/// \param side side of the TPC
297
+
/// \param stoppingConvergence stopping criterion used in the poisson solver (relative change of the potential between two multi grid V cycles, see MGParameters::relativeConvergence).
298
+
/// It is set with PoissonSolver<DataT>::setConvergenceError() and is also used by later calls
299
+
voidpoissonSolver(const Side side, const DataT stoppingConvergence);
300
+
301
+
/// step 1: use the O2TPCPoissonSolver class to numerically calculate the potential with set space charge density and boundary conditions from potential for the A and C side (one after the other).
302
+
/// The stopping criterion of the poisson solver is the one set with PoissonSolver<DataT>::setConvergenceError() (default: PoissonSolver<DataT>::getDefaultConvergenceError())
303
+
voidpoissonSolver();
304
+
305
+
/// step 1: use the O2TPCPoissonSolver class to numerically calculate the potential with set space charge density and boundary conditions from potential for the A and C side (one after the other)
306
+
/// \param stoppingConvergence stopping criterion used in the poisson solver (relative change of the potential between two multi grid V cycles, see MGParameters::relativeConvergence).
307
+
/// It is set with PoissonSolver<DataT>::setConvergenceError() and is also used by later calls
Copy file name to clipboardExpand all lines: Detectors/TPC/spacecharge/include/TPCSpaceCharge/TriCubic.h
+4-1Lines changed: 4 additions & 1 deletion
Original file line number
Diff line number
Diff line change
@@ -83,7 +83,7 @@ class TriCubicInterpolator
83
83
84
84
/// interpolate the data of this and of two other interpolators at the same coordinate (e.g. the three components of a vector field).
85
85
/// The position in the grid and the weights of the interpolation are calculated only once, which is faster than three separate interpolations.
86
-
/// The other interpolators have to use the same grid. The extrapolation settings of this interpolator are used for all three.
86
+
/// If the data layout, the grid or the extrapolation settings of the other interpolators are different, the three values are interpolated separately.
87
87
/// \param z z coordinate
88
88
/// \param r r coordinate
89
89
/// \param phi phi coordinate
@@ -117,6 +117,9 @@ class TriCubicInterpolator
117
117
ExtrapolationType mExtrapolationType = ExtrapolationType::Parabola; ///< sets which type of extrapolation for missing points at boundary is used
118
118
boolmExtraPolateValues = true; ///< extrapolating values outside of the grid or restricting query to inside of the grid
119
119
120
+
/// \return returns whether the other interpolator (with a different grid object) has the same data layout and grid (number of vertices, minimum and spacing)
0 commit comments