Skip to content
Snippets Groups Projects
Commit 6b2ff6d4 authored by François Henrotte's avatar François Henrotte
Browse files

added explanation in .pro file

parent 62dca156
No related branches found
No related tags found
No related merge requests found
/* -------------------------------------------------------------------
File "microstrip.pro"
Tutorial 1 : electrostatic field of a microstrip
This file defines the problem dependent data structures for the
microstrip problem.
To compute the solution:
To compute the solution in a terminal:
getdp microstrip -solve EleSta_v
To compute post-processing results:
getdp microstrip -pos Map
getdp microstrip -pos Cut
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 {
Group { // CG utiliser le mot-cle Regions {} ?
/* One starts by giving explicit meaningful names to
the Physical regions defined in the "microstrip.msh" mesh file.
There are 2 volume regions and 3 surface regions in this model. */
......@@ -26,12 +28,12 @@ Group {
Electrode = Region[121];
SurfInf = Region[130];
/* We now define abstract regions to be used in the definition
of the scalar electric potential formulation:
/* We now define abstract regions to be used below
in the definition of the scalar electric potential formulation:
Vol_Dielectric_Ele : dielectric volume regions where "div epsr[] grad v = 0" is solved
Vol_Dielectric_Ele : dielectric volume regions where "Div ( epsr[] Grad v)" is solved
Sur_Dir_Ele : Dirichlet boundary condition (v imposed)
Sur_Neu_Ele : Neumann bondary condition ( epsr[] n.grad v = 0 )
Sur_Neu_Ele : Neumann bondary condition ( epsr[] n.Grad v = 0 )
Vol_xxx groups contain only volume elements of the mesh (triangles here).
Sur_xxx groups contain only surface elements of the mesh (lines here).
......@@ -43,14 +45,15 @@ Group {
}
Function {
/* Material laws (here the relative permittivity)
are defined piecewise in terms of the above defined physical regions */
/* The relative permittivity is defined piecewise using the above defined Groups */
epsr[Air] = 1.;
epsr[Diel1] = 9.8;
}
Constraint {
/* As for material laws, the Dirichlet boundary conditions can be defined piecewise
/* As for material laws, the Dirichlet boundary condition is defined piecewise.
The constraint "Dirichlet_Ele" is invoked in the FunctionSpace below */
{ Name Dirichlet_Ele; Type Assign;
......@@ -61,22 +64,32 @@ Constraint {
}
}
Group{
/* The domain of definition lists all regions on which the unknown field "v" is defined.
/* The function space of the electric scalar potential v
is definied by
- a domain of definition, which is the Group "Dom_Hgrad_v_Ele",
- a type (Form0 means scalar field)
- a set of scalar basis functions (here nodal basis functions "BF_Node")
- constraints (here the Dirichlet boundary conditions)
Contrary to the above defined groups, which contain either volure or surface elements,
Dom_xxx groups contain both volume and surface elements */
Contrary Vol_xxx and Sur_xxx regions, which contain only volume or surface regions, resp.,
Dom_xxx regions may contain both volume and surface regions.
Hence the use of the prefixes Vol_, Sur_ and Dom_ to avoid confusions.*/
Group{
Dom_Hgrad_v_Ele = Region[ {Vol_Dielectric_Ele, Sur_Dir_Ele, Sur_Neu_Ele} ];
}
FunctionSpace {
/* The function space in which we shall pick the electric scalar potential v solution
is definied by
- a domain of definition ("Dom_Hgrad_v_Ele")
- a type ("Form0" means scalar field)
- 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 fields now reads
v = Sum_k sn_k vn_k
where the "vn_k" are the nodal values and "sn_k" the basis functions.
Not all nodal values are unknowns of the problem, due to the "Constraint"
which assigns particular values to the nodes of the region "Sur_Dir_Ele". */
{ Name Hgrad_v_Ele; Type Form0;
BasisFunction {
{ Name sn; NameOfCoef vn; Function BF_Node;
......@@ -108,6 +121,41 @@ Integration {
}
Formulation {
/* The syntax of the Formulation{} section is a hard nut to crack.
So let's deal with it carefully.
The weak form of the electrostatic problem is particularly simple in this model,
as it has only one term:
(epsr[] Grad v, Grad vp)_Vol_Dielectric_Ele = 0 for all vp in S0(Vol_Dielectric_Ele).
The corresponding Euler-Lagrange are:
* Div ( epsr[] Grad v) = 0 on Vol_Dielectric_Ele
* epsr[] n.Grad v = 0 on Sur_Neu_Ele
The Galerkin{} statement is a symbolic representation of this weak formulation term.
The first argument represents the density to be integrated over finite elements,
the result of which being the bilinear form, i.e., the matrix
that will be assembled in the linear system of FE equations.
What is a bit confusing is that the two comma-separated arguments of the bracket []
are not exactly interpreted the same way.
The symbol "d" represents the exterior derivative and is a synonym of "Grad".
The expression "d v" stands then for the gradient of the v field,
i.e., for the electric field up to a minus sign.
The Dof{} operator, finally, returns the array of basis functions of its argument.
As the Galerkin method uses as trial (test) functions the basis functions
of the unknown field, the density of the bilinear form should be
[ epsr[] * Dof{d v} , Dof{d v} ]
to highlight the symmetry.
But the second argument of the density is always reserved for the trial functions
and the the Dof operator should therefore be there systematically.
It can thus be made implicit and, according to the GetDP syntax, it must be omitted.
Hence the final expression of the density below.
*/
{ Name Electrostatics_v; Type FemEquation;
Quantity {
{ Name v; Type Local; NameOfSpace Hgrad_v_Ele; }
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment