# [Getdp] Concatenation of field operators

Thomas_W e0208228 at student.tuwien.ac.at
Sun Oct 28 20:04:24 CET 2012

```Hi Michael,

I have tried your suggestion with the new quantity. I attach the
files. The result is that the computed result is just zero everywhere.

Are there some other ideas how to compute the term?  - Probably, there
is still a mistake in the files.

Thanks, in any event!
Thomas

Zitat von michael.asam at infineon.com:

> Hi Thomas,
>
> I think this cannot be done in the post-processing.
> Probably one solution would be to calculate the first gradient
> in the Formulation{} part by adding a new quantity (e.g. dv)
> and do a kind of protection:
>
> Formulation{
> ...
>   Galerkin { [ -Norm[{v}], {dv} ]; In All; Integration I1; Jacobian JVol;  }
>   Galerkin { [ Dof{dv}, {dv} ]; In All; Integration I1; Jacobian JVol;  }
> }
>
> You have to make sure that v has at least 2nd order basis functions
> in order to be twice differentiable.
>
> I hope this works :-)
>
> Best regards
> Michael
>
>
>
>
>
> -----Original Message-----
> From: getdp-bounces at ace20.montefiore.ulg.ac.be
> [mailto:getdp-bounces at ace20.montefiore.ulg.ac.be] On Behalf Of
> Thomas_W
> Sent: Monday, October 15, 2012 3:03 PM
> To: getdp at geuz.org
> Subject: [Getdp] Concatenation of field operators
>
> Dear all,
>
> I try to evaluate the following quantity in the Post-processing part,
> where v is the discretized quantity:
>
> I tried variants such as:
>
> PostProcessing {
>    { Name The; NameOfFormulation f_v;
>      Quantity {
>        { Name jj;
>          Value{
>          Local{ [ {d  {d v}*{d v} } ] ; In All;Jacobian JVol; }
>          }
>        }
>      }
>    }
> }
>
> or also
> Local{ [ {d  {Norm[{v}]} } ] ; In All;Jacobian JVol; }.
>
> All these variants yield an error, namely 'syntax error ({)' in the
> corresponding line. It probably has to do with the fact that only v,
> and not Norm[{v}] is the original discretized quantity.
>
> Do you have a solution how to evaluate such a function with GetDP?
>
> Thank you,
> Thomas
>
>
> _______________________________________________
> getdp mailing list
> getdp at geuz.org
> http://www.geuz.org/mailman/listinfo/getdp
>
>

-------------- next part --------------
// File "LaplacianNeumann.geo".

// We include the file containing the numbering of the geometry.
// This is usefull at the end of this file, and used to "synchronise" GMSH and GetDP
Include "param.geo";

//Caracteristic length of the finite elements (reffinement is also possible after the mesh is built):
lc = 0.05;
// This parameter could be placed for instance in "param.geo", to separate more easyly the geometry
// 	and the discretization parameters.

// The parameters of the border of the domain :
x_max = 1;
x_min = 0;
y_max= 1;
y_min = 0;

//Creation of the 4 angle points of the domain Omega (=square)
p1 = newp; Point(p1) = {x_min,y_min,0,lc};
p2 = newp; Point(p2) = {x_min,y_max,0,lc};
p3 = newp; Point(p3) = {x_max,y_max,0,lc};
p4 = newp; Point(p4) = {x_max,y_min,0,lc};
// Remarks:
// -"newp" is a GMSH function that give the first available number for describing a point.
// 	For any other entity, like Line, Surface, etc. We recommand the use of "newreg" (see below).
// - By default, GMSH create a 3D domain. The z-coordinate must always be precised.

//The four edges of the square
L1 = newreg; Line(L1) = {p1,p2};
L2 = newreg; Line(L2) = {p2,p3};
L3 = newreg; Line(L3) = {p3,p4};
L4 = newreg; Line(L4) = {p4,p1};

// Line Loop (= boundary of the square)
Bound = newreg; Line Loop(Bound) = {L1,L2,L3,L4};

//Surface of the square
SurfaceOmega = newreg; Plane Surface(SurfaceOmega) = {Bound};

// To conclude, we define the physical entities, that is "what GetDP could see/use".
// "Omega" is a number imported from the file "param.geo".
Physical Surface(Omega) = {SurfaceOmega};
// Do not forget to let a blank line at the end, this could make GMSH crash...
Physical Line(Left)={L1};
Physical Line(Top) ={L2};
Physical Line(Right)={L3};
Physical Line(Bottom)={L4};

-------------- next part --------------

Include "param.geo";

Group{
Omega = Region[{Omega}];
Top =#Top;
Bottom = #Bottom;
Left = #Left;
Right = #Right;
Boundlined = Region[{Left,Right,Top,Bottom}];
}

Function{
Coef = Pi;
f[] = (1+2*Coef*Coef)*Cos[Coef*X[]]*Cos[Coef*Y[]];
}

Jacobian {
{ Name JVol ;
Case {
{ Region All ; Jacobian Vol ; }
}
}
{ Name JSur ;
Case {
{ Region All ; Jacobian Sur ; }
}
}
{ Name JLin ;
Case {
{ Region All ; Jacobian Lin ; }
}
}
}

Integration {
{ Name I1 ;
Case {
{ Type Gauss ;
Case {
{ GeoElement Point ; NumberOfPoints  1 ; }
{ GeoElement Line ; NumberOfPoints  4 ; }
{ GeoElement Triangle ; NumberOfPoints  6 ; }
{ GeoElement Quadrangle ; NumberOfPoints 7 ; }
{ GeoElement Tetrahedron ; NumberOfPoints 15 ; }
{ GeoElement Hexahedron ; NumberOfPoints 34 ; }
}
}
}
}
}

Constraint {
{ Name Sta_V ; Type Assign;
Case {
{ Region Bottom; Value 0; }
{ Region Top; Value 1;}
}
Name Sta_V ; Type Assign;
Case {
{ Region Boundlined; Value 0; }
}
}
}

FunctionSpace{ // 2nd order basis function
{ Name Vh; Type Form0;
BasisFunction{
{Name wn; NameOfCoef vn; Function BF_Node;    Support Omega; Entity NodesOf[All];}
{Name wn2; NameOfCoef vn2; Function BF_Node_2E;    Support Omega; Entity EdgesOf[All];}
}

Constraint {
{ NameOfCoef vn2; EntityType EdgesOf ; NameOfConstraint Sta_V; }
}
}

{ Name Vh2; Type Form0;
BasisFunction{
{Name wn; NameOfCoef vnr; Function BF_Node;    Support Omega; Entity NodesOf[All];}
}
}
}

Formulation{
{Name Laplacian; Type FemEquation;
Quantity{
{Name u; Type Local; NameOfSpace Vh;}
{Name du; Type Local; NameOfSpace Vh2;} // this is for the evaluation of the change of current density
}
Equation{
In Omega; Jacobian JVol; Integration I1;}

Galerkin{ [  Dof{u}, {u}];
In Omega; Jacobian JVol; Integration I1;}

Galerkin{ [-f[], {u}];
In Omega; Jacobian JVol; Integration I1;}

Galerkin{ [  Dof{du}, {du}];
In Omega; Jacobian JVol; Integration I1;}

In Omega; Jacobian JVol; Integration I1;}

}
}
}

Resolution{
{Name Laplacian;
System{
{Name Syst; NameOfFormulation Laplacian;}
}
Operation{
Generate[Syst]; Solve[Syst]; SaveSolution[Syst];
}
}
}

PostProcessing{
{Name Laplacian; NameOfFormulation Laplacian;
Quantity{
{Name u; Value {Local{[{u}];In Omega;Jacobian JVol;}}}
// this gives the current density vector grad u
{Name CDe1; Value {Local{[{Grad u}]; In Omega; Jacobian JVol;}}}
// this does not give the absolute value |grad u| of the current density vector, but zero
{Name CDe2 ; Value { Term { [{du}]   ; In Omega ; Jacobian JVol;} } }
{Name CDe3 ; Value { Term { [Norm[{d du}]]   ; In Omega ; Jacobian JVol;} } }

}
}
}

PostOperation{
{Name Map_u; NameOfPostProcessing Laplacian;
Operation{
Print[u, OnElementsOf Omega, File "u.pos"];
Print[CDe1, OnElementsOf Omega, File "CDe1.pos"];
Print[CDe2, OnElementsOf Omega, File "CDe2.pos"];
Print[CDe3, OnElementsOf Omega, File "CDe3.pos"];
}
}
}

-------------- next part --------------
// File "param.geo"

//Numbers that caracterise the interior of the square (Omega) and its boundary (Gama):
Omega = 1000;
Top = 111;
Left = 112;
Bottom = 113;
Right = 114;

// Three remarks on these numbers :
// - They are arbitrary choosen.
// - They are placed in a separated file to be readable by both GMSH and GetDP.
// - "Gamma" is a special word used by GMSH/GetDP, that is why the boundary is named "Gama", with one "m"...
// Do not forget to let a blank line at the end, this could make GMSH crash...
```