From 5f52022f3ef45cbebe351df0e10e7e36998f5c89 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Fran=C3=A7ois=20Henrotte?= <francois.henrotte@uclouvain.be> Date: Thu, 5 Mar 2020 16:11:55 +0100 Subject: [PATCH] reorganization --- Conductors3D/w3d.geo | 11 +- Conductors3D/w3d.pro | 490 +++++++++++++++++++++++------------- Conductors3D/w3d_common.pro | 51 ++-- 3 files changed, 354 insertions(+), 198 deletions(-) diff --git a/Conductors3D/w3d.geo b/Conductors3D/w3d.geo index d103979..58a8812 100644 --- a/Conductors3D/w3d.geo +++ b/Conductors3D/w3d.geo @@ -6,19 +6,19 @@ Mesh.VolumeEdges = 0; // hide volume edges Geometry.ExactExtrusion = 0; // to allow rotation of extruded shapes Solver.AutoMesh = 2; // always remesh if necessary (don't reuse mesh on disk) -lc1 = (Flag_Thin) ? MeshSizeThinWire*mm : MeshSizeThinWire*mm ; +lc1 = (Flag_Thin) ? MeshSizeThinWire*mm : MeshSizeWire*mm ; lc2 = box/10; // mesh size at outer surface llWires[] = {}; centerWires[] = {}; surfWires[] = {}; -If( !Flag_Thin || (Flag_Thin && Flag_Corr==2)) +If( !Flag_Thin || (Flag_Thin && Flag_RegSleeve)) // then circles are in the geometry // their radius is R=rw if Flag_Thin != 0 // or R=rs in the "regular sleeve case" (Flag_Corr == 2) - R = ( !Flag_Thin ) ? rw : rs; + R = ( Flag_Thin ) ? rs : rw; For i In {1:NumWires} pC = newp; Point(pC) = {WX~{i} , WY~{i} , 0, lc1}; @@ -96,7 +96,7 @@ If( !Flag_Thin ) // Wires with finite radius e[] = Extrude {0,0,Lz} { Surface{ s1 } ; /*Layers{1}; Recombine;*/ } ; Physical Volume ("AIR", 1) = { e[1] }; -ElseIf( Flag_Corr == 2 ) // regular sleeve +ElseIf( Flag_RegSleeve ) // regular sleeve Wires[]= {}; For i In {1:NumWires} @@ -122,19 +122,18 @@ Else e[] = Extrude {0,0,Lz} { Surface{ s1 } ; /*Layers{1}; Recombine;*/ } ; s2 = e[0]; v1 = e[1]; + Physical Volume ("AIR", 1) = { v1 }; For i In {1:NumWires} e[] = Extrude {0,0,Lz} { Point{ centerWires[i-1] } ; /*Layers{1};*/ } ; //Printf("e=", e[]); Physical Line (Sprintf("LWIRE_%g",i), 50+i) = { e[1] }; Physical Point (Sprintf("PANODE_%g",i), 60+i) = { centerWires[i-1] }; Physical Point (Sprintf("PCATHODE_%g",i), 70+i) = { e[0] }; - Physical Volume ("AIR", 1) = { v1 }; // embed line conductors in mesh Point { centerWires[i-1] } In Surface { s1 }; Point { e[0] } In Surface { s2 }; Curve { e[1] } In Volume { v1 }; EndFor - EndIf Physical Surface("INF", 2) = { CombinedBoundary{ Volume{ : }; } }; diff --git a/Conductors3D/w3d.pro b/Conductors3D/w3d.pro index 1441768..c4c5a47 100644 --- a/Conductors3D/w3d.pro +++ b/Conductors3D/w3d.pro @@ -22,7 +22,7 @@ DefineConstant[ DefineConstant[ R_ = {"Dynamic", Name "GetDP/1ResolutionChoices", Visible 0}, - C_ = {"-sol -pos", Name "GetDP/9ComputeCommand", Visible 0}, + C_ = {"-solve -v2", Name "GetDP/9ComputeCommand", Visible 0}, P_ = {"", Name "GetDP/2PostOperationChoices", Visible 0} ResDir = "Res/" ]; @@ -63,43 +63,33 @@ Group{ // Abstract regions - // Cette formulation peut contenir des conducteurs volumiques Vol_C - // et lin�iques lVol_C. If( !Flag_Thin ) Vol_C = Region[ VWIRES ]; Vol_CC = Region[ AIR ]; - //lVol_C = Region[ {} ]; Else Vol_C = Region[ {LWIRES} ]; Vol_CC = Region[ {VWIRES, AIR} ]; - //lVol_C = Region[ LWIRES ]; EndIf Vol_nu = Region[ {Vol_C, Vol_CC} ]; Sur_Dirichlet_a = Region[ INF ]; Electrodes = Region[ {ANODES, CATHODES} ]; - Vol_Tree = Region[ { Vol_nu } ]; - Sur_Tree = Region[ { Sur_Dirichlet_a /*, SKIN*/ } ]; - // If( Flag_Corr!= 0 ) - // Lin_Tree = Region[ { LWIRES /*, LINTREE*/ } ]; - // Else - Lin_Tree = Region[ {} ]; - - Dom_Hcurl_a = Region[ Vol_nu ]; Dom_Hgrad_u = Region[ Vol_C ]; - //Dom_Hregion_i = Region[ Vol_nu ]; - - // additional Groups for the local correction method - //Dom_Hthin_a = ElementsOf[ Vol_nu, OnOneSideOf LWIRES ]; - - Dom_Hthin_a = Region[ Vol_nu ]; - - // For integration on the volume elements (excluding line elements) - // of the sleeve - Vol_Sleeve = ElementsOf[ {Vol_CC}, OnOneSideOf LWIRES ]; + Dom_Hregion_i = Region[ Vol_C ]; + + // additional Groups for the semi_analytic approach + Dom_Hthin_a = ElementsOf[ Vol_nu, OnOneSideOf LWIRES ]; + //Vol_Tree = ElementsOf[ Vol_nu, Not Dom_Hthin_a ]; + Vol_Tree = Region[ { Vol_nu } ]; + Sur_Tree = Region[ { Sur_Dirichlet_a /*, SKIN*/ } ]; + If( !Flag_SemiAnalytic ) + Lin_Tree = Region[ {} ]; + Else + Lin_Tree = Region[ { LWIRES } ]; + EndIf } @@ -109,13 +99,13 @@ Function{ mur = 1.; sigma = 5.96e7; - i[] = Complex[0,1]; - mu[] = mu0; nu[] = 1/mu[]; sigma[] = sigma; skin_depth = Sqrt[ 2/(omega*sigma*mu0*mur) ]; - + + // Functions for the semi-analytic approach + i[] = Complex[0,1]; tau[] = Complex[-1,1]/skin_depth*rw; J0[] = JnComplex[0,$1]; J1[] = JnComplex[1,$1]; @@ -131,7 +121,7 @@ Function{ EndFor // shape function of the current, wI[r] - wI[] = tau[] * J0[tau[]*$1/rw] / J1[tau[]] /(2*A_c); + 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 @@ -140,8 +130,8 @@ Function{ 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) * J1[tau[]*$1/rw] / J1[tau[]] ); + (($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); @@ -155,11 +145,11 @@ Function{ If ( NumWires == 2 ) Exact_B[] = Analytic_B[R_1[]] + Analytic_B[R_2[]] ; Correction_B[] = ((R_1[]<rs) ? Analytic_B[R_1[]] : - ((R_2[]<rs) ? Analytic_B[R_2[]] : 0 )); + ((R_2[]<rs) ? Analytic_B[R_2[]] : 0 )); EndIf If ( NumWires == 3 ) Exact_B[] = Analytic_B[R_1[]] + Analytic_B[R_2[]] + Analytic_B[R_3[]] ; - Correction_A[] = ((R_1[]<rs) ? Analytic_B[R_1[]] : + 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 @@ -180,14 +170,14 @@ Integration { Case { { Type Gauss ; Case { - { GeoElement Point ; NumberOfPoints 1 ; } - { GeoElement Line ; NumberOfPoints 3 ; } - { GeoElement Triangle ; NumberOfPoints 3 ; } - { GeoElement Quadrangle ; NumberOfPoints 4 ; } - { GeoElement Prism ; NumberOfPoints 21 ; } - { GeoElement Tetrahedron ; NumberOfPoints 4 ; } - { GeoElement Hexahedron ; NumberOfPoints 6 ; } - { GeoElement Pyramid ; NumberOfPoints 8 ; } + { GeoElement Point ; NumberOfPoints 1 ; } + { GeoElement Line ; NumberOfPoints 3 ; } + { GeoElement Triangle ; NumberOfPoints 3 ; } + { GeoElement Quadrangle ; NumberOfPoints 4 ; } + { GeoElement Prism ; NumberOfPoints 21 ; } + { GeoElement Tetrahedron ; NumberOfPoints 4 ; } + { GeoElement Hexahedron ; NumberOfPoints 6 ; } + { GeoElement Pyramid ; NumberOfPoints 8 ; } } } } @@ -197,19 +187,19 @@ Integration { Constraint{ { Name Impose_U ; Case { - For i In {1:NumWires} - If( Flag_Form == 1) // massive + For i In {1:NumWires} + If( !Flag_Stranded ) // massive { Region ANODE~{i} ; Value 0 ; } EndIf - If( Flag_U ) - - If( Flag_Form == 1) // massive + + 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 } } @@ -217,30 +207,48 @@ Constraint{ Case { For i In {1:NumWires} If( !Flag_U ) - - If( Flag_Form == 1) // massive + If( !Flag_Stranded ) // massive { Region ANODE~{i} ; Value WI~{i} ; } { Region CATHODE~{i} ; Value WI~{i} ; } Else // stranded { Region VWIRE~{i} ; Value WI~{i} ; } { Region LWIRE~{i} ; Value WI~{i} ; } - EndIf - + EndIf EndIf EndFor } } - { Name Hcurl_a_3D_ac; Type Assign; + { Name Hcurl_a_3D; Type Assign; Case { { Region Region[ Sur_Dirichlet_a ]; Value 0.; } } } - + { Name GaugeCondition_a; Type Assign; Case { { Region Vol_Tree; SubRegion Region[ { Sur_Tree, Lin_Tree } ]; - Value 0.; } + Value 0; } + } + } + + // Semi-analytic approach + { Name Impose_corr ; + Case { + For i In {1:NumWires} + { Region LWIRE~{i} ; Value 1. ; } + EndFor + } + } + + { Name Hcurl_acorr_3D; Type Assign; + Case { + { Region Region[ Sur_Dirichlet_a ]; Value 0.; } + If( Flag_Thin ) + For i In {1:NumWires} + { Region Region[ LWIRE~{i} ]; Value 1 ; } + EndFor + EndIf } } } @@ -256,7 +264,7 @@ FunctionSpace { } Constraint { { NameOfCoef ac; - EntityType Auto; NameOfConstraint Hcurl_a_3D_ac; } + EntityType Auto; NameOfConstraint Hcurl_a_3D; } { NameOfCoef ac; EntityType EdgesOfTreeIn ; EntitySubType StartingOn ; NameOfConstraint GaugeCondition_a ; } } @@ -280,10 +288,11 @@ FunctionSpace { } } + // Electric current per regionl for stranded conductors { Name Hregion_i_3D; Type Vector; BasisFunction { { Name sr; NameOfCoef ir; Function BF_RegionZ; - Support Dom_Hthin_a; Entity Vol_C; } + Support Dom_Hregion_i; Entity Vol_C; } } GlobalQuantity { { Name Is; Type AliasOf ; NameOfCoef ir; } @@ -295,52 +304,93 @@ FunctionSpace { } } - // For the local correction method - { Name Hcurl_acorr_3D; Type Form1; + // Function spaces for the semi-analytical approach + + // { Name Hcurl_acorr_3D; Type Form1; + // BasisFunction { + // { Name sw; NameOfCoef aw; Function BF_Edge; + // Support Dom_Hthin_a; Entity EdgesOf[ Vol_nu, ConnectedTo LWIRES ]; } + // } + // Constraint { + // { NameOfCoef aw; + // EntityType Auto; NameOfConstraint Hcurl_a_3D; } + // { NameOfCoef aw; EntityType EdgesOfTreeIn ; EntitySubType StartingOn ; + // NameOfConstraint GaugeCondition_a ; } + // } + // } + + { Name Hcurl_athin_3D; Type Form1; BasisFunction { - { Name sw; NameOfCoef aw; Function BF_Edge; - Support Vol_nu; Entity EdgesOf[ Vol_nu, ConnectedTo LWIRES ]; } + { Name si; NameOfCoef ai; Function BF_GroupOfEdges; + Support Dom_Hcurl_a; Entity GroupsOfEdgesOf[ LWIRES ]; } + { Name sj; NameOfCoef aj; Function BF_Edge; + Support Dom_Hcurl_a; Entity EdgesOf[ Vol_nu, Not LWIRES ]; } + } + GlobalQuantity { + { Name F; Type AliasOf ; NameOfCoef ai; } + { Name I; Type AssociatedWith; NameOfCoef ai; } } Constraint { - // { NameOfCoef aw; EntityType EdgesOfTreeIn ; EntitySubType StartingOn ; - // NameOfConstraint GaugeCondition_a ; } + { NameOfCoef I; EntityType Region; NameOfConstraint Impose_corr; } + { NameOfCoef aj; EntityType Auto; NameOfConstraint Hcurl_a_3D; } + { NameOfCoef aj; EntityType EdgesOfTreeIn ; EntitySubType StartingOn ; + NameOfConstraint GaugeCondition_a ; } + } + } + { Name Hcurl_asleeve_3D; Type Form1; + BasisFunction { + { Name si; NameOfCoef ai; Function BF_GroupOfEdges; + Support Dom_Hthin_a; Entity GroupsOfEdgesOf[ LWIRES ]; } + { Name sj; NameOfCoef aj; Function BF_Edge; + Support Dom_Hthin_a; Entity EdgesOf[ Vol_nu, ConnectedTo LWIRES ]; } + } + GlobalQuantity { + { Name F; Type AliasOf ; NameOfCoef ai; } + { Name I; Type AssociatedWith; NameOfCoef ai; } + } + Constraint { + { NameOfCoef I; EntityType Region; NameOfConstraint Impose_corr; } + { NameOfCoef aj; EntityType Auto; NameOfConstraint Hcurl_a_3D; } } } } 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; } - If( Flag_Form == 1) // Massive conductor - Quantity { - { Name a; Type Local; NameOfSpace Hcurl_a_3D; } - - { Name v; Type Local; NameOfSpace Hgrad_u_3D; } - { Name U; Type Global; NameOfSpace Hgrad_u_3D [U]; } - { Name I; Type Global; NameOfSpace Hgrad_u_3D [I]; } - } + { Name v; Type Local; NameOfSpace Hgrad_u_3D; } + { Name U; Type Global; NameOfSpace Hgrad_u_3D [U]; } + { Name I; Type Global; NameOfSpace Hgrad_u_3D [I]; } + } - Equation { - Integral { [ nu[] * Dof{d a} , {d a} ]; - In Vol_nu; Jacobian Vol; Integration I1; } + Equation { + Integral { [ nu[] * Dof{d a} , {d a} ]; + In Vol_nu; Jacobian Vol; Integration I1; } - Integral { DtDof[ sigma[] * Dof{a} , {a} ]; - In Vol_C; Jacobian Vol; Integration I1; } - 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; } - Integral { [ sigma[] * Dof{d v} , {d v} ]; - In Vol_C; Jacobian Vol; Integration I1; } + Integral { DtDof[ sigma[] * Dof{a} , {a} ]; + In Vol_C; Jacobian Vol; Integration I1; } + 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; } + Integral { [ sigma[] * Dof{d v} , {d v} ]; + In Vol_C; Jacobian Vol; Integration I1; } - GlobalTerm { [Dof{I} , {U} ]; In CATHODES; } - } + GlobalTerm { [Dof{I} , {U} ]; In CATHODES; } + } - Else // stranded conductor - If( Flag_Corr==0 ) // naive thin wire model + 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]; } @@ -351,7 +401,7 @@ Formulation { 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; } @@ -361,44 +411,53 @@ Formulation { In Vol_C; Jacobian Vol; Integration I1;} Integral { [ Dof{i}/(sigma[]*A_c) , {i} ]; In Vol_C; Jacobian Vol; Integration I1;} - GlobalTerm { [ Dof{Us}*A_c , {Is} ]; In Vol_C; } - } - Else // use local correction method - Quantity { - { Name a; Type Local; NameOfSpace Hcurl_a_3D; } - { Name as; Type Local; NameOfSpace Hcurl_acorr_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]; } + GlobalTerm { [ Dof{Us}*A_c , {Is} ]; In Vol_C; } } + EndIf + + EndIf - 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 { [ nu[] * Dof{d as} , {d as} ]; - In Vol_nu; Jacobian Vol; Integration I1; } - Integral { [ -Dof{i}/A_c , {as} ]; - In Vol_C; Jacobian Vol; Integration I1; } + If( Flag_SemiAnalytic ) + // 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]; } - /* - 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} ]; - In Vol_C;} - GlobalTerm { [ Dof{Us} , {Is} ]; - In Vol_C; } - */ - } - EndIf + { 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}*1e1 , {F} ]; + In Vol_nu; } + + Integral { [ nu[] * Dof{d as} , {d as} ]; + In Vol_nu; Jacobian Vol; Integration I1; } + GlobalTerm { [ -Dof{Is}*1e1 , {Fs} ]; + In Vol_nu; } + + // Integral { [ -Dof{i}/A_c , {as} ]; + // In Vol_C; Jacobian Vol; Integration I1; } + + /* + 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} ]; + In Vol_C;} + GlobalTerm { [ Dof{Us} , {Is} ]; + In Vol_C; } + */ + } EndIf } } @@ -443,41 +502,65 @@ PostProcessing { { Name b; Value { Local { [ {d a}]; In Dom_Hcurl_a; Jacobian Vol; }}} { Name by; - Value { Local { [ Re[Cart2Pol[CompY [ {d a}]]]]; + Value { Local { [ Cart2Pol[ Norm[ {d a}]] ]; In Dom_Hcurl_a; Jacobian Vol; }}} - - If( Flag_Form == 1 ) // massive - { Name J; - Value { Local { [ -sigma[] * ( Dt[{a}] + {d v} ) ]; - In Vol_C; Jacobian Vol; }}} - { Name U; - Value { Term { [ {U} ]; In Electrodes; }}} + + If( Flag_SemiAnalytic ) + { Name bs; + Value { Local { [ {d as} ]; In Vol_nu; Jacobian Vol; }}} + { Name bsy; + Value { Local { [ Norm[ {d as} ] ]; In Vol_nu; Jacobian Vol; }}} + { 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[] ]; + In Vol_nu; Jacobian Vol; }}} + EndIf + + + + If( !Flag_SemiAnalytic ) + If( !Flag_Stranded ) // massive + { Name J; + Value { Local { [ -sigma[] * ( Dt[{a}] + {d v} ) ]; + In Vol_C; Jacobian Vol; }}} + { Name U; + Value { Term { [ {U} ]; In Electrodes; }}} + { Name I; + Value { Term { [ {I} ]; In Electrodes; }}} + { Name R; + Value { Term { [ Re[ -{U}/{I}/Lz ] ]; In Electrodes; }}} + { Name L; // actually total flux/I + Value { Term { [ Im[ -{U}/{I}/omega/Lz ] ]; In Electrodes; }}} + Else // stranded + { Name J; + Value { Local { [ {Is}/A_c ]; + In Vol_C; Jacobian Vol; }}} + { Name U; + Value { Term { [ {Us} ]; In Vol_C; }}} + { Name I; + Value { Term { [ {Is} ]; In Vol_C; }}} + { Name R; + Value { Term { [ Re[ -{Us}/{Is}/Lz ] ]; In Vol_C; }}} + { Name L; // actually total flux/I + Value { Term { [ Im[ -{Us}/{Is}/omega/Lz ] ]; In Vol_C; }}} + EndIf + EndIf + + If( Flag_SemiAnalytic ) + { Name F; + Value { Term { [ {F} ]; In Vol_C; }}} { Name I; - Value { Term { [ {I} ]; In Electrodes; }}} - { Name R; - Value { Term { [ Re[ -{U}/{I}/Lz ] ]; In Electrodes; }}} - { Name L; // actually total flux/I - Value { Term { [ Im[ -{U}/{I}/omega/Lz ] ]; In Electrodes; }}} - Else // stranded + 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 { [ {Us} ]; In Vol_C; }}} - { Name I; - Value { Term { [ {Is} ]; In Vol_C; }}} + Value { Term { [ i[]*omega*{F}*10 + Analytic_R[]*{I} ]; In Vol_C; }}} { Name R; - Value { Term { [ Re[ -{Us}/{Is}/Lz ] ]; In Vol_C; }}} + Value { Term { [ Analytic_R[] + {I}*0 ]; In Vol_C; }}} { Name L; // actually total flux/I - Value { Term { [ Im[ -{Us}/{Is}/omega/Lz ] ]; In Vol_C; }}} - EndIf - - If( (Flag_Form == 2) && Flag_Corr ) - { Name bs; - Value { Local { [ {d as}]; In Dom_Hthin_a; Jacobian Vol; }}} - { Name bcorry; - Value { Local { [ Re[Cart2Pol[CompY [{d a}-{d as}]]] + Correction_B[] ]; - In Dom_Hthin_a; Jacobian Vol; }}} + Value { Term { [ {F}/{I}/Lz ]; In Vol_C; }}} EndIf } @@ -485,28 +568,48 @@ PostProcessing { } PostOperation map UsingPost MagnetoDynamics { - Print[ J, OnElementsOf Vol_C, File "j.pos"]; + Print[ b, OnElementsOf Vol_nu, File "b.pos"]; - If( (Flag_Form == 2) && Flag_Corr ) + + If( !Flag_SemiAnalytic ) + Print[ J, OnElementsOf Vol_C, File "j.pos"]; + Else Print[ bs, OnElementsOf Vol_nu, File "bs.pos"]; PrintGroup[ EdgesOf[ Vol_nu, ConnectedTo LWIRES ], In Vol_nu, File "group.pos"]; EndIf + PrintGroup[ EdgesOfTreeIn[ { Vol_Tree }, StartingOn { Sur_Tree, Lin_Tree } ], + In Vol_Tree, File "Tree.pos"]; + //PrintGroup[ _CO_Entity_44, In Vol_nu, File "Tree.pos"]; + } PostOperation integaz UsingPost MagnetoDynamics { For i In {1:NumWires} - If( Flag_Form == 1 ) // massive - - Print[ U, OnRegion CATHODES, File > "U.dat", Format Table, - SendToServer Sprintf["Results/1Voltage/Wire %g",i]]; - Print[ I, OnRegion CATHODES, File > "I.dat", Format Table, - SendToServer Sprintf["Results/2Current/Wire %g",i]]; - Print[ R, OnRegion CATHODES, File > "Z.dat", Format Table, - SendToServer Sprintf["Results/3Resistance p.u.l./Wire %g",i]]; - Print[ L, OnRegion CATHODES, File > "Z.dat", Format Table, - SendToServer Sprintf["Results/4Inductance p.u.l./Wire %g",i]]; - Else + If( !Flag_SemiAnalytic ) + If( !Flag_Stranded ) // massive + Print[ U, OnRegion CATHODES, File > "U.dat", Format Table, + SendToServer Sprintf["Results/1Voltage/Wire %g",i]]; + Print[ I, OnRegion CATHODES, File > "I.dat", Format Table, + SendToServer Sprintf["Results/2Current/Wire %g",i]]; + Print[ R, OnRegion CATHODES, File > "Z.dat", Format Table, + SendToServer Sprintf["Results/3Resistance p.u.l./Wire %g",i]]; + Print[ L, OnRegion CATHODES, File > "Z.dat", Format Table, + SendToServer Sprintf["Results/4Inductance p.u.l./Wire %g",i]]; + Else + Print[ U, OnRegion Vol_C, File > "U.dat", Format Table, + SendToServer Sprintf["Results/1Voltage/Wire %g",i]]; + Print[ I, OnRegion Vol_C, File > "I.dat", Format Table, + SendToServer Sprintf["Results/2Current/Wire %g",i]]; + Print[ R, OnRegion Vol_C, File > "Z.dat", Format Table, + SendToServer Sprintf["Results/3Resistance p.u.l./Wire %g",i]]; + Print[ L, OnRegion Vol_C, File > "Z.dat", Format Table, + SendToServer Sprintf["Results/4Inductance p.u.l./Wire %g",i]]; + EndIf + EndIf + If( Flag_SemiAnalytic ) + Print[ F, OnRegion Vol_C, File > "Flux.dat", Format Table, + SendToServer Sprintf["Results/0Flux/Wire %g",i]]; Print[ U, OnRegion Vol_C, File > "U.dat", Format Table, SendToServer Sprintf["Results/1Voltage/Wire %g",i]]; Print[ I, OnRegion Vol_C, File > "I.dat", Format Table, @@ -519,27 +622,66 @@ PostOperation integaz UsingPost MagnetoDynamics { EndFor } -NbPoints = 400; -Xmax = 0.15*box; +NbPoints = 1000; +Xcut = 0.1*box; +Zcut = Lz/2; PostOperation cut UsingPost MagnetoDynamics { - Print [ by, OnLine { {-Xmax,0,0} {Xmax,0,0} }{NbPoints}, - Format Gmsh, File "by.pos" ]; + Print [ by, OnLine { {-Xcut,0,0} {Xcut,0,Zcut} }{NbPoints}, + Format Gmsh, File "by.pos" ]; + Echo[ StrCat["l=PostProcessing.NbViews-1;", + "View[l].Name = 'cut by';", + "View[l].Axes = 3;", + "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" ]; Echo[ StrCat["l=PostProcessing.NbViews-1;", - "View[l].Name = 'cut by';", - "View[l].Axes = 3;", - "View[l].LineWidth = 3;", - "View[l].Type = 2;"], - File "tmp.geo", LastTimeStepOnly]; - If( (Flag_Form == 2) && Flag_Corr ) - Print [ bcorry, OnLine { {-Xmax,0,0} {Xmax,0,0} }{NbPoints}, - Format Gmsh, File "by.pos" ]; + "View[l].Name = 'cut bsy ';", + "View[l].Axes = 3;", + "View[l].LineWidth = 3;", + "View[l].Type = 2;"], + File "tmp.geo", LastTimeStepOnly]; + + Print [ bbsy, OnLine { {-Xcut,0,0} {Xcut,0,Zcut} }{NbPoints}, + Format Gmsh, File "bcorry.pos" ]; Echo[ StrCat["l=PostProcessing.NbViews-1;", - "View[l].Name = 'cut by';", + "View[l].Name = 'cut by-bsy ';", "View[l].Axes = 3;", "View[l].LineWidth = 3;", "View[l].Type = 2;"], File "tmp.geo", LastTimeStepOnly]; + + Print [ bcorry, OnLine { {-Xcut,0,0} {Xcut,0,Zcut} }{NbPoints}, + Format Gmsh, File "bcorry.pos" ]; + Echo[ StrCat["l=PostProcessing.NbViews-1;", + "View[l].Name = 'cut by corrected';", + "View[l].Axes = 3;", + "View[l].LineWidth = 3;", + "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" ]; + Else + Print [ by, OnLine { {0,0,0} {Xcut,0,0} } {NbPoints}, + Format SimpleTable, File "Cut_by.txt" ]; + Print [ bsy, OnLine { {0,0,0} {Xcut,0,0} } {NbPoints}, + Format SimpleTable, File "Cut_bsy.txt" ]; + Print [ bbsy, OnLine { {0,0,0} {Xcut,0,0} } {NbPoints}, + Format SimpleTable, File "Cut_bbsy.txt" ]; + Print [ bcorry, OnLine { {0,0,0} {Xcut,0,0} } {NbPoints}, + Format SimpleTable, File "Cut_bcorry.txt" ]; EndIf } diff --git a/Conductors3D/w3d_common.pro b/Conductors3D/w3d_common.pro index e8dc467..3f34eed 100644 --- a/Conductors3D/w3d_common.pro +++ b/Conductors3D/w3d_common.pro @@ -4,32 +4,43 @@ deg = 180/Pi; /* Expected inductance of a rectilinear 1mm radius wire: 4 nH/cm + + git/documentation/models/Inductor/magstadyn_av_js0_3d.pro */ DefineConstant[ - NumWires = {1, Name "Parameters/0Number of wires", + NumWires = {1, Name "Parameters/00Number of wires", Choices {1="1", 2="2", 3="3"}, Visible 1} - Flag_Form = {2, Name "Parameters/1Conductor type", Visible 1, - Choices {1="Massive", 2="Stranded"}} - Flag_Thin = {1, Name "Parameters/2Use thin wires", Choices {0,1}, Visible 1} - Flag_Corr = {1, Name "Parameters/3Correction method", - Visible (Flag_Form==2) && Flag_Thin, - Choices {0="None", 1="Raw sleeve", 2="regular sleeve"}} - Flag_U = {0, Name "Parameters/4Impose U", Choices {0,1}, Visible 1} + + Flag_Thin = {1, Name "Parameters/01Thin wires", + Choices {0,1}, Visible 1} + Flag_RegSleeve = {0, Name "Parameters/02Regular sleeves", + Choices {0,1}, Visible Flag_Thin} + + Flag_SemiAnalytic = {1, Name "Parameters/03Semi-analytic approach", + Choices {0,1}, Visible Flag_Thin} + + Flag_Stranded = {0, Name "Parameters/04Stranded conductors", + Choices {0,1}, Visible !Flag_SemiAnalytic || !Flag_Thin} + Flag_Dual = {0, Name "Parameters/05Dual approach", + Choices {0,1}, Visible Flag_SemiAnalytic && Flag_Thin} + Flag_U = {0, Name "Parameters/06Impose U", Choices {0,1}, Visible 1} + + WireRadius = {1, Name "Parameters/10Wire radius [mm]"} + MeshSizeWire = {0.25, Name "Parameters/11Mesh size in wire [mm]", + Visible !Flag_Thin} + MeshSizeThinWire = {2, Name "Parameters/12Mesh size on thin wire [mm]", + Visible Flag_Thin} + //SleeveRadius = {2, Name "Parameters/13Sleeve radius [mm]", Visible Flag_RegSleeve} + box = {100*mm, Name "Parameters/5box half-width [m]"} + Flag_Maps = {1, Name "Input/2Display maps", Choices {0,1}, Visible 1} LogFreq = {4, Min 1, Max 8, Step 0.25, Name "Input/4Log of frequency"} Freq = 10^LogFreq - WireRadius = {1, Name "Parameters/5Wire radius [mm]"} - MeshSizeWire = {1, Name "Parameters/6Mesh size in wire [mm]", - Visible !Flag_Thin} - MeshSizeThinWire = {1, Name "Parameters/7Mesh size on thin wire [mm]", - Visible Flag_Thin} - SleeveRadius = {2, Name "Parameters/8Sleeve radius [mm]"} - box = {100*mm, Name "Parameters/9box half-width [m]"} - + rw = WireRadius*mm - rs = SleeveRadius*mm + rs = MeshSizeThinWire*mm A_c = Pi*rw^2 // wire cross section rout = box Lz = 10*mm @@ -41,6 +52,10 @@ DefineConstant[ */ ]; +If( Flag_SemiAnalytic ) + SetNumber("Parameters/01Thin wires", 1); + Flag_Thin=1; +EndIf Xlocs() = { 0, 5*mm, -5*mm}; @@ -50,7 +65,7 @@ For i In {1:NumWires} WX~{i} = { Xlocs(i-1), // place wires symmetrically around x=0 Name Sprintf("Parameters/Wire %g/0X position [m]", i), Closed } WY~{i} = { Ylocs(i-1), Name Sprintf("Parameters/Wire %g/0Y position [m]", i) } - WI~{i} = { 1.23456, Name Sprintf("Parameters/Wire %g/2Current [A]", i) } + WI~{i} = { 1, Name Sprintf("Parameters/Wire %g/2Current [A]", i) } WP~{i} = { 0*NumWires*(i-1), Name Sprintf("Parameters/Wire %g/3Phase [deg]", i) } ]; -- GitLab