Modification de treillis_cours.m pour traiter le problème du TP

Work/tp1_treillis_sujet.m

@@ -0,0 +1,146 @@
+close all ; clc;
+%  Script de calcul statique du treillis traite le TP1
+% 
+% Fonctions utilisees
+%   statiqueUR      : calcul de la reponse statique [U,R]
+%   plotstr         : trace du maillage avec ne noeuds et elements
+%   plotdef         : trace de la deformee
+%   resultante      : calcul de la resultante en x, y, z d'un vecteur nodal
+%   barre_stress    : calcul de la contrainte dans un element barre
+%
+% Initialisation des variables globales pour un treillis
+%   nddln : nb de ddl par noeud
+% 	nnod  : nb de noeuds 
+%  	nddlt : nb de ddl total(=ndln*nnod)   
+%	  nelt  : nb d'elements  
+% 	nnode : nb de noeuds par element (2)
+% 	ndim  : dimension du probleme (1D,2D ou 3D)
+%   ncld  : nb de conditions de champ impose (dirichlet)
+%
+% 	Coord(nnod,ndim): coordonnees des noeuds
+%   Connec(nelt,2)	: connectivites	des elements
+%   Typel(nelt)     : Type des elements (barre_ke)
+% 	Nprop(nelt)		  : Ne de caracteristique pour chaque element
+%   Prop(nprop,ncar): Tableau des caracteristiques mecaniques (ES, f)     
+%	  Ncl(1,nddlt)	  : vaut 1 si le ddl est impose (deplacements imposes)
+%	  Vcl(1,nddlt)	  : valeur du deplacement impose 
+%	  F (nddlt,1)		  : vecteur des charges nodales donnees
+%
+% L'objectif de scripts de donnees est d'initialiser ses variables globales,
+% avant de lancer les calculs et d'exploiter les resultats (post-traitement).
+% 
+global nddln nnod nddlt nelt nnode ndim ncld
+global Coord Connec Typel Nprop Prop Ncl Vcl F
+disp(' ');
+disp('Structure etudiee : treillis traite dans TP1');
+disp('Travail de Jovian');
+disp('Unites : mm, MPa, kN'); % 1N / mm^2 = 1MPa
+disp('==================');
+
+% Definition du maillage (X, Y)
+Coord = [ 0   , 0 ; ...
+          500 , 0 ; ...
+          300 , -150 ];
+
+% Nombre de noeuds (3) et dimension du problème (3D)
+[nnod, ndim] = size(Coord);
+
+% Deux degrés de libertés par noeud : u (deplacement selon x) et v (selon y)
+nddln = 2;
+
+% Nb de ddl total
+nddlt = nddln * nnod;
+
+% Definition de la matrice de connectivite i , j
+Connec=[ 1 , 2 ; ...
+         1 , 3 ; ...
+         2 , 3 ];
+
+[nelt, nnode] = size(Connec);
+
+% Definition du modele EF : type des elements
+Typel = 'barre_ke';
+
+for i = 1 : nelt
+    Typel = char('barre_ke', Typel);
+end
+
+% Definition des caracteristiques mecaniques elementaires (ES fx fy)
+E = 210000 %MPa
+S1 = 4900 %mm^2
+S2 = 2500 %mm^2
+
+Nprop = [1; 2; 2];
+
+% Tableau : ES fx fy
+Prop = [  E*S2 0 0; ...
+          E*S1 0 0];
+      
+% definition des CL en deplacement
+% Ne du noeud, type sur u et v (1 ddl impose ,0 ddl libre)
+CL = [  1 , 0 , 1 ; ...
+        2 , 1 , 1 ];
+     
+Ncl = zeros(1, nddlt);
+ncld = 0;
+Vcl = zeros(1, nddlt);  % Valeurs des deplacements imposes
+%Vcl(2)=1; % à utiliser pour imposer une valeur non nulle sur un ddl i
+for i = 1 : size(CL, 1)
+   for j = 1 : nddln 
+       if CL(i, 1+j) == 1 
+           Ncl(1, (CL(i,1)-1) * nddln + j) = 1;
+           ncld = ncld + 1; 
+       end
+   end
+end
+
+% Definition des charges nodales (Ne du noeud , Fx , Fy)
+Charg = [ 3 0 -950000 ];
+
+% Vecteur sollicitation
+F = zeros(nddlt, 1);
+for iclf = 1:size(Charg, 1)           
+	noeud = Charg(iclf,1);  
+	for i = 1 : nddln
+       F((noeud-1)*nddln+i)=F((noeud-1)*nddln+i) + Charg(iclf,i+1);
+    end
+end
+
+[Fx, Fy, Fz] = feval('resultante', F); % Resultante des charges nodales
+
+disp('Les variables globales sont initialisees');
+disp('Fin de lecture des donnees');
+
+% Trace du maillage pour validation des donnees 
+plotstr                     
+reponse = input('Voulez-vous continuer? O/N [O]: ','s');
+if isempty(reponse) | reponse =='O'
+    U = zeros(nddlt, 1);
+    R = zeros(nddlt, 1);
+    [U(:, 1), R(:, 1)] = statiqueUR; % Resolution du probleme
+% Format d'impression des vecteurs
+    form =' %8.3e   %8.3e   %8.3e  '; format = [form(1:8*nddln),' \n']; 
+    disp(' ');disp('------- deplacements nodaux sur (x,y,z) ----------');
+    fprintf(format,U)
+    plotdef(U)
+% Post-traitement
+    disp(' ');disp('------- Efforts aux appuis  ----------');
+    fprintf(format,R(:,1));
+    [Rx,Ry,Rz] = feval('resultante', R); % Resultantes et reactions
+    disp(' ');
+    fprintf('La resultante des charges nodales    en (x,y,z) est : %8.3e   %8.3e   %8.3e \n',Fx,Fy,Fz);                    
+    fprintf('La resultante des charges reparties  en (x,y,z) est : %8.3e   %8.3e   %8.3e \n',-Rx-Fx,-Ry-Fy,-Rz-Fz);
+    fprintf('La resultante des efforts aux appuis en (x,y,z) est : %8.3e   %8.3e   %8.3e \n',Rx,Ry,Rz);
+    disp(' ');disp('------- Contraintes sur les elements ----------');
+    for iel=1:nelt          %----- boucle sur les elements
+        loce=[]; for i=1:nnode loce=[loce,(Connec(iel,i)-1)*nddln+[1:nddln]];end
+        Ue=U(loce);
+        Ne = feval('barre_stress',iel,Ue);
+        fprintf('Dans l''element %3i l''effort normal est %8.3e\n',iel,Ne)
+    end 
+    clear all
+    return
+else
+disp(' ');disp('---------------- arret du calcul----------------');
+clear all
+end