GetDP
Samuel Kvasnica
kvasnica at iaee.tuwien.ac.at
Tue Oct 5 20:01:29 CEST 1999
Hi Christophe,
Thank you very much for answering my questions.
First let me congratulate you ! I was playing with getDP last week and found it very
useful - I like it's
structure, really a good idea. And it's pretty stable. Well, gmsh is another beast
... but GUI programing
was never a grateful task. After playing with manual examples and some easy
electrostatic problems like
dipole etc. I tried to simulate magnetic field of our magnetron consisting of iron
core with 5 inductor coils
(each other has oposite direction of current flow). See attached picture with z
component of vector potential Az.
The form of the field is almost correct, but not completely. The symetry should be
better in my opinion.
I think its a problem of boundary conditions in infinity. I'd like to ask you
following questions:
A assumed, in cylindrical coordinates x corresponds to r and y to z coordinate, x
always >=0. Correct ?
Do I have to set any constraints on axis r=0 when using cylindrical coordinates ?
Z-component of vector potential A seems to be useful when displaying 2D solenoidal
fields. But in case
of permanent magnets, I need to use also scalar potential (Bperm=grad phi). Correct
?
What is the meaning of Spherical shell transformation ? Is it a projection of
<Rint,Rext> on <Rint,infinity>
to set the value of potential or induction to 0 in infinity ? What should I use in
case of axial symetry ? I've just
assigned induction to 0 somewhere in a distance.
In your magnetostatic example MagSta_a_2D, what is the meaning of SurfaceGh0 and
SurfaceGe0 ?
I just noticed A is set to 0 at SurfaceGe0. How getDP recognises the symetry ? Does
setting A to 0
at SurfaceGe0 mean mirror symetry ? What about Gh0 ? Well, I'm really confused...
I'd like to generate a table of B on grid. I don't understand the syntax of OnGrid
Plot[ b, OnGrid {{-0.07,eps,0}{0.09, eps, 0}} {500}, File "b_phi.cut"] ;
Could you explain me please, what are 2 vectors following keyword 'OnGrid'
and what means {500} ?
What is the difference between OnCut and OnGrid ?
What is the meaning of Depth keyword ?
What means '>' in File <'>'> "string" ?
And finally: do you have any experience with use of GID or SolidMesh in combination with getDP ?
Sorry again for having so many questions, but maybe it could inspirate you once you'll write
a FAQ manual :-)
Best regards,
Sam
--
--------------------------------------------------------------
Mag.Samuel Kvasnica
IAEE TU-WIEN
Gusshausstrase 27-29 / 359 email: kvasnica at iaee.tuwien.ac.at
1040 Wien tel: +43-1-58801-35914
Austria fax: +43-1-58801-35997
--------------------------------------------------------------
Content-Type: text/plain; charset=iso-8859-2;
name="electrode.geo"
Content-Transfer-Encoding: 7bit
Content-Disposition: inline;
filename="electrode.geo"
x0 = 3.5;
y0 = 72.5;
tx = 2.5;
ty = 5;
dx = 9.25;
dy = 25;
xInt = 200;
xExt = 220;
yInt = 300;
yExt = 320;
pCore = 4;
pInt = 7;
pExt = 10;
pAxis = 7;
Point(1) = {x0,y0,0,pCore};
Point(2) = {x0+tx+dx,y0,0,pCore};
Point(3) = {x0+tx+dx,y0+ty,0,pCore};
Point(4) = {x0+tx,y0+ty,0,pCore};
Point(5) = {x0+tx,y0+ty+dy,0,pCore};
Point(6) = {x0+tx+dx,y0+ty+dy,0,pCore};
Point(7) = {x0+tx+dx,y0+2*ty+dy,0,pCore};
Point(8) = {x0+tx,y0+2*ty+dy,0,pCore};
Point(9) = {x0+tx,y0+2*ty+2*dy,0,pCore};
Point(10) = {x0+tx+dx,y0+2*ty+2*dy,0,pCore};
Point(11) = {x0+tx+dx,y0+3*ty+2*dy,0,pCore};
Point(12) = {x0+tx,y0+3*ty+2*dy,0,pCore};
Point(13) = {x0+tx,y0+3*ty+3*dy,0,pCore};
Point(14) = {x0+tx+dx,y0+3*ty+3*dy,0,pCore};
Point(15) = {x0+tx+dx,y0+4*ty+3*dy,0,pCore};
Point(16) = {x0+tx,y0+4*ty+3*dy,0,pCore};
Point(17) = {x0+tx,y0+4*ty+4*dy,0,pCore};
Point(18) = {x0+tx+dx,y0+4*ty+4*dy,0,pCore};
Point(19) = {x0+tx+dx,y0+5*ty+4*dy,0,pCore};
Point(21) = {0,0,0,pAxis};
Point(22) = {xInt,0,0,pInt};
Point(23) = {xExt,-20,0,pExt};
Point(24) = {xExt,yExt,0,pExt};
Point(25) = {0,yExt,0,pExt};
Point(26) = {0,yInt,0,pAxis};
Point(27) = {xInt,yInt,0,pInt};
Point(28) = {0,-20,0,pExt};
Point(29) = {0,y0,0,pAxis};
Point(31) = {x0+tx,y0+5*ty+4*dy,0,pCore};
Point(32) = {x0+tx,y0+5*ty+5*dy,0,pCore};
Point(33) = {x0+tx+dx,y0+5*ty+5*dy,0,pCore};
Point(34) = {x0+tx+dx,y0+6*ty+5*dy,0,pCore};
Point(35) = {x0,y0+6*ty+5*dy,0,pCore};
Point(36) = {0,y0+6*ty+5*dy,0,pAxis};
Line(1) = {16,15};
Line(2) = {15,18};
Line(3) = {18,17};
Line(4) = {17,16};
Line Loop(5) = {2,3,4,1};
Plane Surface(6) = {5}; /* ind 1 */
Line(7) = {12,11};
Line(8) = {11,14};
Line(9) = {14,13};
Line(10) = {13,12};
Line Loop(11) = {8,9,10,7};
Plane Surface(12) = {11}; /* ind 2 */
Line(13) = {8,7};
Line(14) = {7,10};
Line(15) = {10,9};
Line(16) = {9,8};
Line Loop(17) = {14,15,16,13}; /* ind 3 */
Plane Surface(18) = {17};
Line(19) = {4,3};
Line(20) = {3,6};
Line(21) = {6,5};
Line(22) = {5,4};
Line Loop(23) = {20,21,22,19}; /*ind 4*/
Plane Surface(24) = {23};
Line(25) = {1,2};
Line(26) = {2,3};
Line(27) = {6,7};
Line(28) = {10,11};
Line(29) = {14,15};
Line(30) = {18,19};
Line(1105) = {29,1};
Line(1107) = {25,26};
Line(1108) = {25,24};
Line(1109) = {24,23};
Line(1110) = {23,28};
Line(1111) = {28,21};
Line(1112) = {21,29};
Line(1113) = {21,22};
Line(1114) = {22,27};
Line(1115) = {27,26};
Line Loop(1120) = {-1107,1108,1109,1110,1111,1113,1114,1115};
Plane Surface(1121) = {1120}; /* infinity */
Line(1122) = {19,31};
Line(1123) = {31,32};
Line(1124) = {32,33};
Line(1125) = {33,34};
Line(1126) = {19,33};
Line(1127) = {34,35};
Line(1128) = {35,1};
Line(1129) = {29,36};
Line(1130) = {35,36};
Line(1131) = {36,26};
Line Loop(1132) = {-1128,1130,-1129,1105};
Plane Surface(1133) = {1132}; /* Air 1 */
Line Loop(1134) = {-1124,-1123,-1122,1126};
Plane Surface(1135) = {1134}; /* ind 5 */
Line Loop(1136) = {25,26,-19,-22,-21,27,-13,-16,-15,28,-7,-10,-9,29,-1,-4,-3,30,1122,1123,1124,1125,1127,1128};
Plane Surface(1137) = {1136}; /* Core */
Line Loop(1138) = {25,26,20,27,14,28,8,29,2,30,1126,1125,1127,1130,1131,-1115,-1114,-1113,1112,1105};
Plane Surface(1139) = {1138}; /* Air 2 */
Point(199)={0,150,0,pAxis};
Physical Surface(103) = {1137}; /* core */
/*Physical Surface(102) = {6,12,18,24};*/ /* inductor */
Physical Surface(101) = {1133,1139}; /* air */
Physical Surface(111) = {1121}; /* infinity air */
Physical Surface(105) = {6}; /* inductor 1 */
Physical Surface(106) = {12}; /* inductor 2 */
Physical Surface(107) = {18}; /* inductor 3 */
Physical Surface(108) = {24}; /* inductor 4 */
Physical Surface(109) = {1135}; /* inductor 5 */
Physical Line(1102) = {28,23,24,25}; /* SurfaceGh0 */
-------------- next part --------------
/* -------------------------------------------------------------------
File "electrode.pro"
This file defines the problem dependant data structures for the
static inductor problem.
To compute the solution: getdp CoreSta -msh Core.msh -solve MagSta_a_2D
To compute post-results: getdp CoreSta -msh Core.msh -pos Map_a
or getdp CoreSta -msh Core.msh -pos Cutb_a
------------------------------------------------------------------- */
Group {
Air = # 101 ;
Ind1 = # 105 ;
Ind2 = # 106 ;
Ind3 = # 107 ;
Ind4 = # 108 ;
Ind5 = # 109 ;
/* Ind = # 102 ;*/
AirInf = # 111 ; Core = #103;
/* SurfaceGh0 = # 1100 ; SurfaceGe0 = # 1101 ;*/
SurfaceGInf = # 1102 ;
DomainCC_Mag = Region[ {Air, AirInf, Ind1,Ind2,Ind3,Ind4,Ind5, Core} ] ;
DomainC_Mag = Region[ {} ] ;
DomainS_Mag = Region[ {Ind1,Ind2,Ind3,Ind4,Ind5} ] ; // Stranded inductor
DomainInf = Region[ {AirInf} ] ;
Val_Rint = 200 ; Val_Rext = 400 ;
Domain_Mag = Region[ {DomainCC_Mag, DomainC_Mag} ] ;
}
Function {
mu0 = 4.e-7 * Pi ;
murCore = 1000.;
nu [ #{Air, Ind1,Ind2,Ind3,Ind4,Ind5, AirInf} ] = 1. / mu0 ;
nu [ Core ] = 1./(murCore*mu0);
Sc[ #{Ind1,Ind2,Ind3,Ind4,Ind5} ] = 2.5e-2 * 5.e-2 ;
}
Constraint {
{ Name MagneticVectorPotential_2D ;
Case {
/* { Region SurfaceGe0 ; Value 0. ; }*/
/* { Region SurfaceGh0 ; Value 0. ; }*/
{ Region SurfaceGInf ; Value 0. ; }
}
}
Val_I_1_ = 0.01 * 1000. ;
{ Name SourceCurrentDensityZ ;
Case {
{ Region Ind1 ; Value Val_I_1_/Sc[] ; }
{ Region Ind2 ; Value -Val_I_1_/Sc[] ; }
{ Region Ind3 ; Value Val_I_1_/Sc[] ; }
{ Region Ind4 ; Value -Val_I_1_/Sc[] ; }
{ Region Ind5 ; Value -Val_I_1_/Sc[] ; }
}
}
}
Include "Jacobian2_Lib.pro"
Include "Integration_Lib.pro"
Include "MagSta_a_2D.pro"
PostOperation Map_a UsingPost MagSta_a_2D {
Plot[ az, OnRegion Domain_Mag, File "Electrode_az.pos"] ;
Plot[ a, OnRegion Domain_Mag, File "Electrode_a.pos"] ;
Plot[ b, OnRegion Domain_Mag, File "Electrode_b.pos", Depth 0] ; Plot[ h, OnRegion Domain_Mag, File "Electrode_h.pos", Depth 0] ;
}
e = 1.e-5 ;
PostOperation Cutb_a UsingPost MagSta_a_2D {
Plot[ b, OnLine{{e,e,0}{0.12,e,0}} {1000}, File "k_a"] ;
}
-------------- next part --------------
/* -------------------------------------------------------------------
File "Jacobian_Lib.pro"
Definition of a jacobian method
-------------------------------------------------------------------
I N P U T
---------
GlobalGroup :
-----------
DomainInf Regions with Spherical Shell Transformation
Parameters :
----------
Val_Rint, Val_Rext Inner and outer radius of the Spherical Shell
of DomainInf
*/
Group {
DefineGroup[ DomainInf ] ;
DefineVariable[ Val_Rint, Val_Rext ] ;
}
JacobianMethod {
{ Name Vol ;
Case { { Region DomainInf ;
Jacobian VolAxi ; }
{ Region All ; Jacobian VolAxi ; }
}
}
}
-------------- next part --------------
/* -------------------------------------------------------------------
File "MagSta_a_2D.pro"
Magnetostatics - Magnetic vector potential a formulation (2D)
-------------------------------------------------------------------
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_2D
Fixed magnetic vector potential (2D)
(classical boundary condition)
SourceCurrentDensityZ Fixed source current density (in Z direction)
*/
Group {
DefineGroup[ Domain_Mag, DomainS_Mag ] ;
}
Function {
DefineFunction[ nu ] ;
}
FunctionSpace {
// Magnetic vector potential a (b = curl a)
{ Name Hcurl_a_Mag_2D ; Type Form1P ;
BasisFunction {
// a = a s
// e e
{ Name se ; NameOfCoef ae ; Function BF_PerpendicularEdge ;
Support Domain_Mag ; Entity NodesOf[ All ] ; }
}
Constraint {
{ NameOfCoef ae ;
Entity NodesOf UsingConstraint MagneticVectorPotential_2D ; }
}
}
// Source current density js (fully fixed space)
{ Name Hregion_j_Mag_2D ; Type Vector ;
BasisFunction {
{ Name sr ; NameOfCoef jsr ; Function BF_RegionZ ;
Support DomainS_Mag ; Entity DomainS_Mag ; }
}
Constraint {
{ NameOfCoef jsr ;
Entity Region UsingConstraint SourceCurrentDensityZ ; }
}
}
}
Formulation {
{ Name Magnetostatics_a_2D ; Type FemEquation ;
Quantity {
{ Name a ; Type Local ; NameOfSpace Hcurl_a_Mag_2D ; }
{ Name js ; Type Local ; NameOfSpace Hregion_j_Mag_2D ; }
}
Equation {
Galerkin { [ nu[] * Dof{d a} , {d a} ] ; In Domain_Mag ;
Jacobian Vol ; Integration CurlCurl ; }
Galerkin { [ - Dof{js} , {a} ] ; In DomainS_Mag ;
Jacobian Vol ; Integration CurlCurl ; }
}
}
}
Resolution {
{ Name MagSta_a_2D ;
System {
{ Name Sys_Mag ; NameOfFormulation Magnetostatics_a_2D ; }
}
Operation {
Generate Sys_Mag ; Solve Sys_Mag ; SaveSolution Sys_Mag ;
}
}
}
PostProcessing {
{ Name MagSta_a_2D ; NameOfFormulation Magnetostatics_a_2D ;
PostQuantity {
{ Name a ; Value { Term { [ {a} ] ; In Domain_Mag ; } } }
{ Name az ; Value { Term { [ CompZ[{a}] ] ; In Domain_Mag ; } } }
{ Name b ; Value { Term { [ {d a} ] ; In Domain_Mag ; } } }
{ Name h ; Value { Term { [ nu[] * {d a} ] ; In Domain_Mag ; } } }
}
}
}
-------------- next part --------------
/* -------------------------------------------------------------------
File "Integration_Lib.pro"
Definition of integration methods
------------------------------------------------------------------- */
Integration {
{ Name GradGrad ;
Case { {Type Gauss ;
Case { { GeoElement Triangle ; NumberOfPoints 4 ; }
{ GeoElement Quadrangle ; NumberOfPoints 4 ; }
{ GeoElement Tetrahedron ; NumberOfPoints 4 ; }
{ GeoElement Hexahedron ; NumberOfPoints 6 ; }
{ GeoElement Prism ; NumberOfPoints 9 ; } }
}
}
}
{ Name CurlCurl ;
Case { {Type Gauss ;
Case { { GeoElement Triangle ; NumberOfPoints 4 ; }
{ GeoElement Quadrangle ; NumberOfPoints 4 ; }
{ GeoElement Tetrahedron ; NumberOfPoints 4 ; }
{ GeoElement Hexahedron ; NumberOfPoints 6 ; }
{ GeoElement Prism ; NumberOfPoints 9 ; } }
}
}
}
}