Line data Source code
1 : // Copyright Amazon.com, Inc. or its affiliates. All Rights Reserved.
2 : // SPDX-License-Identifier: Apache-2.0
3 :
4 : #ifndef PALACE_MODELS_POST_OPERATOR_CSV_HPP
5 : #define PALACE_MODELS_POST_OPERATOR_CSV_HPP
6 :
7 : #include <memory>
8 : #include <optional>
9 : #include "fem/errorindicator.hpp"
10 : #include "models/curlcurloperator.hpp"
11 : #include "models/laplaceoperator.hpp"
12 : #include "models/spaceoperator.hpp"
13 : #include "utils/configfile.hpp"
14 : #include "utils/filesystem.hpp"
15 : #include "utils/tablecsv.hpp"
16 : #include "utils/units.hpp"
17 :
18 : namespace palace
19 : {
20 : class IoData;
21 : class DomainPostOperator;
22 : class SurfacePostOperator;
23 : class InterpolationOperator;
24 : class SurfaceCurrentOperator;
25 : class LumpedPortOperator;
26 :
27 : // Advance declaration.
28 : template <ProblemType solver_t>
29 : class PostOperator;
30 :
31 : // Statically map solver (ProblemType) to finite element operator.
32 :
33 : template <ProblemType solver_t>
34 : struct fem_op_map_type
35 : {
36 : using type = SpaceOperator;
37 : };
38 : template <>
39 : struct fem_op_map_type<ProblemType::ELECTROSTATIC>
40 : {
41 : using type = LaplaceOperator;
42 : };
43 : template <>
44 : struct fem_op_map_type<ProblemType::MAGNETOSTATIC>
45 : {
46 : using type = CurlCurlOperator;
47 : };
48 :
49 : template <ProblemType solver_t>
50 : using fem_op_t = typename fem_op_map_type<solver_t>::type;
51 :
52 : // Results of measurements on fields. Not all measurements are sensible to define for all
53 : // solvers.
54 : struct Measurement
55 : {
56 : // Mini storage structs for data measurements.
57 : struct DomainData
58 : {
59 : int idx;
60 : double energy;
61 : double participation_ratio;
62 : };
63 :
64 : struct FluxData
65 : {
66 : int idx; // Surface index
67 : std::complex<double> Phi; // Integrated flux
68 : SurfaceFlux type;
69 : };
70 :
71 : struct InterfaceData
72 : {
73 : int idx; // Interface index
74 : double energy; // Surface Electric Field Energy
75 : double tandelta; // Dissipation tangent tan(δ)
76 : double energy_participation; // ratio of interface energy / total_energy
77 : double quality_factor; // 1 / (energy_participation * tan δ)
78 : };
79 :
80 : struct FarFieldData
81 : {
82 : // Theta: polar angle (0 to pi radians).
83 : // Phi: azimuthal angle (0 to 2pi radians).
84 : std::vector<std::pair<double, double>> thetaphis;
85 :
86 : // Components of the electric field.
87 : std::vector<std::array<std::complex<double>, 3>> E_field;
88 : };
89 :
90 : // Data for both lumped and wave port.
91 34 : struct PortPostData
92 : {
93 : std::complex<double> P = 0.0;
94 : std::complex<double> V = 0.0;
95 : std::complex<double> I = 0.0;
96 : // Separate R, L, and C branches for current via Z.
97 : std::array<std::complex<double>, 3> I_RLC = {0.0, 0.0, 0.0};
98 :
99 : // S-Parameter.
100 : std::complex<double> S = 0.0;
101 :
102 : // Energies (currently only for lumped port).
103 : double inductor_energy = 0.0; // E_ind = ∑_j 1/2 L_j I_mj².
104 : double capacitor_energy = 0.0; // E_cap = ∑_j 1/2 C_j V_mj².
105 :
106 : // Resistive lumped port (only eigenmode).
107 : double mode_port_kappa = 0.0;
108 : double quality_factor = mfem::infinity();
109 :
110 : // Inductive lumped port (only eigenmode).
111 : double inductive_energy_participation = 0.0;
112 : };
113 :
114 : // "Pseudo-measurements": input required during measurement or data which is stored here
115 : // in order to pass it along to the printers.
116 :
117 : int ex_idx = 0; // driven
118 :
119 : std::complex<double> freq = {0.0, 0.0}; // driven || eigenvalue.
120 :
121 : // Modulation factor for input excitation:
122 : // - I_inc(t) = J(t) I_in, for transient
123 : // - I_inc(omega) = I_in, for driven so that Jcoeff_excitation = 1.0
124 : double Jcoeff_excitation = 1.0; // transient || driven
125 :
126 : // Eigenmode data including error from solver.
127 : double eigenmode_Q = 0.0;
128 : double error_bkwd = 0.0;
129 : double error_abs = 0.0;
130 :
131 : // "Actual measurements".
132 :
133 : double domain_E_field_energy_all = 0.0;
134 : double domain_H_field_energy_all = 0.0;
135 :
136 : std::vector<DomainData> domain_E_field_energy_i;
137 : std::vector<DomainData> domain_H_field_energy_i;
138 :
139 : double lumped_port_capacitor_energy = 0.0;
140 : double lumped_port_inductor_energy = 0.0;
141 :
142 : std::map<int, PortPostData> lumped_port_vi;
143 : std::map<int, PortPostData> wave_port_vi;
144 :
145 : // Probe data is ordered as [Fx1, Fy1, Fz1, Fx2, Fy2, Fz2, ...].
146 : // TODO: Replace with proper matrix: mdspan (C++23) / Eigen.
147 : std::vector<std::complex<double>> probe_E_field;
148 : std::vector<std::complex<double>> probe_B_field;
149 :
150 : std::vector<FluxData> surface_flux_i;
151 : std::vector<InterfaceData> interface_eps_i;
152 : FarFieldData farfield;
153 :
154 : // Dimensionalize and nondimensionalize a set of measurements
155 : static Measurement Dimensionalize(const Units &units,
156 : const Measurement &nondim_measurement_cache);
157 : static Measurement Nondimensionalize(const Units &units,
158 : const Measurement &dim_measurement_cache);
159 : // Helpers for converting complex variable to magnitude in dB and phase.
160 5 : static double Magnitude(std::complex<double> x) { return 20.0 * std::log10(std::abs(x)); }
161 5 : static double Phase(std::complex<double> x) { return std::arg(x) * 180.0 / M_PI; }
162 : };
163 :
164 : namespace _impl
165 : {
166 : // Filling pattern of rows of column groups — needed to validate reload position of
167 : // previous data. Make it public in an _impl namespace for testing.
168 : std::vector<std::size_t> table_expected_filling(std::size_t m_idx_row, std::size_t ex_idx_i,
169 : std::size_t nr_rows,
170 : std::size_t nr_col_blocks);
171 :
172 : } // namespace _impl
173 :
174 : // Helper class to PostOperator to collect csv tables and printers for measurement that will
175 : // be saved to file. This class contains a pointer to the corresponding PostOperator class
176 : // and is a friend to a PostOperator class; this is equivalent to having these members
177 : // and methods in PostOperator. It exists for code clarity.
178 : template <ProblemType solver_t>
179 : class PostOperatorCSV
180 : {
181 : protected:
182 : // Copy savepath from PostOperator for simpler dependencies.
183 : fs::path post_dir;
184 : bool reload_table = false; // True only for driven simulation with non-default restart
185 :
186 : // Dimensionalized measurement cache. Converted from the PostOperator member variable.
187 : Measurement measurement_cache;
188 :
189 : // Cursor location & cursor value.
190 :
191 : std::size_t row_i = 0; // Plain count of current row (measurement index)
192 : std::size_t ex_idx_i = 0; // Plain count of current column group (excitation)
193 :
194 : double row_idx_v; // Value of row index (time, freq..); must be dimensionful
195 : std::size_t m_ex_idx = 0; // ex_idx_v: Excitation index value (= ex_idx_v_all[ex_idx_i])
196 :
197 : // Required in validation of re-loaded table (driven), otherwise just to reserve space.
198 : // Transient (adaptive time-stepping) or eigenvalue (converged eigenvalues) solver output
199 : // may differ from expectation.
200 : std::size_t nr_expected_measurement_rows = 1;
201 :
202 : // Stored column groups (excitations). Default single "0" for solvers without excitations.
203 : std::vector<std::size_t> ex_idx_v_all = {std::size_t(0)};
204 3 : bool HasSingleExIdx() const { return ex_idx_v_all.size() == 1; }
205 :
206 : void MoveTableValidateReload(TableWithCSVFile &t_csv_base, Table &&t_ref);
207 :
208 : // Data tables.
209 : //
210 : // These are all std::optional since: (a) should only be instantiated on the root mpi
211 : // process, (b) they should only be written if the data is non-empty.
212 :
213 : // Initialize and print methods for various output quantities: The initialize methods
214 : // prepare the tables for data insertion, whilst the print methods insert data
215 : // appropriately. Methods are only enabled when valid given the problem type.
216 :
217 : // Base (all solvers).
218 : std::optional<TableWithCSVFile> domain_E;
219 : void InitializeDomainE(const DomainPostOperator &dom_post_op);
220 : void PrintDomainE();
221 :
222 : std::optional<TableWithCSVFile> surface_F;
223 : void InitializeSurfaceF(const SurfacePostOperator &surf_post_op);
224 : void PrintSurfaceF();
225 :
226 : std::optional<TableWithCSVFile> surface_Q;
227 : void InitializeSurfaceQ(const SurfacePostOperator &surf_post_op);
228 : void PrintSurfaceQ();
229 :
230 : std::optional<TableWithCSVFile> probe_E;
231 : void InitializeProbeE(const InterpolationOperator &interp_op);
232 : void PrintProbeE(const InterpolationOperator &interp_op);
233 :
234 : std::optional<TableWithCSVFile> probe_B;
235 : void InitializeProbeB(const InterpolationOperator &interp_op);
236 : void PrintProbeB(const InterpolationOperator &interp_op);
237 :
238 : // TODO(C++20): Upgrade SFINAE to C++20 concepts to simplify static selection since we can
239 : // just use `void Function(...) requires (solver_t == Type::A);`.
240 :
241 : // Driven + Transient.
242 : std::optional<TableWithCSVFile> surface_I;
243 : template <ProblemType U = solver_t>
244 : auto InitializeSurfaceI(const SurfaceCurrentOperator &surf_j_op)
245 : -> std::enable_if_t<U == ProblemType::DRIVEN || U == ProblemType::TRANSIENT, void>;
246 : template <ProblemType U = solver_t>
247 : auto PrintSurfaceI(const SurfaceCurrentOperator &surf_j_op, const Units &units)
248 : -> std::enable_if_t<U == ProblemType::DRIVEN || U == ProblemType::TRANSIENT, void>;
249 :
250 : // Eigenmode + Driven + Transient.
251 : std::optional<TableWithCSVFile> port_V;
252 : std::optional<TableWithCSVFile> port_I;
253 : template <ProblemType U = solver_t>
254 : auto InitializePortVI(const SpaceOperator &fem_op)
255 : -> std::enable_if_t<U == ProblemType::EIGENMODE || U == ProblemType::DRIVEN ||
256 : U == ProblemType::TRANSIENT,
257 : void>;
258 : template <ProblemType U = solver_t>
259 : auto PrintPortVI(const LumpedPortOperator &lumped_port_op, const Units &units)
260 : -> std::enable_if_t<U == ProblemType::EIGENMODE || U == ProblemType::DRIVEN ||
261 : U == ProblemType::TRANSIENT,
262 : void>;
263 :
264 : // Driven.
265 : std::optional<TableWithCSVFile> port_S;
266 : template <ProblemType U = solver_t>
267 : auto InitializePortS(const SpaceOperator &fem_op)
268 : -> std::enable_if_t<U == ProblemType::DRIVEN, void>;
269 : template <ProblemType U = solver_t>
270 : auto PrintPortS() -> std::enable_if_t<U == ProblemType::DRIVEN, void>;
271 :
272 : // Driven + Eigenmode.
273 : std::optional<TableWithCSVFile> farfield_E;
274 : template <ProblemType U = solver_t>
275 : auto InitializeFarFieldE(const SurfacePostOperator &surf_post_op)
276 : -> std::enable_if_t<U == ProblemType::DRIVEN || U == ProblemType::EIGENMODE, void>;
277 : template <ProblemType U = solver_t>
278 : auto PrintFarFieldE(const SurfacePostOperator &surf_post_op)
279 : -> std::enable_if_t<U == ProblemType::DRIVEN || U == ProblemType::EIGENMODE, void>;
280 :
281 : // Eigenmode.
282 : std::optional<TableWithCSVFile> eig;
283 : template <ProblemType U = solver_t>
284 : auto InitializeEig() -> std::enable_if_t<U == ProblemType::EIGENMODE, void>;
285 : template <ProblemType U = solver_t>
286 : auto PrintEig() -> std::enable_if_t<U == ProblemType::EIGENMODE, void>;
287 :
288 : std::vector<int> ports_with_L;
289 : std::vector<int> ports_with_R;
290 : std::optional<TableWithCSVFile> port_EPR;
291 : template <ProblemType U = solver_t>
292 : auto InitializeEigPortEPR(const LumpedPortOperator &lumped_port_op)
293 : -> std::enable_if_t<U == ProblemType::EIGENMODE, void>;
294 : template <ProblemType U = solver_t>
295 : auto PrintEigPortEPR() -> std::enable_if_t<U == ProblemType::EIGENMODE, void>;
296 :
297 : std::optional<TableWithCSVFile> port_Q;
298 : template <ProblemType U = solver_t>
299 : auto InitializeEigPortQ(const LumpedPortOperator &lumped_port_op)
300 : -> std::enable_if_t<U == ProblemType::EIGENMODE, void>;
301 : template <ProblemType U = solver_t>
302 : auto PrintEigPortQ() -> std::enable_if_t<U == ProblemType::EIGENMODE, void>;
303 :
304 : public:
305 : // Print all data from nondim_measurement_cache.
306 : void PrintAllCSVData(const PostOperator<solver_t> &post_op,
307 : const Measurement &nondim_measurement_cache,
308 : double idx_value_dimensionful, int step);
309 :
310 : // Driven specific overload for specifying excitation index.
311 : template <ProblemType U = solver_t>
312 : auto PrintAllCSVData(const PostOperator<solver_t> &post_op,
313 : const Measurement &nondim_measurement_cache,
314 : double idx_value_dimensionful, int step, int ex_idx)
315 : -> std::enable_if_t<U == ProblemType::DRIVEN, void>
316 : {
317 3 : m_ex_idx = ex_idx;
318 3 : PrintAllCSVData(post_op, nondim_measurement_cache, idx_value_dimensionful, step);
319 3 : }
320 :
321 : // Special case of global indicator — init and print all at once.
322 : void PrintErrorIndicator(bool is_root,
323 : const ErrorIndicator::SummaryStatistics &indicator_stats);
324 :
325 : // "Delayed ctor" so that PostOperator can call it once it is fully constructed.
326 : // Set-up all files to be called from post_op.
327 : void InitializeCSVDataCollection(const PostOperator<solver_t> &post_op);
328 :
329 : explicit PostOperatorCSV(const IoData &iodata, const fem_op_t<solver_t> &fem_op);
330 : };
331 :
332 : } // namespace palace
333 :
334 : #endif // PALACE_MODELS_POST_OPERATOR_CSV_HPP
|
