Skip to content
Snippets Groups Projects
Commit 87ce5349 authored by François Henrotte's avatar François Henrotte
Browse files

Merge branch 'master' of http://gitlab.onelab.info/doc/tutorials

parents d7ce1b2d b5c4660d
No related branches found
No related tags found
No related merge requests found
Pipeline #5862 passed
......@@ -15,7 +15,8 @@ xBox = w/2. * 6. ; yBox = h * 12. ;
/* Definition of parameters for local mesh dimensions */
s = 1. ;
DefineConstant[ s = {1., Name "Parameters/Global mesh size factor"} ] ;
p0 = h / 10. * s ;
pLine0 = w/2. / 10. * s ; pLine1 = w/2. / 50. * s ;
pxBox = xBox / 10. * s ; pyBox = yBox / 8. * s ;
......
......@@ -4,14 +4,18 @@
Features:
- Use of a template formulation library
- Identical to Tutorial 2 for a static current source
- Frequency-domain solution (phasor) for a dynamic current source
- Frequency-domain or time-domain solution for a time-dependent
current source
To compute the static solution in a terminal:
getdp electromagnet -solve Magnetostatics2D_a -pos Map_a
To compute the time-harmonic dynamic solution in a terminal:
To compute the frequency-domain solution in a terminal:
getdp electromagnet -solve Magnetodynamics2D_av -pos Map_a
To compute the time-dependent dynamic solution in a terminal:
getdp electromagnet -setnumber TimeDomain 1 -solve Magnetodynamics2D_av -pos Map_a
To compute the solution interactively from the Gmsh GUI:
File > Open > electromagnet.pro
You may choose the Resolution in the left panel:
......@@ -42,10 +46,27 @@ Function {
DefineConstant[
murCore = {100, Name "Model parameters/Mur core"},
Current = {0.01, Name "Model parameters/Current"},
frequency = {1, Name "Model parameters/Frequency"}
TimeDomain = {0, Choices{0 = "Frequency-domain", 1 = "Time-domain"},
Name "Model parameters/03Analysis type"}
frequency = {50, Visible !TimeDomain,
Name "Model parameters/Frequency" }
];
If(TimeDomain)
// Fix parameters from the "Lib_Magnetodynamics2D_av_Cir.pro" template:
Flag_FrequencyDomain = 0;
TimeInit = 0; // start simulation at time = 0s
TimeFinal = 20e-3; // stop simulation at time = 20 ms
DeltaTime = 1e-3; // use time steps equal to 1 ms
// Define the time modulation of the current source, i.e. a linear ramp from
// 0 to 1 until 10 ms, then a constant value of 1:
myModulation[] = ($Time < TimeFinal / 2) ? (2 / TimeFinal * $Time) : 1;
Else
// Fix parameters from the "Lib_Magnetodynamics2D_av_Cir.pro" template:
Flag_FrequencyDomain = 1;
Freq = frequency;
EndIf
mu0 = 4.e-7 * Pi;
nu[ Region[{Air, Ind, AirInf}] ] = 1. / mu0;
nu[ Core ] = 1. / (murCore * mu0);
......@@ -55,10 +76,13 @@ Function {
Ns[ Ind ] = 1000 ; // number of turns in coil
Sc[ Ind ] = SurfaceArea[] ; // area of coil cross section
// Current density in each coil portion for a unit current (will be multiplied
// by the actual total current in the coil)
js0[ Ind ] = Ns[] / Sc[] * Vector[0, 0, -1];
CoefGeos[] = 1;
// For a correct definition of the voltage:
CoefGeos[] = 1; // planar model, 1 meter thick
}
Constraint {
......@@ -70,8 +94,13 @@ Constraint {
}
{ Name Current_2D;
Case {
// represents the phasor amplitude (peak to peak value) for a dynamic analysis
If(Flag_FrequencyDomain)
// Amplitude of the phasor is set to "Current"
{ Region Ind; Value Current; }
Else
// Time-dependent value is set to "Current * myModulation[]"
{ Region Ind; Value Current; TimeFunction myModulation[]; }
EndIf
}
}
{ Name Voltage_2D;
......@@ -87,8 +116,8 @@ PostOperation {
{ Name Map_a; NameOfPostProcessing Magnetodynamics2D_av;
Operation {
Print[ a, OnElementsOf Vol_Mag, File "a.pos" ];
Print[ b, OnElementsOf Vol_Mag, File "b.pos" , HarmonicToTime 20];
Print[ j, OnElementsOf Vol_Mag, File "j.pos", HarmonicToTime 20];
Print[ b, OnElementsOf Vol_Mag, File "b.pos" ];
Print[ j, OnElementsOf Vol_Mag, File "j.pos" ];
}
}
}
......@@ -3,7 +3,6 @@
Features:
- Use of a generic template formulation library
- Frequency- and time-domain dynamic solutions
- Circuit coupling used as a black-box (see Tutorial 8 for details)
To compute the solution in a terminal:
......@@ -17,13 +16,12 @@
Include "transfo_common.pro";
DefineConstant[
type_Conds = {2, Choices{1 = "Massive", 2 = "Coil"}, Highlight "Blue",
// The "Massive" option is non-physical, but it's interesting to highlight the
// fact that both massive and stranded conductors can be handled in the same
// way as far as circuit-coupling is concerned
ConductorType = {2, Choices{1 = "Massive", 2 = "Coil"}, Highlight "Blue",
Name "Parameters/01Conductor type"}
type_Source = {2, Choices{1 = "Current", 2 = "Voltage"}, Highlight "Blue",
Name "Parameters/02Source type"}
type_Analysis = {1, Choices{1 = "Frequency-domain", 2 = "Time-domain"}, Highlight "Blue",
Name "Parameters/03Analysis type"}
Freq = {50, Min 0, Max 1e3, Step 1,
Freq = {1, Min 0, Max 1e3, Step 1,
Name "Parameters/Frequency"}
];
......@@ -43,81 +41,71 @@ Group {
// Abstract regions that will be used in the "Lib_Magnetodynamics2D_av_Cir.pro"
// template file included below;
Vol_Mag = Region[{Air, Core, Coils}]; // full magnetic domain
If (type_Conds == 1)
If (ConductorType == 1)
Vol_C_Mag = Region[{Coils}]; // massive conductors
ElseIf (type_Conds == 2)
ElseIf (ConductorType == 2)
Vol_S_Mag = Region[{Coils}]; // stranded conductors (coils)
EndIf
}
Function {
mu0 = 4e-7*Pi;
mur_Core = 1000;
mu[Air] = 1 * mu0;
mur_Core = 100;
mu[Core] = mur_Core * mu0;
mu[Coils] = 1 * mu0;
sigma[Coils] = 1e7;
// For a correct definition of the voltage
CoefGeo = thickness_Core;
nu[] = 1 / mu[];
// To be defined separately for each coil portion
Sc[Coil_1_P] = SurfaceArea[];
SignBranch[Coil_1_P] = 1; // To fix the convention of positive current (1:
// along Oz, -1: along -Oz)
sigma[Coils] = 1e7;
Sc[Coil_1_M] = SurfaceArea[];
// To be defined separately for each coil portion, to fix the convention of
// positive current (1: along Oz, -1: along -Oz)
SignBranch[Coil_1_P] = 1;
SignBranch[Coil_1_M] = -1;
Sc[Coil_2_P] = SurfaceArea[];
SignBranch[Coil_2_P] = 1;
Sc[Coil_2_M] = SurfaceArea[];
SignBranch[Coil_2_M] = -1;
// Number of turns (same for PLUS and MINUS portions) (half values because
// half coils are defined)
Ns[Coil_1] = 1;
Ns[Coil_2] = 1;
If(ConductorType == 2)
// Number of turns in the coils (same for PLUS and MINUS portions) - half
// values because half coils are defined geometrically due to symmetry!
Ns[Coil_1] = 100;
Ns[Coil_2] = 200;
// Global definitions (nothing to change):
// To be defined separately for each coil portion:
Sc[Coil_1_P] = SurfaceArea[];
Sc[Coil_1_M] = SurfaceArea[];
Sc[Coil_2_P] = SurfaceArea[];
Sc[Coil_2_M] = SurfaceArea[];
// Current density in each coil portion for a unit current (will be multiplied
// by the actual total current in the coil)
// Current density in each coil portion for a unit current (will be
// multiplied by the actual total current in the coil), in the case of
// stranded conductors
js0[Coils] = Ns[] / Sc[] * Vector[0, 0, SignBranch[]];
CoefGeos[Coils] = SignBranch[] * CoefGeo;
// The reluctivity will be used
nu[] = 1/mu[];
}
If(type_Analysis == 1)
Flag_FrequencyDomain = 1;
Else
Flag_FrequencyDomain = 0;
EndIf
If (type_Source == 1) // current
Flag_CircuitCoupling = 0;
ElseIf (type_Source == 2) // voltage
// For a correct definition of the voltage
CoefGeos[Coils] = SignBranch[] * thickness_Core;
}
// PLUS and MINUS coil portions to be connected in series, with applied
// voltage on the resulting branch
// We will use a circuit coupling to connect the PLUS and MINUS portions of the
// coil in series, for both the primary and secondary. We will then apply a
// voltage source (with a small resistance in series) on the resulting branch on
// the primary, and connect a resistive load on the resulting branch on the
// secondary.
Flag_CircuitCoupling = 1;
// Note that the voltage will not be equally distributed in the PLUS and MINUS
// parts, which is the reason why we must apply the total voltage through a
// circuit -- and we cannot simply use a current source like in Tutorial 7a.
// Here is the definition of the circuits on primary and secondary sides:
Group {
// Empty Groups to be filled
Resistance_Cir = Region[{}];
Inductance_Cir = Region[{}] ;
Capacitance_Cir = Region[{}] ;
SourceV_Cir = Region[{}]; // Voltage sources
SourceI_Cir = Region[{}]; // Current sources
Resistance_Cir = Region[{}]; // all resistances
Inductance_Cir = Region[{}] ; // all inductances
Capacitance_Cir = Region[{}] ; // all capacitances
SourceV_Cir = Region[{}]; // all voltage sources
SourceI_Cir = Region[{}]; // all current sources
// Primary side
E_in = Region[10001]; // arbitrary region number (not linked to the mesh)
......@@ -132,8 +120,8 @@ ElseIf (type_Source == 2) // voltage
Function {
deg = Pi/180;
// Input RMS voltage (half of the voltage because of symmetry; half coils
// are defined)
// Input RMS voltage (half of the voltage because of symmetry; half coils are
// defined!)
val_E_in = 1.;
phase_E_in = 90 *deg; // Phase in radian (from phase in degree)
......@@ -143,27 +131,41 @@ ElseIf (type_Source == 2) // voltage
// End-winding primary winding resistance for more realistic primary coil
// model
Resistance[R_in] = 1;
Resistance[R_in] = 1e-3;
}
Constraint {
{ Name MagneticVectorPotential_2D;
Case {
{ Region Sur_Air_Ext; Value 0; }
}
}
{ Name Current_2D;
Case {
}
}
{ Name Voltage_2D;
Case {
}
}
{ Name Current_Cir ;
Case {
}
}
{ Name Voltage_Cir ;
Case {
{ Region E_in; Value val_E_in; TimeFunction F_Cos_wt_p[]{2*Pi*Freq, phase_E_in}; }
{ Region E_in; Value val_E_in;
// F_Cos_wt_p[] is a built-in function with two parameters (w and p),
// that can be used to evaluate cos(w * t + p) in both frequency- and
// time-domain
TimeFunction F_Cos_wt_p[]{2*Pi*Freq, phase_E_in}; }
}
}
{ Name ElectricalCircuit ; Type Network ;
Case Circuit_1 {
// PLUS and MINUS coil portions to be connected in series, together with
// E_in (an additional resistor should be defined to represent the
// Coil_1 end-winding (not considered in the 2D model))
// E_in; an additional resistor is defined in series to represent the
// Coil_1 end-winding, which is not considered in the 2D model.
{ Region E_in; Branch {1,4}; }
{ Region R_in; Branch {4,2}; }
......@@ -172,38 +174,14 @@ ElseIf (type_Source == 2) // voltage
}
Case Circuit_2 {
// PLUS and MINUS coil portions to be connected in series, together with
// R_out (an additional resistor should be defined to represent the
// Coil_2 end-winding (not considered in the 2D model))
// R_out (an additional resistor could be defined to represent the Coil_2
// end-winding - but we can directly add it to R_out as well)
{ Region R_out; Branch {1,2}; }
{ Region Coil_2_P; Branch {2,3} ; }
{ Region Coil_2_M; Branch {3,1} ; }
}
}
}
EndIf
Constraint {
{ Name MagneticVectorPotential_2D;
Case {
{ Region Sur_Air_Ext; Value 0; }
}
}
{ Name Current_2D;
Case {
If (type_Source == 1)
// Current in each coil (same for PLUS and MINUS portions)
{ Region Coil_1; Value 1; TimeFunction F_Sin_wt_p[]{2*Pi*Freq, 0}; }
{ Region Coil_2; Value 0; }
EndIf
}
}
{ Name Voltage_2D;
Case {
}
}
}
Include "Lib_Magnetodynamics2D_av_Cir.pro";
......@@ -214,40 +192,19 @@ PostOperation {
Print[ j, OnElementsOf Region[{Vol_C_Mag, Vol_S_Mag}], Format Gmsh, File "j.pos" ];
Print[ b, OnElementsOf Vol_Mag, Format Gmsh, File "b.pos" ];
Print[ az, OnElementsOf Vol_Mag, Format Gmsh, File "az.pos" ];
If (type_Analysis == 1) // frequency domain
If (type_Source == 1) // current
// In text file UI.txt: voltage and current for each coil portion (note
// that the voltage is not equally distributed in PLUS and MINUS
// portions, which is the reason why we must apply the total voltage
// through a circuit -> type_Source == 2)
Echo[ "Coil_1_P", Format Table, File "UI.txt" ];
Print[ U, OnRegion Coil_1_P, Format FrequencyTable, File > "UI.txt" ];
Print[ I, OnRegion Coil_1_P, Format FrequencyTable, File > "UI.txt"];
Echo[ "Coil_1_M", Format Table, File > "UI.txt" ];
Print[ U, OnRegion Coil_1_M, Format FrequencyTable, File > "UI.txt" ];
Print[ I, OnRegion Coil_1_M, Format FrequencyTable, File > "UI.txt"];
Echo[ "Coil_2_P", Format Table, File > "UI.txt" ];
Print[ U, OnRegion Coil_2_P, Format FrequencyTable, File > "UI.txt" ];
Print[ I, OnRegion Coil_2_P, Format FrequencyTable, File > "UI.txt"];
Echo[ "Coil_2_M", Format Table, File > "UI.txt" ];
Print[ U, OnRegion Coil_2_M, Format FrequencyTable, File > "UI.txt" ];
Print[ I, OnRegion Coil_2_M, Format FrequencyTable, File > "UI.txt"];
ElseIf (type_Source == 2)
// In text file UI.txt: voltage and current of the primary coil (from E_in)
// (real and imaginary parts!)
If (Flag_FrequencyDomain)
// In text file UI.txt: voltage and current of the primary coil (from
// E_in) (real and imaginary parts!)
Echo[ "E_in", Format Table, File "UI.txt" ];
Print[ U, OnRegion E_in, Format FrequencyTable, File > "UI.txt" ];
Print[ I, OnRegion E_in, Format FrequencyTable, File > "UI.txt"];
// In text file UI.txt: voltage and current of the secondary coil (from R_out)
// In text file UI.txt: voltage and current of the secondary coil (from
// R_out)
Echo[ "R_out", Format Table, File > "UI.txt" ];
Print[ U, OnRegion R_out, Format FrequencyTable, File > "UI.txt" ];
Print[ I, OnRegion R_out, Format FrequencyTable, File > "UI.txt"];
EndIf
EndIf
}
}
}
......@@ -218,10 +218,10 @@ Integration {
/* The function space "Hregion_j_Mag_2D" provides one basis function, and hence
one degree of freedom, per physical region in the abstract region
"Vol_S_Mag". The constraint "SourceCurrentDensityZ" fixes all these dofs, so
the FunctionSpace "Hregion_j_Mag_2D" is fully fixed and has no FE unknowns.
One could thus have replaced it by a simple function and the Integral term
below in the weak formulation could have been written as
"Vol_S_Mag". The constraint "SourceCurrentDensityZ" fixes all these dofs (a
single one!), so the FunctionSpace "Hregion_j_Mag_2D" is fully fixed and has
no FE unknowns. One could thus have replaced it by a simple function and the
Integral term below in the weak formulation could have been written as
Integral { [ Vector[ 0,0,-js_fct[] ] , {a} ];
In Vol_S_Mag; Jacobian Vol; Integration Int; }
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment