Séance 3 : Exemple de problème numérique

Formation Python du département Mécatronique

14h00 - 15h30

Exemple de problème mécanique simple nécessitant de mettre en œuvre un algorithme numérique (pivot de Gauss pour résoudre « Ax=b »).

Le problème choisi est l'analyse statique d'un treillis pour calculer les actions mécaniques s'exerçant sur les barres

  • Description du problème.
  • Construction du problème et 1er tracés du treillis (mise en œuvre de Matplotlib)
  • Construction du système d'équations linéaires (matrices « A et b ») de façon automatisée (en empruntant un peu de théorie des graphs)
  • Résolution de l'équation « Ax=b » par pivot de Gauss. Comparaison avec numpy.linalg
  • Représentation graphique des tractions/compressions.

A) Description du problème du treillis

Le treillis est un ensemble de barres connectées entre elles "en triangles". Ce type de structure est par exemple utilisé dans les grues et certains ponts :

Pont en treillis

Le contact entre les barres est souvent modélisé par des liaisons pivots. Par conséquent, Lorsque la structure est soumis à un chargement (par ex. une masse accrochée au bout de la grue) les barres travaillent en traction/compression.

Pont en treillis

L'objectif de du calcul est de quantifier ces actions mécaniques.

Pourquoi le problème du treillis ?

  • Il s'agit d'un problème inspiré de la mécanique,
  • aboutissant à un système d'équations linéaires ("Ax=b") que l'on peut résoudre par pivot de Gauss
  • Ces équations seront issues du principe de la statique.

Remarque : dans un domaine tout autre (électricité), on obtiendrait des équations similaires en cherchant les courants dans un réseau de résistances.

B) Construction du problème

Avant de résoudre Ax=b, il faut construire la matrice A et le vecteur b, mais pour construire ces matrices, il faut d'abord créer les données qui décrivent le problème.

À propos de la généricité du programme :

  • On va ici s'intéresser à un treillis particulier qu'on pourrait sans doute traiter analytiquement...
  • ...mais tout l'intérêt d'un tel programme est de résoudre un treillis de taille arbitraire
  • ceci explique que le programme doit être parfois écrit de façon plus complexe que nécessaire en apparence.
  • le script S3_Statique_Treillis/treillis.py permet de fabriquer un treillis avec un nombre arbitraire de pivot et une forme y=f(x) ajustable.

B1) Construction des données

On construit d'abord la liste contenant la position (dans un plan) des pivots. On peut choisir (par exemple) d'utiliser une liste de tuples :

In [1]:
pivots = [(0, 0), (1, 2), (2, 0), (3, 2), (4, 0)]

N_piv = len(pivots)
print('treillis à {} pivots (dont 2 pivots de fixation)'.format(N_piv))
treillis à 5 pivots (dont 2 pivots de fixation)

