[Getdp] GetDP-Postprocessing Question

[Paul Fulmek]

Dear Christophe, dear GetDP-List,

first, thank you for your work! I managed to do some simple 
magnetostatics FEM-calculations very easily, great work.
I am sure I do not understand anything of most of the program and am 
looking forward to see some more examples to learn from.

My question regards the postprocessing abilities of GetDP:
I have calculated the solution of a (2D) scalar potential formulation 
and want to do some postprocessing.

1) How can I calculate the mean value of flux- and field-vector in some 
domain? I managed to calculate the area of my domain and the 
area-integral of field and flux; I don't know how to teach GetDP to 
print the simple quotient (integral field)/area.

2) I have defined a closed path around the magnet in the mesh, and want 
to calculate the force by doing a line-integral of Maxwells Stress 
Tensor on this path. I have no idea how to teach GetDP to do this 
calculation, how to find the normal-vector, and how to integrate!

By now I have solved both problems by postprocessing the GetDP 
output-files by hand or C-programs. It should be possible to have it 
done with GetDP!

Thank you very much,

sincerely

Paul

// simple-2d.pro
Group {
MSU = Region[ 100 ];  // Magnet-Surface
ISU = Region[ 101 ];  // Air-Surface
ASU = Region[ 102 ];  // Air->Infinity Surface
BLL = Region[ 110 ];  // Border-Line, Infinity!
MXP = Region[ 1001, 1002 ]; // Closed line around magnet MSU
ALL = Region[ { 100, 101, 102 }];
Val_Rint = 0.05;      // Radii for VolSphShell
Val_Rext = 0.08;
}

Function {
   mu0  = 4.e-7 * Pi;  // Permeability air
   mur  = 1.0;         // Relative permeability magnet
   hcx  = 160000:      // Coercivity magnet

   hc[ MSU ] = Vector[ hcx, 0., 0. ];

   mu[ Region[{ ASU, ISU }]] = mu0;
   mu[ MSU ]                 = mur * mu0;
}

Constraint {
   { Name MyConstraint;
   Case {
     { Region BLL; Type Assign; Value 0.; }}  // Infinite air borderline
   }
}

Jacobian {
   { Name MyJacobian;
   Case {
     { Region ASU ; Jacobian VolSphShell{ Val_Rint, Val_Rext }; }
     { Region All ; Jacobian Vol; } }
   }
}

Integration {
   { Name MyIntegration;
     Case {
       { Type Gauss ;
         Case {
           { GeoElement Line       ; NumberOfPoints 4; }
           { GeoElement Triangle   ; NumberOfPoints 4; }
           { GeoElement Quadrangle ; NumberOfPoints 4; } }
       }
     }
   }
}

FunctionSpace {
   { Name MagSta_phi; Type Form0;
   BasisFunction {
     { Name sn ; NameOfCoef phin ; Function BF_Node ;
       Support ALL ; Entity NodesOf[ All ] ; } }
   Constraint {
     { NameOfCoef phin ; EntityType NodesOf ;
       NameOfConstraint MyConstraint ; } }
   }
}

Formulation {
   { Name MyFormulation ; Type FemEquation ;
   Quantity {
     { Name phi ; Type Local ; NameOfSpace MagSta_phi ; } }
   Equation {
     Galerkin { [ - mu[] * Dof{d phi} , {d phi} ] ;
         In ALL; Jacobian MyJacobian; Integration MyIntegration;}
     Galerkin { [ - mu[] * hc[] , {d phi} ] ;
         In MSU; Jacobian MyJacobian; Integration MyIntegration; } }
   }
}

Resolution {
   { Name MyResolution ;
     System {
       { Name A ; NameOfFormulation MyFormulation ; } }
     Operation {
       Generate[A] ; Solve[A] ; SaveSolution[A] ; } }
}

PostProcessing {
   { Name MyPostProc ; NameOfFormulation MyFormulation ;
     Quantity {
       { Name b   ;
         Value {Local{[- mu[]*{d phi}]; In ALL; Jacobian MyJacobian;}
                Local{[- mu[]*hc[] ]  ; In MSU; Jacobian MyJacobian;}}}
       { Name h   ;
         Value{Local{[-{d phi}]       ; In ALL; Jacobian MyJacobian; }}}
/** |B|-mean **/
       { Name bm  ;
         Value { Integral { [Norm[ - mu[] * {d phi} ]] ; In ALL;
                   Integration MyIntegration; Jacobian MyJacobian; }
                 Integral { [Norm[ - mu[] * hc[] ]]    ; In MSU;
                   Integration MyIntegration; Jacobian MyJacobian; }}}
/** |H|-mean **/
       { Name hm  ;
         Value { Integral { [Norm[ - {d phi} ]]        ; In ALL;
                   Integration MyIntegration; Jacobian MyJacobian; }}}
     }
   }
}

PostOperation {

   { Name phi ; NameOfPostProcessing MyPostProc;
     Operation {
// that does not work ??:
//      Print[ b,   OnElementsOf MPX, File "b_mxp.pos", Depth 0 ] ;
//      Print[ h,   OnElementsOf MPX, File "h_mxp.pos", Depth 0 ] ;
// that does not work ??:
//      Print[ b[MPX],   OnGlobal, File "b_mxp.pos", Depth 0 ] ;
//      Print[ h[MPX],   OnGlobal, File "h_mxp.pos", Depth 0 ] ;
// this works, not very elegant:
       Print[ b, OnLine {{ 0.0115,-0.0040,0}{ 0.0115, 0.0040,0}}{50}, 

         File  "b_mxp.dat", Format Table ] ;
       Print[ b, OnLine {{ 0.0115, 0.0040,0}{-0.0115, 0.0040,0}}{100},
         File >"b_mxp.dat", Format Table ] ;
       Print[ b,   OnLine {{-0.0115, 0.0040,0}{-0.0115,-0.0040,0}}{50},
          File >"b_mxp.dat", Format Table ] ;
       Print[ b,   OnLine {{-0.0115,-0.0040,0}{ 0.0115,-0.0040,0}}{100},
          File >"b_mxp.dat", Format Table ] ;
// this works: for H-field
       Print[ h, OnLine {{ 0.0115,-0.0040,0}{ 0.0115, 0.0040,0}}{50}, 

         File  "h_mxp.dat", Format Table ] ;
       Print[ h, OnLine {{ 0.0115, 0.0040,0}{-0.0115, 0.0040,0}}{100},
         File >"h_mxp.dat", Format Table ] ;
       Print[ h,   OnLine {{-0.0115, 0.0040,0}{-0.0115,-0.0040,0}}{50},
         File >"h_mxp.dat", Format Table ] ;
       Print[ h,   OnLine {{-0.0115,-0.0040,0}{ 0.0115,-0.0040,0}}{100},
          File >"h_mxp.dat", Format Table ] ;
       Print[ bm[MSU],  OnGlobal ];  // integral of flux over magnet area
       Print[ hm[MSU],  OnGlobal ];
       Print[ mar[MSU], OnGlobal ];

     }
   }
}

**********************************************************************
Paul Fulmek
Technische Universitaet Wien
Austria Europe
Email: fulmek at tuwien.ac.at
Phone: +43-1-58801-36665
Fax:   +43-1-58801-36699
**********************************************************************