diff --git a/HelicoidalSymmetricCable/helicoidal_power_cable.pro b/HelicoidalSymmetricCable/helicoidal_power_cable.pro new file mode 100644 index 0000000000000000000000000000000000000000..315142169707ecd10594125ea5775a743313d762 --- /dev/null +++ b/HelicoidalSymmetricCable/helicoidal_power_cable.pro @@ -0,0 +1,707 @@ +Printf("GetDP solver code for helicoidally symmetric power cables"); + +/* + +Documentation of the model: +original IEEE version: https://ieeexplore.ieee.org/document/10006756 +arXiv version: https://arxiv.org/abs/2301.03370 + + +The 10 basic GetDP steps are: +01. Group → Generate interface to mesh-file (.msh) via number of physical groups, + s.t. e.g. GetDP routines like FunctionSpace or Formulation can use those regions. +02. Function → Definition of "magnetoquasistatic" material mappings σ (conductivity) and + μ (permeablility) for different regions in the model. +03. Jacobian → Jacobian methods can be referred to in Formulation and PostProcessing objects + to be used in the computation of integral terms and for changes of coordinates. +04. Integration → Integration methods can be referred to in Formulation and PostProcessing + objects to be used in the computation of integral terms ... with options. +05. Constraint → Constraints can be referred to in FunctionSpace objects to be used for + boundary conditions ... (HERE): Fixing the cohomology conditions of the BVP! +06. FunctionSpace → is characterized by the type of its interpolated fields, + one or several basis functions and optional constraints (in space and time). + → (HERE): V_uvw(Ω) (where magnetic field H_uvw is searched, i.e., H_uvw ∈ V_uvw(Ω)) + = V_uv(Ω) ⊕ V_w(Ω) + = V_uv(Ω_i) ⊕ V_uv(Ω_c) ⊕ V_w(Ω_i) ⊕ V_w(Ω_c) + = [grad(nodal fcns.)(Ω_i)_uv & ℋ¹(Ω_i)_uv & edge_fcns.(Ω_c)_uv] + ⊕ [const.(Ω_i)_w & perpendicular_edge_fcns.(Ω_c)_w] +07. Formulation → The Formulation tool permits to deal with volume, surface and line + integrals with many kinds of densities to integrate, written in a form that + is similar to their symbolic expressions ... +08. Resolution → "∂ₜ → jω" & Solving systems of equations. +09. PostProcessing → The PostProcessing object is based on the quantities defined in a + Formulation and permits the construction (thanks to the expression syntax) + of any useful piecewise defined quantity of interest. +10. PostOperation → The PostOperation is the bridge between results obtained with GetDP + and the external world. +*/ + +/////////////////////////////////////////////// +//////////// step 01. Group /////////////////// +/////////////////////////////////////////////// + +Group { + + Omega_c = Region[{3}]; // (whole) conducting domain, i.e. contains its boundary + Omega_i = Region[{5}]; // (whole) insulating domain, i.e. contains its boundary + Omega = Region[{Omega_c, Omega_i}]; // (whole) computational domain, i.e. contains its boundary + + BdOmega_c = Region[{4}]; // boundary of conducting domain (consists of several bean-shaped outlines of the conductors) + BdOmega_i = Region[{6}]; // boundary of insulating domain (consists of several bean-shaped outlines of the conductors & an outer circle) + BdOmega = Region[{2}]; // boundary of computational domain (is an outer circle) + + // NOTE: here comes now N thick cuts for N cohomology basis functions for N conductors: + + Cut1TO = Region[{33}]; // thick cut for the 1. cohomology basis function + Cut2TO = Region[{34}]; // thick cut for the 2. cohomology basis function + Cut3TO = Region[{35}]; // thick cut for the 3. cohomology basis function + Cut4TO = Region[{36}]; // thick cut for the 4. cohomology basis function + Cut5TO = Region[{37}]; // thick cut for the 5. cohomology basis function + Cut6TO = Region[{38}]; // thick cut for the 6. cohomology basis function + Cut7TO = Region[{39}]; // thick cut for the 7. cohomology basis function + Cut8TO = Region[{40}]; // thick cut for the 8. cohomology basis function + Cut9TO = Region[{41}]; // thick cut for the 9. cohomology basis function + Cut10TO = Region[{42}]; // thick cut for the 10. cohomology basis function + Cut11TO = Region[{43}]; // thick cut for the 11. cohomology basis function + Cut12TO = Region[{44}]; // thick cut for the 12. cohomology basis function + Cut13TO = Region[{45}]; // thick cut for the 13. cohomology basis function + } // end of Group block + +/////////////////////////////////////////////// +//////////// step 02. Function //////////////// +/////////////////////////////////////////////// + +Function { + + // Definition of "magnetoquasistatic" material mappings σ (conductivity) and μ (permeablility) in the cartesian chart: + + // values: + mu_0 = 4*Pi*1e-7; // vacuum permeablility μ₀ [VsA⁻¹m⁻¹] + sigma_Cu = 58.1e6; // electric conductivity of copper σ_Cu [AV⁻¹m⁻¹] + + // permeability μ [VsA⁻¹m⁻¹]: + mu[Omega_i] = mu_0; // ∀ x ∈ Ω_i : μ(x) = μ₀ + mu[Omega_c] = mu_0; // ∀ x ∈ Ω_c : μ(x) = μ₀ + mu[BdOmega_c] = mu_0; // ∀ x ∈ ∂Ω_c: μ(x) = μ₀ + + // conductivity σ [AV⁻¹m⁻¹] and resistivity ρ [VmA¯¹]: + sigma[Omega_i] = 0.; // ∀ x ∈ Ω_i: σ(x) = 0 + sigma[Omega_c] = sigma_Cu; // ∀ x ∈ Ω_c: σ(x) = σ_Cu + rho[Omega_c] = 1./sigma_Cu; // ∀ x ∈ Ω_c: ρ(x) = 1/σ_Cu + + // excitation parameters: + Itot[] = Sqrt[2]/13; // amplitude of the applied AC current through each of the 13 conductors + Freq = 50; // frequency of the applied AC current + + // parameter used in the helicoidal transformation (coordinate transformation between (x, y, z) and (u, v, w) coords): + alpha = 1.0; + beta = 0.03183098861837907; + + xyz_to_uvw[] = Vector[(+X[]*Cos[alpha/beta*Z[]]+Y[]*Sin[alpha/beta*Z[]]), (-X[]*Sin[alpha/beta*Z[]]+Y[]*Cos[alpha/beta*Z[]]), Z[]]; // coordinate transformation φ: p_xyz → p_uvw + uvw_to_xyz[] = Vector[(+X[]*Cos[alpha/beta*Z[]]-Y[]*Sin[alpha/beta*Z[]]), (+X[]*Sin[alpha/beta*Z[]]+Y[]*Cos[alpha/beta*Z[]]), Z[]]; // coordinate transformation φ⁻¹: p_uvw → p_xyz + + // transformation matrices: + + // full J: + // cos(w*α/β) -sin(w*α/β) -α*(u*sin(w*α/β) + v*cos(w*α/β))/β + // sin(w*α/β) cos(w*α/β) α*(u*cos(w*α/β) - v*sin(w*α/β))/β + // 0 0 1 + J[] = Tensor[1, 0, -alpha/beta*Y[], 0, 1, alpha/beta*X[], 0, 0, 1]; // Jacobian of φ¯¹ at w = 0 + + // full Jtrans: + // +cos(w*α/β) sin(w*α/β) 0 + // -sin(w*α/β) cos(w*α/β) 0 + // -α*(u*sin(w*α/β) + v*cos(w*α/β))/β α*(u*cos(w*α/β) - v*sin(w*α/β))/β 1 + Jtrans[] = Tensor[1, 0, 0, 0, 1, 0, -alpha/beta*Y[], alpha/beta*X[], 1]; // Transposed Jacobian of φ¯¹ at w = 0 + + // full Jinv: + // cos(w*α/β) sin(w*α/β) v*α/β + // -sin(w*α/β) cos(w*α/β) -u*α/β + // 0 0 1 + Jinv[] = Tensor[1, 0, alpha/beta*Y[], 0, 1, -alpha/beta*X[], 0, 0, 1]; // inverse Jacobian of φ¯¹ at w = 0 + + // full Jinvtrans: + // cos(w*α/β) -sin(w*α/β) 0 + // sin(w*α/β) cos(w*α/β) 0 + // v*α/β -u*α/β 1 + Jinvtrans[] = Tensor[1, 0, 0, 0, 1, 0, alpha/beta*Y[], -alpha/beta*X[], 1]; // Transposed inverse Jacobian of φ¯¹ at w = 0 + + // full T: T = Jtrans * J + // 1 0 -v*α/β + // 0 1 u*α/β + // -v*α/β u*α/β (u^2*α^2 + v^2*α^2 + β^2)/β^2 + T[] = Tensor[1, 0, -alpha/beta*Y[], 0, 1, alpha/beta*X[], -alpha/beta*Y[], alpha/beta*X[], (X[]^2*alpha^2 + Y[]^2*alpha^2 + beta^2)/beta^2]; + + // Tinv at w = 0: (Jinv * Jinvtrans)(w=0) + // v^2*α^2/β^2 + 1 -u*v*α^2/β^2 v*α/β + // -u*v*α^2/β^2 u^2*α^2/β^2 + 1 -u*α/β + // v*α/β -u*α/β 1 + // Tinv[] = Tensor[1, 0, -alpha/beta*Y[], 0, 1, alpha/beta*X[], -alpha/beta*Y[], alpha/beta*X[], (X[]^2*alpha^2 + Y[]^2*alpha^2 + beta^2)/beta^2]; + Tinv[] = Tensor[(alpha/beta)^2 * Y[]^2 + 1, -(alpha/beta)^2 * X[] * Y[], alpha/beta*Y[], -(alpha/beta)^2 * X[] * Y[], (alpha/beta)^2 * X[]^2 + 1, -alpha/beta*X[], alpha/beta*Y[], -alpha/beta*X[], 1]; + + // material transformation: + mu_uvw[] = mu[] * Tinv[]; // transformation μ_xyz → μ_uvw, see Section III.A in the paper + rho_uvw[] = rho[] * T[]; // transformation ρ_xyz → ρ_uvw, see Section III.A in the paper + +} // end of Function block + +/////////////////////////////////////////////// +//////////// step 03. Jacobian //////////////// +/////////////////////////////////////////////// + +Jacobian { + { Name Vol ; + Case { { Region All ; Jacobian Vol ; } } + } + { Name Sur ; + Case { { Region All ; Jacobian Sur ; } } + } + { Name Lin ; + Case { { Region All ; Jacobian Lin ; } } + } +} // end of Jacobian block + + +/////////////////////////////////////////////// +//////////// step 04. Integration ///////////// +/////////////////////////////////////////////// + +Integration { + { Name Int ; + Case { + { Type Gauss ; + Case { + { GeoElement Point ; NumberOfPoints 1 ; } + { GeoElement Line ; NumberOfPoints 3 ; } + { GeoElement Triangle ; NumberOfPoints 4 ; } + { GeoElement Quadrangle ; NumberOfPoints 4 ; } + { GeoElement Tetrahedron ; NumberOfPoints 15 ; } // valid choices are 1, 4, 5, 15, 16, ... + { GeoElement Hexahedron ; NumberOfPoints 6 ; } + { GeoElement Prism ; NumberOfPoints 6 ; } + } // end of inner Case + } // end of "Type Gauss" + } // end of outer Case + } // end of "Int" +} // end of Integration block + +/////////////////////////////////////////////// +//////////// step 05. Constraint ////////////// +/////////////////////////////////////////////// + +Constraint { + + // constraint for magnetic scalar potential: + { Name phi_constraint ; + // not needed here + } // end Constraint phi_constraint + + // constraint for w-component of magnetic field: + { Name H_w_constraint ; + Case { + {Region Omega_i ; Value Sqrt[2] * 1.0/(0.03183098861837907*2*Pi) -2.335961833704558 ;} + {Region BdOmega_c ; Value Sqrt[2] * 1.0/(0.03183098861837907*2*Pi) -2.335961833704558 ;} + } // end of Case + } // end Constraint H_w_constraint + + { Name VoltageTO ; + Case { + } + } // end of Constraint VoltageTO + + // constraint for cohomomology basis functions, which are part of the function space for h_uv: + { Name CurrentTO ; + Case { + { Region Cut1TO; Value Itot[] ; } + + { Region Cut2TO; Value Itot[] ; } + + { Region Cut3TO; Value Itot[] ; } + + { Region Cut4TO; Value Itot[] ; } + + { Region Cut5TO; Value Itot[] ; } + + { Region Cut6TO; Value Itot[] ; } + + { Region Cut7TO; Value Itot[] ; } + + { Region Cut8TO; Value Itot[] ; } + + { Region Cut9TO; Value Itot[] ; } + + { Region Cut10TO; Value Itot[] ; } + + { Region Cut11TO; Value Itot[] ; } + + { Region Cut12TO; Value Itot[] ; } + + { Region Cut13TO; Value Itot[] ; } + + } // end Case + } // end of Constraint CurrentTO +} // end of Constraint block + +/////////////////////////////////////////////// +////////// step 06. FunctionSpace ///////////// +/////////////////////////////////////////////// + +FunctionSpace { + + + // function space for h_w: + { Name V_w; Type Form1P; + BasisFunction { + { Name sn; NameOfCoef hn; Function BF_PerpendicularEdge; + Support Omega_c; Entity NodesOf[All, Not BdOmega_c]; } + { Name rn; NameOfCoef occi; Function BF_RegionZ; + Support Region[Omega_i]; Entity Region[Omega_i]; } + { Name rn; NameOfCoef occi2; Function BF_GroupOfPerpendicularEdges; + Support Omega_c; Entity GroupsOfNodesOf[BdOmega_c]; } + } + Constraint { + { NameOfCoef occi ; EntityType Region ; NameOfConstraint H_w_constraint ; } + { NameOfCoef occi2 ; EntityType GroupsOfNodesOf ; NameOfConstraint H_w_constraint ; } + } + } // end of FunctionSpace V_w + + // function space for h_uv: + { Name V_uv; Type Form1; + BasisFunction { + { Name sn; NameOfCoef phin; Function BF_GradNode; + Support Omega_i; Entity NodesOf[Omega_i]; } + + { Name sn; NameOfCoef phin2; Function BF_GroupOfEdges; + Support Omega_c; Entity GroupsOfEdgesOnNodesOf[BdOmega_c]; } + + { Name se; NameOfCoef t; Function BF_Edge; + Support Omega_c; Entity EdgesOf[All, Not BdOmega_c]; } + + { Name sc1; NameOfCoef I1; Function BF_GroupOfEdges; + Support Omega; Entity GroupsOfEdgesOf[Cut1TO]; } + + { Name sc2; NameOfCoef I2; Function BF_GroupOfEdges; + Support Omega; Entity GroupsOfEdgesOf[Cut2TO]; } + + { Name sc3; NameOfCoef I3; Function BF_GroupOfEdges; + Support Omega; Entity GroupsOfEdgesOf[Cut3TO]; } + + { Name sc4; NameOfCoef I4; Function BF_GroupOfEdges; + Support Omega; Entity GroupsOfEdgesOf[Cut4TO]; } + + { Name sc5; NameOfCoef I5; Function BF_GroupOfEdges; + Support Omega; Entity GroupsOfEdgesOf[Cut5TO]; } + + { Name sc6; NameOfCoef I6; Function BF_GroupOfEdges; + Support Omega; Entity GroupsOfEdgesOf[Cut6TO]; } + + { Name sc7; NameOfCoef I7; Function BF_GroupOfEdges; + Support Omega; Entity GroupsOfEdgesOf[Cut7TO]; } + + { Name sc8; NameOfCoef I8; Function BF_GroupOfEdges; + Support Omega; Entity GroupsOfEdgesOf[Cut8TO]; } + + { Name sc9; NameOfCoef I9; Function BF_GroupOfEdges; + Support Omega; Entity GroupsOfEdgesOf[Cut9TO]; } + + { Name sc10; NameOfCoef I10; Function BF_GroupOfEdges; + Support Omega; Entity GroupsOfEdgesOf[Cut10TO]; } + + { Name sc11; NameOfCoef I11; Function BF_GroupOfEdges; + Support Omega; Entity GroupsOfEdgesOf[Cut11TO]; } + + { Name sc12; NameOfCoef I12; Function BF_GroupOfEdges; + Support Omega; Entity GroupsOfEdgesOf[Cut12TO]; } + + { Name sc13; NameOfCoef I13; Function BF_GroupOfEdges; + Support Omega; Entity GroupsOfEdgesOf[Cut13TO]; } + +} // end of BasisFunction + +GlobalQuantity { + +{ Name Current1 ; Type AliasOf ; NameOfCoef I1 ; } +{ Name Voltage1 ; Type AssociatedWith ; NameOfCoef I1 ; } +{ Name Current2 ; Type AliasOf ; NameOfCoef I2 ; } +{ Name Voltage2 ; Type AssociatedWith ; NameOfCoef I2 ; } +{ Name Current3 ; Type AliasOf ; NameOfCoef I3 ; } +{ Name Voltage3 ; Type AssociatedWith ; NameOfCoef I3 ; } +{ Name Current4 ; Type AliasOf ; NameOfCoef I4 ; } +{ Name Voltage4 ; Type AssociatedWith ; NameOfCoef I4 ; } +{ Name Current5 ; Type AliasOf ; NameOfCoef I5 ; } +{ Name Voltage5 ; Type AssociatedWith ; NameOfCoef I5 ; } +{ Name Current6 ; Type AliasOf ; NameOfCoef I6 ; } +{ Name Voltage6 ; Type AssociatedWith ; NameOfCoef I6 ; } +{ Name Current7 ; Type AliasOf ; NameOfCoef I7 ; } +{ Name Voltage7 ; Type AssociatedWith ; NameOfCoef I7 ; } +{ Name Current8 ; Type AliasOf ; NameOfCoef I8 ; } +{ Name Voltage8 ; Type AssociatedWith ; NameOfCoef I8 ; } +{ Name Current9 ; Type AliasOf ; NameOfCoef I9 ; } +{ Name Voltage9 ; Type AssociatedWith ; NameOfCoef I9 ; } +{ Name Current10 ; Type AliasOf ; NameOfCoef I10 ; } +{ Name Voltage10 ; Type AssociatedWith ; NameOfCoef I10 ; } +{ Name Current11 ; Type AliasOf ; NameOfCoef I11 ; } +{ Name Voltage11 ; Type AssociatedWith ; NameOfCoef I11 ; } +{ Name Current12 ; Type AliasOf ; NameOfCoef I12 ; } +{ Name Voltage12 ; Type AssociatedWith ; NameOfCoef I12 ; } +{ Name Current13 ; Type AliasOf ; NameOfCoef I13 ; } +{ Name Voltage13 ; Type AssociatedWith ; NameOfCoef I13 ; } +} // end of GlobalQuantity + +Constraint { +{ NameOfCoef Current1 ; EntityType GroupsOfEdgesOf ; NameOfConstraint CurrentTO ; } +{ NameOfCoef Voltage1 ; EntityType GroupsOfEdgesOf ; NameOfConstraint VoltageTO ; } + + +{ NameOfCoef Current2 ; EntityType GroupsOfEdgesOf ; NameOfConstraint CurrentTO ; } +{ NameOfCoef Voltage2 ; EntityType GroupsOfEdgesOf ; NameOfConstraint VoltageTO ; } + + +{ NameOfCoef Current3 ; EntityType GroupsOfEdgesOf ; NameOfConstraint CurrentTO ; } +{ NameOfCoef Voltage3 ; EntityType GroupsOfEdgesOf ; NameOfConstraint VoltageTO ; } + + +{ NameOfCoef Current4 ; EntityType GroupsOfEdgesOf ; NameOfConstraint CurrentTO ; } +{ NameOfCoef Voltage4 ; EntityType GroupsOfEdgesOf ; NameOfConstraint VoltageTO ; } + + +{ NameOfCoef Current5 ; EntityType GroupsOfEdgesOf ; NameOfConstraint CurrentTO ; } +{ NameOfCoef Voltage5 ; EntityType GroupsOfEdgesOf ; NameOfConstraint VoltageTO ; } + + +{ NameOfCoef Current6 ; EntityType GroupsOfEdgesOf ; NameOfConstraint CurrentTO ; } +{ NameOfCoef Voltage6 ; EntityType GroupsOfEdgesOf ; NameOfConstraint VoltageTO ; } + + +{ NameOfCoef Current7 ; EntityType GroupsOfEdgesOf ; NameOfConstraint CurrentTO ; } +{ NameOfCoef Voltage7 ; EntityType GroupsOfEdgesOf ; NameOfConstraint VoltageTO ; } + + +{ NameOfCoef Current8 ; EntityType GroupsOfEdgesOf ; NameOfConstraint CurrentTO ; } +{ NameOfCoef Voltage8 ; EntityType GroupsOfEdgesOf ; NameOfConstraint VoltageTO ; } + + +{ NameOfCoef Current9 ; EntityType GroupsOfEdgesOf ; NameOfConstraint CurrentTO ; } +{ NameOfCoef Voltage9 ; EntityType GroupsOfEdgesOf ; NameOfConstraint VoltageTO ; } + + +{ NameOfCoef Current10 ; EntityType GroupsOfEdgesOf ; NameOfConstraint CurrentTO ; } +{ NameOfCoef Voltage10 ; EntityType GroupsOfEdgesOf ; NameOfConstraint VoltageTO ; } + + +{ NameOfCoef Current11 ; EntityType GroupsOfEdgesOf ; NameOfConstraint CurrentTO ; } +{ NameOfCoef Voltage11 ; EntityType GroupsOfEdgesOf ; NameOfConstraint VoltageTO ; } + + +{ NameOfCoef Current12 ; EntityType GroupsOfEdgesOf ; NameOfConstraint CurrentTO ; } +{ NameOfCoef Voltage12 ; EntityType GroupsOfEdgesOf ; NameOfConstraint VoltageTO ; } + + +{ NameOfCoef Current13 ; EntityType GroupsOfEdgesOf ; NameOfConstraint CurrentTO ; } +{ NameOfCoef Voltage13 ; EntityType GroupsOfEdgesOf ; NameOfConstraint VoltageTO ; } + + + { NameOfCoef phin ; EntityType NodesOf; NameOfConstraint phi_constraint; } + { NameOfCoef phin2; EntityType NodesOf; NameOfConstraint phi_constraint; } + + } // end of Constraint block (in FunctionSpace block) + } // end of FunctionSpace V_uv +} // end of FunctionSpace + + +/////////////////////////////////////////////// +////////// step 07. Formulation /////////////// +/////////////////////////////////////////////// + +Formulation { + + { Name MagDyn_htot_full; Type FemEquation; + + Quantity { + { Name h_uv; Type Local; NameOfSpace V_uv; } + { Name h_w; Type Local; NameOfSpace V_w; } + + +{ Name I1; Type Global; NameOfSpace V_uv[Current1]; } +{ Name V1; Type Global; NameOfSpace V_uv[Voltage1]; } + + +{ Name I2; Type Global; NameOfSpace V_uv[Current2]; } +{ Name V2; Type Global; NameOfSpace V_uv[Voltage2]; } + + +{ Name I3; Type Global; NameOfSpace V_uv[Current3]; } +{ Name V3; Type Global; NameOfSpace V_uv[Voltage3]; } + + +{ Name I4; Type Global; NameOfSpace V_uv[Current4]; } +{ Name V4; Type Global; NameOfSpace V_uv[Voltage4]; } + + +{ Name I5; Type Global; NameOfSpace V_uv[Current5]; } +{ Name V5; Type Global; NameOfSpace V_uv[Voltage5]; } + + +{ Name I6; Type Global; NameOfSpace V_uv[Current6]; } +{ Name V6; Type Global; NameOfSpace V_uv[Voltage6]; } + + +{ Name I7; Type Global; NameOfSpace V_uv[Current7]; } +{ Name V7; Type Global; NameOfSpace V_uv[Voltage7]; } + + +{ Name I8; Type Global; NameOfSpace V_uv[Current8]; } +{ Name V8; Type Global; NameOfSpace V_uv[Voltage8]; } + + +{ Name I9; Type Global; NameOfSpace V_uv[Current9]; } +{ Name V9; Type Global; NameOfSpace V_uv[Voltage9]; } + + +{ Name I10; Type Global; NameOfSpace V_uv[Current10]; } +{ Name V10; Type Global; NameOfSpace V_uv[Voltage10]; } + + +{ Name I11; Type Global; NameOfSpace V_uv[Current11]; } +{ Name V11; Type Global; NameOfSpace V_uv[Voltage11]; } + + +{ Name I12; Type Global; NameOfSpace V_uv[Current12]; } +{ Name V12; Type Global; NameOfSpace V_uv[Voltage12]; } + + +{ Name I13; Type Global; NameOfSpace V_uv[Current13]; } +{ Name V13; Type Global; NameOfSpace V_uv[Voltage13]; } + + } // end of Quantity + +Equation { + // Weak formumation reads: + // Seek H_uvw ∈ V_uvw(Ω), s.t. ∀ H_uvw' ∈ V₀_uvw(Ω): + // ∫_Ω_uvw ∂ₜ μ_uvw H_uvw ⋅ H_uvw' dV + ∫_Ω_c_uvw ρ_uvw curl_xyz H_uvw ⋅ curl_xyz H_uvw' dV = 0 + + // ∫_Ω_uvw ∂ₜ μ_uvw H_uvw ⋅ H_uvw' dV: + Galerkin { DtDof [ mu_uvw[] * Dof{h_uv} , {h_uv} ]; + In Omega; Integration Int; Jacobian Vol; } // combination 1: uv-components tested by uv-components of test functions + Galerkin { DtDof [ mu_uvw[] * Dof{h_w} , {h_w} ] ; + In Omega; Integration Int; Jacobian Vol; } // combination 2: w-component tested by w-component of test functions + Galerkin { DtDof [ mu_uvw[] * Dof{h_w} , {h_uv} ] ; + In Omega; Integration Int; Jacobian Vol; } // combination 3: w-components tested by uv-components of test functions + Galerkin { DtDof [ mu_uvw[] * Dof{h_uv} , {h_w} ] ; + In Omega; Integration Int; Jacobian Vol; } // combination 4: uv-components tested by w-component of test functions + + // ∫_Ω_c_uvw ρ_uvw curl_xyz H_uvw ⋅ curl_xyz H_uvw' dV: + Galerkin { [ rho_uvw[] * Dof{d h_uv} , {d h_uv} ]; + In Omega_c; Integration Int; Jacobian Vol; } // combination 1: uv-components tested by uv-components of test functions + Galerkin { [ rho_uvw[] * Dof{d h_w} , {d h_w} ]; + In Omega_c; Integration Int; Jacobian Vol; } // combination 2: w-component tested by w-component of test functions + Galerkin { [ rho_uvw[] * Dof{d h_w} , {d h_uv} ]; + In Omega_c; Integration Int; Jacobian Vol; } // combination 3: w-components tested by uv-components of test functions + Galerkin { [ rho_uvw[] * Dof{d h_uv} , {d h_w} ]; + In Omega_c; Integration Int; Jacobian Vol; } // combination 4: uv-components tested by w-component of test functions + + GlobalTerm { [ Dof{V1} , {I1} ] ; In Cut1TO ; } + GlobalTerm { [ Dof{V2} , {I2} ] ; In Cut2TO ; } + GlobalTerm { [ Dof{V3} , {I3} ] ; In Cut3TO ; } + GlobalTerm { [ Dof{V4} , {I4} ] ; In Cut4TO ; } + GlobalTerm { [ Dof{V5} , {I5} ] ; In Cut5TO ; } + GlobalTerm { [ Dof{V6} , {I6} ] ; In Cut6TO ; } + GlobalTerm { [ Dof{V7} , {I7} ] ; In Cut7TO ; } + GlobalTerm { [ Dof{V8} , {I8} ] ; In Cut8TO ; } + GlobalTerm { [ Dof{V9} , {I9} ] ; In Cut9TO ; } + GlobalTerm { [ Dof{V10} , {I10} ] ; In Cut10TO ; } + GlobalTerm { [ Dof{V11} , {I11} ] ; In Cut11TO ; } + GlobalTerm { [ Dof{V12} , {I12} ] ; In Cut12TO ; } + GlobalTerm { [ Dof{V13} , {I13} ] ; In Cut13TO ; } +} // end of Equation +} // end of MagDyn_htot_full +} // end of Formulation + +/////////////////////////////////////////////// +////////// step 08. Resolution //////////////// +/////////////////////////////////////////////// + +Resolution { + { Name MagDynTOComplex; + System { + { Name A; NameOfFormulation MagDyn_htot_full; + Type ComplexValue; Frequency Freq;} // "∂ₜ → jω" + } + Operation { + Generate[A]; Solve[A]; SaveSolution[A]; + } + } // end of MagDynTOComplex block +} // end of Resolution block + + /////////////////////////////////////////////// + ////////// step 09. PostProcessing //////////// + /////////////////////////////////////////////// + + PostProcessing { + + // h-formulation, total field + { Name MagDyn_htot_full; NameOfFormulation MagDyn_htot_full; + Quantity { + + // magnetic scalar potential: + { Name phi; Value{ Local{ [ {dInv h_uv} ] ; In Omega_i; Jacobian Vol; } } } + + // magnetic field in helicoidal and cartesian coords: + { Name h_uvw; Value{ Local{ [ ({h_uv} + {h_w}) ] ; In Omega; Jacobian Vol; } } } + { Name h_xyz; Value{ Local{ [ Jinvtrans[]*({h_uv} + {h_w}) ] ; In Omega; Jacobian Vol; } } } + { Name h_xyz_abs; Value{ Local{ [Norm[ Jinvtrans[]*({h_uv} + {h_w}) ]] ; In Omega; Jacobian Vol; } } } + + // current density in helicoidal and cartesian coords: + { Name j_uvw; Value{ Local{ [ ({d h_uv} + {d h_w}) ] ; In Omega_c; Jacobian Vol; } } } + { Name j_xyz; Value{ Local{ [J[] * ({d h_uv} + {d h_w}) ] ; In Omega_c; Jacobian Vol; } } } + { Name j_xyz_abs; Value{ Local{ [Norm[ J[] * ({d h_uv} + {d h_w}) ]] ; In Omega_c; Jacobian Vol; } } } + + // loss density: + { Name power_losses_xyz; Value{ Local{ [ rho[] * SquNorm [ J[] * ({d h_uv} + {d h_w}) ] ]; In Omega_c; Jacobian Vol; } } } + + // global quantities: + { Name dissPower; + Value{ + Integral{ + [rho[]/2 * SquNorm[ J[] * ({d h_uv} + {d h_w}) ]]; In Omega_c; Integration Int; Jacobian Vol; + } // end of Integral + } // end of Value + } + + + { Name I1 ; Value { Term { [ {I1} ] ; In Cut1TO ; } } } + { Name V1 ; Value { Term { [ {V1} ] ; In Cut1TO ; } } } + + + { Name I2 ; Value { Term { [ {I2} ] ; In Cut2TO ; } } } + { Name V2 ; Value { Term { [ {V2} ] ; In Cut2TO ; } } } + + + { Name I3 ; Value { Term { [ {I3} ] ; In Cut3TO ; } } } + { Name V3 ; Value { Term { [ {V3} ] ; In Cut3TO ; } } } + + + { Name I4 ; Value { Term { [ {I4} ] ; In Cut4TO ; } } } + { Name V4 ; Value { Term { [ {V4} ] ; In Cut4TO ; } } } + + + { Name I5 ; Value { Term { [ {I5} ] ; In Cut5TO ; } } } + { Name V5 ; Value { Term { [ {V5} ] ; In Cut5TO ; } } } + + + { Name I6 ; Value { Term { [ {I6} ] ; In Cut6TO ; } } } + { Name V6 ; Value { Term { [ {V6} ] ; In Cut6TO ; } } } + + + { Name I7 ; Value { Term { [ {I7} ] ; In Cut7TO ; } } } + { Name V7 ; Value { Term { [ {V7} ] ; In Cut7TO ; } } } + + + { Name I8 ; Value { Term { [ {I8} ] ; In Cut8TO ; } } } + { Name V8 ; Value { Term { [ {V8} ] ; In Cut8TO ; } } } + + + { Name I9 ; Value { Term { [ {I9} ] ; In Cut9TO ; } } } + { Name V9 ; Value { Term { [ {V9} ] ; In Cut9TO ; } } } + + + { Name I10 ; Value { Term { [ {I10} ] ; In Cut10TO ; } } } + { Name V10 ; Value { Term { [ {V10} ] ; In Cut10TO ; } } } + + + { Name I11 ; Value { Term { [ {I11} ] ; In Cut11TO ; } } } + { Name V11 ; Value { Term { [ {V11} ] ; In Cut11TO ; } } } + + + { Name I12 ; Value { Term { [ {I12} ] ; In Cut12TO ; } } } + { Name V12 ; Value { Term { [ {V12} ] ; In Cut12TO ; } } } + + + { Name I13 ; Value { Term { [ {I13} ] ; In Cut13TO ; } } } + { Name V13 ; Value { Term { [ {V13} ] ; In Cut13TO ; } } } + + } // end of Quantity +} // end of MagDyn_htot_full PostProcessing +} // end of PostProcessing + + +/////////////////////////////////////////////// +////////// step 10. PostOperation ///////////// +/////////////////////////////////////////////// + +PostOperation { + + { Name MagDyn; NameOfPostProcessing MagDyn_htot_full ; + Operation{ + Print[ phi, OnElementsOf Omega , File "phi.pos", Name "φ_magnetic [A]" ]; + Print[ h_xyz, OnElementsOf Omega , File "h_xyz.pos", Name "H(xyz) [A/m]" ]; + Print[ h_uvw, OnElementsOf Omega , File "h_uvw.pos", Name "H(uvw) [A/m]" ]; + Print[ h_xyz_abs, OnElementsOf Omega, File "h_xyz_abs.pos", Name "|H|(xyz) [A/m]" ]; + Print[ j_xyz, OnElementsOf Omega_c , File "j_xyz.pos", Name "J(xyz) [A/m²]" ]; + Print[ j_uvw, OnElementsOf Omega_c , File "j_uvw.pos", Name "J(uvw) [A/m²]" ]; + Print[ j_xyz_abs, OnElementsOf Omega_c , File "j_xyz_abs.pos", Name "|J|(xyz) [A/m²]" ]; + Print[ power_losses_xyz, OnElementsOf Omega_c , File "p_loss_xyz.pos", Name "p_loss(xyz) [W/m]" ]; + + + // Print[ phi, OnLine {{-0.045, 0, 0}{0.045, 0, 0}} {100}, File "φ_table.txt"]; + // Print[ h_xyz, OnLine {{-0.045, 0, 0}{0.045, 0, 0}} {1100}, File "h_xyz_table.txt"]; + // Print[ h_uvw, OnLine {{-0.045, 0, 0}{0.045, 0, 0}} {1100}, File "h_uvw_table.txt"]; + // Print[ j_xyz, OnLine {{-0.045, 0, 0}{0.045, 0, 0}} {1100}, File "j_xyz_table.txt"]; + // Print[ j_uvw, OnLine {{-0.045, 0, 0}{0.045, 0, 0}} {1100}, File "j_uvw_table.txt"]; + + + + // print "global" results: + Print[ dissPower[Omega_c], OnGlobal, Format Table, File "jouleLosses.txt" ]; + + + + + Print[I1, OnRegion Cut1TO, Format FrequencyTable, File >> "currentsAndVoltages.txt"]; // real and imaginary part of each current + Print[V1, OnRegion Cut1TO, Format FrequencyTable, File >> "currentsAndVoltages.txt"]; // real and imaginary part of each voltage (drop) + + Print[I2, OnRegion Cut2TO, Format FrequencyTable, File >> "currentsAndVoltages.txt"]; // real and imaginary part of each current + Print[V2, OnRegion Cut2TO, Format FrequencyTable, File >> "currentsAndVoltages.txt"]; // real and imaginary part of each voltage (drop) + + Print[I3, OnRegion Cut3TO, Format FrequencyTable, File >> "currentsAndVoltages.txt"]; // real and imaginary part of each current + Print[V3, OnRegion Cut3TO, Format FrequencyTable, File >> "currentsAndVoltages.txt"]; // real and imaginary part of each voltage (drop) + + Print[I4, OnRegion Cut4TO, Format FrequencyTable, File >> "currentsAndVoltages.txt"]; // real and imaginary part of each current + Print[V4, OnRegion Cut4TO, Format FrequencyTable, File >> "currentsAndVoltages.txt"]; // real and imaginary part of each voltage (drop) + + Print[I5, OnRegion Cut5TO, Format FrequencyTable, File >> "currentsAndVoltages.txt"]; // real and imaginary part of each current + Print[V5, OnRegion Cut5TO, Format FrequencyTable, File >> "currentsAndVoltages.txt"]; // real and imaginary part of each voltage (drop) + + Print[I6, OnRegion Cut6TO, Format FrequencyTable, File >> "currentsAndVoltages.txt"]; // real and imaginary part of each current + Print[V6, OnRegion Cut6TO, Format FrequencyTable, File >> "currentsAndVoltages.txt"]; // real and imaginary part of each voltage (drop) + + Print[I7, OnRegion Cut7TO, Format FrequencyTable, File >> "currentsAndVoltages.txt"]; // real and imaginary part of each current + Print[V7, OnRegion Cut7TO, Format FrequencyTable, File >> "currentsAndVoltages.txt"]; // real and imaginary part of each voltage (drop) + + Print[I8, OnRegion Cut8TO, Format FrequencyTable, File >> "currentsAndVoltages.txt"]; // real and imaginary part of each current + Print[V8, OnRegion Cut8TO, Format FrequencyTable, File >> "currentsAndVoltages.txt"]; // real and imaginary part of each voltage (drop) + + Print[I9, OnRegion Cut9TO, Format FrequencyTable, File >> "currentsAndVoltages.txt"]; // real and imaginary part of each current + Print[V9, OnRegion Cut9TO, Format FrequencyTable, File >> "currentsAndVoltages.txt"]; // real and imaginary part of each voltage (drop) + + Print[I10, OnRegion Cut10TO, Format FrequencyTable, File >> "currentsAndVoltages.txt"]; // real and imaginary part of each current + Print[V10, OnRegion Cut10TO, Format FrequencyTable, File >> "currentsAndVoltages.txt"]; // real and imaginary part of each voltage (drop) + + Print[I11, OnRegion Cut11TO, Format FrequencyTable, File >> "currentsAndVoltages.txt"]; // real and imaginary part of each current + Print[V11, OnRegion Cut11TO, Format FrequencyTable, File >> "currentsAndVoltages.txt"]; // real and imaginary part of each voltage (drop) + + Print[I12, OnRegion Cut12TO, Format FrequencyTable, File >> "currentsAndVoltages.txt"]; // real and imaginary part of each current + Print[V12, OnRegion Cut12TO, Format FrequencyTable, File >> "currentsAndVoltages.txt"]; // real and imaginary part of each voltage (drop) + + Print[I13, OnRegion Cut13TO, Format FrequencyTable, File >> "currentsAndVoltages.txt"]; // real and imaginary part of each current + Print[V13, OnRegion Cut13TO, Format FrequencyTable, File >> "currentsAndVoltages.txt"]; // real and imaginary part of each voltage (drop) + + } // end of Operation + } // end of PostProcessing "MagDyn" + } // end of PostOperation + + // some GetDP constants: + DefineConstant[ + R_ = {"MagDynTOComplex" , Name "GetDP/1ResolutionChoices" , Visible 1}, + C_ = {"-solve -v2 -pos" , Name "GetDP/9ComputeCommand" , Visible 1}, + P_ = {"MagDyn" , Name "GetDP/2PostOperationChoices", Visible 1} + ];