Représentation graphique (ou l'intérêt de travailler interactivement)

In [2]:
# Conversion en tablea NumPy pour bénéficier des facilités du slicing
pivots_arr = np.array(pivots)

# Tracés en boucle:
for piv in pivots:
    plt.plot(..., ..., 'o')

# Ajustements des limites :
xlim(-.5, 4.5)
ylim(-.5, 2.5);

Choix des points d'accroche du treillis :

In [3]:
# indice des points d'accroche:
iP_A = 0 # 1er point du treillis
iP_B = 1 # 2ème point

P_A = pivots[iP_A]
P_B = pivots[iP_B]

On créer ensuite la liste des barres. Ici on choisit une liste de paires contenant les points. On aurait pu aussi utiliser une liste de paires contenant les indices de ces points dans la liste pivots.

In [4]:
barres = [
 ((0, 0), (2, 0)), # barre entre le pt (0,0) et le pt (2,0)
 ((1, 2), (2, 0)),
 ((1, 2), (3, 2)),
 ((2, 0), (3, 2)),
 ((2, 0), (4, 0)),
 ((3, 2), (4, 0))]

N_bar = len(barres)
print('treillis à {} barres'.format(N_bar))
treillis à 6 barres

Représentation graphique :

In [5]:
color=(0.6, 0.6, 0.8)

# Tracé des barres, en boucle :
for j, bj in enumerate(barres):
    # Coordonnées des 2 pivots d'accroche:
    (x1,y1), (x2, y2) = bj
    color = (.5, .5, .8)
    plt.plot((x1, x2), (y1, y2), '-', color = color, lw=4, zorder=1)

# Tracés des pivots :
for piv in pivots:
    plt.plot(piv[0], piv[1], 'ro')

# Ajustements des limites :
xlim(-.5, 4.5)
ylim(-.5, 2.5);

Description du chargement des pivots

Pour décrire une force extérieure (2 composantes) s'appliquant au niveau de chaque liaison, on utilise un tableau NumPy de taille $(N_{piv}, 2)$.

Question 1: Créer F_ext, un tableau de zéros de dimension N_piv * 2

In [6]:
F_ext = np.zeros(...)
F_ext.T # .T signifie "transposée". (Utilisé ici pour l'esthétique)
Out[6]:
array([[ 0.,  0.,  0.,  0.,  0.],
       [ 0.,  0.,  0.,  0.,  0.]])

Question 2 : Modifier les éléments F_ext pour décrire une force de 1 N vers le bas, s'excerçant sur le dernier pivot

In [8]:
# Appui sur le dernier pivot:

...

F_ext.T 
Out[8]:
array([[ 0.,  0.,  0.,  0.,  0.],
       [ 0.,  0.,  0.,  0., -1.]])

B2) Construction des équations (matrices A et b)

Les inconnues sont la force de traction qui s'exerce sur chaque barre (6). "Traction" indique le choix (arbitraire) de convention de signe :

  • $f_j>0$ si la barre j est "étirée"
  • $f_j<0$ si la barre j est "comprimée"

le vecteur des inconnues est $x=(f_1, \dots, f_{Nbar})$. On veut construire les matrices A et b telles que $A.x = b$

Les équations se basent de l'application du principe de la statique à chaque pivot.

$$\sum_j \vec{F}(\text{barre}_j\rightarrow \text{pivot}_i) = \vec{0}, \quad \text{pour chaque pivot $i$}$$

Équation d'équilibre d'un pivot

Remarque importante : la complexité de ce qui suit eu égard à la simplicité du problème (5 pivots, 6 barres) est dû à la volonté d'avoir un programme qui puisse résoudre le problème générique ($N_{piv}$ pivots arbitrairement grand)

