quelques questions sur l'Axi
trophime christophe
trophime at labs.polycnrs-gre.fr
Wed Mar 21 13:58:03 CET 2001
Christophe Geuzaine wrote:
>
> trophime christophe wrote:
> >
> > Il y a des points que je comprends pas bien. Est-ce que pourrai eclaire
> > ma
> > lanterne?
> >
> > 1) Pour la formulation MagDyn_Axi j'ai declare la formulation de la
> > facon
> > suivante :
> >
> > Formulation {
> > { Name Magnetodynamics_av_Axi ; Type FemEquation ;
> > Quantity {
> > { Name a ; Type Local ; NameOfSpace Hcurl_a_Mag_Axi ; }
> > { Name ur ; Type Local ; NameOfSpace Hregion_u_Mag_Axi ; }
> > { Name I ; Type Global ; NameOfSpace Hregion_u_Mag_Axi [I] ; }
> > { Name U ; Type Global ; NameOfSpace Hregion_u_Mag_Axi [U] ; }
> > { Name js ; Type Local ; NameOfSpace Hregion_j_Mag_Axi ; }
> > }
> > Equation {
> > Galerkin { [ nu[] * Dof{d a} , {d a} ] ; In Domain_Mag ;
> > Jacobian VolAxi ; Integration CurlCurl ; }
> >
> > Galerkin { DtDof [ sigma[] * Dof{a} , {a} ] ; In DomainC_Mag ;
> > Jacobian VolAxi ; Integration CurlCurl ; }
> > Galerkin { [ sigma[] * Factor[] * Dof{ur} , {a} ] ; In
> > DomainC_Mag ;
> > Jacobian VolAxi ; Integration CurlCurl ; }
> >
> > Galerkin { [ - sigma[] * (Velocity[] *^ Dof{d a}) , {a} ] ;
> > In DomainV_Mag ;
> > Jacobian VolAxi ; Integration CurlCurl ; }
> >
> > Galerkin { [ - Dof{js} , {a} ] ; In DomainS_Mag ;
> > Jacobian VolAxi ;
> > Integration CurlCurl ; }
> >
> > Galerkin { DtDof [ sigma[] * Dof{a} , {ur} ] ; In DomainC_Mag ;
> > Jacobian VolAxi ; Integration CurlCurl ; }
> > Galerkin { [ sigma[] * Factor[] * Dof{ur} , {ur} ] ; In
> > DomainC_Mag ;
> > Jacobian VolAxi ; Integration CurlCurl ; }
> > GlobalTerm { [ Dof{I} , {U} ] ; In DomainC_Mag ; }
> > }
> > }
> > }
> >
> > sans rien changer aux definitions des espaces fonctionnels fournis dans
> > l'exemple
> > du manuel.
> >
> > La fonction Factor[] (qui vaut 1/(2*Pi) dans le cas Axi) est juste la
> > pour que
> > je garde la meme formulation que le 2D.
> >
> > Dans le post-processing je definis les grandeurs suivantes :
> >
> > { Name V ; Value { Local { [ CompZ[{ur}] ] ; In
> > Domain_Mag ; Jacobian VolAxi ;} } }
> > { Name jc ;
> > Value {
> > Local { [ - sigma[]*CompZ[Factor[] *{ur}] ] ; In DomainC_Mag
> > ; Jacobian VolAxi ;}
> > }
> > }
> >
> > En toute logique je me dis que si je trace V (que je multiplie ensuite
> > par sigma[]/(2*Pi)) et jc
> > par exemple dans Gnuplot je devrais trouver la meme chose... Et bien
> > non! La je comprends pas ce
> > qui ce passe.
>
> C'est la fonction Factor[] qui pose probleme (aussi bien dans la
> formulation que dans le post-pro), et qui doit etre completement
> supprimee. A partir du moment ou tu utilises un jacobien de
> transformation axisymetrique, tu ne dois pas rajouter ce facteur
> manuellement. La formulation est valable aussi bien dans le cas 2D plan
> que 2D axi. La chose qu'il faut bien saisir, c'est que les inconnues en
> potentiel vecteur sont des _circulations_ de ce potentiel le long d'un
> arc de cercle unitaire. Tu devrais introduire des corrections si tu
> pre-supposais connue une grandeur qui n'est pas tansformee, alors
> qu'elle devrait (e.g. une densite de courant imposee par l'intermediaire
> d'une fonction, sans interpolation).
>
Je viens de supprimer Factor[] et j'obtiens les resulats de MagSta_a_Axi
a un facteur \frac{1}{2 \pi} pres. Si je ne dois pas introduire facteur
Factor[] dans la formulation c'est donc qu'il y a un probleme ailleurs.
Je te mets en attachment mon fichier geo et mes fichiers pro pour
que tu vois si il ya un probleme avec les donnees ou avec getdp...
Par contre effectivement si je trace jz et -sigma*ur ca colle
maintenant.
> >
> > Par ailleurs j'ai fait un modele "complet" de deux inducteurs
> > concentriques
> > et un modele avec symetrie suivant l'axe Or pour avoir simplement une
> > moitie du
> > probleme precedent. J'obtiens bien les memes resultats en terme de
> > potentiel a,v et
> > de densite de courant mais par contre pour le courant total I j'obtiens
> > des resultats
> > totalement incoherent. Mais cette fois je pense qu'il y a un "probleme"
> > dans la definition
> > des contraintes... J'ai simplement ajoute une "ligne" dans la definition
> > des Group sans
> > preciser de contrainte pour a et ur. Au niveau de a c'est bien un
> > neumann homogene qui vient
> > mais pour ur qu'est-ce qu'il faudrait mettre?
>
> Tu dois imposer soit le courant global, soit la tension globale.
>
C'est pourtant ce que j'ai fait.
> >
> > 3) Enfin je me suis dit que le resultat du probleme magnetodynamique
> > devrait...
> > >
> Si je te suis, tu imposes un "step" de tension ? je ne sais pas quels
> parametres tu as utilises pour la resolution temporelle. Il faut garder
> a l'esprit que l'integration temporelle est dissipative (numeriquement).
> Si tu as utilise Euler implicite, le schema est seulement d'ordre 1. Tu
> pourrais essayer Cranck-Nicholson (i.e. theta=0.5, et qui est d'ordre 2,
> mais pas inconditionnellement stable). Il faudrait surtout qu'on
> implemente des schemas temporels d'ordre plus eleve dans GetDP.
>
Tout a fait d'accord avec toi. Mais j'esperai voir un changement en
passant
d'un schema Euler-Implicite a un Cranck-Nicholson... Ca veut dire que
pour ameliorer
la convergence il faudrait regarder des methodes de type time-splitting
par exemple...
Au fait c'est quoi le schema Newmark.
> >
> > 4) Par ailleurs, il y a toujours le meme probleme avec la transformation
> > pour AirInf
> > si je remets des valeurs de Rint et Rext identiques dans le fichier de
> > donnees et le
> > fichier equation.
>
> Les differences sont-elles petites ? Comme le maillage ne respecte pas
> exactement les frontieres courbes, ca pourrait expliquer que certains
> points soient hors de l'intervalle. Sinon, envoie-moi ton fichier de
> donnees, que je jette un oeil.
>
Non ca vient de moi. Je viens de verifier. Desole de t'avoir fait douter
de tes modifs.
-------------- next part --------------
/* -------------------------------------------------------------------
File "MagDyn_av_Axi.pro"
Magnetodynamics - Magnetic vector potential and electric scalar
potential a-v formulation (Axi)
-------------------------------------------------------------------
I N P U T
---------
GlobalGroup : (Extension '_Mag' is for Magnetic problem)
-----------
Domain_Mag Whole magnetic domain
DomainCC_Mag Nonconducting regions (not used)
DomainC_Mag Conducting regions
DomainS_Mag Inductor regions (Source)
DomainV_Mag All regions in movement (for speed term)
Function :
--------
nu[] Magnetic reluctivity
sigma[] Electric conductivity
Velocity[] Velocity of regions
Constraint :
----------
MagneticVectorPotential_Axi
Fixed magnetic vector potential (Axi)
(classical boundary condition)
SourceCurrentDensityZ Fixed source current density (in Z direction)
VolAxitage_Axi Fixed voltage
Current_Axi Fixed Current
Parameters :
----------
Freq Frequency (Hz)
Parameters for time loop with theta scheme :
Mag_Time0, Mag_TimeMax, Mag_DTime
Initial time, Maximum time, Time step (s)
Mag_Theta Theta (e.g. 1. : Implicit Euler,
0.5 : Cranck Nicholson)
*/
Group {
DefineGroup[ Domain_Mag, DomainCC_Mag, DomainC_Mag,
DomainS_Mag, DomainV_Mag ] ;
}
Function {
DefineFunction[ nu, sigma ] ;
DefineFunction[ Velocity ] ;
DefineVariable[ Freq ] ;
DefineVariable[ Mag_Time0, Mag_TimeMax, Mag_DTime, Mag_Theta ] ;
DefineFunction[ Factor ];
}
FunctionSpace {
// Magnetic vector potential a (b = curl a)
{ Name Hcurl_a_Mag_Axi ; Type Form1P ;
BasisFunction {
// a = a s
// e e
{ Name se ; NameOfCoef ae ; Function BF_PerpendicularEdge ;
Support Domain_Mag ; Entity NodesOf[ All ] ; }
}
Constraint {
{ NameOfCoef ae ; EntityType NodesOf ;
NameOfConstraint MagneticVectorPotential_Axi ; }
}
}
// Gradient of Electric scalar potential (Axi)
{ Name Hregion_u_Mag_Axi ; Type Form1P ;
BasisFunction {
{ Name sr ; NameOfCoef ur ; Function BF_RegionZ ;
Support DomainC_Mag ; Entity DomainC_Mag ; }
}
GlobalQuantity {
{ Name U ; Type AliasOf ; NameOfCoef ur ; }
{ Name I ; Type AssociatedWith ; NameOfCoef ur ; }
}
Constraint {
{ NameOfCoef U ; EntityType Region ;
NameOfConstraint Voltage_Axi ; }
{ NameOfCoef I ; EntityType Region ;
NameOfConstraint Current_Axi ; }
}
}
// Source current density js (fully fixed space)
{ Name Hregion_j_Mag_Axi ; Type Vector ;
BasisFunction {
{ Name sr ; NameOfCoef jsr ; Function BF_RegionZ ;
Support DomainS_Mag ; Entity DomainS_Mag ; }
}
Constraint {
{ NameOfCoef jsr ; EntityType Region ;
NameOfConstraint SourceCurrentDensityZ ; }
}
}
}
Formulation {
{ Name Magnetodynamics_av_Axi ; Type FemEquation ;
Quantity {
{ Name a ; Type Local ; NameOfSpace Hcurl_a_Mag_Axi ; }
{ Name ur ; Type Local ; NameOfSpace Hregion_u_Mag_Axi ; }
{ Name I ; Type Global ; NameOfSpace Hregion_u_Mag_Axi [I] ; }
{ Name U ; Type Global ; NameOfSpace Hregion_u_Mag_Axi [U] ; }
{ Name js ; Type Local ; NameOfSpace Hregion_j_Mag_Axi ; }
}
Equation {
Galerkin { [ nu[] * Dof{d a} , {d a} ] ; In Domain_Mag ;
Jacobian VolAxi ; Integration CurlCurl ; }
Galerkin { DtDof [ sigma[] * Dof{a} , {a} ] ; In DomainC_Mag ;
Jacobian VolAxi ; Integration CurlCurl ; }
Galerkin { [ sigma[] * Dof{ur} , {a} ] ; In DomainC_Mag ;
Jacobian VolAxi ; Integration CurlCurl ; }
Galerkin { [ - sigma[] * (Velocity[] *^ Dof{d a}) , {a} ] ;
In DomainV_Mag ;
Jacobian VolAxi ; Integration CurlCurl ; }
Galerkin { [ - Dof{js} , {a} ] ; In DomainS_Mag ;
Jacobian VolAxi ;
Integration CurlCurl ; }
Galerkin { DtDof [ sigma[] * Dof{a} , {ur} ] ; In DomainC_Mag ;
Jacobian VolAxi ; Integration CurlCurl ; }
Galerkin { [ sigma[] * Dof{ur} , {ur} ] ; In DomainC_Mag ;
Jacobian VolAxi ; Integration CurlCurl ; }
GlobalTerm { [ Dof{I} , {U} ] ; In DomainC_Mag ; }
}
}
}
Resolution {
{ Name MagDyn_av_Axi ;
System {
{ Name Sys_Mag ; NameOfFormulation Magnetodynamics_av_Axi ;
Type ComplexValue ; Frequency Freq ; }
}
Operation {
Generate Sys_Mag ; Solve Sys_Mag ; SaveSolution Sys_Mag ;
}
}
{ Name MagDyn_t_av_Axi ;
System {
{ Name Sys_Mag ; NameOfFormulation Magnetodynamics_av_Axi ; }
}
Operation {
InitSolution Sys_Mag ; SaveSolution Sys_Mag ;
TimeLoopTheta{
Time0 Mag_Time0 ; TimeMax Mag_TimeMax ;
DTime Mag_DTime[] ; Theta Mag_Theta[] ;
Operation {
Generate Sys_Mag ; Solve Sys_Mag ; SaveSolution Sys_Mag ;
}
}
}
}
}
PostProcessing {
{ Name MagDyn_av_Axi ; NameOfFormulation Magnetodynamics_av_Axi ;
PostQuantity {
{ Name a ; Value { Local { [ {a} ] ; In Domain_Mag ; Jacobian VolAxi ;} } }
{ Name az ; Value { Local { [ CompZ[{a}] ] ; In Domain_Mag ; Jacobian VolAxi ;} } }
{ Name b ; Value { Local { [ {d a} ] ; In Domain_Mag ; Jacobian VolAxi ;} } }
{ Name h ; Value { Local { [ nu[] * {d a} ] ; In Domain_Mag ; Jacobian VolAxi ;} } }
{ Name j ;
Value {
Local { [ - sigma[]*(Dt[{a}]+{ur}) ] ; In DomainC_Mag ; Jacobian VolAxi ;}
}
}
{ Name jz ;
Value {
Local { [ - sigma[]*CompZ[Dt[{a}]+{ur}] ] ; In DomainC_Mag ; Jacobian VolAxi ;}
}
}
{ Name jt ;
Value {
Local { [ - sigma[]*(Dt[{a}]+{ur}) + {js} ] ; In DomainC_Mag ; Jacobian VolAxi ;}
}
}
{ Name jtz ;
Value {
Local { [ - sigma[]*CompZ[Dt[{a}]+{ur}] + CompZ[{js}] ] ; In DomainC_Mag ; Jacobian VolAxi ;}
}
}
{ Name roj2 ;
Value {
Local { [ sigma[]*SquNorm[Dt[{a}]+{ur}] ] ; In DomainC_Mag ; Jacobian VolAxi ;}
}
}
{ Name U ; Value { Local { [ {U} ] ; In DomainC_Mag ; Jacobian VolAxi ;} } }
{ Name I ; Value { Local { [ {I} ] ; In DomainC_Mag ; Jacobian VolAxi ;} } }
{ Name Z ; Value { Local { [ {U}/{I} ] ; In DomainC_Mag ; Jacobian VolAxi ;} } }
}
}
{ Name MagDyn_t_av_Axi ; NameOfFormulation Magnetodynamics_av_Axi ;
PostQuantity {
{ Name a ; Value { Local { [ {a} ] ; In Domain_Mag ; Jacobian VolAxi ;} } }
{ Name V ; Value { Local { [ sigma[]*CompZ[{ur}] ] ; In Domain_Mag ; Jacobian VolAxi ;} } }
{ Name az ; Value { Local { [ CompZ[{a}] ] ; In Domain_Mag ; Jacobian VolAxi ;} } }
{ Name b ; Value { Local { [ {d a} ] ; In Domain_Mag ; Jacobian VolAxi ;} } }
{ Name h ; Value { Local { [ nu[] * {d a} ] ; In Domain_Mag ; Jacobian VolAxi ;} } }
{ Name j ;
Value {
Local { [ - sigma[]*(Dt[{a}]+{ur}) ] ; In DomainC_Mag ; Jacobian VolAxi ;}
}
}
{ Name jz ;
Value {
Local { [ - sigma[]*CompZ[Dt[{a}]+{ur}] ] ; In DomainC_Mag ; Jacobian VolAxi ;}
}
}
{ Name jt ;
Value {
Local { [ - sigma[]*(Dt[{a}]+{ur}) + {js} ] ; In DomainC_Mag ; Jacobian VolAxi ;}
}
}
{ Name jtz ;
Value {
Local { [ - sigma[]*CompZ[Dt[{a}]+{ur}] + CompZ[{js}] ] ; In DomainC_Mag ; Jacobian VolAxi ;}
}
}
{ Name jc ;
Value {
Local { [ - sigma[]*CompZ[{ur}] ] ; In DomainC_Mag ; Jacobian VolAxi ;}
}
}
{ Name jf ;
Value {
Local { [ - sigma[]*CompZ[Dt[{a}]] ] ; In DomainC_Mag ; Jacobian VolAxi ;}
}
}
{ Name roj2 ;
Value {
Local { [ sigma[]*SquNorm[Dt[{a}]+{ur}] ] ; In DomainC_Mag ; Jacobian VolAxi ;}
}
}
{ Name U ; Value { Local { [ {U} ] ; In DomainC_Mag ; Jacobian VolAxi ;} } }
{ Name I ; Value { Local { [ {I} ] ; In DomainC_Mag ; Jacobian VolAxi ;} } }
{ Name R ; Value { Local { [ Fabs[{U}/{I}] ] ; In DomainC_Mag ; Jacobian VolAxi ;} } }
}
}
}
-------------- next part --------------
/* -------------------------------------------------------------------
File "CoreMassive.pro"
This file defines the problem dependent data structures for the
dynamic core-inductor problem.
To compute the solution: getdp CoreMassive -msh core.msh -solve MagDyn_av_Axi
To compute post-results: getdp CoreMassive -msh core.msh -pos Map_a
or getdp CoreMassive -msh core.msh -pos U_av
------------------------------------------------------------------- */
Group {
Ind = Region[ 240 ];
Ind1 = Region[ 241 ];
Air = Region[ 242 ];
AirInf = Region[ 243 ];
SurfaceAxe = Region[ 244 ];
SurfaceGInf = Region[ 245 ];
IntSurfaceInd = Region[ 246 ];
BottomSurfaceInd = Region[ 247 ];
ExtSurfaceInd = Region[ 248 ];
TopSurfaceInd = Region[ 249 ];
IntSurfaceInd1 = Region[ 250 ];
BottomSurfaceInd1 = Region[ 251 ];
ExtSurfaceInd1 = Region[ 252 ];
TopSurfaceInd1 = Region[ 253 ];
// Domains for Thermal Problems
SurfaceHeatFlux = Region[ {}];
SurfaceConvection = Region[ {IntSurfaceInd,ExtSurfaceInd,TopSurfaceInd,BottomSurfaceInd,
IntSurfaceInd1,ExtSurfaceInd1,TopSurfaceInd1,BottomSurfaceInd1}];
//SurfaceHeatFlux = Region[{}];
DomainC_Th = Region[ {}];
DomainS_Th = Region[ {Ind, Ind1}];
Domain_Th = Region[ {DomainS_Th, DomainC_Th}];
// Domains for Electromagnetics Problems
DomainCC_Mag = Region[ {Air, AirInf} ] ;
/* Massive inductor (Ind), in DomainC_Mag */
DomainB_Mag = Region[ {} ] ;
DomainS_Mag = Region[ {} ] ;
// if Static analysis
//DomainS_Mag = Region[ {DomainS_Th} ] ;
//DomainC_Mag = Region[ {} ] ;
// if transient analysis
DomainS_Mag = Region[ {} ] ;
DomainC_Mag = Region[ {DomainS_Th} ] ;
DomainInf = Region[ {AirInf} ] ;
Val_Rint = 399.999e-3 ;
Val_Rext = 500.001e-3 ;
Domain_Mag = Region[ {DomainCC_Mag, DomainC_Mag, DomainS_Mag} ] ;
}
Function {
Factor [ Region[{DomainC_Mag} ] ] = 1.0/(2.0*Pi); //in case of Transient Axi
mu0 = 4.e-7 * Pi ;
nu [ Region[{Air, Ind, Ind1, AirInf} ] ] = 1. / mu0 ;
T_ref = 20; // Temperature de reference 20?C
T_kelvin = 273; // 0?C = 273 K
Wiedemann[ Ind ] = 2.4e-8;
Wiedemann[ Ind1 ] = 2.4e-8;
Alpha[ Ind ] = 3.9e-3;
Alpha[ Ind1 ] = 3.9e-3;
sigma [ Ind ] = 5.9e7;
sigma [ Ind1 ] = 5.9e7;
rho[ Region[{Ind, Ind1}] ]= 1;
Cp[ Region[{Ind, Ind1}] ] = 1;
k[ Region[{Ind, Ind1}] ] = Wiedemann[] * (20+T_kelvin) * 5.9e7 / (1 + Alpha[] * (20-T_ref)) ;
T_Top = 20 + T_kelvin;
T_Bottom = 40 + T_kelvin;
T_ambiant[ Region[{IntSurfaceInd,IntSurfaceInd1}] ] = (T_Bottom-T_Top)*(Y[]+25.e-3)/(50.e-3);
T_ambiant[ Region[{ExtSurfaceInd,ExtSurfaceInd1}] ] = (T_Bottom-T_Top)*(Y[]+25.e-3)/(50.e-3);
T_ambiant[ Region[{BottomSurfaceInd,BottomSurfaceInd1}] ] = T_Bottom;
T_ambiant[ Region[{TopSurfaceInd,TopSurfaceInd1}] ] = T_Top;
h[ Region[{BottomSurfaceInd,TopSurfaceInd}] ] = 1000;
h[ Region[{BottomSurfaceInd1,TopSurfaceInd1}] ] = 1000;
h[ Region[{IntSurfaceInd,ExtSurfaceInd}] ] = 1000;
h[ Region[{IntSurfaceInd1,ExtSurfaceInd1}] ] = 1000;
//
Surface[ Region[{Ind}] ] = 25.e-3*50e-3;
Surface[ Region[{Ind1}] ] = 12.5e-3*25e-3;
//Current[ Region[{Ind}] ] = -135068.76;
//Current[ Region[{Ind1}] ] = -67534.382;
// Aubert's Value
Current[ Region[{Ind}] ] = -134145.81;
Current[ Region[{Ind1}] ] = -67073;
Distribution[ Region[{Ind}] ] = Current[]*75.e-3 / ((75.e-3 * Log[100./75.] * 50.e-3) * X[]);
Distribution[ Region[{Ind1}] ] = Current[]*37.5e-3 / ((37.5e-3 * Log[50./37.5] * 25.e-3) * X[]);
//Parameters for Transient Analysis
Nbr_Time_Steps = 100;
TC = 1.0;
Mag_time0 = 0;
Mag_TimeMax = 0.5*TC;
Mag_DTime[] = (Mag_TimeMax-Mag_time0)/Nbr_Time_Steps;
Mag_Theta[] = 0.5;
Tau = 0.1*TC;
Time_Fct_Ramp[] = ($Time < Tau) ? $Time : Tau;
Time_Fct_Step[] = ($Time > 0) ? 1 : 0;
//Parameters for NonLinear Thermal problems
NL_Eps = 1.e-4; NL_Relax = 1.; NL_NbrMax = 1600;
}
Constraint {
{ Name MagneticVectorPotential_Axi ;
Case {
{ Region SurfaceAxe ; Value 0. ; }
{ Region SurfaceGInf ; Value 0. ; }
}
}
// Set to Zero if transient analysis else to bitter distribution
{ Name SourceCurrentDensityZ ;
Case {
// To use if Static analysis
{ Region DomainS_Mag ; Value Distribution[] ; }
}
}
Val_I_ = 0.01 * 1000. ;
Val_U_ = 1 ;
{ Name Current_Axi ;
Case {
//{ Region Ind ; Value Val_I_; }
//{ Region Ind1 ; Value Val_I_; }
}
}
{ Name Voltage_Axi ;
Case {
// To use if transient analysis with U(t)
//{ Region Ind ; Value Val_U_; TimeFunction Time_Fct_Ramp[];}
// To use if transient analysis with U step function
{ Region Ind ; Value Val_U_; TimeFunction Time_Fct_Step[];}
{ Region Ind1 ; Value Val_U_; TimeFunction Time_Fct_Step[];}
}
}
{ Name FixedTemperature ; Type Assign;
Case {
//{ Region Region[ {BottomSurfaceInd,TopSurfaceInd}] ; Value 0. ; }
}
}
{ Name SourceHeatDensity ; Type Assign;
Case {
//{ Region Ind ; Value 100 ; }
//{ Region Ind1 ; Value 100 ; }
}
}
}
Include "Jacobian_Lib.pro"
Include "Integration_Lib.pro"
Include "MagDyn_av_Axi.pro"
Include "MagSta_a_Axi.pro"
Include "Heat_Axi.pro"
PostOperation Map_av UsingPost MagDyn_t_av_Axi {
Print[ az, OnElementsOf Domain_Mag, File "serie_MagDyn_a.pos", TimeStep {1:Nbr_Time_Steps}] ;
Print[ jz, OnElementsOf DomainC_Mag, File "serie_MagDyn_j.pos", TimeStep {1:Nbr_Time_Steps}] ;
Print[ roj2, OnElementsOf DomainS_Th, File "serie_MagDyn_Joules.pos", TimeStep {1:Nbr_Time_Steps}] ;
Print[V, OnLine {{37.5e-3,0,0}{50.e-3,0,0}}{10}, File "serie_MagDyn_V1.pos", Format Gnuplot, TimeStep {Nbr_Time_Steps:Nbr_Time_Steps}];
Print[V, OnLine {{75.e-3,0,0}{100.e-3,0,0}}{10}, File "serie_MagDyn_V2.pos", Format Gnuplot, TimeStep {Nbr_Time_Steps:Nbr_Time_Steps}];
}
PostOperation Map_T UsingPost HeatDyn_Axi {
Print[ T, OnElementsOf Ind, File "serie_HeatDyn_ind_T.pos"] ;
Print[ T, OnElementsOf Ind1, File "serie_HeatDyn_ind1_T.pos"] ;
Print[ Source_Th, OnElementsOf DomainS_Th, File "serie_HeatDyn_Qth.pos", TimeStep {1:Nbr_Time_Steps}] ;
Print[T, OnLine {{37.5e-3,0,0}{50.e-3,0,0}}{10}, File "serie_HeatDyn_t1.pos", Format Table, TimeStep {1:Nbr_Time_Steps}];
Print[T, OnLine {{75.e-3,0,0}{100.e-3,0,0}}{10}, File "serie_HeatDyn_t2.pos", Format Table, TimeStep {1:Nbr_Time_Steps}];
Print[ az, OnElementsOf Domain_Mag, File "serie_HeatDyn_a.pos", TimeStep {1:Nbr_Time_Steps}] ;
Print[ jz, OnElementsOf DomainC_Mag, File "serie_HeatDyn_j.pos", TimeStep {1:Nbr_Time_Steps}] ;
Print[ roj2, OnElementsOf DomainC_Mag, File "serie_HeatDyn_Joules.pos", TimeStep {1:Nbr_Time_Steps}] ;
}
PostOperation U_av UsingPost MagDyn_t_av_Axi {
Print[ R, OnRegion #{Ind, Ind1},File "serie_MagDyn_R.pos" , Format Table, TimeStep {Nbr_Time_Steps:Nbr_Time_Steps}] ;
Print[ U, OnRegion #{Ind, Ind1},File "serie_MagDyn_U.pos" , Format TimeTable] ;
Print[ I, OnRegion #{Ind, Ind1},File "serie_MagDyn_I.pos" , Format TimeTable] ;
Print[ b, OnLine {{1.e-06,0,0}{250.e-3, 0, 0}} {50}, File "serie_MagDyn_b.pos" , Format Gnuplot, TimeStep {Nbr_Time_Steps:Nbr_Time_Steps}] ;
Print[jz, OnLine {{37.5e-3,0,0}{50.e-3,0,0}}{10}, File "serie_MagDyn_j1.pos", Format Gnuplot, TimeStep {Nbr_Time_Steps:Nbr_Time_Steps}];
Print[jz, OnLine {{75.e-3,0,0}{100.e-3,0,0}}{10}, File "serie_MagDyn_j2.pos", Format Gnuplot, TimeStep {Nbr_Time_Steps:Nbr_Time_Steps}];
}
-------------- next part --------------
/* -------------------------------------------------------------------
File "Core.geo"
This file is the geometrical description used by GMSH to produce
the file "Core.msh".
------------------------------------------------------------------- */
r_int_Ind = 75.e-3 ;
r_ext_Ind = 100.e-3 ;
h_Ind = 25.e-3;
r_int_Ind1 = 37.5e-3 ;
r_ext_Ind1 = 50.e-3 ;
h_Ind1 = 12.5e-3;
rInt = 400.e-3 ;
rExt = 500.0001e-3;
s = 1. ;
p0 = 6.e-3 *s;
pCorex = 2.e-3 *s;
pCorey0 = 4.e-3 *s;
pCorey = 2.e-3 *s;
pIndx = 1.25e-3 *s;
pIndy = 1.25e-3 *s;
pInd1x = 1.25e-3 *s;
pInd1y = 1.25e-3 *s;
pInt = 6.25e-3 *s;
pExt = 6.25e-3 *s;
Point(1) = {0,0,0,pCorey};
Point(2) = {0,h_Ind,0,pCorey};
Point(3) = {0,-h_Ind,0,pCorey};
// First inducteur
Point(4) = {r_int_Ind,0,0,pIndx};
Point(5) = {r_ext_Ind,0,0,pIndx};
Point(6) = {r_ext_Ind,h_Ind,0,pIndx};
Point(7) = {r_int_Ind,h_Ind,0,pIndx};
Point(8) = {r_int_Ind,-h_Ind,0,pIndx};
Point(9) = {r_ext_Ind,-h_Ind,0,pIndx};
Point(10) = {rInt,0,0,pInt};
Point(11) = {rExt,0,0,pInt};
Point(12) = {0,rInt,0,pInt};
Point(13) = {0,rExt,0,pInt};
Point(14) = {0,-rInt,0,pInt};
Point(15) = {0,-rExt,0,pInt};
// First inducteur
Point(16) = {0,h_Ind1,0,pCorey};
Point(17) = {0,-h_Ind1,0,pCorey};
Point(18) = {r_int_Ind1,0,0,pInd1x};
Point(19) = {r_ext_Ind1,0,0,pInd1x};
Point(20) = {r_ext_Ind1,h_Ind1,0,pInd1x};
Point(21) = {r_int_Ind1,h_Ind1,0,pInd1x};
Point(22) = {r_int_Ind1,-h_Ind1,0,pInd1x};
Point(23) = {r_ext_Ind1,-h_Ind1,0,pInd1x};
Line(1) = {1,4};
Line(2) = {4,5};
Line(3) = {5,10};
Line(4) = {10,11};
Line(5) = {5,6};
Line(6) = {6,7};
Line(7) = {7,4};
Line(8) = {4,8};
Line(9) = {8,9};
Line(10) = {9,5};
Line(11) = {7,2};
Line(12) = {2,1};
Line(13) = {1,3};
Line(14) = {3,8};
Line(15) = {2,12};
Line(16) = {12,13};
Line(17) = {3,14};
Line(18) = {14,15};
Circle(19) = {14,1,10};
Circle(20) = {15,1,11};
Circle(21) = {10,1,12};
Circle(22) = {11,1,13};
Line(23) = {1,18};
Line(24) = {18,19};
Line(25) = {19,4};
Line(26) = {7,20};
Line(27) = {19,20};
Line(28) = {20,21};
Line(29) = {21,16};
Line(30) = {21,18};
Line(31) = {16,1};
Line(32) = {8,23};
Line(33) = {23,22};
Line(34) = {19,23};
Line(35) = {22,18};
Line(36) = {22,17};
Line(37) = {17,1};
Line(38) = {17,3};
Line(39) = {2,16};
Line Loop(100) = {8,9,10,5,6,7};
Plane Surface(130) = {100}; /* Inducteur */
Line Loop(101) = {-35,-33,-34,27,28,30};
Plane Surface(131) = {101}; /* Inducteur 1*/
Line Loop(102) = {-32,-8,-7,26,-27,34};
Line Loop(103) = {38,14,32,33,36};
Line Loop(104) = {-37,-36,35,-23};
Line Loop(105) = {31,23,-30,29};
Line Loop(106) = {39,-29,-28,-26,11};
Line Loop(107) = {17,19,-3,-10,-9,-14};
Line Loop(108) = {-15,-11,-6,-5,3,21};
Plane Surface(132) = {102}; // Air
Plane Surface(133) = {103}; // Air
Plane Surface(134) = {104}; // Air
Plane Surface(135) = {105}; // Air
Plane Surface(136) = {106}; // Air
Plane Surface(137) = {107}; // Air
Plane Surface(138) = {108}; // Air
Line Loop(109) = {18,20,-4,-19};
Line Loop(110) = {-16,-21,4,22};
Plane Surface(139) = {109}; // Infini
Plane Surface(140) = {110}; // Infini
// To get quadrangles
//Recombine Surface {139,140};
Physical Surface(240) = {130};
Physical Surface(241) = {131};
Physical Surface(242) = {132,133,134,135,136,137,138};
Physical Surface(243) = {139,140};
Physical Line (244) = {-16,-15,39,31,-37,38,17,18};
Physical Line (245) = {20,22};
Physical Line (246) = {7,8}; //IntSurfaceInd
Physical Line (247) = {9}; //BottomSurfaceInd
Physical Line (248) = {10,5}; //ExtSurfaceInd
Physical Line (249) = {6}; //TopSurfaceInd
Physical Line (250) = {30,-35}; //IntSurfaceInd1
Physical Line (251) = {-33}; //BottomSurfaceInd1
Physical Line (252) = {-34,27}; //ExtSurfaceInd1
Physical Line (253) = {28}; //TopSurfaceInd1
-------------- next part --------------
/* -------------------------------------------------------------------
File "CoreMassive.pro"
This file defines the problem dependent data structures for the
dynamic core-inductor problem.
To compute the solution: getdp CoreMassive -msh core.msh -solve MagDyn_av_Axi
To compute post-results: getdp CoreMassive -msh core.msh -pos Map_a
or getdp CoreMassive -msh core.msh -pos U_av
------------------------------------------------------------------- */
Group {
Ind = Region[ 240 ];
Ind1 = Region[ 241 ];
Air = Region[ 242 ];
AirInf = Region[ 243 ];
SurfaceAxe = Region[ 244 ];
SurfaceGInf = Region[ 245 ];
IntSurfaceInd = Region[ 246 ];
BottomSurfaceInd = Region[ 247 ];
ExtSurfaceInd = Region[ 248 ];
TopSurfaceInd = Region[ 249 ];
IntSurfaceInd1 = Region[ 250 ];
BottomSurfaceInd1 = Region[ 251 ];
ExtSurfaceInd1 = Region[ 252 ];
TopSurfaceInd1 = Region[ 253 ];
// Domains for Thermal Problems
SurfaceHeatFlux = Region[ {}];
SurfaceConvection = Region[ {IntSurfaceInd,ExtSurfaceInd,TopSurfaceInd,BottomSurfaceInd,
IntSurfaceInd1,ExtSurfaceInd1,TopSurfaceInd1,BottomSurfaceInd1}];
//SurfaceHeatFlux = Region[{}];
DomainC_Th = Region[ {}];
DomainS_Th = Region[ {Ind, Ind1}];
Domain_Th = Region[ {DomainS_Th, DomainC_Th}];
// Domains for Electromagnetics Problems
DomainCC_Mag = Region[ {Air, AirInf} ] ;
/* Massive inductor (Ind), in DomainC_Mag */
DomainB_Mag = Region[ {} ] ;
DomainS_Mag = Region[ {} ] ;
// if Static analysis
DomainS_Mag = Region[ {DomainS_Th} ] ;
DomainC_Mag = Region[ {} ] ;
// if transient analysis
//DomainS_Mag = Region[ {} ] ;
//DomainC_Mag = Region[ {DomainS_Th} ] ;
DomainInf = Region[ {AirInf} ] ;
Val_Rint = 399.999e-3 ;
Val_Rext = 500.001e-3 ;
Domain_Mag = Region[ {DomainCC_Mag, DomainC_Mag, DomainS_Mag} ] ;
}
Function {
Factor [ Region[{DomainC_Mag} ] ] = 1.0/(2.0*Pi); //in case of Transient Axi
mu0 = 4.e-7 * Pi ;
nu [ Region[{Air, Ind, Ind1, AirInf} ] ] = 1. / mu0 ;
T_ref = 20; // Temperature de reference 20?C
T_kelvin = 273; // 0?C = 273 K
Wiedemann[ Ind ] = 2.4e-8;
Wiedemann[ Ind1 ] = 2.4e-8;
Alpha[ Ind ] = 3.9e-3;
Alpha[ Ind1 ] = 3.9e-3;
sigma [ Ind ] = 5.9e7;
sigma [ Ind1 ] = 5.9e7;
rho[ Region[{Ind, Ind1}] ]= 1;
Cp[ Region[{Ind, Ind1}] ] = 1;
k[ Region[{Ind, Ind1}] ] = Wiedemann[] * (20+T_kelvin) * 5.9e7 / (1 + Alpha[] * (20-T_ref)) ;
T_Top = 20 + T_kelvin;
T_Bottom = 40 + T_kelvin;
T_ambiant[ Region[{IntSurfaceInd,IntSurfaceInd1}] ] = (T_Bottom-T_Top)*(Y[]+25.e-3)/(50.e-3);
T_ambiant[ Region[{ExtSurfaceInd,ExtSurfaceInd1}] ] = (T_Bottom-T_Top)*(Y[]+25.e-3)/(50.e-3);
T_ambiant[ Region[{BottomSurfaceInd,BottomSurfaceInd1}] ] = T_Bottom;
T_ambiant[ Region[{TopSurfaceInd,TopSurfaceInd1}] ] = T_Top;
h[ Region[{BottomSurfaceInd,TopSurfaceInd}] ] = 1000;
h[ Region[{BottomSurfaceInd1,TopSurfaceInd1}] ] = 1000;
h[ Region[{IntSurfaceInd,ExtSurfaceInd}] ] = 1000;
h[ Region[{IntSurfaceInd1,ExtSurfaceInd1}] ] = 1000;
//
Surface[ Region[{Ind}] ] = 25.e-3*50e-3;
Surface[ Region[{Ind1}] ] = 12.5e-3*25e-3;
//Current[ Region[{Ind}] ] = -135068.76;
//Current[ Region[{Ind1}] ] = -67534.382;
// Aubert's Value
Current[ Region[{Ind}] ] = -134145.81;
Current[ Region[{Ind1}] ] = -67073;
Distribution[ Region[{Ind}] ] = Current[]*75.e-3 / ((75.e-3 * Log[100./75.] * 50.e-3) * X[]);
Distribution[ Region[{Ind1}] ] = Current[]*37.5e-3 / ((37.5e-3 * Log[50./37.5] * 25.e-3) * X[]);
//Parameters for Transient Analysis
Nbr_Time_Steps = 100;
TC = 1.0;
Mag_time0 = 0;
Mag_TimeMax = 0.5*TC;
Mag_DTime[] = (Mag_TimeMax-Mag_time0)/Nbr_Time_Steps;
Mag_Theta[] = 1;
Tau = 0.1*TC;
Time_Fct_Ramp[] = ($Time < Tau) ? $Time : Tau;
//Parameters for NonLinear Thermal problems
NL_Eps = 1.e-4; NL_Relax = 1.; NL_NbrMax = 1600;
}
Constraint {
{ Name MagneticVectorPotential_Axi ;
Case {
{ Region SurfaceAxe ; Value 0. ; }
{ Region SurfaceGInf ; Value 0. ; }
}
}
// Set to Zero if transient analysis else to bitter distribution
{ Name SourceCurrentDensityZ ;
Case {
// To use if Static analysis
{ Region DomainS_Mag ; Value Distribution[] ; }
//{ Region DomainS_Mag ; Value 0.0 ; }
}
}
Val_I_ = 0.01 * 1000. ;
Val_U_ = 1 ;
{ Name Current_Axi ;
Case {
//{ Region Ind ; Value Val_I_; }
//{ Region Ind1 ; Value Val_I_; }
}
}
{ Name Voltage_Axi ;
Case {
// To use if transient analysis with U(t)
//{ Region Ind ; Value Val_U_; TimeFunction Time_Fct_Ramp[];}
// To use if transient analysis with U step function
//{ Region Ind ; Value Val_U_;}
//{ Region Ind1 ; Value Val_U_;}
}
}
{ Name FixedTemperature ; Type Assign;
Case {
//{ Region Region[ {BottomSurfaceInd,TopSurfaceInd}] ; Value 0. ; }
}
}
{ Name SourceHeatDensity ; Type Assign;
Case {
//{ Region Ind ; Value 100 ; }
//{ Region Ind1 ; Value 100 ; }
{ Region DomainS_Mag ; Value SquNorm[Distribution[]]/sigma[] ; }
}
}
}
Include "Jacobian_Lib.pro"
Include "Integration_Lib.pro"
Include "MagDyn_av_Axi.pro"
Include "MagSta_a_Axi.pro"
Include "Heat_Axi.pro"
PostOperation Map_a UsingPost MagSta_a_Axi {
Print[ a_theta, OnElementsOf Domain_Mag, File "serie_magsta_a.pos"] ;
Print[ j_theta, OnElementsOf DomainS_Mag, File "serie_magsta_j.pos"] ;
Print[ j_source, OnElementsOf DomainS_Mag, File "serie_magsta_js.pos"] ;
Print[ roj2, OnElementsOf DomainS_Mag, File "serie_magsta_Joules.pos"] ;
Print[j_source, OnLine {{37.5e-3,0,0}{50.e-3,0,0}}{10}, File "serie_magsta_j1.pos", Format Gnuplot];
Print[j_source, OnLine {{75.e-3,0,0}{100.e-3,0,0}}{10}, File "serie_magsta_j2.pos", Format Gnuplot];
Print[j_theta, OnLine {{37.5e-3,0,0}{50.e-3,0,0}}{10}, File "serie_magsta_jt1.pos", Format Gnuplot];
Print[j_theta, OnLine {{75.e-3,0,0}{100.e-3,0,0}}{10}, File "serie_magsta_jt2.pos", Format Gnuplot];
Print[ b, OnLine {{1.e-06,0,0}{250.e-3, 0, 0}} {50}, File "serie_magsta_b.pos" , Format Gnuplot ] ;
}
PostOperation Map_T UsingPost HeatSta_Axi {
Print[ T, OnElementsOf Ind, File "serie_HeatSta_ind_T.pos"] ;
Print[ T, OnElementsOf Ind1, File "serie_HeatSta_ind1_T.pos"] ;
Print[ Source_Th, OnElementsOf DomainS_Th, File "serie_t_Qth.pos"] ;
Print[T, OnLine {{37.5e-3,0,0}{50.e-3,0,0}}{10}, File "serie_HeatSta_t1.pos", Format Table];
Print[T, OnLine {{75.e-3,0,0}{100.e-3,0,0}}{10}, File "serie_HeatSta_t2.pos", Format Table];
//Print[ roj2, OnElementsOf DomainC_Mag, File "serie_t_Joules.pos"] ;
}
-------------- next part --------------
/* -------------------------------------------------------------------
File "MagSta_a_Axi.pro"
Magnetostatics - Magnetic vector potential a formulation (Axi)
-------------------------------------------------------------------
I N P U T
---------
GlobalGroup : (Extension '_Mag' is for Magnetic problem)
-----------
Domain_Mag Whole magnetic domain
DomainS_Mag Inductor regions (Source)
Function :
--------
nu[] Magnetic reluctivity
Constraint :
----------
MagneticVectorPotential_Axi
Fixed magnetic vector potential (Axi)
(classical boundary condition)
SourceCurrentDensityY Fixed source current density (in Y direction)
*/
Group {
DefineGroup[ Domain_Mag, DomainS_Mag, DomainC_Mag ] ;
}
Function {
DefineFunction[ nu ] ;
DefineFunction[ Distribution ] ;
}
FunctionSpace {
// Magnetic vector potential a (b = curl a)
{ Name Hcurl_a_MagSta_Axi ; Type Form1P ;
BasisFunction {
// a = a s
// e e
{ Name se ; NameOfCoef ae ; Function BF_PerpendicularEdge ;
Support Domain_Mag ; Entity NodesOf[ All ] ; }
}
Constraint {
{ NameOfCoef ae ; EntityType NodesOf ;
NameOfConstraint MagneticVectorPotential_Axi ; }
}
}
// Source current density js (fully fixed space)
{ Name Hregion_j_MagSta_Axi ; Type Vector ;
BasisFunction {
{ Name sr ; NameOfCoef jsr ; Function BF_NodeZ ;
Support DomainS_Mag ; Entity NodesOf[DomainS_Mag ]; }
}
Constraint {
{ NameOfCoef jsr ; EntityType NodesOf ;
NameOfConstraint SourceCurrentDensityZ ; }
}
}
// // for constant Current Densities
// { Name Hregion_j_MagSta_Axi ; Type Vector ;
// BasisFunction {
// { Name sr ; NameOfCoef jsr ; Function BF_RegionZ ;
// Support DomainS_Mag ; Entity DomainS_Mag; }
// }
// Constraint {
// { NameOfCoef jsr ; EntityType Region ;
// NameOfConstraint SourceCurrentDensityZ ; }
// }
// }
}
Formulation {
{ Name Magnetostatics_a_Axi ; Type FemEquation ;
Quantity {
{ Name a ; Type Local ; NameOfSpace Hcurl_a_MagSta_Axi ; }
{ Name js ; Type Local ; NameOfSpace Hregion_j_MagSta_Axi ; }
}
Equation {
Galerkin { [ nu[] * Dof{d a} , {d a} ] ; In Domain_Mag ;
Jacobian VolAxi ; Integration CurlCurl ; }
Galerkin { [ - Dof{js} , {a} ] ; In DomainS_Mag ;
Jacobian VolAxi ; Integration CurlCurl ; }
//Galerkin { [ - Vector[0,0,Distribution[]] , {a} ] ; In DomainS_Mag ; //just to check
// Jacobian VolAxi ; Integration CurlCurl ; }
}
}
}
Resolution {
{ Name MagSta_a_Axi ;
System {
{ Name Sys_Mag ; NameOfFormulation Magnetostatics_a_Axi ; }
}
Operation {
Generate Sys_Mag ; Solve Sys_Mag ; SaveSolution Sys_Mag ;
}
}
}
PostProcessing {
{ Name MagSta_a_Axi ; NameOfFormulation Magnetostatics_a_Axi ;
PostQuantity {
{ Name a ; Value { Local { [ {a} ] ; In Domain_Mag ; Jacobian VolAxi ; } } }
{ Name a_theta ; Value { Local { [ CompZ[{a}] ] ; In Domain_Mag ; Jacobian VolAxi ; } } }
{ Name b ; Value { Local { [ {d a} ] ; In Domain_Mag ; Jacobian VolAxi ; } } }
{ Name h ; Value { Local { [ nu[] * {d a} ] ; In Domain_Mag ; Jacobian VolAxi ;} } }
{ Name j_theta ; Value { Local { [ CompZ[{js}] ] ; In DomainS_Mag ; Jacobian VolAxi ; } } }
{ Name j_source ; Value { Local { [ Distribution[] ] ; In DomainS_Mag ; Jacobian VolAxi ; } } }
{ Name roj2 ; Value { Local { [ SquNorm[CompZ[{js}]] / sigma[] ] ; In DomainS_Mag ; Jacobian VolAxi ; } } }
}
}
}