-
Notifications
You must be signed in to change notification settings - Fork 0
/
iteration_zg.m
110 lines (94 loc) · 2.3 KB
/
iteration_zg.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
close all;
clc;
clear all;
N=20;
ITER=350;
XHIST = zeros(ITER,N+1);
TEHIST = zeros(ITER,N);
RESHIST = zeros(ITER,N);
THIST = zeros(ITER,N);
SOLERRHIST = zeros(ITER,1);
WHIST = zeros(ITER,N+1);
ERRHIST = zeros(ITER,N);
% Anfangsgitter
%XI = linspace(0,1,N+1);
XMIN=0.0;
XMAX=1.0;
ALPHA=0.5;
XI = zeros(1,N+1);
for I=1:N+1
XI(I) = XMIN + (ALPHA^(I-1)-1)/(ALPHA^N-1)*(XMAX-XMIN);
end
%XI=[0,sort(rand(1,N-1)),1];
for A=1:ITER
XHIST(A,:) = XI;
% TE ausrechnen
[TERRI, RESI, TI,SERRI,ERRI] = dif1d_orth_it(N,XI);
XFINE = zeros(1,2*N+1);
for I=1:N
XFINE(2*I-1) = XI(I);
XFINE(2*I) = (XI(I+1)+XI(I))/2;
end
XFINE(2*N+1)=XI(N+1);
%finer solution
[TERRI2, RESI2, TI2,SERRI2,ERRI2] = dif1d_orth_it(N*2,XFINE);
% %fehler wzischen feinem und groebn gitter
TR = zeros(1,N);
ERR=0;
for I=1:N
TR(I) = (ERRI2(2*I)+ERRI2(2*I-1))-ERRI(I);
ERR=ERR+TR(I)^2;
end
SERR(A)=sqrt(ERR/(N));
% fprintf('Summierter Fehler %16.10e I=%g\n', SERR(A), A);
% Gitter anpassen
[XI,WI] = rref_te(N,XI,TR.-TI');
% Werte speichern
TEHIST(A,:) = TR;
RESHIST(A,:) =RESI;
THIST(A,:) = TI;
SOLERRHIST(A,:) = SERRI;
WHIST(A,:) = WI;
ERRHIST(A,:) = ERRI;
end
% Suche bestes Ergebnis im Vergleich zu 1
figure(1);
hold on;
X1 = XHIST(1,:);
XC1 = (X1(1:N)+X1(2:N+1))/2;
TE1 = TEHIST(1,:);
RES1 = RESHIST(1,:);
plot(XC1,TE1,'kx-');
plot(XC1,RES1,'gx-');
IDXMIN=1;
ERRMIN=SERR(1);
for A=1:ITER
if SERR(A) <= ERRMIN
IDXMIN=A;
ERRMIN=SERR(A);
end
end
X2 = XHIST(IDXMIN,:);
XC2 = (X2(1:N)+X2(2:N+1))/2;
TE2 = TEHIST(IDXMIN,:);
plot(XC2,RESHIST(IDXMIN,:),'rx-');
plot(XC2,TE2,'bx-');
legend('TE1','RES1','TE2','RES2')
figure(2);
hold on;
semilogy(1:ITER, SERR)
semilogy(1:ITER, ERRMIN*ones(1,ITER),'g-')
semilogy(1:ITER, SERR(1)*ones(1,ITER),'r-')
legend('Summierter Abbruchfehler','Anfangsfehler','Kleinster Fehler');
fprintf('Verbesserung relativ: %4.2f%%\n', abs(ERRMIN-SERR(1))/SERR(1)*100);
figure(3);
hold on;
semilogy(1:ITER, SOLERRHIST)
semilogy(1:ITER, SOLERRHIST(1)*ones(1,ITER),'g-')
semilogy(1:ITER, SOLERRHIST(IDXMIN)*ones(1,ITER),'r-')
legend('Summierter Fehler','Anfangsfehler','Kleinster Fehler');
fprintf('Verbesserung absoluter Fehler relativ: %4.2f%%\n', abs(SOLERRHIST(1)-SOLERRHIST(IDXMIN))/SOLERRHIST(1)*100);
figure(4);
hold on;
plot(XC1, ERRHIST(1,:),'g-')
plot(XC2, ERRHIST(IDXMIN,:),'r-')