Coverage for fiqus/post_processors/PostProcessCCT.py: 92%
290 statements
« prev ^ index » next coverage.py v6.4.4, created at 2024-05-20 03:24 +0200
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_dict = json.load(open(os.path.join(self.geom_folder, f'{w_name}.normals')))
324 winding_obj = Winding(self.cctdm, w_i, post_proc=True)
325 ww2 = self.cctdm.geometry.windings.wwws[w_i] / 2 # half of channel width
326 wh2 = self.cctdm.geometry.windings.wwhs[w_i] / 2 # half of channel height
327 wwns = self.cctdm.postproc.windings_wwns[w_i] # number of wires in width direction
328 whns = self.cctdm.postproc.windings_whns[w_i] # number of wires in height direction
329 wsw = self.cctdm.geometry.windings.wwws[w_i] / wwns # wire size in width
330 wsh = self.cctdm.geometry.windings.wwhs[w_i] / whns # wire size in height
331 for i_w in range(wwns):
332 for i_h in range(whns):
333 corner_i = 0
334 corners_lists = {'x': [], 'y': [], 'z': []}
335 for far in [False, True]:
336 for ii_w, ii_h in zip([0, 0, 1, 1], [0, 1, 1, 0]):
337 wt = (-ww2 + (ii_w + i_w) * wsw) / ww2
338 ht = (-wh2 + (ii_h + i_h) * wsh) / wh2
339 corner_i += 1
340 xs, ys, zs, turns_from_rotation = winding_obj.calc_turns_corner_coords(wt, ht, far=far)
341 zs = [z - winding_obj.z_corr for z in zs]
342 corners_lists['x'].append(xs)
343 corners_lists['y'].append(ys)
344 corners_lists['z'].append(zs)
345 nodeToPhysicalTurn = [(global_channel_pos - 1) * self.cctdm.geometry.windings.n_turnss[w_i] + t for t in turns_from_rotation]
346 nodeToHalfTurn = [winding_order.index(global_channel_pos) * self.cctdm.geometry.windings.n_turnss[w_i] + t for t in turns_from_rotation]
347 coord_center, ds = self._calc_coord_center_and_ds(corners_lists)
348 self._check_normals_match_corners(coord_center, normals_dict, w_name)
349 # if global_channel_pos == 1:
350 # s_ave_elem_len = np.copy(ds)
351 # s_ave_elem_len[0] = s_ave_elem_len[0] / 2
352 # s_ave = np.cumsum(s_ave_elem_len)
353 # else:
354 # s_ave = np.cumsum(ds) + s_max
355 # s_max = np.max(s_ave)
356 current = [abs(self.cctdm.solve.windings.currents[w_i]) / (wwns * whns)] * len(nodeToHalfTurn)
357 Bxs, Bys, Bzs, Bhs, Bws, Bls = self._get_fields(coord_center, normals_dict)
358 if gui:
359 self._plot_fields_in_views(w_name, global_channel_pos, coord_center, Bxs, Bys, Bzs)
360 # self._plot_turns_in_view(data_turns, coord_center, nodeToHalfTurn)
361 csv_data_dict[winding_order.index(global_channel_pos) + 1] = ds, corners_lists, nodeToHalfTurn, nodeToPhysicalTurn, current, Bxs, Bys, Bzs, Bls, Bhs, Bws
362 global_channel_pos += 1
363 for f_i, f_name in enumerate(self.cctdm.geometry.fqpls.names):
364 # ----- calculate centers of coordinates for hexes of fqpl
365 fqpl_obj = FQPL(self.cctdm, f_i)
366 corners_lists = {'x': [], 'y': [], 'z': []}
367 for corner_i in range(8):
368 for cord in ['x', 'y', 'z']:
369 corners_lists[cord].append([v_dict[corner_i][cord] for v_dict in fqpl_obj.hexes.values()])
370 coord_center, ds = self._calc_coord_center_and_ds(corners_lists)
371 normals_dict = json.load(open(os.path.join(self.geom_folder, f'{f_name}.normals')))
372 self._check_normals_match_corners(coord_center, normals_dict, f_name)
373 Bxs, Bys, Bzs, Bhs, Bws, Bls = self._get_fields(coord_center, normals_dict)
374 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
375 for z_close, z_far in zip(corners_lists['z'][0], corners_lists['z'][-1]):
376 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]:
377 self.masks_fqpls[f_i].append(True)
378 else:
379 self.masks_fqpls[f_i].append(False)
380 self.masks_fqpls[f_i] = np.array(self.masks_fqpls[f_i], dtype=bool)
381 for cord in ['x', 'y', 'z']:
382 for corner_i, corner in enumerate(corners_lists[cord]):
383 corners_lists[cord][corner_i] = list(np.array(corner)[self.masks_fqpls[f_i]])
384 lists_to_trim = [Bxs, Bys, Bzs, Bhs, Bws, Bls, ds]
385 for idx, array in enumerate(lists_to_trim):
386 lists_to_trim[idx] = self._trim_array(array, self.masks_fqpls[f_i])
387 [Bxs, Bys, Bzs, Bhs, Bws, Bls, ds] = lists_to_trim
388 for cord in ['x', 'y', 'z']:
389 coord_center[cord] = self._trim_array(coord_center[cord], self.masks_fqpls[f_i])
390 # s_ave = np.cumsum(ds) + s_max
391 # s_max = np.max(s_ave)
392 go_dict = {'from': 0, 'to': len(Bxs)//2} # fqpl 'go' part
393 ret_dict = {'from': len(Bxs)//2, 'to': len(Bxs)}
394 for go_ret in [go_dict, ret_dict]:
395 total_n_turns += 1
396 ds_half = ds[go_ret['from']:go_ret['to']]
397 corners_lists_half = {}
398 for cord in ['x', 'y', 'z']:
399 corners_lists_half[cord] = []
400 for corner_i, corner in enumerate(corners_lists[cord]):
401 corners_lists_half[cord].append(corner[go_ret['from']:go_ret['to']])
402 nodeToHalfTurn = [total_n_turns] * len(Bxs[go_ret['from']:go_ret['to']]) # total expected number of turns + one + fqpl number
403 nodeToPhysicalTurn = [total_n_turns] * len(Bxs[go_ret['from']:go_ret['to']])
404 current = [self.cctdm.solve.fqpls.currents[f_i]] * len(Bxs[go_ret['from']:go_ret['to']])
405 Bxs_half = Bxs[go_ret['from']:go_ret['to']]
406 Bys_half = Bys[go_ret['from']:go_ret['to']]
407 Bzs_half = Bzs[go_ret['from']:go_ret['to']]
408 Bls_half = Bls[go_ret['from']:go_ret['to']]
409 Bhs_half = Bhs[go_ret['from']:go_ret['to']]
410 Bws_half = Bws[go_ret['from']:go_ret['to']]
411 #s_ave_half = s_ave[go_ret['from']:go_ret['to']]
412 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
413 if gui:
414 self._plot_fields_in_views(f_name, global_channel_pos, coord_center, Bxs, Bys, Bzs)
415 global_channel_pos += 1
417 ds_global = []
418 csv_data_sorded = dict(sorted(csv_data_dict.items())).values()# sorted list of tuples with the csv data
419 for tuple_with_values in csv_data_sorded:
420 ds_global.append(tuple_with_values[0]) # 0 index in tuple is for ds, changing to list to easily extend.
422 s_ave_global = []
423 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
424 for ds_i, ds_np in enumerate(ds_global):
425 s_ave_np = np.cumsum(ds_np) + s_max
426 s_ave_global.append(s_ave_np)
427 s_max = np.max(s_ave_np)
429 for sorted_v_tuple, s_ave_np in zip(csv_data_sorded, s_ave_global):
430 self._load_to_fields_csv_ds(*sorted_v_tuple, s_ave_np) # this is to avoid rearranging s_ave by electrical order
432 #el_order_half_turns = [int(t) for t in list(pd.unique(self.data_for_csv['nodeToPhysicalTurn [-]']))]
433 #json.dump({'el_order_half_turns': el_order_half_turns}, open(os.path.join(self.model_folder, "el_order_half_turns.json"), 'w'))
434 self._save_fields_csv()
435 if gui:
436 self.gu.launch_interactive_GUI()
437 else:
438 gmsh.clear()
439 gmsh.finalize()