Coverage for fiqus/post_processors/PostProcessCCT.py: 92%
294 statements
« prev ^ index » next coverage.py v7.4.4, created at 2025-11-03 01:32 +0000
1import os
2import math
3import csv
5import gmsh
6import json
7import numpy as np
8import pandas as pd
9from pathlib import Path
11from fiqus.geom_generators.GeometryCCT import Winding, FQPL
12from fiqus.data.DataWindingsCCT import WindingsInformation
13from fiqus.parsers.ParserDAT import ParserDAT
14from fiqus.utils.Utils import GmshUtils, FilesAndFolders
17class Post_Process:
18 def __init__(self, fdm, verbose=True):
19 """
20 Class to cct models postprocessing
21 :param fdm: FiQuS data model
22 :param verbose: If True more information is printed in python console.
23 """
24 self.cctdm = fdm.magnet
25 self.model_folder = os.path.join(os.getcwd())
26 self.magnet_name = fdm.general.magnet_name
28 self.verbose = verbose
29 self.gu = GmshUtils(self.model_folder, self.verbose)
30 self.gu.initialize()
31 self.masks_fqpls = {}
32 self.pos_names = []
33 for variable, volume, file_ext in zip(self.cctdm.postproc.variables, self.cctdm.postproc.volumes, self.cctdm.postproc.file_exts):
34 self.pos_names.append(f'{variable}_{volume}.{file_ext}')
35 self.pos_name = self.pos_names[0]
36 self.field_map_3D = os.path.join(self.model_folder, 'field_map_3D.csv')
37 self.model_file = self.field_map_3D
38 self.geom_folder = Path(self.model_folder).parent.parent
39 # csv output definition for fields
40 self.csv_ds = {'ds [m]': None} # length of each hexahedron
41 self.csv_sAve = {'sAvePositions [m]': None} # cumulative positions along the electrical order of the center of each hexahedron
42 self.csv_3Dsurf = { # basically this is the mapping of coordinates of corners to output. This dictionary defines which corner of hexahedra gets
43 # used and which coordinate for each name like 'x3Dsurf 1'. Look into Hexahedron.py base class for drawing of corner numbering
44 'x3Dsurf 1 [m]': {5: 'x'},
45 'x3Dsurf 2 [m]': {6: 'x'},
46 'y3Dsurf 1 [m]': {5: 'y'},
47 'y3Dsurf 2 [m]': {6: 'y'},
48 'z3Dsurf 1 [m]': {5: 'z'},
49 'z3Dsurf 2 [m]': {6: 'z'},
50 'x3Dsurf_2 1 [m]': {7: 'x'},
51 'x3Dsurf_2 2 [m]': {6: 'x'},
52 'y3Dsurf_2 1 [m]': {7: 'y'},
53 'y3Dsurf_2 2 [m]': {6: 'y'},
54 'z3Dsurf_2 1 [m]': {7: 'z'},
55 'z3Dsurf_2 2 [m]': {6: 'z'}}
56 self.csv_nodeToTurns = {'nodeToHalfTurn [-]': None, # hexahedron (node) mapping to "Electrical" turns (electrical order)
57 'nodeToPhysicalTurn [-]': None} # hexahedron (node) mapping to "Physical" turns (this is used for defining thermal connections between turns)
58 self.csv_TransportCurrent = {'TransportCurrent [A]': None} # Transport current in the turn (not channel) for which the magnetic field was calculated.
59 # The above does not need to be changed for running models at different current.
60 self.csv_Bs = {
61 # lists with magnetic field flux density components along cartesian coordinates (x, y, z) and l, h, w i.e. length, height and width of the channel (or wires if aligned with channel).
62 # x, y, z is used for display, l, h, w are used in LEDET for Ic scaling (well only h and w for perpendicular field).
63 'Bx [T]': 'Vx',
64 'By [T]': 'Vy',
65 'Bz [T]': 'Vz',
66 'Bl [T]': 'Vl',
67 'Bh [T]': 'Vh',
68 'Bw [T]': 'Vw',
69 }
70 self.headers_dict = {**self.csv_ds, **self.csv_sAve, **self.csv_3Dsurf, **self.csv_nodeToTurns, **self.csv_TransportCurrent, **self.csv_Bs}
71 self.data_for_csv = {}
72 self.precision_for_csv = {}
73 # -------- make keys in data for csv ----
74 for key in list(self.headers_dict.keys()):
75 self.precision_for_csv[key] = '%.6f' # data will be written to csv output file with 6 decimal places, i.e. down to um and uT level.
76 for key in list(self.csv_ds.keys()):
77 self.data_for_csv[key] = []
78 for key in list(self.csv_sAve.keys()):
79 self.data_for_csv[key] = []
80 for key in list(self.csv_3Dsurf.keys()):
81 self.data_for_csv[key] = []
82 for key in list(self.csv_nodeToTurns.keys()):
83 self.data_for_csv[key] = []
84 for key in list(self.csv_TransportCurrent.keys()):
85 self.data_for_csv[key] = []
86 for key in list(self.csv_Bs.keys()):
87 self.data_for_csv[key] = []
88 self.distance_key = list(self.csv_ds.keys())[0]
89 self.sAve_key = list(self.csv_sAve.keys())[0]
91 @staticmethod
92 def _get_fields(coord_center, normals):
93 """
94 Helper funciton to probe magnetic field solution form gmsh view and calculate magnetic field along the height, width and length of the channel section (hex element)
95 :param coord_center: dictionary with 'x', 'y' and 'z' coordinates stored as lists with items for each hexahedra. Magnetic field is probed at these locations.
96 :param normals: dictionary with normals along the height, width and length stored per coordinate direction 'n_x', 'n_y' and 'n_z' and lists with items for each hexahedra.
97 :return: tuple with list of Bxs, Bys, Bzs, Bhs, Bws, Bls
98 """
99 Bxs = []
100 Bys = []
101 Bzs = []
102 Bhs = []
103 Bws = []
104 Bls = []
105 view_tag = gmsh.view.getTags()[0]
106 for v, _ in enumerate(coord_center['x']):
107 field = gmsh.view.probe(view_tag, coord_center['x'][v], coord_center['y'][v], coord_center['z'][v])[0]
108 Bx = field[0]
109 By = field[1]
110 Bz = field[2]
111 Bmod = math.sqrt(Bx ** 2 + By ** 2 + Bz ** 2)
112 Bh = abs(Bx * normals['normals_h']['n_x'][v] + By * normals['normals_h']['n_y'][v] + Bz * normals['normals_h']['n_z'][v])
113 Bw = abs(Bx * normals['normals_w']['n_x'][v] + By * normals['normals_w']['n_y'][v] + Bz * normals['normals_w']['n_z'][v])
114 Bl = math.sqrt(abs(Bmod ** 2 - Bh ** 2 - Bw ** 2)) # sometimes it is very small and negative, so force it to be positive before taking the root.
115 for arr, val in zip([Bxs, Bys, Bzs, Bhs, Bws, Bls], [Bx, By, Bz, Bh, Bw, Bl]):
116 arr.append(val)
117 return Bxs, Bys, Bzs, Bhs, Bws, Bls
119 @staticmethod
120 def _plot_fields_in_views(p_name, global_channel_pos, coord_center, Bxs, Bys, Bzs):
121 """
122 Add gmsh list data view with name 'B cartesian {p_name}, {global_channel_pos}' at coordinates stored in coord_center and magnetic field values from list Bxs, Bys, Bzs
123 :param p_name: string with powered region name
124 :param global_channel_pos: integer with wire position in the channel
125 :param coord_center: dictionary with 'x', 'y' and 'z' coordinates stored as lists with items for each hexahedra. Magnetic field is plotted at these locations.
126 :param Bxs: list of magnetic field along the x axis to use in the view
127 :param Bys: list of magnetic field along the y axis to use in the view
128 :param Bzs: list of magnetic field along the z axis to use in the view
129 :return: none, adds view to currently initialized gmsh model and synchronizes
130 """
131 data_cartesian = []
132 for v, _ in enumerate(coord_center['x']):
133 data_cartesian.append(coord_center['x'][v])
134 data_cartesian.append(coord_center['y'][v])
135 data_cartesian.append(coord_center['z'][v])
136 data_cartesian.append(Bxs[v])
137 data_cartesian.append(Bys[v])
138 data_cartesian.append(Bzs[v])
139 gmsh.view.addListData(gmsh.view.add(f'B cartesian {p_name}, {global_channel_pos}'), "VP", len(data_cartesian) // 6, data_cartesian)
140 gmsh.model.occ.synchronize()
142 @staticmethod
143 def _plot_turns_in_view(data_turns, coord_center, turns_for_ch_pos):
144 """
145 Add gmsh list data view with name 'Turn numbering' at coordinates stored in coord_center and values from turns_for_ch_pos
146 :param data_turns: This is list to which the data gets appended so at the end one view with all the turns labeled is created with looping through turns in the channels.
147 :param coord_center: dictionary with 'x', 'y' and 'z' coordinates stored as lists with items for each hexahedra. Magnetic field is plotted at these locations.
148 :param turns_for_ch_pos: Lists with turn numbers to plot in the view.
149 :return: none, adds view to currently initialized gmsh model and synchronizes
150 """
151 for v, _ in enumerate(coord_center['x']):
152 data_turns.append(coord_center['x'][v])
153 data_turns.append(coord_center['y'][v])
154 data_turns.append(coord_center['z'][v])
155 data_turns.append(turns_for_ch_pos[v])
156 gmsh.view.addListData(gmsh.view.add('Turn numbering'), "SP", len(data_turns) // 4, data_turns)
157 gmsh.model.occ.synchronize()
159 @staticmethod
160 def _check_normals_match_corners(coord_center, normals_dict, p_name):
161 """
162 This is error checking function. Coord_centers are calculated in this class 'on the 'fly' form cctdm file. Normals are generated form geom_generators before meshing, then saved to file and now loaded here.
163 There is a chance that thses could get 'out of sync'. A basic check of the length (i.e. number of hexahedra involved is performed) of these lists is performed.
164 This will only detect if number of turns or discretization of each powered geom_generators has changed.
165 :param coord_center: list with coordinate centers. The content is not used only list length is queried.
166 :param normals_dict: list with normals read form json file. The content is not used only list length is queried.
167 :param p_name: string with powered volume name. Only used for error message display to say which region is not matching.
168 :return: None, prints error to the python console, if any.
169 """
170 for cord in ['x', 'y', 'z']:
171 if coord_center[cord].size != len(normals_dict['normals_h'][cord]) or coord_center[cord].size != len(normals_dict['normals_w'][cord]):
172 raise ValueError(f'Number of volumes in normals does not match the magnet definition for {p_name} winding!')
174 @staticmethod
175 def _calc_coord_center_and_ds(corners_lists):
176 """
177 Calculates a geometrical centre of each hexahedron base on its corners coordinates
178 :param corners_lists: lists of lists corresponding to a list of hexahedrons with each list containing 8 corner coordinates for each of the 8 points in the corners of hexahedron
179 :return: coord_center: dictionary with 'x', 'y' and 'z' coordinates stored as lists with items for each hexahedra. and ds
180 """
181 coord_center = {}
182 for cord in ['x', 'y', 'z']:
183 coord_center[cord] = np.mean(corners_lists[cord], axis=0) # coordinates center
184 surf_close = {}
185 surf_far = {}
186 for cord in ['x', 'y', 'z']:
187 surf_close[cord] = np.mean(corners_lists[cord][0:4], axis=0) # coordinates surface close
188 surf_far[cord] = np.mean(corners_lists[cord][4:8], axis=0) # coordinates surface far
189 ds = np.sqrt(np.square(surf_close['x'] - surf_far['x']) + np.square(surf_close['y'] - surf_far['y']) + np.square(surf_close['z'] - surf_far['z']))
190 return coord_center, ds
192 def _load_to_fields_csv_ds(self, ds, h3Dsurf, nodeToHalfTurn, nodeToPhysicalTurn, current, Bxs, Bys, Bzs, Bls, Bhs, Bws, s_ave):
193 """
194 Method to load data to csv dictionary to be dumped to a csv file.
195 :param ds: length of each hexahedron
196 :param h3Dsurf: dict with coordinates to export
197 :param nodeToHalfTurn: hexahedron (node) mapping to "Electrical" turns (electrical order)
198 :param nodeToPhysicalTurn: hexahedron (node) mapping to "Physical" turns (this is used for defining thermal connections between turns)
199 :param current: Transport current in the turn (not channel) for which the magnetic field was calculated. It does not need to be changed for running models at different current.
200 :param Bxs: list of magnetic flux density component along x direction
201 :param Bys: list of magnetic flux density component along y direction
202 :param Bzs: list of magnetic flux density component along z direction
203 :param Bls: list of magnetic flux density component along l (length) of wire direction
204 :param Bhs: list of magnetic flux density component along h (height) of wire direction
205 :param Bws: list of magnetic flux density component along w (wight) of wire direction
206 :param s_ave: cumulative positions along the electrical order of the center of each hexahedron
207 :return: nothing, appends data to this class attribute data_for_csv
208 """
209 for arr, key in zip([ds, s_ave], list(self.csv_ds.keys()) + list(self.csv_sAve.keys())):
210 self.data_for_csv[key].extend(arr)
211 for key, dict_what in self.csv_3Dsurf.items():
212 for p_n, p_c in dict_what.items():
213 self.data_for_csv[key].extend(h3Dsurf[p_c][p_n])
214 for arr, key in zip([Bxs, Bys, Bzs, Bls, Bhs, Bws, current, nodeToHalfTurn, nodeToPhysicalTurn], list(self.csv_Bs.keys()) + list(self.csv_TransportCurrent.keys()) + list(self.csv_nodeToTurns.keys())):
215 self.data_for_csv[key].extend(arr)
217 def _save_fields_csv(self):
218 """
219 Helper method to save csv file from data_for_csv class attribute.
220 :return: nothing, saves csv to file with name magnet_name(from yaml input file)_total_n_turns.csv
221 """
222 # csv_file_path = os.path.join(self.model_folder, f'{self.magnet_name[:self.magnet_name.index("_")]}_{int(total_n_turns)}.csv')
224 df = pd.DataFrame(self.data_for_csv)
225 for column in df:
226 df[column] = df[column].map(lambda x: self.precision_for_csv[column] % x)
227 # put NaN at the end of some columns to match what LEDET needs.
228 df.loc[df.index[-1], self.distance_key] = 'NaN'
229 df.loc[df.index[-1], self.sAve_key] = 'NaN'
230 # df.loc[df.index[-1], list(self.csv_nodeToTurns.keys())[0]] = 'NaN'
231 for h in list(self.csv_nodeToTurns.keys()):
232 df.loc[df.index[-1], h] = 'NaN'
233 df.loc[df.index[-1], list(self.csv_TransportCurrent.keys())[0]] = 'NaN'
234 for h in list(self.csv_Bs.keys()):
235 df.loc[df.index[-1], h] = 'NaN'
236 df.to_csv(self.field_map_3D, index=False, sep=',') # , float_format='%.7f')
238 @staticmethod
239 def _trim_array(array, mask):
240 return list(np.array(array)[mask])
242 @staticmethod
243 def _all_equal(iterator):
244 iterator = iter(iterator)
245 try:
246 first = next(iterator)
247 except StopIteration:
248 return True
249 return all(first == x for x in iterator)
251 def postporcess_inductance(self):
252 """
253 Function to postprocess .dat file with inductance saved by GetDP to selfMutualInductanceMatrix.csv required by LEDET for CCT magnet solve
254 :return: Nothing, writes selfMutualInductanceMatrix.csv file on disc.
255 :rtype: None
256 """
257 inductance_file_from_getdp = os.path.join(self.model_folder, 'Inductance.dat')
258 channels_inductance = ParserDAT(inductance_file_from_getdp).pqv # postprocessed quantity value (pqv) for Inductance.dat contains magnet self-inductance (channel, not wire turns)
260 total_n_turns_channel = 0 # including fqpls
261 for n_turns in self.cctdm.geometry.windings.n_turnss:
262 total_n_turns_channel = total_n_turns_channel + n_turns
263 if self._all_equal(self.cctdm.postproc.windings_wwns) and self._all_equal(self.cctdm.postproc.windings_whns):
264 number_of_turns_in_channel = self.cctdm.postproc.windings_wwns[0] * self.cctdm.postproc.windings_whns[0]
265 total_n_turns_channel = int(total_n_turns_channel)
266 total_n_turns_channel_and_fqpls = total_n_turns_channel
268 #windings_inductance = channels_inductance * total_n_turns_windings**2
269 windings_inductance = channels_inductance * number_of_turns_in_channel**3 # scaling with square for number of turns in the channel as the model has current in the channel but not number of turns
272 geometry_folder = Path(self.model_folder).parent.parent
273 winding_information_file = os.path.join(geometry_folder, f"{self.magnet_name}.wi")
274 cctwi = FilesAndFolders.read_data_from_yaml(winding_information_file, WindingsInformation)
275 windings_avg_length = cctwi.windings_avg_length
276 windings_inductance_per_m = windings_inductance / windings_avg_length
278 print(f"Channels self-inductance: {channels_inductance} H")
279 print(f"Number of turns in channel: {number_of_turns_in_channel} H")
280 print(f"Total number of channel turns: {total_n_turns_channel} turns")
281 print(f"Total number of inductance blocks: {total_n_turns_channel}")
282 print(f"Windings self-inductance: {windings_inductance} H")
283 print(f"Windings self-inductance per meter: {windings_inductance_per_m} H/m")
284 print(f"Magnetic length: {windings_avg_length} m")
286 fqpls_dummy_inductance_block = 1e-9 # H/m (H per meter of magnet!)
287 for _ in self.cctdm.geometry.fqpls.names:
288 total_n_turns_channel_and_fqpls = total_n_turns_channel_and_fqpls + 1
289 windings_inductance_per_m = windings_inductance_per_m - fqpls_dummy_inductance_block # subtracting fqpls as they will be added later to the matrix
291 win_ind_for_matrix = windings_inductance_per_m / total_n_turns_channel**2
292 M_matrix = np.ones(shape=(total_n_turns_channel_and_fqpls, total_n_turns_channel_and_fqpls)) * win_ind_for_matrix
293 M_matrix[:, total_n_turns_channel:total_n_turns_channel_and_fqpls] = 0.0 # mutual inductance of windings to fqpl set to zero for the columns
294 M_matrix[total_n_turns_channel:total_n_turns_channel_and_fqpls, :] = 0.0 # mutual inductance of windings to fqpl set to zero for the rows
295 for ij in range(total_n_turns_channel, total_n_turns_channel_and_fqpls, 1):
296 M_matrix[ij, ij] = fqpls_dummy_inductance_block # self-inductance of fqpls inductance block
297 # below code is the same as in BuilderLEDET in steam_sdk
298 csv_write_path = os.path.join(self.model_folder, 'selfMutualInductanceMatrix.csv')
299 print(f'Saving square matrix of size {total_n_turns_channel_and_fqpls}x{total_n_turns_channel_and_fqpls} into: {csv_write_path}')
300 with open(csv_write_path, 'w', newline='') as file:
301 reader = csv.writer(file)
302 reader.writerow(["# Extended self mutual inductance matrix [H/m]"])
303 for i in range(M_matrix.shape[0]):
304 reader.writerow(M_matrix[i])
306 def postprocess_fields(self, gui=False):
307 """
308 Methods to calculated 'virtual' turn positions in the channels and output information for about them for LEDET as a csv file.
309 :param gui: if True, graphical user interface (gui) of gmsh is shown with views created
310 :return: nothing, writes csv to file
311 """
312 winding_order = self.cctdm.postproc.winding_order
313 total_n_turns = 0
314 for n_turns, wwn, whn in zip(self.cctdm.geometry.windings.n_turnss, self.cctdm.postproc.windings_wwns, self.cctdm.postproc.windings_whns):
315 total_n_turns = total_n_turns + n_turns * wwn * whn
317 gmsh.open(self.pos_name)
318 global_channel_pos = 1
319 # data_turns = [] # used in self._plot_turns_in_view() - do not delete
320 csv_data_dict = {}
321 s_max = 0
322 for w_i, w_name in enumerate(self.cctdm.geometry.windings.names): # + self.cctwi.f_names:
323 normals_path = os.path.join(self.geom_folder, f"{w_name}.normals")
324 with open(normals_path, "r", encoding="cp1252") as f:
325 normals_dict = json.load(f)
326 winding_obj = Winding(self.cctdm, w_i, post_proc=True)
327 ww2 = self.cctdm.geometry.windings.wwws[w_i] / 2 # half of channel width
328 wh2 = self.cctdm.geometry.windings.wwhs[w_i] / 2 # half of channel height
329 wwns = self.cctdm.postproc.windings_wwns[w_i] # number of wires in width direction
330 whns = self.cctdm.postproc.windings_whns[w_i] # number of wires in height direction
331 wsw = self.cctdm.geometry.windings.wwws[w_i] / wwns # wire size in width
332 wsh = self.cctdm.geometry.windings.wwhs[w_i] / whns # wire size in height
333 for i_w in range(wwns):
334 for i_h in range(whns):
335 corner_i = 0
336 corners_lists = {'x': [], 'y': [], 'z': []}
337 for far in [False, True]:
338 for ii_w, ii_h in zip([0, 0, 1, 1], [0, 1, 1, 0]):
339 wt = (-ww2 + (ii_w + i_w) * wsw) / ww2
340 ht = (-wh2 + (ii_h + i_h) * wsh) / wh2
341 corner_i += 1
342 xs, ys, zs, turns_from_rotation = winding_obj.calc_turns_corner_coords(wt, ht, far=far)
343 zs = [z - winding_obj.z_corr for z in zs]
344 corners_lists['x'].append(xs)
345 corners_lists['y'].append(ys)
346 corners_lists['z'].append(zs)
347 nodeToPhysicalTurn = [(global_channel_pos - 1) * self.cctdm.geometry.windings.n_turnss[w_i] + t for t in turns_from_rotation]
348 nodeToHalfTurn = [winding_order.index(global_channel_pos) * self.cctdm.geometry.windings.n_turnss[w_i] + t for t in turns_from_rotation]
349 coord_center, ds = self._calc_coord_center_and_ds(corners_lists)
350 self._check_normals_match_corners(coord_center, normals_dict, w_name)
351 # if global_channel_pos == 1:
352 # s_ave_elem_len = np.copy(ds)
353 # s_ave_elem_len[0] = s_ave_elem_len[0] / 2
354 # s_ave = np.cumsum(s_ave_elem_len)
355 # else:
356 # s_ave = np.cumsum(ds) + s_max
357 # s_max = np.max(s_ave)
358 current = [abs(self.cctdm.solve.windings.currents[w_i]) / (wwns * whns)] * len(nodeToHalfTurn)
359 Bxs, Bys, Bzs, Bhs, Bws, Bls = self._get_fields(coord_center, normals_dict)
360 if gui:
361 self._plot_fields_in_views(w_name, global_channel_pos, coord_center, Bxs, Bys, Bzs)
362 # self._plot_turns_in_view(data_turns, coord_center, nodeToHalfTurn)
363 csv_data_dict[winding_order.index(global_channel_pos) + 1] = ds, corners_lists, nodeToHalfTurn, nodeToPhysicalTurn, current, Bxs, Bys, Bzs, Bls, Bhs, Bws
364 global_channel_pos += 1
365 for f_i, f_name in enumerate(self.cctdm.geometry.fqpls.names):
366 # ----- calculate centers of coordinates for hexes of fqpl
367 fqpl_obj = FQPL(self.cctdm, f_i)
368 corners_lists = {'x': [], 'y': [], 'z': []}
369 for corner_i in range(8):
370 for cord in ['x', 'y', 'z']:
371 corners_lists[cord].append([v_dict[corner_i][cord] for v_dict in fqpl_obj.hexes.values()])
372 coord_center, ds = self._calc_coord_center_and_ds(corners_lists)
373 normals_path = os.path.join(self.geom_folder, f"{f_name}.normals")
374 with open(normals_path, "r", encoding="cp1252") as f:
375 normals_dict = json.load(f)
376 self._check_normals_match_corners(coord_center, normals_dict, f_name)
377 Bxs, Bys, Bzs, Bhs, Bws, Bls = self._get_fields(coord_center, normals_dict)
378 self.masks_fqpls[f_i] = [] # dictionary to keep masks for trimming one end of fqpls that is long to extend to the end of air region, but not needed that long in simulations
379 for z_close, z_far in zip(corners_lists['z'][0], corners_lists['z'][-1]):
380 if z_close > -self.cctdm.geometry.fqpls.z_ends[f_i] * self.cctdm.postproc.fqpl_export_trim_tol[f_i] and z_far > -self.cctdm.geometry.fqpls.z_ends[f_i] * self.cctdm.postproc.fqpl_export_trim_tol[f_i]:
381 self.masks_fqpls[f_i].append(True)
382 else:
383 self.masks_fqpls[f_i].append(False)
384 self.masks_fqpls[f_i] = np.array(self.masks_fqpls[f_i], dtype=bool)
385 for cord in ['x', 'y', 'z']:
386 for corner_i, corner in enumerate(corners_lists[cord]):
387 corners_lists[cord][corner_i] = list(np.array(corner)[self.masks_fqpls[f_i]])
388 lists_to_trim = [Bxs, Bys, Bzs, Bhs, Bws, Bls, ds]
389 for idx, array in enumerate(lists_to_trim):
390 lists_to_trim[idx] = self._trim_array(array, self.masks_fqpls[f_i])
391 [Bxs, Bys, Bzs, Bhs, Bws, Bls, ds] = lists_to_trim
392 for cord in ['x', 'y', 'z']:
393 coord_center[cord] = self._trim_array(coord_center[cord], self.masks_fqpls[f_i])
394 # s_ave = np.cumsum(ds) + s_max
395 # s_max = np.max(s_ave)
396 go_dict = {'from': 0, 'to': len(Bxs)//2} # fqpl 'go' part
397 ret_dict = {'from': len(Bxs)//2, 'to': len(Bxs)}
398 for go_ret in [go_dict, ret_dict]:
399 total_n_turns += 1
400 ds_half = ds[go_ret['from']:go_ret['to']]
401 corners_lists_half = {}
402 for cord in ['x', 'y', 'z']:
403 corners_lists_half[cord] = []
404 for corner_i, corner in enumerate(corners_lists[cord]):
405 corners_lists_half[cord].append(corner[go_ret['from']:go_ret['to']])
406 nodeToHalfTurn = [total_n_turns] * len(Bxs[go_ret['from']:go_ret['to']]) # total expected number of turns + one + fqpl number
407 nodeToPhysicalTurn = [total_n_turns] * len(Bxs[go_ret['from']:go_ret['to']])
408 current = [self.cctdm.solve.fqpls.currents[f_i]] * len(Bxs[go_ret['from']:go_ret['to']])
409 Bxs_half = Bxs[go_ret['from']:go_ret['to']]
410 Bys_half = Bys[go_ret['from']:go_ret['to']]
411 Bzs_half = Bzs[go_ret['from']:go_ret['to']]
412 Bls_half = Bls[go_ret['from']:go_ret['to']]
413 Bhs_half = Bhs[go_ret['from']:go_ret['to']]
414 Bws_half = Bws[go_ret['from']:go_ret['to']]
415 #s_ave_half = s_ave[go_ret['from']:go_ret['to']]
416 csv_data_dict[total_n_turns] = ds_half, corners_lists_half, nodeToHalfTurn, nodeToPhysicalTurn, current, Bxs_half, Bys_half, Bzs_half, Bls_half, Bhs_half, Bws_half
417 if gui:
418 self._plot_fields_in_views(f_name, global_channel_pos, coord_center, Bxs, Bys, Bzs)
419 global_channel_pos += 1
421 ds_global = []
422 csv_data_sorded = dict(sorted(csv_data_dict.items())).values()# sorted list of tuples with the csv data
423 for tuple_with_values in csv_data_sorded:
424 ds_global.append(tuple_with_values[0]) # 0 index in tuple is for ds, changing to list to easily extend.
426 s_ave_global = []
427 s_max = -ds_global[0][0]/2 # make a s_max to a negative half of the length of the first hexahedra. This is to set the first s_ave max to a half of the first ds element
428 for ds_i, ds_np in enumerate(ds_global):
429 s_ave_np = np.cumsum(ds_np) + s_max
430 s_ave_global.append(s_ave_np)
431 s_max = np.max(s_ave_np)
433 for sorted_v_tuple, s_ave_np in zip(csv_data_sorded, s_ave_global):
434 self._load_to_fields_csv_ds(*sorted_v_tuple, s_ave_np) # this is to avoid rearranging s_ave by electrical order
436 #el_order_half_turns = [int(t) for t in list(pd.unique(self.data_for_csv['nodeToPhysicalTurn [-]']))]
437 #json.dump({'el_order_half_turns': el_order_half_turns}, open(os.path.join(self.model_folder, "el_order_half_turns.json"), 'w'))
438 self._save_fields_csv()
439 if gui:
440 self.gu.launch_interactive_GUI()
441 else:
442 gmsh.clear()
443 gmsh.finalize()