FINAL PROJECT
>> % Here is my program:
>> type cheby.m
function cheby(A,b,v)
format long;
[n o]=size(A);
lmax=input('Please enter upper bound for the eigenvalues of A: ');
lmin=input('Please enter lower bound for the eigenvalues of A: ');
m=input('Please enter the degree of the Chebyshev polynomial being used: ');
iter=input('Please enter desired number of cycles :');
for i=1:m
theta=((2*(i-1)+1)/(2*m))*pi;
tauinverse = ( ( lmax - lmin ) / 2.0 )*cos(theta) + ( ( lmax + lmin ) / 2.0 );
tau(i) = 1.0 / tauinverse;
end
%fid=fopen('cheby.txt','w');
for j=1:iter
% fprintf(fid,'Solution vector after %d cycle(s): \n',j);
for k=1:m
u=v-tau(k)*(A*v-b);
v=u;
end
% for l=1:n
% fprintf(fid,'u(%d) = %12.8f\n',l,u(l));
% end
end
%fclose(fid);
v
>> % Now for some test cases. The first test case is the model problem given
>> % in class.
>> load cheby1
>> cheby(A,b,v)
Please enter upper bound for the eigenvalues of A: 100
Please enter lower bound for the eigenvalues of A: 1
Please enter the degree of the Chebyshev polynomial being used: 10
Please enter desired number of cycles :10
v =
1.00000000000000
0.50000000017162
0.33333405299766
0.25000064334448
0.20000001051919
0.16666666666671
0.14285714285909
0.12500001239379
0.11111146984138
0.10000131085979
0.09091035278820
0.08333371344827
0.07692310714712
0.07142857171569
0.06666666666667
0.06250000000001
0.05882352961066
0.05555557286799
0.05263178613746
0.05000082214169
0.04762054019876
0.04545595680242
0.04347897541501
0.04166684897424
0.04000001937698
0.03846153898709
0.03703703703776
0.03571428571429
0.03448275862107
0.03333333366316
0.03225807701384
0.03125012024974
0.03030353323762
0.02941291222112
0.02857301653006
0.02777916827018
0.02702779501855
0.02631604379964
0.02564107034788
0.02500000319689
0.02439024395029
0.02380952380954
0.02325581395349
0.02272727272888
0.02222222267729
0.02173914262174
0.02127669941010
0.02083375860328
0.02040917706494
0.02000153649482
0.01960938024692
0.01923178424891
0.01886835065186
0.01851862247607
0.01818183041309
0.01785714331424
0.01754385965074
0.01724137931034
0.01694915254238
0.01666666671490
0.01639344584806
0.01612907740136
0.01587327292662
0.01562577699143
0.01538602360197
0.01515312504719
0.01492653777667
0.01470639337702
0.01449287595298
0.01428572701552
0.01408450737865
0.01388888888928
0.01369863013699
0.01351351351426
0.01333333387262
0.01315791465561
0.01298720075044
0.01282125026806
0.01265968768713
0.01250154761861
0.01234653374054
0.01219533798481
0.01204821088106
0.01190476211357
0.01176470588235
0.01162790697674
0.01149425317923
0.01136366873428
0.01123636506780
0.01111248375516
0.01099045149425
0.01086996440230
0.01075270218407
0.01063829787458
0.01052631578957
0.01041667967874
0.01031012729997
0.01020515011383
0.01010101044078
0.01000163355733
>> % Some other test cases:
>> A=[2 -1 0;-1 2 -1;0 -1 2]
A =
2 -1 0
-1 2 -1
0 -1 2
>> b=[3 -3 5]'
b =
3
-3
5
>> v=[1 1 1]'
v =
1
1
1
>> cheby(A,b,v)
Please enter upper bound for the eigenvalues of A: 3.42
Please enter lower bound for the eigenvalues of A: 0.5
Please enter the degree of the Chebyshev polynomial being used: 3
Please enter desired number of cycles :7
v =
2.00000313031781
0.99999546693576
3.00000313031767
>> A=[2 -1 1;-1 2 -1;1 -1 2]
A =
2 -1 1
-1 2 -1
1 -1 2
>> b=[9 -13 10]'
b =
9
-13
10
>> v=[1 1 1]'
v =
1
1
1
>> cheby(A,b,v)
Please enter upper bound for the eigenvalues of A: 4
Please enter lower bound for the eigenvalues of A: 1
Please enter the degree of the Chebyshev polynomial being used: 3
Please enter desired number of cycles :7
v =
1.00000005655918
-4.99999998384024
2.00000004443935
>> % It looks like the program works, and the ordering of the tau's does
>> % not matter. I believe the ordering of the tau's matters when we chose
>> % C to be a matrix other than the identity matrix.
>> % The concern is when the matrix has a large condition number. As an experiment I used
>> % the following matrix:
>> A
A =
1 0 0 0
0 10 0 0
0 0 20 0
0 0 0 5000
>> % With b=[1 1 1 1]' and v=[1 1 1 1]'. I used p=4 and found that the order of the tau's
>> % did not matter, but the program needed a large amount of cycles to converge. In fact
>> % the number of cycles was equal to 5000 (the largest eigenvalue):
>> cheby(A,b,v)
Please enter upper bound for the eigenvalues of A: 5000
Please enter lower bound for the eigenvalues of A: 1
Please enter the degree of the Chebyshev polynomial being used: 4
Please enter desired number of cycles :5000
v =
1.00000000000000
0.10000000000000
0.05000000000000
0.00020000000001
>> % This was also true for the matrix with 10,000 instead of 5,000.
>> A(4,4)=10000;
>> cheby(A,b,v)
Please enter upper bound for the eigenvalues of A: 10000
Please enter lower bound for the eigenvalues of A: 1
Please enter the degree of the Chebyshev polynomial being used: 4
Please enter desired number of cycles :10000
v =
1.00000000000000
0.10000000000000
0.05000000000000
0.00010000000001
>>% Actually with further experimentation I found that it converged with p=400 and only
>>% one cycle:
>> cheby(A,b,v)
Please enter upper bound for the eigenvalues of A: 10000
Please enter lower bound for the eigenvalues of A: 1
Please enter the degree of the Chebyshev polynomial being used: 400
Please enter desired number of cycles :1
v =
1.00000000000000
0.10025868704727
0.04939644815039
0.00010000000000
>> diary off
The last experiment involved using different values for tau. Instead of theta = ( (2l + 1)/2p )*pi, I let theta be 2l+1/2p. This implies that 0.5403 < cos(theta) < 1 which gives ( 1.5403*lmax + .4597*lmin ) / 2 < 1/tau < lmax. This means that 1/lmax < tau < 2 / ( 1.5403*lmax +.4597*lmin ). Just a guess: our tau's are closer to 2 / ( lmax + lmin ) which was the optimal tau. Below, you will find the graphs ( the y-axis being the values of the vector x; the x-axis being the indices of the entries in x ) of the solutions of our model problem, the first with the tau's from the book, and the second with the new tau's. Both are with p=10 and 3 cycles. You can see that the method with the new tau's is converging faster. The third graph is the method with the new tau's applied to Randy's matrix ( diagonal entries of our model matrix squared ). Here p=100 with 3 cycles and convergence to the solution is very fast.
