From 7eda25d5025c2624cb74b5f782500d45765b5160 Mon Sep 17 00:00:00 2001 From: Philipp Windischhofer Date: Mon, 24 Jun 2024 21:28:27 -0500 Subject: [PATCH 1/6] try to run clang-format in CI --- .github/workflows/build-ci.yml | 13 +++++++++++++ eisvogel/core/Current.hh | 7 ++++--- 2 files changed, 17 insertions(+), 3 deletions(-) diff --git a/.github/workflows/build-ci.yml b/.github/workflows/build-ci.yml index bcc441647..ca125569e 100644 --- a/.github/workflows/build-ci.yml +++ b/.github/workflows/build-ci.yml @@ -3,6 +3,19 @@ name: tests on: [push, pull_request, workflow_dispatch] jobs: + + code_format: + name: "Check code formatting" + runs-on: ubuntu-22.04 + steps: + + - name: Checkout Eisvogel repository + uses: actions/checkout@v4 + + - name: Run clang-format + run: | + find . -iname '*.hxx' -o -iname '*.hh' -o -iname '*.cpp' -o -iname '*.h' | xargs clang-format --dry-run --Werror --style=file + build_and_test: name: "Run build and tests" runs-on: ubuntu-22.04 diff --git a/eisvogel/core/Current.hh b/eisvogel/core/Current.hh index 058b735ba..5fc498f36 100644 --- a/eisvogel/core/Current.hh +++ b/eisvogel/core/Current.hh @@ -3,14 +3,15 @@ #include "Vector.hh" struct LineCurrentSegment { - - LineCurrentSegment(XYZCoordVector&& start_pos, XYZCoordVector&& end_pos, scalar_t start_time, scalar_t end_time, scalar_t charge); + + LineCurrentSegment(XYZCoordVector&& start_pos, XYZCoordVector&& end_pos, + scalar_t start_time, scalar_t end_time, scalar_t charge); XYZCoordVector start_pos; XYZCoordVector end_pos; scalar_t start_time; scalar_t end_time; - scalar_t charge; + scalar_t charge; }; #include "Current.hxx" From 585d7a4e2c76814abe2c075739078fc6aa02fd8b Mon Sep 17 00:00:00 2001 From: Philipp Windischhofer Date: Mon, 24 Jun 2024 21:48:10 -0500 Subject: [PATCH 2/6] add .clang-format --- .clang-format | 93 +++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 93 insertions(+) create mode 100644 .clang-format diff --git a/.clang-format b/.clang-format new file mode 100644 index 000000000..8ce06de2e --- /dev/null +++ b/.clang-format @@ -0,0 +1,93 @@ +--- +Language: Cpp + +ColumnLimit: 80 # line width, 0: respect programmer's choice +TabWidth: 8 +UseTab: Never + +AccessModifierOffset: -2 +AlignAfterOpenBracket: Align +AlignConsecutiveAssignments: false +AlignConsecutiveDeclarations: false +AlignEscapedNewlinesLeft: true +AlignOperands: true +AlignTrailingComments: true +AllowAllParametersOfDeclarationOnNextLine: true +AllowShortBlocksOnASingleLine: true +AllowShortCaseLabelsOnASingleLine: false +AllowShortFunctionsOnASingleLine: All +AllowShortIfStatementsOnASingleLine: false +AllowShortLoopsOnASingleLine: true +AlwaysBreakAfterDefinitionReturnType: None +AlwaysBreakAfterReturnType: None +AlwaysBreakBeforeMultilineStrings: true +BinPackArguments: true +BinPackParameters: true +BraceWrapping: + AfterClass: false + AfterControlStatement: false + AfterEnum: false + AfterFunction: false + AfterNamespace: true + AfterObjCDeclaration: false + AfterStruct: false + AfterUnion: false + BeforeCatch: false + BeforeElse: false + IndentBraces: false +BreakBeforeBinaryOperators: None +BreakBeforeBraces: Stroustrup +BreakBeforeTernaryOperators: true +BreakConstructorInitializersBeforeComma: true +# BreakInheritanceListBeforeComma: true + +CommentPragmas: '^ IWYU pragma:' +ConstructorInitializerAllOnOneLineOrOnePerLine: false +ConstructorInitializerIndentWidth: 4 +ContinuationIndentWidth: 4 +Cpp11BracedListStyle: true +DerivePointerAlignment: false +DisableFormat: false +ExperimentalAutoDetectBinPacking: false +ForEachMacros: [ foreach, Q_FOREACH, BOOST_FOREACH ] +IncludeCategories: + - Regex: '^<.*\.h>' + Priority: 1 + - Regex: '^<.*' + Priority: 2 + - Regex: '.*' + Priority: 3 +SortIncludes: false +IndentCaseLabels: true +IndentWidth: 2 +IndentWrappedFunctionNames: false +KeepEmptyLinesAtTheStartOfBlocks: true +MacroBlockBegin: '' +MacroBlockEnd: '' +MaxEmptyLinesToKeep: 1 +NamespaceIndentation: All +ObjCBlockIndentWidth: 2 +ObjCSpaceAfterProperty: false +ObjCSpaceBeforeProtocolList: false +PenaltyBreakBeforeFirstCallParameter: 1 +PenaltyBreakComment: 300 +PenaltyBreakFirstLessLess: 120 +PenaltyBreakString: 1000 +PenaltyExcessCharacter: 1000000 +PenaltyReturnTypeOnItsOwnLine: 200 +PointerAlignment: Left +ReflowComments: true +SpaceAfterCStyleCast: false +SpaceBeforeAssignmentOperators: true +SpaceBeforeParens: ControlStatements +SpaceInEmptyParentheses: false +SpacesBeforeTrailingComments: 1 +SpacesInAngles: false +SpacesInContainerLiterals: true +SpacesInCStyleCastParentheses: false +SpacesInParentheses: false +SpacesInSquareBrackets: false +Standard: Auto + +# If true, always break after the template<...> of a template declaration. +AlwaysBreakTemplateDeclarations: true From 6795e8e77809152f16e25b96cb97ffb8a96ce75f Mon Sep 17 00:00:00 2001 From: Philipp Windischhofer Date: Mon, 24 Jun 2024 22:03:02 -0500 Subject: [PATCH 3/6] local run with clang-format --- .../core/AnalyticGreensFunctionCalculator.hh | 10 +- eisvogel/core/Antenna.hh | 25 +- eisvogel/core/Geometry.hh | 20 +- eisvogel/core/SignalCalculator.hh | 10 +- eisvogel/core/SignalExport.hh | 5 +- .../impl/AnalyticGreensFunctionCalculator.hxx | 164 ++- eisvogel/core/impl/Antenna.hxx | 53 +- eisvogel/core/impl/Cache.hh | 58 +- eisvogel/core/impl/Cache.hxx | 175 ++- eisvogel/core/impl/Current.hxx | 18 +- .../core/impl/DefaultSerializationTraits.hxx | 192 +-- eisvogel/core/impl/DistributedNDVecArray.hh | 347 +++-- eisvogel/core/impl/DistributedNDVecArray.hxx | 1306 ++++++++++------- eisvogel/core/impl/Geometry.hxx | 24 +- eisvogel/core/impl/GreensFunction.hh | 115 +- eisvogel/core/impl/GreensFunction.hxx | 550 ++++--- eisvogel/core/impl/Interpolation.hh | 20 +- eisvogel/core/impl/Interpolation.hxx | 187 +-- eisvogel/core/impl/IteratorUtils.hh | 143 +- eisvogel/core/impl/MathUtils.hh | 13 +- eisvogel/core/impl/MathUtils.hxx | 37 +- eisvogel/core/impl/MemoryUtils.hh | 28 +- eisvogel/core/impl/NDVecArray.hh | 247 ++-- eisvogel/core/impl/NDVecArray.hxx | 426 +++--- eisvogel/core/impl/NDVecArrayOperations.hh | 19 +- eisvogel/core/impl/NDVecArrayOperations.hxx | 81 +- eisvogel/core/impl/NDVecArraySerialization.hh | 277 ++-- eisvogel/core/impl/NDVecArrayStreamer.hh | 57 +- eisvogel/core/impl/NDVecArrayStreamer.hxx | 505 ++++--- eisvogel/core/impl/Quadrature.hh | 36 +- eisvogel/core/impl/Serialization.hh | 32 +- eisvogel/core/impl/Serialization.hxx | 8 +- eisvogel/core/impl/SignalCalculator.hxx | 25 +- eisvogel/core/impl/SignalExport.hxx | 14 +- eisvogel/core/impl/Symmetry.hh | 16 +- eisvogel/core/impl/Symmetry.hxx | 35 +- eisvogel/core/impl/Vector.hh | 414 ++++-- eisvogel/core/impl/Vector.hxx | 114 +- eisvogel/core/impl/VectorView.hh | 65 +- ...MEEPCylindricalGreensFunctionCalculator.hh | 49 +- ...EEPCylindricalGreensFunctionCalculator.hxx | 827 ++++++----- .../01generate_weighting_field.cpp | 14 +- .../02calculate_shower_emission.cpp | 155 +- .../charge_excess_profile.cpp | 96 +- .../environment/ice_profile.cpp | 16 +- .../shower_profile/environment/ice_profile.h | 15 +- extern/shower_profile/main.cpp | 57 +- extern/shower_profile/shower_1D/shower_1D.cpp | 251 ++-- extern/shower_profile/shower_1D/shower_1D.h | 56 +- extern/shower_profile/shower_2D/shower_2D.cpp | 180 +-- extern/shower_profile/shower_2D/shower_2D.h | 54 +- .../shower_creator/shower_creator.cpp | 194 ++- .../shower_creator/shower_creator.h | 82 +- extern/utilities/constants.h | 14 +- extern/utilities/units.h | 173 ++- tests/test_utils/CSVReader.hh | 4 +- tests/test_utils/CSVReader.hxx | 44 +- tests/test_utils/testUtils.hh | 44 +- 58 files changed, 4746 insertions(+), 3450 deletions(-) diff --git a/eisvogel/core/AnalyticGreensFunctionCalculator.hh b/eisvogel/core/AnalyticGreensFunctionCalculator.hh index ebd3507c8..b32851ad2 100644 --- a/eisvogel/core/AnalyticGreensFunctionCalculator.hh +++ b/eisvogel/core/AnalyticGreensFunctionCalculator.hh @@ -5,9 +5,13 @@ namespace GreensFunctionCalculator::Analytic { - void ElectricDipole(std::filesystem::path gf_path, const RZCoordVector& start_coords, const RZCoordVector& end_coords, scalar_t t_end, scalar_t ior, - scalar_t filter_t_peak, unsigned int filter_order, scalar_t r_min, - scalar_t os_factor = 10, std::size_t max_pts_in_chunk = 400); + void ElectricDipole(std::filesystem::path gf_path, + const RZCoordVector& start_coords, + const RZCoordVector& end_coords, scalar_t t_end, + scalar_t ior, scalar_t filter_t_peak, + unsigned int filter_order, scalar_t r_min, + scalar_t os_factor = 10, + std::size_t max_pts_in_chunk = 400); } diff --git a/eisvogel/core/Antenna.hh b/eisvogel/core/Antenna.hh index 422ce99c6..6f5349e02 100644 --- a/eisvogel/core/Antenna.hh +++ b/eisvogel/core/Antenna.hh @@ -8,37 +8,34 @@ class Antenna : public meep::src_time { public: Antenna(scalar_t start_time, scalar_t end_time, scalar_t z_pos, - std::function impulse_response_func); - virtual ~Antenna() { }; - + std::function impulse_response_func); + virtual ~Antenna(){}; + virtual void AddToGeometry(meep::fields& f, CylinderGeometry& geom) const = 0; public: virtual std::complex current(double time, double dt) const; virtual std::complex dipole(double time) const; - - virtual double last_time() const {return m_end_time;} + + virtual double last_time() const { return m_end_time; } protected: scalar_t m_start_time, m_end_time, m_z_pos; - + private: std::function m_impulse_response_func; - }; class InfEDipoleAntenna : public Antenna { - + public: InfEDipoleAntenna(scalar_t start_time, scalar_t end_time, scalar_t z_pos, - std::function impulse_response_func); - virtual ~InfEDipoleAntenna() { }; - + std::function impulse_response_func); + virtual ~InfEDipoleAntenna(){}; + virtual void AddToGeometry(meep::fields& f, CylinderGeometry& geom) const; - virtual meep::src_time* clone() const { - return new InfEDipoleAntenna(*this); - } + virtual meep::src_time* clone() const { return new InfEDipoleAntenna(*this); } }; #include "Antenna.hxx" diff --git a/eisvogel/core/Geometry.hh b/eisvogel/core/Geometry.hh index ac025a65f..64238bec0 100644 --- a/eisvogel/core/Geometry.hh +++ b/eisvogel/core/Geometry.hh @@ -9,32 +9,30 @@ struct RZCoordVector; class CylinderGeometry : public meep::material_function { public: - CylinderGeometry(scalar_t r_max, scalar_t z_min, scalar_t z_max, - std::function eps_func); - - scalar_t GetRMax() {return m_r_max;}; - scalar_t GetZMin() {return m_z_min;}; - scalar_t GetZMax() {return m_z_max;}; + std::function eps_func); + + scalar_t GetRMax() { return m_r_max; }; + scalar_t GetZMin() { return m_z_min; }; + scalar_t GetZMax() { return m_z_max; }; meep::vec toMeepCoords(const RZCoordVector& coords); RZCoordVector toEisvogelCoords(const meep::vec& meep_coords); - + public: - virtual double eps(const meep::vec& pos); - virtual double chi1p1(meep::field_type ft, const meep::vec &r) { + virtual double chi1p1(meep::field_type ft, const meep::vec& r) + { (void)ft; return eps(r); } - + private: scalar_t m_r_max, m_z_min, m_z_max; public: std::function m_eps_func; - }; #include "Geometry.hxx" diff --git a/eisvogel/core/SignalCalculator.hh b/eisvogel/core/SignalCalculator.hh index e8621e135..7bb0c3a0c 100644 --- a/eisvogel/core/SignalCalculator.hh +++ b/eisvogel/core/SignalCalculator.hh @@ -11,15 +11,15 @@ class LineCurrentSegment; class SignalCalculator { public: - - SignalCalculator(const std::filesystem::path& geometry_path, std::size_t cache_depth = 5); - void AccumulateSignal(const LineCurrentSegment& seg, scalar_t t_sig_start, scalar_t t_sig_samp, std::size_t num_samples, std::vector& signal); + SignalCalculator(const std::filesystem::path& geometry_path, + std::size_t cache_depth = 5); + void AccumulateSignal(const LineCurrentSegment& seg, scalar_t t_sig_start, + scalar_t t_sig_samp, std::size_t num_samples, + std::vector& signal); private: - std::filesystem::path m_geometry_path; std::shared_ptr m_gf; - }; #include "SignalCalculator.hxx" diff --git a/eisvogel/core/SignalExport.hh b/eisvogel/core/SignalExport.hh index 8baa4964f..40e132a43 100644 --- a/eisvogel/core/SignalExport.hh +++ b/eisvogel/core/SignalExport.hh @@ -4,7 +4,8 @@ #include #include "Common.hh" -void ExportSignal(std::vector signal_times, std::vector signal_values, - std::string outpath, int time_prec = 3, int value_prec = 15); +void ExportSignal(std::vector signal_times, + std::vector signal_values, std::string outpath, + int time_prec = 3, int value_prec = 15); #include "SignalExport.hxx" diff --git a/eisvogel/core/impl/AnalyticGreensFunctionCalculator.hxx b/eisvogel/core/impl/AnalyticGreensFunctionCalculator.hxx index d8eaa97de..9191b5439 100644 --- a/eisvogel/core/impl/AnalyticGreensFunctionCalculator.hxx +++ b/eisvogel/core/impl/AnalyticGreensFunctionCalculator.hxx @@ -7,60 +7,68 @@ namespace { - void EvaluateElectricDipoleGreensFunction(const RZTCoordVector& start_coords, const RZTCoordVector& stepsize, scalar_t ior, - scalar_t filter_t_peak, unsigned int filter_order, scalar_t r_min, - SpatialSymmetry::Cylindrical::chunk_t& buffer) { + void EvaluateElectricDipoleGreensFunction( + const RZTCoordVector& start_coords, const RZTCoordVector& stepsize, + scalar_t ior, scalar_t filter_t_peak, unsigned int filter_order, + scalar_t r_min, SpatialSymmetry::Cylindrical::chunk_t& buffer) + { scalar_t Qw = 1.0; - scalar_t eps0 = 1.0; // vacuum dielectric constant + scalar_t eps0 = 1.0; // vacuum dielectric constant scalar_t ds = 1.0; - scalar_t c = 1.0; // speed of light in vacuum - + scalar_t c = 1.0; // speed of light in vacuum + auto filtered_theta = [&](scalar_t t) -> scalar_t { - if(t <= 0) { - return 0.0; + if (t <= 0) { + return 0.0; } - return 1.0 - MathUtils::incomplete_gamma(1 + filter_order, filter_order * t / filter_t_peak) / std::exp(std::lgamma(filter_order + 1)); + return 1.0 - MathUtils::incomplete_gamma( + 1 + filter_order, filter_order * t / filter_t_peak) / + std::exp(std::lgamma(filter_order + 1)); }; - + auto filtered_delta = [&](scalar_t t) -> scalar_t { - if(t <= 0) { - return 0.0; + if (t <= 0) { + return 0.0; } - return std::pow(t / filter_t_peak * filter_order, filter_order) * std::exp(-t / filter_t_peak * filter_order) / (filter_t_peak * std::exp(std::lgamma(filter_order))); + return std::pow(t / filter_t_peak * filter_order, filter_order) * + std::exp(-t / filter_t_peak * filter_order) / + (filter_t_peak * std::exp(std::lgamma(filter_order))); }; auto filtered_delta_prime = [&](scalar_t t) -> scalar_t { - if(t <= 0) { - return 0.0; + if (t <= 0) { + return 0.0; } - return filtered_delta(t) * (filter_t_peak - t) * filter_order / (filter_t_peak * t); - }; + return filtered_delta(t) * (filter_t_peak - t) * filter_order / + (filter_t_peak * t); + }; // Weighting field in spherical coordinates auto E_r = [&](scalar_t t, scalar_t r_xy, scalar_t z) -> scalar_t { r_xy = std::fabs(r_xy); scalar_t r = std::sqrt(r_xy * r_xy + z * z); - if(r < r_min) { - return std::nan(""); + if (r < r_min) { + return std::nan(""); } scalar_t t_prop = r * ior / c, t_del = t - t_prop; scalar_t cos_theta = z / r; - return -2.0 * Qw * ds / (eps0 * 4 * M_PI) * cos_theta / std::pow(r, 3) * (filtered_theta(t_del) + - t_prop * filtered_delta(t_del)); + return -2.0 * Qw * ds / (eps0 * 4 * M_PI) * cos_theta / std::pow(r, 3) * + (filtered_theta(t_del) + t_prop * filtered_delta(t_del)); }; auto E_theta = [&](scalar_t t, scalar_t r_xy, scalar_t z) -> scalar_t { r_xy = std::fabs(r_xy); scalar_t r = std::sqrt(r_xy * r_xy + z * z); - if(r < r_min) { - return std::nan(""); + if (r < r_min) { + return std::nan(""); } scalar_t t_prop = r * ior / c, t_del = t - t_prop; scalar_t sin_theta = r_xy / r; - return -Qw * ds / (eps0 * 4 * M_PI) * sin_theta / std::pow(r, 3) * (filtered_theta(t_del) + t_prop * filtered_delta(t_del) - + std::pow(t_prop, 2) * filtered_delta_prime(t_del)); + return -Qw * ds / (eps0 * 4 * M_PI) * sin_theta / std::pow(r, 3) * + (filtered_theta(t_del) + t_prop * filtered_delta(t_del) + + std::pow(t_prop, 2) * filtered_delta_prime(t_del)); }; // Weighting field in cylindrical coordinates @@ -70,19 +78,19 @@ namespace { scalar_t cos_theta = z / r, sin_theta = r_xy / r; return E_r(t, r_xy, z) * sin_theta + E_theta(t, r_xy, z) * cos_theta; }; - + auto E_z = [&](scalar_t t, scalar_t r_xy, scalar_t z) -> scalar_t { r_xy = std::fabs(r_xy); scalar_t r = std::sqrt(r_xy * r_xy + z * z); scalar_t cos_theta = z / r, sin_theta = r_xy / r; return E_r(t, r_xy, z) * cos_theta - E_theta(t, r_xy, z) * sin_theta; }; - - using view_t = typename SpatialSymmetry::Cylindrical::view_t; - auto evaluator = [&](const RZTIndexVector& ind, view_t cur_element) { + using view_t = typename SpatialSymmetry::Cylindrical::view_t; + auto evaluator = [&](const RZTIndexVector& ind, view_t cur_element) { // Coordinates of current point - RZTCoordVector cur_pt = start_coords + stepsize * (ind.template as_type()); + RZTCoordVector cur_pt = + start_coords + stepsize * (ind.template as_type()); // Fields at this location scalar_t cur_E_rxy = E_rxy(cur_pt.t(), cur_pt.r(), cur_pt.z()); @@ -92,46 +100,58 @@ namespace { cur_element[0] = cur_E_rxy; cur_element[1] = cur_E_z; }; - + // evaluate the Green's function and fill the buffer - buffer.index_loop_over_elements(evaluator); + buffer.index_loop_over_elements(evaluator); } -} +} // namespace namespace GreensFunctionCalculator::Analytic { - - void ElectricDipole(std::filesystem::path gf_path, const RZCoordVector& start_coords, const RZCoordVector& end_coords, scalar_t t_end, scalar_t ior, - scalar_t filter_t_peak, unsigned int filter_order, scalar_t r_min, - scalar_t os_factor, std::size_t max_pts_in_chunk) { - + + void ElectricDipole(std::filesystem::path gf_path, + const RZCoordVector& start_coords, + const RZCoordVector& end_coords, scalar_t t_end, + scalar_t ior, scalar_t filter_t_peak, + unsigned int filter_order, scalar_t r_min, + scalar_t os_factor, std::size_t max_pts_in_chunk) + { + // make sure to start from scratch std::filesystem::remove_all(gf_path); - + // compute required step size for sampling of weighting field - scalar_t c = 1.0; // speed of light in vacuum - scalar_t fmax = (scalar_t)(filter_order) / (2 * M_PI * filter_t_peak) * std::sqrt(std::pow(2.0, 1.0 / (filter_order + 1)) - 1); + scalar_t c = 1.0; // speed of light in vacuum + scalar_t fmax = (scalar_t)(filter_order) / (2 * M_PI * filter_t_peak) * + std::sqrt(std::pow(2.0, 1.0 / (filter_order + 1)) - 1); scalar_t lambda_min = c / (fmax * ior); scalar_t delta_t = 1.0 / (2 * fmax * os_factor); scalar_t delta_pos = lambda_min / (2.0 * os_factor); scalar_t t_start = 0.0; - + std::cout << "---------------------------" << std::endl; std::cout << "Using oversampling factor = " << os_factor << std::endl; std::cout << "delta_t = " << delta_t << std::endl; std::cout << "delta_r = delta_z = " << delta_pos << std::endl; - std::cout << "start_coords: t = " << t_start << ", r = " << start_coords.r() << ", z = " << start_coords.z() << std::endl; - std::cout << "end_coords: t = " << t_end << ", r = " << end_coords.r() << ", z = " << end_coords.z() << std::endl; + std::cout << "start_coords: t = " << t_start << ", r = " << start_coords.r() + << ", z = " << start_coords.z() << std::endl; + std::cout << "end_coords: t = " << t_end << ", r = " << end_coords.r() + << ", z = " << end_coords.z() << std::endl; std::cout << "---------------------------" << std::endl; - RZTCoordVector start_coords_rzt{start_coords.r(), start_coords.z(), t_start}; + RZTCoordVector start_coords_rzt{start_coords.r(), start_coords.z(), + t_start}; RZTCoordVector end_coords_rzt{end_coords.r(), end_coords.z(), t_end}; - // Determine step size so that an integer number of samples fits into the domain of the Green's function + // Determine step size so that an integer number of samples fits into the + // domain of the Green's function RZTVector stepsize_requested{delta_pos, delta_pos, delta_t}; - RZTVector number_pts = ((end_coords_rzt - start_coords_rzt) / stepsize_requested).template as_type(); - RZTVector stepsize = (end_coords_rzt - start_coords_rzt) / number_pts.template as_type(); - + RZTVector number_pts = + ((end_coords_rzt - start_coords_rzt) / stepsize_requested) + .template as_type(); + RZTVector stepsize = (end_coords_rzt - start_coords_rzt) / + number_pts.template as_type(); + // Resulting start- and end indices RZTIndexVector start_inds(0); RZTIndexVector end_inds = number_pts + 1; @@ -139,40 +159,52 @@ namespace GreensFunctionCalculator::Analytic { // Prepare chunk buffer: need 3-dim array storing 2-dim vectors using chunk_t = typename SpatialSymmetry::Cylindrical::chunk_t; using darr_t = typename SpatialSymmetry::Cylindrical::darr_t; - + RZTVector chunk_size(max_pts_in_chunk); chunk_t chunk_buffer(chunk_size); - + // Prepare distributed array - std::size_t cache_depth = 5; // all chunks are prepared as final, no need for a large cache here + std::size_t cache_depth = + 5; // all chunks are prepared as final, no need for a large cache here RZTVector init_cache_el_shape = chunk_size; - RZTVector streamer_chunk_shape(stor::INFTY); streamer_chunk_shape[0] = 1; // serialize one outermost slice at a time - darr_t darr(gf_path, cache_depth, init_cache_el_shape, streamer_chunk_shape); + RZTVector streamer_chunk_shape(stor::INFTY); + streamer_chunk_shape[0] = 1; // serialize one outermost slice at a time + darr_t darr(gf_path, cache_depth, init_cache_el_shape, + streamer_chunk_shape); // Loop over chunks, fill them, and register them in the distributed array - auto fill_and_register_chunk = [&](const RZTIndexVector& chunk_start_ind, const RZTIndexVector& chunk_end_ind){ - + auto fill_and_register_chunk = [&](const RZTIndexVector& chunk_start_ind, + const RZTIndexVector& chunk_end_ind) { // Prepare the start- and end coordinates of this chunk ... RZTVector cur_chunk_size = chunk_end_ind - chunk_start_ind; - RZTCoordVector chunk_start_coords = start_coords_rzt + chunk_start_ind.template as_type() * stepsize; + RZTCoordVector chunk_start_coords = + start_coords_rzt + + chunk_start_ind.template as_type() * stepsize; // ... and make sure the chunk buffer matches the required shape chunk_buffer.resize(cur_chunk_size); // Fill the buffer ... - EvaluateElectricDipoleGreensFunction(chunk_start_coords, stepsize, ior, filter_t_peak, filter_order, r_min, chunk_buffer); + EvaluateElectricDipoleGreensFunction(chunk_start_coords, stepsize, ior, + filter_t_peak, filter_order, r_min, + chunk_buffer); // ... and register it darr.RegisterChunk(chunk_buffer, chunk_start_ind); - }; - IteratorUtils::index_loop_over_chunks(start_inds, end_inds, chunk_size, fill_and_register_chunk); - - // Reshape the chunks to make sure the overlap required for the interpolation is respected + }; + IteratorUtils::index_loop_over_chunks(start_inds, end_inds, chunk_size, + fill_and_register_chunk); + + // Reshape the chunks to make sure the overlap required for the + // interpolation is respected std::size_t overlap = 2; std::filesystem::path workdir_tmp = "./darr_test_tmp"; - darr.RebuildChunks(chunk_size, workdir_tmp, overlap, SpatialSymmetry::Cylindrical::boundary_evaluator); - + darr.RebuildChunks( + chunk_size, workdir_tmp, overlap, + SpatialSymmetry::Cylindrical::boundary_evaluator); + // Create the actual Green's function from the sampled data - CylindricalGreensFunction(start_coords_rzt, end_coords_rzt, stepsize, std::move(darr)); + CylindricalGreensFunction(start_coords_rzt, end_coords_rzt, stepsize, + std::move(darr)); } -} +} // namespace GreensFunctionCalculator::Analytic diff --git a/eisvogel/core/impl/Antenna.hxx b/eisvogel/core/impl/Antenna.hxx index 10c9e93f5..6a2cd8ef2 100644 --- a/eisvogel/core/impl/Antenna.hxx +++ b/eisvogel/core/impl/Antenna.hxx @@ -1,41 +1,54 @@ #include "Vector.hh" Antenna::Antenna(scalar_t start_time, scalar_t end_time, scalar_t z_pos, - std::function impulse_response_func) : - m_start_time(start_time), m_end_time(end_time), m_z_pos(z_pos), m_impulse_response_func(impulse_response_func) { + std::function impulse_response_func) + : m_start_time(start_time) + , m_end_time(end_time) + , m_z_pos(z_pos) + , m_impulse_response_func(impulse_response_func) +{ - // This is important: `is_integrated = false` specifies that MEEP expects us to - // provide the feed current and not the dipole moment of the antenna + // This is important: `is_integrated = false` specifies that MEEP expects us + // to provide the feed current and not the dipole moment of the antenna is_integrated = false; }; -std::complex Antenna::current(double time, [[maybe_unused]] double dt) const { - +std::complex Antenna::current(double time, + [[maybe_unused]] double dt) const +{ + float rtime = float(time); - if(rtime >= m_start_time && rtime <= m_end_time) { + if (rtime >= m_start_time && rtime <= m_end_time) { return m_impulse_response_func(time); - } else { + } + else { return 0.0; - } + } }; - std::complex Antenna::dipole([[maybe_unused]] double time) const { - // Because we set `is_integrated = false`, no need to provide anything here. - // Return NaN instead of 0.0 (or any other value) to make it obvious if this makes its way into any downstream calculations - return std::numeric_limits::quiet_NaN(); +std::complex Antenna::dipole([[maybe_unused]] double time) const +{ + // Because we set `is_integrated = false`, no need to provide anything here. + // Return NaN instead of 0.0 (or any other value) to make it obvious if this + // makes its way into any downstream calculations + return std::numeric_limits::quiet_NaN(); }; -InfEDipoleAntenna::InfEDipoleAntenna(scalar_t start_time, scalar_t end_time, scalar_t z_pos, - std::function impulse_response_func) : - Antenna(start_time, end_time, z_pos, impulse_response_func) { }; +InfEDipoleAntenna::InfEDipoleAntenna( + scalar_t start_time, scalar_t end_time, scalar_t z_pos, + std::function impulse_response_func) + : Antenna(start_time, end_time, z_pos, impulse_response_func){}; -void InfEDipoleAntenna::AddToGeometry(meep::fields& f, CylinderGeometry& geom) const { +void InfEDipoleAntenna::AddToGeometry(meep::fields& f, + CylinderGeometry& geom) const +{ - // According to workaround described in https://github.com/NanoComp/meep/issues/2704, need to place source at r = 1.5 * delta_R, - // where delta_R is the radial voxel extent + // According to workaround described in + // https://github.com/NanoComp/meep/issues/2704, need to place source at r + // = 1.5 * delta_R, where delta_R is the radial voxel extent scalar_t r_pos = 1.5 / f.a; RZCoordVector antenna_pos{r_pos, m_z_pos}; - + // orientation hard-coded for now ... need to change f.add_point_source(meep::Ez, *this, geom.toMeepCoords(antenna_pos)); } diff --git a/eisvogel/core/impl/Cache.hh b/eisvogel/core/impl/Cache.hh index 11c94100f..96960696e 100644 --- a/eisvogel/core/impl/Cache.hh +++ b/eisvogel/core/impl/Cache.hh @@ -5,17 +5,16 @@ template class Cache { - -public: - - template - Cache(std::size_t depth, PayloadArgs&& ... args); + +public: + template + Cache(std::size_t depth, PayloadArgs&&... args); bool has_free_slot(); bool contains(const IndexT& elem); std::vector contained_elements(); - + // fast cache entry lookup PayloadT& get(const IndexT& elem); @@ -23,41 +22,47 @@ public: template void insert_no_overwrite(const IndexT& index, const DataT& payload); - // reserves a new element in the cache and returns a reference to that clients can directly write to it + // reserves a new element in the cache and returns a reference to that clients + // can directly write to it PayloadT& insert_ref_no_overwrite(const IndexT& index); - // removes this element from the cache, creating a free slot and returning a reference to the element at this location - // cache slot can be overwritten at the next insert call - // is now the responsibility of the caller to do anything that needs to be done to not lose the information + // removes this element from the cache, creating a free slot and returning a + // reference to the element at this location cache slot can be overwritten at + // the next insert call is now the responsibility of the caller to do anything + // that needs to be done to not lose the information PayloadT& evict(const IndexT& elem); - - // assumes that the cache is full, evicts the oldest entry and returns it so that - // it can be properly descoped - // cache slot can be overwritten at the next insert call - // is now the responsibility of the caller to do anything that needs to be done to not lose the information + + // assumes that the cache is full, evicts the oldest entry and returns it so + // that it can be properly descoped cache slot can be overwritten at the next + // insert call is now the responsibility of the caller to do anything that + // needs to be done to not lose the information PayloadT& evict_oldest_from_full_cache(); - + template void loop_over_elements_old_to_new(CallableT&& worker); template void loop_over_elements_new_to_old(CallableT&& worker); - -private: +private: struct CacheElement { - - template - CacheElement(PayloadConstructorArgs&& ... args) : - payload(std::forward(args)...), occupied(false), next(nullptr), prev(nullptr) { } - + + template + CacheElement(PayloadConstructorArgs&&... args) + : payload(std::forward(args)...) + , occupied(false) + , next(nullptr) + , prev(nullptr) + { + } + PayloadT payload; IndexT index; bool occupied; - + CacheElement* next; CacheElement* prev; - }; + }; // moves to the "oldestmost" position of the list void mark_as_oldest(CacheElement* element); @@ -66,9 +71,8 @@ private: void mark_as_newest(CacheElement* element); private: - std::size_t m_depth; - + std::unordered_map m_accessor; std::vector m_storage; CacheElement* m_oldest; diff --git a/eisvogel/core/impl/Cache.hxx b/eisvogel/core/impl/Cache.hxx index dffd00005..eb42ca398 100644 --- a/eisvogel/core/impl/Cache.hxx +++ b/eisvogel/core/impl/Cache.hxx @@ -1,27 +1,30 @@ #include template -template -Cache::Cache(std::size_t depth, PayloadConstructorArgs&& ... args) : m_depth(depth) { - +template +Cache::Cache(std::size_t depth, + PayloadConstructorArgs&&... args) + : m_depth(depth) +{ + std::cout << "building cache with depth " << depth << std::endl; // allocate all required cache elements - m_storage.reserve(depth); - for(std::size_t ind = 0; ind < depth; ind++) { + m_storage.reserve(depth); + for (std::size_t ind = 0; ind < depth; ind++) { m_storage.emplace_back(std::forward(args)...); } // link them together CacheElement* prev_element = nullptr; - for(CacheElement& cur_element : m_storage) { - + for (CacheElement& cur_element : m_storage) { + // link current entry to previous one ... cur_element.prev = prev_element; // ... and previous one to this one - if(prev_element != nullptr) { - prev_element -> next = &cur_element; + if (prev_element != nullptr) { + prev_element->next = &cur_element; } prev_element = &cur_element; @@ -34,185 +37,197 @@ Cache::Cache(std::size_t depth, PayloadConstructorArgs&& ... a template template -void Cache::insert_no_overwrite(const IndexT& index, const DataT& payload) { +void Cache::insert_no_overwrite(const IndexT& index, + const DataT& payload) +{ PayloadT& insert_location = insert_ref_no_overwrite(index); insert_location = payload; - } template -PayloadT& Cache::insert_ref_no_overwrite(const IndexT& index) { +PayloadT& Cache::insert_ref_no_overwrite(const IndexT& index) +{ // an element with the new index must not already exist assert(!contains(index)); - + // try to insert new element into the oldest slot CacheElement* insert_location = m_oldest; - assert(!(insert_location -> occupied)); + assert(!(insert_location->occupied)); // prepare the new payload ... - insert_location -> index = index; - insert_location -> occupied = true; + insert_location->index = index; + insert_location->occupied = true; // .. which now becomes the most-recently changed element in the cache mark_as_newest(insert_location); // update the index mapping m_accessor[index] = insert_location; - - return insert_location -> payload; + + return insert_location->payload; } template -bool Cache::has_free_slot() { +bool Cache::has_free_slot() +{ return m_accessor.size() < m_depth; } template -bool Cache::contains(const IndexT& elem) { +bool Cache::contains(const IndexT& elem) +{ return m_accessor.contains(elem); } template -std::vector Cache::contained_elements() { - +std::vector Cache::contained_elements() +{ + std::vector elements; elements.reserve(m_accessor.size()); - - for (auto& [index, elem]: m_accessor) { - elements.push_back(index); - } + + for (auto& [index, elem] : m_accessor) { elements.push_back(index); } return elements; } template -PayloadT& Cache::get(const IndexT& elem) { - +PayloadT& Cache::get(const IndexT& elem) +{ + CacheElement* accessed_element = m_accessor.at(elem); // mark this as the most-recently accessed element mark_as_newest(accessed_element); - - return accessed_element -> payload; + + return accessed_element->payload; } template -PayloadT& Cache::evict_oldest_from_full_cache() { +PayloadT& Cache::evict_oldest_from_full_cache() +{ - // this function assumes that the cache is full, i.e. that the oldest cache slot is actually occupied + // this function assumes that the cache is full, i.e. that the oldest cache + // slot is actually occupied assert(!has_free_slot()); - - // (don't need to move it, is already at the location of the oldest cache element) + + // (don't need to move it, is already at the location of the oldest cache + // element) CacheElement* evicted_entry = m_oldest; - evicted_entry -> occupied = false; + evicted_entry->occupied = false; // remove its entry from the index mapping - m_accessor.erase(evicted_entry -> index); - - return evicted_entry -> payload; -} + m_accessor.erase(evicted_entry->index); + + return evicted_entry->payload; +} template -PayloadT& Cache::evict(const IndexT& elem) { +PayloadT& Cache::evict(const IndexT& elem) +{ // this function assumes that the element actually exists in the cache assert(contains(elem)); - + CacheElement* evicted_entry = m_accessor.at(elem); - evicted_entry -> occupied = false; + evicted_entry->occupied = false; // move it to the insert location at the "oldest" side of the cache mark_as_oldest(evicted_entry); // remove its entry from the index mapping - m_accessor.erase(evicted_entry -> index); + m_accessor.erase(evicted_entry->index); - return evicted_entry -> payload; + return evicted_entry->payload; } // ------ template -void Cache::mark_as_oldest(CacheElement* element) { +void Cache::mark_as_oldest(CacheElement* element) +{ // this is already the oldest element, nothing needs to be done - if(element -> prev == nullptr) { + if (element->prev == nullptr) { return; } // connect the elements to the left and to the right of this one - element -> prev -> next = element -> next; - - if(element -> next != nullptr) { - element -> next -> prev = element -> prev; + element->prev->next = element->next; + + if (element->next != nullptr) { + element->next->prev = element->prev; } else { // `m_newest` points to `element`, make sure to update - m_newest = element -> prev; + m_newest = element->prev; } // move the element to the front - element -> prev = nullptr; - element -> next = m_oldest; - m_oldest -> prev = element; + element->prev = nullptr; + element->next = m_oldest; + m_oldest->prev = element; // reseat pointer to the new element at the front of the queue m_oldest = element; } template -void Cache::mark_as_newest(CacheElement* element) { +void Cache::mark_as_newest(CacheElement* element) +{ // this is already the newest element, nothing needs to be done - if(element -> next == nullptr) { + if (element->next == nullptr) { return; } // connect the elements to the left and to the right of this one - element -> next -> prev = element -> prev; + element->next->prev = element->prev; - if(element -> prev != nullptr) { - element -> prev -> next = element -> next; + if (element->prev != nullptr) { + element->prev->next = element->next; } - else { + else { // `m_oldest` points to `element`, make sure to update - m_oldest = element -> next; + m_oldest = element->next; } // move the element to the back - element -> next = nullptr; - element -> prev = m_newest; - m_newest -> next = element; + element->next = nullptr; + element->prev = m_newest; + m_newest->next = element; - // reseat pointer to the new element at the back of the queue - m_newest = element; + // reseat pointer to the new element at the back of the queue + m_newest = element; } template template -void Cache::loop_over_elements_old_to_new(CallableT&& worker) { +void Cache::loop_over_elements_old_to_new(CallableT&& worker) +{ CacheElement* cur_element = m_oldest; - while(cur_element != nullptr) { - - if(cur_element -> occupied) { - worker(cur_element -> payload); + while (cur_element != nullptr) { + + if (cur_element->occupied) { + worker(cur_element->payload); } - cur_element = cur_element -> next; - } + cur_element = cur_element->next; + } } template template -void Cache::loop_over_elements_new_to_old(CallableT&& worker) { +void Cache::loop_over_elements_new_to_old(CallableT&& worker) +{ CacheElement* cur_element = m_newest; - while(cur_element != nullptr) { + while (cur_element != nullptr) { - if(cur_element -> occupied) { - worker(cur_element -> payload); - } - cur_element = cur_element -> prev; - } + if (cur_element->occupied) { + worker(cur_element->payload); + } + cur_element = cur_element->prev; + } } diff --git a/eisvogel/core/impl/Current.hxx b/eisvogel/core/impl/Current.hxx index 4deadcddb..14fa8477d 100644 --- a/eisvogel/core/impl/Current.hxx +++ b/eisvogel/core/impl/Current.hxx @@ -1,7 +1,17 @@ -LineCurrentSegment::LineCurrentSegment(XYZCoordVector&& start_pos, XYZCoordVector&& end_pos, scalar_t start_time, scalar_t end_time, scalar_t charge) : - start_pos(start_pos), end_pos(end_pos), start_time(start_time), end_time(end_time), charge(charge) { +LineCurrentSegment::LineCurrentSegment(XYZCoordVector&& start_pos, + XYZCoordVector&& end_pos, + scalar_t start_time, scalar_t end_time, + scalar_t charge) + : start_pos(start_pos) + , end_pos(end_pos) + , start_time(start_time) + , end_time(end_time) + , charge(charge) +{ - if(start_time >= end_time) { - throw std::runtime_error("Error: have line current segment that spans zero or negative time period!"); + if (start_time >= end_time) { + throw std::runtime_error( + "Error: have line current segment that spans zero or negative time " + "period!"); } } diff --git a/eisvogel/core/impl/DefaultSerializationTraits.hxx b/eisvogel/core/impl/DefaultSerializationTraits.hxx index 64c3f035a..d116df011 100644 --- a/eisvogel/core/impl/DefaultSerializationTraits.hxx +++ b/eisvogel/core/impl/DefaultSerializationTraits.hxx @@ -8,154 +8,173 @@ namespace stor { struct Traits { using type = uint32_t; using ser_type = type; - - static void serialize(std::iostream& stream, const type& val) { + + static void serialize(std::iostream& stream, const type& val) + { ser_type outval = htonl(val); stream.write((char*)&outval, sizeof(outval)); } - - static type deserialize(std::iostream& stream) { + + static type deserialize(std::iostream& stream) + { ser_type inval; stream.read((char*)&inval, sizeof(inval)); return ntohl(inval); } }; - + template <> struct Traits { using type = uint64_t; using ser_type = type; - - static void serialize(std::iostream& stream, const type& val) { + + static void serialize(std::iostream& stream, const type& val) + { ser_type outval = htonll(val); stream.write((char*)&outval, sizeof(outval)); } - - static type deserialize(std::iostream& stream) { + + static type deserialize(std::iostream& stream) + { ser_type inval; stream.read((char*)&inval, sizeof(inval)); return ntohll(inval); } }; - + template <> struct Traits { using type = float; using ser_type = uint32_t; - - static void serialize(std::iostream& stream, const type& val) { + + static void serialize(std::iostream& stream, const type& val) + { const ser_type ser_val = reinterpret_cast(val); Traits::serialize(stream, ser_val); } - - static type deserialize(std::iostream& stream) { + + static type deserialize(std::iostream& stream) + { ser_type ser_val = Traits::deserialize(stream); type retval = 0.0f; std::memcpy(&retval, &ser_val, sizeof(ser_val)); return retval; } }; - - template<> + + template <> struct Traits { using type = double; using ser_type = uint64_t; - - static void serialize(std::iostream& stream, const type& val) { + + static void serialize(std::iostream& stream, const type& val) + { const ser_type ser_val = reinterpret_cast(val); Traits::serialize(stream, ser_val); } - - static type deserialize(std::iostream& stream) { + + static type deserialize(std::iostream& stream) + { ser_type ser_val = Traits::deserialize(stream); type retval = 0.0; std::memcpy(&retval, &ser_val, sizeof(ser_val)); return retval; } }; - - template<> + + template <> struct Traits> { using type = std::vector; using ser_type = uint32_t; - - static void serialize(std::iostream& stream, const type& val, std::size_t block_size = 10000) { + + static void serialize(std::iostream& stream, const type& val, + std::size_t block_size = 10000) + { std::size_t vec_size = val.size(); Traits::serialize(stream, vec_size); std::size_t vec_ind = 0; - while(vec_ind < vec_size) { - std::vector outbuf(std::min(vec_size - vec_ind, block_size)); - for(std::size_t buf_ind = 0; (buf_ind < block_size) && (vec_ind < vec_size); buf_ind++, vec_ind++) { - ser_type ser_val = reinterpret_cast(val[vec_ind]); - outbuf[buf_ind] = htonl(ser_val); - } - stream.write((char*)outbuf.data(), sizeof(outbuf[0]) * outbuf.size()); + while (vec_ind < vec_size) { + std::vector outbuf(std::min(vec_size - vec_ind, block_size)); + for (std::size_t buf_ind = 0; + (buf_ind < block_size) && (vec_ind < vec_size); + buf_ind++, vec_ind++) { + ser_type ser_val = reinterpret_cast(val[vec_ind]); + outbuf[buf_ind] = htonl(ser_val); + } + stream.write((char*)outbuf.data(), sizeof(outbuf[0]) * outbuf.size()); } } - - static type deserialize(std::iostream& stream, std::size_t block_size = 10000) { + + static type deserialize(std::iostream& stream, + std::size_t block_size = 10000) + { std::size_t vec_size = Traits::deserialize(stream); std::size_t vec_ind = 0; type val; - while(vec_ind < vec_size) { - std::vector inbuf(std::min(vec_size - vec_ind, block_size)); - stream.read((char*)inbuf.data(), sizeof(inbuf[0]) * inbuf.size()); - for(std::size_t buf_ind = 0; (buf_ind < block_size) && (vec_ind < vec_size); buf_ind++, vec_ind++) { - ser_type ser_val = ntohl(inbuf[buf_ind]); - float cur_val = 0.0f; - std::memcpy(&cur_val, &ser_val, sizeof(ser_val)); - val.push_back(cur_val); - } + while (vec_ind < vec_size) { + std::vector inbuf(std::min(vec_size - vec_ind, block_size)); + stream.read((char*)inbuf.data(), sizeof(inbuf[0]) * inbuf.size()); + for (std::size_t buf_ind = 0; + (buf_ind < block_size) && (vec_ind < vec_size); + buf_ind++, vec_ind++) { + ser_type ser_val = ntohl(inbuf[buf_ind]); + float cur_val = 0.0f; + std::memcpy(&cur_val, &ser_val, sizeof(ser_val)); + val.push_back(cur_val); + } } return val; } }; - + // template // struct Traits> { - + // }; - + // For general vectors template struct Traits> { using type = std::vector; - - static void serialize(std::iostream& stream, const type& val) { + + static void serialize(std::iostream& stream, const type& val) + { Traits::serialize(stream, val.size()); - for(const T& cur : val) { - Traits::serialize(stream, cur); - } + for (const T& cur : val) { Traits::serialize(stream, cur); } } - - static type deserialize(std::iostream& stream) { + + static type deserialize(std::iostream& stream) + { type val; std::size_t size = Traits::deserialize(stream); - // TODO: speed this up by deserializing the entire array instead of values one by one - for(std::size_t ind = 0; ind < size; ind++) { - val.push_back(Traits::deserialize(stream)); + // TODO: speed this up by deserializing the entire array instead of values + // one by one + for (std::size_t ind = 0; ind < size; ind++) { + val.push_back(Traits::deserialize(stream)); } return val; } }; - + // For general arrays template struct Traits> { using type = std::array; - - static void serialize(std::iostream& stream, const type& val) { - // TODO: speed this up by serializing the entire array instead of values one by one - for(const T& cur : val) { - Traits::serialize(stream, cur); - } + + static void serialize(std::iostream& stream, const type& val) + { + // TODO: speed this up by serializing the entire array instead of values + // one by one + for (const T& cur : val) { Traits::serialize(stream, cur); } } - - static type deserialize(std::iostream& stream) { + + static type deserialize(std::iostream& stream) + { type val; - // TODO: speed this up by serializing the entire array instead of values one by one - for(std::size_t ind = 0; ind < n; ind++) { - val[ind] = Traits::deserialize(stream); + // TODO: speed this up by serializing the entire array instead of values + // one by one + for (std::size_t ind = 0; ind < n; ind++) { + val[ind] = Traits::deserialize(stream); } return val; } @@ -167,53 +186,58 @@ namespace stor { struct Traits> { using type = std::map; - static void serialize(std::iostream& stream, const type& val) { + static void serialize(std::iostream& stream, const type& val) + { std::vector keys; std::vector values; - + for (auto const& [key, val] : val) { - keys.push_back(key); - values.push_back(val); + keys.push_back(key); + values.push_back(val); } Traits>::serialize(stream, keys); Traits>::serialize(stream, values); } - static type deserialize(std::iostream& stream) { + static type deserialize(std::iostream& stream) + { std::vector keys = Traits>::deserialize(stream); std::vector values = Traits>::deserialize(stream); - if(keys.size() != values.size()) { - throw std::runtime_error("Error: trying to deserialize malformed std::map!"); + if (keys.size() != values.size()) { + throw std::runtime_error( + "Error: trying to deserialize malformed std::map!"); } type retval; - for(std::size_t ind = 0; ind < keys.size(); ind++) { - retval.insert(retval.end(), std::pair{keys[ind], values[ind]}); + for (std::size_t ind = 0; ind < keys.size(); ind++) { + retval.insert(retval.end(), std::pair{keys[ind], values[ind]}); } - + return retval; } }; - + template <> struct Traits { using type = std::string; - - static void serialize(std::iostream& stream, const type& val) { + + static void serialize(std::iostream& stream, const type& val) + { std::size_t num_chars = val.size(); Traits::serialize(stream, num_chars); stream.write(val.data(), num_chars); } - - static type deserialize(std::iostream& stream) { + + static type deserialize(std::iostream& stream) + { std::size_t num_chars = Traits::deserialize(stream); std::vector string_data(num_chars); stream.read(string_data.data(), num_chars); return std::string(string_data.begin(), string_data.end()); } - }; -} + }; +} // namespace stor diff --git a/eisvogel/core/impl/DistributedNDVecArray.hh b/eisvogel/core/impl/DistributedNDVecArray.hh index 040d4c9ea..a98f15079 100644 --- a/eisvogel/core/impl/DistributedNDVecArray.hh +++ b/eisvogel/core/impl/DistributedNDVecArray.hh @@ -10,10 +10,7 @@ #include "Vector.hh" #include "NDVecArrayStreamer.hh" -enum class ChunkType : std::size_t { - specified = 0, - all_null = 1 -}; +enum class ChunkType : std::size_t { specified = 0, all_null = 1 }; // forward declarations (needed for `operator<<` decleared below) template @@ -31,12 +28,15 @@ struct ChunkMetadata { ChunkMetadata(); // fully-specified constructor - ChunkMetadata(const ChunkType& chunk_type, const std::filesystem::path& filepath, const id_t& chunk_id, - const Vector start_ind, const Vector shape, std::size_t overlap); + ChunkMetadata(const ChunkType& chunk_type, + const std::filesystem::path& filepath, const id_t& chunk_id, + const Vector start_ind, + const Vector shape, std::size_t overlap); // pretty printing - friend std::ostream& operator<< (std::ostream& stream, const ChunkMetadata& meta); - + friend std::ostream& operator<<(std::ostream& stream, + const ChunkMetadata& meta); + // Accessor that grows the recorded shape of a chunk template void GrowChunk(std::size_t shape_growth); @@ -44,32 +44,35 @@ struct ChunkMetadata { // Accessor to swap two axes template void SwapAxes(); - + // data members - + // type of this chunk ChunkType chunk_type; - + // path to the file where this chunk lives std::filesystem::path filepath; // unique id: does not change when slices are appended id_t chunk_id; - + // offset within the file where this chunk lives std::streampos start_pos; // the shape of this chunk ... Vector shape; - - // ... as well as the indices of the elements marking the start and end element in the chunk + + // ... as well as the indices of the elements marking the start and end + // element in the chunk Vector start_ind; Vector end_ind; - // The number of samples (in each coordinate direction) this chunk overlaps with the neighbouring ones + // The number of samples (in each coordinate direction) this chunk overlaps + // with the neighbouring ones std::size_t overlap; - // subtract this from a global index to get the corresponding index within the chunk array data structure + // subtract this from a global index to get the corresponding index within the + // chunk array data structure Vector loc_ind_offset; }; @@ -80,27 +83,29 @@ template class ChunkIndex { public: - using metadata_t = ChunkMetadata; using shape_t = Vector; using ind_t = Vector; - -public: +public: // build from stored index file ChunkIndex(std::filesystem::path index_path); ~ChunkIndex(); - metadata_t& RegisterChunk(const Vector& start_ind, const Vector& shape, std::size_t overlap); - + metadata_t& RegisterChunk(const Vector& start_ind, + const Vector& shape, + std::size_t overlap); + // get chunk that contains the element with index `ind` const metadata_t& GetChunk(const Vector& ind) const; metadata_t& GetChunk(const Vector& ind); - // get chunks that, taken together, cover the rectangular region between `start_ind` and `end_ind` - std::vector> GetChunks(const Vector& start_ind, - const Vector& end_ind); - + // get chunks that, taken together, cover the rectangular region between + // `start_ind` and `end_ind` + std::vector> GetChunks( + const Vector& start_ind, + const Vector& end_ind); + shape_t GetShape(); // housekeeping and moving operations @@ -108,43 +113,47 @@ public: void MoveIndex(std::filesystem::path new_index_path); void ClearIndex(); void ImportIndex(std::filesystem::path index_path); - + // iterators over chunk metadata sets auto begin(); auto end(); public: - - // to check whether a specific index, or index range, is contained in or overlaps with a chunk - static bool is_in_chunk(const metadata_t& chunk, const Vector& ind); - static bool is_in_region(const Vector& start_ind, const Vector& shape, - const Vector& ind); - static bool chunk_overlaps_region(const metadata_t& chunk, const Vector& start_ind, const Vector& end_ind); - - // calculate the actual start / end index of the region where `chunk` overlaps with a specified index range - static Vector get_overlap_start_ind(const metadata_t& chunk, const Vector& start_ind); - static Vector get_overlap_end_ind(const metadata_t& chunk, const Vector& end_ind); - -public: + // to check whether a specific index, or index range, is contained in or + // overlaps with a chunk + static bool is_in_chunk(const metadata_t& chunk, + const Vector& ind); + static bool is_in_region(const Vector& start_ind, + const Vector& shape, + const Vector& ind); + static bool chunk_overlaps_region(const metadata_t& chunk, + const Vector& start_ind, + const Vector& end_ind); + + // calculate the actual start / end index of the region where `chunk` overlaps + // with a specified index range + static Vector get_overlap_start_ind( + const metadata_t& chunk, const Vector& start_ind); + static Vector get_overlap_end_ind( + const metadata_t& chunk, const Vector& end_ind); +public: Vector get_start_ind(); Vector get_end_ind(); - -private: +private: void invalidate_cached_index_metadata(); void calculate_and_cache_index_metadata(); Vector calculate_start_ind(); Vector calculate_end_ind(); - - metadata_t& find_chunk_by_index(const Vector& ind); + + metadata_t& find_chunk_by_index(const Vector& ind); id_t get_next_chunk_id(); - + void load_and_rebuild_index(); std::vector load_index_entries(std::filesystem::path index_path); - -private: +private: id_t m_next_chunk_id; std::filesystem::path m_index_path; @@ -153,7 +162,7 @@ private: shape_t m_shape; ind_t m_start_ind; ind_t m_end_ind; - + // TODO: use R-tree to quickly find the chunks in this list std::vector m_chunk_list; std::size_t m_last_accessed_ind; @@ -163,77 +172,90 @@ private: namespace CacheStatus { - struct Nothing { }; - struct Serialize { }; - + struct Nothing { + }; + struct Serialize { + }; + struct Append { - Append(std::size_t axis) : axis(axis) { }; + Append(std::size_t axis) + : axis(axis){}; std::size_t axis; - }; -} + }; +} // namespace CacheStatus // The elements stored in the cache -template class ArrayT, - typename T, std::size_t dims, std::size_t vec_dims> -struct CacheEntry { +template