From: Peter Schaefer Date: Mon, 2 Apr 2012 18:43:31 +0000 (+0200) Subject: [src] optionale Konditionierung eingebaut X-Git-Url: https://git.leopard-lacewing.eu/?a=commitdiff_plain;h=3f76ae0e889d3ee9e155c7142acad9166529b421;p=bacc.git [src] optionale Konditionierung eingebaut [src] test_sol Skript weiter ausgebaut --- diff --git a/src/A_plots.m b/src/A_plots.m index 12188c3..8288ae5 100644 --- a/src/A_plots.m +++ b/src/A_plots.m @@ -59,7 +59,7 @@ for i = 1:length(files) }'; leg3 = {leg3{:}... - ['cond(A_{fine})' l0 l1{i}]... + ['cond(A_{h/2})' l0 l1{i}]... }'; sym = {sym{:} type2sym{data(1,[2+(i-1)*rows])}}'; end @@ -157,14 +157,14 @@ loglog(repmat(X(:,i+1),1,3),[... G_D(:,2+6+rows*i)... ],type2sym{i+1}); end -loglog(X(:,1),[7*X(:,1).^(-1/2),3*X(:,1).^(-1/4),2*X(:,1).^(-3/4)],'-.') +% loglog(X(:,1),[7*X(:,1).^(-1/2),3*X(:,1).^(-1/4),2*X(:,1).^(-3/4)],'-.') hold off title('hmin hmax') xlabel('Elemente'); ylabel('Schaetzer'); legend({leg2{:} ... - 'N^{-1/2}' 'N^{-1/4}' 'N^{-3/4}'... +% 'N^{-1/2}' 'N^{-1/4}' 'N^{-3/4}'... } ,'location','southwest','box', 'off'); print('-r600','-depsc',[printt '_hminmax.eps']) @@ -182,15 +182,15 @@ loglog(repmat(X(:,i+1),1,1),[... G_D(:,2+7+rows*i)... ],type2sym{i+1}); end -loglog(X(:,1),[7*X(:,1).^(-1/2),3*X(:,1).^(-1/4),2*X(:,1).^(-3/4)],'-.') +% loglog(X(:,1),[7*X(:,1).^(-1/2),3*X(:,1).^(-1/4),2*X(:,1).^(-3/4)],'-.') hold off -title('hmin hmax') +title('KonditionierungsZahlen von V') xlabel('Elemente'); -ylabel('Schaetzer'); +ylabel('Condition'); legend({leg3{:} ... - 'N^{-1/2}' 'N^{-1/4}' 'N^{-3/4}'... - } ,'location','southwest','box', 'off'); +% 'N^{-1/2}' 'N^{-1/4}' 'N^{-3/4}'... + } ,'location','northwest','box', 'off'); print('-r600','-depsc',[printt '_hminmax.eps']) diff --git a/src/A_run.m b/src/A_run.m index 2b68d73..063936e 100644 --- a/src/A_run.m +++ b/src/A_run.m @@ -1,11 +1,12 @@ -function A_run(file,times,zeta,type,theta,nu,out) -% A_AnIso(file,times,zeta,type,theta,nu,out) +function A_run(file,times,zeta,type,theta,nu,vcon,out) +% A_AnIso(file,times,zeta,type,theta,nu,vcon,out) % file - starting mesh % times - how often % zeta - Zulaessigkeitsbed % type - art des Tests % theta - alle verfeinern oder nur wichtige? % nu - wie verfeinern iso oder aniso +% vcon - soll Vorkonditioniert werden? 0 1 % out - Dateiname der Ausgabe % % P. Schaefer @@ -36,7 +37,7 @@ for i = 1:times end % Ein kompletter Verfeinerungschritt - A_step(zeta,type,theta,nu,out); + A_step(zeta,type,theta,nu,vcon,out); %Zeit Speichern usedTime = G_T(end,:); diff --git a/src/A_step.m b/src/A_step.m index 6f3fa5c..ed89d9d 100644 --- a/src/A_step.m +++ b/src/A_step.m @@ -1,10 +1,11 @@ -function [data time er] = A_step(zeta,type,theta,nu,varargin) -% [data time er] = A_step(zeta,type,eta,eps [, out]) +function [data time er] = A_step(zeta,type,theta,nu,vcon,varargin) +% [data time er] = A_step(zeta,type,theta,nu,vcon [, out]) % Führt einen Verfeinerungsschritt aus % % zeta & type - Bestimmen die Art der Matrix Berechnung -% eta - adaptiv? -% eps - isotrop? +% theta - adaptiv? +% nu - isotrop? +% vcon - Vorkonditionierung der Matrix? 1 oder 0 % out - (optional) Dateizusatz um Netz, Lösung & co zu speichern % % P. Schaefer @@ -66,18 +67,21 @@ time = zeros(1,3); plotMark([r';c'],coo_fine,ele_fine) title('Fehlerhafte Elemente') end - + + if(~vcon) %Lösung Berechnen - x_fine = V_fine\b_fine; - + x_fine = V_fine\b_fine; + con = cond(V_fine); + else %Vorkonditionierte Lösung! -% D = diag(diag(V_fine.^(-1/2))); -% -% A = D * V_fine * D; -% c = D*b_fine; -% y = A\c; -% x_fine = D*y; - + D = diag(diag(V_fine.^(-1/2))); + + A = D * V_fine * D; + c = D*b_fine; + y = A\c; + x_fine = D*y; + con = cond(A); + end % \tilde \mu ( \Pi h -h + L_2 ) tmu = hmin.* b .* sum((x_fine(f2s)'-repmat(sum(x_fine(f2s)',1)/4,4,1)).^2)' /4; @@ -85,15 +89,16 @@ time = zeros(1,3); %Fehlerschätzer 2 aufbauen V = mex_build_V(G_C,G_E,zeta,type(i)); - - x = V\b; - -% D = diag(diag(V.^(-1/2))); -% -% A = D * V * D; -% c = D*b; -% y = A\c; -% x = D*y; + if(~vcon) + x = V\b; + else + D = diag(diag(V.^(-1/2))); + + A = D * V * D; + c = D*b; + y = A\c; + x = D*y; + end xo_fine(f2s) = repmat(x,1,4); xd_fine = xo_fine'-x_fine; @@ -123,7 +128,7 @@ time = zeros(1,3); data = [data type(i) sqrt(sum(tmu)) sqrt(eta) xe sqrt(sum(mu))... - min(hmin)/max(hmax) min(hmax)/max(hmax) min(hmin./hmax) cond(V_fine)]; + min(hmin)/max(hmax) min(hmax)/max(hmax) min(hmin./hmax) con]; end time(2) = toc; diff --git a/src/test_sol.m b/src/test_sol.m index 8c7858f..f2d7470 100644 --- a/src/test_sol.m +++ b/src/test_sol.m @@ -9,21 +9,25 @@ mex mex_build_V.cpp slpRectangle.cpp %Anzahl der Schritte steps = 33; -%LShape adaptiv anisotrop -A_run('exmpl_2DQuad', steps, 0.7, [1], 0.5, 0.5, 'testAA_') - -%LShape adaptiv isotrop -% A_run('exmpl_2DLShape', steps, 0, 0.5, 0, 0, 'testAI_') - -%LShape uniform anisotrop -% A_run('exmpl_2DLShape', steps, 0, 1, 0.5, 0, 'testUA_') +%Art der Berechnungen +type = [1]; -%LShape uniform isotrop -% A_run('exmpl_2DLShape', steps, 0, 1, 0, 0, 'testUI_') - - -%Ausgabe: -% data = AnzElement Typ Schätzer Energienorm^2 +%LShape adaptiv anisotrop +A_run('exmpl_2DQuad', steps, 0.7, type, 0.5, 0.5, 0, 'testAA_') +A_run('exmpl_2DQuad', steps, 0.7, type, 0.5, 0.5, 1, 'testAAvcon_') + +% %LShape adaptiv isotrop +% % A_run('exmpl_2DLShape', steps, 0, 0.5, 0, 0, 'testAI_') +% +% %LShape uniform anisotrop +% % A_run('exmpl_2DLShape', steps, 0, 1, 0.5, 0, 'testUA_') +% +% %LShape uniform isotrop +% % A_run('exmpl_2DLShape', steps, 0, 1, 0, 0, 'testUI_') +% +% +% %Ausgabe: +% % data = AnzElement Typ Schätzer Energienorm^2 %meshSave Dateien % test* - enthält die Ergebnisse (wie Ausgabe) sowie das aktuelle Netz @@ -31,4 +35,7 @@ A_run('exmpl_2DQuad', steps, 0.7, [1], 0.5, 0.5, 'testAA_') %Ergebnis Plotten und Teil der Daten ausgeben sowie abspeichern der figure -A_plots({['meshSave/testAA_1_' num2str(steps)]},'exmplAA_2DQuad') +typeN = int2str(type); +typeN = typeN(typeN ~= ' '); +A_plots({['meshSave/testAA_' typeN '_' num2str(steps)]... + ['meshSave/testAAvcon_' typeN '_' num2str(steps)]},'exmplAA_2DQuad')