Coverage for fiqus/data/DataFiQuSPancake3D.py: 81%
1148 statements
« prev ^ index » next coverage.py v7.4.4, created at 2025-05-04 03:30 +0200
1from typing import Literal, Optional, Annotated, Union
2from contextvars import ContextVar
3import logging
4import math
5import pathlib
6import scipy.integrate
7from functools import cached_property
8from pydantic import (
9 BaseModel,
10 PositiveFloat,
11 NonNegativeFloat,
12 PositiveInt,
13 Field,
14 field_validator,
15 model_validator,
16 computed_field,
17 ValidationInfo,
18)
19import pandas as pd
20import numpy as np
21import matplotlib
22from annotated_types import Len
24logger = logging.getLogger(__name__)
26# ======================================================================================
27# Available materials: =================================================================
28NormalMaterialName = Literal[
29 "Copper", "Hastelloy", "Silver", "Indium", "Stainless Steel", "Kapton", "G10"
30]
31SuperconductingMaterialName = Literal["HTSSuperPower", "HTSFujikura", "HTSSucci"]
32# ======================================================================================
33# ======================================================================================
35# ======================================================================================
36# Material information: ================================================================
37resistivityMacroNames = {
38 "Copper": "MATERIAL_Resistivity_Copper_T_B",
39 "Hastelloy": "MATERIAL_Resistivity_Hastelloy_T",
40 "Silver": "MATERIAL_Resistivity_Silver_T_B",
41 "Indium": "MATERIAL_Resistivity_Indium_T",
42 "Stainless Steel": "MATERIAL_Resistivity_SSteel_T",
43}
44thermalConductivityMacroNames = {
45 "Copper": "MATERIAL_ThermalConductivity_Copper_T_B",
46 "Hastelloy": "MATERIAL_ThermalConductivity_Hastelloy_T",
47 "Silver": "MATERIAL_ThermalConductivity_Silver_T",
48 "Indium": "MATERIAL_ThermalConductivity_Indium_T",
49 "Stainless Steel": "MATERIAL_ThermalConductivity_SSteel_T",
50 "Kapton": "MATERIAL_ThermalConductivity_Kapton_T",
51 "G10": "MATERIAL_ThermalConductivity_G10_T",
52}
53heatCapacityMacroNames = {
54 "Copper": "MATERIAL_SpecificHeatCapacity_Copper_T",
55 "Hastelloy": "MATERIAL_SpecificHeatCapacity_Hastelloy_T",
56 "Silver": "MATERIAL_SpecificHeatCapacity_Silver_T",
57 "Indium": "MATERIAL_SpecificHeatCapacity_Indium_T",
58 "Stainless Steel": "MATERIAL_SpecificHeatCapacity_SSteel_T",
59 "Kapton": "MATERIAL_SpecificHeatCapacity_Kapton_T",
60 "G10": "MATERIAL_SpecificHeatCapacity_G10_T",
61}
62getdpTSAStiffnessThermalConductivityMacroNames = {
63 "Indium": "MATERIAL_ThermalConductivity_Indium_TSAStiffness_T",
64 "Stainless Steel": "MATERIAL_ThermalConductivity_SSteel_TSAStiffness_T",
65 "Kapton": "MATERIAL_ThermalConductivity_Kapton_TSAStiffness_T",
66 "G10": "MATERIAL_ThermalConductivity_G10_TSAStiffness_T",
67 "Copper": "MATERIAL_ThermalConductivity_Copper_TSAStiffness_T",
68}
69getdpTSAMassThermalConductivityMacroNames = {
70 "Indium": "MATERIAL_ThermalConductivity_Indium_TSAMass_T",
71 "Stainless Steel": "MATERIAL_ThermalConductivity_SSteel_TSAMass_T",
72 "Kapton": "MATERIAL_ThermalConductivity_Kapton_TSAMass_T",
73 "G10": "MATERIAL_ThermalConductivity_G10_TSAMass_T",
74 "Copper": "MATERIAL_ThermalConductivity_Copper_TSAMass_T",
75}
76getdpTSAMassHeatCapacityMacroNames = {
77 "Indium": "MATERIAL_SpecificHeatCapacity_Indium_TSAMass_T",
78 "Stainless Steel": "MATERIAL_SpecificHeatCapacity_SSteel_TSAMass_T",
79 "Kapton": "MATERIAL_SpecificHeatCapacity_Kapton_TSAMass_T",
80 "G10": "MATERIAL_SpecificHeatCapacity_G10_TSAMass_T",
81 "Copper": "MATERIAL_SpecificHeatCapacity_Copper_TSAMass_T",
82}
83getdpTSARHSFunctions = {
84 "Indium": "TSA_CFUN_rhoIn_T_constantThickness_rhs",
85 "Stainless Steel": None,
86}
87getdpTSATripleFunctions = {
88 "Indium": "TSA_CFUN_rhoIn_T_constantThickness_triple",
89 "Stainless Steel": None,
90}
91getdpTSAOnlyResistivityFunctions = {
92 "Indium": "TSA_CFUN_rhoIn_T_constantThickness_fct_only",
93 "Stainless Steel": None,
94}
95getdpTSAMassResistivityMacroNames = {
96 "Indium": "MATERIAL_Resistivity_Indium_TSAMass_T",
97 "Stainless Steel": None,
98 "Copper": "MATERIAL_Resistivity_Copper_TSAMass_T",
99}
100getdpTSAStiffnessResistivityMacroNames = {
101 "Indium": "MATERIAL_Resistivity_Indium_TSAStiffness_T",
102 "Stainless Steel": None,
103 "Copper": "MATERIAL_Resistivity_Copper_TSAStiffness_T",
104}
105getdpCriticalCurrentDensityFunctions = {
106 "HTSSuperPower": "CFUN_HTS_JcFit_SUPERPOWER_T_B_theta",
107 "HTSFujikura": "CFUN_HTS_JcFit_Fujikura_T_B_theta",
108 "HTSSucci": "CFUN_HTS_JcFit_Succi_T_B",
109}
110getdpNormalMaterialNames = {
111 "Copper": "Copper",
112 "Hastelloy": "Hastelloy",
113 "Silver": "Silver",
114 "Indium": "Indium",
115 "Stainless Steel": "StainlessSteel",
116 "Kapton": "Kapton",
117 "G10": "G10",
118}
119# ======================================================================================
120# ======================================================================================
122# ======================================================================================
123# Available quantities: ================================================================
124PositionRequiredQuantityName = Literal[
125 "magneticField",
126 "magnitudeOfMagneticField",
127 "currentDensity",
128 "magnitudeOfCurrentDensity",
129 "resistiveHeating",
130 "temperature",
131 "criticalCurrentDensity",
132 "heatFlux",
133 "resistivity",
134 "thermalConductivity",
135 "specificHeatCapacity",
136 "jHTSOverjCritical",
137 "criticalCurrent",
138 "axialComponentOfTheMagneticField",
139 "debug",
140 "jHTS",
141 "currentSharingIndex",
142 "arcLength",
143 "turnNumber"
144]
145PositionNotRequiredQuantityName = Literal[
146 "currentThroughCoil",
147 "voltageBetweenTerminals",
148 "inductance",
149 "timeConstant",
150 "totalResistiveHeating",
151 "magneticEnergy",
152 "maximumTemperature",
153 "cryocoolerAveragePower",
154 "cryocoolerAverageTemperature"
155]
156# ======================================================================================
157# ======================================================================================
159# ======================================================================================
160# Quantity information: ================================================================
161EMQuantities = [
162 "magneticField",
163 "magnitudeOfMagneticField",
164 "currentDensity",
165 "magnitudeOfCurrentDensity",
166 "resistiveHeating",
167 "criticalCurrentDensity",
168 "resistivity",
169 "jHTSOverjCritical",
170 "criticalCurrent",
171 "debug",
172 "inductance",
173 "timeConstant",
174 "currentThroughCoil",
175 "voltageBetweenTerminals",
176 "totalResistiveHeating",
177 "magneticEnergy",
178 "axialComponentOfTheMagneticField",
179 "jHTS",
180 "currentSharingIndex",
181 "arcLength",
182 "turnNumber"
183]
184ThermalQuantities = [
185 "temperature",
186 "heatFlux",
187 "thermalConductivity",
188 "specificHeatCapacity",
189 "maximumTemperature",
190 "debug",
191 "cryocoolerAveragePower",
192 "cryocoolerAverageTemperature"
193]
194quantityProperNames = {
195 "magneticField": "Magnetic Field",
196 "magneticEnergy": "Magnetic Energy",
197 "magnitudeOfMagenticField": "Magnitude of Magnetic Field",
198 "currentDensity": "Current Density",
199 "magnitudeOfCurrentDensity": "Magnitude of Current Density",
200 "resistiveHeating": "Resistive Heating",
201 "totalResistiveHeating": "Total Resistive Heating",
202 "temperature": "Temperature",
203 "currentThroughCoil": "Current Through Coil",
204 "voltageBetweenTerminals": "Voltage Between Terminals",
205 "criticalCurrentDensity": "Critical Current Density",
206 "heatFlux": "Heat Flux",
207 "resistivity": "Resistivity",
208 "thermalConductivity": "Thermal Conductivity",
209 "specificHeatCapacity": "Specific Heat Capacity",
210 "jHTSOverjCritical": "jHTS/jCritical",
211 "criticalCurrent": "Critical Current",
212 "debug": "Debug",
213 "inductance": "Inductance",
214 "timeConstant": "Time Constant",
215 "axialComponentOfTheMagneticField": "Axial Component of the Magnetic Field",
216 "maximumTemperature": "Maximum Temperature",
217 "jHTS": "Current Density in HTS Layer",
218 "currentSharingIndex": "Current Sharing Index",
219 "cryocoolerAveragePower": "Cryocooler Average Power",
220 "arcLength": "Arc Length",
221 "turnNumber": "Turn Number",
222 "cryocoolerAverageTemperature": "Cryocooler Average Temperature"
223}
225quantityUnits = {
226 "magneticField": "T",
227 "magneticEnergy": "J",
228 "magnitudeOfMagneticField": "T",
229 "currentDensity": "A/m^2",
230 "magnitudeOfCurrentDensity": "A/m^2",
231 "resistiveHeating": "W",
232 "totalResistiveHeating": "W",
233 "temperature": "K",
234 "currentThroughCoil": "A",
235 "voltageBetweenTerminals": "V",
236 "criticalCurrentDensity": "A/m^2",
237 "heatFlux": "W/m^2",
238 "resistivity": "Ohm*m",
239 "thermalConductivity": "W/m*K",
240 "specificHeatCapacity": "J/kg*K",
241 "jHTSOverjCritical": "-",
242 "criticalCurrent": "A",
243 "debug": "1",
244 "inductance": "H",
245 "timeConstant": "s",
246 "axialComponentOfTheMagneticField": "T",
247 "maximumTemperature": "K",
248 "jHTS": "A/m^2",
249 "currentSharingIndex": "-",
250 "cryocoolerAveragePower": "W",
251 "arcLength": "m",
252 "turnNumber": "-",
253 "cryocoolerAverageTemperature": "K"
254}
256getdpQuantityNames = {
257 "magneticField": "RESULT_magneticField",
258 "magneticEnergy": "RESULT_magneticEnergy",
259 "magnitudeOfMagneticField": "RESULT_magnitudeOfMagneticField",
260 "currentDensity": "RESULT_currentDensity",
261 "magnitudeOfCurrentDensity": "RESULT_magnitudeOfCurrentDensity",
262 "resistiveHeating": "RESULT_resistiveHeating",
263 "totalResistiveHeating": "RESULT_totalResistiveHeating",
264 "temperature": "RESULT_temperature",
265 "currentThroughCoil": "RESULT_currentThroughCoil",
266 "voltageBetweenTerminals": "RESULT_voltageBetweenTerminals",
267 "criticalCurrentDensity": "RESULT_criticalCurrentDensity",
268 "heatFlux": "RESULT_heatFlux",
269 "resistivity": "RESULT_resistivity",
270 "thermalConductivity": "RESULT_thermalConductivity",
271 "specificHeatCapacity": "RESULT_specificHeatCapacity",
272 "jHTSOverjCritical": "RESULT_jHTSOverjCritical",
273 "criticalCurrent": "RESULT_criticalCurrent",
274 "debug": "RESULT_debug",
275 "inductance": "RESULT_inductance",
276 "timeConstant": "RESULT_timeConstant",
277 "axialComponentOfTheMagneticField": "RESULT_axialComponentOfTheMagneticField",
278 "maximumTemperature": "RESULT_maximumTemperature",
279 "jHTS": "RESULT_jHTS",
280 "currentSharingIndex": "RESULT_currentSharingIndex",
281 "cryocoolerAveragePower": "RESULT_cryocoolerAveragePower",
282 "arcLength": "RESULT_arcLength",
283 "turnNumber": "RESULT_turnNumber",
284 "cryocoolerAverageTemperature": "RESULT_cryocoolerAverageTemperature"
285}
287getdpPostOperationNames = {
288 "magneticField": "POSTOP_magneticField",
289 "magneticEnergy": "RESULT_magneticEnergy",
290 "magnitudeOfMagneticField": "POSTOP_magnitudeOfMagneticField",
291 "currentDensity": "POSTOP_currentDensity",
292 "magnitudeOfCurrentDensity": "POSTOP_magnitudeOfCurrentDensity",
293 "resistiveHeating": "POSTOP_resistiveHeating",
294 "totalResistiveHeating": "POSTOP_totalResistiveHeating",
295 "temperature": "POSTOP_temperature",
296 "currentThroughCoil": "POSTOP_currentThroughCoil",
297 "voltageBetweenTerminals": "POSTOP_voltageBetweenTerminals",
298 "criticalCurrentDensity": "POSTOP_criticalCurrentDensity",
299 "heatFlux": "POSTOP_heatFlux",
300 "resistivity": "POSTOP_resistivity",
301 "thermalConductivity": "POSTOP_thermalConductivity",
302 "specificHeatCapacity": "POSTOP_specificHeatCapacity",
303 "jHTSOverjCritical": "POSTOP_jHTSOverjCritical",
304 "criticalCurrent": "POSTOP_criticalCurrent",
305 "debug": "POSTOP_debug",
306 "inductance": "POSTOP_inductance",
307 "timeConstant": "POSTOP_timeConstant",
308 "axialComponentOfTheMagneticField": "POSTOP_axialComponentOfTheMagneticField",
309 "maximumTemperature": "POSTOP_maximumTemperature",
310 "jHTS": "POSTOP_jHTS",
311 "currentSharingIndex": "POSTOP_currentSharingIndex",
312 "cryocoolerAveragePower": "POSTOP_cryocoolerAveragePower",
313 "arcLength": "POSTOP_arcLength",
314 "turnNumber": "POSTOP_turnNumber",
315 "cryocoolerAverageTemperature": "POSTOP_cryocoolerAverageTemperature"
316}
318# ======================================================================================
319# ======================================================================================
321# Global variables
322geometry_input = ContextVar("geometry")
323mesh_input = ContextVar("mesh")
324solve_input = ContextVar("solve")
325input_file_path = ContextVar("input_file_path")
326all_break_points = []
329def getWindingOuterRadius():
330 """Return outer radius of the winding."""
331 geometry = geometry_input.get()
332 return (
333 geometry["winding"]["innerRadius"]
334 + geometry["winding"]["thickness"]
335 + geometry["winding"]["numberOfTurns"]
336 * (geometry["winding"]["thickness"] + geometry["contactLayer"]["thickness"])
337 )
340def getAirHeight():
341 """Return the height of the air."""
342 geometry = geometry_input.get()
343 h = (
344 geometry["numberOfPancakes"]
345 * (geometry["winding"]["height"] + geometry["gapBetweenPancakes"])
346 - geometry["gapBetweenPancakes"]
347 + 2 * geometry["air"]["axialMargin"]
348 )
349 return h
352def getTransitionNotchAngle():
353 """Return transition notch angle of the winding."""
354 mesh = mesh_input.get()
356 azimuthalNumberOfElementsPerTurn = max(
357 mesh["winding"]["azimuthalNumberOfElementsPerTurn"]
358 )
360 transitionNotchAngle = 2 * math.pi / azimuthalNumberOfElementsPerTurn
362 return transitionNotchAngle
365def checkIfAirOrTerminalMeshIsStructured():
366 geometry = geometry_input.get()
367 mesh = mesh_input.get()
369 structuredAirMesh = False
370 structuredTerminalMesh = False
371 if "air" in mesh:
372 structuredAirMesh = mesh["air"]["structured"]
373 if "terminals" in mesh:
374 structuredTerminalMesh = mesh["terminals"]["structured"]
375 structuredMesh = structuredAirMesh or structuredTerminalMesh
377 return structuredMesh
380# ======================================================================================
381# FUNDAMENTAL CLASSES STARTS ===========================================================
382# ======================================================================================
383class Pancake3DPositionInCoordinates(BaseModel):
384 x: float = Field(
385 title="x coordinate",
386 description="x coordinate of the position.",
387 )
388 y: float = Field(
389 title="y coordinate",
390 description="y coordinate of the position.",
391 )
392 z: float = Field(
393 title="z coordinate",
394 description="z coordinate of the position.",
395 )
398class Pancake3DPositionInTurnNumbers(BaseModel):
399 turnNumber: float = Field(
400 title="Turn Number",
401 description=(
402 "Winding turn number as a position input. It starts from 0 and it can be a"
403 " float."
404 ),
405 )
406 whichPancakeCoil: Optional[PositiveInt] = Field(
407 default=None,
408 title="Pancake Coil Number",
409 description="The first pancake coil is 1, the second is 2, etc.",
410 )
412 @field_validator("turnNumber")
413 @classmethod
414 def check_turnNumber(cls, turnNumber):
415 geometry = geometry_input.get()
417 if turnNumber < 0:
418 raise ValueError("Turn number cannot be less than 0.")
419 elif turnNumber > geometry["winding"]["numberOfTurns"]:
420 raise ValueError(
421 "Turn number cannot be greater than the number of turns of the winding"
422 f" ({geometry['numberOfPancakes']})."
423 )
425 return turnNumber
427 @field_validator("whichPancakeCoil")
428 @classmethod
429 def check_whichPancakeCoil(cls, whichPancakeCoil):
430 geometry = geometry_input.get()
432 if whichPancakeCoil is not None:
433 if whichPancakeCoil < 1:
434 raise ValueError(
435 "Pancake coil numbers start from 1. Therefore, it cannot be less"
436 " than 1."
437 )
438 elif whichPancakeCoil > geometry["numberOfPancakes"]:
439 raise ValueError(
440 "Pancake coil number cannot be greater than the number of pancakes"
441 f" ({geometry['numberOfPancakes']})."
442 )
443 else:
444 return 1
446 return whichPancakeCoil
448 def compute_coordinates(self):
449 geometry = geometry_input.get()
450 mesh = mesh_input.get()
452 if geometry["contactLayer"]["thinShellApproximation"]:
453 windingThickness = (
454 geometry["winding"]["thickness"]
455 + geometry["contactLayer"]["thickness"]
456 * (geometry["winding"]["numberOfTurns"] - 1)
457 / geometry["winding"]["numberOfTurns"]
458 )
459 gapThickness = 0
460 else:
461 windingThickness = geometry["winding"]["thickness"]
462 gapThickness = geometry["contactLayer"]["thickness"]
464 innerRadius = geometry["winding"]["innerRadius"]
465 initialTheta = 0.0
466 if isinstance(mesh["winding"]["azimuthalNumberOfElementsPerTurn"], list):
467 ane = mesh["winding"]["azimuthalNumberOfElementsPerTurn"][0]
468 elif isinstance(mesh["winding"]["azimuthalNumberOfElementsPerTurn"], int):
469 ane = mesh["winding"]["azimuthalNumberOfElementsPerTurn"]
470 else:
471 raise ValueError(
472 "The azimuthal number of elements per turn must be either an integer"
473 " or a list of integers."
474 )
476 numberOfPancakes = geometry["numberOfPancakes"]
477 gapBetweenPancakes = geometry["gapBetweenPancakes"]
478 windingHeight = geometry["winding"]["height"]
480 turnNumber = self.turnNumber
481 whichPancake = self.whichPancakeCoil
483 elementStartTurnNumber = math.floor(turnNumber / (1 / ane)) * (1 / ane)
484 elementEndTurnNumber = elementStartTurnNumber + 1 / ane
486 class point:
487 def __init__(self, x, y, z):
488 self.x = x
489 self.y = y
490 self.z = z
492 def __add__(self, other):
493 return point(self.x + other.x, self.y + other.y, self.z + other.z)
495 def __sub__(self, other):
496 return point(self.x - other.x, self.y - other.y, self.z - other.z)
498 def __mul__(self, scalar):
499 return point(self.x * scalar, self.y * scalar, self.z * scalar)
501 def __truediv__(self, scalar):
502 return point(self.x / scalar, self.y / scalar, self.z / scalar)
504 def rotate(self, degrees):
505 return point(
506 self.x * math.cos(degrees) - self.y * math.sin(degrees),
507 self.x * math.sin(degrees) + self.y * math.cos(degrees),
508 self.z,
509 )
511 def normalize(self):
512 return self / math.sqrt(self.x**2 + self.y**2 + self.z**2)
514 if whichPancake % 2 == 1:
515 # If the spiral is counter-clockwise, the initial theta angle decreases,
516 # and r increases as the theta angle decreases.
517 multiplier = 1
518 elif whichPancake % 2 == 0:
519 # If the spiral is clockwise, the initial theta angle increases, and r
520 # increases as the theta angle increases.
521 multiplier = -1
523 # Mesh element's starting point:
524 elementStartTheta = 2 * math.pi * elementStartTurnNumber * multiplier
525 elementStartRadius = (
526 innerRadius
527 + elementStartTheta
528 / (2 * math.pi)
529 * (gapThickness + windingThickness)
530 * multiplier
531 )
532 elementStartPointX = elementStartRadius * math.cos(
533 initialTheta + elementStartTheta
534 )
535 elementStartPointY = elementStartRadius * math.sin(
536 initialTheta + elementStartTheta
537 )
538 elementStartPointZ = (
539 -(
540 numberOfPancakes * windingHeight
541 + (numberOfPancakes - 1) * gapBetweenPancakes
542 )
543 / 2
544 + windingHeight / 2
545 + (whichPancake - 1) * (windingHeight + gapBetweenPancakes)
546 )
547 elementStartPoint = point(
548 elementStartPointX, elementStartPointY, elementStartPointZ
549 )
551 # Mesh element's ending point:
552 elementEndTheta = 2 * math.pi * elementEndTurnNumber * multiplier
553 elementEndRadius = (
554 innerRadius
555 + elementEndTheta
556 / (2 * math.pi)
557 * (gapThickness + windingThickness)
558 * multiplier
559 )
560 elementEndPointX = elementEndRadius * math.cos(initialTheta + elementEndTheta)
561 elementEndPointY = elementEndRadius * math.sin(initialTheta + elementEndTheta)
562 elementEndPointZ = elementStartPointZ
563 elementEndPoint = point(elementEndPointX, elementEndPointY, elementEndPointZ)
565 turnNumberFraction = (turnNumber - elementStartTurnNumber) / (
566 elementEndTurnNumber - elementStartTurnNumber
567 )
568 location = (
569 elementStartPoint
570 + (elementEndPoint - elementStartPoint) * turnNumberFraction
571 ) + (elementEndPoint - elementStartPoint).rotate(
572 -math.pi / 2
573 ).normalize() * windingThickness / 2 * multiplier
575 return location.x, location.y, location.z
577 @computed_field
578 @cached_property
579 def x(self) -> float:
580 return self.compute_coordinates()[0]
582 @computed_field
583 @cached_property
584 def y(self) -> float:
585 return self.compute_coordinates()[1]
587 @computed_field
588 @cached_property
589 def z(self) -> float:
590 return self.compute_coordinates()[2]
593Pancake3DPosition = Pancake3DPositionInCoordinates | Pancake3DPositionInTurnNumbers
595# ======================================================================================
596# FUNDAMENTAL CLASSES ENDS =============================================================
597# ======================================================================================
600# ======================================================================================
601# GEOMETRY CLASSES STARTS ==============================================================
602# ======================================================================================
603class Pancake3DGeometryWinding(BaseModel):
604 # Mandatory:
605 innerRadius: PositiveFloat = Field(
606 title="Inner Radius",
607 description="Inner radius of the winding.",
608 )
609 thickness: PositiveFloat = Field(
610 title="Winding Thickness",
611 description="Thickness of the winding.",
612 )
613 numberOfTurns: float = Field(
614 ge=3,
615 title="Number of Turns",
616 description="Number of turns of the winding.",
617 )
618 height: PositiveFloat = Field(
619 title="Winding Height",
620 description="Height/width of the winding.",
621 )
623 # Optionals:
624 name: str = Field(
625 default="winding",
626 title="Winding Name",
627 description="The The name to be used in the mesh..",
628 examples=["winding", "myWinding"],
629 )
630 numberOfVolumesPerTurn: int = Field(
631 default=2,
632 validate_default=True,
633 ge=2,
634 title="Number of Volumes Per Turn (Advanced Input)",
635 description="The number of volumes per turn (CAD related, not physical).",
636 )
638 @field_validator("numberOfVolumesPerTurn")
639 @classmethod
640 def check_numberOfVolumesPerTurn(cls, numberOfVolumesPerTurn):
641 geometry = geometry_input.get()
642 mesh = mesh_input.get()
644 # Check if the numberOfVolumesPerTurn is compatible swith the azimuthal number of
645 # elements per turn:
646 for i, ane in enumerate(mesh["winding"]["azimuthalNumberOfElementsPerTurn"]):
647 if ane % numberOfVolumesPerTurn != 0:
648 raise ValueError(
649 "The azimuthal number of elements per turn for the pancake coil"
650 f" number {i+1} is ({ane}), but it must be divisible by the number"
651 f" of volumes per turn ({geometry['winding']['numberOfVolumesPerTurn']})!"
652 " So it needs to be rounded to"
653 f" {math.ceil(ane/numberOfVolumesPerTurn)*numberOfVolumesPerTurn:.5f} or"
654 f" {math.floor(ane/numberOfVolumesPerTurn)*numberOfVolumesPerTurn:.5f}."
655 )
657 structured = checkIfAirOrTerminalMeshIsStructured()
659 if structured:
660 # If the mesh is structured, the number of volumes per turn must be 4:
661 numberOfVolumesPerTurn = 4
663 return numberOfVolumesPerTurn
665 @computed_field
666 @cached_property
667 def theta_i(self) -> float:
668 """Return start angle of the winding."""
669 return 0.0
671 @computed_field
672 @cached_property
673 def r_o(self) -> float:
674 """Return outer radius of the winding."""
675 return getWindingOuterRadius()
677 @computed_field
678 @cached_property
679 def turnTol(self) -> float:
680 """Return turn tolerance of the winding."""
681 geometry: Pancake3DGeometry = geometry_input.get()
682 mesh: Pancake3DMesh = mesh_input.get()
684 # Calculate the turn tolerance required due to the geometrymetry input:
685 # Turn tolerance is the smallest turn angle (in turns) that is allowed.
686 if "dimTol" in geometry:
687 dimTol = geometry["dimTol"]
688 else:
689 dimTol = 1e-8
690 turnTol = geometry["winding"]["numberOfTurns"] % 1
691 if math.isclose(turnTol, 0, abs_tol=dimTol):
692 turnTol = 0.5
694 turnTolDueToTransition = getTransitionNotchAngle() / (2 * math.pi)
696 # Calculate the minimum turn tolerance possible due to the mesh input:
697 minimumTurnTol = 1 / min(mesh["winding"]["azimuthalNumberOfElementsPerTurn"])
699 if turnTol < minimumTurnTol:
700 numberOfTurns = geometry["winding"]["numberOfTurns"]
702 raise ValueError(
703 "The azimuthal number of elements per turn for one of the pancakes is"
704 f" {min(mesh['winding']['azimuthalNumberOfElementsPerTurn'])}, and the"
705 " number of turns is"
706 f" {numberOfTurns:.5f}."
707 " The number of turns must always be divisible by the (1/(the"
708 " azimuthal number of elements per turn)) to ensure conformality."
709 " Please change the number of turns or the azimuthal number of"
710 " elemenets per turn. The closest possible number of turns value is"
711 f" {round(numberOfTurns * min(mesh['winding']['azimuthalNumberOfElementsPerTurn']))/min(mesh['winding']['azimuthalNumberOfElementsPerTurn']):.5f}"
712 )
713 else:
714 # Minimum possible sections per turn is 16 (otherwise splines might collide
715 # into each other). But it should be greater than the number of volumes per
716 # turn and it should be divisible by both 1/turnTol and the number of
717 # volumes per turn.
718 sectionsPerTurn = 16
719 if "numberOfVolumesPerTurn" in geometry["winding"]:
720 numberOfVolumesPerTurn = geometry["winding"]["numberOfVolumesPerTurn"]
721 else:
722 numberOfVolumesPerTurn = 2
724 while (
725 (
726 math.fmod(sectionsPerTurn, (1 / turnTol)) > 1e-8
727 and math.fmod(sectionsPerTurn, (1 / turnTol)) - (1 / turnTol)
728 < -1e-8
729 )
730 or (
731 math.fmod(sectionsPerTurn, (1 / turnTolDueToTransition)) > 1e-8
732 and math.fmod(sectionsPerTurn, (1 / turnTolDueToTransition))
733 - (1 / turnTolDueToTransition)
734 < -1e-8
735 )
736 or sectionsPerTurn % numberOfVolumesPerTurn != 0
737 or sectionsPerTurn < numberOfVolumesPerTurn
738 ):
739 sectionsPerTurn += 1
741 # Sections per turn will set the turn tolerance value as well.
742 return 1.0 / sectionsPerTurn
744 @computed_field
745 @cached_property
746 def spt(self) -> float:
747 """Return sections per turn of the winding."""
748 return int(1.0 / self.turnTol)
750 @computed_field
751 @cached_property
752 def totalTapeLength(self) -> float:
753 """Return total tape length of the winding."""
754 geometry: Pancake3DGeometry = geometry_input.get()
756 # Calculate the total tape length of the coil:
758 # The same angle can be subtracted from both theta_1 and theta_2 to simplify the
759 # calculations:
760 theta2 = geometry["winding"]["numberOfTurns"] * 2 * math.pi
761 theta1 = 0
763 # Since r = a * theta + b, r_1 = b since theta_1 = 0:
764 b = geometry["winding"]["innerRadius"]
766 # Since r = a * theta + b, r_2 = a * theta2 + b:
767 a = (getWindingOuterRadius() - b) / theta2
769 def integrand(t):
770 return math.sqrt(a**2 + (a * t + b) ** 2)
772 totalTapeLength = abs(scipy.integrate.quad(integrand, theta1, theta2)[0])
774 return totalTapeLength
777class Pancake3DGeometryContactLayer(BaseModel):
778 # Mandatory:
779 thinShellApproximation: bool = Field(
780 title="Use Thin Shell Approximation",
781 description=(
782 "If True, the contact layer will be modeled with 2D shell elements (thin"
783 " shell approximation), and if False, the contact layer will be modeled"
784 " with 3D elements."
785 ),
786 )
787 thickness: PositiveFloat = Field(
788 title="Contact Layer Thickness",
789 description=("Thickness of the contact layer."
790 "It is the total thickness of the contact or insulation layer."
791 "In particular, for perfect insulation this would be the sum of the insulation layer of the two adjacent CC with an insulation layer of "
792 "thickness t/2 on each side."
793 ),
794 )
796 # Optionals:
797 name: str = Field(
798 default="contactLayer",
799 title="Contact Layer Name",
800 description="The name to be used in the mesh.",
801 examples=["myContactLayer"],
802 )
805class Pancake3DGeometryTerminalBase(BaseModel):
806 # Mandatory:
807 thickness: PositiveFloat = Field(
808 title="Terminal Thickness",
809 description="Thickness of the terminal's tube.",
810 ) # thickness
812 @field_validator("thickness")
813 @classmethod
814 def check_t(cls, thickness):
815 geometry = geometry_input.get()
817 if thickness < geometry["winding"]["thickness"] / 2:
818 raise ValueError(
819 "Terminal's thickness is smaller than half of the winding's thickness!"
820 " Please increase the terminal's thickness."
821 )
823 return thickness
826class Pancake3DGeometryInnerTerminal(Pancake3DGeometryTerminalBase):
827 name: str = Field(
828 default="innerTerminal",
829 title="Terminal Name",
830 description="The name to be used in the mesh.",
831 examples=["innerTerminal", "outerTeminal"],
832 )
834 @computed_field
835 @cached_property
836 def r(self) -> float:
837 """Return inner radius of the inner terminal."""
838 geometry = geometry_input.get()
840 innerRadius = geometry["winding"]["innerRadius"] - 2 * self.thickness
841 if innerRadius < 0:
842 raise ValueError(
843 "Inner terminal's radius is smaller than 0! Please decrease the inner"
844 " terminal's thickness or increase the winding's inner radius."
845 )
847 return innerRadius
850class Pancake3DGeometryOuterTerminal(Pancake3DGeometryTerminalBase):
851 name: str = Field(
852 default="outerTerminal",
853 title="Terminal Name",
854 description="The name to be used in the mesh.",
855 examples=["innerTerminal", "outerTeminal"],
856 )
858 @computed_field
859 @cached_property
860 def r(self) -> float:
861 """Return outer radius of the outer terminal."""
862 outerRadius = getWindingOuterRadius() + 2 * self.thickness
864 return outerRadius
867class Pancake3DGeometryTerminals(BaseModel):
868 # 1) User inputs:
869 inner: Pancake3DGeometryInnerTerminal = Field()
870 outer: Pancake3DGeometryOuterTerminal = Field()
872 # Optionals:
873 firstName: str = Field(
874 default="firstTerminal", description="name of the first terminal"
875 )
876 lastName: str = Field(
877 default="lastTerminal", description="name of the last terminal"
878 )
880 @computed_field
881 @cached_property
882 def transitionNotchAngle(self) -> float:
883 """Return transition notch angle of the terminals."""
884 return getTransitionNotchAngle()
887class Pancake3DGeometryAirBase(BaseModel):
888 # Mandatory:
889 axialMargin: PositiveFloat = Field(
890 title="Axial Margin of the Air",
891 description=(
892 "Axial margin between the ends of the air and first/last pancake coils."
893 ),
894 ) # axial margin
896 # Optionals:
897 name: str = Field(
898 default="air",
899 title="Air Name",
900 description="The name to be used in the mesh.",
901 examples=["air", "myAir"],
902 )
903 shellTransformation: bool = Field(
904 default=False,
905 title="Use Shell Transformation",
906 description=(
907 "Generate outer shell air to apply shell transformation if True (GetDP"
908 " related, not physical)"
909 ),
910 )
911 shellTransformationMultiplier: float = Field(
912 default=1.2,
913 gt=1.1,
914 title="Shell Transformation Multiplier (Advanced Input)",
915 description=(
916 "multiply the air's outer dimension by this value to get the shell's outer"
917 " dimension"
918 ),
919 )
920 cutName: str = Field(
921 default="Air-Cut",
922 title="Air Cut Name",
923 description="name of the cut (cochain) to be used in the mesh",
924 examples=["Air-Cut", "myAirCut"],
925 )
926 shellVolumeName: str = Field(
927 default="air-Shell",
928 title="Air Shell Volume Name",
929 description="name of the shell volume to be used in the mesh",
930 examples=["air-Shell", "myAirShell"],
931 )
932 generateGapAirWithFragment: bool = Field(
933 default=False,
934 title="Generate Gap Air with Fragment (Advanced Input)",
935 description=(
936 "generate the gap air with gmsh/model/occ/fragment if true (CAD related,"
937 " not physical)"
938 ),
939 )
941 @field_validator("axialMargin")
942 @classmethod
943 def check_axialMargin(cls, axialMargin):
944 geometry = geometry_input.get()
945 windingHeight = geometry["winding"]["height"]
947 if axialMargin < windingHeight / 2:
948 raise ValueError(
949 "Axial margin is smaller than half of the winding's height! Please"
950 " increase the axial margin."
951 )
953 return axialMargin
955 @computed_field
956 @cached_property
957 def h(self) -> float:
958 """Return total height of the air."""
959 h = getAirHeight()
961 return h
964class Pancake3DGeometryAirCylinder(Pancake3DGeometryAirBase):
965 type: Literal["cylinder"] = Field(default="cylinder", title="Air Type")
966 radius: PositiveFloat = Field(
967 default=None,
968 title="Air Radius",
969 description="Radius of the air (for cylinder type air).",
970 )
972 @field_validator("radius")
973 @classmethod
974 def check_r(cls, radius):
975 geometry = geometry_input.get()
976 outerTerminalOuterRadius = (
977 getWindingOuterRadius() + 2 * geometry["terminals"]["outer"]["thickness"]
978 )
980 if radius < outerTerminalOuterRadius * 1.5:
981 raise ValueError(
982 "Radius of the air must be at least 1.5 times the outer radius of the"
983 " winding! Please increase the radius of the air."
984 )
986 return radius
988 @computed_field
989 @cached_property
990 def shellOuterRadius(self) -> float:
991 """Return outer radius of the air."""
992 shellOuterRadius = self.shellTransformationMultiplier * self.radius
994 return shellOuterRadius
997class Pancake3DGeometryAirCuboid(Pancake3DGeometryAirBase):
998 type: Literal["cuboid"] = Field(default="cuboid", title="Air Type")
999 sideLength: PositiveFloat = Field(
1000 default=None,
1001 title="Air Side Length",
1002 description="Side length of the air (for cuboid type air).",
1003 )
1005 @field_validator("sideLength")
1006 @classmethod
1007 def check_a(cls, sideLength):
1008 geometry = geometry_input.get()
1009 outerTerminalOuterRadius = (
1010 getWindingOuterRadius() + 2 * geometry["terminals"]["outer"]["thickness"]
1011 )
1013 if sideLength / 2 < outerTerminalOuterRadius * 1.5:
1014 raise ValueError(
1015 "Half of the side length of the air must be at least 1.5 times the"
1016 " outer radius of the winding! Please increase the side length of the"
1017 " air."
1018 )
1020 return sideLength
1022 @computed_field
1023 @cached_property
1024 def shellSideLength(self) -> float:
1025 """Return outer radius of the air."""
1026 shellSideLength = self.shellTransformationMultiplier * self.sideLength
1028 return shellSideLength
1031Pancake3DGeometryAir = Annotated[
1032 Pancake3DGeometryAirCylinder | Pancake3DGeometryAirCuboid,
1033 Field(discriminator="type"),
1034]
1035# ======================================================================================
1036# GEOMETRY CLASSES ENDS ================================================================
1037# ======================================================================================
1040# ======================================================================================
1041# MESH CLASSES STARTS ==================================================================
1042# ======================================================================================
1043class Pancake3DMeshWinding(BaseModel):
1044 # Mandatory:
1045 axialNumberOfElements: list[PositiveInt] | PositiveInt = Field(
1046 title="Axial Number of Elements",
1047 description=(
1048 "The number of axial elements for the whole height of the coil. It can be"
1049 " either a list of integers to specify the value for each pancake coil"
1050 " separately or an integer to use the same setting for each pancake coil."
1051 ),
1052 )
1054 azimuthalNumberOfElementsPerTurn: list[PositiveInt] | PositiveInt = Field(
1055 title="Azimuthal Number of Elements Per Turn",
1056 description=(
1057 "The number of azimuthal elements per turn of the coil. It can be either a"
1058 " list of integers to specify the value for each pancake coil separately or"
1059 " an integer to use the same setting for each pancake coil."
1060 ),
1061 )
1063 radialNumberOfElementsPerTurn: list[PositiveInt] | PositiveInt = Field(
1064 title="Winding Radial Number of Elements Per Turn",
1065 description=(
1066 "The number of radial elements per tape of the winding. It can be either a"
1067 " list of integers to specify the value for each pancake coil separately or"
1068 " an integer to use the same setting for each pancake coil."
1069 ),
1070 )
1072 # Optionals:
1073 axialDistributionCoefficient: list[PositiveFloat] | PositiveFloat = Field(
1074 default=[1],
1075 title="Axial Bump Coefficients",
1076 description=(
1077 "If 1, it won't affect anything. If smaller than 1, elements will get finer"
1078 " in the axial direction at the ends of the coil. If greater than 1,"
1079 " elements will get coarser in the axial direction at the ends of the coil."
1080 " It can be either a list of floats to specify the value for each pancake"
1081 " coil separately or a float to use the same setting for each pancake coil."
1082 ),
1083 )
1085 elementType: (
1086 list[Literal["tetrahedron", "hexahedron", "prism"]]
1087 | Literal["tetrahedron", "hexahedron", "prism"]
1088 ) = Field(
1089 default=["tetrahedron"],
1090 title="Element Type",
1091 description=(
1092 "The element type of windings and contact layers. It can be either a"
1093 " tetrahedron, hexahedron, or a prism. It can be either a list of strings"
1094 " to specify the value for each pancake coil separately or a string to use"
1095 " the same setting for each pancake coil."
1096 ),
1097 )
1099 @field_validator("axialNumberOfElements", "azimuthalNumberOfElementsPerTurn", "radialNumberOfElementsPerTurn", "axialDistributionCoefficient", "elementType")
1100 @classmethod
1101 def check_inputs(cls, value, info: ValidationInfo):
1102 geometry = geometry_input.get()
1104 numberOfPancakes = geometry["numberOfPancakes"]
1106 structuredMesh = checkIfAirOrTerminalMeshIsStructured()
1108 if not isinstance(value, list):
1109 value = [value] * numberOfPancakes
1110 elif len(value) == 1:
1111 value = value * numberOfPancakes
1112 else:
1113 if len(value) != numberOfPancakes:
1114 raise ValueError(
1115 "The length of the input list must be equal to the number of"
1116 " pancake coils!"
1117 )
1118 if info.field_name == "ane":
1119 if value[0] < 7:
1120 raise ValueError(
1121 "Azimuthal number of elements per turn must be greater than or"
1122 " equal to 7!"
1123 )
1125 if structuredMesh:
1126 if len(set(value)) != 1:
1127 raise ValueError(
1128 "If structured mesh is used, the same mesh setting must be used for"
1129 " all pancake coils!"
1130 )
1132 if info.field_name == "elementType":
1133 if value[0] != "tetrahedron":
1134 raise ValueError(
1135 "If structured air or terminal mesh is used, the element type"
1136 " must be tetrahedron!"
1137 )
1139 if info.field_name == "ane":
1140 if value[0] % 4 != 0:
1141 raise ValueError(
1142 "If structured mesh is used, the number of azimuthal elements"
1143 " per turn must be divisible by 4!"
1144 )
1146 return value
1149class Pancake3DMeshContactLayer(BaseModel):
1150 # Mandatory:
1151 radialNumberOfElementsPerTurn: list[PositiveInt] = Field(
1152 title="Contact Layer Radial Number of Elements Per Turn",
1153 description=(
1154 "The number of radial elements per tape of the contact layer. It can be"
1155 " either a list of integers to specify the value for each pancake coil"
1156 " separately or an integer to use the same setting for each pancake coil."
1157 ),
1158 )
1160 @field_validator("radialNumberOfElementsPerTurn")
1161 @classmethod
1162 def check_inputs(cls, value):
1163 geometry = geometry_input.get()
1165 structuredMesh = checkIfAirOrTerminalMeshIsStructured()
1167 numberOfPancakeCoils = geometry["numberOfPancakes"]
1169 if not isinstance(value, list):
1170 value = [value] * numberOfPancakeCoils
1171 elif len(value) == 1:
1172 value = value * numberOfPancakeCoils
1173 else:
1174 if len(value) != numberOfPancakeCoils:
1175 raise ValueError(
1176 "The length of the input list must be equal to the number of"
1177 " pancake coils!"
1178 )
1180 if structuredMesh:
1181 if len(set(value)) != 1:
1182 raise ValueError(
1183 "If structured mesh is used, the same mesh setting must be used for"
1184 " all pancake coils!"
1185 )
1187 return value
1190class Pancake3DMeshAirAndTerminals(BaseModel):
1191 # Optionals:
1192 structured: bool = Field(
1193 default=False,
1194 title="Structure Mesh",
1195 description=(
1196 "If True, the mesh will be structured. If False, the mesh will be"
1197 " unstructured."
1198 ),
1199 )
1200 radialElementSize: PositiveFloat = Field(
1201 default=1,
1202 title="Radial Element Size",
1203 description=(
1204 "If structured mesh is used, the radial element size can be set. It is the"
1205 " radial element size in terms of the winding's radial element size."
1206 ),
1207 )
1210# ======================================================================================
1211# MESH CLASSES ENDS ====================================================================
1212# ======================================================================================
1215# ======================================================================================
1216# SOLVE CLASSES STARTS =================================================================
1217# ======================================================================================
1218class Pancake3DSolveAir(BaseModel):
1219 # 1) User inputs:
1221 # Mandatory:
1222 permeability: PositiveFloat = Field(
1223 title="Permeability of Air",
1224 description="Permeability of air.",
1225 )
1228class Pancake3DSolveIcVsLengthList(BaseModel):
1229 lengthValues: list[float] = Field(
1230 title="Tape Length Values",
1231 description="Tape length values that corresponds to criticalCurrentValues.",
1232 )
1233 criticalCurrentValues: list[float] = Field(
1234 title="Critical Current Values",
1235 description="Critical current values that corresponds to lengthValues.",
1236 )
1237 lengthUnit: str = Field(
1238 title="Unit",
1239 description=(
1240 "Unit of the critical current values. "
1241 "It can be either the arc length in meter or "
1242 "the number of turns."
1243 ),
1244 examples=["meter", "turnNumber"],
1245 )
1247class Pancake3DSolveIcVsLengthCSV(BaseModel):
1248 csvFile: str = Field(
1249 title="CSV File",
1250 description="The path of the CSV file that contains the critical current values.",
1251 )
1253 lengthUnit: str = Field(
1254 title="Unit",
1255 description=(
1256 "Unit of the critical current values. "
1257 "It can be either the arc length in meter or "
1258 "the number of turns."
1259 ),
1260 examples=["meter", "turnNumber"],
1261 )
1264class Pancake3DSolveMaterialBase(BaseModel):
1265 name: str
1267 # Optionals:
1268 RRR: PositiveFloat = Field(
1269 default=100,
1270 title="Residual Resistance Ratio",
1271 description=(
1272 "Residual-resistivity ratio (also known as Residual-resistance ratio or"
1273 " just RRR) is the ratio of the resistivity of a material at reference"
1274 " temperature and at 0 K."
1275 ),
1276 )
1277 RRRRefTemp: PositiveFloat = Field(
1278 default=295,
1279 title="Residual Resistance Ratio Reference Temperature",
1280 description="Reference temperature for residual resistance ratio",
1281 )
1283 @computed_field
1284 @cached_property
1285 def resistivityMacroName(self) -> str:
1286 """Return the resistivity macro name of the material."""
1287 if self.name not in resistivityMacroNames:
1288 return "NOT_DEFINED_IN_DATA_FIQUS_PANCAKE3D"
1290 return resistivityMacroNames[self.name]
1292 @computed_field
1293 @cached_property
1294 def thermalConductivityMacroName(self) -> str:
1295 """Return the thermal conductivity macro name of the material."""
1296 if self.name not in thermalConductivityMacroNames:
1297 return "NOT_DEFINED_IN_DATA_FIQUS_PANCAKE3D"
1299 return thermalConductivityMacroNames[self.name]
1301 @computed_field
1302 @cached_property
1303 def heatCapacityMacroName(self) -> str:
1304 """Return the heat capacity macro name of the material."""
1305 if self.name not in heatCapacityMacroNames:
1306 return "NOT_DEFINED_IN_DATA_FIQUS_PANCAKE3D"
1308 return heatCapacityMacroNames[self.name]
1310 @computed_field
1311 @cached_property
1312 def getdpTSAOnlyResistivityFunction(self) -> str:
1313 """Return the GetDP function name of the material's resistivity."""
1314 if self.name not in getdpTSAOnlyResistivityFunctions:
1315 return "NOT_DEFINED_IN_DATA_FIQUS_PANCAKE3D"
1317 return getdpTSAOnlyResistivityFunctions[self.name]
1319 @computed_field
1320 @cached_property
1321 def getdpTSAMassResistivityMacroName(self) -> str:
1322 """Return the GetDP function name of the material's mass resistivity."""
1323 if self.name not in getdpTSAMassResistivityMacroNames:
1324 return "NOT_DEFINED_IN_DATA_FIQUS_PANCAKE3D"
1326 return getdpTSAMassResistivityMacroNames[self.name]
1328 @computed_field
1329 @cached_property
1330 def getdpTSAStiffnessResistivityMacroName(self) -> str:
1331 """Return the GetDP function name of the material's stiffness resistivity."""
1332 if self.name not in getdpTSAStiffnessResistivityMacroNames:
1333 return "NOT_DEFINED_IN_DATA_FIQUS_PANCAKE3D"
1335 return getdpTSAStiffnessResistivityMacroNames[self.name]
1337 @computed_field
1338 @cached_property
1339 def getdpTSAMassThermalConductivityMacroName(self) -> str:
1340 """Return the GetDP function name of the material's mass thermal conductivity."""
1341 if self.name not in getdpTSAMassThermalConductivityMacroNames:
1342 return "NOT_DEFINED_IN_DATA_FIQUS_PANCAKE3D"
1344 return getdpTSAMassThermalConductivityMacroNames[self.name]
1346 @computed_field
1347 @cached_property
1348 def getdpTSAStiffnessThermalConductivityMacroName(self) -> str:
1349 """Return the GetDP function name of the material's stiffness thermal conductivity."""
1350 if self.name not in getdpTSAStiffnessThermalConductivityMacroNames:
1351 return "NOT_DEFINED_IN_DATA_FIQUS_PANCAKE3D"
1353 return getdpTSAStiffnessThermalConductivityMacroNames[self.name]
1355 @computed_field
1356 @cached_property
1357 def getdpTSAMassHeatCapacityMacroName(self) -> str:
1358 """Return the GetDP function name of the material's mass heat capacity."""
1359 if self.name not in getdpTSAMassHeatCapacityMacroNames:
1360 return "NOT_DEFINED_IN_DATA_FIQUS_PANCAKE3D"
1362 return getdpTSAMassHeatCapacityMacroNames[self.name]
1364 @computed_field
1365 @cached_property
1366 def getdpTSARHSFunction(self) -> str:
1367 """Return the GetDP function name of the material's RHS."""
1368 if self.name not in getdpTSARHSFunctions:
1369 return "NOT_DEFINED_IN_DATA_FIQUS_PANCAKE3D"
1371 return getdpTSARHSFunctions[self.name]
1373 @computed_field
1374 @cached_property
1375 def getdpTSATripleFunction(self) -> str:
1376 """Return the GetDP function name of the material's triple."""
1377 if self.name not in getdpTSATripleFunctions:
1378 return "NOT_DEFINED_IN_DATA_FIQUS_PANCAKE3D"
1380 return getdpTSATripleFunctions[self.name]
1383class Pancake3DSolveNormalMaterial(Pancake3DSolveMaterialBase):
1384 # Mandatory:
1385 name: NormalMaterialName = Field(
1386 title="Material Name",
1387 )
1389 @computed_field
1390 @cached_property
1391 def getdpNormalMaterialGetDPName(self) -> str:
1392 """Return the GetDP material name of the normal conducting material. To make names with spaces like 'Stainless Steel' more robust."""
1393 if self.name not in getdpNormalMaterialNames:
1394 raise ValueError(
1395 f"Normal conducting material '{self.name}' is not defined"
1396 " in FiQuS!"
1397 )
1399 return getdpNormalMaterialNames[self.name]
1402class Pancake3DSolveSuperconductingMaterial(Pancake3DSolveMaterialBase):
1403 # Mandatory:
1404 name: SuperconductingMaterialName = Field(
1405 title="Superconducting Material Name",
1406 )
1407 nValue: PositiveFloat = Field(
1408 default=30,
1409 description="N-value for E-J power law.",
1410 )
1411 IcAtTAndBref: PositiveFloat | Pancake3DSolveIcVsLengthCSV | Pancake3DSolveIcVsLengthList = Field(
1412 title="Critical Current at Reference Temperature and Field in A",
1413 description=(
1414 "Critical current in A at reference temperature and magnetic field."
1415 "The critical current value will"
1416 " change with temperature depending on the superconductor material.Either"
1417 " the same critical current for the whole tape or the critical current with"
1418 " respect to the tape length can be specified. To specify the same critical"
1419 " current for the entire tape, just use a scalar. To specify critical"
1420 " current with respect to the tape length: a CSV file can be used, or"
1421 " lengthValues and criticalCurrentValues can be given as lists. The data"
1422 " will be linearly interpolated.If a CSV file is to be used, the input"
1423 " should be the name of a CSV file (which is in the same folder as the"
1424 " input file) instead of a scalar. The first column of the CSV file will be"
1425 " the tape length in m, and the second column will be the critical current in A. "
1426 ),
1427 examples=[230, "IcVSlength.csv"],
1428 )
1430 # Optionals:
1431 electricFieldCriterion: PositiveFloat = Field(
1432 default=1e-4,
1433 title="Electric Field Criterion",
1434 description=(
1435 "The electric field that defines the critical current density, i.e., the"
1436 " electric field at which the current density reaches the critical current"
1437 " density."
1438 ),
1439 )
1440 jCriticalScalingNormalToWinding: PositiveFloat = Field(
1441 default=1,
1442 title="Critical Current Scaling Normal to Winding",
1443 description=(
1444 "Critical current scaling normal to winding, i.e., along the c_axis. "
1445 " We have Jc_cAxis = scalingFactor * Jc_abPlane."
1446 " A factor of 1 means no scaling such that the HTS layer is isotropic."
1447 ),
1448 )
1450 IcReferenceTemperature: PositiveFloat = Field(
1451 default=77,
1452 title="Critical Current Reference Temperature",
1453 description="Critical current reference temperature in Kelvin.",
1454 )
1456 IcReferenceBmagnitude: NonNegativeFloat = Field(
1457 default=0.0,
1458 title="Critical Current Reference Magnetic Field Magnitude",
1459 description="Critical current reference magnetic field magnitude in Tesla.",
1460 )
1462 IcReferenceBangle: NonNegativeFloat = Field(
1463 default=90.0,
1464 title="Critical Current Reference Magnetic Field Angle",
1465 description= (
1466 "Critical current reference magnetic field angle in degrees."
1467 "0 degrees means the magnetic field is normal to the tape's wide surface"
1468 "and 90 degrees means the magnetic field is parallel to the tape's wide"
1469 "surface."
1470 ),
1471 )
1473 @computed_field
1474 @cached_property
1475 def IcValues(self) -> list[float]:
1476 """Return the critical current values of the material."""
1477 if hasattr(self.IcAtTAndBref, "criticalCurrentValues"):
1478 return self.IcAtTAndBref.criticalCurrentValues
1479 elif hasattr(self.IcAtTAndBref, "csvFile"):
1480 csv_file_path = pathlib.Path(input_file_path.get()).parent / self.IcAtTAndBref.csvFile
1481 # return the second column:
1482 IcValues = list(pd.read_csv(csv_file_path, header=None).iloc[:, 1])
1483 for Ic in IcValues:
1484 if Ic < 0:
1485 raise ValueError(
1486 "Critical current values in the CSV file should be positive!"
1487 )
1488 return IcValues
1489 else:
1490 return [self.IcAtTAndBref]
1492 @computed_field
1493 @cached_property
1494 def lengthValues(self) -> list[float]:
1495 """Return the length values of the material."""
1496 if hasattr(self.IcAtTAndBref, "lengthValues"):
1497 return self.IcAtTAndBref.lengthValues
1498 elif hasattr(self.IcAtTAndBref, "csvFile"):
1499 csv_file_path = pathlib.Path(input_file_path.get()).parent / self.IcAtTAndBref.csvFile
1500 # return the first column:
1501 lengthValues = list(pd.read_csv(csv_file_path, header=None).iloc[:, 0])
1502 for length in lengthValues:
1503 if length < 0:
1504 raise ValueError("Tape lengths in the CSV file should be positive!")
1505 return lengthValues
1506 else:
1507 return [1]
1509 @computed_field
1510 @cached_property
1511 def IcInArcLength(self) -> int:
1512 """Return 1 if Ic is given as a function of arc length in meters.
1513 Return 0 if Ic is given as a function of number of turns.
1514 Return -1 if Ic is given as a scalar value.
1515 """
1516 if hasattr(self.IcAtTAndBref, "lengthUnit"):
1517 if self.IcAtTAndBref.lengthUnit == "meter":
1518 return 1
1519 else:
1520 return 0
1521 else:
1522 return -1
1524 @computed_field
1525 @cached_property
1526 def getdpCriticalCurrentDensityFunction(self) -> str:
1527 """Return the GetDP function name of the material's critical current density."""
1528 if self.name not in getdpCriticalCurrentDensityFunctions:
1529 raise ValueError(
1530 f"Critical current density of the material '{self.name}' is not defined"
1531 " in FiQuS!"
1532 )
1534 return getdpCriticalCurrentDensityFunctions[self.name]
1536class Pancake3DSolveHTSMaterialBase(BaseModel):
1537 relativeThickness: float = Field(
1538 le=1,
1539 title="Relative Thickness (only for winding)",
1540 description=(
1541 "Winding tapes generally consist of more than one material. Therefore, when"
1542 " materials are given as a list in winding, their relative thickness,"
1543 " (thickness of the material) / (thickness of the bare conductor), should be"
1544 " specified."
1545 ),
1546 )
1549class Pancake3DSolveHTSNormalMaterial(
1550 Pancake3DSolveHTSMaterialBase, Pancake3DSolveNormalMaterial
1551):
1552 pass
1555class Pancake3DSolveHTSSuperconductingMaterial(
1556 Pancake3DSolveHTSMaterialBase, Pancake3DSolveSuperconductingMaterial
1557):
1558 pass
1561class Pancake3DSolveHTSShuntLayerMaterial(Pancake3DSolveNormalMaterial):
1562 name: NormalMaterialName = Field(
1563 default="Copper",
1564 title="Material Name",
1565 )
1566 relativeHeight: float = Field(
1567 default=0.0,
1568 ge=0,
1569 le=1,
1570 title="Relative Height of the Shunt Layer",
1571 description=(
1572 "HTS 2G coated conductor are typically plated, usually "
1573 " using copper. The relative height of the shunt layer is the "
1574 " width of the shunt layer divided by the width of the tape. "
1575 " 0 means no shunt layer."
1576 ),
1577 )
1580class Pancake3DSolveMaterial(BaseModel):
1581 # 1) User inputs:
1583 # Mandatory:
1585 # Optionals:
1586 resistivity: Optional[PositiveFloat] = Field(
1587 default=None,
1588 title="Resistivity",
1589 description=(
1590 "A scalar value. If this is given, material properties won't be used for"
1591 " resistivity."
1592 ),
1593 )
1594 thermalConductivity: Optional[PositiveFloat] = Field(
1595 default=None,
1596 title="Thermal Conductivity",
1597 description=(
1598 "A scalar value. If this is given, material properties won't be used for"
1599 " thermal conductivity."
1600 ),
1601 )
1602 specificHeatCapacity: Optional[PositiveFloat] = Field(
1603 default=None,
1604 title="Specific Heat Capacity",
1605 description=(
1606 "A scalar value. If this is given, material properties won't be used for"
1607 " specific heat capacity."
1608 ),
1609 )
1610 material: Optional[Pancake3DSolveNormalMaterial] = Field(
1611 default=None,
1612 title="Material",
1613 description="Material from STEAM material library.",
1614 )
1616 @model_validator(mode="after")
1617 @classmethod
1618 def check_material(cls, model: "Pancake3DSolveMaterial"):
1619 if model.material is None:
1620 if model.resistivity is None:
1621 raise ValueError(
1622 "Resistivity of the material is not given, and no material is"
1623 " provided!"
1624 )
1625 if model.thermalConductivity is None:
1626 raise ValueError(
1627 "Thermal conductivity of the material is not given, and no material"
1628 " is provided!"
1629 )
1630 if model.specificHeatCapacity is None:
1631 raise ValueError(
1632 "Specific heat capacity of the material is not given, and no"
1633 " material is provided!"
1634 )
1636 return model
1639class Pancake3DSolveShuntLayerMaterial(Pancake3DSolveMaterial):
1640 material: Optional[Pancake3DSolveHTSShuntLayerMaterial] = Field(
1641 default=Pancake3DSolveHTSShuntLayerMaterial(),
1642 title="Material",
1643 description="Material from STEAM material library.",
1644 )
1647class Pancake3DSolveContactLayerMaterial(Pancake3DSolveMaterial):
1648 resistivity: PositiveFloat | Literal["perfectlyInsulating"] = Field(
1649 default=None,
1650 title="Resistivity",
1651 description=(
1652 'A scalar value or "perfectlyInsulating". If "perfectlyInsulating" is'
1653 " given, the contact layer will be perfectly insulating. If this value is"
1654 " given, material properties won't be used for resistivity."
1655 ),
1656 )
1657 numberOfThinShellElements: PositiveInt = Field(
1658 default=1,
1659 title="Number of Thin Shell Elements (Advanced Input)",
1660 description=(
1661 "Number of thin shell elements in the FE formulation (GetDP related, not"
1662 " physical and only used when TSA is set to True)"
1663 ),
1664 )
1666 @field_validator("resistivity")
1667 @classmethod
1668 def checkPerfectlyInsulatingCase(cls, resistivity):
1669 if resistivity == "perfectlyInsulating":
1670 geometry = geometry_input.get()
1671 if geometry["numberOfPancakes"] > 1:
1672 raise ValueError(
1673 "Contact layer can't be perfectly insulating for multi-pancake"
1674 " coils!"
1675 )
1677 return resistivity
1680Pancake3DHTSMaterial = Annotated[
1681 Pancake3DSolveHTSNormalMaterial | Pancake3DSolveHTSSuperconductingMaterial,
1682 Field(discriminator="name"),
1683]
1686class Pancake3DSolveWindingMaterial(Pancake3DSolveMaterial):
1687 material: Optional[list[Pancake3DHTSMaterial]] = Field(
1688 default=None,
1689 title="Materials of HTS CC",
1690 description="List of materials of HTS CC.",
1691 )
1693 shuntLayer: Pancake3DSolveShuntLayerMaterial = Field(
1694 default=Pancake3DSolveShuntLayerMaterial(),
1695 title="Shunt Layer Properties",
1696 description="Material properties of the shunt layer.",
1697 )
1699 isotropic: Optional[bool] = Field(
1700 default=False,
1701 title="Isotropic Material",
1702 description=(
1703 "If True, resistivity and thermal conductivity are isotropic. If False, they are anisotropic. "
1704 "The default is anisotropic material."
1705 ),
1706 )
1708 minimumPossibleResistivity: NonNegativeFloat = Field(
1709 default=1e-20,
1710 title="Minimum Possible Resistivity",
1711 description=(
1712 "The resistivity of the winding won't be lower than this value, no matter"
1713 " what."
1714 ),
1715 )
1717 maximumPossibleResistivity: PositiveFloat = Field(
1718 default=0.01,
1719 title="Maximum Possible Resistivity",
1720 description=(
1721 "The resistivity of the winding won't be higher than this value, no matter"
1722 " what."
1723 ),
1724 )
1726 @field_validator("material")
1727 @classmethod
1728 def checkIfRelativeThicknessesSumToOne(cls, material):
1729 if material is not None:
1730 totalRelativeThickness = sum(
1731 material.relativeThickness for material in material
1732 )
1733 if not math.isclose(totalRelativeThickness, 1, rel_tol=1e-3):
1734 raise ValueError(
1735 "Sum of relative thicknesses of HTS materials should be 1!"
1736 )
1738 return material
1740 @computed_field
1741 @cached_property
1742 def relativeThicknessOfNormalConductor(self) -> PositiveFloat:
1743 """Return the relative thickness of the normal conductor."""
1744 if self.material is None:
1745 return 0
1746 else:
1747 # look at normal materials inside self.material and sum their relativeThickness
1748 return sum(
1749 material.relativeThickness
1750 for material in self.material
1751 if isinstance(material, Pancake3DSolveHTSNormalMaterial)
1752 )
1754 @computed_field
1755 @cached_property
1756 def relativeThicknessOfSuperConductor(self) -> PositiveFloat:
1757 """Return the relative thickness of the super conductor."""
1758 if self.material is None:
1759 return 0
1760 else:
1761 # look at superconducting materials inside self.material and sum their relativeThickness
1762 return sum(
1763 material.relativeThickness
1764 for material in self.material
1765 if isinstance(material, Pancake3DSolveHTSSuperconductingMaterial)
1766 )
1768 @computed_field
1769 @cached_property
1770 def normalConductors(self) -> list[Pancake3DSolveHTSNormalMaterial]:
1771 """Return the normal conductors of the winding."""
1772 if self.material is None:
1773 return []
1774 else:
1775 return [
1776 material
1777 for material in self.material
1778 if isinstance(material, Pancake3DSolveHTSNormalMaterial)
1779 ]
1781 @computed_field
1782 @cached_property
1783 def superConductor(self) -> Pancake3DSolveHTSSuperconductingMaterial:
1784 """Return the super conductor of the winding."""
1785 if self.material is None:
1786 return None
1787 else:
1788 superConductors = [
1789 material
1790 for material in self.material
1791 if isinstance(material, Pancake3DSolveHTSSuperconductingMaterial)
1792 ]
1793 if len(superConductors) == 0:
1794 return None
1795 elif len(superConductors) == 1:
1796 return superConductors[0]
1797 else:
1798 raise ValueError(
1799 "There should be only one superconductor in the winding!"
1800 )
1802class Pancake3DTerminalCryocoolerLumpedMass(Pancake3DSolveMaterial):
1804 material: Optional[Pancake3DSolveNormalMaterial] = Field(
1805 default=Pancake3DSolveNormalMaterial(name="Copper", RRR=295),
1806 title="Material",
1807 description="Material from STEAM material library.",
1808 )
1810 volume: NonNegativeFloat = Field(
1811 default=0,
1812 title="Cryocooler Lumped Block Volume",
1813 description=(
1814 "Volume of the lumped thermal mass between second stage of the cryocooler and pancake coil in m^3. "
1815 "A zero value effectively disables the lumped thermal mass between second stage of the cryocooler and pancake coil."
1816 )
1817 )
1819 numberOfThinShellElements: PositiveInt = Field(
1820 default=1,
1821 title="Number of Thin Shell Elements for Cryocooler Lumped Mass",
1822 description=(
1823 "Number of thin shell elements in the FE formulation (GetDP related, not"
1824 " physical and only used when TSA is set to True)"
1825 ),
1826 )
1828class Pancake3DTerminalCryocoolerBoundaryCondition(BaseModel):
1830 coolingPowerMultiplier: NonNegativeFloat = Field(
1831 default=1,
1832 title="Cooling Power Multiplier",
1833 description=(
1834 "Multiplier for the cooling power. It can be used to scale"
1835 " the cooling power given by the coldhead capacity map by a non-negative float factor."
1836 ),
1837 )
1839 staticHeatLoadPower: NonNegativeFloat = Field(
1840 default=0,
1841 title="Static Heat Load Power",
1842 description=(
1843 "Static heat load power in W. It can be used to add a static heat load"
1844 " to the cryocooler, i.e., decrease the power available for cooling. "
1845 " The actual cooling power is P(t) = P_cryocooler(T) - P_staticLoad."
1846 ),
1847 )
1849 lumpedMass: Pancake3DTerminalCryocoolerLumpedMass = Field(
1850 default = Pancake3DTerminalCryocoolerLumpedMass(),
1851 title="Cryocooler Lumped Mass",
1852 description="Thermal lumped mass between second stage of the cryocooler and pancake coil modeled via TSA.",
1853 )
1855class Pancake3DSolveTerminalMaterialAndBoundaryCondition(Pancake3DSolveMaterial):
1856 cooling: Literal["adiabatic", "fixedTemperature", "cryocooler"] = Field(
1857 default="fixedTemperature",
1858 title="Cooling condition",
1859 description=(
1860 "Cooling condition of the terminal. It can be either adiabatic, fixed"
1861 " temperature, or cryocooler."
1862 ),
1863 )
1865 cryocoolerOptions: Optional[Pancake3DTerminalCryocoolerBoundaryCondition] = Field(
1866 default=Pancake3DTerminalCryocoolerBoundaryCondition(),
1867 title="Cryocooler Boundary Condition",
1868 description="Additional inputs for the cryocooler boundary condition.",
1869 )
1871 transitionNotch: Pancake3DSolveMaterial = Field(
1872 title="Transition Notch Properties",
1873 description="Material properties of the transition notch volume.",
1874 )
1875 terminalContactLayer: Pancake3DSolveMaterial = Field(
1876 title="Transition Layer Properties",
1877 description=(
1878 "Material properties of the transition layer between terminals and"
1879 " windings."
1880 ),
1881 )
1884class Pancake3DSolveToleranceBase(BaseModel):
1885 # Mandatory:
1886 quantity: str
1887 relative: NonNegativeFloat = Field(
1888 title="Relative Tolerance",
1889 description="Relative tolerance for the quantity.",
1890 )
1891 absolute: NonNegativeFloat = Field(
1892 title="Absolute Tolerance", description="Absolute tolerance for the quantity"
1893 )
1895 # Optionals:
1896 normType: Literal["L1Norm", "MeanL1Norm", "L2Norm", "MeanL2Norm", "LinfNorm"] = (
1897 Field(
1898 default="L2Norm",
1899 title="Norm Type",
1900 description=(
1901 "Sometimes, tolerances return a vector instead of a scalar (ex,"
1902 " solutionVector). Then, the magnitude of the tolerance should be"
1903 " calculated with a method. Norm type selects this method."
1904 ),
1905 )
1906 )
1909class Pancake3DSolvePositionRequiredTolerance(Pancake3DSolveToleranceBase):
1910 # Mandatory:
1911 quantity: PositionRequiredQuantityName = Field(
1912 title="Quantity", description="Name of the quantity for tolerance."
1913 )
1914 position: Pancake3DPosition = Field(
1915 title="Probing Position of the Quantity",
1916 description="Probing position of the quantity for tolerance.",
1917 )
1920class Pancake3DSolvePositionNotRequiredTolerance(Pancake3DSolveToleranceBase):
1921 # Mandatory:
1922 quantity: (
1923 Literal[
1924 "electromagneticSolutionVector",
1925 "thermalSolutionVector",
1926 "coupledSolutionVector",
1927 ]
1928 | PositionNotRequiredQuantityName
1929 ) = Field(
1930 title="Quantity",
1931 description="Name of the quantity for tolerance.",
1932 )
1935Pancake3DSolveTolerance = Annotated[
1936 Pancake3DSolvePositionRequiredTolerance
1937 | Pancake3DSolvePositionNotRequiredTolerance,
1938 Field(discriminator="quantity"),
1939]
1942class Pancake3DSolveSettingsWithTolerances(BaseModel):
1943 tolerances: list[Pancake3DSolveTolerance] = Field(
1944 title="Tolerances for Adaptive Time Stepping",
1945 description=(
1946 "Time steps or nonlinear iterations will be refined until the tolerances"
1947 " are satisfied."
1948 ),
1949 )
1951 @computed_field
1952 @cached_property
1953 def postOperationTolerances(self) -> list[Pancake3DSolveTolerance]:
1954 """Return the post operation tolerances."""
1955 tolerances = [
1956 tolerance
1957 for tolerance in self.tolerances
1958 if "SolutionVector" not in tolerance.quantity
1959 ]
1960 return tolerances
1963 @computed_field
1964 @cached_property
1965 def systemTolerances(self) -> list[Pancake3DSolveTolerance]:
1966 """Return the system tolerances."""
1967 tolerances = [
1968 tolerance
1969 for tolerance in self.tolerances
1970 if "SolutionVector" in tolerance.quantity
1971 ]
1972 return tolerances
1975class Pancake3DSolveAdaptiveTimeLoopSettings(Pancake3DSolveSettingsWithTolerances):
1976 # Mandatory:
1977 initialStep: PositiveFloat = Field(
1978 title="Initial Step for Adaptive Time Stepping",
1979 description="Initial step for adaptive time stepping",
1980 )
1981 minimumStep: PositiveFloat = Field(
1982 title="Minimum Step for Adaptive Time Stepping",
1983 description=(
1984 "The simulation will be aborted if a finer time step is required than this"
1985 " minimum step value."
1986 ),
1987 )
1988 maximumStep: PositiveFloat = Field(
1989 description="Bigger steps than this won't be allowed",
1990 )
1992 # Optionals:
1993 integrationMethod: Literal[
1994 "Euler", "Gear_2", "Gear_3", "Gear_4", "Gear_5", "Gear_6"
1995 ] = Field(
1996 default="Euler",
1997 title="Integration Method",
1998 description="Integration method for transient analysis",
1999 )
2000 breakPoints_input: list[float] = Field(
2001 default=[0],
2002 title="Break Points for Adaptive Time Stepping",
2003 description="Make sure to solve the system for these times.",
2004 )
2006 @field_validator("breakPoints_input")
2007 @classmethod
2008 def updateGlobalBreakPointsList(cls, breakPoints_input):
2009 all_break_points.extend(breakPoints_input)
2010 return breakPoints_input
2012 @model_validator(mode="after")
2013 @classmethod
2014 def check_time_steps(cls, model: "Pancake3DSolveAdaptiveTimeLoopSettings"):
2015 """
2016 Checks if the time steps are consistent.
2018 :param values: dictionary of time steps
2019 :type values: dict
2020 :return: dictionary of time steps
2021 :rtype: dict
2022 """
2024 if model.initialStep < model.minimumStep:
2025 raise ValueError(
2026 "Initial time step cannot be smaller than the minimum time step!"
2027 )
2029 if model.initialStep > model.maximumStep:
2030 raise ValueError(
2031 "Initial time step cannot be greater than the maximum time step!"
2032 )
2034 if model.minimumStep > model.maximumStep:
2035 raise ValueError(
2036 "Minimum time step cannot be greater than the maximum time step!"
2037 )
2039 return model
2041 @computed_field
2042 @cached_property
2043 def breakPoints(self) -> list[float]:
2044 """Return the break points for adaptive time stepping."""
2045 breakPoints = list(set(all_break_points))
2046 breakPoints.sort()
2047 return breakPoints
2050class Pancake3DSolveFixedTimeLoopSettings(BaseModel):
2051 # Mandatory:
2052 step: PositiveFloat = Field(
2053 title="Step for Fixed Time Stepping",
2054 description="Time step for fixed time stepping.",
2055 )
2058class Pancake3DSolveFixedLoopInterval(BaseModel):
2059 # Mandatory:
2060 startTime: NonNegativeFloat = Field(
2061 title="Start Time of the Interval",
2062 description="Start time of the interval.",
2063 )
2064 endTime: NonNegativeFloat = Field(
2065 title="End Time of the Interval",
2066 description="End time of the interval.",
2067 )
2068 step: PositiveFloat = Field(
2069 title="Step for the Interval",
2070 description="Time step for the interval",
2071 )
2073 @model_validator(mode="after")
2074 @classmethod
2075 def check_time_steps(cls, model: "Pancake3DSolveFixedLoopInterval"):
2076 """
2077 Checks if the time steps are consistent.
2079 :param values: dictionary of time steps
2080 :type values: dict
2081 :return: dictionary of time steps
2082 :rtype: dict
2083 """
2085 if model.startTime > model.endTime:
2086 raise ValueError("Start time cannot be greater than the end time!")
2088 interval = model.endTime - model.startTime
2089 if (
2090 math.fmod(interval, model.step) > 1e-8
2091 and math.fmod(interval, model.step) - model.step < -1e-8
2092 ):
2093 raise ValueError("Time interval must be a multiple of the time step!")
2095 return model
2098class Pancake3DSolveTimeBase(BaseModel):
2099 # Mandatory:
2100 start: float = Field(
2101 title="Start Time", description="Start time of the simulation."
2102 )
2103 end: float = Field(title="End Time", description="End time of the simulation.")
2105 # Optionals:
2106 extrapolationOrder: Literal[0, 1, 2, 3] = Field(
2107 default=1,
2108 title="Extrapolation Order",
2109 description=(
2110 "Before solving for the next time steps, the previous solutions can be"
2111 " extrapolated for better convergence."
2112 ),
2113 )
2115 @model_validator(mode="after")
2116 @classmethod
2117 def check_time_steps(cls, model: "Pancake3DSolveTimeBase"):
2118 """
2119 Checks if the time steps are consistent.
2120 """
2122 if model.start > model.end:
2123 raise ValueError("Start time cannot be greater than the end time!")
2125 return model
2128class Pancake3DSolveTimeAdaptive(Pancake3DSolveTimeBase):
2129 timeSteppingType: Literal["adaptive"] = "adaptive"
2130 adaptiveSteppingSettings: Pancake3DSolveAdaptiveTimeLoopSettings = Field(
2131 title="Adaptive Time Loop Settings",
2132 description=(
2133 "Adaptive time loop settings (only used if stepping type is adaptive)."
2134 ),
2135 )
2137 @model_validator(mode="after")
2138 @classmethod
2139 def check_time_steps(cls, model: "Pancake3DSolveTimeAdaptive"):
2140 """
2141 Checks if the time steps are consistent.
2142 """
2143 if model.adaptiveSteppingSettings.initialStep > model.end:
2144 raise ValueError("Initial time step cannot be greater than the end time!")
2146 return model
2149class Pancake3DSolveTimeFixed(Pancake3DSolveTimeBase):
2150 timeSteppingType: Literal["fixed"] = "fixed"
2151 fixedSteppingSettings: (
2152 list[Pancake3DSolveFixedLoopInterval] | Pancake3DSolveFixedTimeLoopSettings
2153 ) = Field(
2154 title="Fixed Time Loop Settings",
2155 description="Fixed time loop settings (only used if stepping type is fixed).",
2156 )
2158 @model_validator(mode="after")
2159 @classmethod
2160 def check_time_steps(cls, model: "Pancake3DSolveTimeFixed"):
2161 if isinstance(model.fixedSteppingSettings, list):
2162 for i in range(len(model.fixedSteppingSettings) - 1):
2163 if model.fixedSteppingSettings[i].endTime != model.fixedSteppingSettings[i + 1].startTime:
2164 raise ValueError(
2165 "End time of the previous interval must be equal to the start"
2166 " time of the next interval!"
2167 )
2169 if model.fixedSteppingSettings[0].startTime != model.start:
2170 raise ValueError(
2171 "Start time of the first interval must be equal to the start time"
2172 " of the simulation!"
2173 )
2175 else:
2176 if model.fixedSteppingSettings.step > model.end:
2177 raise ValueError("Time step cannot be greater than the end time!")
2179 if not (
2180 math.isclose(
2181 (model.end - model.start) % model.fixedSteppingSettings.step, 0, abs_tol=1e-8
2182 )
2183 ):
2184 raise ValueError("Time interval must be a multiple of the time step!")
2187Pancake3DSolveTime = Annotated[
2188 Pancake3DSolveTimeAdaptive | Pancake3DSolveTimeFixed,
2189 Field(discriminator="timeSteppingType"),
2190]
2193class Pancake3DSolveNonlinearSolverSettings(Pancake3DSolveSettingsWithTolerances):
2194 # Optionals:
2195 maximumNumberOfIterations: PositiveInt = Field(
2196 default=100,
2197 title="Maximum Number of Iterations",
2198 description="Maximum number of iterations allowed for the nonlinear solver.",
2199 )
2200 relaxationFactor: float = Field(
2201 default=0.7,
2202 gt=0,
2203 title="Relaxation Factor",
2204 description=(
2205 "Calculated step changes of the solution vector will be multiplied with"
2206 " this value for better convergence."
2207 ),
2208 )
2211class Pancake3DSolveInitialConditions(BaseModel):
2212 # 1) User inputs:
2214 # Mandatory:
2215 temperature: PositiveFloat = Field(
2216 title="Initial Temperature",
2217 description="Initial temperature of the pancake coils.",
2218 )
2220class Pancake3DSolveImposedField(BaseModel):
2222 imposedAxialField: float = Field(
2223 title="Imposed Axial Magnetic Field",
2224 description="Imposed axial magnetic field in Tesla. Only constant, purely axial magnetic fields are supported at the moment.",
2225)
2227class Pancake3DSolveLocalDefect(BaseModel):
2228 # Mandatory:
2229 value: NonNegativeFloat = Field(
2230 title="Value",
2231 description="Value of the local defect.",
2232 )
2233 startTurn: NonNegativeFloat = Field(
2234 title="Start Turn",
2235 description="Start turn of the local defect.",
2236 )
2237 endTurn: PositiveFloat = Field(
2238 title="End Turn",
2239 description="End turn of the local defect.",
2240 )
2242 startTime: NonNegativeFloat = Field(
2243 title="Start Time",
2244 description="Start time of the local defect.",
2245 )
2247 # Optionals:
2248 transitionDuration: NonNegativeFloat = Field(
2249 default=0,
2250 title="Transition Duration",
2251 description=(
2252 "Transition duration of the local defect. If not given, the transition will"
2253 " be instantly."
2254 ),
2255 )
2256 whichPancakeCoil: Optional[PositiveInt] = Field(
2257 default=None,
2258 title="Pancake Coil Number",
2259 description="The first pancake coil is 1, the second is 2, etc.",
2260 )
2262 @field_validator("startTime")
2263 @classmethod
2264 def updateGlobalBreakPointsList(cls, startTime):
2265 all_break_points.append(startTime)
2266 return startTime
2268 @field_validator("startTime")
2269 @classmethod
2270 def check_start_time(cls, startTime):
2271 solve = solve_input.get()
2272 start_time = solve["time"]["start"]
2273 end_time = solve["time"]["end"]
2275 if startTime < start_time or startTime > end_time:
2276 raise ValueError(
2277 "Start time of the local defect should be between the start and end"
2278 " times!"
2279 )
2281 return start_time
2283 @field_validator("endTurn")
2284 @classmethod
2285 def check_end_turn(cls, endTurn, info: ValidationInfo):
2286 geometry = geometry_input.get()
2288 if endTurn > geometry["winding"]["numberOfTurns"]:
2289 raise ValueError(
2290 "End turn of the local defect can't be greater than the number of"
2291 " turns!"
2292 )
2294 if endTurn < info.data["startTurn"]:
2295 raise ValueError(
2296 "End turn of the local defect can't be smaller than the start turn!"
2297 )
2299 return endTurn
2301 @field_validator("whichPancakeCoil")
2302 @classmethod
2303 def check_which_pancake_coil(cls, whichPancakeCoil):
2304 geometry = geometry_input.get()
2306 if whichPancakeCoil is None:
2307 if geometry["numberOfPancakes"] == 1:
2308 whichPancakeCoil = 1
2309 else:
2310 raise ValueError(
2311 "whickPancakeCoil (pancake coil number) should be given if there"
2312 " are more than one pancake coil!"
2313 )
2315 return whichPancakeCoil
2317 @computed_field
2318 @cached_property
2319 def zTop(self) -> float:
2320 """Return the z-coordinate of the top of the pancake coil."""
2321 geometry = geometry_input.get()
2323 zTop = self.zBottom + geometry["winding"]["height"]
2325 return zTop
2327 @computed_field
2328 @cached_property
2329 def zBottom(self) -> float:
2330 """Return the z-coordinate of the bottom of the pancake coil."""
2331 geometry = geometry_input.get()
2333 zBottom = -(
2334 geometry["numberOfPancakes"] * geometry["winding"]["height"]
2335 + (geometry["numberOfPancakes"] - 1)
2336 * geometry["gapBetweenPancakes"]
2337 ) / 2 + (self.whichPancakeCoil - 1) * (
2338 geometry["winding"]["height"] + geometry["gapBetweenPancakes"]
2339 )
2341 return zBottom
2344class Pancake3DSolveLocalDefects(BaseModel):
2345 # 1) User inputs:
2347 criticalCurrentDensity: Optional[Pancake3DSolveLocalDefect] = Field(
2348 default=None,
2349 title="Local Defect for Critical Current Density",
2350 description="Set critical current density locally.",
2351 )
2353 @field_validator("criticalCurrentDensity")
2354 @classmethod
2355 def check_criticalCurrentDensity_local_defect(cls, criticalCurrentDensity):
2356 if criticalCurrentDensity is not None:
2357 solve = solve_input.get()
2359 if "material" in solve["winding"]:
2360 windingMaterials = [
2361 material["name"] for material in solve["winding"]["material"]
2362 ]
2363 else:
2364 windingMaterials = []
2366 superconducting_material_is_used = "HTSSuperPower" in windingMaterials or "HTSFujikura" in windingMaterials or "HTSSucci" in windingMaterials
2368 if not superconducting_material_is_used:
2369 raise ValueError(
2370 "Local defects can only be set if superconducting material is used!"
2371 )
2373 return criticalCurrentDensity
2376class Pancake3DSolveConvectiveCooling(BaseModel):
2378 heatTransferCoefficient: Optional[Union[NonNegativeFloat, Literal["nitrogenBath"]]] = Field(
2379 default=0,
2380 title="Heat Transfer Coefficient",
2381 description=(
2382 "The heat transfer coefficient for the heat transfer between the winding and the air. "
2383 "If zero, no heat transfer to the air is considered."
2384 "This feature is only implemented for the thin shell approximation."
2385 "At the moment, only constant values are supported."
2386 ),
2387 )
2389 exteriorBathTemperature: Optional[NonNegativeFloat] = Field(
2390 default=4.2,
2391 title="Exterior Bath Temperature",
2392 description=(
2393 "The temperature of the exterior bath for convective cooling boundary condition. "
2394 ),
2395 )
2397class Pancake3DSolveEECircuit(BaseModel):
2399 inductanceInSeriesWithPancakeCoil: Optional[NonNegativeFloat] = Field(
2400 default=0,
2401 title="Inductance in Series with Pancake Coil",
2402 description=(
2403 "A lumped inductance in series with the pancake coil to model a bigger coil. "
2404 ),
2405 )
2407 enable: bool = Field(
2408 default=False,
2409 title="Enable Detection Circuit",
2410 description=(
2411 "Enable the detection circuit for the pancake coil. "
2412 ),
2413 )
2415 ResistanceEnergyExtractionOpenSwitch: Optional[NonNegativeFloat] = Field(
2416 default=1E6,
2417 title="Resistance of Energy Extraction Open Switch",
2418 description=(
2419 "The resistance of the energy extraction switch when modeled as open. "
2420 ),
2421 )
2423 ResistanceEnergyExtractionClosedSwitch: Optional[NonNegativeFloat] = Field(
2424 default=1E-6,
2425 title="Resistance of Energy Extraction Closed Switch",
2426 description=(
2427 "The resistance of the energy extraction switch when modeled as closed. "
2428 ),
2429 )
2431 ResistanceCrowbarOpenSwitch: Optional[NonNegativeFloat] = Field(
2432 default=1E6,
2433 title="Resistance of Crowbar Open Switch",
2434 description=(
2435 "The resistance of the crowbar switch when modeled as open. "
2436 ),
2437 )
2439 ResistanceCrowbarClosedSwitch: Optional[NonNegativeFloat] = Field(
2440 default=1E-6,
2441 title="Resistance of Crowbar Closed Switch",
2442 description=(
2443 "The resistance of the crowbar switch when modeled as closed. "
2444 ),
2445 )
2447 stopSimulationAtCurrent: Optional[NonNegativeFloat] = Field(
2448 default=0,
2449 title="Stop Simulation at Current",
2450 description=(
2451 "If a quench is detected and the current reaches this value, the simulation will be stopped after. "
2452 "stopSimulationWaitingTime seconds."
2453 ),
2454 )
2456 stopSimulationWaitingTime: Optional[NonNegativeFloat] = Field(
2457 default=0,
2458 title="Stop Simulation Waiting Time",
2459 description=(
2460 "The time to wait after a quench is detected and the current reaches stopSimulationAtCurrent before stopping the simulation."
2461 ),
2462 )
2464 # or use t_off from power supply? I don't like it since it's used as time value to turn off and not time interval.
2465 TurnOffDeltaTimePowerSupply: Optional[NonNegativeFloat] = Field(
2466 default=0,
2467 title="Time Interval Between Quench Detection and Power Supply Turn Off",
2468 description=(
2469 "The time it takes for the power supply to be turned off after quench detection. "
2470 "A linear ramp-down is assumed between the time of quench detection and the time of power supply turn off."
2471 ),
2472 )
2474class Pancake3DSolvePowerDensity(BaseModel):
2475 power: Optional[NonNegativeFloat] = Field(
2476 default=0,
2477 title="Power Density",
2478 description=(
2479 "The power in W for an imposed power density in the winding. "
2480 "'startTime', 'endTime', 'startTurn', and 'endTurn' "
2481 "are also required to be set."
2482 ),
2483 )
2485 startTime: Optional[NonNegativeFloat] = Field(
2486 default=0,
2487 title="Power Density Start Time",
2488 description=(
2489 "The start time for the imposed power density in the winding. "
2490 "'power', 'endTime', 'startTurn', and 'endTurn' "
2491 "are also required to be set."
2492 ),
2493 )
2495 endTime: Optional[NonNegativeFloat] = Field(
2496 default=0,
2497 title="Power Density End Time",
2498 description=(
2499 "The end time for the imposed power density in the winding. "
2500 "'power', 'startTime', 'startTurn', and 'endTurn' "
2501 "are also required to be set."
2502 ),
2503 )
2505 startArcLength: Optional[NonNegativeFloat] = Field(
2506 default=0,
2507 title="Power Density Start Arc Length",
2508 description=(
2509 "The start arc length in m for the imposed power density in the winding. "
2510 "'power', 'startTime', 'endTime', and 'endArcLength' "
2511 "are also required to be set."
2512 ),
2513 )
2515 endArcLength: Optional[NonNegativeFloat] = Field(
2516 default=0,
2517 title="Power Density End Arc Length",
2518 description=(
2519 "The end arc length in m for the imposed power density in the winding. "
2520 "'power', 'startTime', 'endTime', and 'startArcLength' "
2521 "are also required to be set."
2522 ),
2523 )
2525class Pancake3DSolveQuantityBase(BaseModel):
2526 # Mandatory:
2527 quantity: PositionNotRequiredQuantityName | PositionRequiredQuantityName = Field(
2528 title="Quantity",
2529 description="Name of the quantity to be saved.",
2530 )
2532 @computed_field
2533 @cached_property
2534 def getdpQuantityName(self) -> str:
2535 """Return the GetDP name of the quantity."""
2536 if self.quantity not in getdpQuantityNames:
2537 raise ValueError(f"Quantity '{self.quantity}' is not defined in FiQuS!")
2539 return getdpQuantityNames[self.quantity]
2541 @computed_field
2542 @cached_property
2543 def getdpPostOperationName(self) -> str:
2544 """Return the GetDP name of the post operation."""
2545 if self.quantity not in getdpPostOperationNames:
2546 raise ValueError(f"Quantity '{self.quantity}' is not defined in FiQuS!")
2548 return getdpPostOperationNames[self.quantity]
2551class Pancake3DSolveSaveQuantity(Pancake3DSolveQuantityBase):
2552 # Optionals:
2553 timesToBeSaved: Union[list[float], None] = Field(
2554 default=None,
2555 title="Times to be Saved",
2556 description=(
2557 "List of times that wanted to be saved. If not given, all the time steps"
2558 " will be saved."
2559 ),
2560 )
2562 @field_validator("timesToBeSaved")
2563 @classmethod
2564 def updateGlobalBreakPointsList(cls, timesToBeSaved):
2565 if timesToBeSaved is not None:
2566 all_break_points.extend(timesToBeSaved)
2567 return timesToBeSaved
2569 @field_validator("timesToBeSaved")
2570 @classmethod
2571 def check_times_to_be_saved(cls, timesToBeSaved):
2572 solve = solve_input.get()
2573 start_time = solve["time"]["start"]
2574 end_time = solve["time"]["end"]
2576 if timesToBeSaved is not None:
2577 for time in timesToBeSaved:
2578 if time < start_time or time > end_time:
2579 raise ValueError(
2580 "Times to be saved should be between the start and end times!"
2581 )
2584# ======================================================================================
2585# SOLVE CLASSES ENDS ===================================================================
2586# ======================================================================================
2588# ======================================================================================
2589# POSTPROCESS CLASSES STARTS ===========================================================
2590# ======================================================================================
2593class Pancake3DPostprocessTimeSeriesPlotBase(Pancake3DSolveQuantityBase):
2594 # Mandatory:
2595 quantity: str
2597 @computed_field
2598 @cached_property
2599 def fileName(self) -> str:
2600 """Return the name of the file to be saved."""
2601 return f"{self.quantity}(TimeSeriesPlot)"
2603 @computed_field
2604 @cached_property
2605 def quantityProperName(self) -> str:
2606 """Return the proper name of the quantity."""
2607 if self.quantity not in quantityProperNames:
2608 raise ValueError(
2609 f"Quantity '{self.quantity}'s proper name is not defined! Please"
2610 ' update the dictionary "quantityProperNames" in the file'
2611 ' "fiqus/data/DataFiQuSPancake3D.py".'
2612 )
2614 return quantityProperNames[self.quantity]
2616 @computed_field
2617 @cached_property
2618 def units(self) -> str:
2619 """Return the units of the quantity."""
2620 if self.quantity not in quantityUnits:
2621 raise ValueError(
2622 f"Quantity '{self.quantity}'s units are not defined! Please update"
2623 ' the dictionary "quantityUnits" in the file'
2624 ' "fiqus/data/DataFiQuSPancake3D.py".'
2625 )
2627 return quantityUnits[self.quantity]
2630class Pancake3DPostprocessTimeSeriesPlotPositionRequired(
2631 Pancake3DPostprocessTimeSeriesPlotBase
2632):
2633 # Mandatory:
2634 quantity: PositionRequiredQuantityName = Field(
2635 title="Quantity",
2636 description="Name of the quantity to be plotted.",
2637 )
2639 position: Pancake3DPosition = Field(
2640 title="Probing Position",
2641 description="Probing position of the quantity for time series plot.",
2642 )
2645class Pancake3DPostprocessTimeSeriesPlotPositionNotRequired(
2646 Pancake3DPostprocessTimeSeriesPlotBase
2647):
2648 # Mandatory:
2649 quantity: PositionNotRequiredQuantityName = Field(
2650 title="Quantity",
2651 description="Name of the quantity to be plotted.",
2652 )
2655Pancake3DPostprocessTimeSeriesPlot = Optional[Annotated[
2656 Pancake3DPostprocessTimeSeriesPlotPositionRequired
2657 | Pancake3DPostprocessTimeSeriesPlotPositionNotRequired,
2658 Field(discriminator="quantity"),
2659]]
2662class Pancake3DPostprocessMagneticFieldOnPlane(BaseModel):
2663 # Optional:
2664 colormap: str = Field(
2665 default="viridis",
2666 title="Colormap",
2667 description="Colormap for the plot.",
2668 )
2669 streamLines: bool = Field(
2670 default=True,
2671 title="Stream Lines",
2672 description=(
2673 "If True, streamlines will be plotted. Note that magnetic field vectors may"
2674 " have components perpendicular to the plane, and streamlines will be drawn"
2675 " depending on the vectors' projection onto the plane."
2676 ),
2677 )
2678 interpolationMethod: Literal["nearest", "linear", "cubic"] = Field(
2679 default="linear",
2680 title="Interpolation Method",
2681 description=(
2682 "Interpolation type for the plot.Because of the FEM basis function"
2683 " selections of FiQuS, each mesh element has a constant magnetic field"
2684 " vector. Therefore, for smooth 2D plots, interpolation can be"
2685 " used.Types:nearest: it will plot the nearest magnetic field value to"
2686 " the plotting point.linear: it will do linear interpolation to the"
2687 " magnetic field values.cubic: it will do cubic interpolation to the"
2688 " magnetic field values."
2689 ),
2690 )
2691 timesToBePlotted: Optional[list[float]] = Field(
2692 default=None,
2693 title="Times to be Plotted",
2694 description=(
2695 "List of times that wanted to be plotted. If not given, all the time steps"
2696 " will be plotted."
2697 ),
2698 )
2699 planeNormal: Annotated[list[float], Len(min_length=3, max_length=3)] = Field(
2700 default=[1, 0, 0],
2701 title="Plane Normal",
2702 description="Normal vector of the plane. The default is YZ-plane (1, 0, 0).",
2703 )
2704 planeXAxisUnitVector: Annotated[list[float], Len(min_length=3, max_length=3)] = (
2705 Field(
2706 default=[0, 1, 0],
2707 title="Plane X Axis",
2708 description=(
2709 "If an arbitrary plane is wanted to be plotted, the arbitrary plane's X"
2710 " axis unit vector must be specified. The dot product of the plane's X"
2711 " axis and the plane's normal vector must be zero."
2712 ),
2713 )
2714 )
2716 @field_validator("timesToBePlotted")
2717 @classmethod
2718 def updateGlobalBreakPointsList(cls, timesToBePlotted):
2719 all_break_points.extend(timesToBePlotted)
2720 return timesToBePlotted
2722 @field_validator("colormap")
2723 @classmethod
2724 def check_color_map(cls, colorMap):
2725 """
2726 Check if the colormap exists.
2727 """
2728 if colorMap not in matplotlib.pyplot.colormaps():
2729 raise ValueError(
2730 f"{colorMap} is not a valid colormap! Please see"
2731 " https://matplotlib.org/stable/gallery/color/colormap_reference.html"
2732 " for available colormaps."
2733 )
2735 return colorMap
2737 @computed_field
2738 @cached_property
2739 def onSection(self) -> list[list[float]]:
2740 """Return the three corner points of the plane."""
2742 class unitVector:
2743 def __init__(self, u, v, w) -> None:
2744 length = math.sqrt(u**2 + v**2 + w**2)
2745 self.u = u / length
2746 self.v = v / length
2747 self.w = w / length
2749 def rotate(self, theta, withRespectTo):
2750 # Rotate with respect to the withRespectTo vector by theta degrees:
2751 # https://en.wikipedia.org/wiki/Rotation_matrix#Rotation_matrix_from_axis_and_angle
2752 a = withRespectTo.u
2753 b = withRespectTo.v
2754 c = withRespectTo.w
2756 rotationMatrix = np.array(
2757 [
2758 [
2759 math.cos(theta) + a**2 * (1 - math.cos(theta)),
2760 a * b * (1 - math.cos(theta)) - c * math.sin(theta),
2761 a * c * (1 - math.cos(theta)) + b * math.sin(theta),
2762 ],
2763 [
2764 b * a * (1 - math.cos(theta)) + c * math.sin(theta),
2765 math.cos(theta) + b**2 * (1 - math.cos(theta)),
2766 b * c * (1 - math.cos(theta)) - a * math.sin(theta),
2767 ],
2768 [
2769 c * a * (1 - math.cos(theta)) - b * math.sin(theta),
2770 c * b * (1 - math.cos(theta)) + a * math.sin(theta),
2771 math.cos(theta) + c**2 * (1 - math.cos(theta)),
2772 ],
2773 ]
2774 )
2775 vector = np.array([[self.u], [self.v], [self.w]])
2776 rotatedVector = rotationMatrix @ vector
2777 return unitVector(
2778 rotatedVector[0][0],
2779 rotatedVector[1][0],
2780 rotatedVector[2][0],
2781 )
2783 def __pow__(self, otherUnitVector):
2784 # Cross product:
2785 u = self.v * otherUnitVector.w - self.w * otherUnitVector.v
2786 v = self.w * otherUnitVector.u - self.u * otherUnitVector.w
2787 w = self.u * otherUnitVector.v - self.v * otherUnitVector.u
2788 return unitVector(u, v, w)
2790 def __mul__(self, otherUnitVector) -> float:
2791 # Dot product:
2792 return (
2793 self.u * otherUnitVector.u
2794 + self.v * otherUnitVector.v
2795 + self.w * otherUnitVector.w
2796 )
2798 planeNormal = unitVector(
2799 self.planeNormal[0], self.planeNormal[1], self.planeNormal[2]
2800 )
2801 planeXAxis = unitVector(
2802 self.planeXAxis[0], self.planeXAxis[1], self.planeXAxis[2]
2803 )
2805 if not math.isclose(planeNormal * planeXAxis, 0, abs_tol=1e-8):
2806 raise ValueError(
2807 "planeNormal and planeXAxis must be perpendicular to each"
2808 " other! If planeNormal is chosen arbitrarily, planeXAxis must"
2809 " be specified."
2810 )
2812 # A plane that passes through the origin with the normal vector planeNormal
2813 # can be defined as:
2814 # a*x + b*y + c*z = 0
2815 a = planeNormal.u
2816 b = planeNormal.v
2817 c = planeNormal.w
2819 # Pick three points on the plane to define a rectangle. The points will
2820 # be the corners of the rectangle.
2822 # To do that, change coordinate system:
2824 # Find a vector that is perpendicular to planeNormal:
2825 randomVector = unitVector(c, a, b)
2826 perpendicularVector1 = planeNormal**randomVector
2828 # Rotate perperndicular vector with respect to the plane's normal vector
2829 # by 90 degrees to find the second perpendicular vector:
2830 perpendicularVector2 = perpendicularVector1.rotate(math.pi / 2, planeNormal)
2832 # Build the transformation matrix to change from the plane's coordinate
2833 # system to the global coordinate system:
2834 transformationMatrix = np.array(
2835 [
2836 [
2837 perpendicularVector1.u,
2838 perpendicularVector1.v,
2839 perpendicularVector1.w,
2840 ],
2841 [
2842 perpendicularVector2.u,
2843 perpendicularVector2.v,
2844 perpendicularVector2.w,
2845 ],
2846 [planeNormal.u, planeNormal.v, planeNormal.w],
2847 ]
2848 )
2849 transformationMatrix = np.linalg.inv(transformationMatrix)
2851 # Select three points to define the rectangle. Pick large points because
2852 # only the intersection of the rectangle and the mesh will be used.
2853 geometry = geometry_input.get()
2854 if geometry["air"]["type"] == "cuboid":
2855 sideLength = geometry["air"]["sideLength"]
2856 airMaxWidth = 2 * math.sqrt((sideLength / 2) ** 2 + (sideLength / 2) ** 2)
2857 if geometry["air"]["type"] == "cylinder":
2858 airMaxWidth = geometry["air"]["radius"]
2860 airHeight = getAirHeight()
2862 point1InPlaneCoordinates = np.array(
2863 [5 * max(airMaxWidth, airHeight), 5 * max(airMaxWidth, airHeight), 0]
2864 )
2865 point1 = transformationMatrix @ point1InPlaneCoordinates
2867 point2InPlaneCoordinates = np.array(
2868 [-5 * max(airMaxWidth, airHeight), 5 * max(airMaxWidth, airHeight), 0]
2869 )
2870 point2 = transformationMatrix @ point2InPlaneCoordinates
2872 point3InPlaneCoordinates = np.array(
2873 [
2874 -5 * max(airMaxWidth, airHeight),
2875 -5 * max(airMaxWidth, airHeight),
2876 0,
2877 ]
2878 )
2879 point3 = transformationMatrix @ point3InPlaneCoordinates
2881 # Round the point coordinates to the nearest multiple of the dimTol:
2882 if "dimTol" in geometry:
2883 dimTol = geometry["dimTol"]
2884 else:
2885 dimTol = 1e-8
2887 point1[0] = round(point1[0] / dimTol) * dimTol
2888 point1[1] = round(point1[1] / dimTol) * dimTol
2889 point1[2] = round(point1[2] / dimTol) * dimTol
2890 point2[0] = round(point2[0] / dimTol) * dimTol
2891 point2[1] = round(point2[1] / dimTol) * dimTol
2892 point2[2] = round(point2[2] / dimTol) * dimTol
2893 point3[0] = round(point3[0] / dimTol) * dimTol
2894 point3[1] = round(point3[1] / dimTol) * dimTol
2895 point3[2] = round(point3[2] / dimTol) * dimTol
2897 onSection = [
2898 [float(point1[0]), float(point1[1]), float(point1[2])],
2899 [float(point2[0]), float(point2[1]), float(point2[2])],
2900 [float(point3[0]), float(point3[1]), float(point3[2])],
2901 ]
2903 return onSection
2906# ======================================================================================
2907# POSTPROCESS CLASSES ENDS =============================================================
2908# ======================================================================================
2911class Pancake3DGeometry(BaseModel):
2912 conductorWrite: bool = Field(
2913 default=False,
2914 title="Flag:to Write the Conductor File",
2915 description="To Write the Conductor File"
2916 )
2918 # Mandatory:
2919 numberOfPancakes: PositiveInt = Field(
2920 ge=1,
2921 title="Number of Pancakes",
2922 description="Number of pancake coils stacked on top of each other.",
2923 )
2925 gapBetweenPancakes: PositiveFloat = Field(
2926 title="Gap Between Pancakes",
2927 description="Gap distance between the pancake coils.",
2928 )
2930 winding: Pancake3DGeometryWinding = Field(
2931 title="Winding Geometry",
2932 description="This dictionary contains the winding geometry information.",
2933 )
2935 contactLayer: Pancake3DGeometryContactLayer = Field(
2936 title="Contact Layer Geometry",
2937 description="This dictionary contains the contact layer geometry information.",
2938 )
2940 terminals: Pancake3DGeometryTerminals = Field(
2941 title="Terminals Geometry",
2942 description="This dictionary contains the terminals geometry information.",
2943 )
2945 air: Pancake3DGeometryAir = Field(
2946 title="Air Geometry",
2947 description="This dictionary contains the air geometry information.",
2948 )
2950 # Optionals:
2951 dimensionTolerance: PositiveFloat = Field(
2952 default=1e-8,
2953 description="dimension tolerance (CAD related, not physical)",
2954 )
2955 pancakeBoundaryName: str = Field(
2956 default="PancakeBoundary",
2957 description=(
2958 "name of the pancake's curves that touches the air to be used in the mesh"
2959 ),
2960 )
2961 contactLayerBoundaryName: str = Field(
2962 default="contactLayerBoundary",
2963 description=(
2964 "name of the contact layers's curves that touches the air to be used in the"
2965 " mesh (only for TSA)"
2966 ),
2967 )
2970class Pancake3DMesh(BaseModel):
2971 # Mandatory:
2972 winding: Pancake3DMeshWinding = Field(
2973 title="Winding Mesh",
2974 description="This dictionary contains the winding mesh information.",
2975 )
2976 contactLayer: Pancake3DMeshContactLayer = Field(
2977 title="Contact Layer Mesh",
2978 description="This dictionary contains the contact layer mesh information.",
2979 )
2981 # Optionals:
2982 terminals: Pancake3DMeshAirAndTerminals = Field(
2983 default=Pancake3DMeshAirAndTerminals(),
2984 title="Terminal Mesh",
2985 description="This dictionary contains the terminal mesh information.",
2986 )
2987 air: Pancake3DMeshAirAndTerminals = Field(
2988 default=Pancake3DMeshAirAndTerminals(),
2989 title="Air Mesh",
2990 description="This dictionary contains the air mesh information.",
2991 )
2993 computeCohomologyForInsulating: bool = Field(
2994 default=True,
2995 title="Compute Cohomology for Insulating",
2996 description=(
2997 "Expert option only. "
2998 "If False, the cohomology regions needed for simulating an insulating coil"
2999 "will not be computed. This will reduce the time spent for the meshing "
3000 "or more accurately the cohomology computing phase. BEWARE: The simulation "
3001 "will fail if set to False and a perfectlyInsulating coil is simulated."
3002 ),
3003 )
3005 # Mandatory:
3006 minimumElementSize: PositiveFloat = Field(
3007 title="Minimum Element Size",
3008 description=(
3009 "The minimum mesh element size in terms of the largest mesh size in the"
3010 " winding. This mesh size will be used in the regions close the the"
3011 " winding, and then the mesh size will increate to maximum mesh element"
3012 " size as it gets away from the winding."
3013 ),
3014 )
3015 maximumElementSize: PositiveFloat = Field(
3016 title="Maximum Element Size",
3017 description=(
3018 "The maximum mesh element size in terms of the largest mesh size in the"
3019 " winding. This mesh size will be used in the regions close the the"
3020 " winding, and then the mesh size will increate to maximum mesh element"
3021 " size as it gets away from the winding."
3022 ),
3023 )
3025 @field_validator("maximumElementSize")
3026 @classmethod
3027 def check_rel_size_max(cls, maximumElementSize, info: ValidationInfo):
3028 if maximumElementSize < info.data["minimumElementSize"]:
3029 raise ValueError(
3030 "Maximum mesh element size (maximumElementSize) cannot be smaller than"
3031 " the minimum mesh element size (minimumElementSize)!"
3032 )
3034 return maximumElementSize
3036 @computed_field
3037 @cached_property
3038 def sizeMin(self) -> float:
3039 """Return the minimum mesh element size in real dimensions."""
3040 geometry = geometry_input.get()
3042 meshLengthsPerElement = []
3044 # azimuthal elements:
3045 for numberOfElements in self.winding.azimuthalNumberOfElementsPerTurn:
3046 terminalOuterRadius = (
3047 getWindingOuterRadius()
3048 + 2 * geometry["terminals"]["outer"]["thickness"]
3049 )
3050 meshLengthsPerElement.append(
3051 2 * math.pi * terminalOuterRadius / numberOfElements
3052 )
3054 # radial elements:
3055 for numberOfElements in self.winding.radialNumberOfElementsPerTurn:
3056 meshLengthsPerElement.append(
3057 geometry["winding"]["thickness"] / numberOfElements
3058 )
3060 if not geometry["contactLayer"]["thinShellApproximation"]:
3061 # radial elements:
3062 for numberOfElements in self.contactLayer.radialNumberOfElementsPerTurn:
3063 meshLengthsPerElement.append(
3064 geometry["contactLayer"]["thickness"] / numberOfElements
3065 )
3067 # axial elements:
3068 for numberOfElements in self.winding.axialNumberOfElements:
3069 meshLengthsPerElement.append(
3070 geometry["winding"]["height"] / numberOfElements
3071 )
3073 return max(meshLengthsPerElement) * self.minimumElementSize
3075 @computed_field
3076 @cached_property
3077 def sizeMax(self) -> float:
3078 """Return the minimum mesh element size in real dimensions."""
3079 geometry = geometry_input.get()
3081 meshLengthsPerElement = []
3083 # azimuthal elements:
3084 for numberOfElements in self.winding.azimuthalNumberOfElementsPerTurn:
3085 terminalOuterRadius = (
3086 getWindingOuterRadius()
3087 + 2 * geometry["terminals"]["outer"]["thickness"]
3088 )
3089 meshLengthsPerElement.append(
3090 2 * math.pi * terminalOuterRadius / numberOfElements
3091 )
3093 # radial elements:
3094 for numberOfElements in self.winding.radialNumberOfElementsPerTurn:
3095 meshLengthsPerElement.append(
3096 geometry["winding"]["thickness"] / numberOfElements
3097 )
3099 if not geometry["contactLayer"]["thinShellApproximation"]:
3100 # radial elements:
3101 for numberOfElements in self.contactLayer.radialNumberOfElementsPerTurn:
3102 meshLengthsPerElement.append(
3103 geometry["contactLayer"]["thickness"] / numberOfElements
3104 )
3106 # axial elements:
3107 for numberOfElements in self.winding.axialNumberOfElements:
3108 meshLengthsPerElement.append(
3109 geometry["winding"]["height"] / numberOfElements
3110 )
3112 return max(meshLengthsPerElement) * self.maximumElementSize
3114 @computed_field
3115 @cached_property
3116 def stopGrowingDistance(self) -> float:
3117 """Return the distance from the pancake coils to start growing the mesh elements."""
3118 geometry = geometry_input.get()
3119 terminalOuterRadius = (
3120 getWindingOuterRadius() + 2 * geometry["terminals"]["outer"]["thickness"]
3121 )
3123 if geometry["air"]["type"] == "cuboid":
3124 sideLength = geometry["air"]["sideLength"]
3125 stopGrowingDistance = sideLength / 2 - terminalOuterRadius
3126 if geometry["air"]["type"] == "cylinder":
3127 stopGrowingDistance = geometry["air"]["radius"] - terminalOuterRadius
3129 return stopGrowingDistance
3131 @computed_field
3132 @cached_property
3133 def startGrowingDistance(self) -> float:
3134 geometry = geometry_input.get()
3135 terminalOuterRadius = (
3136 getWindingOuterRadius() + 2 * geometry["terminals"]["outer"]["thickness"]
3137 )
3138 terminalInnerRadius = (
3139 geometry["winding"]["innerRadius"]
3140 - 2 * geometry["terminals"]["inner"]["thickness"]
3141 )
3143 startGrowingDistance = (terminalOuterRadius - terminalInnerRadius) / 2
3145 return startGrowingDistance
3148class Pancake3DSolve(BaseModel):
3149 # 1) User inputs:
3150 time: Pancake3DSolveTime = Field(
3151 title="Time Settings",
3152 description="All the time related settings for transient analysis.",
3153 )
3155 nonlinearSolver: Optional[Pancake3DSolveNonlinearSolverSettings] = Field(
3156 title="Nonlinear Solver Settings",
3157 description="All the nonlinear solver related settings.",
3158 )
3160 winding: Pancake3DSolveWindingMaterial = Field(
3161 title="Winding Properties",
3162 description="This dictionary contains the winding material properties.",
3163 )
3164 contactLayer: Pancake3DSolveContactLayerMaterial = Field(
3165 title="Contact Layer Properties",
3166 description="This dictionary contains the contact layer material properties.",
3167 )
3168 terminals: Pancake3DSolveTerminalMaterialAndBoundaryCondition = Field(
3169 title="Terminals Properties",
3170 description=(
3171 "This dictionary contains the terminals material properties and cooling"
3172 " condition."
3173 ),
3174 )
3175 air: Pancake3DSolveAir = Field(
3176 title="Air Properties",
3177 description="This dictionary contains the air material properties.",
3178 )
3180 initialConditions: Pancake3DSolveInitialConditions = Field(
3181 title="Initial Conditions",
3182 description="Initial conditions of the problem.",
3183 )
3185 boundaryConditions: Union[Literal["vanishingTangentialElectricField"], Pancake3DSolveImposedField] = Field(
3186 default="vanishingTangentialElectricField",
3187 title="Boundary Conditions",
3188 description="Boundary conditions of the problem.",
3189 )
3191 quantitiesToBeSaved: list[Pancake3DSolveSaveQuantity] = Field(
3192 default=None,
3193 title="Quantities to be Saved",
3194 description="List of quantities to be saved.",
3195 )
3197 # Mandatory:
3198 type: Literal["electromagnetic", "thermal", "weaklyCoupled", "stronglyCoupled"] = (
3199 Field(
3200 title="Simulation Type",
3201 description=(
3202 "FiQuS/Pancake3D can solve only electromagnetic and thermal or"
3203 " electromagnetic and thermal coupled. In the weaklyCoupled setting,"
3204 " thermal and electromagnetics systems will be put into different"
3205 " matrices, whereas in the stronglyCoupled setting, they all will be"
3206 " combined into the same matrix. The solution should remain the same."
3207 ),
3208 )
3209 )
3211 # Optionals:
3212 proTemplate: str = Field(
3213 default="Pancake3D_template.pro",
3214 description="file name of the .pro template file",
3215 )
3217 localDefects: Pancake3DSolveLocalDefects = Field(
3218 default=Pancake3DSolveLocalDefects(),
3219 title="Local Defects",
3220 description=(
3221 "Local defects (like making a small part of the winding normal conductor at"
3222 " some time) can be introduced."
3223 ),
3224 )
3226 initFromPrevious: str = Field(
3227 default="",
3228 title="Full path to res file to continue from",
3229 description=(
3230 "The simulation is continued from an existing .res file. The .res file is"
3231 " from a previous computation on the same geometry and mesh. The .res file"
3232 " is taken from the folder Solution_<<initFromPrevious>>."
3233 " IMPORTANT: When the option is used, the start time should be identical to the last "
3234 " time value for the <<initFromPrevious>> simulation."
3235 ),
3236 )
3238 isothermalInAxialDirection: bool = Field(
3239 default=False,
3240 title="Equate DoF along axial direction",
3241 description=(
3242 "If True, the DoF along the axial direction will be equated. This means"
3243 " that the temperature will be the same along the axial direction reducing"
3244 " the number of DoF. This is only valid for the thermal analysis."
3245 ),
3246 )
3248 voltageTapPositions: Optional[list[Pancake3DPosition]] = Field(
3249 default=[],
3250 title="Voltage Tap Positions",
3251 description=(
3252 "List of voltage tap positions. The "
3253 "position can be given in the form of a list of [x, y, z] coordinates or "
3254 "as turnNumber and number of pancake coil."
3255 ),
3256 )
3258 EECircuit: Optional[Pancake3DSolveEECircuit] = Field(
3259 default = Pancake3DSolveEECircuit(),
3260 title="Detection Circuit",
3261 description=(
3262 "This dictionary contains the detection circuit settings."
3263 ),
3264 )
3266 noOfMPITasks: Optional[Union[bool, int]] = Field(
3267 default=False,
3268 title="No. of tasks for MPI parallel run of GetDP",
3269 description=(
3270 "If integer, GetDP will be run in parallel using MPI. This is only valid"
3271 " if MPI is installed on the system and an MPI-enabled GetDP is used."
3272 " If False, GetDP will be run in serial without invoking mpiexec."
3273 ),
3274 )
3276 resistiveHeatingTerminals: bool = Field(
3277 default=True,
3278 title="Resistive Heating in Terminals",
3279 description=(
3280 "If True, terminals are subject to Joule heating. If False, terminal regions are"
3281 " not subject to Joule heating. In both cases, heat conduction through the terminal is "
3282 " considered."
3283 ),
3284 )
3286 solveHeatEquationTerminalsTransitionNotch: bool = Field(
3287 default=True,
3288 title="Solve Heat Equation in Terminals",
3289 description=(
3290 "If True, the heat equation is solved in the terminals and transition notch."
3291 "If False, the heat equation is not solved in the terminals and transition notches."
3292 "In the latter case, neither heat conduction nor generation are considered."
3293 "In other words, the temperature is not an unknown of the problem in the terminals."
3294 ),
3295 )
3297 heatFlowBetweenTurns: bool = Field(
3298 default=True,
3299 title="Heat Equation Between Turns",
3300 description=(
3301 "If True, heat flow between turns is considered. If False, it is not considered. "
3302 "In the latter case, heat conduction is only considered to the middle of the winding in the thin shell approximation "
3303 "in order to keep the thermal mass of the insulation included. In the middle between the turns, an adiabatic condition is applied. "
3304 "Between the turns refers to the region between the winding turns, NOT to the region between terminals "
3305 "and the first and last turn. "
3306 "This feature is only implemented for the thin shell approximation."
3307 ),
3308 )
3310 convectiveCooling: Optional[Pancake3DSolveConvectiveCooling] = Field(
3311 default=Pancake3DSolveConvectiveCooling(),
3312 title="Convective Cooling",
3313 description=(
3314 "This dictionary contains the convective cooling settings."
3315 ),
3316 )
3318 imposedPowerDensity: Optional[Pancake3DSolvePowerDensity] = Field(
3319 default=None,
3320 title="Power Density",
3321 description=(
3322 "The power density for an imposed power density in the winding."
3323 ),
3324 )
3326 materialParametersUseCoilField: bool = Field(
3327 default=True,
3328 title="Use Coil Field for Critical Current",
3329 description=(
3330 "If True, the total field (i.e., coil field plus potentially imposed field)"
3331 "will be used for the material (default)."
3332 "If False, only the imposed field (can be zero) will be used."
3333 ),
3334 )
3336 stopWhenTemperatureReaches: Optional[float] = Field(
3337 default=0,
3338 title="Stop When Temperature Reaches",
3339 description=(
3340 "If the maximum temperature reaches this value, the simulation will"
3341 " be stopped."
3342 ),
3343 )
3345 @computed_field
3346 @cached_property
3347 def systemsOfEquationsType(self) -> str:
3349 windingMaterialIsGiven = self.winding.material is not None
3350 contactLayerMaterialIsGiven = self.contactLayer.material is not None
3351 terminalMaterialIsGiven = self.terminals.material is not None
3352 notchMaterialIsGiven = self.terminals.transitionNotch.material is not None
3353 terminalContactLayerMaterialIsGiven = self.terminals.terminalContactLayer.material is not None
3355 if (
3356 windingMaterialIsGiven
3357 or contactLayerMaterialIsGiven
3358 or terminalMaterialIsGiven
3359 or notchMaterialIsGiven
3360 or terminalContactLayerMaterialIsGiven
3361 ):
3362 systemsOfEquationsType = "nonlinear"
3363 else:
3364 systemsOfEquationsType = "linear"
3366 return systemsOfEquationsType
3368 @model_validator(mode="before")
3369 @classmethod
3370 def check_nls_system_tolerances(cls, model: "Pancake3DSolve"):
3372 if not "nonlinearSolver" in model or not "tolerances" in model["nonlinearSolver"]:
3373 return model
3375 all_tolerances = model["nonlinearSolver"]["tolerances"]
3377 for tolerance in all_tolerances:
3378 if tolerance["quantity"] == "electromagneticSolutionVector":
3379 if model["type"] == "thermal" or model["type"] == "stronglyCoupled":
3380 raise ValueError(
3381 "Nonlinear iteration:"
3382 "The 'electromagneticSolutionVector' tolerance can be used only"
3383 " in 'electromagnetic' or 'weaklyCoupled' simulations."
3384 )
3386 if tolerance["quantity"] == "thermalSolutionVector":
3387 if model["type"] == "electromagnetic" or model["type"] == "stronglyCoupled":
3388 raise ValueError(
3389 "Nonlinear iteration:"
3390 "The 'thermalSolutionVector' tolerance can be used only"
3391 " in 'thermal' or 'weaklyCoupled' simulations."
3392 )
3394 if tolerance["quantity"] == "coupledSolutionVector":
3395 if model["type"] == "electromagnetic" or model["type"] == "thermal" or model["type"] == "weaklyCoupled":
3396 raise ValueError(
3397 "Nonlinear iteration:"
3398 "The 'coupledSolutionVector' tolerance can be used only"
3399 " in 'stronglyCoupled' simulations."
3400 )
3401 return model
3403 @model_validator(mode="before")
3404 @classmethod
3405 def check_adaptive_system_tolerances(cls, model: "Pancake3DSolve"):
3407 if model["time"]["timeSteppingType"] == "fixed":
3408 return model
3410 all_tolerances = model["time"]["adaptiveSteppingSettings"]["tolerances"]
3412 for tolerance in all_tolerances:
3413 if tolerance["quantity"] == "electromagneticSolutionVector":
3414 if model["type"] == "thermal" or model["type"] == "stronglyCoupled":
3415 raise ValueError(
3416 "Adaptive time stepping:"
3417 "The 'electromagneticSolutionVector' tolerance can be used only"
3418 " in 'electromagnetic' or 'weaklyCoupled' simulations."
3419 )
3421 if tolerance["quantity"] == "thermalSolutionVector":
3422 if model["type"] == "electromagnetic" or model["type"] == "stronglyCoupled":
3423 raise ValueError(
3424 "Adaptive time stepping:"
3425 "The 'thermalSolutionVector' tolerance can be used only"
3426 " in 'thermal' or 'weaklyCoupled' simulations."
3427 )
3429 if tolerance["quantity"] == "coupledSolutionVector":
3430 if model["type"] == "electromagnetic" or model["type"] == "thermal" or model["type"] == "weaklyCoupled":
3431 raise ValueError(
3432 "Adaptive time stepping:"
3433 "The 'coupledSolutionVector' tolerance can be used only"
3434 " in 'stronglyCoupled' simulations."
3435 )
3436 return model
3438 @model_validator(mode="before")
3439 @classmethod
3440 def check_convective_cooling(cls, model: "Pancake3DSolve"):
3441 if "convectiveCooling" in model:
3443 convectiveCooling = model["convectiveCooling"]
3445 if not "heatTransferCoefficient" in convectiveCooling or not "exteriorBathTemperature" in convectiveCooling:
3446 raise ValueError(
3447 "If convective cooling is turned on, the heatTransferCoefficient and exteriorBathTemperature must be provided."
3448 )
3450 return model
3452 @model_validator(mode="before")
3453 @classmethod
3454 def check_power_density(cls, model: "Pancake3DSolve"):
3455 if "imposedPowerDensity" in model:
3457 powerDensity = model["imposedPowerDensity"]
3459 if powerDensity:
3460 if not "startTime" in powerDensity or not "endTime" in powerDensity or not "startArcLength" in powerDensity or not "endArcLength" in powerDensity or not "power" in powerDensity:
3461 raise ValueError(
3462 "If the power density provided, the power, startTime, endTime, startTurn, and endTurn must be provided."
3463 )
3465 return model
3467 @model_validator(mode="before")
3468 @classmethod
3469 def check_heat_equation_between_turns(cls, model: "Pancake3DSolve"):
3470 if "heatFlowBetweenTurns" in model and not model["heatFlowBetweenTurns"]:
3472 geometry = geometry_input.get()
3474 if not geometry["contactLayer"]["thinShellApproximation"]:
3475 raise ValueError(
3476 "The heat equation between turns can only be turned off in the thin shell approximation."
3477 )
3479 if not model["contactLayer"]["resistivity"] == "perfectlyInsulating":
3480 raise ValueError(
3481 "The heat equation between turns can only be turned off for perfectly insulated coils. That is, coils with contactLayer resistivity set to 'perfectlyInsulating'."
3482 )
3484 return model
3486 @model_validator(mode="before")
3487 @classmethod
3488 def check_voltage_taps(cls, model: "Pancake3DSolve"):
3490 if "voltageTapPositions" in model:
3491 if len(model['voltageTapPositions'])>0:
3492 if model["type"] == "thermal":
3493 raise ValueError(
3494 "Voltage taps can only be used in electromagnetic simulations."
3495 )
3497 if not model["contactLayer"]["resistivity"] == "perfectlyInsulating":
3498 raise ValueError(
3499 "For now, voltage taps can only be used for perfectly insulated coils. That is, coils with contactLayer resistivity set to 'perfectlyInsulating'."
3500 )
3502 geometry = geometry_input.get()
3504 if not geometry["contactLayer"]["thinShellApproximation"]:
3505 raise ValueError(
3506 "Voltage taps are only available in the thin shell approximation for now."
3507 )
3510 return model
3512 # TODO: add model_validator to check if MPI is installed on the system and an MPI-enabled GetDP is used if useMPI is True
3514 # TODO: add model_validator to check postprocess quantities are available for this solve type (e.g. thermal quantities cannot be chosen for electromagnetic solve)
3516 # TODO: add model_validator to check convergence quantities are available for this solve type (e.g. thermal quantities cannot be chosen for electromagnetic solve)
3518class Pancake3DPostprocess(BaseModel):
3519 """
3520 TO BE UPDATED
3521 """
3523 # 1) User inputs:
3524 timeSeriesPlots: Optional[list[Pancake3DPostprocessTimeSeriesPlot]] = Field(
3525 default=None,
3526 title="Time Series Plots",
3527 description="Values can be plotted with respect to time.",
3528 )
3530 magneticFieldOnCutPlane: Optional[Pancake3DPostprocessMagneticFieldOnPlane] = Field(
3531 default=None,
3532 title="Magnetic Field on a Cut Plane",
3533 description=(
3534 "Color map of the magnetic field on the YZ plane can be plotted with"
3535 " streamlines."
3536 ),
3537 )
3540class Pancake3D(BaseModel):
3541 """
3542 Level 1: Class for FiQuS Pancake3D
3543 """
3545 type: Literal["Pancake3D"]
3546 geometry: Pancake3DGeometry = Field(
3547 default=None,
3548 title="Geometry",
3549 description="This dictionary contains the geometry information.",
3550 )
3551 mesh: Pancake3DMesh = Field(
3552 default=None,
3553 title="Mesh",
3554 description="This dictionary contains the mesh information.",
3555 )
3556 solve: Pancake3DSolve = Field(
3557 default=None,
3558 title="Solve",
3559 description="This dictionary contains the solve information.",
3560 )
3561 postproc: Pancake3DPostprocess = Field(
3562 default=Pancake3DPostprocess(),
3563 title="Postprocess",
3564 description="This dictionary contains the postprocess information.",
3565 )
3566 input_file_path: str = Field(
3567 default=None,
3568 description="path of the input file (calculated by FiQuS)",
3569 exclude=True,
3570 )
3572 @model_validator(mode="before")
3573 @classmethod
3574 def set_context_variables(cls, model: "Pancake3D"):
3575 """Set the context variables (geometry data model, mesh data model, solve data
3576 model) of the Pancake3D model. They will be used in the submodel validators.
3577 """
3579 if isinstance(
3580 model["mesh"]["winding"]["azimuthalNumberOfElementsPerTurn"], int
3581 ):
3582 model["mesh"]["winding"]["azimuthalNumberOfElementsPerTurn"] = [
3583 model["mesh"]["winding"]["azimuthalNumberOfElementsPerTurn"]
3584 ] * model["geometry"]["numberOfPancakes"]
3586 if isinstance(model["mesh"]["winding"]["radialNumberOfElementsPerTurn"], int):
3587 model["mesh"]["winding"]["radialNumberOfElementsPerTurn"] = [
3588 model["mesh"]["winding"]["radialNumberOfElementsPerTurn"]
3589 ] * model["geometry"]["numberOfPancakes"]
3591 if isinstance(model["mesh"]["winding"]["axialNumberOfElements"], int):
3592 model["mesh"]["winding"]["axialNumberOfElements"] = [
3593 model["mesh"]["winding"]["axialNumberOfElements"]
3594 ] * model["geometry"]["numberOfPancakes"]
3596 if isinstance(model["mesh"]["winding"]["elementType"], str):
3597 model["mesh"]["winding"]["elementType"] = [
3598 model["mesh"]["winding"]["elementType"]
3599 ] * model["geometry"]["numberOfPancakes"]
3601 if isinstance(
3602 model["mesh"]["contactLayer"]["radialNumberOfElementsPerTurn"], int
3603 ):
3604 model["mesh"]["contactLayer"]["radialNumberOfElementsPerTurn"] = [
3605 model["mesh"]["contactLayer"]["radialNumberOfElementsPerTurn"]
3606 ] * model["geometry"]["numberOfPancakes"]
3608 geometry_input.set(model["geometry"])
3609 mesh_input.set(model["mesh"])
3610 solve_input.set(model["solve"])
3611 input_file_path.set(model["input_file_path"])
3613 return model