La somme des forces des équations d'équilibre ne doit se faire que sur les barres incidentes au pivot considéré. Pour cela on introduit la notion de matrice d'incidence (outil issu de la théorie des graphs, cf. fr.wikipedia.org/wiki/Matrice_d'incidence)

Matrice d'incidence : définition

Il s'agit de considérer le treillis comme un graph orienté :

  • les sommets sont les pivots
  • les arcs sont les barres. Chaque barre (Piv1, Piv2) sera considéré orienté de Piv1 vers Piv2. Cette orientation dirigera leur vecteur directeur.

On va construire la matrice d'incidence $M$ qui contient l'information :

  • $m_{ij} = +1$ si la barre $j$ "arrive" pivot $i$
  • $m_{ij} = -1$ si la barre $j$ "quitte" pivot $i$

Équation d'équilibre d'un pivot

Grâce à la matrice d'incidence, on peut connaître la force excercée par une barre sur n'importe quelle pivot.

$$\vec{F} (\text{barre}_j\rightarrow \text{pivot}_i ) = - M[i,j] . f_j . \vec{u}_j$$

Grâce à cette formule, on peut construire les équations d'équilibre génériques. Il reste à :

  • construire la matrice d'indidence
  • calculer les vecteur directeurs $\vec{u}_j$

Matrice d'incidence : construction

In [10]:
Inc_mat = np.zeros((N_piv, N_bar), dtype=int)
Inc_mat
Out[10]:
array([[0, 0, 0, 0, 0, 0],
       [0, 0, 0, 0, 0, 0],
       [0, 0, 0, 0, 0, 0],
       [0, 0, 0, 0, 0, 0],
       [0, 0, 0, 0, 0, 0]])

On va remplir cette matrice colonne par colonne (i.e barre par barre).

Pour trouver dans quel ligne placer les "+1" et les "-1", on a besoin d'une méthode pour retrouver l'indice d'un élement dans une liste

In [14]:
print(pivots)
pivots.index((1,2))
# pivots.index((9,9)) # la recherche d'un élement inexistant génère une ValueError
[(0, 0), (1, 2), (2, 0), (3, 2), (4, 0)]
Out[14]:
1
In [15]:
for j, bar_j in enumerate(barres):
    P1, P2 = bar_j
    # Remarquons la convention de signe:
    i1 = ...
    Inc_mat[...] = -1 # la barre j "quitte" P1
    i2 = ...
    Inc_mat[...] = +1 # la barre j "arrive à" P2

print('Matrice d\'incidence:')
print(str(Inc_mat).replace('0','.')) # Méthode d'affichage cosmétique
Matrice d'incidence:
[[-1  .  .  .  .  .]
 [ . -1 -1  .  .  .]
 [ 1  1  . -1 -1  .]
 [ .  .  1  1  . -1]
 [ .  .  .  .  1  1]]

Chaque ligne peut se lire comme l'équation d'équilibre de chaque pivot, indiquant quelle barre lui est incidente.

Retirer les pivots d'accroche

Les pivots d'accroches $P_A$ et $P_B$ ne sont pas à l'équilibre (c'est une force d'accroche extérieure supplémentaire qui permet l'équilibre)

On enlève de la matrice d'incidence les lignes correspondants aux deux pivots d'accroche à l'aide de l'indexation par liste

In [17]:
# Étape 1 : construire une liste d'indice contenant les indices de tous les pivots sauf P_A et P_B
piv_ind = range(N_piv)
piv_ind.remove(iP_A)
piv_ind.remove(iP_B)
piv_ind
Out[17]:
[2, 3, 4]
In [18]:
# Étape 2 : utiliser la liste `piv_ind` pour sélectionner les lignes de la matrice d'incidence :
Inc_mat_red = Inc_mat[piv_ind, :]

print(str(Inc_mat_red).replace('0','.'))
[[ 1  1  . -1 -1  .]
 [ .  .  1  1  . -1]
 [ .  .  .  .  1  1]]

Observation des équations

La matrice d'incidence réduite contient 3 lignes correspondant aux 3 pivots auxquels on veut appliquer le principe de la statique.

$$\sum_j \vec{F}(\text{barre}_j\rightarrow \text{pivot}_i) = \vec{0}, \quad \text{pour chaque pivot $i$}$$

Ce principe d'équilibre va donner 2 équations scalaires par pivot, soit 6 eq. au total.

Comme on a 6 barres en traction/compression, on a bien 6 inconnues $f_j$ et le système d'équation devrait être soluble : le système mécanique est isostatique (à moins que la position des pivots mène à une matrice A de rang < 6).

Vecteur directeur des barres

Pour projeter les équations de la statique selon les axes (x,y), on a besoin du vecteur directeur de chaque barre.

In [43]:
# Conversion des barres en trableau NumPy
barres_arr = np.array(barres)
print(barres_arr.shape) # On obtient un tableau 3D
barres_arr
(6, 2, 2)
Out[43]:
array([[[0, 0],
        [2, 0]],

       [[1, 2],
        [2, 0]],

       [[1, 2],
        [3, 2]],

       [[2, 0],
        [3, 2]],

       [[2, 0],
        [4, 0]],

       [[3, 2],
        [4, 0]]])

Question de réflexion : à quoi correspond chacune des 3 dimensions de barres_arr (regarder pour cela la structure de barres)

Pour calculer le vecteur directeur, on calcule d'abord $P_1P_2 = OP_2 - OP_1$. Il faut donc d'abord récupérer $OP_1$ et $OP_2$ pour chaque barre.

On procède de manière vectorielle :

In [44]:
barres_OP1 = barres_arr[:,0,:]
barres_OP2 = barres_arr[:,1,:]

# À quoi correspond ce barres_OP1 ?
barres_OP1
Out[44]:
array([[0, 0],
       [1, 2],
       [1, 2],
       [2, 0],
       [2, 0],
       [3, 2]])

De même l'opération $P_1P_2 = OP_2 - OP_1$ est vectorisée :

In [30]:
barres_P1P2 = barres_OP2 - barres_OP1
barres_P1P2
Out[30]:
array([[ 2,  0],
       [ 1, -2],
       [ 2,  0],
       [ 1,  2],
       [ 2,  0],
       [ 1, -2]])

Exercice de calcul vectorisé : Il faut à présent normaliser les vecteurs P1P2 pour obtenir les vecteurs directeurs $\vec{u}$

$$ \vec{u} = \frac{P_1P_2}{|P_1P_2|} $$

In [35]:
# Étape 1 : calcul de la longueur des barres :
barre_l = np.sqrt(....)
barre_l
Out[35]:
array([ 2.        ,  2.23606798,  2.        ,  2.23606798,  2.        ,
        2.23606798])
In [37]:
# Étape 2 : Transformation en vecteur colonne : 
barre_l = ....
barre_l
Out[37]:
array([[ 2.        ],
       [ 2.23606798],
       [ 2.        ],
       [ 2.23606798],
       [ 2.        ],
       [ 2.23606798]])
In [39]:
# Étape 3 : normalisation vectorisée
barre_dir = barres_P1P2 / barre_l

barre_dir
Out[39]:
array([[ 1.        ,  0.        ],
       [ 0.4472136 , -0.89442719],
       [ 1.        ,  0.        ],
       [ 0.4472136 ,  0.89442719],
       [ 1.        ,  0.        ],
       [ 0.4472136 , -0.89442719]])

Construction des matrices A et b

La matrice $A$ est fabriquée par superposition de 2 blocs :

  • $A_x$ contient les équations d'équilibre projetée sur l'axe horizontal $A_x = \text{IncMatRed} \times \text{diag}(u_{j}^{(x)})$
  • $A_y$ contient les équations d'équilibre projetée sur l'axe vertical $A_y = \text{IncMatRed} \times \text{diag}(u_{j}^{(y)})$
In [45]:
# 1) Matrice A
Ax = Inc_mat_red*barre_dir[:,0]
Ay = Inc_mat_red*barre_dir[:,1]
# ou bien: Ax = np.dot(Inc_mat_red, np.diag(barre_dir[:,0]))

Ax.round(2)
Out[45]:
array([[ 1.  ,  0.45,  0.  , -0.45, -1.  ,  0.  ],
       [ 0.  ,  0.  ,  1.  ,  0.45,  0.  , -0.45],
       [ 0.  ,  0.  ,  0.  ,  0.  ,  1.  ,  0.45]])

La superposition (concaténation) des blocs se fait avec np.vstack (il existe np.hstack et np.concatenate)

In [42]:
A = np.vstack((Ax, Ay))
# Image de la matrice:
# plt.imshow(A, interpolation='nearest')

print(A.round(2))
[[ 1.    0.45  0.   -0.45 -1.    0.  ]
 [ 0.    0.    1.    0.45  0.   -0.45]
 [ 0.    0.    0.    0.    1.    0.45]
 [ 0.   -0.89  0.   -0.89 -0.   -0.  ]
 [ 0.   -0.    0.    0.89  0.    0.89]
 [ 0.   -0.    0.    0.    0.   -0.89]]

Le vecteur $b$ est construit selon le même principe de concaténation

In [41]:
# 2) Vecteur b : force extérieure appliquée à chaque pivot:
# Empilement des composantes selon x et y:
b_ext = np.concatenate((F_ext[piv_ind,0], 
                        F_ext[piv_ind,1]))
b_ext
Out[41]:
array([ 0.,  0.,  0.,  0.,  0., -1.])

C) Résolution du problème, pivot de Gauss

C1) Pivot de Gauss

