La matrice de raideur lorsque la solution est un champ vectoriel

General Discussion
1

comment recuperer la matrice de raideur lorsque la solution est un champ vectoriel ?
voila mon code que je developpe
// Définir la forme variationnelle pour la matrice de raideur
varf a([u, v], [w, s]) = int2d(th)(
lambda * div(w, s) * div(u, v)
+ 2. * mu * ( epsilon(w, s)’ * epsilon(u, v) )
)
+ on(4, u = 0, v = 0);

// Définir la forme variationnelle pour le vecteur de charges
varf L([uu, vv], [w, s]) = int1d(th, 2)(P * w);

2

Hello,
You have an example here:

3

// Paramètres
int n = 100;
real rayon = 5;
real E = 210; // en GPa
real sigma = 0.2;
real P = 1.; // traction en MPa

// Bords du carré (sens trigonométrique)
border a(t = -20, 20){ x = t; y = -20; label = 1; } // bas
border b(t = -20, 20){ x = 20; y = t; label = 2; } // droite
border c(t = 20, -20){ x = t; y = 20; label = 3; } // haut
border d(t = 20, -20){ x = -20; y = t; label = 4; } // gauche

// Trou circulaire au centre (sens horaire pour le trou intérieur)
border cercle(t = 0, 2*pi){ x = rayon * cos(t); y = rayon * sin(t); label = 5; }

mesh th = buildmesh(b(n) + c(n) + d(n) + a(n) + cercle(-n));

// Espace
fespace Vh(th, [P1, P1]);

// Déformation & divergence
real sqrt2 = sqrt(2.);
macro epsilon(u1, u2) [dx(u1), dy(u2), (dy(u1)+dx(u2))/sqrt2] //
macro div(u,v) (dx(u) + dy(v)) //

// Paramètres matériaux
real mu = E/(2*(1+sigma));
real lambda = Esigma/((1+sigma)(1-2*sigma));

// Définir la forme variationnelle pour la matrice de raideur
varf a([u, v], [w, s]) = int2d(th)(
lambda * div(w, s) * div(u, v)

      • mu * ( epsilon(w, s)’ * epsilon(u, v) )
        )
  • on(4, u = 0, v = 0);

// Définir la forme variationnelle pour le vecteur de charges
varf L([uu, vv], [w, s]) = int1d(th, 2)(P * w);

// Assemblage
matrix K = avar(Vh2, Vh2, tgv = -1); // Matrice du système
real[int] b = Lrhs(0, Vh2); // Vecteur second membre

// Export K vers fichier CSV
ofstream fk(“kk.csv”);
for (int i = 0; i < K.n; ++i)
for (int j = 0; j < K.n; ++j)
if (abs(K(i,j)) > 1e-14) // ignorer les zéros
fk << i << " " << j << " " << K(i,j) << endl;

// Export rhs vers fichier CSV
ofstream fp(“pp.csv”);
for (int i = 0; i < rhs.n; ++i)
fp << b[i] << endl;

// Exporter les coordonnées des sommets
ofstream fxy(“xy.csv”);
for (int i = 0; i < th.nv; ++i) {
fxy << th(i).x << " " << th(i).y << endl;
}
il me donne une erreur de compilation à la ligne 36 a cause d,un point virgule !
aidez moi s’il vous plais

4

the error message says “The identifier a exists”.
Your varf needs another name instead of “a” because “a” is a border in your code.
→ call it “avar” “varf avar([u, v], [w, s]) =…”
Then the name “Vh2” has to be replaced by “Vh”
“Lrhs” has to be replaced by “L”
“b” needs also another name that you need to choose, different from the variables you already have in your code (“b” is a border in your code).

5

Merciiiiiii merciiiiiii beaucoup :folded_hands:t2::white_heart:

6

Bonjour Mr
j’espère que vous allez bien et en bon santé
Si je veux récupérer la matrice de raideur, quelle est la différence entre :

varf avar([u, v], [w, s]) = int2d(th)(
lambda * div(w, s) * div(u, v)

      • mu * (epsilon(w, s)’ * epsilon(u, v))
        )
  • on(1, u = 0, v = -1)
  • on(5, u = 0, v = 0);

// Assemblage
matrix K = avar(Vh, Vh, tgv = -1); // Matrice du système

(((((( et les cas où tgv = 1 et tgv = 0 ? )))))))

7

The effect of the value of tgv for Dirichlet boundary condition is explained in the documentation

in the “Note” just before the section “Numerical Integration”.
tgv=-1 is the most accurate, it replaces the lines of the linear system corresponding to degrees of freedom of the test function on the considered boundary, by the exact equation u=Dirichlet value.
tgv=large value (1e30 is the default) uses a penalty term similar to tgv*\int_\Gamma (u-u_D)*w, where w is the test function.