diff --git a/Conductors3D/w3d.pro b/Conductors3D/w3d.pro index 95e4db75a56e269a6c270bfe041573acfe4e025b..10eb5fa4b38546101ea214f5950694dccbf3e5a4 100644 --- a/Conductors3D/w3d.pro +++ b/Conductors3D/w3d.pro @@ -3,14 +3,14 @@ Include "w3d_common.pro"; Struct S::SOLVER [ Enum, mumps, gmres, fgmres, bcgs ]; DefineConstant[ - FileId = "" - DataDir = StrCat["data/" , FileId] + FileId = "" + DataDir = StrCat["data/" , FileId] ElekinDir = StrCat["Elekin_Line_Const/" , FileId] MagDynDir = StrCat["MagDyn_Line_Const/" , FileId] FullWaveDir = StrCat["FullWave_Line_Const/", FileId] type_solver = { - S::SOLVER.mumps, - Choices { + S::SOLVER.mumps, + Choices { S::SOLVER.mumps ="mumps", S::SOLVER.gmres="gmres", S::SOLVER.fgmres= "fgmres", @@ -24,7 +24,7 @@ DefineConstant[ R_ = {"Dynamic", Name "GetDP/1ResolutionChoices", Visible 0}, C_ = {"-solve -v2", Name "GetDP/9ComputeCommand", Visible 0}, P_ = {"", Name "GetDP/2PostOperationChoices", Visible 0} - ResDir = "Res/" + ResDir = "Res/" ]; @@ -43,7 +43,7 @@ Group{ For i In {1:NumWires} LWIRE~{i} = Region[ (50+i) ]; LWIRES += Region[ LWIRE~{i} ]; - + If( !Flag_Thin ) VWIRE~{i} = Region[ (10+i) ]; VWIRES += Region[ VWIRE~{i} ]; @@ -53,10 +53,10 @@ Group{ Else VWIRE~{i} = Region[ {} ]; ANODE~{i} = Region[ (60+i) ]; - CATHODE~{i} = Region[ (70+i) ]; + CATHODE~{i} = Region[ (70+i) ]; SWIRE~{i} = ElementsOf[ AIR, OnOneSideOf LWIRE~{i} ]; EndIf - + ANODES += Region[ ANODE~{i} ]; CATHODES += Region[ CATHODE~{i} ]; EndFor @@ -73,8 +73,8 @@ Group{ Vol_nu = Region[ {Vol_C, Vol_CC} ]; Sur_Dirichlet_a = Region[ INF ]; - Electrodes = Region[ {ANODES, CATHODES} ]; - + Electrodes = Region[ {ANODES, CATHODES} ]; + Dom_Hcurl_a = Region[ Vol_nu ]; Dom_Hgrad_u = Region[ Vol_C ]; Dom_Hregion_i = Region[ Vol_C ]; @@ -94,7 +94,7 @@ Function{ omega = 2*Pi*Freq; mu0 = Pi*4e-7; mur = 1.; - sigma = 5.96e7; + sigma = 5.96e7; mu[] = mu0; nu[] = 1/mu[]; @@ -110,26 +110,26 @@ Function{ For i In {1:NumWires} // Current density - J[ Region[ {VWIRE~{i}, LWIRE~{i}} ] ] = + J[ Region[ {VWIRE~{i}, LWIRE~{i}} ] ] = Vector[ 0, 0, WI~{i}/A_c * Exp[ i[]*WP~{i}*deg ] ]; // local radial coordinate with origin on wire i R~{i}[] = Hypot[ X[]-WX~{i} , Y[]-WY~{i} ]; EndFor - // shape function of the current, wI[r] - wI[] = tau[] * J0[tau[]*$1/rw]/J1[tau[]]/(2*A_c); + // shape function of the current, wI[r] + wI[] = tau[] * J0[tau[]*$1/rw]/J1[tau[]]/(2*A_c); //radial analytical solutions with a zero flux condition imposed at $1(=r)=rs Analytic_A[] = // per Amp (($1>rw) ? mu0/(2*Pi) * Log[rs/$1] : - mu0/(2*Pi) * Log[rs/rw] + i[]*( wI[$1]-wI[rw] )/(sigma*omega) ); + mu0/(2*Pi) * Log[rs/rw] + i[]*( wI[$1]-wI[rw] )/(sigma*omega) ); AnalyticStatic_A[] = mu0/(2*Pi) * // per Amp (($1>rw) ? Log[rs/$1] : Log[rs/rw] + mur*(1-($1/rw)^2)/2 ); Analytic_B[] = // per Amp (($1>=rw) ? mu0/(2*Pi*$1) : mu0*mur/(2*Pi*rw) * J1[tau[]*$1/rw] / J1[tau[]] ); - + // Impedance of thin wire p.u. length R_DC = 1./(sigma*A_c); Analytic_R[] = Re[ wI[rw]/sigma ]; @@ -146,8 +146,8 @@ Function{ EndIf If ( NumWires == 3 ) Exact_B[] = Analytic_B[R_1[]] + Analytic_B[R_2[]] + Analytic_B[R_3[]] ; - Correction_B[] = ((R_1[]<rs) ? Analytic_B[R_1[]] : - ((R_2[]<rs) ? Analytic_B[R_2[]] : + Correction_B[] = ((R_1[]<rs) ? Analytic_B[R_1[]] : + ((R_2[]<rs) ? Analytic_B[R_2[]] : ((R_3[]<rs) ? Analytic_B[R_3[]] : 0 ))); EndIf } @@ -184,19 +184,19 @@ Integration { Constraint{ { Name Impose_U ; Case { - For i In {1:NumWires} + For i In {1:NumWires} If( !Flag_Stranded ) // massive { Region ANODE~{i} ; Value 0 ; } EndIf - - If( Flag_U ) + + If( Flag_U ) If( !Flag_Stranded ) // massive { Region CATHODE~{i} ; Value -R_DC*Lz*WI~{i} ; } Else // stranded { Region VWIRE~{i} ; Value -R_DC*Lz*WI~{i} ; } EndIf EndIf - + EndFor } } @@ -210,7 +210,7 @@ Constraint{ Else // stranded { Region VWIRE~{i} ; Value WI~{i} ; } { Region LWIRE~{i} ; Value WI~{i} ; } - EndIf + EndIf EndIf EndFor } @@ -221,14 +221,14 @@ Constraint{ { Region Region[ Sur_Dirichlet_a ]; Value 0.; } } } - + { Name GaugeCondition_a; Type Assign; Case { - { Region Vol_Tree; SubRegion Region[ { Sur_Tree, Lin_Tree } ]; + { Region Vol_Tree; SubRegion Region[ { Sur_Tree, Lin_Tree } ]; Value 0; } } } - + // Semi-analytic approach { Name Impose_corr ; Case { @@ -266,20 +266,20 @@ FunctionSpace { NameOfConstraint GaugeCondition_a ; } } } - + // Electric scalar potential for massive conductors - { Name Hgrad_u_3D; Type Form0; + { Name Hgrad_u_3D; Type Form0; BasisFunction { { Name sn; NameOfCoef un; Function BF_Node; Support Dom_Hgrad_u; Entity NodesOf[ Dom_Hgrad_u, Not Electrodes]; } { Name sn2; NameOfCoef un2; Function BF_GroupOfNodes; - Support Dom_Hgrad_u; Entity GroupsOfNodesOf[ Electrodes ]; } + Support Dom_Hgrad_u; Entity GroupsOfNodesOf[ Electrodes ]; } } GlobalQuantity { { Name U; Type AliasOf ; NameOfCoef un2; } { Name I; Type AssociatedWith; NameOfCoef un2; } } - Constraint { + Constraint { { NameOfCoef U; EntityType GroupsOfNodesOf; NameOfConstraint Impose_U; } { NameOfCoef I; EntityType GroupsOfNodesOf; NameOfConstraint Impose_I; } } @@ -302,7 +302,7 @@ FunctionSpace { } // Function spaces for the semi-analytical approach - + // { Name Hcurl_acorr_3D; Type Form1; // BasisFunction { // { Name sw; NameOfCoef aw; Function BF_Edge; @@ -315,7 +315,7 @@ FunctionSpace { // NameOfConstraint GaugeCondition_a ; } // } // } - + { Name Hcurl_athin_3D; Type Form1; BasisFunction { { Name si; NameOfCoef ai; Function BF_GroupOfEdges; @@ -354,11 +354,11 @@ FunctionSpace { Formulation { { Name MagnetoDynamics; Type FemEquation; - + If( !Flag_SemiAnalytic ) // Conventional massive and stranded conductor formulation // If Flag_Thin is true, naive thine wire formulations - + If( !Flag_Stranded ) // Massive conductor Quantity { { Name a; Type Local; NameOfSpace Hcurl_a_3D; } @@ -368,7 +368,7 @@ Formulation { { Name I; Type Global; NameOfSpace Hgrad_u_3D [I]; } } - Equation { + Equation { Integral { [ nu[] * Dof{d a} , {d a} ]; In Vol_nu; Jacobian Vol; Integration I1; } @@ -377,31 +377,31 @@ Formulation { Integral { [ sigma[] * Dof{d v} , {a} ]; In Vol_C; Jacobian Vol; Integration I1; } Integral { DtDof[ sigma[] * Dof{a} , {d v} ]; - In Vol_C; Jacobian Vol; Integration I1; } + In Vol_C; Jacobian Vol; Integration I1; } Integral { [ sigma[] * Dof{d v} , {d v} ]; - In Vol_C; Jacobian Vol; Integration I1; } - - GlobalTerm { [Dof{I} , {U} ]; In CATHODES; } + In Vol_C; Jacobian Vol; Integration I1; } + + GlobalTerm { [Dof{I} , {U} ]; In CATHODES; } } - + Else // stranded conductor Quantity { { Name a; Type Local; NameOfSpace Hcurl_a_3D; } - + { Name i; Type Local; NameOfSpace Hregion_i_3D; } { Name Is; Type Global; NameOfSpace Hregion_i_3D [Is]; } { Name Us; Type Global; NameOfSpace Hregion_i_3D [Us]; } } - + Equation { Integral { [ nu[] * Dof{d a} , {d a} ]; In Vol_nu; Jacobian Vol; Integration I1; } Integral { [ -Dof{i}/A_c , {a} ]; In Vol_C; Jacobian Vol; Integration I1; } - + // Integral { [ -J[] , {a} ]; // In Vol_C; Jacobian Vol; Integration I1; } - + // Us is the voltage drop // grad v = -Dt a - j/sigma Integral { DtDof [ Dof{a} , {i} ]; @@ -411,47 +411,47 @@ Formulation { GlobalTerm { [ Dof{Us}*A_c , {Is} ]; In Vol_C; } } EndIf - + EndIf - + If( Flag_SemiAnalytic ) - // Semi-analytic correction methods - + // Semi-analytic correction methods + Quantity { { Name a; Type Local; NameOfSpace Hcurl_athin_3D; } { Name F; Type Global; NameOfSpace Hcurl_athin_3D [F]; } { Name I; Type Global; NameOfSpace Hcurl_athin_3D [I]; } - + { Name as; Type Local; NameOfSpace Hcurl_asleeve_3D; } { Name Fs; Type Global; NameOfSpace Hcurl_asleeve_3D [F]; } { Name Is; Type Global; NameOfSpace Hcurl_asleeve_3D [I]; } } - + Equation { - + Integral { [ nu[] * Dof{d a} , {d a} ]; In Vol_nu; Jacobian Vol; Integration I1; } - GlobalTerm { [ -Dof{I}*NbDivision , {F} ]; - In Vol_nu; } - + GlobalTerm { [ -Dof{I}*NbDivision , {F} ]; + In LWIRES ; } + Integral { [ nu[] * Dof{d as} , {d as} ]; In Vol_nu; Jacobian Vol; Integration I1; } - GlobalTerm { [ -Dof{Is}*NbDivision , {Fs} ]; - In Vol_nu; } - + GlobalTerm { [ -Dof{Is}*NbDivision , {Fs} ]; + In LWIRES ; } + // Integral { [ -Dof{i}/A_c , {as} ]; // In Vol_C; Jacobian Vol; Integration I1; } - - /* - GlobalTerm { [ Analytic_R[] * Dof{Is} , {Is} ]; + + /* + GlobalTerm { [ Analytic_R[] * Dof{Is} , {Is} ]; In Vol_C;} Integral { DtDof [ Dof{a}/A_c , {i} ]; In Vol_C; Jacobian Vol; Integration I1; } Integral { DtDof [ -Dof{as}/A_c , {i} ]; In Vol_C; Jacobian Vol; Integration I1;} - GlobalTerm { DtDof [ Analytic_L[] * Dof{Is} , {Is} ]; + GlobalTerm { DtDof [ Analytic_L[] * Dof{Is} , {Is} ]; In Vol_C;} - GlobalTerm { [ Dof{Us} , {Is} ]; + GlobalTerm { [ Dof{Us} , {Is} ]; In Vol_C; } */ } @@ -468,23 +468,23 @@ Resolution { } Operation { CreateDir[DataDir]; - - Generate[S]; Solve[S]; SaveSolution[S]; - + + Generate[S]; Solve[S]; SaveSolution[S]; + If( Flag_Maps == 1 ) // DeleteFile["U.dat"]; // DeleteFile["I.dat"]; // DeleteFile["Z.dat"]; - + PostOperation[map]; EndIf PostOperation[integaz]; PostOperation[cut]; For i In {1:NumWires} - Print[ {i, rw, rs, A_c, skin_depth, R_DC, WI~{i}/A_c, mu0*WI~{i}/(2*Pi*rw)}, + Print[ {i, rw, rs, A_c, skin_depth, R_DC, WI~{i}/A_c, mu0*WI~{i}/(2*Pi*rw)}, Format "WIRE %2g: rw = %7.3e rs = %7.3e A_c = %7.3e delta = %7.3e R_DC = %7.3e J = %7.3e Bmax = %7.3e"] ; - // Print[ {i, rw, Sqrt[ $sarea~{i}/Pi], rs, $sarea~{i}, Pi*rs^2, skin_depth, $intby~{i}/$sarea~{i}}, + // Print[ {i, rw, Sqrt[ $sarea~{i}/Pi], rs, $sarea~{i}, Pi*rs^2, skin_depth, $intby~{i}/$sarea~{i}}, // Format "WIRE %2g: rw = %7.3e rs = %7.3e (%7.3e) as = %7.3e (%7.3e) delta = %7.3e bave = %7.3e" ] ; EndFor } @@ -499,9 +499,9 @@ PostProcessing { { Name b; Value { Local { [ {d a}]; In Dom_Hcurl_a; Jacobian Vol; }}} { Name by; - Value { Local { [ Cart2Pol[ Norm[ {d a}]] ]; + Value { Local { [ Cart2Pol[ Norm[ {d a}]] ]; In Dom_Hcurl_a; Jacobian Vol; }}} - + If( Flag_SemiAnalytic ) { Name bs; Value { Local { [ {d as} ]; In Vol_nu; Jacobian Vol; }}} @@ -510,16 +510,16 @@ PostProcessing { { Name bbsy; Value { Local { [ Norm[ {d a} - {d as} ] ]; In Vol_nu; Jacobian Vol; }}} { Name bcorry; - Value { Local { [ Norm[ {d a} - {d as} ] + Correction_B[] ]; + Value { Local { [ Norm[ {d a} - {d as} ] + Correction_B[] ]; In Vol_nu; Jacobian Vol; }}} EndIf - - - + + + If( !Flag_SemiAnalytic ) If( !Flag_Stranded ) // massive { Name J; - Value { Local { [ -sigma[] * ( Dt[{a}] + {d v} ) ]; + Value { Local { [ -sigma[] * ( Dt[{a}] + {d v} ) ]; In Vol_C; Jacobian Vol; }}} { Name U; Value { Term { [ {U} ]; In Electrodes; }}} @@ -531,7 +531,7 @@ PostProcessing { Value { Term { [ Im[ -{U}/{I}/omega/Lz ] ]; In Electrodes; }}} Else // stranded { Name J; - Value { Local { [ {Is}/A_c ]; + Value { Local { [ {Is}/A_c ]; In Vol_C; Jacobian Vol; }}} { Name U; Value { Term { [ {Us} ]; In Vol_C; }}} @@ -543,14 +543,14 @@ PostProcessing { Value { Term { [ Im[ -{Us}/{Is}/omega/Lz ] ]; In Vol_C; }}} EndIf EndIf - - If( Flag_SemiAnalytic ) + + If( Flag_SemiAnalytic ) { Name F; Value { Term { [ {F} ]; In Vol_C; }}} { Name I; Value { Term { [ {I} ]; In Vol_C; }}} { Name J; - Value { Local { [ {I}/A_c ]; + Value { Local { [ {I}/A_c ]; In Vol_C; Jacobian Vol; }}} { Name U; Value { Term { [ i[]*omega*{F}*10 + Analytic_R[]*{I} ]; In Vol_C; }}} @@ -559,15 +559,15 @@ PostProcessing { { Name L; // actually total flux/I Value { Term { [ {F}/{I}/Lz ]; In Vol_C; }}} EndIf - + } } } PostOperation map UsingPost MagnetoDynamics { - + Print[ b, OnElementsOf Vol_Tree, File "b.pos"]; - + If( !Flag_SemiAnalytic ) Print[ J, OnElementsOf Vol_C, File "j.pos"]; Else @@ -632,7 +632,7 @@ PostOperation cut UsingPost MagnetoDynamics { "View[l].LineWidth = 3;", "View[l].Type = 2;"], File "tmp.geo", LastTimeStepOnly]; - + If( Flag_SemiAnalytic ) Print [ bsy, OnLine { {-Xcut,0,0} {Xcut,0,Zcut} }{NbPoints}, Format Gmsh, File "bcorry.pos" ]; @@ -661,12 +661,12 @@ PostOperation cut UsingPost MagnetoDynamics { "View[l].Type = 2;"], File "tmp.geo", LastTimeStepOnly]; EndIf - + // cuts in txt format for Gnuplot // Print [ exact, OnLine { {0,0,0} {Xcut,0,0} } {NbPoints}, // Format SimpleTable, File "Cut_analytic.txt" ]; - + If( !Flag_SemiAnalytic ) Print [ by, OnLine { {0,0,0} {Xcut,0,0} } {NbPoints}, Format SimpleTable, File "Cut_byref.txt" ]; @@ -681,5 +681,3 @@ PostOperation cut UsingPost MagnetoDynamics { Format SimpleTable, File "Cut_bcorry.txt" ]; EndIf } - -