Méthode d'élimination itérative pour résoudre un système d'équations linéaires. Cf Wikipedia pour le détail de la méthode.

In [178]:
# On assemble A et b (on "borde" la matrice)
Ab = np.hstack((A,b_ext.reshape(-1,1))) # attention au reshape
print(np.round(Ab,2))
[[ 1.    0.45  0.   -0.45 -1.    0.    0.  ]
 [ 0.    0.    1.    0.45  0.   -0.45  0.  ]
 [ 0.    0.    0.    0.    1.    0.45  0.  ]
 [ 0.   -0.89  0.   -0.89 -0.   -0.    0.  ]
 [ 0.   -0.    0.    0.89  0.    0.89  0.  ]
 [ 0.   -0.    0.    0.    0.   -0.89 -1.  ]]
In [101]:
N = A.shape[0]

for k in range(N):
    # Trouver le pivot
    piv_candidats = Ab[k:,k]
    k_piv = np.argmax(np.abs(piv_candidats)) + k
    piv = Ab[k_piv,k]
    if piv == 0:
        print(np.round(Ab,2))
        raise ValueError('Systeme non inversible')
    
    # échange des lignes k et k_piv:
    Ab[[k,k_piv], :] = Ab[[k_piv,k], :]
    
    # Normalisation de la ligne k
    Ab[k,k:] /= piv
    # Mise à zéros de la colonne k
    for i in range(N):
        if i == k:
            continue
        Ab[i,:] -= Ab[k,:]*Ab[i,k]

