NUM_exo1.m 3.4 KB

1234567891011121314151617181920212223242526272829303132333435363738394041424344454647484950515253545556575859606162636465666768697071727374757677787980818283848586
  1. close all ; clc;
  2. % Exercice 1 du chapitre de cours NUM - elements T3
  3. % Maillage en 2 element d'une plaque rectangulaire
  4. %
  5. % Fonctions utilisees
  6. % statique : calcul de la reponse statique [U,R]
  7. % plotstr : trace du maillage avec numero noeuds et elements
  8. % plotdef : trace de la deformee
  9. % resultante : calcul de la resultante en x, y, z d'un vecteur nodal
  10. % T3_ep_stress : calcul de la contrainte dans un element T3
  11. % plot_sig : visualisation de l'etat de contrainte
  12. %
  13. global nddln nnod nddlt nelt nnode ndim ncld
  14. global Coord Connec Typel Nprop Prop Ncl Vcl F
  15. disp(' ');
  16. disp('Maillage en deux T3 de plaque de l''exercice 1 du chapitre NUM');
  17. disp('===================');
  18. % definition du maillage
  19. L = 2 ; h = 1;
  20. Coord=[ 0 0; L 0; L h;0 h ]; % definition des coordonnees des noeuds
  21. [nnod,ndim]=size(Coord);
  22. nddln=2; nddlt=nddln*nnod;
  23. Connec=[ 1 2 4; ... % Tableau de connectivite i , j
  24. 2 3 4 ];
  25. [nelt,nnode]=size(Connec);
  26. Typel = 'T3_ep'; % definition du type des elements
  27. for i=1:nelt Typel = char('T3_ep',Typel); end
  28. % definition des caracteristiques mecaniques
  29. Nprop=[1;1;1;1];
  30. Prop=[ 210000 0.3 2 0 0 ]; % valeurs (E nu e fx fy)
  31. % CL en deplacement
  32. CL=[ 1 , 1 , 1; ...
  33. 2 , 0 , 1; ...
  34. 4 , 1 , 1;];
  35. Ncl=zeros(1,nddlt);ncld=0;
  36. Vcl=zeros(1,nddlt);
  37. for i=1:size(CL,1)
  38. for j=1:nddln
  39. if CL(i,1+j)==1 Ncl(1,(CL(i,1)-1)*nddln+j)=1;ncld=ncld+1; end
  40. end
  41. end
  42. % charges nodales: numero du noeud , Fx,Fy
  43. Charg=[ 3 5 0 ];
  44. F=zeros(nddlt,1);
  45. for iclf=1:size(Charg,1)
  46. noeud=Charg(iclf,1);
  47. for i=1:nddln
  48. F((noeud-1)*nddln+i)=F((noeud-1)*nddln+i) + Charg(iclf,i+1);
  49. end
  50. end
  51. [Fx,Fy,Fz] = feval('resultante',F); %----- resultante des charges nodales
  52. plotstr % trace du maillage
  53. U = zeros(nddlt,1);
  54. R = zeros(nddlt,1);
  55. [U(:,1),R(:,1)] = statique; % ----- resolution du probleme (pivot = 1)
  56. %----- post-traitement
  57. plotdef(U)
  58. form =' %8.3e %8.3e %8.3e '; format = [form(1:8*nddln),' \n'];
  59. disp(' ');disp('------- deplacements nodaux sur (x,y,z) ----------');
  60. fprintf(format,U)
  61. disp(' ');disp('------- Efforts aux appuis ----------');
  62. fprintf(format,R(:,1));
  63. [Rx,Ry,Rz] = feval('resultante',R); %----- resultantes et reactions
  64. disp(' ');
  65. fprintf('La resultante des charges nodales en (x,y,z) est : %8.3e %8.3e %8.3e \n',Fx,Fy,Fz);
  66. fprintf('La resultante des charges reparties en (x,y,z) est : %8.3e %8.3e %8.3e \n',-Rx-Fx,-Ry-Fy,-Rz-Fz);
  67. fprintf('La resultante des efforts aux appuis en (x,y,z) est : %8.3e %8.3e %8.3e \n',Rx,Ry,Rz);
  68. disp(' ');disp('------- Contraintes sur les elements ----------');
  69. for iel=1:nelt %----- boucle sur les elements
  70. loce=[];
  71. for i=1:nnode
  72. if Connec(iel,i) > 0 loce=[loce,(Connec(iel,i)-1)*nddln+[1:nddln]]; end
  73. end;
  74. Ue=U(loce);
  75. sig = feval([deblank(Typel(iel,:)),'_stress'],iel,Ue);
  76. % contrainte moyenne de Von Mises sur l'element
  77. sVM = sqrt(sig(1)^2+sig(2)^2-sig(1)*sig(2)+3*sig(3)^2);
  78. VM(iel)=sVM ;
  79. fprintf('Valeur Moyenne des contraintes dans l''element %3i\n',iel)
  80. fprintf('sigxx = %8.3e sigyy = %8.3e sigxy = %8.3e sigVM = %8.3e \n',sig(1),sig(2),sig(3),sVM)
  81. end
  82. plot_sig(VM)
  83. return