# [Getdp] OneTooth Example with non linear material

Matthias Lang lang_matthias at gmx.de
Tue Oct 22 23:55:37 CEST 2013

```Hello ererybody,

I took the 3D magnetostatic example OneTooth.pro
and included nonlinear material. (found in Material.pro)
But unfortunatly I get the following Error Message:

.
.
Info    : GetDP - E n d   P r e - P r o c e s s i n g
Info    : GetDP - P r o c e s s i n g . . .
Info    : GetDP - IterativeLoop ...
Info    : GetDP - Non Linear Iteration Relaxation 1
Info    : GetDP - GenerateJac[Sys_Mag]
Error   : GetDP - Function _Exp_Fct_17 called with too few arguments
Info    : Done running 'GetDP'

Where is the problem? - inside or in front of the computer?

Regards

Matthias Lang
-------------- next part --------------
/* -------------------------------------------------------------------
File "MagSta_a_3D_nl.pro"
Magnetostatics - Magnetic vector potential a formulation (3D)
-------------------------------------------------------------------
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 :
----------
*/
Group {
DefineGroup[ Domain_Mag, DomainS_Mag, Domain_NL];
}

Function {
DefineFunction[ nu, dhdb_NL ];
DefineConstant[
Nb_max_iter        = 30,
relaxation_factor  = 1,
stop_criterion     = 1e-5,
NbT        	   = 1,
Steel_1020_dhdb_NL
];
}

FunctionSpace {
// Magnetic vector potential a (b = curl a)
{ Name Hcurl_a_Mag_3D; Type Form1;
BasisFunction {
// a = a s
//  e   e
{ Name se; NameOfCoef ae; Function BF_Edge;
Support Domain_Mag; Entity EdgesOf[ All ]; }
}
Constraint {
{ NameOfCoef ae; EntityType EdgesOf;
NameOfConstraint MagneticVectorPotential_3D; }
{ NameOfCoef ae; EntityType EdgesOfTreeIn; EntitySubType StartingOn;
NameOfConstraint Gauge; }
}
}

}

Formulation {
{ Name Magnetostatics_a_3D; Type FemEquation;
Quantity {
{ Name a ; Type Local; NameOfSpace Hcurl_a_Mag_3D; }
}
Equation {
Galerkin { [ nu[] * Dof{Curl a} , {Curl a} ]; In Domain_Mag;
Jacobian Vol; Integration CurlCurl; }
Galerkin { JacNL [ dhdb_NL[{d a}] * Dof{d a} , {d a} ] ;
In Domain_NL ; Jacobian Vol ; Integration CurlCurl ; }

Galerkin { [ - js[] , {a} ]; In DomainS_Mag;
Jacobian Vol; Integration CurlCurl; }

}
}
}

Resolution {
{   Name MagSta_a_3D;
System {
{ Name Sys_Mag; NameOfFormulation Magnetostatics_a_3D; }
}
Operation {
//             Generate[Sys_Mag]; Solve[Sys_Mag]; //linear case
IterativeLoop[Nb_max_iter, stop_criterion, relaxation_factor]{
GenerateJac[Sys_Mag] ; SolveJac[Sys_Mag] ; }
SaveSolution[Sys_Mag];
}
}
}

PostProcessing {
{   Name MagSta_a_3D; NameOfFormulation Magnetostatics_a_3D;
Quantity {
{ Name a;
Value {
Local { [ {a} ]; In Domain_Mag; Jacobian Vol; }
}
}
{ Name a_norm ;
Value {
Local { [ Norm[{a}] ] ; In Domain_Mag ; Jacobian Vol ; }
}
}
{ Name b;
Value {
Local { [ {d a} ]; In Domain_Mag; Jacobian Vol; }
}
}
{ Name h;
Value {
Local { [ nu[] * {d a} ]; In Domain_Mag; Jacobian Vol; }
}
}
{ Name js;
Value {
Local { [js[]] ; In DomainS_Mag; Jacobian Vol; }
}
}

}
}
}
-------------- next part --------------
Function{
/*--------------------------------------------------------------------------------------*/
// nu = 100. + 10. * exp ( 1.8 * b * b )
// analytical
nu_1a[] = 100. + 10. * Exp[1.8*SquNorm[\$1]] ;
dnudb2_1a[] = 18. * Exp[1.8*SquNorm[\$1]] ;
h_1a[] = nu_1a[\$1]*\$1 ;
dhdb_1a[] = TensorDiag[1,1,1] * nu_1a[\$1#1] + 2*dnudb2_1a[#1] * SquDyadicProduct[#1]  ;
dhdb_1a_NL[] = 2*dnudb2_1a[\$1#1] * SquDyadicProduct[#1]  ;
/*--------------------------------------------------------------------------------------*/

/*--------------------------------------------------------------------------------------*/

/*--------------------------------------------------------------------------------------*/

// BH data EIcore

MatEIcore_h = {0.,
2.1827e+02,3.8272e+02,5.2036e+02,7.1167e+02,
8.4921e+02,1.0405e+03,1.4741e+03,2.0409e+03,
3.0109e+03,5.0572e+03,7.5335e+03,1.0037e+04,
1.2486e+04,1.5015e+04,1.7464e+04,2.0021e+04,
2.2040e+04,2.5000e+04} ;

MatEIcore_b = {0.,
2.0329e-01,4.0287e-01,6.0986e-01,8.0575e-01,
1.0053e+00,1.1975e+00,1.4008e+00,1.5154e+00,
1.5967e+00,1.6706e+00,1.7076e+00,1.7335e+00,
1.7446e+00,1.7593e+00,1.7667e+00,1.7815e+00,
1.7889e+00,1.7963e+00} ;

MatEIcore_b2 = List[MatEIcore_b]^2 ;
MatEIcore_nu = List[MatEIcore_h]/List[MatEIcore_b] ;
MatEIcore_nu(0) = MatEIcore_nu(1);

MatEIcore_nu_b2  = ListAlt[MatEIcore_b2, MatEIcore_nu] ;
//nu_EIcore[] = InterpolationLinear[SquNorm[\$1]]{List[MatEIcore_nu_b2]} ; // used later in *.pro file
MatEIcore_interp_nu[] = InterpolationLinear[SquNorm[\$1]]{List[MatEIcore_nu_b2]} ; // used later in *.pro file
//dnudb2_EIcore[] = dInterpolationLinear[SquNorm[\$1]]{List[MatEIcore_nu_b2]} ;
MatEIcore_interp_dnudb2[] = dInterpolationLinear[SquNorm[\$1]]{List[MatEIcore_nu_b2]} ;
//h_EIcore[] = nu_EIcore[\$1] * \$1 ;
MatEIcore_interp_h[] = MatEIcore_interp_nu[\$1] * \$1 ;
//dhdb_EIcore[] = TensorDiag[1,1,1]*nu_EIcore[\$1#1] + 2*dnudb2_EIcore[#1] * SquDyadicProduct[#1] ;
MatEIcore_dhdb[] = TensorDiag[1,1,1]*MatEIcore_interp_nu[\$1#1] + 2*MatEIcore_interp_dnudb2[#1] * SquDyadicProduct[#1] ;
//dhdb_EIcore_NL[] = 2*dnudb2_EIcore[\$1] * SquDyadicProduct[\$1] ; // used later in *.pro file
MatEIcore_dhdb_NL[] = 2*MatEIcore_interp_dnudb2[\$1] * SquDyadicProduct[\$1] ; // used later in *.pro file
/*--------------------------------------------------------------------------------------*/
/*--------------------------------------------------------------------------------------*/

// BH data for: 1020 Steel

Steel_1020_h = {0.000000,
79.577472,    100.182101,    126.121793,    158.777930,
199.889571,    251.646061,    316.803620,    398.832128,
502.099901,    632.106325,    795.774715,   1001.821011,
1261.217929,   1587.779301,   1998.895710,   2516.460605,
3168.036204,   3988.321282,   5020.999013,   6321.063250,
7957.747155,  10018.210114,  12612.179293,  15877.793010,
19988.957103,  25164.606052,  31680.362037,  39883.212823,
50209.990127,  63210.632497,  79577.471546, 100182.101136,
126121.792926, 158777.930096, 199889.571030, 251646.060522,
316803.620370};

Steel_1020_b = {0.000000,
0.166252, 0.208725, 0.261635, 0.327146,
0.407495, 0.504604, 0.619366, 0.750502,
0.893195, 1.038258, 1.173276, 1.286731,
1.373533, 1.437658, 1.489189, 1.538344,
1.591321, 1.649723, 1.711927, 1.774712,
1.834976, 1.891442, 1.945262, 1.998705,
2.052852, 2.106184, 2.155149, 2.196375,
2.229008, 2.255450, 2.280053, 2.307136,
2.339975, 2.381044, 2.432708, 2.497748,
2.579627};

Steel_1020_b2 = List[Steel_1020_b]^2 ;
Steel_1020_nu = List[Steel_1020_h]/List[Steel_1020_b] ;
Steel_1020_nu(0) = Steel_1020_nu(1);

Steel_1020_nu_b2  = ListAlt[Steel_1020_b2, Steel_1020_nu] ;
//nu_EIcore[] = InterpolationLinear[SquNorm[\$1]]{List[Steel_1020_nu_b2]} ; // used later in *.pro file
Steel_1020_interp_nu[] = InterpolationLinear[SquNorm[\$1]]{List[Steel_1020_nu_b2]} ; // used later in *.pro file
//dnudb2_EIcore[] = dInterpolationLinear[SquNorm[\$1]]{List[Steel_1020_nu_b2]} ;
Steel_1020_interp_dnudb2[] = dInterpolationLinear[SquNorm[\$1]]{List[Steel_1020_nu_b2]} ;
//h_EIcore[] = nu_EIcore[\$1] * \$1 ;
Steel_1020_interp_h[] = Steel_1020_interp_nu[\$1] * \$1 ;
//dhdb_EIcore[] = TensorDiag[1,1,1]*nu_EIcore[\$1#1] + 2*dnudb2_EIcore[#1] * SquDyadicProduct[#1] ;
Steel_1020_dhdb[] = TensorDiag[1,1,1]*Steel_1020_interp_nu[\$1#1] + 2*Steel_1020_interp_dnudb2[#1] * SquDyadicProduct[#1] ;
//dhdb_EIcore_NL[] = 2*dnudb2_EIcore[\$1] * SquDyadicProduct[\$1] ; // used later in *.pro file
Steel_1020_dhdb_NL[] = 2*Steel_1020_interp_dnudb2[\$1] * SquDyadicProduct[\$1] ; // used later in *.pro file
/*--------------------------------------------------------------------------------------*/
}

-------------- next part --------------
A non-text attachment was scrubbed...
Name: OneTooth.geo
Type: application/vnd.dynageo
Size: 7369 bytes
Desc: not available
URL: <http://www.geuz.org/pipermail/getdp/attachments/20131022/186a4a61/attachment.geo>
-------------- next part --------------
// Fichier "OneTooth.pro"
// Magnétostatique

Include "OneTooth-dat.pro";
Include "Materials.pro"; // nonlinear BH caracteristic of magnetic material

Group {
Air          = #AIR;
Core         = #CORE;
Ind          = #COIL;
Plate        = #PLATE;

SurfaceGe0   = #CUT_SYMMETRY ;
SurfaceGInf  = #AIR_INF;
SkinCore     = #SKIN_CORE;
SkinPlate    = #SKIN_PLATE;

DomainS_Mag = Region[{Ind}];
Domain_Mag  = Region[{Core, Air, Ind, Plate}];
Domain_NL  = Region[{Core, Plate}];
}

Function {
mu0 = 4.e-7 * Pi;
murSteel = 100.;
nu [ Region[{Air, Ind}] ] = 1. / mu0;
//linear Case
//    nu [ Region[{Core}]] = 1. / (murSteel*mu0);
//    nu [ Region[{Plate}]] = 1. / (murSteel*mu0);
//Nonlinear Case
//	nu [ Domain_NL ] = MatEIcore_interp_nu[\$1] ;
//	dhdb_NL [ Domain_NL ] = MatEIcore_dhdb_NL[\$1];
nu [ Domain_NL ] = Steel_1020_interp_nu[\$1] ;
dhdb_NL [ Domain_NL ] = Steel_1020_dhdb_NL[\$1];

Sc          = th_coil*ly_coil;  //cross section of coil
Ns          = 100;		//Number of turns of the inductor Ind
Val_Current = 1.0;		//current in the inductor

//calculation of current vectors
xc = 0.;  // center of the coil
zc = 0.;   // center of the coil

js[ Ind ] = Ns * Val_Current * Vector[-Sin[Atan2[Z[]-zc,X[]-xc]],
0.,
Cos[Atan2[Z[]-zc,X[]-xc]]] ;
}

Constraint {
{ Name MagneticVectorPotential_3D; Type Assign;
Case {
{ Region SurfaceGInf ; Value 0.; }
{ Region SurfaceGe0 ; Value 0.; }
}
}

{ Name Gauge; Type Assign;
Case {
{ Region Domain_Mag; SubRegion #{SurfaceGInf, SurfaceGe0,
SkinPlate, SkinCore}; Value 0.; }
}

}
}
Jacobian {
{ Name Vol ;
Case {
{ Region All ; Jacobian Vol ; }
}
}
}

Integration {
{ Name CurlCurl ;
Case { {Type Gauss ;
Case {
{ GeoElement Triangle ; NumberOfPoints 4 ; }
{ GeoElement Tetrahedron ; NumberOfPoints 4 ; }
}
}
}
}
}

Include "Magsta_a_3D_nl.pro"

PostOperation {
{ Name Map_a; NameOfPostProcessing MagSta_a_3D;
Operation {
Print[ a, OnElementsOf Domain_Mag, File "Magsta_a_3D_a.pos" ];
Print[ a_norm, OnElementsOf Domain_Mag, File "Magsta_a_3D_a_norm.pos" ];
Print[ b,  OnElementsOf Domain_Mag, File "Magsta_a_3D_b.pos" ];
Print[ js, OnElementsOf DomainS_Mag, File "Magsta_a_3D_js.pos" ];
Print[ js, OnElementsOf DomainS_Mag, Format Table,File  "js.txt" ];
}
}
}

// EOF
-------------- next part --------------
// Fichier "OneTooth-dat.pro"

mm = 1.e-3; // mm -> m

/*
x -> longueur		-- gauche, droite
y -> hauteur		-- haut, bas
z -> ?paisseur	-- avant, arri?re
*/

// core
lx_core = 30. * mm;
ly_core = 100. * mm;
lz_core = 30. * mm;

// coil
th_coil = 10. * mm;
ri_coil = (lx_core^2. + lz_core^2.)^0.5 / 2. + 1. * mm;
re_coil = ri_coil + th_coil;
ly_coil = 50. * mm;

// plate
airgap = 2. * mm;
th_plate = 20. * mm;
ri_plate = ((lx_core/2.)^2. + (ly_core/2.)^2.)^0.5 + airgap;
re_plate = ri_plate+th_plate;
lz_plate = lz_core + 10. * mm;

// la bo?te ? air
bordure = 20. * mm;
axl = -re_plate - bordure; //left of airbox
axr =  re_plate + bordure; //right of airbox
ayb = -re_plate - bordure; //bottom of airbox
ayt =  re_plate + bordure; //top of airbox
azb =  0.; //back of airbox
azf =  re_coil +  bordure; //front of airbox

// les entit?s physiques
AIR    = 9001;
CORE   = 9002;
COIL   = 9003;
PLATE  = 9004;

AIR_INF      = 9011;
SKIN_PLATE   = 9012;
SKIN_CORE    = 9013;
CUT_SYMMETRY = 9014;

// EOF
```