print(np.round(Ab,2))
[[ 1.    0.    0.    0.    0.    0.   -1.5 ]
 [ 0.    1.    0.    0.    0.    0.    1.12]
 [ 0.    0.    1.    0.    0.    0.    1.  ]
 [ 0.    0.    0.    1.    0.    0.   -1.12]
 [ 0.    0.    0.    0.    1.    0.   -0.5 ]
 [ 0.    0.    0.    0.    0.    1.    1.12]]

Une fois la matrice transformée en $(I_N|v)$ la solution du système $Ax=b$ est le vecteur $v$

In [104]:
# Extraction du vecteur (colonne N+1) :
trac_barres_gauss = Ab[:,N]
trac_barres_gauss.round(2)
Out[104]:
array([-1.5 ,  1.12,  1.  , -1.12, -0.5 ,  1.12])

C2) Résolution des force de traction avec numpy.linalg

Au delà de l'intérêt de comprendre l'algorithme du pivot de Gauss, il existe déjà une méthode de résolution de système linéaire : numpy.linalg.linsolv

(fonction du sous-module linalg qui contient aussi l'inversion de matrices, la diagonalisation, la décomposition en valeurs singulières, ...)

In [ ]:
np.linalg. # taper TAB

Exercice : Résoudre d'une manière ou d'une autre l'équation $A.x = b_{ext}$

In [56]:
trac_barres = ...
In [59]:
print('Effort de traction sur chaque barre:')
print(trac_barres.round(2))

# Recherche de l'effort maximal
trac_max = np.max(np.abs(trac_barres))
print(' -> effort max (en val. absolue) : {:.1f}'.format(trac_max))
Effort de traction sur chaque barre:
[-1.5   1.12  1.   -1.12 -0.5   1.12]
 -> effort max (en val. absolue) : 1.5

Connaissant les efforts sur chaque barrre, on peut en déduire les forces s'exerçant sur les points d'accroche

In [60]:
# Efforts sur les points d'accroche (de la part du treillis et de l'extérieur):
resul_A = -np.inner(Inc_mat[iP_A,:]*trac_barres, barre_dir.T) + F_ext[iP_A]
resul_B = -np.inner(Inc_mat[iP_B,:]*trac_barres, barre_dir.T) + F_ext[iP_B]

print('action du treillis sur pt A: {!s} ({:.2f})'.format(
      resul_A, np.linalg.norm(resul_A,2)) )
print('action du treillis sur pt B: {!s} ({:.2f})'.format(
      resul_B, np.linalg.norm(resul_B,2)) )
action du treillis sur pt A: [-1.5  0. ] (1.50)
action du treillis sur pt B: [ 1.5 -1. ] (1.80)

Remarque : on vérifie que la somme des actions du treills sur les pts d'accroche est égale à la somme des forces extérieur sur le treillis F_ext

In [61]:
F_ext.sum(axis=0)
Out[61]:
array([ 0., -1.])

D) Représentation graphique de la solution

