clear all

% Berechnungsmodell (Torsionsproblem fuer Rechteckquerschnitt):
b  =  1 ; h  =  2 ;                % Breite und Hoehe des Rechtecks
nb = 30 ; nh = 60 ;                % nb*nh*8 Dreieckselemente
[ne , nk , xy , km , kr] = TorsnetD6 (b , h , nb , nh) ;
ep = [] ;
bk = zeros (nk,1) ;

% FEM-Gleichungssystem Ks*x=fq:
[succ Ks fq]   = syswbc_m (xy , km , ep , kr , bk) ;

Ks = SymBand2Sparse (Ks) ;         % Umwandlung in "Sparse matrix"
disp ('Gleichungssystem für Torsionsproblem') ;
[m,n]=size(Ks) ;
disp (['Anzahl der Gleichungen:              ' , num2str(m)]) ;

t1 = cputime ;
xE = Ks \ fq ;                     % Matlab-Standard-Solver
disp (['Zeit für Standard-Solver:            ' num2str(cputime-t1)]) ;

t1 = cputime   ;
R  = chol (Ks) ;                   % Cholesky-Zerlegung
xC = R\(R'\fq) ;                   % Vorwaert-Rueckwaertseinsetzen
disp (['Zeit für chol + Vorw.-Rückw.-Eins.:  ' num2str(cputime-t1)]) ;

t1  = cputime ;
p   = symamd(Ks) ;                 % Fill-in-Optimierer
Ksf = Ks(p,p) ;                    % Optimierte Koeffizientenmatrix
fqf = fq(p) ;                      % Tauschoperationen f. rechte Seite
Rf  = chol(Ksf) ;
xf(p,1) = Rf\(Rf'\fqf) ;           % Ruecktausch nach Vorw.-Rueckw.-E.
disp (['Zeit Fill-in-optimiert:              ' num2str(cputime-t1)]) ;