Coverage for steam_pysigma\utils\Utils.py: 91%
243 statements
« prev ^ index » next coverage.py v7.4.3, created at 2024-12-16 17:09 +0100
1import getpass
2import logging
3import os
4import subprocess
5from pathlib import Path
7import numpy as np
8import pandas as pd
9import ruamel.yaml
10from matplotlib import pyplot as plt
11from steam_pysigma.data.DataSettings import DataSettings
13logger = logging.getLogger(__name__)
16def read_data_from_yaml(full_file_path, data_class):
17 with open(full_file_path, 'r') as stream:
18 yaml = ruamel.yaml.YAML(typ='safe', pure=True)
19 yaml_str = yaml.load(stream)
20 return data_class(**yaml_str)
23def displayWaitAndClose(waitTimeBeforeMessage: float, waitTimeAfterMessage: float = 0, flag_show_text: bool = True):
24 """
26 **Function useful in Pycharm tests; it allows closing plots after some time that they are displayed **
28 Wait a certain time, display a message, and wait a certain time
30 :param waitTimeBeforeMessage: Time to wait before the message [s]
31 :type waitTimeBeforeMessage: float
32 :param waitTimeAfterMessage: Time to wait after the message [s]
33 :type waitTimeAfterMessage: float
34 :return:
35 """
36 # Show plots in Pycharm, wait a certain time, alert time is up, and close the window
37 plt.ion()
38 plt.tight_layout()
39 plt.show()
40 plt.draw()
41 plt.pause(waitTimeBeforeMessage)
42 if flag_show_text == True:
43 plt.title('Figure will close in {} seconds...'.format(waitTimeAfterMessage))
44 plt.pause(waitTimeAfterMessage)
47def create_coordinate_file(path_map2d, coordinate_file_path):
48 """
49 Creates a csv file with same coordinates as the map2d.
51 :param path_map2d: Map2d file to read coordinates from
52 :param coordinate_file_path: Path to csv filw to be created
53 :return:
54 """
55 df = pd.read_csv(path_map2d, delim_whitespace=True)
56 df_new = pd.DataFrame()
57 df_new["X-POS/MM"] = df["X-POS/MM"].apply(lambda x: x / 1000)
58 df_new["Y-POS/MM"] = df["Y-POS/MM"].apply(lambda x: x / 1000)
59 df_new.to_csv(coordinate_file_path, header=None, index=False)
62def arcCenter(C, iH, oH, iL, oL, diff_radius=None):
63 inner_radius = (np.sqrt(np.square(iH.x - C.x) + np.square(iH.y - C.y)) +
64 np.sqrt(np.square(iL.x - C.x) + np.square(iL.y - C.y))) / 2
65 if diff_radius:
66 outer_radius = inner_radius + diff_radius
67 else:
68 outer_radius = (np.sqrt(np.square(oH.x - C.x) + np.square(oH.y - C.y)) +
69 np.sqrt(np.square(oL.x - C.x) + np.square(oL.y - C.y))) / 2
70 d_inner = [0.5 * abs((iL.x - iH.x)), 0.5 * abs((iH.y - iL.y))]
71 d_outer = [0.5 * abs((oL.x - oH.x)), 0.5 * abs((oH.y - oL.y))]
72 aa = [np.sqrt(np.square(d_inner[0]) + np.square(d_inner[1])),
73 np.sqrt(np.square(d_outer[0]) + np.square(d_outer[1]))]
74 bb = [np.sqrt(np.square(inner_radius) - np.square(aa[0])), np.sqrt(np.square(outer_radius) - np.square(aa[1]))]
75 if iL.y < iH.y:
76 M_inner = [iH.x + d_inner[0], iL.y + d_inner[1]]
77 M_outer = [oH.x + d_outer[0], oL.y + d_outer[1]]
78 if iL.y >= 0.:
79 sign = [-1, -1]
80 else:
81 sign = [1, 1]
82 else:
83 M_inner = [iH.x + d_inner[0], iH.y + d_inner[1]]
84 M_outer = [oH.x + d_outer[0], oH.y + d_outer[1]]
85 if iL.y >= 0.:
86 sign = [1, -1]
87 else:
88 sign = [-1, 1]
89 inner = [M_inner[0] + sign[0] * bb[0] * d_inner[1] / aa[0], M_inner[1] + sign[1] * bb[0] * d_inner[0] / aa[0]]
90 outer = [M_outer[0] + sign[0] * bb[1] * d_outer[1] / aa[1], M_outer[1] + sign[1] * bb[1] * d_outer[0] / aa[1]]
91 return inner, outer
94def get_user_settings(settings_folder):
95 """
96 Helper function for getting user settings.
97 :return: dictionary with user settings
98 :rtype: dict
99 """
100 user_name = getpass.getuser()
101 if user_name == 'root':
102 user_name = 'SYSTEM'
103 elif user_name == 'MP-WIN-02$':
104 user_name = 'MP_WIN_02'
105 if not settings_folder:
106 settings_path = Path.joinpath(Path(__file__).parent.parent.parent, 'tests', f"settings.{user_name}.yaml")
107 else:
108 settings_path = Path(settings_folder).joinpath(f"settings.{user_name}.yaml")
109 if Path.exists(settings_path):
110 return read_data_from_yaml(settings_path, DataSettings)
111 else:
112 raise Exception(f'Settings file: {settings_path} does not exist!')
115def make_folder_if_not_existing(folder: str, verbose: bool = False):
116 if not os.path.isdir(folder):
117 Path(folder).mkdir(parents=True, exist_ok=True)
118 if verbose:
119 print("Folder {} does not exist. Making it now".format(folder))
122def run(magnet_name, working_dir_path, model_java_file_path, model_class_file_path, output_path, model_name,
123 split_java_file_path, COMSOL_compile_path, COMSOL_batch_path, compile_batch_file_path, java_jdk_path):
124 os.chdir(working_dir_path)
125 with open(model_java_file_path) as java_file:
126 print("Creating java files initialized.")
127 max_no_lines = 6e4
128 public_indexes = []
129 files = []
130 content = java_file.readlines()
131 for i, line in enumerate(content):
132 if "public static" in line:
133 public_indexes += [i]
135 no_lines = public_indexes[-1] - public_indexes[0]
136 no_files = int(np.ceil(no_lines / max_no_lines))
137 max_no_lines = round(no_lines / no_files)
138 real_indexes = [public_indexes[i] - public_indexes[0] for i in range(len(public_indexes))]
139 closest = [min(real_indexes, key=lambda x: abs(x - max_no_lines * (i + 1))) + public_indexes[0]
140 for i in range(no_files)]
141 no_run = [int(content[i][content[i].index('run') + 3:content[i].index('(')]) for i in closest[0:-1]]
142 no_run += [int(content[public_indexes[-2]][content[public_indexes[-2]].index('run')
143 + 3:content[public_indexes[-2]].index('(')]) + 1]
144 additional_lines = {}
145 for i in range(no_files):
146 file_path = os.path.join(output_path, model_name)
147 files += [open(file_path + '_%d.java' % i, 'w')]
148 name = model_name + '_%d' % i
149 split_java_file_path += [f"{os.path.join(output_path, model_name + '_%d' % i)}"]
150 files[i].writelines(content[0:public_indexes[0] - 2] + ['public class ' + name + ' {\n', '\n'])
151 if i == 0:
152 files[i].writelines('\tpublic static Model run1(Model mph) {\n')
153 files[i].writelines(content[public_indexes[0] + 2:closest[i]] + ['}\n'])
154 additional_lines[name] = {'start': 2, 'end': no_run[i] - 1}
155 elif i + 1 == no_files:
156 files[i].writelines(content[closest[i - 1]:public_indexes[-1]] + ['}\n'])
157 additional_lines[name] = {'start': no_run[i - 1], 'end': len(public_indexes) - 1}
158 else:
159 files[i].writelines(content[closest[i - 1]:closest[i]] + ['}\n'])
160 additional_lines[name] = {'start': no_run[i - 1], 'end': no_run[i] - 1}
161 files[i].close()
163 with open(model_java_file_path, 'w') as java_file:
164 content = content[0:public_indexes[0] + 2] + ['\n'] + \
165 content[public_indexes[-1]:public_indexes[-1] + 2] + ['\t}\n'] + ['}\n']
166 content.insert(public_indexes[0] + 2, '\t\tmph = ' + model_name + '_0.run1(mph);\n')
167 ll = 1
168 for name, item in additional_lines.items():
169 for j in range(item['end'] - item['start'] + 1):
170 content.insert(public_indexes[0] + 2 + ll + j,
171 '\t\tmph = ' + name + '.run' + str(item['start'] + j) + '(mph);\n')
172 ll += j + 1
173 content.insert(public_indexes[0] + 2 + ll, '\t\treturn mph;\n')
174 content.insert(public_indexes[0] + 3 + ll, '\t}\n')
175 java_file.writelines(content)
176 print("Creating java files DONE.")
177 print("Excecuting bat file")
178 script_lines = []
179 class_paths = ''
180 for file in split_java_file_path:
181 script_lines += [f'"{COMSOL_compile_path}" -jdkroot "{java_jdk_path}" -XX:+DisableExplicitGC "{file}.java"',
182 f'"{java_jdk_path}\\bin\\jar.exe" cf "{file}.jar" "{file}.class"']
183 class_paths += f'-classpathadd "{file}.jar" '
185 script_lines += [
186 f'"{COMSOL_compile_path}" -jdkroot "{java_jdk_path}" {class_paths}"{model_java_file_path}" -XX:+DisableExplicitGC',
187 f'"{COMSOL_batch_path}" -inputfile "{model_class_file_path}" '
188 f'-outputfile "{os.path.join(output_path, magnet_name)}.mph" ']
190 with open(compile_batch_file_path, "w") as outfile:
191 outfile.write("\n".join(str(line) for line in script_lines))
193 # Excecute process without DriverSIGMA
194 proc = subprocess.Popen([compile_batch_file_path], stdout=subprocess.PIPE, stderr=subprocess.PIPE,
195 universal_newlines=True)
196 (stdout, stderr) = proc.communicate()
198 if proc.returncode != 0:
199 print(stderr)
200 else:
201 print(stdout)
204def build_global_variables(g, model_data):
205 map = g.gateway.jvm.java.util.HashMap()
206 constants = g.constants
207 end_sim_time = model_data.Options_SIGMA.time_vector_solution.time_step[-1][-1]
208 # Cliq variables:
209 R_crow = model_data.Circuit.R_circuit
210 L_circuit = model_data.Circuit.L_circuit
211 C_cliq = model_data.Quench_Protection.CLIQ.C
212 L_cliq = model_data.Quench_Protection.CLIQ.L
214 V_cliq_0 = model_data.Quench_Protection.CLIQ.U0
215 I_cliq_0 = model_data.Quench_Protection.CLIQ.I0
217 sym_factor = model_data.Quench_Protection.CLIQ.sym_factor
218 cliq_time = model_data.Quench_Protection.CLIQ.t_trigger
219 if cliq_time > end_sim_time:
220 with_cliq = 0
221 else:
222 with_cliq = 1
224 if V_cliq_0 == None:
225 V_cliq_0 = 0
226 if I_cliq_0 == None:
227 I_cliq_0 = 0
228 if sym_factor == None:
229 sym_factor = 1
230 if L_circuit == None:
231 L_circuit = "1e-6"
232 if L_cliq == None:
233 L_cliq = "1e-6"
234 map.put(constants.LABEL_CLIQ_RCROW, str(R_crow))
235 map.put(constants.LABEL_CLIQ_LCIR, str(L_circuit))
236 map.put(constants.LABEL_CLIQ_CAPASITOR, str(C_cliq))
237 map.put(constants.LABEL_CLIQ_INDUCTANCE, str(L_cliq))
238 map.put(constants.LABEL_CLIQ_VOLTAGE_INITIAL, str(V_cliq_0))
239 map.put(constants.LABEL_CLIQ_CURRENT_INITIAL, str(I_cliq_0))
240 map.put(constants.LABEL_CLIQ_SYMFACTOR, str(sym_factor))
241 map.put(constants.LABEL_CLIQ_SWITCH, str(with_cliq))
243 FLAG_M_pers = model_data.Options_SIGMA.physics.FLAG_M_pers
244 FLAG_M_pers = "0" if FLAG_M_pers is None else FLAG_M_pers
246 FLAG_ifcc = model_data.Options_SIGMA.physics.FLAG_ifcc
247 FLAG_ifcc = "0" if FLAG_ifcc is None else FLAG_ifcc
249 FLAG_iscc_crossover = model_data.Options_SIGMA.physics.FLAG_iscc_crossover
250 FLAG_iscc_crossover = "0" if FLAG_iscc_crossover is None else FLAG_iscc_crossover
252 FLAG_iscc_adjw = model_data.Options_SIGMA.physics.FLAG_iscc_adjw
253 FLAG_iscc_adjw = "0" if FLAG_iscc_adjw is None else FLAG_iscc_adjw
255 FLAG_iscc_adjn = model_data.Options_SIGMA.physics.FLAG_iscc_adjn
256 FLAG_iscc_adjn = "0" if FLAG_iscc_adjn is None else FLAG_iscc_adjn
258 FLAG_quench_all = model_data.Options_SIGMA.quench_initialization.FLAG_quench_all
259 FLAG_quench_all = "0" if FLAG_quench_all is None else FLAG_quench_all
261 FLAG_quench_off = model_data.Options_SIGMA.quench_initialization.FLAG_quench_off
262 FLAG_quench_off = "0" if FLAG_quench_off is None else FLAG_quench_off
264 PARAM_time_quench = model_data.Options_SIGMA.quench_initialization.PARAM_time_quench
265 PARAM_time_quench = "0" if PARAM_time_quench is None else PARAM_time_quench
267 magnetic_length = model_data.GeneralParameters.magnetic_length
268 T_initial = model_data.GeneralParameters.T_initial
270 map.put(constants.LABEL_FLAG_IFCC, FLAG_ifcc)
271 map.put(constants.LABEL_FLAG_ISCC_CROSSOVER, FLAG_iscc_crossover)
272 map.put(constants.LABEL_FLAG_ISCC_ADJW, FLAG_iscc_adjw)
273 map.put(constants.LABEL_FLAG_ISCC_ADJN, FLAG_iscc_adjn)
274 map.put(constants.LABEL_FLAG_MPERS, FLAG_M_pers)
275 map.put(constants.LABEL_FLAG_QUENCH_ALL, FLAG_quench_all)
276 map.put(constants.LABEL_FLAG_QUENCH_OFF, FLAG_quench_off)
277 map.put(constants.LABEL_PARAM_QUENCH_TIME, PARAM_time_quench)
278 map.put(constants.LABEL_MAGNETIC_LENGTH, magnetic_length)
279 map.put(constants.LABEL_OPERATIONAL_TEMPERATUR, str(T_initial))
280 map.put(constants.LABEL_INIT_QUENCH_HEAT, str(50000))
281 map.put(constants.LABEL_QUENCH_TEMP, str(10))
283 num_qh = model_data.Quench_Protection.Quench_Heaters.N_strips
284 ins_list = model_data.Quench_Protection.Quench_Heaters.s_ins
285 w_list = model_data.Quench_Protection.Quench_Heaters.w
286 qh_to_bath_list = model_data.Quench_Protection.Quench_Heaters.s_ins_He
287 qh_steel_strip = model_data.Quench_Protection.Quench_Heaters.h
288 tau = [round(a * b, 4) for a, b in zip(model_data.Quench_Protection.Quench_Heaters.R_warm,
289 model_data.Quench_Protection.Quench_Heaters.C)]
290 num_qh_div = model_data.Options_SIGMA.quench_initialization.num_qh_div
291 u_init = model_data.Quench_Protection.Quench_Heaters.U0
292 frac_heater = [round(a / b, 4) for a, b in
293 zip(model_data.Quench_Protection.Quench_Heaters.l_stainless_steel, [sum(x) for x in zip(
294 model_data.Quench_Protection.Quench_Heaters.l_stainless_steel,
295 model_data.Quench_Protection.Quench_Heaters.l_copper)])]
296 trigger_time = model_data.Quench_Protection.Quench_Heaters.t_trigger
297 ins_thick_to_coil = model_data.Options_SIGMA.quench_initialization.th_coils
298 lengths_qh = model_data.Quench_Protection.Quench_Heaters.l
300 for i in range(num_qh):
301 if model_data.Options_SIGMA.time_vector_solution.time_step[-1][-1] < trigger_time[i]:
302 trigger_time[i] = end_sim_time + model_data.Options_SIGMA.time_vector_solution.time_step[-1][-2]
303 map.put(constants.LABEL_INSULATION_THICKNESS_QH_TO_COIL + str(i + 1), str(ins_list[i]))
304 map.put(constants.LABEL_WIDTH_QH + str(i + 1), str(w_list[i]))
305 map.put(constants.LABEL_INSULATION_THICKNESS_QH_TO_BATH + str(i + 1), str(qh_to_bath_list[i]))
306 map.put(constants.LABEL_INSULATION_THICKNESS_QH_STRIP + str(i + 1), str(qh_steel_strip[i]))
307 map.put(constants.LABEL_EXPONENTIAL_TIME_CONSTANT_DECAY + str(i + 1), str(tau[i]))
308 map.put(constants.LABEL_NUMBER_OF_QH_SUBDIVISIONS + str(i + 1), str(num_qh_div[i]))
309 map.put(constants.LABEL_QH + constants.LABEL_L + str(i + 1), str(lengths_qh[i]))
310 map.put(constants.LABEL_QH + str(i + 1) + constants.LABEL_FRACTION_OF_QH_STATION, str(frac_heater[i]))
311 map.put(constants.LABEL_TRIGGER_TIME_QH + str(i + 1), str(trigger_time[i]))
312 map.put(constants.LABEL_INSULATION_THICKNESS_TO_COIL + str(i + 1), str(ins_thick_to_coil[i]))
313 map.put(constants.LABEL_INITIAL_QH_VOLTAGE + str(i + 1), str(u_init[i]))
315 return map
318def build_study(study_type, srv, cfg, study_API, timeRange):
319 if study_type is not None:
320 if study_type == "Stationary":
321 # Add code to create and run study
322 study_API.setNewBasicStationaryStudy(srv, cfg, "sol1")
323 elif (study_type == "Transient"):
324 study_API.setNewMonolithicStudy(srv, cfg, "Default_study", timeRange)
327def create_results(srv, cfg, variables1DvsTime, time2DConverted, variables2DConverted, time1DConverted,
328 variables1DvsTimeVector,
329 result_api, input_coordinates_path, path_to_results):
330 for i in range(len(variables2DConverted)):
331 if len(time2DConverted) > 1:
332 time_vector_2D = ', '.join(str(x) for x in time2DConverted[i])
334 else:
335 time_vector_2D = time2DConverted[0]
336 result_api.create2DResultNode(srv, cfg, variables2DConverted[i], time_vector_2D, f"data {i}",
337 input_coordinates_path, path_to_results)
338 for j in range(len(time1DConverted)):
339 if len(time2DConverted) > 1:
340 time_vector_1D = ', '.join(str(x) for x in time1DConverted[j])
342 else:
343 time_vector_1D = time1DConverted[0]
345 if cfg.getStudyType() == "Transient":
346 for j in range(len(variables1DvsTime)):
347 result_api.create1DResultNodeAllTimes(srv, cfg, variables1DvsTime[j], f"1DExport_{j}", path_to_results)
349 for k in range(len(variables1DvsTimeVector)):
350 result_api.create1DResultNodeTimeVector(srv, cfg, variables1DvsTimeVector[k], time_vector_1D,
351 f"data {i + j + 1 + k}", path_to_results)