| 1 |
- function U = statiqueU
% Calcul, assemblage et resolution d'un probleme statique KU = F
% resolution apres elimination des CL (ligne - colonne)
% Nous suivons la demarche proposee dans le cours "MEF" Chapitre II
%================================================================
% appel : U = statiqueU
% en sortie valeurs nodale du champ inconnu
%
% H.Oudin
% Vos commentaires pour ameliorer ce script seront bienvenus.
% mailto : herve.oudin@ec-nantes.fr
%================================================================
global nddln nddlt nelt nnode
global Connec Typel Ncl Vcl F
t0=cputime; %----- initialisation du temps de calcul
K=zeros(nddlt); %----- initialisation de la matrice rigidite
Fg=zeros(nddlt,1); %----- initialisation du vecteur force generalisee
for iel=1:nelt %----- boucle sur les elements
[Ke,Fe] = feval(Typel(iel,:),iel);
loce=[]; %----- localisation des ddl de l'element
for i=1:nnode
if Connec(iel,i) > 0
loce=[loce,(Connec(iel,i)-1)*nddln+[1:nddln]];
end
end;
K(loce,loce)=K(loce,loce) + Ke; %----- assemblage
Fg(loce)=Fg(loce) + Fe;
end
F = F + Fg; %----- vecteur du chargement
ir=0;
for i=1:nddlt %----- prise en compte des conditions aux limites
if ( Ncl(i) == 1 ) %----- valeurs imposes dans F
F = F - K(:,i)*Vcl(i); ir=ir+1;
end
end
for i=nddlt:-1:1
if ( Ncl(i) == 1 )
K(i,:) = []; K(:,i) = []; %----- supression ligne colonne dans K
F(i)=[]; %----- supression ligne dans F
end
end
%disp(K);disp(F);
Sol = K\ F; %----- resolution
%disp(Sol);
ir=0;
for i=1:nddlt %----- on replace les valeurs imposees
if ( Ncl(i) == 1 ) U(i,1) = Vcl(i);
else ir=ir+1; U(i,1)=Sol(ir);
end
end
%(' ');fprintf(1,' temps de calcul %10.4f secondes\n',cputime-t0)
|
