new tuto: electrostatics with global quantities

Electrostatics/microstrip.pro

@@ -14,8 +14,6 @@
    To compute the solution interactively from the Gmsh GUI:
        File > Open > microstrip.pro
        Run (button at the bottom of the left panel)
-
-   The microstrip.pro file describes the finite element (FE) model.
    ------------------------------------------------------------------- */
 
 Group { 
@@ -90,7 +88,7 @@ FunctionSpace {
      - a set of scalar basis functions ("BF_Node" means nodal basis functions)
      - a constraint (here the Dirichlet boundary conditions) 
 
-     The FE representation of the unknown field "v" reads
+     The FE expansion of the unknown field "v" reads
      
      v = Sum_k vn_k sn_k 
 

ElectrostaticsFloating/floating.geo

+/* -------------------------------------------------------------------
+   File "microstrip.geo"
+
+   This file is the geometrical description used by GMSH to produce
+   the file "microstrip.msh".
+   ------------------------------------------------------------------- */
+
+/* Definition of some parameters for geometrical dimensions, i.e.
+   h (height of 'Diel1'), w (width of 'Line'), t (thickness of 'Line')
+   xBox (width of the air box) and yBox (height of the air box) */
+
+h = 1.e-3 ; w = 4.72e-3 ;  t = 0.035e-3 ;
+xBox = w/2. * 6. ;  yBox = h * 12. ;
+
+/* Definition of parameters for local mesh dimensions */
+
+s = 1. ;
+p0 = h / 10. * s ;
+pLine0 = w/2. / 10. * s ;  pLine1 = w/2. / 50. * s ;
+pxBox = xBox / 10. * s ;  pyBox = yBox / 8. * s ;
+
+/* Definition of gemetrical points */
+
+Point(1) = { 0   , 0, 0, p0} ;
+Point(2) = { xBox, 0, 0, pxBox} ;
+Point(3) = { xBox, h, 0, pxBox} ;
+Point(4) = { 0   , h, 0, pLine0} ;
+Point(5) = { w/2., h, 0, pLine1} ;
+Point(6) = { 0   , h+t, 0, pLine0} ;
+Point(7) = { w/2., h+t, 0, pLine1} ;
+Point(8) = { 0   , yBox, 0, pyBox} ;
+Point(9) = { xBox, yBox, 0, pyBox} ;
+
+/* Definition of gemetrical lines */
+
+Line(1) = {1,2};   Line(2) = {2,3};  Line(3) = {3,9};
+Line(4) = {9,8};   Line(5) = {8,6};  Line(7) = {4,1};
+Line(8) = {5,3};   Line(9) = {4,5};  Line(10) = {6,7};
+Line(11) = {5,7};
+
+/* Definition of geometrical surfaces */
+
+Line Loop(12) = {1, 2, -8, -9, 7};   Plane Surface(13) = {12};
+Line Loop(14) = {10,-11,8,3,4,5}; Plane Surface(15) = {14};
+
+/* Definition of Physical entities (surfaces, lines). 
+   The definition of Physical entities (Surfaces and Lines) 
+   tells GMSH the elements and associated region numbers
+   that have to be saved in the mesh file 'microstrip.msh'. 
+   For example, Region 111 is made of the triangle elements of the geometric surface 13, 
+   whereas Region 121 is made of line elements of the geometric lines 9, 10 and 11. */
+
+Physical Surface ("Air", 101) = {15};
+Physical Surface ("Dielectric", 111) = {13};
+
+Physical Line ("Ground", 120) = {1} ;
+Physical Line ("Electrode", 121) = {9,10,11} ;
+Physical Line ("Surface infinity", 130) = {2,3,4} ;

ElectrostaticsFloating/floating.pro

+/* -------------------------------------------------------------------
+   Tutorial 4 : floating potential of a microstrip electrode
+
+   Features:
+   - Global quantities and their special shape functions
+   - Duality 
+   - More on ONELAB parameters (flags, model options, check boxes, menus, ...)
+
+   To compute the solution interactively from the Gmsh GUI:
+       File > Open > floating.pro
+       Run (button at the bottom of the left panel)
+
+   ------------------------------------------------------------------- */
+
+/*
+  A Thing GetDP is pretty good at is the management of global (non-local) basis functions.
+  Finite element expansions typically associate basis functions to individual nodes 
+  or edges in the mesh. But consider the situation where a scalar field is set to be uniform 
+  over a region of the problem (a floating potential electrode in an Electrostatics problem,
+  to fix the idea).
+  By factorizing the identical nodal value "v_electrode", 
+  a global (non-local) basis function "BF_electrode" is obtained as factor 
+  which is the sum of the shape functions of all the nodes in the electrode region.
+  BF_electrode 
+  - is a continuous function
+  - is equal to 1 at the nodes of the electrode region, and to 0 at all other nodes
+  - decreases from 1 to 0 over the one element thick layer of outside finite elements
+    immediately in contact with the electrode region
+  One such glabal basis function can be associated with each electrode in the system,
+  so that the finite element expansion of the electric scalar potential reads:
+
+   v = Sum_k sn_k vn_k + Sum_electrode v_electrode BF_electrode
+
+   We show in this tutorial how GetDP takes advantage of global quantities
+   and the associated global basis functions
+   - to reduce the number of unknowns
+   - to deal with floating potentials
+   - to compute efficiently the electrode charges "Q_electrode", 
+     which are precisely the energy duals of the global "v_electrode" quantities
+   - to provide output quantities (charges, armature voltages, capacitances, ...)
+     that can be immediately used in a external circuit.
+ */
+
+
+Group { 
+  /* Geometrical regions: */ 
+
+  Air = Region[101]; 
+  Diel1 = Region[111];
+
+  Ground = Region[120];
+  Microstrip = Region[121];
+  SurfInf = Region[130];
+
+  /* Abstract regions:
+     Vol_Dielectric_Ele : dielectric volume regions where 
+                          "Div ( epsr[] Grad v)" is solved
+     Sur_Neu_Ele        : Neumann bondary condition ( epsr[] n.Grad v = 0 )
+     Electrodes_Ele     : electrode regions
+     
+
+     No prefix (Vol_ or Sur_) for the region "Electrodes_Ele", 
+     which can contain both surface or volume regions.
+     There are two electrodes in this model: Ground and Microstrip
+  */
+ 
+  Vol_Dielectric_Ele = Region[ {Air, Diel1} ];
+  Sur_Neu_Ele = Region[ {SurfInf} ];
+  Electrodes_Ele = Region [ {Ground, Microstrip} ]; 
+}
+
+/* A number of ONELAB parameters are defined to define model parameters 
+   or model options interactively. */
+MicrostripTypeBC = 
+  DefineNumber[0, Name "1Microstrip excitation/Type", 
+	       Choices{ 0="Fixed voltage", 1="Fixed charge"} ] ;
+MicrostripValueBC = 
+  DefineNumber[1e-3, Name "1Microstrip excitation/Value"] ;
+EpsilonRelDiel = 
+  DefineNumber[9.8, Name "2Dielectric/Relative permittivity"] ;
+DisplayGlobalBF = 
+  DefineNumber[0, Name "3Options/Display global basis functions", Choices {0,1} ] ;
+OverwriteOutput = 
+  DefineNumber[1, Name "3Options/Overwrite output.txt file", Choices {0,1} ] ;
+
+Function {
+  epsr[Air] = 1.;
+  epsr[Diel1] = EpsilonRelDiel;
+  eps0 = 8.854187818e-12;  // permittivity of empty space
+}
+
+
+Constraint {
+  /* Dirichlet boundary condition is no longer used. 
+     The microstrip and the ground are now treated as electrodes   */
+    { Name Dirichlet_Ele; Type Assign;
+      Case {}
+    }
+
+  { Name SetGlobalPotential; Type Assign;
+    Case {
+      /* Define the imposed potential regionwise on the different parts of "Electrodes_Ele".
+	 No voltage imposed to the Microstrip electrode
+	 when the "Fixed charge" option is enabled ( MicrostripTypeBC = true ). */
+      { Region Ground; Value 0; }
+      If( !MicrostripTypeBC )
+	{ Region Microstrip; Value MicrostripValueBC; }
+      EndIf
+    }
+  }
+  { Name SetArmatureCharge; Type Assign;
+    Case {
+      If( MicrostripTypeBC )
+	{ Region Microstrip; Value MicrostripValueBC; }
+      EndIf
+    }
+  }
+}
+
+Group{
+  /* The domain of definition of lists all regions 
+     on which the field "v" is defined.*/
+  Dom_Hgrad_v_Ele =  Region[ {Vol_Dielectric_Ele,
+			      Sur_Neu_Ele,
+			      Electrodes_Ele } ];
+}
+
+FunctionSpace {
+  /* The magic in the treatment of global quantitities by GetDP is in the fact 
+     that nearly all is done at the level of the FunctionSpace definition.
+
+     The finite element expansion is 
+
+     v = Sum_k sn_k vn_k + Sum_electrode v_electrode BF_electrode
+
+     The sum_k runs over all nodes except those of the electrode regions.
+     "sf" stands for "BF_electrode"
+     "vf" stands for "v_electrode"
+  */
+
+  { Name Hgrad_v_Ele; Type Form0;
+    BasisFunction {
+      { Name sn; NameOfCoef vn; Function BF_Node;
+        Support Dom_Hgrad_v_Ele; Entity NodesOf[ All, Not Electrodes_Ele ]; }
+      { Name sf; NameOfCoef vf; Function BF_GroupOfNodes; 
+        Support Dom_Hgrad_v_Ele; Entity GroupsOfNodesOf[ Electrodes_Ele ]; }
+    }
+    GlobalQuantity {
+      { Name GlobalPotential; Type AliasOf       ; NameOfCoef vf; }
+      { Name ArmatureCharge ; Type AssociatedWith; NameOfCoef vf; }
+    }
+    Constraint {
+      { NameOfCoef vn; EntityType NodesOf; 
+        NameOfConstraint Dirichlet_Ele; } // unused in this model
+      { NameOfCoef GlobalPotential; EntityType GroupsOfNodesOf;
+	NameOfConstraint SetGlobalPotential; }
+      { NameOfCoef ArmatureCharge; EntityType GroupsOfNodesOf; 
+	NameOfConstraint SetArmatureCharge; }
+    }
+    // Subspace definition only needed to display "vf" in PostProcessing
+    SubSpace { 
+      { Name vf; NameOfBasisFunction sf; }
+    }
+  }
+}
+
+Jacobian {
+  { Name Vol ;
+    Case { 
+      { Region All ; Jacobian Vol ; }
+    }
+  }
+}
+
+Integration {
+  { Name Int ;
+    Case { {Type Gauss ;
+            Case { { GeoElement Triangle    ; NumberOfPoints  4 ; }
+                   { GeoElement Quadrangle  ; NumberOfPoints  4 ; } }
+      }
+    }
+  }
+}
+
+Formulation {
+  /* Minor changes in the formulation.
+     Global quantities are declared in the "Quantity{}" section.
+     The Global term triggers the creation of the addition equation
+     in the system that computes the charge Q_electrode carried by each electrode
+
+     Q_electrode = (-epsr[] Grad v, Grad BF_electrode)_Vol_Dielectric_Ele
+   */
+  { Name Electrostatics_v; Type FemEquation;
+    Quantity {
+      { Name v; Type Local; NameOfSpace Hgrad_v_Ele; }
+      { Name U; Type Global; NameOfSpace Hgrad_v_Ele [GlobalPotential]; }
+      { Name Q; Type Global; NameOfSpace Hgrad_v_Ele [ArmatureCharge]; }
+      // next line only needed to display the global BF in PostProcessing
+      { Name vf; Type Local; NameOfSpace Hgrad_v_Ele [vf]; }
+    }
+    Equation {
+      Galerkin { [ epsr[] * Dof{d v} , {d v} ]; In Vol_Dielectric_Ele; 
+	Jacobian Vol; Integration Int; }
+      GlobalTerm { [ -Dof{Q}/eps0 , {U} ]; In Electrodes_Ele; }
+    }
+  }
+}
+
+Resolution {
+  { Name EleSta_v;
+    System {
+      { Name Sys_Ele; NameOfFormulation Electrostatics_v; }
+    }
+    Operation { 
+      Generate[Sys_Ele]; Solve[Sys_Ele]; SaveSolution[Sys_Ele];
+      If( OverwriteOutput )
+	DeleteFile[ "output.txt" ];
+      EndIf
+    }
+  }
+}
+
+PostProcessing {
+ { Name EleSta_v; NameOfFormulation Electrostatics_v;
+    Quantity {
+      { Name v; Value { Term { [ {v} ]; In Dom_Hgrad_v_Ele; Jacobian Vol; } } }
+      { Name e; Value { Term { [ -{d v} ]; In Dom_Hgrad_v_Ele; Jacobian Vol; } } }
+      { Name d; Value { Term { [ -eps0*epsr[] * {d v} ]; 
+	    In Dom_Hgrad_v_Ele; Jacobian Vol; } } }
+      { Name Q; Value { Term { [ {Q} ]; In Electrodes_Ele; } } }
+      { Name U; Value { Term { [ {U} ]; In Electrodes_Ele; } } }
+      { Name C; Value { Term { [ {Q}/{U} ]; 
+	    In Electrodes_Ele; } } }
+      { Name energy;
+	Value { Integral { Type Global;
+	    [ eps0*epsr[] / 2. * SquNorm[{d v}] ];
+	    In Vol_Dielectric_Ele; Jacobian Vol; Integration Int;
+	  }
+	}
+      }
+      // next lines only needed to display global basis functions in PostProcessing
+      { Name BFGround; Value { Term { [ BF {vf} ]; In Dom_Hgrad_v_Ele; 
+	    SubRegion Ground; Jacobian Vol; } } }
+      { Name BFMicrostrip; Value { Term { [ BF {vf} ]; In Dom_Hgrad_v_Ele; 
+	    SubRegion Microstrip; Jacobian Vol; } } }
+
+    }
+  }
+}
+
+/* Various output results associated with the global quantities are generated.
+   They are both displayed in the graphical user interface, and stored in disk files.
+   In particular, all global quantities related results 
+   are stored in the "output.txt" file. 
+   There is a user option to display the global Basis functions of the two electrodes
+   in the system. Another user option allow the user to chose to not overwrite
+   the "output.txt" file when running a new simulation. */
+
+PostOperation {
+  { Name Map; NameOfPostProcessing EleSta_v;
+     Operation {
+       If( DisplayGlobalBF )
+	 Print[ BFGround, OnElementsOf Dom_Hgrad_v_Ele, File "BFGround.pos" ];
+         Echo[ StrCat["l=PostProcessing.NbViews-1;", 
+		      "View[l].IntervalsType = 1;",
+		      "View[l].NbIso = 40;",
+		      "View[l].ShowElement = 1;"],
+	       File "BFGround.opt", LastTimeStepOnly] ;
+
+	 Print[ BFMicrostrip, OnElementsOf Dom_Hgrad_v_Ele, File "BFMicrostrip.pos" ];
+	 Echo[ StrCat["l=PostProcessing.NbViews-1;", 
+		      "View[l].IntervalsType = 1;",
+		      "View[l].NbIso = 40;",
+		      "View[l].ShowElement = 1;"],
+	       File "BFMicrostrip.opt", LastTimeStepOnly] ;
+       EndIf
+
+       Print[ v, OnElementsOf Dom_Hgrad_v_Ele, File "v.pos" ];
+       Echo[ StrCat["l=PostProcessing.NbViews-1;", 
+                    "View[l].IntervalsType = 3;",
+                    "View[l].NbIso = 40;"],
+             File "mStrip_v.geo", LastTimeStepOnly] ;
+
+       Echo[ "Microstrip charge [C]:", Format Table, File > "output.txt"] ;
+       Print[ Q, OnRegion Microstrip, File > "output.txt", Color "AliceBlue",
+	      Format Table, SendToServer "Output/Microstrip/Charge [C]" ];
+       Echo[ "Microstrip potential [V]:", Format Table, File > "output.txt"] ;
+       Print[ U, OnRegion Microstrip, File > "output.txt", Color "AliceBlue",
+	      Format Table, SendToServer "Output/Microstrip/Potential [V]" ];
+       Echo[ "Ground charge [C]:", Format Table, File > "output.txt"] ;
+       Print[ Q, OnRegion Ground, File > "output.txt", Color "AliceBlue",
+	      Format Table, SendToServer "Output/Ground/Charge [C]" ];
+       Echo[ "Microstrip capacitance [F]:", Format Table, File > "output.txt"] ;
+       Print[ C, OnRegion Microstrip, File > "output.txt", Color "AliceBlue",
+	      Format Table, SendToServer "Output/Global/Capacitance [F]" ];
+       Echo[ "Electrostatic energy [J]:", Format Table, File > "output.txt"] ;
+       Print[ energy, OnRegion Vol_Dielectric_Ele, File > "output.txt", Color "AliceBlue",
+	      Format Table, SendToServer "Output/Global/Energy [J]" ];
+     }
+  }
+}
+
+