Représentation graphique des tractions/compression avec Matplotlib

D1) Colorer les barres selon l'effort (signe et intensité)

Outil : les colormaps de Matplotlib.

Celles ci servent le plus souvent avec imshow pour coloriser des images en niveau de gris (ou simplement "représenter" les valeurs d'une matrice)

In [62]:
cmap = plt.cm.gray
#cmap = plt.cm.winter
#cmap = plt.cm.jet # couleurs par défaut, comme dans Matlab
plt.imshow(A, interpolation='none', cmap=cmap)
plt.title('matrice A');

Les colormaps sont en fait des fonctions qui transforme un nombre entre [0,1] (un "niveau de gris") en un quadruplet (R,G,B, alpha) (chaque composante dans [0,1])

In [63]:
col = plt.cm.jet(0.9)
col # On reconnait du rouge
#plot(0,0, 'D', color=col)
Out[63]:
(0.94563279857397531, 0.029774872912127992, 0.0, 1.0)

On peut fabriquer sa propre colormap les préexistantes ne conviennent pas. (chercher "colormaps demo" pour avoir un aperçu)

In [64]:
# Colormap pour colorer les barres selon l'effort subit: rouge-bleu
col_list = [(0.9,0,0.0), (0.7,0.7,0.7), (0,0,0.9)]
cmap = mpl.colors.LinearSegmentedColormap.from_list('red-blue', col_list)

On peut utiliser une colormap pour attribuer une couleur à une force $F$

$$cmap(\frac{F}{2 F_{max}} + 1/2) $$

In [65]:
# Tracé des barres et des efforts
for j, bj in enumerate(barres):
    # Coordonnées des 2 pivots d'accroche:
    (x1,y1), (x2, y2) = bj
    # direction :
    uj = barre_dir[j]
    # effort de traction sur la barre bj:
    trac = trac_barres[j]
    color = cmap(trac/(2*trac_max)+0.5)
    plt.plot((x1, x2), (y1, y2), '-', color = color, lw=4, zorder=1)


# Tracés des pivots :
for piv in pivots:
    plt.plot(piv[0], piv[1], 'wo')

# Ajustements des limites :
xlim(-.5, 4.5)
ylim(-.5, 2.5);

D2) Représenter les forces... avec des flèches !

In [172]:
color=(0.8, 0.8, 0.8)
# Échelle pour le tracé des forces
F_scale = 0.3*barre_l.mean()/trac_max
# Couleur des efforts:
F_color = (1., 0, 0) # rouge
dx, dy  = 1,2

# Tracé des barres et des efforts
for j, bj in enumerate(barres):
    # Coordonnées des 2 pivots d'accroche:
    (x1,y1), (x2, y2) = bj
    # direction :
    uj = barre_dir[j]
    trac = trac_barres[j]
    
    plt.plot((x1, x2), (y1, y2), '-', color = color, lw=4, zorder=1)
    
    # Tracé des efforts barre -> pivot1 et pivot 2
    plt.arrow(x1, y1, +trac*uj[0]*F_scale, +trac*uj[1]*F_scale,
              zorder=2, head_width=0.05*dx, lw=0, width=0.02*dx, color=F_color)
    plt.arrow(x2, y2, -trac*uj[0]*F_scale, -trac*uj[1]*F_scale,
              zorder=2, head_width=0.05*dx, lw=0, width=0.02*dx, color=F_color)

# Tracés des pivots :
for piv in pivots:
    plt.plot(piv[0], piv[1], 'wo', zorder=3)

# Ajustements des limites :
xlim(-.5, 4.5)
ylim(-.5, 2.5)

title(u'Forces exercées par les barres sur les liasons')

