From 2c1a7a86e143adf7edde26c3e8de15595870e880 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Fran=C3=A7ois=20Henrotte?= <francois.henrotte@uclouvain.be>
Date: Mon, 30 Mar 2020 16:25:01 +0200
Subject: [PATCH] update
---
Conductors3D/w2d.pro | 475 ++++++++++++++----------------------
Conductors3D/w3d.geo | 8 +-
Conductors3D/w3d.pro | 25 +-
Conductors3D/w3d_common.pro | 4 +-
4 files changed, 206 insertions(+), 306 deletions(-)
diff --git a/Conductors3D/w2d.pro b/Conductors3D/w2d.pro
index f43b73e..b09e70d 100644
--- a/Conductors3D/w2d.pro
+++ b/Conductors3D/w2d.pro
@@ -36,13 +36,13 @@ Group{
LWIRES = Region[ {} ];
For i In {1:NumWires}
LWIRE~{i} = Region[ (50+i) ];
- LWIRES += Region[ LWIRE~{i} ];
- If( Flag_Thin == 0 )
+ LWIRES += Region[ LWIRE~{i} ];
+ If( Flag_Thin == 0 )
VWIRE~{i} = Region[ (10+i) ];
VWIRES += Region[ VWIRE~{i} ];
- SWIRE~{i} = Region[ {} ];
- Else
- VWIRE~{i} = Region[ {} ];
+ SWIRE~{i} = Region[ {} ];
+ Else
+ VWIRE~{i} = Region[ {} ];
SWIRE~{i} = ElementsOf[ AIR, OnOneSideOf LWIRE~{i} ];
EndIf
EndFor
@@ -100,9 +100,9 @@ Function{
For i In {1:NumWires}
// Current density
J[ Region[ {VWIRE~{i}, LWIRE~{i}} ] ] =
- Vector[ 0, 0, WI~{i}/A_c * Exp[ i[]*WP~{i}*deg ] ];
+ Vector[ 0, 0, WI~{i}/A_c * Exp[ i[]*WP~{i}*deg ] ];
- // local radial coordinate with origin on wire i
+ // local radial coordinate with origin on wire i
R~{i}[] = Hypot[ X[]-WX~{i} , Y[]-WY~{i} ];
EndFor
@@ -111,10 +111,10 @@ Function{
//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) );
+ (($1>rw) ? mu0/(2*Pi) * Log[rs/$1] :
+ 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 );
+ (($1>rw) ? Log[rs/$1] : Log[rs/rw] + mur*(1-($1/rw)^2)/2 );
// Impedance of thin wire p.u. length
R_DC = 1./(sigma*A_c);
@@ -136,7 +136,7 @@ Function{
+ WI_3*Analytic_A[R_3[]] + (WI_1+WI_2+WI_3)*mu0/(2*Pi) * NumWires * Log[rout/rs] ;
Correction_A[] = ((R_1[]<rs) ? WI_1*Analytic_A[R_1[]] :
((R_2[]<rs) ? WI_2*Analytic_A[R_2[]] :
- ((R_3[]<rs) ? WI_3*Analytic_A[R_3[]] : 0 )));
+ ((R_3[]<rs) ? WI_3*Analytic_A[R_3[]] : 0 )));
EndIf
}
@@ -155,14 +155,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 ; }
}
}
}
@@ -182,7 +182,7 @@ Constraint{
Case {
For i In {1:NumWires}
{ Region VWIRE~{i} ; Value WI~{i} ; }
- { Region LWIRE~{i} ; Value WI~{i} ; }
+ { Region LWIRE~{i} ; Value WI~{i} ; }
EndFor
}
}
@@ -206,9 +206,9 @@ FunctionSpace {
{ Name Hcurl_a_2D; Type Form1P;
BasisFunction {
{ Name sc; NameOfCoef ac; Function BF_PerpendicularEdge;
- Support Dom_Hcurl_a; Entity NodesOf[All]; }
+ Support Dom_Hcurl_a; Entity NodesOf[All]; }
{ Name sw; NameOfCoef aw; Function BF_PerpendicularEdge;
- Support Dom_Hthin_a; Entity NodesOf[lVol_C]; }
+ Support Dom_Hthin_a; Entity NodesOf[lVol_C]; }
}
SubSpace {
{ Name ac ; NameOfBasisFunction { sc } ; }
@@ -216,7 +216,7 @@ FunctionSpace {
}
Constraint {
{ NameOfCoef ac;
- EntityType Auto; NameOfConstraint Hcurl_a_2D_ac; }
+ EntityType Auto; NameOfConstraint Hcurl_a_2D_ac; }
}
}
@@ -253,9 +253,9 @@ FunctionSpace {
{ Name Hcurl_a_3D; Type Form1;
BasisFunction {
{ Name sc; NameOfCoef ac; Function BF_Edge;
- Support Dom_Hcurl_a; Entity EdgesOf[All]; }
+ Support Dom_Hcurl_a; Entity EdgesOf[All]; }
{ Name sw; NameOfCoef aw; Function BF_Edge;
- Support Dom_Hthin_a; Entity EdgesOf[lVol_C]; }
+ Support Dom_Hthin_a; Entity EdgesOf[lVol_C]; }
}
SubSpace {
{ Name ac ; NameOfBasisFunction { sc } ; }
@@ -263,16 +263,16 @@ FunctionSpace {
}
Constraint {
{ NameOfCoef ac;
- EntityType Auto; NameOfConstraint Hcurl_a_3D_ac; }
+ EntityType Auto; NameOfConstraint Hcurl_a_3D_ac; }
}
}
{ Name Hgrad_u_3D; Type Form0;
BasisFunction {
{ Name sn; NameOfCoef un; Function BF_Node;
- Support Dom_Hgrad_u; Entity NodesOf[ Vol_C, Not Electrodes]; }
+ Support Dom_Hgrad_u; Entity NodesOf[ Vol_C, 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; }
@@ -287,7 +287,6 @@ FunctionSpace {
Formulation {
{ Name MagnetoDynamics; Type FemEquation;
- If( !Flag_3D )
Quantity {
{ Name a; Type Local; NameOfSpace Hcurl_a_2D; }
{ Name ac; Type Local; NameOfSpace Hcurl_a_2D[ac]; }
@@ -331,56 +330,19 @@ Formulation {
GlobalTerm { [ Dof{I} , {U} ]; In Vol_C; }
GlobalTerm { [ Dof{Us} , {Is} ];
- In lVol_C; }
- GlobalTerm { [ Analytic_R[] * Dof{Is} , {Is} ];
- In lVol_C;}
+ In lVol_C; }
+ GlobalTerm { [ Analytic_R[] * Dof{Is} , {Is} ];
+ In lVol_C;}
Integral { DtDof [ Dof{ac}/A_c , {i} ];
In lVol_C; Jacobian Vol; Integration I1; }
+
If( Flag_Thin != 3)
Integral { DtDof [ -Dof{aw}/A_c , {i} ];
- In lVol_C; Jacobian Vol; Integration I1;}
- GlobalTerm { DtDof [ Analytic_L[] * Dof{Is} , {Is} ];
- In lVol_C;}
+ In lVol_C; Jacobian Vol; Integration I1;}
+ GlobalTerm { DtDof [ Analytic_L[] * Dof{Is} , {Is} ];
+ In lVol_C;}
EndIf
}
- Else // ############### 3D
- Quantity {
- { Name a; Type Local; NameOfSpace Hcurl_a_3D; }
- { Name ac; Type Local; NameOfSpace Hcurl_a_3D[ac]; }
- { Name aw; Type Local; NameOfSpace Hcurl_a_3D[aw]; }
- // { 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 ac} , {d ac} ];
- In Vol_nu; Jacobian Vol; Integration I1; }
- Integral { [ nu[] * Dof{d aw} , {d aw} ];
- In Vol_nu; Jacobian Vol; Integration I1; }
-
- Integral { [ -Vector[ 0, 0, I0/A_c] , {ac} ];
- In lVol_C; Jacobian Vol; Integration I1; }
- Integral { [ -Vector[ 0, 0, I0/A_c] , {aw} ];
- In lVol_C; Jacobian Vol; Integration I1; }
-
- Integral { [ -Vector[ 0, 0, I0/A_c] , {a} ];
- 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{U} , {I} ]; In Electrodes; }
- */
-
- }
- EndIf
}
}
@@ -388,126 +350,125 @@ Resolution {
{ Name Dynamic;
System {
{ Name S; NameOfFormulation MagnetoDynamics;
- Type ComplexValue; Frequency Freq; }
+ Type ComplexValue; Frequency Freq; }
}
Operation {
CreateDir[DataDir];
// SOLVER OPTION
If(type_solver == S::SOLVER.gmres)
- SetGlobalSolverOptions[
- "-pc_type ilu -ksp_type gmres -ksp_max_it 1000 -ksp_rtol 1e-6 -ksp_atol 1e-6"];
- ElseIf(type_solver == S::SOLVER.fgmres)
- SetGlobalSolverOptions[
- "-pc_type ilu -ksp_type fgmres -ksp_max_it 1000 -ksp_rtol 1e-6 -ksp_atol 1e-6"];
- ElseIf(type_solver == S::SOLVER.bcgs)
- SetGlobalSolverOptions[
- "-pc_type ilu -ksp_type bcgs -ksp_max_it 1000 -ksp_rtol 1e-6 -ksp_atol 1e-6"];
+ SetGlobalSolverOptions[
+ "-pc_type ilu -ksp_type gmres -ksp_max_it 1000 -ksp_rtol 1e-6 -ksp_atol 1e-6"];
+ ElseIf(type_solver == S::SOLVER.fgmres)
+ SetGlobalSolverOptions[
+ "-pc_type ilu -ksp_type fgmres -ksp_max_it 1000 -ksp_rtol 1e-6 -ksp_atol 1e-6"];
+ ElseIf(type_solver == S::SOLVER.bcgs)
+ SetGlobalSolverOptions[
+ "-pc_type ilu -ksp_type bcgs -ksp_max_it 1000 -ksp_rtol 1e-6 -ksp_atol 1e-6"];
EndIf
- SetGlobalSolverOptions["-petsc_prealloc 900"];
-
+ SetGlobalSolverOptions["-petsc_prealloc 900"];
+
Generate[S]; Solve[S]; SaveSolution[S];
- PostOperation[integaz];
-
- If( Flag_Maps == 1 )
- PostOperation[mapaz];
- PostOperation[cutaz];
-
- DeleteFile["U_full.dat"];
- DeleteFile["U_raw.dat"];
- DeleteFile["U_reg.dat"];
- DeleteFile["U_naive.dat"];
- DeleteFile["joule_full.dat"];
- DeleteFile["joule_raw.dat"];
- DeleteFile["joule_reg.dat"];
- DeleteFile["joule_naive.dat"];
- EndIf
-
- If( Flag_Thin != 0 )
- For i In {1:NumWires}
- 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
- EndIf
+ PostOperation[integaz];
+
+ If( Flag_Maps == 1 )
+ PostOperation[mapaz];
+ PostOperation[cutaz];
+
+ DeleteFile["U_full.dat"];
+ DeleteFile["U_raw.dat"];
+ DeleteFile["U_reg.dat"];
+ DeleteFile["U_naive.dat"];
+ DeleteFile["joule_full.dat"];
+ DeleteFile["joule_raw.dat"];
+ DeleteFile["joule_reg.dat"];
+ DeleteFile["joule_naive.dat"];
+ EndIf
+
+ If( Flag_Thin != 0 )
+ For i In {1:NumWires}
+ 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" ] ;
+ //Print[ {$sarea~{i}}, Format "DEBUG %7.3e" ] ;
+ EndFor
+ EndIf
}
- }
+ }
}
PostProcessing {
{ Name PostProcessings; NameOfFormulation MagnetoDynamics;
PostQuantity {
{ Name az;
- Value { Local { [ CompZ[{ac}] - (Flag_Thin!=0)*(CompZ[{aw}]-Correction_A[])];
- In Dom_Hcurl_a; Jacobian Vol; }}}
+ Value { Local { [ CompZ[{ac}] - (Flag_Thin!=0)*(CompZ[{aw}]-Correction_A[])];
+ In Dom_Hcurl_a; Jacobian Vol; }}}
{ Name acz;
- Value { Local { [ CompZ[ {ac}] ];
- In Dom_Hcurl_a; Jacobian Vol; }}}
+ Value { Local { [ CompZ[ {ac}] ];
+ In Dom_Hcurl_a; Jacobian Vol; }}}
{ Name awz;
- Value { Local { [ CompZ[ {aw}] ];
- In Dom_Hthin_a; Jacobian Vol; }}}
+ Value { Local { [ CompZ[ {aw}] ];
+ In Dom_Hthin_a; Jacobian Vol; }}}
{ Name acwz;
- Value { Local { [ CompZ[ {ac} - {aw} ] ];
- In Dom_Hcurl_a; Jacobian Vol; }}}
+ Value { Local { [ CompZ[ {ac} - {aw} ] ];
+ In Dom_Hcurl_a; Jacobian Vol; }}}
{ Name exact;
- Value { Local { [ Exact_A[] ] ;
- In Dom_Hcurl_a; Jacobian Vol; }}}
+ Value { Local { [ Exact_A[] ] ;
+ In Dom_Hcurl_a; Jacobian Vol; }}}
{ Name corr;
- Value { Local { [ Correction_A[] ] ;
- In Dom_Hcurl_a; Jacobian Vol; }}}
+ Value { Local { [ Correction_A[] ] ;
+ In Dom_Hcurl_a; Jacobian Vol; }}}
{ Name bcw;
- Value { Local { [ {d ac} - {d aw} ];
- In Dom_Hcurl_a; Jacobian Vol; }}}
+ Value { Local { [ {d ac} - {d aw} ];
+ In Dom_Hcurl_a; Jacobian Vol; }}}
{ Name j;
- Value { Local { [ -sigma[] * ( Dt[{a}] + {dv} ) ];
- In Vol_C; Jacobian Vol; }}}
-
- { Name un;
- Value{ Local { [ 1 ] ;
- In Vol_nu; Integration I1 ; Jacobian Vol; }}}
- { Name surf;
- Value{ Integral{ [ 1 ] ;
- In Vol_nu; Integration I1 ; Jacobian Vol; }}}
-
-
- If( Flag_Thin == 0 )
- { Name flux;
- Value { Integral { [ CompZ[ {a} ] / A_c ];
- In Vol_C; Integration I1 ; Jacobian Vol; }}}
+ Value { Local { [ -sigma[] * ( Dt[{a}] + {dv} ) ];
+ In Vol_C; Jacobian Vol; }}}
+
+ { Name un;
+ Value{ Local { [ 1 ] ;
+ In Vol_nu; Integration I1 ; Jacobian Vol; }}}
+ { Name surf;
+ Value{ Integral{ [ 1 ] ;
+ In Vol_nu; Integration I1 ; Jacobian Vol; }}}
+
+ If( Flag_Thin == 0 )
+ { Name flux;
+ Value { Integral { [ CompZ[ {a} ] / A_c ];
+ In Vol_C; Integration I1 ; Jacobian Vol; }}}
{ Name JouleLosses;
- Value { Integral { [ sigma[]*SquNorm[ Dt[{a}] + {dv} ] ] ;
- In Vol_C ; Integration I1 ; Jacobian Vol ; }}}
- { Name U;
+ Value { Integral { [ sigma[]*SquNorm[ Dt[{a}] + {dv} ] ] ;
+ In Vol_C ; Integration I1 ; Jacobian Vol ; }}}
+ { Name U;
Value { Term { [ Cart2Pol [ {U} ] ]; In Vol_C; }}}
- { Name L;
- Value { Term { [ Im [ {U}/{I}/omega ] ]; In Vol_C; }}}
+ { Name L;
+ Value { Term { [ Im [ {U}/{I}/omega ] ]; In Vol_C; }}}
- Else
- { Name flux; // naive flux
- Value { Integral { [ CompZ[ {ac} ] ];
- In lVol_C; Integration I1 ; Jacobian Vol; }}}
+ Else
+ { Name flux; // naive flux
+ Value { Integral { [ CompZ[ {ac} ] ];
+ In lVol_C; Integration I1 ; Jacobian Vol; }}}
- { Name JouleLosses;
- Value { Term { [ Analytic_R[]*{Is}*{Is} ] ; In lVol_C; }}}
- { Name U;
+ { Name JouleLosses;
+ Value { Term { [ Analytic_R[]*{Is}*{Is} ] ; In lVol_C; }}}
+ { Name U;
Value { Term { [ Cart2Pol [ {Us} ] ] ; In lVol_C; }}}
{ Name L; // only valid with one current
Value { Term { [ Im [ {Us}/{Is}/omega ] ] ; In lVol_C; }}}
-
- For i In {1:NumWires}
- { Name area~{i};
- Value{ Integral{ [ 1 ] ;
- In SWIRE~{i}; Integration I1 ; Jacobian Vol; }}}
- { Name intbx~{i}; // integral of b, 1st step for computing b average
- Value { Integral { [ Norm [ CompX[ {d ac} - {d aw} ] ] ] ;
- In SWIRE~{i}; Integration I1 ; Jacobian Vol; }}}
- { Name intby~{i}; // integral of b, 1st step for computing b average
- Value { Integral { [ Norm [ CompY[ {d ac} - {d aw} ] ] ] ;
- In SWIRE~{i}; Integration I1 ; Jacobian Vol; }}}
- EndFor
-
- EndIf
+
+ For i In {1:NumWires}
+ { Name area~{i};
+ Value{ Integral{ [ 1 ] ;
+ In SWIRE~{i}; Integration I1 ; Jacobian Vol; }}}
+ { Name intbx~{i}; // integral of b, 1st step for computing b average
+ Value { Integral { [ Norm [ CompX[ {d ac} - {d aw} ] ] ] ;
+ In SWIRE~{i}; Integration I1 ; Jacobian Vol; }}}
+ { Name intby~{i}; // integral of b, 1st step for computing b average
+ Value { Integral { [ Norm [ CompY[ {d ac} - {d aw} ] ] ] ;
+ In SWIRE~{i}; Integration I1 ; Jacobian Vol; }}}
+ EndFor
+ EndIf
}
}
}
@@ -515,22 +476,22 @@ PostProcessing {
PostOperation mapaz UsingPost PostProcessings {
Print [az, OnElementsOf Vol_nu, File "az.pos"];
Echo[ StrCat["l=PostProcessing.NbViews-1;",
- "View[l].IntervalsType = 3;",
- "View[l].NbIso = 30;",
- "View[l].RaiseZ = 15000;",
- "View[l].NormalRaise = 0;"],
- File "tmp.geo", LastTimeStepOnly];
+ "View[l].IntervalsType = 3;",
+ "View[l].NbIso = 30;",
+ "View[l].RaiseZ = 15000;",
+ "View[l].NormalRaise = 0;"],
+ File "tmp.geo", LastTimeStepOnly];
If(Flag_Thin != 0)
Print [bcw, OnElementsOf Vol_Sleeve, File "b.pos"];
Print [acwz, OnElementsOf Vol_nu, File "acz.pos"];
- Echo[ StrCat["l=PostProcessing.NbViews-1;",
- "View[l].IntervalsType = 3;",
- "View[l].NbIso = 30;",
- "View[l].RaiseZ = 15000;",
- "View[l].NormalRaise = 0;"],
- File "tmp.geo", LastTimeStepOnly] ;
+ Echo[ StrCat["l=PostProcessing.NbViews-1;",
+ "View[l].IntervalsType = 3;",
+ "View[l].NbIso = 30;",
+ "View[l].RaiseZ = 15000;",
+ "View[l].NormalRaise = 0;"],
+ File "tmp.geo", LastTimeStepOnly] ;
EndIf
}
@@ -542,53 +503,53 @@ PostOperation cutaz UsingPost PostProcessings {
// cut of az displayed in GUI
If( Flag_Thin == 0 )
Print [ az, OnLine { {-Xmax,0,0} {Xmax,0,0} }{NbPoints},
- Format Gmsh, File "Cut_az.pos" ];
+ Format Gmsh, File "Cut_az.pos" ];
Echo[ StrCat["l=PostProcessing.NbViews-1;",
- "View[l].Name = 'cut az full';",
- "View[l].Axes = 3;",
- "View[l].LineWidth = 3;",
- "View[l].Type = 2;"],
- File "tmp.geo", LastTimeStepOnly];
+ "View[l].Name = 'cut az full';",
+ "View[l].Axes = 3;",
+ "View[l].LineWidth = 3;",
+ "View[l].Type = 2;"],
+ File "tmp.geo", LastTimeStepOnly];
- Else
+ Else
- Print [ acz, OnLine { {-Xmax,0,0} {Xmax,0,0} }{NbPoints},
- Format Gmsh, File "Cut_az.pos" ];
+ Print [ acz, OnLine { {-Xmax,0,0} {Xmax,0,0} }{NbPoints},
+ Format Gmsh, File "Cut_az.pos" ];
Echo[ StrCat["l=PostProcessing.NbViews-1;",
- "View[l].Name = 'cut az not corrected';",
- "View[l].Axes = 3;",
- "View[l].LineWidth = 3;",
- "View[l].Type = 2;"],
- File "tmp.geo", LastTimeStepOnly];
+ "View[l].Name = 'cut az not corrected';",
+ "View[l].Axes = 3;",
+ "View[l].LineWidth = 3;",
+ "View[l].Type = 2;"],
+ File "tmp.geo", LastTimeStepOnly];
Print [ az, OnLine { {-Xmax,0,0} {Xmax,0,0} }{NbPoints},
- Format Gmsh, File "Cut_az.pos" ];
+ Format Gmsh, File "Cut_az.pos" ];
Echo[ StrCat["l=PostProcessing.NbViews-1;",
- "View[l].Name = 'cut az corrected';",
- "View[l].Axes = 3;",
- "View[l].LineWidth = 3;",
- "View[l].Type = 2;"],
- File "tmp.geo", LastTimeStepOnly];
+ "View[l].Name = 'cut az 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} {Xmax,0,0} } {NbPoints},
- Format SimpleTable, File "Cut_analytic.txt" ];
+ Format SimpleTable, File "Cut_analytic.txt" ];
If( Flag_Thin == 0 )
Print [ az, OnLine { {0,0,0} {Xmax,0,0} } {NbPoints},
- Format SimpleTable, File "Cut_fullfem.txt" ];
+ Format SimpleTable, File "Cut_fullfem.txt" ];
Else
Print [ az, OnLine { {0,0,0} {Xmax,0,0} } {NbPoints},
- Format SimpleTable, File "Cut_thinfem.txt" ];
+ Format SimpleTable, File "Cut_thinfem.txt" ];
Print [ acz, OnLine { {0,0,0} {Xmax,0,0} } {NbPoints},
- Format SimpleTable, File "Cut_ac.txt" ];
+ Format SimpleTable, File "Cut_ac.txt" ];
Print [ awz, OnLine { {0,0,0} {Xmax,0,0} } {NbPoints},
- Format SimpleTable, File "Cut_aw.txt" ];
+ Format SimpleTable, File "Cut_aw.txt" ];
Print [ acwz, OnLine { {0,0,0} {Xmax,0,0} } {NbPoints},
- Format SimpleTable, File "Cut_acw.txt" ];
+ Format SimpleTable, File "Cut_acw.txt" ];
Print [ corr, OnLine { {0,0,0} {Xmax,0,0} } {NbPoints},
- Format SimpleTable, File "Cut_correction.txt" ];
+ Format SimpleTable, File "Cut_correction.txt" ];
EndIf
}
@@ -597,111 +558,53 @@ PostOperation integaz UsingPost PostProcessings {
For i In {1:NumWires}
If( Flag_Thin == 0 )
- Print[ U, OnRegion Vol_C, File > "U_full.dat", Format Table,
- SendToServer Sprintf["}Voltage/Wire %g/1Full",i]];
+ Print[ U, OnRegion Vol_C, File > "U_full.dat", Format Table,
+ SendToServer Sprintf["}Voltage/Wire %g/1Full",i]];
Print[ L, OnRegion Vol_C, File > "L_full.dat", Format Table,
SendToServer Sprintf["}Inductance/Wire %g/1Full",i]];
- Print[ JouleLosses[ VWIRE~{i}], OnGlobal,
- Format Table, File > "joule_full.dat",
- SendToServer Sprintf["}Joule Losses/Wire %g/1Full",i]];
+ Print[ JouleLosses[ VWIRE~{i}], OnGlobal,
+ Format Table, File > "joule_full.dat",
+ SendToServer Sprintf["}Joule Losses/Wire %g/1Full",i]];
EndIf
- If( Flag_Thin == 1 )
- Print[ U, OnRegion lVol_C, File > "U_raw.dat", Format Table,
+
+ If( Flag_Thin == 1 )
+ Print[ U, OnRegion lVol_C, File > "U_raw.dat", Format Table,
SendToServer Sprintf["}Voltage/Wire %g/2Raw",i]];
Print[ L, OnRegion lVol_C, File > "L_raw.dat", Format Table,
SendToServer Sprintf["}Inductance/Wire %g/2Raw",i]];
- Print[ JouleLosses, OnRegion LWIRE~{i},
+ Print[ JouleLosses, OnRegion LWIRE~{i},
Format Table, File > "joule_raw.dat",
SendToServer Sprintf["}Joule Losses/Wire %g/2Raw",i]];
EndIf
- If( Flag_Thin == 2 )
- Print[ U, OnRegion lVol_C, File > "U_reg.dat", Format Table,
- SendToServer Sprintf["}Voltage/Wire %g/3Reg",i]];
+
+ If( Flag_Thin == 2 )
+ Print[ U, OnRegion lVol_C, File > "U_reg.dat", Format Table,
+ SendToServer Sprintf["}Voltage/Wire %g/3Reg",i]];
Print[ L, OnRegion lVol_C, File > "L_reg.dat", Format Table,
SendToServer Sprintf["}Inductance/Wire %g/3Reg",i]];
- Print[ JouleLosses, OnRegion LWIRE~{i},
- Format Table, File > "joule_reg.dat",
- SendToServer Sprintf["}Joule Losses/Wire %g/3Reg",i]];
+ Print[ JouleLosses, OnRegion LWIRE~{i},
+ Format Table, File > "joule_reg.dat",
+ SendToServer Sprintf["}Joule Losses/Wire %g/3Reg",i]];
EndIf
- If( Flag_Thin == 3 )
- Print[ U, OnRegion lVol_C, File > "U_naive.dat", Format Table,
- SendToServer Sprintf["}Voltage/Wire %g/4Naive",i]];
+
+ If( Flag_Thin == 3 )
+ Print[ U, OnRegion lVol_C, File > "U_naive.dat", Format Table,
+ SendToServer Sprintf["}Voltage/Wire %g/4Naive",i]];
Print[ L, OnRegion lVol_C, File > "L_naive.dat", Format Table,
SendToServer Sprintf["}Inductance/Wire %g/4naive",i]];
- Print[ JouleLosses, OnRegion LWIRE~{i},
- Format Table, File > "joule_naive.dat",
- SendToServer Sprintf["}Joule Losses/Wire %g/4Naiv",i]];
+ Print[ JouleLosses, OnRegion LWIRE~{i},
+ Format Table, File > "joule_naive.dat",
+ SendToServer Sprintf["}Joule Losses/Wire %g/4Naiv",i]];
EndIf
- If( Flag_Thin != 0 )
- Print[ area~{i}[SWIRE~{i}], OnGlobal,
- Format Table, File >"zzz",
- StoreInVariable $sarea~{i}];
+ If( Flag_Thin != 0 )
+ Print[ surf[SWIRE~{i}], OnGlobal,
+ Format Table, File >"zzz", StoreInVariable $sarea~{i}];
Print[ intby~{i}[SWIRE~{i}], OnGlobal,
- Format Table, File > "zzz",
- StoreInVariable $intby~{i}];
- EndIf
+ Format Table, File > "zzz", StoreInVariable $intby~{i}];
+ EndIf
EndFor
}
-
-
- /*{ Name flux_corr; // corrected flux
- Value {
- Integral { [ CompZ[ {ac} - {aw} ] ];
- In lVol_C; Integration I1 ; Jacobian Vol; }
- Integral { [ Analytic_L[]*{Is} ];
- In lVol_C; Integration I1 ; Jacobian Vol; }
- }
- }*/
-
- /* Print [awz, OnElementsOf Vol_nu, File "awz.pos"];
- Echo[ StrCat["l=PostProcessing.NbViews-1;",
- "View[l].IntervalsType = 3;",
- "View[l].NbIso = 30;",
- "View[l].RaiseZ = 15000;",
- "View[l].NormalRaise = 0;"],
- File "tmp.geo", LastTimeStepOnly] ;
-
- Print [acz, OnElementsOf Vol_nu, File "acz.pos"];
- Echo[ StrCat["l=PostProcessing.NbViews-1;",
- "View[l].IntervalsType = 3;",
- "View[l].NbIso = 30;",
- "View[l].RaiseZ = 15000;",
- "View[l].NormalRaise = 0;"],
- File "tmp.geo", LastTimeStepOnly] ;
-
-
-
- Print[ flux[ VWIRE~{i} ], OnGlobal,
- Format Table, File > "flux_full.dat",
- SendToServer Sprintf["}Flux/Wire %g/1Full",i]];
- Print[ flux[ LWIRE~{i} ], OnGlobal,
- Format Table, File > "flux_raw.dat",
- SendToServer Sprintf["}Flux/Wire %g/2Raw",i]];
- Print[ flux[ LWIRE~{i} ], OnGlobal,
- Format Table, File > "flux_reg.dat",
- SendToServer Sprintf["}Flux/Wire %g/4Reg",i]];
-*/
- /*
- X = rw/skin_depth;
- pI[] = (X/2)*Re[(1+i[])*(J_0[(1+i[])*X]/J_1[(1+i[])*X])];
- Resist[] = Re[pI[]*R_dc];
- // radial analytical solution with a zero flux condition imposed at r=rs
- Analytic_a~{i}[] = mu0/(2*Pi) * WI~{i} *
- ( (R~{i}[]>WR~{i}) ? Log[rs/R~{i}[]] :
- Log[rs/WR~{i}] + (J0[tau[]*R~{i}[]/rw] - J0[tau[]]) / J1[tau[]] / tau[] );
-
- Analytic_L~{i}[] = mu0 / (2*Pi) *
- ( 2/tau[]/tau[]- J0[tau[]] / J1[tau[]] / tau[] - Log[rs/WR~{i}] ) ;
-
- // Idem in static case (thus particular case of the above)
- AnalyticStatic_a~{i}[] = mu0 * WI~{i} / (2*Pi) *
- ( (R~{i}[]>WR~{i}) ? Log[rs/R~{i}[]] : Log[rs/WR~{i}] + (1-(R~{i}[]/WR~{i})^2)/2. );
- Analytic_e~{i}[] = mu0 * WI~{i} / (2*Pi) * i[] * omega
- * J0[tau[]*R~{i}[]/WR~{i}] / J1[tau[]] / tau[] ;
- Analytic_flux~{i}[] = mu0 * WI~{i} / (2*Pi) *
- ( i[] * (skin_depth/WR~{i})^2 + Log[rs/WR~{i}] - J0[tau[]] / J1[tau[]] / tau[] ) ;
- */
diff --git a/Conductors3D/w3d.geo b/Conductors3D/w3d.geo
index cdc2a4c..521752c 100644
--- a/Conductors3D/w3d.geo
+++ b/Conductors3D/w3d.geo
@@ -116,7 +116,10 @@ ElseIf( Flag_RegSleeve ) // regular sleeve
s1 = news; Plane Surface(s1) = {ll1,-llWires[]};
e[] = Extrude {0,0,Lz} { Surface{ s1 } ; /*Layers{1}; Recombine;*/ } ;
Physical Volume ("AIR", 1) = { e[1], Wires[] };
- Physical Volume ("BLA", 5) = { e[1] };
+ // auxiliary Physical Volume to test the tree.
+ // Note that two Physical Volumes defined on a same region
+ // yield duplicate finite element in GetDP.
+ //Physical Volume ("BLA", 5) = { e[1] }; this
Else
@@ -148,4 +151,5 @@ If( !Flag_Thin ) // Wires with finite radius
Electrodes[] += {20+i, 30+i};
EndFor
EndIf
-Physical Line("LINTREE", 4) = { Boundary { Physical Surface { Electrodes[] }; } };
+
+Physical Line("LINTREE", 4) = { CombinedBoundary { Physical Surface { Electrodes[] }; } };
diff --git a/Conductors3D/w3d.pro b/Conductors3D/w3d.pro
index 69d4578..95e4db7 100644
--- a/Conductors3D/w3d.pro
+++ b/Conductors3D/w3d.pro
@@ -33,8 +33,7 @@ Group{
AIR = Region[ 1 ];
INF = Region[ 2 ];
SKIN = Region[ 3 ];
- LINTREE = Region[ 4 ]; // not used
- BLA = Region[ 5 ];
+ LINTREE = Region[ 4 ];
VWIRES = Region[ {} ];
LWIRES = Region[ {} ];
@@ -62,7 +61,6 @@ Group{
CATHODES += Region[ CATHODE~{i} ];
EndFor
-
// Abstract regions
If( !Flag_Thin )
@@ -76,7 +74,7 @@ Group{
Vol_nu = Region[ {Vol_C, Vol_CC} ];
Sur_Dirichlet_a = Region[ INF ];
Electrodes = Region[ {ANODES, CATHODES} ];
-
+
Dom_Hcurl_a = Region[ Vol_nu ];
Dom_Hgrad_u = Region[ Vol_C ];
Dom_Hregion_i = Region[ Vol_C ];
@@ -84,7 +82,6 @@ Group{
// additional Groups for the semi_analytic approach
Dom_Hthin_a = ElementsOf[ Vol_nu, OnOneSideOf LWIRES ];
Vol_Tree = ElementsOf[ Vol_nu, DisjointOf LWIRES ];
- //Vol_Tree = Region[ { BLA } ];
Sur_Tree = Region[ { Sur_Dirichlet_a /*, SKIN*/ } ];
If( !Flag_SemiAnalytic )
Lin_Tree = Region[ {} ];
@@ -92,7 +89,6 @@ Group{
Lin_Tree = Region[ { LWIRES } ];
EndIf
}
-
Function{
omega = 2*Pi*Freq;
@@ -267,7 +263,7 @@ FunctionSpace {
{ NameOfCoef ac;
EntityType Auto; NameOfConstraint Hcurl_a_3D; }
{ NameOfCoef ac; EntityType EdgesOfTreeIn ; EntitySubType StartingOn ;
- NameOfConstraint GaugeCondition_a ; }
+ NameOfConstraint GaugeCondition_a ; }
}
}
@@ -576,16 +572,13 @@ PostOperation map UsingPost MagnetoDynamics {
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"];
- PrintGroup[ _CO_Entity_71, In Vol_nu, File "Tree.pos"];
+ // PrintGroup[ EdgesOf[ Vol_nu, ConnectedTo LWIRES ],
+ // In Vol_nu, File "group.pos"];
+ //PrintGroup[ _CO_Entity_71, In Vol_nu, File "Tree.pos"];
EndIf
-
- // PrintGroup[ EdgesOfTreeIn[ { Vol_Tree }, StartingOn { Sur_Tree, Lin_Tree } ],
- // In Vol_nu, File "Tree.pos"];
-
-
-
+ //PrintGroup[ _CO_Entity_30, In Vol_nu, File "Tree.pos"];
+ PrintGroup[ EdgesOfTreeIn[ Vol_Tree , StartingOn { Sur_Tree, Lin_Tree } ],
+ In Vol_nu, File "Tree.pos"];
}
PostOperation integaz UsingPost MagnetoDynamics {
diff --git a/Conductors3D/w3d_common.pro b/Conductors3D/w3d_common.pro
index 8bb6be1..3f3af2a 100644
--- a/Conductors3D/w3d_common.pro
+++ b/Conductors3D/w3d_common.pro
@@ -22,8 +22,8 @@ DefineConstant[
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_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]"}
--
GitLab