Coverage for fiqus/data/DataFiQuSPancake3D.py: 79%
1010 statements
« prev ^ index » next coverage.py v6.4.4, created at 2024-05-20 03:24 +0200
1from typing import Literal, Optional, Annotated
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"
30]
31SuperconductingMaterialName = Literal["HTSSuperPower"]
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}
51heatCapacityMacroNames = {
52 "Copper": "MATERIAL_SpecificHeatCapacity_Copper_T",
53 "Hastelloy": "MATERIAL_SpecificHeatCapacity_Hastelloy_T",
54 "Silver": "MATERIAL_SpecificHeatCapacity_Silver_T",
55 "Indium": "MATERIAL_SpecificHeatCapacity_Indium_T",
56 "Stainless Steel": "MATERIAL_SpecificHeatCapacity_SSteel_T",
57}
58getdpTSAOnlyResistivityFunctions = {
59 "Indium": "TSA_CFUN_rhoIn_T_constantThickness_fct_only",
60 "Stainless Steel": None,
61}
62getdpTSAMassResistivityFunctions = {
63 "Indium": "TSA_CFUN_rhoIn_T_constantThickness_mass",
64 "Stainless Steel": None,
65}
66getdpTSAStiffnessResistivityFunctions = {
67 "Indium": "TSA_CFUN_rhoIn_T_constantThickness_stiffness",
68 "Stainless Steel": None,
69}
70getdpTSAMassThermalConductivityFunctions = {
71 "Indium": "TSA_CFUN_kIn_constantThickness_mass",
72 "Stainless Steel": "TSA_CFUN_kSteel_T_constantThickness_mass",
73}
74getdpTSAStiffnessThermalConductivityFunctions = {
75 "Indium": "TSA_CFUN_kIn_constantThickness_stiffness",
76 "Stainless Steel": "TSA_CFUN_kSteel_T_constantThickness_stiffness",
77}
78getdpTSAMassHeatCapacityFunctions = {
79 "Indium": "TSA_CFUN_CvIn_constantThickness_mass",
80 "Stainless Steel": "TSA_CFUN_CvSteel_T_constantThickness_mass",
81}
82getdpTSARHSFunctions = {
83 "Indium": "TSA_CFUN_rhoIn_T_constantThickness_rhs",
84 "Stainless Steel": None,
85}
86getdpTSATripleFunctions = {
87 "Indium": "TSA_CFUN_rhoIn_T_constantThickness_triple",
88 "Stainless Steel": None,
89}
90getdpCriticalCurrentDensityFunctions = {
91 "HTSSuperPower": "CFUN_HTS_JcFit_SUPERPOWER_T_B_theta",
92}
93# ======================================================================================
94# ======================================================================================
96# ======================================================================================
97# Available quantities: ================================================================
98PositionRequiredQuantityName = Literal[
99 "magneticField",
100 "axialComponentOfTheMagneticField",
101 "magnitudeOfMagneticField",
102 "currentDensity",
103 "magnitudeOfCurrentDensity",
104 "resistiveHeating",
105 "temperature",
106 "criticalCurrentDensity",
107 "heatFlux",
108 "resistivity",
109 "thermalConductivity",
110 "specificHeatCapacity",
111 "jHTSOverjCritical",
112 "criticalCurrent",
113 "debug",
114]
115PositionNotRequiredQuantityName = Literal[
116 "currentThroughCoil",
117 "voltageBetweenTerminals",
118 "inductance",
119 "timeConstant",
120 "totalResistiveHeating",
121 "magneticEnergy",
122]
123# ======================================================================================
124# ======================================================================================
126# ======================================================================================
127# Quantity information: ================================================================
128EMQuantities = [
129 "magneticField",
130 "magnitudeOfMagneticField",
131 "currentDensity",
132 "magnitudeOfCurrentDensity",
133 "resistiveHeating",
134 "criticalCurrentDensity",
135 "resistivity",
136 "jHTSOverjCritical",
137 "criticalCurrent",
138 "debug",
139 "inductance",
140 "timeConstant",
141 "currentThroughCoil",
142 "voltageBetweenTerminals",
143 "totalResistiveHeating",
144 "magneticEnergy",
145]
146ThermalQuantities = [
147 "temperature",
148 "heatFlux",
149 "thermalConductivity",
150 "specificHeatCapacity",
151 "debug",
152]
153quantityProperNames = {
154 "axialComponentOfTheMagneticField": "Axila Component of the Magnetic Field",
155 "magneticField": "Magnetic Field",
156 "magnitudeOfMagenticField": "Magnitude of Magnetic Field",
157 "currentDensity": "Current Density",
158 "magnitudeOfCurrentDensity": "Magnitude of Current Density",
159 "resistiveHeating": "Resistive Heating",
160 "temperature": "Temperature",
161 "currentThroughCoil": "Current Through Coil",
162 "voltageBetweenTerminals": "Voltage Between Terminals",
163 "criticalCurrentDensity": "Critical Current Density",
164 "heatFlux": "Heat Flux",
165 "resistivity": "Resistivity",
166 "thermalConductivity": "Thermal Conductivity",
167 "specificHeatCapacity": "Specific Heat Capacity",
168 "jHTSOverjCritical": "jHTS/jCritical",
169 "criticalCurrent": "Critical Current",
170 "debug": "Debug",
171 "inductance": "Inductance",
172 "timeConstant": "Time Constant",
173}
175quantityUnits = {
176 "axialComponentOfTheMagneticField": "T",
177 "magneticField": "T",
178 "magnitudeOfMagneticField": "T",
179 "currentDensity": "$A/m^2$",
180 "magnitudeOfCurrentDensity": "$A/m^2$",
181 "resistiveHeating": "W",
182 "temperature": "K",
183 "currentThroughCoil": "A",
184 "voltageBetweenTerminals": "V",
185 "criticalCurrentDensity": "$A/m^2$",
186 "heatFlux": "$W/m^2$",
187 "resistivity": "Ohm*m",
188 "thermalConductivity": "W/m*K",
189 "specificHeatCapacity": "J/kg*K",
190 "jHTSOverjCritical": "A/A",
191 "criticalCurrent": "A",
192 "debug": "1",
193 "inductance": "H",
194 "timeConstant": "s",
195}
197getdpQuantityNames = {
198 "axialComponentOfTheMagneticField": "RESULT_axialComponentOfTheMagneticField",
199 "magneticField": "RESULT_magneticField",
200 "magnitudeOfMagneticField": "RESULT_magnitudeOfMagneticField",
201 "currentDensity": "RESULT_currentDensity",
202 "magnitudeOfCurrentDensity": "RESULT_magnitudeOfCurrentDensity",
203 "resistiveHeating": "RESULT_resistiveHeating",
204 "temperature": "RESULT_temperature",
205 "currentThroughCoil": "RESULT_currentThroughCoil",
206 "voltageBetweenTerminals": "RESULT_voltageBetweenTerminals",
207 "criticalCurrentDensity": "RESULT_criticalCurrentDensity",
208 "heatFlux": "RESULT_heatFlux",
209 "resistivity": "RESULT_resistivity",
210 "thermalConductivity": "RESULT_thermalConductivity",
211 "specificHeatCapacity": "RESULT_specificHeatCapacity",
212 "jHTSOverjCritical": "RESULT_jHTSOverjCritical",
213 "criticalCurrent": "RESULT_criticalCurrent",
214 "debug": "RESULT_debug",
215 "inductance": "RESULT_inductance",
216 "timeConstant": "RESULT_timeConstant",
217}
219getdpPostOperationNames = {
220 "axialComponentOfTheMagneticField": "POSTOP_axialComponentOfTheMagneticField",
221 "magneticField": "POSTOP_magneticField",
222 "magnitudeOfMagneticField": "POSTOP_magnitudeOfMagneticField",
223 "currentDensity": "POSTOP_currentDensity",
224 "magnitudeOfCurrentDensity": "POSTOP_magnitudeOfCurrentDensity",
225 "resistiveHeating": "POSTOP_resistiveHeating",
226 "temperature": "POSTOP_temperature",
227 "currentThroughCoil": "POSTOP_currentThroughCoil",
228 "voltageBetweenTerminals": "POSTOP_voltageBetweenTerminals",
229 "criticalCurrentDensity": "POSTOP_criticalCurrentDensity",
230 "heatFlux": "POSTOP_heatFlux",
231 "resistivity": "POSTOP_resistivity",
232 "thermalConductivity": "POSTOP_thermalConductivity",
233 "specificHeatCapacity": "POSTOP_specificHeatCapacity",
234 "jHTSOverjCritical": "POSTOP_jHTSOverjCritical",
235 "criticalCurrent": "POSTOP_criticalCurrent",
236 "debug": "POSTOP_debug",
237 "inductance": "POSTOP_inductance",
238 "timeConstant": "POSTOP_timeConstant",
239}
241# ======================================================================================
242# ======================================================================================
244# Global variables
245geometry_input = ContextVar("geometry")
246mesh_input = ContextVar("mesh")
247solve_input = ContextVar("solve")
248input_file_path = ContextVar("input_file_path")
249all_break_points = []
252def getWindingOuterRadius():
253 """Return outer radius of the winding."""
254 geometry = geometry_input.get()
255 return (
256 geometry["winding"]["innerRadius"]
257 + geometry["winding"]["thickness"]
258 + geometry["winding"]["numberOfTurns"]
259 * (geometry["winding"]["thickness"] + geometry["contactLayer"]["thickness"])
260 )
263def getAirHeight():
264 """Return the height of the air."""
265 geometry = geometry_input.get()
266 h = (
267 geometry["numberOfPancakes"]
268 * (geometry["winding"]["height"] + geometry["gapBetweenPancakes"])
269 - geometry["gapBetweenPancakes"]
270 + 2 * geometry["air"]["axialMargin"]
271 )
272 return h
275def getTransitionNotchAngle():
276 """Return transition notch angle of the winding."""
277 mesh = mesh_input.get()
279 azimuthalNumberOfElementsPerTurn = max(
280 mesh["winding"]["azimuthalNumberOfElementsPerTurn"]
281 )
283 transitionNotchAngle = 2 * math.pi / azimuthalNumberOfElementsPerTurn
285 return transitionNotchAngle
288def checkIfAirOrTerminalMeshIsStructured():
289 geometry = geometry_input.get()
290 mesh = mesh_input.get()
292 structuredAirMesh = False
293 structuredTerminalMesh = False
294 if "air" in mesh:
295 structuredAirMesh = mesh["air"]["structured"]
296 if "terminals" in mesh:
297 structuredTerminalMesh = mesh["terminals"]["structured"]
298 structuredMesh = structuredAirMesh or structuredTerminalMesh
300 return structuredMesh
303# ======================================================================================
304# FUNDAMENTAL CLASSES STARTS ===========================================================
305# ======================================================================================
306class Pancake3DPositionInCoordinates(BaseModel):
307 x: float = Field(
308 title="x coordinate",
309 description="x coordinate of the position.",
310 )
311 y: float = Field(
312 title="y coordinate",
313 description="y coordinate of the position.",
314 )
315 z: float = Field(
316 title="z coordinate",
317 description="z coordinate of the position.",
318 )
321class Pancake3DPositionInTurnNumbers(BaseModel):
322 turnNumber: float = Field(
323 title="Turn Number",
324 description=(
325 "Winding turn number as a position input. It starts from 0 and it can be a"
326 " float."
327 ),
328 )
329 whichPancakeCoil: Optional[PositiveInt] = Field(
330 default=None,
331 title="Pancake Coil Number",
332 description="The first pancake coil is 1, the second is 2, etc.",
333 )
335 @field_validator("turnNumber")
336 @classmethod
337 def check_turnNumber(cls, turnNumber):
338 geometry = geometry_input.get()
340 if turnNumber < 0:
341 raise ValueError("Turn number cannot be less than 0.")
342 elif turnNumber > geometry["winding"]["numberOfTurns"]:
343 raise ValueError(
344 "Turn number cannot be greater than the number of turns of the winding"
345 f" ({geometry['numberOfPancakes']})."
346 )
348 return turnNumber
350 @field_validator("whichPancakeCoil")
351 @classmethod
352 def check_whichPancakeCoil(cls, whichPancakeCoil):
353 geometry = geometry_input.get()
355 if whichPancakeCoil is not None:
356 if whichPancakeCoil < 1:
357 raise ValueError(
358 "Pancake coil numbers start from 1. Therefore, it cannot be less"
359 " than 1."
360 )
361 elif whichPancakeCoil > geometry["numberOfPancakes"]:
362 raise ValueError(
363 "Pancake coil number cannot be greater than the number of pancakes"
364 f" ({geometry['numberOfPancakes']})."
365 )
366 else:
367 return 1
369 return whichPancakeCoil
371 def compute_coordinates(self):
372 geometry = geometry_input.get()
373 mesh = mesh_input.get()
375 if geometry["contactLayer"]["thinShellApproximation"]:
376 windingThickness = (
377 geometry["winding"]["thickness"]
378 + geometry["contactLayer"]["thickness"]
379 * (geometry["winding"]["numberOfTurns"] - 1)
380 / geometry["winding"]["numberOfTurns"]
381 )
382 gapThickness = 0
383 else:
384 windingThickness = geometry["winding"]["thickness"]
385 gapThickness = geometry["contactLayer"]["thickness"]
387 innerRadius = geometry["winding"]["innerRadius"]
388 initialTheta = 0.0
389 if isinstance(mesh["winding"]["azimuthalNumberOfElementsPerTurn"], list):
390 ane = mesh["winding"]["azimuthalNumberOfElementsPerTurn"][0]
391 elif isinstance(mesh["winding"]["azimuthalNumberOfElementsPerTurn"], int):
392 ane = mesh["winding"]["azimuthalNumberOfElementsPerTurn"]
393 else:
394 raise ValueError(
395 "The azimuthal number of elements per turn must be either an integer"
396 " or a list of integers."
397 )
399 numberOfPancakes = geometry["numberOfPancakes"]
400 gapBetweenPancakes = geometry["gapBetweenPancakes"]
401 windingHeight = geometry["winding"]["height"]
403 turnNumber = self.turnNumber
404 whichPancake = self.whichPancakeCoil
406 elementStartTurnNumber = math.floor(turnNumber / (1 / ane)) * (1 / ane)
407 elementEndTurnNumber = elementStartTurnNumber + 1 / ane
409 class point:
410 def __init__(self, x, y, z):
411 self.x = x
412 self.y = y
413 self.z = z
415 def __add__(self, other):
416 return point(self.x + other.x, self.y + other.y, self.z + other.z)
418 def __sub__(self, other):
419 return point(self.x - other.x, self.y - other.y, self.z - other.z)
421 def __mul__(self, scalar):
422 return point(self.x * scalar, self.y * scalar, self.z * scalar)
424 def __truediv__(self, scalar):
425 return point(self.x / scalar, self.y / scalar, self.z / scalar)
427 def rotate(self, degrees):
428 return point(
429 self.x * math.cos(degrees) - self.y * math.sin(degrees),
430 self.x * math.sin(degrees) + self.y * math.cos(degrees),
431 self.z,
432 )
434 def normalize(self):
435 return self / math.sqrt(self.x**2 + self.y**2 + self.z**2)
437 if whichPancake % 2 == 1:
438 # If the spiral is counter-clockwise, the initial theta angle decreases,
439 # and r increases as the theta angle decreases.
440 multiplier = 1
441 elif whichPancake % 2 == 0:
442 # If the spiral is clockwise, the initial theta angle increases, and r
443 # increases as the theta angle increases.
444 multiplier = -1
446 # Mesh element's starting point:
447 elementStartTheta = 2 * math.pi * elementStartTurnNumber * multiplier
448 elementStartRadius = (
449 innerRadius
450 + elementStartTheta
451 / (2 * math.pi)
452 * (gapThickness + windingThickness)
453 * multiplier
454 )
455 elementStartPointX = elementStartRadius * math.cos(
456 initialTheta + elementStartTheta
457 )
458 elementStartPointY = elementStartRadius * math.sin(
459 initialTheta + elementStartTheta
460 )
461 elementStartPointZ = (
462 -(
463 numberOfPancakes * windingHeight
464 + (numberOfPancakes - 1) * gapBetweenPancakes
465 )
466 / 2
467 + windingHeight / 2
468 + (whichPancake - 1) * (windingHeight + gapBetweenPancakes)
469 )
470 elementStartPoint = point(
471 elementStartPointX, elementStartPointY, elementStartPointZ
472 )
474 # Mesh element's ending point:
475 elementEndTheta = 2 * math.pi * elementEndTurnNumber * multiplier
476 elementEndRadius = (
477 innerRadius
478 + elementEndTheta
479 / (2 * math.pi)
480 * (gapThickness + windingThickness)
481 * multiplier
482 )
483 elementEndPointX = elementEndRadius * math.cos(initialTheta + elementEndTheta)
484 elementEndPointY = elementEndRadius * math.sin(initialTheta + elementEndTheta)
485 elementEndPointZ = elementStartPointZ
486 elementEndPoint = point(elementEndPointX, elementEndPointY, elementEndPointZ)
488 turnNumberFraction = (turnNumber - elementStartTurnNumber) / (
489 elementEndTurnNumber - elementStartTurnNumber
490 )
491 location = (
492 elementStartPoint
493 + (elementEndPoint - elementStartPoint) * turnNumberFraction
494 ) + (elementEndPoint - elementStartPoint).rotate(
495 -math.pi / 2
496 ).normalize() * windingThickness / 2 * multiplier
498 return location.x, location.y, location.z
500 @computed_field
501 @cached_property
502 def x(self) -> float:
503 return self.compute_coordinates()[0]
505 @computed_field
506 @cached_property
507 def y(self) -> float:
508 return self.compute_coordinates()[1]
510 @computed_field
511 @cached_property
512 def z(self) -> float:
513 return self.compute_coordinates()[2]
516Pancake3DPosition = Pancake3DPositionInCoordinates | Pancake3DPositionInTurnNumbers
518# ======================================================================================
519# FUNDAMENTAL CLASSES ENDS =============================================================
520# ======================================================================================
523# ======================================================================================
524# GEOMETRY CLASSES STARTS ==============================================================
525# ======================================================================================
526class Pancake3DGeometryWinding(BaseModel):
527 # Mandatory:
528 r_i: PositiveFloat = Field(
529 alias="innerRadius",
530 title="Inner Radius",
531 description="Inner radius of the winding.",
532 )
533 t: PositiveFloat = Field(
534 alias="thickness",
535 title="Winding Thickness",
536 description="Thickness of the winding.",
537 )
538 N: float = Field(
539 alias="numberOfTurns",
540 ge=3,
541 title="Number of Turns",
542 description="Number of turns of the winding.",
543 )
544 h: PositiveFloat = Field(
545 alias="height",
546 title="Winding Height",
547 description="Height/width of the winding.",
548 )
550 # Optionals:
551 name: str = Field(
552 default="winding",
553 title="Winding Name",
554 description="The The name to be used in the mesh..",
555 examples=["winding", "myWinding"],
556 )
557 NofVolPerTurn: int = Field(
558 default=2,
559 validate_default=True,
560 alias="numberOfVolumesPerTurn",
561 ge=2,
562 title="Number of Volumes Per Turn (Advanced Input)",
563 description="The number of volumes per turn (CAD related, not physical).",
564 )
566 @field_validator("NofVolPerTurn")
567 @classmethod
568 def check_NofVolPerTurn(cls, NofVolPerTurn):
569 geometry = geometry_input.get()
570 mesh = mesh_input.get()
572 # Check if the NofVolPerTurn is compatible swith the azimuthal number of
573 # elements per turn:
574 for i, ane in enumerate(mesh["winding"]["azimuthalNumberOfElementsPerTurn"]):
575 if ane % NofVolPerTurn != 0:
576 raise ValueError(
577 "The azimuthal number of elements per turn for the pancake coil"
578 f" number {i+1} is ({ane}), but it must be divisible by the number"
579 f" of volumes per turn ({geometry['winding']['NofVolPerTurn']})!"
580 " So it needs to be rounded to"
581 f" {math.ceil(ane/NofVolPerTurn)*NofVolPerTurn:.5f} or"
582 f" {math.floor(ane/NofVolPerTurn)*NofVolPerTurn:.5f}."
583 )
585 structured = checkIfAirOrTerminalMeshIsStructured()
587 if structured:
588 # If the mesh is structured, the number of volumes per turn must be 4:
589 NofVolPerTurn = 4
591 return NofVolPerTurn
593 @computed_field
594 @cached_property
595 def theta_i(self) -> float:
596 """Return start angle of the winding."""
597 return 0.0
599 @computed_field
600 @cached_property
601 def r_o(self) -> float:
602 """Return outer radius of the winding."""
603 return getWindingOuterRadius()
605 @computed_field
606 @cached_property
607 def turnTol(self) -> float:
608 """Return turn tolerance of the winding."""
609 geometry: Pancake3DGeometry = geometry_input.get()
610 mesh: Pancake3DMesh = mesh_input.get()
612 # Calculate the turn tolerance required due to the geometrymetry input:
613 # Turn tolerance is the smallest turn angle (in turns) that is allowed.
614 if "dimTol" in geometry:
615 dimTol = geometry["dimTol"]
616 else:
617 dimTol = 1e-8
618 turnTol = geometry["winding"]["numberOfTurns"] % 1
619 if math.isclose(turnTol, 0, abs_tol=dimTol):
620 turnTol = 0.5
622 turnTolDueToTransition = getTransitionNotchAngle() / (2 * math.pi)
624 # Calculate the minimum turn tolerance possible due to the mesh input:
625 minimumTurnTol = 1 / min(mesh["winding"]["azimuthalNumberOfElementsPerTurn"])
627 if turnTol < minimumTurnTol:
628 numberOfTurns = geometry["winding"]["numberOfTurns"]
630 raise ValueError(
631 "The azimuthal number of elements per turn for one of the pancakes is"
632 f" {min(mesh['winding']['azimuthalNumberOfElementsPerTurn'])}, and the"
633 " number of turns is"
634 f" {numberOfTurns:.5f}."
635 " The number of turns must always be divisible by the (1/(the"
636 " azimuthal number of elements per turn)) to ensure conformality."
637 " Please change the number of turns or the azimuthal number of"
638 " elemenets per turn. The closest possible number of turns value is"
639 f" {round(numberOfTurns * min(mesh['winding']['azimuthalNumberOfElementsPerTurn']))/min(mesh['winding']['azimuthalNumberOfElementsPerTurn']):.5f}"
640 )
641 else:
642 # Minimum possible sections per turn is 16 (otherwise splines might collide
643 # into each other). But it should be greater than the number of volumes per
644 # turn and it should be divisible by both 1/turnTol and the number of
645 # volumes per turn.
646 sectionsPerTurn = 16
647 if "numberOfVolumesPerTurn" in geometry["winding"]:
648 numberOfVolumesPerTurn = geometry["winding"]["numberOfVolumesPerTurn"]
649 else:
650 numberOfVolumesPerTurn = 2
652 while (
653 (
654 math.fmod(sectionsPerTurn, (1 / turnTol)) > 1e-8
655 and math.fmod(sectionsPerTurn, (1 / turnTol)) - (1 / turnTol)
656 < -1e-8
657 )
658 or (
659 math.fmod(sectionsPerTurn, (1 / turnTolDueToTransition)) > 1e-8
660 and math.fmod(sectionsPerTurn, (1 / turnTolDueToTransition))
661 - (1 / turnTolDueToTransition)
662 < -1e-8
663 )
664 or sectionsPerTurn % numberOfVolumesPerTurn != 0
665 or sectionsPerTurn < numberOfVolumesPerTurn
666 ):
667 sectionsPerTurn += 1
669 # Sections per turn will set the turn tolerance value as well.
670 return 1.0 / sectionsPerTurn
672 @computed_field
673 @cached_property
674 def spt(self) -> float:
675 """Return sections per turn of the winding."""
676 return int(1.0 / self.turnTol)
678 @computed_field
679 @cached_property
680 def totalTapeLength(self) -> float:
681 """Return total tape length of the winding."""
682 geometry: Pancake3DGeometry = geometry_input.get()
684 # Calculate the total tape length of the coil:
686 # The same angle can be subtracted from both theta_1 and theta_2 to simplify the
687 # calculations:
688 theta2 = geometry["winding"]["numberOfTurns"] * 2 * math.pi
689 theta1 = 0
691 # Since r = a * theta + b, r_1 = b since theta_1 = 0:
692 b = geometry["winding"]["innerRadius"]
694 # Since r = a * theta + b, r_2 = a * theta2 + b:
695 a = (getWindingOuterRadius() - b) / theta2
697 def integrand(t):
698 return math.sqrt(a**2 + (a * t + b) ** 2)
700 totalTapeLength = abs(scipy.integrate.quad(integrand, theta1, theta2)[0])
702 return totalTapeLength
705class Pancake3DGeometryContactLayer(BaseModel):
706 # Mandatory:
707 tsa: bool = Field(
708 alias="thinShellApproximation",
709 title="Use Thin Shell Approximation",
710 description=(
711 "If True, the contact layer will be modeled with 2D shell elements (thin"
712 " shell approximation), and if False, the contact layer will be modeled"
713 " with 3D elements."
714 ),
715 )
716 t: PositiveFloat = Field(
717 alias="thickness",
718 title="Contact Layer Thickness",
719 description="Thickness of the contact layer.",
720 )
722 # Optionals:
723 name: str = Field(
724 default="contactLayer",
725 title="Contact Layer Name",
726 description="The name to be used in the mesh.",
727 examples=["myContactLayer"],
728 )
731class Pancake3DGeometryTerminalBase(BaseModel):
732 # Mandatory:
733 t: PositiveFloat = Field(
734 alias="thickness",
735 title="Terminal Thickness",
736 description="Thickness of the terminal's tube.",
737 ) # thickness
739 @field_validator("t")
740 @classmethod
741 def check_t(cls, t):
742 geometry = geometry_input.get()
744 if t < geometry["winding"]["thickness"] / 2:
745 raise ValueError(
746 "Terminal's thickness is smaller than half of the winding's thickness!"
747 " Please increase the terminal's thickness."
748 )
750 return t
753class Pancake3DGeometryInnerTerminal(Pancake3DGeometryTerminalBase):
754 name: str = Field(
755 default="innerTerminal",
756 title="Terminal Name",
757 description="The name to be used in the mesh.",
758 examples=["innerTerminal", "outerTeminal"],
759 )
761 @computed_field
762 @cached_property
763 def r(self) -> float:
764 """Return inner radius of the inner terminal."""
765 geometry = geometry_input.get()
767 innerRadius = geometry["winding"]["innerRadius"] - 2 * self.t
768 if innerRadius < 0:
769 raise ValueError(
770 "Inner terminal's radius is smaller than 0! Please decrease the inner"
771 " terminal's thickness or increase the winding's inner radius."
772 )
774 return innerRadius
777class Pancake3DGeometryOuterTerminal(Pancake3DGeometryTerminalBase):
778 name: str = Field(
779 default="outerTerminal",
780 title="Terminal Name",
781 description="The name to be used in the mesh.",
782 examples=["innerTerminal", "outerTeminal"],
783 )
785 @computed_field
786 @cached_property
787 def r(self) -> float:
788 """Return outer radius of the outer terminal."""
789 outerRadius = getWindingOuterRadius() + 2 * self.t
791 return outerRadius
794class Pancake3DGeometryTerminals(BaseModel):
795 # 1) User inputs:
796 i: Pancake3DGeometryInnerTerminal = Field(alias="inner")
797 o: Pancake3DGeometryOuterTerminal = Field(alias="outer")
799 # Optionals:
800 firstName: str = Field(
801 default="firstTerminal", description="name of the first terminal"
802 )
803 lastName: str = Field(
804 default="lastTerminal", description="name of the last terminal"
805 )
807 @computed_field
808 @cached_property
809 def transitionNotchAngle(self) -> float:
810 """Return transition notch angle of the terminals."""
811 return getTransitionNotchAngle()
814class Pancake3DGeometryAirBase(BaseModel):
815 # Mandatory:
816 margin: PositiveFloat = Field(
817 alias="axialMargin",
818 title="Axial Margin of the Air",
819 description=(
820 "Axial margin between the ends of the air and first/last pancake coils."
821 ),
822 ) # axial margin
824 # Optionals:
825 name: str = Field(
826 default="air",
827 title="Air Name",
828 description="The name to be used in the mesh.",
829 examples=["air", "myAir"],
830 )
831 shellTransformation: bool = Field(
832 default=False,
833 alias="shellTransformation",
834 title="Use Shell Transformation",
835 description=(
836 "Generate outer shell air to apply shell transformation if True (GetDP"
837 " related, not physical)"
838 ),
839 )
840 shellTransformationMultiplier: float = Field(
841 default=1.2,
842 gt=1.1,
843 alias="shellTransformationMultiplier",
844 title="Shell Transformation Multiplier (Advanced Input)",
845 description=(
846 "multiply the air's outer dimension by this value to get the shell's outer"
847 " dimension"
848 ),
849 )
850 cutName: str = Field(
851 default="Air-Cut",
852 title="Air Cut Name",
853 description="name of the cut (cochain) to be used in the mesh",
854 examples=["Air-Cut", "myAirCut"],
855 )
856 shellVolumeName: str = Field(
857 default="air-Shell",
858 title="Air Shell Volume Name",
859 description="name of the shell volume to be used in the mesh",
860 examples=["air-Shell", "myAirShell"],
861 )
862 fragment: bool = Field(
863 default=False,
864 alias="generateGapAirWithFragment",
865 title="Generate Gap Air with Fragment (Advanced Input)",
866 description=(
867 "generate the gap air with gmsh/model/occ/fragment if true (CAD related,"
868 " not physical)"
869 ),
870 )
872 @field_validator("margin")
873 @classmethod
874 def check_margin(cls, margin):
875 geometry = geometry_input.get()
876 windingHeight = geometry["winding"]["height"]
878 if margin < windingHeight / 2:
879 raise ValueError(
880 "Axial margin is smaller than half of the winding's height! Please"
881 " increase the axial margin."
882 )
884 return margin
886 @computed_field
887 @cached_property
888 def h(self) -> float:
889 """Return total height of the air."""
890 h = getAirHeight()
892 return h
895class Pancake3DGeometryAirCylinder(Pancake3DGeometryAirBase):
896 type: Literal["cylinder"] = Field(default="cylinder", title="Air Type")
897 r: PositiveFloat = Field(
898 default=None,
899 alias="radius",
900 title="Air Radius",
901 description="Radius of the air (for cylinder type air).",
902 )
904 @field_validator("r")
905 @classmethod
906 def check_r(cls, r):
907 geometry = geometry_input.get()
908 outerTerminalOuterRadius = (
909 getWindingOuterRadius() + 2 * geometry["terminals"]["outer"]["thickness"]
910 )
912 if r < outerTerminalOuterRadius * 1.5:
913 raise ValueError(
914 "Radius of the air must be at least 1.5 times the outer radius of the"
915 " winding! Please increase the radius of the air."
916 )
918 return r
920 @computed_field
921 @cached_property
922 def shellOuterRadius(self) -> float:
923 """Return outer radius of the air."""
924 shellOuterRadius = self.shellTransformationMultiplier * self.r
926 return shellOuterRadius
929class Pancake3DGeometryAirCuboid(Pancake3DGeometryAirBase):
930 type: Literal["cuboid"] = Field(default="cuboid", title="Air Type")
931 a: PositiveFloat = Field(
932 default=None,
933 alias="sideLength",
934 title="Air Side Length",
935 description="Side length of the air (for cuboid type air).",
936 )
938 @field_validator("a")
939 @classmethod
940 def check_a(cls, a):
941 geometry = geometry_input.get()
942 outerTerminalOuterRadius = (
943 getWindingOuterRadius() + 2 * geometry["terminals"]["outer"]["thickness"]
944 )
946 if a / 2 < outerTerminalOuterRadius * 1.5:
947 raise ValueError(
948 "Half of the side length of the air must be at least 1.5 times the"
949 " outer radius of the winding! Please increase the side length of the"
950 " air."
951 )
953 return a
955 @computed_field
956 @cached_property
957 def shellSideLength(self) -> float:
958 """Return outer radius of the air."""
959 shellSideLength = self.shellTransformationMultiplier * self.a
961 return shellSideLength
964Pancake3DGeometryAir = Annotated[
965 Pancake3DGeometryAirCylinder | Pancake3DGeometryAirCuboid,
966 Field(discriminator="type"),
967]
968# ======================================================================================
969# GEOMETRY CLASSES ENDS ================================================================
970# ======================================================================================
973# ======================================================================================
974# MESH CLASSES STARTS ==================================================================
975# ======================================================================================
976class Pancake3DMeshWinding(BaseModel):
977 # Mandatory:
978 axne: list[PositiveInt] | PositiveInt = Field(
979 alias="axialNumberOfElements",
980 title="Axial Number of Elements",
981 description=(
982 "The number of axial elements for the whole height of the coil. It can be"
983 " either a list of integers to specify the value for each pancake coil"
984 " separately or an integer to use the same setting for each pancake coil."
985 ),
986 )
988 ane: list[PositiveInt] | PositiveInt = Field(
989 alias="azimuthalNumberOfElementsPerTurn",
990 title="Azimuthal Number of Elements Per Turn",
991 description=(
992 "The number of azimuthal elements per turn of the coil. It can be either a"
993 " list of integers to specify the value for each pancake coil separately or"
994 " an integer to use the same setting for each pancake coil."
995 ),
996 )
998 rne: list[PositiveInt] | PositiveInt = Field(
999 alias="radialNumberOfElementsPerTurn",
1000 title="Winding Radial Number of Elements Per Turn",
1001 description=(
1002 "The number of radial elements per tape of the winding. It can be either a"
1003 " list of integers to specify the value for each pancake coil separately or"
1004 " an integer to use the same setting for each pancake coil."
1005 ),
1006 )
1008 # Optionals:
1009 axbc: list[PositiveFloat] | PositiveFloat = Field(
1010 default=[1],
1011 alias="axialDistributionCoefficient",
1012 title="Axial Bump Coefficients",
1013 description=(
1014 "If 1, it won't affect anything. If smaller than 1, elements will get finer"
1015 " in the axial direction at the ends of the coil. If greater than 1,"
1016 " elements will get coarser in the axial direction at the ends of the coil."
1017 " It can be either a list of floats to specify the value for each pancake"
1018 " coil separately or a float to use the same setting for each pancake coil."
1019 ),
1020 )
1022 elementType: (
1023 list[Literal["tetrahedron", "hexahedron", "prism"]]
1024 | Literal["tetrahedron", "hexahedron", "prism"]
1025 ) = Field(
1026 default=["tetrahedron"],
1027 title="Element Type",
1028 description=(
1029 "The element type of windings and contact layers. It can be either a"
1030 " tetrahedron, hexahedron, or a prism. It can be either a list of strings"
1031 " to specify the value for each pancake coil separately or a string to use"
1032 " the same setting for each pancake coil."
1033 ),
1034 )
1036 @field_validator("axne", "ane", "rne", "axbc", "elementType")
1037 @classmethod
1038 def check_inputs(cls, value, info: ValidationInfo):
1039 geometry = geometry_input.get()
1041 numberOfPancakes = geometry["numberOfPancakes"]
1043 structuredMesh = checkIfAirOrTerminalMeshIsStructured()
1045 if not isinstance(value, list):
1046 value = [value] * numberOfPancakes
1047 elif len(value) == 1:
1048 value = value * numberOfPancakes
1049 else:
1050 if len(value) != numberOfPancakes:
1051 raise ValueError(
1052 "The length of the input list must be equal to the number of"
1053 " pancake coils!"
1054 )
1055 if info.field_name == "ane":
1056 if value[0] < 7:
1057 raise ValueError(
1058 "Azimuthal number of elements per turn must be greater than or"
1059 " equal to 7!"
1060 )
1062 if structuredMesh:
1063 if len(set(value)) != 1:
1064 raise ValueError(
1065 "If structured mesh is used, the same mesh setting must be used for"
1066 " all pancake coils!"
1067 )
1069 if info.field_name == "elementType":
1070 if value[0] != "tetrahedron":
1071 raise ValueError(
1072 "If structured air or terminal mesh is used, the element type"
1073 " must be tetrahedron!"
1074 )
1076 if info.field_name == "ane":
1077 if value[0] % 4 != 0:
1078 raise ValueError(
1079 "If structured mesh is used, the number of azimuthal elements"
1080 " per turn must be divisible by 4!"
1081 )
1083 return value
1086class Pancake3DMeshContactLayer(BaseModel):
1087 # Mandatory:
1088 rne: list[PositiveInt] = Field(
1089 alias="radialNumberOfElementsPerTurn",
1090 title="Contact Layer Radial Number of Elements Per Turn",
1091 description=(
1092 "The number of radial elements per tape of the contact layer. It can be"
1093 " either a list of integers to specify the value for each pancake coil"
1094 " separately or an integer to use the same setting for each pancake coil."
1095 ),
1096 )
1098 @field_validator("rne")
1099 @classmethod
1100 def check_inputs(cls, value):
1101 geometry = geometry_input.get()
1103 structuredMesh = checkIfAirOrTerminalMeshIsStructured()
1105 numberOfPancakeCoils = geometry["numberOfPancakes"]
1107 if not isinstance(value, list):
1108 value = [value] * numberOfPancakeCoils
1109 elif len(value) == 1:
1110 value = value * numberOfPancakeCoils
1111 else:
1112 if len(value) != numberOfPancakeCoils:
1113 raise ValueError(
1114 "The length of the input list must be equal to the number of"
1115 " pancake coils!"
1116 )
1118 if structuredMesh:
1119 if len(set(value)) != 1:
1120 raise ValueError(
1121 "If structured mesh is used, the same mesh setting must be used for"
1122 " all pancake coils!"
1123 )
1125 return value
1128class Pancake3DMeshAirAndTerminals(BaseModel):
1129 # Optionals:
1130 structured: bool = Field(
1131 default=False,
1132 title="Structure Mesh",
1133 description=(
1134 "If True, the mesh will be structured. If False, the mesh will be"
1135 " unstructured."
1136 ),
1137 )
1138 radialElementSize: PositiveFloat = Field(
1139 default=1,
1140 title="Radial Element Size",
1141 description=(
1142 "If structured mesh is used, the radial element size can be set. It is the"
1143 " radial element size in terms of the winding's radial element size."
1144 ),
1145 )
1148# ======================================================================================
1149# MESH CLASSES ENDS ====================================================================
1150# ======================================================================================
1153# ======================================================================================
1154# SOLVE CLASSES STARTS =================================================================
1155# ======================================================================================
1156class Pancake3DSolveAir(BaseModel):
1157 # 1) User inputs:
1159 # Mandatory:
1160 permeability: PositiveFloat = Field(
1161 title="Permeability of Air",
1162 description="Permeability of air.",
1163 )
1166class Pancake3DSolveIcVsLength(BaseModel):
1167 lengthValues: list[float] = Field(
1168 title="Tape Length Values",
1169 description="Tape length values that corresponds to criticalCurrentValues.",
1170 )
1171 criticalCurrentValues: list[float] = Field(
1172 title="Critical Current Values",
1173 description="Critical current values that corresponds to lengthValues.",
1174 )
1177class Pancake3DSolveMaterialBase(BaseModel):
1178 name: str
1180 # Optionals:
1181 rrr: PositiveFloat = Field(
1182 default=100,
1183 alias="residualResistanceRatio",
1184 title="Residual Resistance Ratio",
1185 description=(
1186 "Residual-resistivity ratio (also known as Residual-resistance ratio or"
1187 " just RRR) is the ratio of the resistivity of a material at reference"
1188 " temperature and at 0 K."
1189 ),
1190 )
1191 rrrRefT: PositiveFloat = Field(
1192 default=295,
1193 alias="residualResistanceRatioReferenceTemperature",
1194 title="Residual Resistance Ratio Reference Temperature",
1195 description="Reference temperature for residual resistance ratio",
1196 )
1198 @computed_field
1199 @cached_property
1200 def resistivityMacroName(self) -> str:
1201 """Return the resistivity macro name of the material."""
1202 if self.name not in resistivityMacroNames:
1203 return "NOT_DEFINED_IN_DATA_FIQUS_PANCAKE3D"
1205 return resistivityMacroNames[self.name]
1207 @computed_field
1208 @cached_property
1209 def thermalConductivityMacroName(self) -> str:
1210 """Return the thermal conductivity macro name of the material."""
1211 if self.name not in thermalConductivityMacroNames:
1212 return "NOT_DEFINED_IN_DATA_FIQUS_PANCAKE3D"
1214 return thermalConductivityMacroNames[self.name]
1216 @computed_field
1217 @cached_property
1218 def heatCapacityMacroName(self) -> str:
1219 """Return the heat capacity macro name of the material."""
1220 if self.name not in heatCapacityMacroNames:
1221 return "NOT_DEFINED_IN_DATA_FIQUS_PANCAKE3D"
1223 return heatCapacityMacroNames[self.name]
1225 @computed_field
1226 @cached_property
1227 def getdpTSAOnlyResistivityFunction(self) -> str:
1228 """Return the GetDP function name of the material's resistivity."""
1229 if self.name not in getdpTSAOnlyResistivityFunctions:
1230 return "NOT_DEFINED_IN_DATA_FIQUS_PANCAKE3D"
1232 return getdpTSAOnlyResistivityFunctions[self.name]
1234 @computed_field
1235 @cached_property
1236 def getdpTSAMassResistivityFunction(self) -> str:
1237 """Return the GetDP function name of the material's mass resistivity."""
1238 if self.name not in getdpTSAMassResistivityFunctions:
1239 return "NOT_DEFINED_IN_DATA_FIQUS_PANCAKE3D"
1241 return getdpTSAMassResistivityFunctions[self.name]
1243 @computed_field
1244 @cached_property
1245 def getdpTSAStiffnessResistivityFunction(self) -> str:
1246 """Return the GetDP function name of the material's stiffness resistivity."""
1247 if self.name not in getdpTSAStiffnessResistivityFunctions:
1248 return "NOT_DEFINED_IN_DATA_FIQUS_PANCAKE3D"
1250 return getdpTSAStiffnessResistivityFunctions[self.name]
1252 @computed_field
1253 @cached_property
1254 def getdpTSAMassThermalConductivityFunction(self) -> str:
1255 """Return the GetDP function name of the material's mass thermal conductivity."""
1256 if self.name not in getdpTSAMassThermalConductivityFunctions:
1257 return "NOT_DEFINED_IN_DATA_FIQUS_PANCAKE3D"
1259 return getdpTSAMassThermalConductivityFunctions[self.name]
1261 @computed_field
1262 @cached_property
1263 def getdpTSAStiffnessThermalConductivityFunction(self) -> str:
1264 """Return the GetDP function name of the material's stiffness thermal conductivity."""
1265 if self.name not in getdpTSAStiffnessThermalConductivityFunctions:
1266 return "NOT_DEFINED_IN_DATA_FIQUS_PANCAKE3D"
1268 return getdpTSAStiffnessThermalConductivityFunctions[self.name]
1270 @computed_field
1271 @cached_property
1272 def getdpTSAMassHeatCapacityFunction(self) -> str:
1273 """Return the GetDP function name of the material's mass heat capacity."""
1274 if self.name not in getdpTSAMassHeatCapacityFunctions:
1275 return "NOT_DEFINED_IN_DATA_FIQUS_PANCAKE3D"
1277 return getdpTSAMassHeatCapacityFunctions[self.name]
1279 @computed_field
1280 @cached_property
1281 def getdpTSARHSFunction(self) -> str:
1282 """Return the GetDP function name of the material's RHS."""
1283 if self.name not in getdpTSARHSFunctions:
1284 return "NOT_DEFINED_IN_DATA_FIQUS_PANCAKE3D"
1286 return getdpTSARHSFunctions[self.name]
1288 @computed_field
1289 @cached_property
1290 def getdpTSATripleFunction(self) -> str:
1291 """Return the GetDP function name of the material's triple."""
1292 if self.name not in getdpTSATripleFunctions:
1293 return "NOT_DEFINED_IN_DATA_FIQUS_PANCAKE3D"
1295 return getdpTSATripleFunctions[self.name]
1298class Pancake3DSolveNormalMaterial(Pancake3DSolveMaterialBase):
1299 # Mandatory:
1300 name: NormalMaterialName = Field(
1301 title="Material Name",
1302 )
1305class Pancake3DSolveSuperconductingMaterial(Pancake3DSolveMaterialBase):
1306 # Mandatory:
1307 name: SuperconductingMaterialName = Field(
1308 title="Superconduncting Material Name",
1309 )
1310 nValue: PositiveFloat = Field(
1311 default=30,
1312 alias="N-Value for E-J Power Law",
1313 description="N-value for E-J power law.",
1314 )
1315 IcAtTinit: PositiveFloat | str | Pancake3DSolveIcVsLength = Field(
1316 alias="criticalCurrentAtInitialTemperature",
1317 title="Critical Current at Initial Temperature",
1318 description=(
1319 "Critical current at initial temperature. The critical current value will"
1320 " change with temperature depending on the superconductor material.\nEither"
1321 " the same critical current for the whole tape or the critical current with"
1322 " respect to the tape length can be specified. To specify the same critical"
1323 " current for the entire tape, just use a scalar. To specify critical"
1324 " current with respect to the tape length: a CSV file can be used, or"
1325 " lengthValues and criticalCurrentValues can be given as lists. The data"
1326 " will be linearly interpolated.\nIf a CSV file is to be used, the input"
1327 " should be the name of a CSV file (which is in the same folder as the"
1328 " input file) instead of a scalar. The first column of the CSV file will be"
1329 " the tape length, and the second column will be the critical current."
1330 ),
1331 examples=[230, "IcVSlength.csv"],
1332 )
1334 # Optionals:
1335 electricFieldCriterion: PositiveFloat = Field(
1336 default=1e-4,
1337 title="Electric Field Criterion",
1338 description=(
1339 "The electric field that defines the critical current density, i.e., the"
1340 " electric field at which the current density reaches the critical current"
1341 " density."
1342 ),
1343 )
1344 jCriticalScalingNormalToWinding: PositiveFloat = Field(
1345 default=1,
1346 title="Critical Current Scaling Normal to Winding",
1347 description=(
1348 "Critical current scaling normal to winding, i.e., along the c_axis. "
1349 " We have Jc_cAxis = scalingFactor * Jc_abPlane."
1350 " A factor of 1 means no scaling such that the HTS layer is isotropic."
1351 ),
1352 )
1353 minimumPossibleResistivity: NonNegativeFloat = Field(
1354 default=0,
1355 title="Minimum Possible Resistivity",
1356 description=(
1357 "The resistivity of the winding won't be lower than this value, no matter"
1358 " what."
1359 ),
1360 )
1361 maximumPossibleResistivity: PositiveFloat = Field(
1362 default=1,
1363 title="Maximum Possible Resistivity",
1364 description=(
1365 "The resistivity of the winding won't be higher than this value, no matter"
1366 " what."
1367 ),
1368 )
1370 @computed_field
1371 @cached_property
1372 def IcValues(self) -> list[float]:
1373 """Return the critical current values of the material."""
1374 if hasattr(self.IcAtTinit, "criticalCurrentValues"):
1375 return self.IcAtTinit.criticalCurrentValues
1376 elif isinstance(self.IcAtTinit, str):
1377 csv_file_path = pathlib.Path(input_file_path.get()).parent / self.IcAtTinit
1378 # return the second column:
1379 IcValues = list(pd.read_csv(csv_file_path, header=None).iloc[:, 1])
1380 for Ic in IcValues:
1381 if Ic < 0:
1382 raise ValueError(
1383 "Critical current values in the CSV file should be positive!"
1384 )
1385 return IcValues
1386 else:
1387 return [self.IcAtTinit]
1389 @computed_field
1390 @cached_property
1391 def lengthValues(self) -> list[float]:
1392 """Return the length values of the material."""
1393 if hasattr(self.IcAtTinit, "lengthValues"):
1394 return self.IcAtTinit.lengthValues
1395 elif isinstance(self.IcAtTinit, str):
1396 csv_file_path = pathlib.Path(input_file_path.get()).parent / self.IcAtTinit
1397 # return the first column:
1398 lengthValues = list(pd.read_csv(csv_file_path, header=None).iloc[:, 0])
1399 for length in lengthValues:
1400 if length < 0:
1401 raise ValueError("Tape lengths in the CSV file should be positive!")
1402 return lengthValues
1403 else:
1404 return [1]
1406 @computed_field
1407 @cached_property
1408 def getdpCriticalCurrentDensityFunction(self) -> str:
1409 """Return the GetDP function name of the material's critical current density."""
1410 if self.name not in getdpCriticalCurrentDensityFunctions:
1411 raise ValueError(
1412 f"Critical current density of the material '{self.name}' is not defined"
1413 " in FiQuS!"
1414 )
1416 return getdpCriticalCurrentDensityFunctions[self.name]
1419class Pancake3DSolveHTSMaterialBase(BaseModel):
1420 relativeThickness: float = Field(
1421 le=1,
1422 title="Relative Thickness (only for winding)",
1423 description=(
1424 "Winding tapes generally consist of more than one material. Therefore, when"
1425 " materials are given as a list in winding, their relative thickness,"
1426 " (thickness of the material) / (thickness of the winding), should be"
1427 " specified."
1428 ),
1429 )
1432class Pancake3DSolveHTSNormalMaterial(
1433 Pancake3DSolveHTSMaterialBase, Pancake3DSolveNormalMaterial
1434):
1435 pass
1438class Pancake3DSolveHTSSuperconductingMaterial(
1439 Pancake3DSolveHTSMaterialBase, Pancake3DSolveSuperconductingMaterial
1440):
1441 pass
1444class Pancake3DSolveHTSShuntLayerMaterial(Pancake3DSolveNormalMaterial):
1445 name: NormalMaterialName = Field(
1446 default="Copper",
1447 title="Material Name",
1448 )
1449 relativeHeight: float = Field(
1450 default=0.0,
1451 ge=0,
1452 le=1,
1453 title="Relative Height of the Shunt Layer",
1454 description=(
1455 "HTS 2G coated conductor are typically plated, usually "
1456 " using copper. The relative height of the shunt layer is the "
1457 " width of the shunt layer divided by the width of the tape. "
1458 " 0 means no shunt layer."
1459 ),
1460 )
1463class Pancake3DSolveMaterial(BaseModel):
1464 # 1) User inputs:
1466 # Mandatory:
1468 # Optionals:
1469 resistivity: Optional[PositiveFloat] = Field(
1470 default=None,
1471 title="Resistivity",
1472 description=(
1473 "A scalar value. If this is given, material properties won't be used for"
1474 " resistivity."
1475 ),
1476 )
1477 thermalConductivity: Optional[PositiveFloat] = Field(
1478 default=None,
1479 title="Thermal Conductivity",
1480 description=(
1481 "A scalar value. If this is given, material properties won't be used for"
1482 " thermal conductivity."
1483 ),
1484 )
1485 specificHeatCapacity: Optional[PositiveFloat] = Field(
1486 default=None,
1487 title="Specific Heat Capacity",
1488 description=(
1489 "A scalar value. If this is given, material properties won't be used for"
1490 " specific heat capacity."
1491 ),
1492 )
1493 material: Optional[Pancake3DSolveNormalMaterial] = Field(
1494 default=None,
1495 title="Material",
1496 description="Material from STEAM material library.",
1497 )
1499 @model_validator(mode="after")
1500 @classmethod
1501 def check_material(cls, model: "Pancake3DSolveMaterial"):
1502 if model.material is None:
1503 if model.resistivity is None:
1504 raise ValueError(
1505 "Resistivity of the material is not given, and no material is"
1506 " provided!"
1507 )
1508 if model.thermalConductivity is None:
1509 raise ValueError(
1510 "Thermal conductivity of the material is not given, and no material"
1511 " is provided!"
1512 )
1513 if model.specificHeatCapacity is None:
1514 raise ValueError(
1515 "Specific heat capacity of the material is not given, and no"
1516 " material is provided!"
1517 )
1519 return model
1522class Pancake3DSolveShuntLayerMaterial(Pancake3DSolveMaterial):
1523 material: Optional[Pancake3DSolveHTSShuntLayerMaterial] = Field(
1524 default=Pancake3DSolveHTSShuntLayerMaterial(),
1525 title="Material",
1526 description="Material from STEAM material library.",
1527 )
1530class Pancake3DSolveContactLayerMaterial(Pancake3DSolveMaterial):
1531 resistivity: PositiveFloat | Literal["perfectlyInsulating"] = Field(
1532 default=None,
1533 title="Resistivity",
1534 description=(
1535 'A scalar value or "perfectlyInsulating". If "perfectlyInsulating" is'
1536 " given, the contact layer will be perfectly insulating. If this value is"
1537 " given, material properties won't be used for resistivity."
1538 ),
1539 )
1540 numberOfThinShellElements: PositiveInt = Field(
1541 default=1,
1542 title="Number of Thin Shell Elements (Advanced Input)",
1543 description=(
1544 "Number of thin shell elements in the FE formulation (GetDP related, not"
1545 " physical and only used when TSA is set to True)"
1546 ),
1547 )
1549 @field_validator("resistivity")
1550 @classmethod
1551 def checkPerfectlyInsulatingCase(cls, resistivity):
1552 if resistivity == "perfectlyInsulating":
1553 geometry = geometry_input.get()
1554 if geometry["numberOfPancakes"] > 1:
1555 raise ValueError(
1556 "Contact layer can't be perfectly insulating for multi-pancake"
1557 " coils!"
1558 )
1560 return resistivity
1563Pancake3DHTSMaterial = Annotated[
1564 Pancake3DSolveHTSNormalMaterial | Pancake3DSolveHTSSuperconductingMaterial,
1565 Field(discriminator="name"),
1566]
1569class Pancake3DSolveWindingMaterial(Pancake3DSolveMaterial):
1570 material: Optional[list[Pancake3DHTSMaterial]] = Field(
1571 default=None,
1572 title="Materials of HTS CC",
1573 description="List of materials of HTS CC.",
1574 )
1576 shuntLayer: Pancake3DSolveShuntLayerMaterial = Field(
1577 default=Pancake3DSolveShuntLayerMaterial(),
1578 title="Shunt Layer Properties",
1579 description="Material properties of the shunt layer.",
1580 )
1582 @field_validator("material")
1583 @classmethod
1584 def checkIfRelativeThicknessesSumToOne(cls, material):
1585 if material is not None:
1586 totalRelativeThickness = sum(
1587 material.relativeThickness for material in material
1588 )
1589 if not math.isclose(totalRelativeThickness, 1, rel_tol=1e-3):
1590 raise ValueError(
1591 "Sum of relative thicknesses of HTS materials should be 1!"
1592 )
1594 return material
1596 @computed_field
1597 @cached_property
1598 def relativeThicknessOfNormalConductor(self) -> PositiveFloat:
1599 """Return the relative thickness of the normal conductor."""
1600 if self.material is None:
1601 return 0
1602 else:
1603 # look at normal materials inside self.material and sum their relativeThickness
1604 return sum(
1605 material.relativeThickness
1606 for material in self.material
1607 if isinstance(material, Pancake3DSolveHTSNormalMaterial)
1608 )
1610 @computed_field
1611 @cached_property
1612 def relativeThicknessOfSuperConductor(self) -> PositiveFloat:
1613 """Return the relative thickness of the super conductor."""
1614 if self.material is None:
1615 return 0
1616 else:
1617 # look at superconducting materials inside self.material and sum their relativeThickness
1618 return sum(
1619 material.relativeThickness
1620 for material in self.material
1621 if isinstance(material, Pancake3DSolveHTSSuperconductingMaterial)
1622 )
1624 @computed_field
1625 @cached_property
1626 def normalConductors(self) -> list[Pancake3DSolveHTSNormalMaterial]:
1627 """Return the normal conductors of the winding."""
1628 if self.material is None:
1629 return []
1630 else:
1631 return [
1632 material
1633 for material in self.material
1634 if isinstance(material, Pancake3DSolveHTSNormalMaterial)
1635 ]
1637 @computed_field
1638 @cached_property
1639 def superConductor(self) -> Pancake3DSolveHTSSuperconductingMaterial:
1640 """Return the super conductor of the winding."""
1641 if self.material is None:
1642 return None
1643 else:
1644 superConductors = [
1645 material
1646 for material in self.material
1647 if isinstance(material, Pancake3DSolveHTSSuperconductingMaterial)
1648 ]
1649 if len(superConductors) == 0:
1650 return None
1651 elif len(superConductors) == 1:
1652 return superConductors[0]
1653 else:
1654 raise ValueError(
1655 "There should be only one superconductor in the winding!"
1656 )
1659class Pancake3DSolveTerminalMaterialAndBoundaryCondition(Pancake3DSolveMaterial):
1660 cooling: Literal["adiabatic", "fixedTemperature", "cryocooler"] = Field(
1661 default="fixedTemperature",
1662 title="Cooling condition",
1663 description=(
1664 "Cooling condition of the terminal. It can be either adiabatic, fixed"
1665 " temperature, or cryocooler."
1666 ),
1667 )
1668 transitionNotch: Pancake3DSolveMaterial = Field(
1669 title="Transition Notch Properties",
1670 description="Material properties of the transition notch volume.",
1671 )
1672 terminalContactLayer: Pancake3DSolveMaterial = Field(
1673 title="Transition Layer Properties",
1674 description=(
1675 "Material properties of the transition layer between terminals and"
1676 " windings."
1677 ),
1678 )
1681class Pancake3DSolveToleranceBase(BaseModel):
1682 # Mandatory:
1683 quantity: str
1684 relative: NonNegativeFloat = Field(
1685 title="Relative Tolerance",
1686 description="Relative tolerance for the quantity.",
1687 )
1688 absolute: NonNegativeFloat = Field(
1689 title="Absolute Tolerance", description="Absolute tolerance for the quantity"
1690 )
1692 # Optionals:
1693 normType: Literal["L1Norm", "MeanL1Norm", "L2Norm", "MeanL2Norm", "LinfNorm"] = (
1694 Field(
1695 default="L2Norm",
1696 alias="normType",
1697 title="Norm Type",
1698 description=(
1699 "Sometimes, tolerances return a vector instead of a scalar (ex,"
1700 " solutionVector). Then, the magnitude of the tolerance should be"
1701 " calculated with a method. Norm type selects this method."
1702 ),
1703 )
1704 )
1707class Pancake3DSolvePositionRequiredTolerance(Pancake3DSolveToleranceBase):
1708 # Mandatory:
1709 quantity: PositionRequiredQuantityName = Field(
1710 title="Quantity", description="Name of the quantity for tolerance."
1711 )
1712 position: Pancake3DPosition = Field(
1713 title="Probing Position of the Quantity",
1714 description="Probing position of the quantity for tolerance.",
1715 )
1718class Pancake3DSolvePositionNotRequiredTolerance(Pancake3DSolveToleranceBase):
1719 # Mandatory:
1720 quantity: (
1721 Literal[
1722 "solutionVector",
1723 ]
1724 | PositionNotRequiredQuantityName
1725 ) = Field(
1726 title="Quantity",
1727 description="Name of the quantity for tolerance.",
1728 )
1731Pancake3DSolveTolerance = Annotated[
1732 Pancake3DSolvePositionRequiredTolerance
1733 | Pancake3DSolvePositionNotRequiredTolerance,
1734 Field(discriminator="quantity"),
1735]
1738class Pancake3DSolveSettingsWithTolerances(BaseModel):
1739 tolerances: list[Pancake3DSolveTolerance] = Field(
1740 title="Tolerances for Adaptive Time Stepping",
1741 description=(
1742 "Time steps or nonlinear iterations will be refined until the tolerances"
1743 " are satisfied."
1744 ),
1745 )
1747 @computed_field
1748 @cached_property
1749 def postOperationTolerances(self) -> list[Pancake3DSolveTolerance]:
1750 """Return the post operation tolerances."""
1751 tolerances = [
1752 tolerance
1753 for tolerance in self.tolerances
1754 if "solutionVector" not in tolerance.quantity
1755 ]
1756 return tolerances
1758 @computed_field
1759 @cached_property
1760 def systemTolerances(self) -> list[Pancake3DSolveTolerance]:
1761 """Return the system tolerances."""
1762 tolerances = [
1763 tolerance
1764 for tolerance in self.tolerances
1765 if "solutionVector" in tolerance.quantity
1766 ]
1767 return tolerances
1770class Pancake3DSolveAdaptiveTimeLoopSettings(Pancake3DSolveSettingsWithTolerances):
1771 # Mandatory:
1772 initialStep: PositiveFloat = Field(
1773 alias="initialStep",
1774 title="Initial Step for Adaptive Time Stepping",
1775 description="Initial step for adaptive time stepping",
1776 )
1777 minimumStep: PositiveFloat = Field(
1778 alias="minimumStep",
1779 title="Minimum Step for Adaptive Time Stepping",
1780 description=(
1781 "The simulation will be aborted if a finer time step is required than this"
1782 " minimum step value."
1783 ),
1784 )
1785 maximumStep: PositiveFloat = Field(
1786 alias="maximumStep",
1787 description="Bigger steps than this won't be allowed",
1788 )
1790 # Optionals:
1791 integrationMethod: Literal[
1792 "Euler", "Gear_2", "Gear_3", "Gear_4", "Gear_5", "Gear_6"
1793 ] = Field(
1794 default="Euler",
1795 alias="integrationMethod",
1796 title="Integration Method",
1797 description="Integration method for transient analysis",
1798 )
1799 breakPoints_input: list[float] = Field(
1800 default=[0],
1801 alias="breakPoints",
1802 title="Break Points for Adaptive Time Stepping",
1803 description="Make sure to solve the system for these times.",
1804 )
1806 @field_validator("breakPoints_input")
1807 @classmethod
1808 def updateGlobalBreakPointsList(cls, breakPoints_input):
1809 all_break_points.extend(breakPoints_input)
1810 return breakPoints_input
1812 @model_validator(mode="after")
1813 @classmethod
1814 def check_time_steps(cls, model: "Pancake3DSolveAdaptiveTimeLoopSettings"):
1815 """
1816 Checks if the time steps are consistent.
1818 :param values: dictionary of time steps
1819 :type values: dict
1820 :return: dictionary of time steps
1821 :rtype: dict
1822 """
1824 if model.initialStep < model.minimumStep:
1825 raise ValueError(
1826 "Initial time step cannot be smaller than the minimum time step!"
1827 )
1829 if model.initialStep > model.maximumStep:
1830 raise ValueError(
1831 "Initial time step cannot be greater than the maximum time step!"
1832 )
1834 if model.minimumStep > model.maximumStep:
1835 raise ValueError(
1836 "Minimum time step cannot be greater than the maximum time step!"
1837 )
1839 return model
1841 @computed_field
1842 @cached_property
1843 def breakPoints(self) -> list[float]:
1844 """Return the break points for adaptive time stepping."""
1845 breakPoints = list(set(all_break_points))
1846 breakPoints.sort()
1847 return breakPoints
1850class Pancake3DSolveFixedTimeLoopSettings(BaseModel):
1851 # Mandatory:
1852 step: PositiveFloat = Field(
1853 title="Step for Fixed Time Stepping",
1854 description="Time step for fixed time stepping.",
1855 )
1858class Pancake3DSolveFixedLoopInterval(BaseModel):
1859 # Mandatory:
1860 startTime: NonNegativeFloat = Field(
1861 title="Start Time of the Interval",
1862 description="Start time of the interval.",
1863 )
1864 endTime: NonNegativeFloat = Field(
1865 title="End Time of the Interval",
1866 description="End time of the interval.",
1867 )
1868 step: PositiveFloat = Field(
1869 title="Step for the Interval",
1870 description="Time step for the interval",
1871 )
1873 @model_validator(mode="after")
1874 @classmethod
1875 def check_time_steps(cls, model: "Pancake3DSolveFixedLoopInterval"):
1876 """
1877 Checks if the time steps are consistent.
1879 :param values: dictionary of time steps
1880 :type values: dict
1881 :return: dictionary of time steps
1882 :rtype: dict
1883 """
1885 if model.startTime > model.endTime:
1886 raise ValueError("Start time cannot be greater than the end time!")
1888 interval = model.endTime - model.startTime
1889 if (
1890 math.fmod(interval, model.step) > 1e-8
1891 and math.fmod(interval, model.step) - model.step < -1e-8
1892 ):
1893 raise ValueError("Time interval must be a multiple of the time step!")
1895 return model
1898class Pancake3DSolveTimeBase(BaseModel):
1899 # Mandatory:
1900 start: float = Field(
1901 title="Start Time", description="Start time of the simulation."
1902 )
1903 end: float = Field(title="End Time", description="End time of the simulation.")
1905 # Optionals:
1906 extrapolationOrder: Literal[0, 1, 2, 3] = Field(
1907 default=1,
1908 alias="extrapolationOrder",
1909 title="Extrapolation Order",
1910 description=(
1911 "Before solving for the next time steps, the previous solutions can be"
1912 " extrapolated for better convergence."
1913 ),
1914 )
1916 @model_validator(mode="after")
1917 @classmethod
1918 def check_time_steps(cls, model: "Pancake3DSolveTimeBase"):
1919 """
1920 Checks if the time steps are consistent.
1921 """
1923 if model.start > model.end:
1924 raise ValueError("Start time cannot be greater than the end time!")
1926 return model
1929class Pancake3DSolveTimeAdaptive(Pancake3DSolveTimeBase):
1930 timeSteppingType: Literal["adaptive"] = "adaptive"
1931 adaptive: Pancake3DSolveAdaptiveTimeLoopSettings = Field(
1932 alias="adaptiveSteppingSettings",
1933 title="Adaptive Time Loop Settings",
1934 description=(
1935 "Adaptive time loop settings (only used if stepping type is adaptive)."
1936 ),
1937 )
1939 @model_validator(mode="after")
1940 @classmethod
1941 def check_time_steps(cls, model: "Pancake3DSolveTimeAdaptive"):
1942 """
1943 Checks if the time steps are consistent.
1944 """
1945 if model.adaptive.initialStep > model.end:
1946 raise ValueError("Initial time step cannot be greater than the end time!")
1948 return model
1951class Pancake3DSolveTimeFixed(Pancake3DSolveTimeBase):
1952 timeSteppingType: Literal["fixed"] = "fixed"
1953 fixed: (
1954 list[Pancake3DSolveFixedLoopInterval] | Pancake3DSolveFixedTimeLoopSettings
1955 ) = Field(
1956 alias="fixedSteppingSettings",
1957 title="Fixed Time Loop Settings",
1958 description="Fixed time loop settings (only used if stepping type is fixed).",
1959 )
1961 @model_validator(mode="after")
1962 @classmethod
1963 def check_time_steps(cls, model: "Pancake3DSolveTimeFixed"):
1964 if isinstance(model.fixed, list):
1965 for i in range(len(model.fixed) - 1):
1966 if model.fixed[i].endTime != model.fixed[i + 1].startTime:
1967 raise ValueError(
1968 "End time of the previous interval must be equal to the start"
1969 " time of the next interval!"
1970 )
1972 if model.fixed[0].startTime != model.start:
1973 raise ValueError(
1974 "Start time of the first interval must be equal to the start time"
1975 " of the simulation!"
1976 )
1978 else:
1979 if model.fixed.step > model.end:
1980 raise ValueError("Time step cannot be greater than the end time!")
1982 if not (
1983 math.isclose(
1984 (model.end - model.start) % model.fixed.step, 0, abs_tol=1e-8
1985 )
1986 ):
1987 raise ValueError("Time interval must be a multiple of the time step!")
1990 return model
1993Pancake3DSolveTime = Annotated[
1994 Pancake3DSolveTimeAdaptive | Pancake3DSolveTimeFixed,
1995 Field(discriminator="timeSteppingType"),
1996]
1999class Pancake3DSolveNonlinearSolverSettings(Pancake3DSolveSettingsWithTolerances):
2000 # Optionals:
2001 maxIter: PositiveInt = Field(
2002 default=100,
2003 alias="maximumNumberOfIterations",
2004 title="Maximum Number of Iterations",
2005 description="Maximum number of iterations allowed for the nonlinear solver.",
2006 )
2007 relaxationFactor: float = Field(
2008 default=0.7,
2009 gt=0,
2010 alias="relaxationFactor",
2011 title="Relaxation Factor",
2012 description=(
2013 "Calculated step changes of the solution vector will be multiplied with"
2014 " this value for better convergence."
2015 ),
2016 )
2019class Pancake3DSolveInitialConditions(BaseModel):
2020 # 1) User inputs:
2022 # Mandatory:
2023 T: PositiveFloat = Field(
2024 alias="temperature",
2025 title="Initial Temperature",
2026 description="Initial temperature of the pancake coils.",
2027 )
2030class Pancake3DSolveLocalDefect(BaseModel):
2031 # Mandatory:
2032 value: NonNegativeFloat = Field(
2033 alias="value",
2034 title="Value",
2035 description="Value of the local defect.",
2036 )
2037 startTurn: NonNegativeFloat = Field(
2038 alias="startTurn",
2039 title="Start Turn",
2040 description="Start turn of the local defect.",
2041 )
2042 endTurn: PositiveFloat = Field(
2043 alias="endTurn",
2044 title="End Turn",
2045 description="End turn of the local defect.",
2046 )
2048 startTime: NonNegativeFloat = Field(
2049 alias="startTime",
2050 title="Start Time",
2051 description="Start time of the local defect.",
2052 )
2054 # Optionals:
2055 transitionDuration: NonNegativeFloat = Field(
2056 default=0,
2057 title="Transition Duration",
2058 description=(
2059 "Transition duration of the local defect. If not given, the transition will"
2060 " be instantly."
2061 ),
2062 )
2063 whichPancakeCoil: Optional[PositiveInt] = Field(
2064 default=None,
2065 title="Pancake Coil Number",
2066 description="The first pancake coil is 1, the second is 2, etc.",
2067 )
2069 @field_validator("startTime")
2070 @classmethod
2071 def updateGlobalBreakPointsList(cls, startTime):
2072 all_break_points.append(startTime)
2073 return startTime
2075 @field_validator("startTime")
2076 @classmethod
2077 def check_start_time(cls, startTime):
2078 solve = solve_input.get()
2079 start_time = solve["time"]["start"]
2080 end_time = solve["time"]["end"]
2082 if startTime < start_time or startTime > end_time:
2083 raise ValueError(
2084 "Start time of the local defect should be between the start and end"
2085 " times!"
2086 )
2088 return startTime
2090 @field_validator("endTurn")
2091 @classmethod
2092 def check_end_turn(cls, endTurn, info: ValidationInfo):
2093 geometry = geometry_input.get()
2095 if endTurn > geometry["winding"]["numberOfTurns"]:
2096 raise ValueError(
2097 "End turn of the local defect can't be greater than the number of"
2098 " turns!"
2099 )
2101 if endTurn < info.data["startTurn"]:
2102 raise ValueError(
2103 "End turn of the local defect can't be smaller than the start turn!"
2104 )
2106 return endTurn
2108 @field_validator("whichPancakeCoil")
2109 @classmethod
2110 def check_which_pancake_coil(cls, whichPancakeCoil):
2111 geometry = geometry_input.get()
2113 if whichPancakeCoil is None:
2114 if geometry["numberOfPancakes"] == 1:
2115 whichPancakeCoil = 1
2116 else:
2117 raise ValueError(
2118 "whickPancakeCoil (pancake coil number) should be given if there"
2119 " are more than one pancake coil!"
2120 )
2122 return whichPancakeCoil
2124 @computed_field
2125 @cached_property
2126 def zTop(self) -> float:
2127 """Return the z-coordinate of the top of the pancake coil."""
2128 geometry = geometry_input.get()
2130 zTop = self.zBottom + geometry["winding"]["height"]
2132 return zTop
2134 @computed_field
2135 @cached_property
2136 def zBottom(self) -> float:
2137 """Return the z-coordinate of the bottom of the pancake coil."""
2138 geometry = geometry_input.get()
2140 if not self.whichPancakeCoil:
2141 self.whichPancakeCoil = 1
2143 zBottom = -(
2144 geometry["numberOfPancakes"] * geometry["winding"]["height"]
2145 + (geometry["numberOfPancakes"] - 1)
2146 * geometry["gapBetweenPancakes"]
2147 ) / 2 + (self.whichPancakeCoil - 1) * (
2148 geometry["winding"]["height"] + geometry["gapBetweenPancakes"]
2149 )
2151 return zBottom
2154class Pancake3DSolveLocalDefects(BaseModel):
2155 # 1) User inputs:
2157 jCritical: Optional[Pancake3DSolveLocalDefect] = Field(
2158 default=None,
2159 alias="criticalCurrentDensity",
2160 title="Local Defect for Critical Current Density",
2161 description="Set critical current density locally.",
2162 )
2164 @field_validator("jCritical")
2165 @classmethod
2166 def check_jCritical_local_defect(cls, jCritical):
2167 if jCritical is not None:
2168 solve = solve_input.get()
2170 if "material" in solve["winding"]:
2171 windingMaterials = [
2172 material["name"] for material in solve["winding"]["material"]
2173 ]
2174 else:
2175 windingMaterials = []
2177 superconducting_material_is_used = "HTSSuperPower" in windingMaterials
2179 if not superconducting_material_is_used:
2180 raise ValueError(
2181 "Local defects can only be set if superconducting material is used!"
2182 )
2184 return jCritical
2187class Pancake3DSolveQuantityBase(BaseModel):
2188 # Mandatory:
2189 quantity: PositionNotRequiredQuantityName | PositionRequiredQuantityName = Field(
2190 title="Quantity",
2191 description="Name of the quantity to be saved.",
2192 )
2194 @computed_field
2195 @cached_property
2196 def getdpQuantityName(self) -> str:
2197 """Return the GetDP name of the quantity."""
2198 if self.quantity not in getdpQuantityNames:
2199 raise ValueError(f"Quantity '{self.quantity}' is not defined in FiQuS!")
2201 return getdpQuantityNames[self.quantity]
2203 @computed_field
2204 @cached_property
2205 def getdpPostOperationName(self) -> str:
2206 """Return the GetDP name of the post operation."""
2207 if self.quantity not in getdpPostOperationNames:
2208 raise ValueError(f"Quantity '{self.quantity}' is not defined in FiQuS!")
2210 return getdpPostOperationNames[self.quantity]
2213class Pancake3DSolveSaveQuantity(Pancake3DSolveQuantityBase):
2214 # Optionals:
2215 timesToBeSaved: list[float] = Field(
2216 default=[],
2217 title="Times to be Saved",
2218 description=(
2219 "List of times that wanted to be saved. If not given, all the time steps"
2220 " will be saved."
2221 ),
2222 )
2224 @field_validator("timesToBeSaved")
2225 @classmethod
2226 def updateGlobalBreakPointsList(cls, timesToBeSaved):
2227 all_break_points.extend(timesToBeSaved)
2228 return timesToBeSaved
2230 @field_validator("timesToBeSaved")
2231 @classmethod
2232 def check_times_to_be_saved(cls, timesToBeSaved):
2233 solve = solve_input.get()
2234 start_time = solve["time"]["start"]
2235 end_time = solve["time"]["end"]
2237 for time in timesToBeSaved:
2238 if time < start_time or time > end_time:
2239 raise ValueError(
2240 "Times to be saved should be between the start and end times!"
2241 )
2243 return timesToBeSaved
2246# ======================================================================================
2247# SOLVE CLASSES ENDS ===================================================================
2248# ======================================================================================
2250# ======================================================================================
2251# POSTPROCESS CLASSES STARTS ===========================================================
2252# ======================================================================================
2255class Pancake3DPostprocessTimeSeriesPlotBase(Pancake3DSolveQuantityBase):
2256 # Mandatory:
2257 quantity: str
2259 @computed_field
2260 @cached_property
2261 def fileName(self) -> str:
2262 """Return the name of the file to be saved."""
2263 return f"{self.quantity}(TimeSeriesPlot)"
2265 @computed_field
2266 @cached_property
2267 def quantityProperName(self) -> str:
2268 """Return the proper name of the quantity."""
2269 if self.quantity not in quantityProperNames:
2270 raise ValueError(
2271 f"Quantity '{self.quantity}'s proper name is not defined! Please"
2272 ' update the dictionary "quantityProperNames" in the file'
2273 ' "fiqus/data/DataFiQuSPancake3D.py".'
2274 )
2276 return quantityProperNames[self.quantity]
2278 @computed_field
2279 @cached_property
2280 def units(self) -> str:
2281 """Return the units of the quantity."""
2282 if self.quantity not in quantityUnits:
2283 raise ValueError(
2284 f"Quantity '{self.quantity}'s units are not defined! Please update"
2285 ' the dictionary "quantityUnits" in the file'
2286 ' "fiqus/data/DataFiQuSPancake3D.py".'
2287 )
2289 return quantityUnits[self.quantity]
2292class Pancake3DPostprocessTimeSeriesPlotPositionRequired(
2293 Pancake3DPostprocessTimeSeriesPlotBase
2294):
2295 # Mandatory:
2296 quantity: PositionRequiredQuantityName = Field(
2297 title="Quantity",
2298 description="Name of the quantity to be plotted.",
2299 )
2301 position: Pancake3DPosition = Field(
2302 title="Probing Position",
2303 description="Probing position of the quantity for time series plot.",
2304 )
2307class Pancake3DPostprocessTimeSeriesPlotPositionNotRequired(
2308 Pancake3DPostprocessTimeSeriesPlotBase
2309):
2310 # Mandatory:
2311 quantity: PositionNotRequiredQuantityName = Field(
2312 title="Quantity",
2313 description="Name of the quantity to be plotted.",
2314 )
2317Pancake3DPostprocessTimeSeriesPlot = Annotated[
2318 Pancake3DPostprocessTimeSeriesPlotPositionRequired
2319 | Pancake3DPostprocessTimeSeriesPlotPositionNotRequired,
2320 Field(discriminator="quantity"),
2321]
2324class Pancake3DPostprocessMagneticFieldOnPlane(BaseModel):
2325 # Optional:
2326 colormap: str = Field(
2327 default="viridis",
2328 title="Colormap",
2329 description="Colormap for the plot.",
2330 )
2331 streamLines: bool = Field(
2332 default=True,
2333 title="Stream Lines",
2334 description=(
2335 "If True, streamlines will be plotted. Note that magnetic field vectors may"
2336 " have components perpendicular to the plane, and streamlines will be drawn"
2337 " depending on the vectors' projection onto the plane."
2338 ),
2339 )
2340 interpolationMethod: Literal["nearest", "linear", "cubic"] = Field(
2341 default="linear",
2342 title="Interpolation Method",
2343 description=(
2344 "Interpolation type for the plot.\nBecause of the FEM basis function"
2345 " selections of FiQuS, each mesh element has a constant magnetic field"
2346 " vector. Therefore, for smooth 2D plots, interpolation can be"
2347 " used.\nTypes:\nnearest: it will plot the nearest magnetic field value to"
2348 " the plotting point.\nlinear: it will do linear interpolation to the"
2349 " magnetic field values.\ncubic: it will do cubic interpolation to the"
2350 " magnetic field values."
2351 ),
2352 )
2353 timesToBePlotted: Optional[list[float]] = Field(
2354 default=None,
2355 title="Times to be Plotted",
2356 description=(
2357 "List of times that wanted to be plotted. If not given, all the time steps"
2358 " will be plotted."
2359 ),
2360 )
2361 planeNormal: Annotated[list[float], Len(min_length=3, max_length=3)] = Field(
2362 default=[1, 0, 0],
2363 title="Plane Normal",
2364 description="Normal vector of the plane. The default is YZ-plane (1, 0, 0).",
2365 )
2366 planeXAxisUnitVector: Annotated[list[float], Len(min_length=3, max_length=3)] = (
2367 Field(
2368 default=[0, 1, 0],
2369 title="Plane X Axis",
2370 description=(
2371 "If an arbitrary plane is wanted to be plotted, the arbitrary plane's X"
2372 " axis unit vector must be specified. The dot product of the plane's X"
2373 " axis and the plane's normal vector must be zero."
2374 ),
2375 )
2376 )
2378 @field_validator("timesToBePlotted")
2379 @classmethod
2380 def updateGlobalBreakPointsList(cls, timesToBePlotted):
2381 all_break_points.extend(timesToBePlotted)
2382 return timesToBePlotted
2384 @field_validator("colormap")
2385 @classmethod
2386 def check_color_map(cls, colorMap):
2387 """
2388 Check if the colormap exists.
2389 """
2390 if colorMap not in matplotlib.pyplot.colormaps():
2391 raise ValueError(
2392 f"{colorMap} is not a valid colormap! Please see"
2393 " https://matplotlib.org/stable/gallery/color/colormap_reference.html"
2394 " for available colormaps."
2395 )
2397 return colorMap
2399 @computed_field
2400 @cached_property
2401 def onSection(self) -> list[list[float]]:
2402 """Return the three corner points of the plane."""
2404 class unitVector:
2405 def __init__(self, u, v, w) -> None:
2406 length = math.sqrt(u**2 + v**2 + w**2)
2407 self.u = u / length
2408 self.v = v / length
2409 self.w = w / length
2411 def rotate(self, theta, withRespectTo):
2412 # Rotate with respect to the withRespectTo vector by theta degrees:
2413 # https://en.wikipedia.org/wiki/Rotation_matrix#Rotation_matrix_from_axis_and_angle
2414 a = withRespectTo.u
2415 b = withRespectTo.v
2416 c = withRespectTo.w
2418 rotationMatrix = np.array(
2419 [
2420 [
2421 math.cos(theta) + a**2 * (1 - math.cos(theta)),
2422 a * b * (1 - math.cos(theta)) - c * math.sin(theta),
2423 a * c * (1 - math.cos(theta)) + b * math.sin(theta),
2424 ],
2425 [
2426 b * a * (1 - math.cos(theta)) + c * math.sin(theta),
2427 math.cos(theta) + b**2 * (1 - math.cos(theta)),
2428 b * c * (1 - math.cos(theta)) - a * math.sin(theta),
2429 ],
2430 [
2431 c * a * (1 - math.cos(theta)) - b * math.sin(theta),
2432 c * b * (1 - math.cos(theta)) + a * math.sin(theta),
2433 math.cos(theta) + c**2 * (1 - math.cos(theta)),
2434 ],
2435 ]
2436 )
2437 vector = np.array([[self.u], [self.v], [self.w]])
2438 rotatedVector = rotationMatrix @ vector
2439 return unitVector(
2440 rotatedVector[0][0],
2441 rotatedVector[1][0],
2442 rotatedVector[2][0],
2443 )
2445 def __pow__(self, otherUnitVector):
2446 # Cross product:
2447 u = self.v * otherUnitVector.w - self.w * otherUnitVector.v
2448 v = self.w * otherUnitVector.u - self.u * otherUnitVector.w
2449 w = self.u * otherUnitVector.v - self.v * otherUnitVector.u
2450 return unitVector(u, v, w)
2452 def __mul__(self, otherUnitVector) -> float:
2453 # Dot product:
2454 return (
2455 self.u * otherUnitVector.u
2456 + self.v * otherUnitVector.v
2457 + self.w * otherUnitVector.w
2458 )
2460 planeNormal = unitVector(
2461 self.planeNormal[0], self.planeNormal[1], self.planeNormal[2]
2462 )
2463 planeXAxis = unitVector(
2464 self.planeXAxis[0], self.planeXAxis[1], self.planeXAxis[2]
2465 )
2467 if not math.isclose(planeNormal * planeXAxis, 0, abs_tol=1e-8):
2468 raise ValueError(
2469 "planeNormal and planeXAxis must be perpendicular to each"
2470 " other! If planeNormal is chosen arbitrarily, planeXAxis must"
2471 " be specified."
2472 )
2474 # A plane that passes through the origin with the normal vector planeNormal
2475 # can be defined as:
2476 # a*x + b*y + c*z = 0
2477 a = planeNormal.u
2478 b = planeNormal.v
2479 c = planeNormal.w
2481 # Pick three points on the plane to define a rectangle. The points will
2482 # be the corners of the rectangle.
2484 # To do that, change coordinate system:
2486 # Find a vector that is perpendicular to planeNormal:
2487 randomVector = unitVector(c, a, b)
2488 perpendicularVector1 = planeNormal**randomVector
2490 # Rotate perperndicular vector with respect to the plane's normal vector
2491 # by 90 degrees to find the second perpendicular vector:
2492 perpendicularVector2 = perpendicularVector1.rotate(math.pi / 2, planeNormal)
2494 # Build the transformation matrix to change from the plane's coordinate
2495 # system to the global coordinate system:
2496 transformationMatrix = np.array(
2497 [
2498 [
2499 perpendicularVector1.u,
2500 perpendicularVector1.v,
2501 perpendicularVector1.w,
2502 ],
2503 [
2504 perpendicularVector2.u,
2505 perpendicularVector2.v,
2506 perpendicularVector2.w,
2507 ],
2508 [planeNormal.u, planeNormal.v, planeNormal.w],
2509 ]
2510 )
2511 transformationMatrix = np.linalg.inv(transformationMatrix)
2513 # Select three points to define the rectangle. Pick large points because
2514 # only the intersection of the rectangle and the mesh will be used.
2515 geometry = geometry_input.get()
2516 if geometry["air"]["type"] == "cuboid":
2517 sideLength = geometry["air"]["sideLength"]
2518 airMaxWidth = 2 * math.sqrt((sideLength / 2) ** 2 + (sideLength / 2) ** 2)
2519 if geometry["air"]["type"] == "cylinder":
2520 airMaxWidth = geometry["air"]["radius"]
2522 airHeight = getAirHeight()
2524 point1InPlaneCoordinates = np.array(
2525 [5 * max(airMaxWidth, airHeight), 5 * max(airMaxWidth, airHeight), 0]
2526 )
2527 point1 = transformationMatrix @ point1InPlaneCoordinates
2529 point2InPlaneCoordinates = np.array(
2530 [-5 * max(airMaxWidth, airHeight), 5 * max(airMaxWidth, airHeight), 0]
2531 )
2532 point2 = transformationMatrix @ point2InPlaneCoordinates
2534 point3InPlaneCoordinates = np.array(
2535 [
2536 -5 * max(airMaxWidth, airHeight),
2537 -5 * max(airMaxWidth, airHeight),
2538 0,
2539 ]
2540 )
2541 point3 = transformationMatrix @ point3InPlaneCoordinates
2543 # Round the point coordinates to the nearest multiple of the dimTol:
2544 if "dimTol" in geometry:
2545 dimTol = geometry["dimTol"]
2546 else:
2547 dimTol = 1e-8
2549 point1[0] = round(point1[0] / dimTol) * dimTol
2550 point1[1] = round(point1[1] / dimTol) * dimTol
2551 point1[2] = round(point1[2] / dimTol) * dimTol
2552 point2[0] = round(point2[0] / dimTol) * dimTol
2553 point2[1] = round(point2[1] / dimTol) * dimTol
2554 point2[2] = round(point2[2] / dimTol) * dimTol
2555 point3[0] = round(point3[0] / dimTol) * dimTol
2556 point3[1] = round(point3[1] / dimTol) * dimTol
2557 point3[2] = round(point3[2] / dimTol) * dimTol
2559 onSection = [
2560 [float(point1[0]), float(point1[1]), float(point1[2])],
2561 [float(point2[0]), float(point2[1]), float(point2[2])],
2562 [float(point3[0]), float(point3[1]), float(point3[2])],
2563 ]
2565 return onSection
2568# ======================================================================================
2569# POSTPROCESS CLASSES ENDS =============================================================
2570# ======================================================================================
2573class Pancake3DGeometry(BaseModel):
2574 # Mandatory:
2575 N: PositiveInt = Field(
2576 ge=1,
2577 alias="numberOfPancakes",
2578 title="Number of Pancakes",
2579 description="Number of pancake coils stacked on top of each other.",
2580 )
2582 gap: PositiveFloat = Field(
2583 alias="gapBetweenPancakes",
2584 title="Gap Between Pancakes",
2585 description="Gap distance between the pancake coils.",
2586 )
2588 wi: Pancake3DGeometryWinding = Field(
2589 alias="winding",
2590 title="Winding Geometry",
2591 description="This dictionary contains the winding geometry information.",
2592 )
2594 ii: Pancake3DGeometryContactLayer = Field(
2595 alias="contactLayer",
2596 title="Contact Layer Geometry",
2597 description="This dictionary contains the contact layer geometry information.",
2598 )
2600 ti: Pancake3DGeometryTerminals = Field(
2601 alias="terminals",
2602 title="Terminals Geometry",
2603 description="This dictionary contains the terminals geometry information.",
2604 )
2606 ai: Pancake3DGeometryAir = Field(
2607 alias="air",
2608 title="Air Geometry",
2609 description="This dictionary contains the air geometry information.",
2610 )
2612 # Optionals:
2613 dimTol: PositiveFloat = Field(
2614 default=1e-8,
2615 alias="dimensionTolerance",
2616 description="dimension tolerance (CAD related, not physical)",
2617 )
2618 pancakeBoundaryName: str = Field(
2619 default="PancakeBoundary",
2620 description=(
2621 "name of the pancake's curves that touches the air to be used in the mesh"
2622 ),
2623 )
2624 contactLayerBoundaryName: str = Field(
2625 default="contactLayerBoundary",
2626 description=(
2627 "name of the contact layers's curves that touches the air to be used in the"
2628 " mesh (only for TSA)"
2629 ),
2630 )
2633class Pancake3DMesh(BaseModel):
2634 # Mandatory:
2635 wi: Pancake3DMeshWinding = Field(
2636 alias="winding",
2637 title="Winding Mesh",
2638 description="This dictionary contains the winding mesh information.",
2639 )
2640 ii: Pancake3DMeshContactLayer = Field(
2641 alias="contactLayer",
2642 title="Contact Layer Mesh",
2643 description="This dictionary contains the contact layer mesh information.",
2644 )
2646 # Optionals:
2647 ti: Pancake3DMeshAirAndTerminals = Field(
2648 default=Pancake3DMeshAirAndTerminals(),
2649 alias="terminals",
2650 title="Terminal Mesh",
2651 description="This dictionary contains the terminal mesh information.",
2652 )
2653 ai: Pancake3DMeshAirAndTerminals = Field(
2654 default=Pancake3DMeshAirAndTerminals(),
2655 alias="air",
2656 title="Air Mesh",
2657 description="This dictionary contains the air mesh information.",
2658 )
2660 # Mandatory:
2661 relSizeMin: PositiveFloat = Field(
2662 alias="minimumElementSize",
2663 title="Minimum Element Size",
2664 description=(
2665 "The minimum mesh element size in terms of the largest mesh size in the"
2666 " winding. This mesh size will be used in the regions close the the"
2667 " winding, and then the mesh size will increate to maximum mesh element"
2668 " size as it gets away from the winding."
2669 ),
2670 )
2671 relSizeMax: PositiveFloat = Field(
2672 alias="maximumElementSize",
2673 title="Maximum Element Size",
2674 description=(
2675 "The maximum mesh element size in terms of the largest mesh size in the"
2676 " winding. This mesh size will be used in the regions close the the"
2677 " winding, and then the mesh size will increate to maximum mesh element"
2678 " size as it gets away from the winding."
2679 ),
2680 )
2682 @field_validator("relSizeMax")
2683 @classmethod
2684 def check_rel_size_max(cls, relSizeMax, info: ValidationInfo):
2685 if relSizeMax < info.data["relSizeMin"]:
2686 raise ValueError(
2687 "Maximum mesh element size (maximumElementSize) cannot be smaller than"
2688 " the minimum mesh element size (minimumElementSize)!"
2689 )
2691 return relSizeMax
2693 @computed_field
2694 @cached_property
2695 def sizeMin(self) -> float:
2696 """Return the minimum mesh element size in real dimensions."""
2697 geometry = geometry_input.get()
2699 meshLengthsPerElement = []
2701 # azimuthal elements:
2702 for numberOfElements in self.wi.ane:
2703 terminalOuterRadius = (
2704 getWindingOuterRadius()
2705 + 2 * geometry["terminals"]["outer"]["thickness"]
2706 )
2707 meshLengthsPerElement.append(
2708 2 * math.pi * terminalOuterRadius / numberOfElements
2709 )
2711 # radial elements:
2712 for numberOfElements in self.wi.rne:
2713 meshLengthsPerElement.append(
2714 geometry["winding"]["thickness"] / numberOfElements
2715 )
2717 if not geometry["contactLayer"]["thinShellApproximation"]:
2718 # radial elements:
2719 for numberOfElements in self.ii.rne:
2720 meshLengthsPerElement.append(
2721 geometry["contactLayer"]["thickness"] / numberOfElements
2722 )
2724 # axial elements:
2725 for numberOfElements in self.wi.axne:
2726 meshLengthsPerElement.append(
2727 geometry["winding"]["height"] / numberOfElements
2728 )
2730 return max(meshLengthsPerElement) * self.relSizeMin
2732 @computed_field
2733 @cached_property
2734 def sizeMax(self) -> float:
2735 """Return the minimum mesh element size in real dimensions."""
2736 geometry = geometry_input.get()
2738 meshLengthsPerElement = []
2740 # azimuthal elements:
2741 for numberOfElements in self.wi.ane:
2742 terminalOuterRadius = (
2743 getWindingOuterRadius()
2744 + 2 * geometry["terminals"]["outer"]["thickness"]
2745 )
2746 meshLengthsPerElement.append(
2747 2 * math.pi * terminalOuterRadius / numberOfElements
2748 )
2750 # radial elements:
2751 for numberOfElements in self.wi.rne:
2752 meshLengthsPerElement.append(
2753 geometry["winding"]["thickness"] / numberOfElements
2754 )
2756 if not geometry["contactLayer"]["thinShellApproximation"]:
2757 # radial elements:
2758 for numberOfElements in self.ii.rne:
2759 meshLengthsPerElement.append(
2760 geometry["contactLayer"]["thickness"] / numberOfElements
2761 )
2763 # axial elements:
2764 for numberOfElements in self.wi.axne:
2765 meshLengthsPerElement.append(
2766 geometry["winding"]["height"] / numberOfElements
2767 )
2769 return max(meshLengthsPerElement) * self.relSizeMax
2771 @computed_field
2772 @cached_property
2773 def stopGrowingDistance(self) -> float:
2774 """Return the distance from the pancake coils to start growing the mesh elements."""
2775 geometry = geometry_input.get()
2776 terminalOuterRadius = (
2777 getWindingOuterRadius() + 2 * geometry["terminals"]["outer"]["thickness"]
2778 )
2780 if geometry["air"]["type"] == "cuboid":
2781 sideLength = geometry["air"]["sideLength"]
2782 stopGrowingDistance = sideLength / 2 - terminalOuterRadius
2783 if geometry["air"]["type"] == "cylinder":
2784 stopGrowingDistance = geometry["air"]["radius"] - terminalOuterRadius
2786 return stopGrowingDistance
2788 @computed_field
2789 @cached_property
2790 def startGrowingDistance(self) -> float:
2791 geometry = geometry_input.get()
2792 terminalOuterRadius = (
2793 getWindingOuterRadius() + 2 * geometry["terminals"]["outer"]["thickness"]
2794 )
2795 terminalInnerRadius = (
2796 geometry["winding"]["innerRadius"]
2797 - 2 * geometry["terminals"]["inner"]["thickness"]
2798 )
2800 startGrowingDistance = (terminalOuterRadius - terminalInnerRadius) / 2
2802 return startGrowingDistance
2805class Pancake3DSolve(BaseModel):
2806 # 1) User inputs:
2807 t: Pancake3DSolveTime = Field(
2808 alias="time",
2809 title="Time Settings",
2810 description="All the time related settings for transient analysis.",
2811 )
2813 nls: Optional[Pancake3DSolveNonlinearSolverSettings] = Field(
2814 alias="nonlinearSolver",
2815 title="Nonlinear Solver Settings",
2816 description="All the nonlinear solver related settings.",
2817 )
2819 wi: Pancake3DSolveWindingMaterial = Field(
2820 alias="winding",
2821 title="Winding Properties",
2822 description="This dictionary contains the winding material properties.",
2823 )
2824 ii: Pancake3DSolveContactLayerMaterial = Field(
2825 alias="contactLayer",
2826 title="Contact Layer Properties",
2827 description="This dictionary contains the contact layer material properties.",
2828 )
2829 ti: Pancake3DSolveTerminalMaterialAndBoundaryCondition = Field(
2830 alias="terminals",
2831 title="Terminals Properties",
2832 description=(
2833 "This dictionary contains the terminals material properties and cooling"
2834 " condition."
2835 ),
2836 )
2837 ai: Pancake3DSolveAir = Field(
2838 alias="air",
2839 title="Air Properties",
2840 description="This dictionary contains the air material properties.",
2841 )
2843 ic: Pancake3DSolveInitialConditions = Field(
2844 alias="initialConditions",
2845 title="Initial Conditions",
2846 description="Initial conditions of the problem.",
2847 )
2849 save: list[Pancake3DSolveSaveQuantity] = Field(
2850 alias="quantitiesToBeSaved",
2851 default=None,
2852 title="Quantities to be Saved",
2853 description="List of quantities to be saved.",
2854 )
2856 # Mandatory:
2857 type: Literal["electromagnetic", "thermal", "weaklyCoupled", "stronglyCoupled"] = (
2858 Field(
2859 title="Simulation Type",
2860 description=(
2861 "FiQuS/Pancake3D can solve only electromagnetics and thermal or"
2862 " electromagnetic and thermal coupled. In the weaklyCoupled setting,"
2863 " thermal and electromagnetics systems will be put into different"
2864 " matrices, whereas in the stronglyCoupled setting, they all will be"
2865 " combined into the same matrix. The solution should remain the same."
2866 ),
2867 )
2868 )
2870 # Optionals:
2871 proTemplate: str = Field(
2872 default="Pancake3D_template.pro",
2873 description="file name of the .pro template file",
2874 )
2876 localDefects: Pancake3DSolveLocalDefects = Field(
2877 default=Pancake3DSolveLocalDefects(),
2878 alias="localDefects",
2879 title="Local Defects",
2880 description=(
2881 "Local defects (like making a small part of the winding normal conductor at"
2882 " some time) can be introduced."
2883 ),
2884 )
2886 initFromPrevious: str = Field(
2887 default="",
2888 title="Full path to res file to continue from",
2889 description=(
2890 "The simulation is continued from an existing .res file. The .res file is"
2891 " from a previous computation on the same geometry and mesh. The .res file"
2892 " is taken from the folder Solution_<<initFromPrevious>>"
2893 ),
2894 )
2896 isothermalInAxialDirection: bool = Field(
2897 default=False,
2898 title="Equate DoF along axial direction",
2899 description=(
2900 "If True, the DoF along the axial direction will be equated. This means"
2901 " that the temperature will be the same along the axial direction reducing"
2902 " the number of DoF. This is only valid for the thermal analysis."
2903 ),
2904 )
2906 @computed_field
2907 @cached_property
2908 def systemsOfEquationsType(self) -> str:
2910 windingMaterialIsGiven = self.wi.material is not None
2911 contactLayerMaterialIsGiven = self.ii.material is not None
2912 terminalMaterialIsGiven = self.ti.material is not None
2913 notchMaterialIsGiven = self.ti.transitionNotch.material is not None
2914 terminalContactLayerMaterialIsGiven = self.ti.terminalContactLayer.material is not None
2916 if (
2917 windingMaterialIsGiven
2918 or contactLayerMaterialIsGiven
2919 or terminalMaterialIsGiven
2920 or notchMaterialIsGiven
2921 or terminalContactLayerMaterialIsGiven
2922 ):
2923 systemsOfEquationsType = "nonlinear"
2924 else:
2925 systemsOfEquationsType = "linear"
2927 return systemsOfEquationsType
2930class Pancake3DPostprocess(BaseModel):
2931 """
2932 TO BE UPDATED
2933 """
2935 # 1) User inputs:
2936 timeSeriesPlots: list[Pancake3DPostprocessTimeSeriesPlot] = Field(
2937 default=None,
2938 title="Time Series Plots",
2939 description="Values can be plotted with respect to time.",
2940 )
2942 magneticFieldOnCutPlane: Optional[Pancake3DPostprocessMagneticFieldOnPlane] = Field(
2943 default=None,
2944 title="Magnetic Field on a Cut Plane",
2945 description=(
2946 "Color map of the magnetic field on the YZ plane can be plotted with"
2947 " streamlines."
2948 ),
2949 )
2952class Pancake3D(BaseModel):
2953 """
2954 Level 1: Class for FiQuS Pancake3D
2955 """
2957 type: Literal["Pancake3D"]
2958 geometry: Pancake3DGeometry = Field(
2959 default=None,
2960 title="Geometry",
2961 description="This dictionary contains the geometry information.",
2962 )
2963 mesh: Pancake3DMesh = Field(
2964 default=None,
2965 title="Mesh",
2966 description="This dictionary contains the mesh information.",
2967 )
2968 solve: Pancake3DSolve = Field(
2969 default=None,
2970 title="Solve",
2971 description="This dictionary contains the solve information.",
2972 )
2973 postproc: Pancake3DPostprocess = Field(
2974 default=None,
2975 title="Postprocess",
2976 description="This dictionary contains the postprocess information.",
2977 )
2978 input_file_path: str = Field(
2979 default=None,
2980 description="path of the input file (calculated by FiQuS)",
2981 exclude=True,
2982 )
2984 @model_validator(mode="before")
2985 @classmethod
2986 def set_context_variables(cls, model: "Pancake3D"):
2987 """Set the context variables (geometry data model, mesh data model, solve data
2988 model) of the Pancake3D model. They will be used in the submodel validators.
2989 """
2991 if isinstance(
2992 model["mesh"]["winding"]["azimuthalNumberOfElementsPerTurn"], int
2993 ):
2994 model["mesh"]["winding"]["azimuthalNumberOfElementsPerTurn"] = [
2995 model["mesh"]["winding"]["azimuthalNumberOfElementsPerTurn"]
2996 ] * model["geometry"]["numberOfPancakes"]
2998 if isinstance(model["mesh"]["winding"]["radialNumberOfElementsPerTurn"], int):
2999 model["mesh"]["winding"]["radialNumberOfElementsPerTurn"] = [
3000 model["mesh"]["winding"]["radialNumberOfElementsPerTurn"]
3001 ] * model["geometry"]["numberOfPancakes"]
3003 if isinstance(model["mesh"]["winding"]["axialNumberOfElements"], int):
3004 model["mesh"]["winding"]["axialNumberOfElements"] = [
3005 model["mesh"]["winding"]["axialNumberOfElements"]
3006 ] * model["geometry"]["numberOfPancakes"]
3008 if isinstance(model["mesh"]["winding"]["elementType"], str):
3009 model["mesh"]["winding"]["elementType"] = [
3010 model["mesh"]["winding"]["elementType"]
3011 ] * model["geometry"]["numberOfPancakes"]
3013 if isinstance(
3014 model["mesh"]["contactLayer"]["radialNumberOfElementsPerTurn"], int
3015 ):
3016 model["mesh"]["contactLayer"]["radialNumberOfElementsPerTurn"] = [
3017 model["mesh"]["contactLayer"]["radialNumberOfElementsPerTurn"]
3018 ] * model["geometry"]["numberOfPancakes"]
3020 geometry_input.set(model["geometry"])
3021 mesh_input.set(model["mesh"])
3022 solve_input.set(model["solve"])
3023 input_file_path.set(model["input_file_path"])
3025 return model