Remarque : on utilise l'argument zorder pour maîtriser la superposition des éléments

D3) Le code complet de visualisation des efforts :

In [118]:
fig = plt.figure('efforts treillis', figsize=(12,5))

ax = fig.add_subplot(111, title='efforts sur le treillis '
                                '({:d} pivots, {:d} barres)'.format(N_piv, N_bar))

# Échelle pour le tracé des forces
F_scale = 0.3*barre_l.mean()/trac_max
# Couleur des efforts:
F_color = (0,0.8,0) # vert
dx, dy  = 1,2

# Colormap pour colorer les barres selon l'effort subit: rouge-bleu
col_list = [(0.9,0,0.0), (0.7,0.7,0.7), (0,0,0.9)]
cm_rb = mpl.colors.LinearSegmentedColormap.from_list('red-blue', col_list)
cm = cm_rb
#cm = plt.cm.coolwarm_r

# Tracé des barres et des efforts
for j, bj in enumerate(barres):
    # Coordonnées des 2 pivots d'accroche:
    (x1,y1), (x2, y2) = bj
    # direction :
    uj = barre_dir[j]
    # effort de traction sur la barre bj:
    trac = trac_barres[j]
    color = cm(trac/(2*trac_max)+0.5)
    plt.plot((x1, x2), (y1, y2), '-', color = color, lw=4, zorder=1)
    
    # Tracé des efforts barre -> pivot1 et pivot 2
    plt.arrow(x1, y1, +trac*uj[0]*F_scale, +trac*uj[1]*F_scale,
              zorder=2, head_width=0.05*dx, lw=0, width=0.02*dx, color=F_color)
    plt.arrow(x2, y2, -trac*uj[0]*F_scale, -trac*uj[1]*F_scale,
              zorder=2, head_width=0.05*dx, lw=0, width=0.02*dx, color=F_color)
# end for each barre

# Couleur des pivots
piv_color = (1.,1.,1.) # blanc
piv_color_AB = (1.,1.,0.5) # jaune clair
piv_alpha = 1 # opaque


# Tracé des pivots
for i, Pi in enumerate(pivots):
    # Tracé du pivot:
    marker = 'D' if Pi in (P_A,P_B) else 'o' # marqueur Diamond 'D' ou disque 'o'
    color = piv_color_AB if Pi in (P_A,P_B) else piv_color
    plt.plot(Pi[0], Pi[1], marker, ms=8, c=color, alpha = piv_alpha, zorder=3)
    # Force extérieur s'appliquant sur le pivot
    Fi = F_ext[i]
    if Fi.any():
        plt.arrow(Pi[0], Pi[1], Fi[0]*F_scale, Fi[1]*F_scale,
                  zorder=2, head_width=0.05*dx, lw=0, width=0.02*dx, color=(1,0,0))
    if Pi in (P_A,P_B):
        F_soutien = -resul_A if Pi==P_A else -resul_B
        plt.arrow(Pi[0], Pi[1], F_soutien[0]*F_scale, F_soutien[1]*F_scale,
                  zorder=2, head_width=0.05*dx, lw=0, width=0.02*dx, color=(1,1,0))
# end for each pivot

# Limites du tracé
plt.xlim(min([x for (x,y) in pivots]) - dx*1,
         max([x for (x,y) in pivots]) + dx*1)
plt.ylim(min([y for (x,y) in pivots]) - dy*.3,
         max([y for (x,y) in pivots]) + dy*.3)

# Couleur de fond:
ax.patch.set_fc((0.9,)*3)
ax.set_aspect('equal')
plt.grid(False)

# "Maximisation" du graphique au sein de la figure :
fig.tight_layout()

fin de la session 3

Licence Creative Commons

Ce <span xmlns:dct="http://purl.org/dc/terms/" href="http://purl.org/dc/dcmitype/InteractiveResource" rel="dct:type">support de formation</span> créé par <span xmlns:cc="http://creativecommons.org/ns#" property="cc:attributionName">Pierre Haessig</span> est mis à disposition selon les termes de la licence Creative Commons Attribution 3.0 France.