WEBVTT 00:00:00.040 --> 00:00:02.470 The following content is provided under a Creative 00:00:02.470 --> 00:00:03.880 Commons license. 00:00:03.880 --> 00:00:06.920 Your support will help MIT OpenCourseWare continue to 00:00:06.920 --> 00:00:10.570 offer high quality educational resources for free. 00:00:10.570 --> 00:00:13.470 To make a donation, or view additional materials from 00:00:13.470 --> 00:00:17.400 hundreds of MIT courses, visit MIT OpenCourseWare at 00:00:17.400 --> 00:00:18.650 ocw.mit.edu. 00:00:20.816 --> 00:00:23.260 PROFESSOR: Ladies and gentlemen, welcome to 00:00:23.260 --> 00:00:24.690 lecture number 12. 00:00:24.690 --> 00:00:27.510 In this lecture, I would like to discuss with you solution 00:00:27.510 --> 00:00:29.900 methods for eigenproblems. 00:00:29.900 --> 00:00:32.500 This is a very large field, and all I can do really in 00:00:32.500 --> 00:00:36.970 this lecture is to give you or discuss this you a few of the 00:00:36.970 --> 00:00:39.670 solution methods that we are using in 00:00:39.670 --> 00:00:41.390 finite element analysis. 00:00:41.390 --> 00:00:44.710 So it will be a very brief introduction, a brief survey 00:00:44.710 --> 00:00:48.390 of the methods of solution that we're using for finite 00:00:48.390 --> 00:00:50.440 element eigenvalue problems. 00:00:50.440 --> 00:00:54.010 The generalized eigenvalue problem, K phi equals lambda M 00:00:54.010 --> 00:00:57.350 phi, is a problem that we encountered in the multiple 00:00:57.350 --> 00:00:59.230 position analysis. 00:00:59.230 --> 00:01:03.480 Please recognize that I'm calling lambda here, the 00:01:03.480 --> 00:01:05.550 eigenvalue, which was our omega 00:01:05.550 --> 00:01:07.940 squared in the last lecture. 00:01:07.940 --> 00:01:09.950 This is really the eigenproblems that we are 00:01:09.950 --> 00:01:14.120 primarily addressing ourselves to finite element analysis. 00:01:14.120 --> 00:01:17.980 The standard eigenproblems problem here, EVP stands for 00:01:17.980 --> 00:01:21.660 eigenproblem, has M equal to the identity matrix. 00:01:21.660 --> 00:01:23.800 So when we're talking about solution methods of this 00:01:23.800 --> 00:01:26.470 problem, we're really also talk about solution methods of 00:01:26.470 --> 00:01:27.460 this problem. 00:01:27.460 --> 00:01:29.530 If we have, in a [? multiple ?] position 00:01:29.530 --> 00:01:33.660 analysis, non-proportional damping, then we would have to 00:01:33.660 --> 00:01:37.040 solve this eigenproblem, which is a quadratic eigenproblem. 00:01:37.040 --> 00:01:39.610 I will not address myself to the solution of this problem, 00:01:39.610 --> 00:01:42.850 so in what follow, we are talking about a solution of 00:01:42.850 --> 00:01:45.650 this problem, and of course, the solution of this problem 00:01:45.650 --> 00:01:49.030 also, when M is equal to the identity matrix. 00:01:49.030 --> 00:01:51.920 In particular, we are interested in the solution of 00:01:51.920 --> 00:01:56.120 large eigenvalue problems when N, is say, greater than 500. 00:01:56.120 --> 00:02:00.120 N, of course, being the order of K and M, and with the half 00:02:00.120 --> 00:02:02.900 bandwidth, I should say, being greater than 60. 00:02:02.900 --> 00:02:05.800 The number of eigenvalues that we will be looking for, that 00:02:05.800 --> 00:02:10.139 we want to calculate, are say, 1 to 1/3 of N. 1/3 of N is 00:02:10.139 --> 00:02:12.460 already large number of eigenvalues. 00:02:12.460 --> 00:02:16.840 Generally, we are talking about 1 to and at most, 100 00:02:16.840 --> 00:02:18.660 for large eigenproblems. 00:02:18.660 --> 00:02:21.710 When we have a small eigenvalue problem, in other 00:02:21.710 --> 00:02:25.530 words, N being small, say 50 or 100, there are a large 00:02:25.530 --> 00:02:27.940 number of different techniques that can be used to solve this 00:02:27.940 --> 00:02:31.970 eigenvalue problem, and it really does not make that much 00:02:31.970 --> 00:02:34.600 difference which technique you are using. 00:02:34.600 --> 00:02:39.610 The cost can still be quite significant, but certainly, 00:02:39.610 --> 00:02:43.550 when we're talking about very large eigenproblems, N being 00:02:43.550 --> 00:02:47.370 even larger than a 1,000 or 2,000, then it is of up most 00:02:47.370 --> 00:02:52.000 concern to select a technique that is very effective. 00:02:52.000 --> 00:02:55.670 And I want to address myself really now towards those 00:02:55.670 --> 00:02:59.590 techniques that we feel are very effectively used in 00:02:59.590 --> 00:03:01.960 finite element analysis. 00:03:01.960 --> 00:03:04.650 In dynamic analysis, of course, with proportional 00:03:04.650 --> 00:03:07.930 damping, this was the eigenvalue problem that we 00:03:07.930 --> 00:03:09.290 wanted to solve. 00:03:09.290 --> 00:03:12.360 Now, if they are 0 frequencies present, in other words, when 00:03:12.360 --> 00:03:15.540 we have rigid body modes in the system, then we can use 00:03:15.540 --> 00:03:18.750 the following procedure to basically getting rid of these 00:03:18.750 --> 00:03:19.630 0 frequencies. 00:03:19.630 --> 00:03:23.040 We simply add, on the left hand side, a mu M phi, and on 00:03:23.040 --> 00:03:25.950 the right hand side, a mu M phi. 00:03:25.950 --> 00:03:28.880 And taking these two together, we are getting this 00:03:28.880 --> 00:03:32.740 coefficient matrix, and a new lambda value. 00:03:32.740 --> 00:03:35.460 That lambda value now is omega squared plus mu, where omega 00:03:35.460 --> 00:03:38.630 squared is equal to lambda minus mu. 00:03:38.630 --> 00:03:41.870 What we have been doing here, is basically, we have been 00:03:41.870 --> 00:03:43.350 applying a shift. 00:03:43.350 --> 00:03:46.210 In other words, if we were to plot the characteristic 00:03:46.210 --> 00:03:50.560 polynomial, giving us the eigenvalues, the 00:03:50.560 --> 00:03:53.600 characteristic polynomial of course, being p of lambda 00:03:53.600 --> 00:03:57.940 being determined K bar minus lambda M, where K bar is equal 00:03:57.940 --> 00:04:00.490 to K plus mu M. 00:04:00.490 --> 00:04:03.270 If we plug this characteristic polynomial, we will see this 00:04:03.270 --> 00:04:06.960 curve here, schematically drawn, of course, where the 00:04:06.960 --> 00:04:10.250 omega squared values would start here, in other words, 00:04:10.250 --> 00:04:14.370 the lowest omega squared value is, in fact, 0. 00:04:14.370 --> 00:04:17.589 Whereas our lambda value now starts right here. 00:04:17.589 --> 00:04:19.500 This is the shift that we [? imposed. ?] 00:04:19.500 --> 00:04:25.755 So, basically, we now recognize that we have all 00:04:25.755 --> 00:04:29.090 eigenvalues being greater than 0, because this is our new 00:04:29.090 --> 00:04:30.100 starting value. 00:04:30.100 --> 00:04:33.150 All eigenvalues measured from here are greater than 0, the 00:04:33.150 --> 00:04:35.780 smallest one being mu, namely when omega 00:04:35.780 --> 00:04:37.200 squared is equal to 0. 00:04:37.200 --> 00:04:39.610 And the larger ones being here. 00:04:39.610 --> 00:04:45.260 It is effective to perform such a shift when we have 0 00:04:45.260 --> 00:04:49.960 eigenvalues, in other words, to operate on this eigenvalue, 00:04:49.960 --> 00:04:54.790 rather than on that one, because our techniques then 00:04:54.790 --> 00:04:58.330 are very stable and effective, and this shift here is, in 00:04:58.330 --> 00:05:00.950 general, therefore performed. 00:05:00.950 --> 00:05:04.540 So all we have to address ourselves to then from now on, 00:05:04.540 --> 00:05:08.390 is really, the solution of eigenvalues and eigenvectors, 00:05:08.390 --> 00:05:12.470 when all of the eigenvalues are greater than 0. 00:05:12.470 --> 00:05:15.720 If we start off with a problem where we have a 0 frequency, 00:05:15.720 --> 00:05:19.160 we apply the shift so that the eigenvalues now of this 00:05:19.160 --> 00:05:23.280 problem here are all greater than 0. 00:05:23.280 --> 00:05:26.680 In buckling analysis, we have this eigenproblem. 00:05:26.680 --> 00:05:28.185 We really did not talk about-- 00:05:28.185 --> 00:05:31.700 we did not really at all about buckling analysis in the set 00:05:31.700 --> 00:05:35.930 of lectures, but you might just be interested in finding 00:05:35.930 --> 00:05:38.140 out how would we solve for these eigenvalues. 00:05:38.140 --> 00:05:40.320 Well, here is the characteristic polynomial. 00:05:40.320 --> 00:05:43.390 The smallest lambda value, of course, gives now a critical 00:05:43.390 --> 00:05:46.270 load, that is the one e would be interested in. 00:05:46.270 --> 00:05:48.640 It is effective here, to rewrite this 00:05:48.640 --> 00:05:50.700 problem in this form. 00:05:50.700 --> 00:05:53.970 All we do, is we are taking the lambda value onto the 00:05:53.970 --> 00:05:54.920 other side. 00:05:54.920 --> 00:05:56.660 In other words, we are dividing this 00:05:56.660 --> 00:05:58.570 equation here by a lambda. 00:05:58.570 --> 00:06:01.320 So that one of the lambda stands on this side, and KG 00:06:01.320 --> 00:06:02.730 phi on that side. 00:06:02.730 --> 00:06:05.940 The resulting problem is this one, where kappa is 1 over 00:06:05.940 --> 00:06:08.540 lambda and we have now here, is a K matrix. 00:06:08.540 --> 00:06:14.630 K is positive definite, and now it is more effective to 00:06:14.630 --> 00:06:17.040 operate on this equation. 00:06:17.040 --> 00:06:19.820 How do we obtain now the largest kappa? 00:06:19.820 --> 00:06:22.130 Because the largest kappa is the one that we are interested 00:06:22.130 --> 00:06:25.010 in, which corresponds to the smallest lambda, the smallest 00:06:25.010 --> 00:06:25.730 buckling load. 00:06:25.730 --> 00:06:28.850 Well, we would now perform a shift to the right, 00:06:28.850 --> 00:06:34.190 subtracting mu K from this side, from this side, and from 00:06:34.190 --> 00:06:37.360 that side, to obtain this eigenvalue problem. 00:06:37.360 --> 00:06:39.980 That means we are really shifting from here up to 00:06:39.980 --> 00:06:45.130 there, and now we are interested in solving kappa N 00:06:45.130 --> 00:06:47.620 by operating on this eigenproblems problem. 00:06:47.620 --> 00:06:50.100 It is this eigenvalue problem that we can operate on 00:06:50.100 --> 00:06:53.540 directly, using the techniques that I would be discussing 00:06:53.540 --> 00:06:55.950 later, to determine search method and the subspace 00:06:55.950 --> 00:06:58.222 iteration method. 00:06:58.222 --> 00:07:00.990 But before we get into those techniques, I would like to 00:07:00.990 --> 00:07:05.810 discuss with you some preliminary considerations 00:07:05.810 --> 00:07:07.710 that are important when we look at 00:07:07.710 --> 00:07:09.440 an eigenvalue solution. 00:07:09.440 --> 00:07:14.030 The important first point that I want to make is that there 00:07:14.030 --> 00:07:17.260 is a traditional approach for the solution of this 00:07:17.260 --> 00:07:18.840 eigenvalue problem here. 00:07:18.840 --> 00:07:20.610 Let's look at this eigenvalue problem. 00:07:20.610 --> 00:07:22.840 Where of course, when we're talking about a buckling 00:07:22.840 --> 00:07:25.700 analysis, we would have different matrices here, but 00:07:25.700 --> 00:07:29.860 basically, we still solve this eigenvalue problem. 00:07:29.860 --> 00:07:34.160 The traditional approach lies in factorizing the M matrix 00:07:34.160 --> 00:07:36.120 into Cholesky factors. 00:07:36.120 --> 00:07:38.880 We talked about the Cholesky factor when we talked about 00:07:38.880 --> 00:07:40.230 solution of equations. 00:07:40.230 --> 00:07:41.750 So here it comes up again. 00:07:41.750 --> 00:07:45.440 We would factorize the M matrix into Cholesky factors, 00:07:45.440 --> 00:07:46.470 as shown here. 00:07:46.470 --> 00:07:50.730 We define now a new eigenvector, phi curl being 00:07:50.730 --> 00:07:54.880 equal as shown here, L curl transpose times 5. 00:07:54.880 --> 00:07:59.970 We substitute from here into there, and we directly obtain 00:07:59.970 --> 00:08:02.640 this eigenvalue problem here. 00:08:02.640 --> 00:08:05.980 Now we have the standard eigenproblem. 00:08:05.980 --> 00:08:09.160 Now once we have a standard eigenproblem, there are a 00:08:09.160 --> 00:08:11.660 large number of techniques that can be used. 00:08:11.660 --> 00:08:14.735 For example, the QR method, the standard 00:08:14.735 --> 00:08:16.790 Jacobi method, et cetera. 00:08:16.790 --> 00:08:18.940 So this is, what we might call. 00:08:18.940 --> 00:08:21.530 a traditional approach to solve the generalized 00:08:21.530 --> 00:08:24.140 eigenproblem, this is the one here that we are really 00:08:24.140 --> 00:08:25.110 interested in. 00:08:25.110 --> 00:08:27.250 But it involves this transformation, and this 00:08:27.250 --> 00:08:32.280 transformation here can involve a very significant 00:08:32.280 --> 00:08:35.070 cost, because we have to find these Cholesky factors and 00:08:35.070 --> 00:08:37.929 then, of course, we have to calculate this K curl. 00:08:37.929 --> 00:08:40.370 Another way of performing a transformation of the 00:08:40.370 --> 00:08:43.549 generalized item problem to a standard form would be to 00:08:43.549 --> 00:08:46.540 solve the eigenproblem of M first. 00:08:46.540 --> 00:08:49.110 W here and W transpose [? store ?] 00:08:49.110 --> 00:08:54.080 the eigenvectors of the M matrix, and D squared is a 00:08:54.080 --> 00:08:57.320 diagonal matrix of the eigenvalues of M. And then 00:08:57.320 --> 00:09:01.020 having factored M in this form, an equation similar to 00:09:01.020 --> 00:09:02.510 this one here can be written. 00:09:02.510 --> 00:09:05.490 Where now, the transformation, of course, involves the 00:09:05.490 --> 00:09:07.400 eigenvectors, W. 00:09:07.400 --> 00:09:12.860 Well, this approach is sometimes still quite 00:09:12.860 --> 00:09:15.990 effective, but in general, it cannot be recommended because 00:09:15.990 --> 00:09:20.420 there's too much cost involved in the factorizing of the M 00:09:20.420 --> 00:09:23.830 matrix, and then particular, in the transformation here, 00:09:23.830 --> 00:09:26.180 and the solution of this problem and the, of course, to 00:09:26.180 --> 00:09:29.670 get from this problem, the phi vector, we have to go back to 00:09:29.670 --> 00:09:30.980 this equation. 00:09:30.980 --> 00:09:36.070 A direct solution of the generalized eigenvalue problem 00:09:36.070 --> 00:09:39.390 that we want to solve is really more effective. 00:09:39.390 --> 00:09:43.680 And now, we are ordering the lambdas as shown here. 00:09:43.680 --> 00:09:48.310 Remember, our smallest lambda now is greater than 0, the 00:09:48.310 --> 00:09:52.490 eigonvectors are stored as shown here, and for most of a 00:09:52.490 --> 00:09:56.150 position analysis, we are really interested in solving 00:09:56.150 --> 00:10:00.550 only for a certain number of eigen pairs, i equals 1 to p, 00:10:00.550 --> 00:10:04.200 as discussed in the last lecture, or i equals r to s. 00:10:04.200 --> 00:10:06.980 If we have a vibration excitation problem, we really 00:10:06.980 --> 00:10:09.800 only want to solve for specific frequencies and 00:10:09.800 --> 00:10:11.810 corresponding mode shapes in a certain 00:10:11.810 --> 00:10:14.300 interval, as shown here. 00:10:14.300 --> 00:10:19.000 The solution procedures on this generalized eigenvalue 00:10:19.000 --> 00:10:23.690 problem can be categorized into the following forms. 00:10:23.690 --> 00:10:29.780 The first solution procedure, or a whole class of solution 00:10:29.780 --> 00:10:33.100 procedures, that I should really say, which we might 00:10:33.100 --> 00:10:36.050 call the vector iteration techniques. 00:10:36.050 --> 00:10:37.620 The idea is the following-- 00:10:37.620 --> 00:10:43.160 if we want to solve this problem, then why not solve it 00:10:43.160 --> 00:10:45.870 in the following way? 00:10:45.870 --> 00:10:55.830 Start off with a vector on this side, that we assume, 00:10:55.830 --> 00:11:00.260 calculate a new vector on this side, then take that new 00:11:00.260 --> 00:11:03.130 vector, plug it back on the right hand side, and cycle 00:11:03.130 --> 00:11:04.570 through this way. 00:11:04.570 --> 00:11:06.640 This is indeed inverse iteration. 00:11:06.640 --> 00:11:07.810 What do we do here? 00:11:07.810 --> 00:11:11.740 We are putting a vector on the right hand side, x1, when K is 00:11:11.740 --> 00:11:15.790 equal to 1, that vector we have to assume. 00:11:15.790 --> 00:11:18.540 Then we are multiplying out the right hand side. 00:11:18.540 --> 00:11:23.380 Now we have a set of simultaneous equations that we 00:11:23.380 --> 00:11:26.460 would solve just in the same way as we solve the equations 00:11:26.460 --> 00:11:27.480 in static analysis. 00:11:27.480 --> 00:11:30.190 We can think of this as a load vector. 00:11:30.190 --> 00:11:34.100 We are calculating this displacement vector, and then 00:11:34.100 --> 00:11:38.780 having got this displacement vector, and not taking care of 00:11:38.780 --> 00:11:42.390 this lambda, or having, basically, taken this lambda 00:11:42.390 --> 00:11:44.270 out of this equation. 00:11:44.270 --> 00:11:49.090 Which means that this vector will grow in lengths, or 00:11:49.090 --> 00:11:50.760 collapse in lengths. 00:11:50.760 --> 00:11:54.160 Therefore, want to scale this vector down to a good length-- 00:11:54.160 --> 00:11:57.100 a length that does not hit overflow, so to say, or 00:11:57.100 --> 00:11:58.550 underflow in the machine. 00:11:58.550 --> 00:12:02.260 And that scaling is achieved as shown here. 00:12:02.260 --> 00:12:05.980 So, the final result then, this is just a scalar here, as 00:12:05.980 --> 00:12:10.520 you can see, is a vector here, which we can put on to the 00:12:10.520 --> 00:12:13.800 right hand side and keep on iterating around. 00:12:13.800 --> 00:12:18.110 In fact, we find as that iteration procedure converges 00:12:18.110 --> 00:12:22.610 to the lowest eigenvector, provided the starting vector 00:12:22.610 --> 00:12:28.210 is not orthogonal, is not M orthogonal to phi 1. 00:12:28.210 --> 00:12:32.330 This is one first iteration technique. 00:12:32.330 --> 00:12:37.130 Maybe I should put a side note down. 00:12:37.130 --> 00:12:42.940 A side note that many of you might have in your mind. 00:12:42.940 --> 00:12:44.730 Many of you might have the following question in your 00:12:44.730 --> 00:12:48.350 mind-- why he iterate, why does he not solve directly for 00:12:48.350 --> 00:12:49.290 the eigenvalues? 00:12:49.290 --> 00:12:53.230 Just like we saw a set of simultaneous equations. 00:12:53.230 --> 00:12:57.890 Well, there's one very important point, and this is 00:12:57.890 --> 00:12:58.570 the point-- 00:12:58.570 --> 00:13:02.800 if we are solving an eigenvalue problem, and shown 00:13:02.800 --> 00:13:05.350 here, then really that means the following-- 00:13:05.350 --> 00:13:10.920 we are saying that determinant, K minus lambda M 00:13:10.920 --> 00:13:12.680 shall be equal to 0. 00:13:12.680 --> 00:13:13.260 Why? 00:13:13.260 --> 00:13:15.740 Because if we write this equation in 00:13:15.740 --> 00:13:17.080 the following form-- 00:13:17.080 --> 00:13:24.760 K minus lambda M times phi equals 0, then we have here a 00:13:24.760 --> 00:13:30.390 set of simultaneous equations with a 0 vector on 00:13:30.390 --> 00:13:31.730 the right hand side. 00:13:31.730 --> 00:13:34.700 And there is an important theorem in linear algebra that 00:13:34.700 --> 00:13:39.960 says, only a solution to this set of equations when this 00:13:39.960 --> 00:13:42.380 coefficient matrix is singular. 00:13:42.380 --> 00:13:44.920 Which means that the determinant of this 00:13:44.920 --> 00:13:47.700 coefficient matrix must be 0. 00:13:47.700 --> 00:13:49.240 So, what do we want? 00:13:49.240 --> 00:13:54.650 We really want to solve p of lambda determinant, K minus 00:13:54.650 --> 00:14:00.480 lambda M. We want to solve for the 0's of this polynomial. 00:14:00.480 --> 00:14:02.270 That's really what we want. 00:14:02.270 --> 00:14:07.700 Now, if we want to solve for the 0's of this polynomial, 00:14:07.700 --> 00:14:13.390 then we also have to understand that when k is of 00:14:13.390 --> 00:14:17.210 order 4, there will be four 0's. 00:14:17.210 --> 00:14:20.080 The polynomial itself will be of order 4. 00:14:20.080 --> 00:14:23.710 If K is of order, N, then there will be N 0's. 00:14:23.710 --> 00:14:29.680 So the polynomial is order N. However, if we want to find 00:14:29.680 --> 00:14:34.320 the solution of a polynomial, or order N, then we very 00:14:34.320 --> 00:14:40.650 quickly find that, in fact, we can only look up the solution 00:14:40.650 --> 00:14:45.080 to such polynomials if they are of low enough order, lower 00:14:45.080 --> 00:14:46.660 than say, 4. 00:14:46.660 --> 00:14:51.080 Above for larger order polynomials, there's no way 00:14:51.080 --> 00:14:53.150 you can look up the solution. 00:14:53.150 --> 00:14:58.160 You have to iterate to obtain the solution of the 0's of 00:14:58.160 --> 00:14:59.620 this polynomial. 00:14:59.620 --> 00:15:04.480 However, since this problem here, is completely equivalent 00:15:04.480 --> 00:15:06.430 to this problem. 00:15:06.430 --> 00:15:09.170 To solve this problem, you have to iterate. 00:15:09.170 --> 00:15:11.480 And this is the important conclusion. 00:15:11.480 --> 00:15:15.640 It is of up most importance that you understand that to 00:15:15.640 --> 00:15:20.570 solve an eigenvalue problem of order larger than say, 4, you 00:15:20.570 --> 00:15:22.330 will have to iterate. 00:15:22.330 --> 00:15:27.290 Because you're solving for the 0's of a polynomial of order 00:15:27.290 --> 00:15:28.600 larger than 4. 00:15:28.600 --> 00:15:32.040 Only for very small order polynomials, we can look up 00:15:32.040 --> 00:15:34.060 the solutions in handbooks. 00:15:34.060 --> 00:15:38.680 So having grasped the fact that we have to iterate to 00:15:38.680 --> 00:15:44.460 solve for the 0's of this polynomial, to solve for the 00:15:44.460 --> 00:15:51.390 eigenvalues of this problem, we are stuck with iteration, 00:15:51.390 --> 00:15:54.320 and since we're stuck with iteration, we now want to come 00:15:54.320 --> 00:15:56.380 up with techniques to iterate. 00:15:56.380 --> 00:15:58.620 And the first technique to iterate is this vector 00:15:58.620 --> 00:16:01.180 iteration technique. 00:16:01.180 --> 00:16:04.470 And the basic idea, once more, is if you want to solve this 00:16:04.470 --> 00:16:09.270 equation, well, then why not iterate in the following way? 00:16:09.270 --> 00:16:13.150 Put some starting value in here, this number is just a 00:16:13.150 --> 00:16:16.580 scalar, let's neglect that effect for the moment. 00:16:16.580 --> 00:16:19.640 And having put a starting value in here, calculate this 00:16:19.640 --> 00:16:22.790 one and keep on cycling around this way. 00:16:22.790 --> 00:16:25.600 And that's what we call inverse iteration. 00:16:25.600 --> 00:16:28.480 Inverse iteration here, a vector inverse iteration, 00:16:28.480 --> 00:16:32.680 because we have to solve a set of equations here, with K as 00:16:32.680 --> 00:16:33.930 the coefficient matrix. 00:16:36.310 --> 00:16:38.290 There are other iteration techniques. 00:16:38.290 --> 00:16:39.900 There's also forward iteration. 00:16:39.900 --> 00:16:42.150 And basically, what we're saying there is-- 00:16:42.150 --> 00:16:45.970 well if we are starting off here and going around in this 00:16:45.970 --> 00:16:49.970 circle, why not start off here and go around in that circle, 00:16:49.970 --> 00:16:51.190 the red line. 00:16:51.190 --> 00:16:54.970 And that, in fact, gives us forward iteration. 00:16:54.970 --> 00:16:57.970 One can prove that that forward iteration converges to 00:16:57.970 --> 00:17:00.140 the largest eigenvalue. 00:17:00.140 --> 00:17:03.190 We are usually interested in the smallest frequencies, the 00:17:03.190 --> 00:17:06.069 smallest eigenvalues, and therefore, we are using 00:17:06.069 --> 00:17:08.810 inverse iteration quite commonly. 00:17:08.810 --> 00:17:11.119 Then there is Rayleigh Quotient iteration. 00:17:11.119 --> 00:17:16.790 The Rayleigh Quotient iteration is again, an inverse 00:17:16.790 --> 00:17:21.160 iteration, however, with an acceleration scheme in it, a 00:17:21.160 --> 00:17:24.190 very effective acceleration scheme. 00:17:24.190 --> 00:17:26.810 This Rayleigh Quotient iteration can be employed to 00:17:26.810 --> 00:17:30.070 calculate one eigenvalue and vector, and then, of course, 00:17:30.070 --> 00:17:34.010 we have to deflate afterwards the matrix, the coefficient 00:17:34.010 --> 00:17:37.170 matrix for the eigenvector that we have calculated 00:17:37.170 --> 00:17:42.750 already, and we then could go on towards the calculation of 00:17:42.750 --> 00:17:44.770 additional eigenvalues and vectors. 00:17:44.770 --> 00:17:47.490 The same deflation also would be used in inverse iteration 00:17:47.490 --> 00:17:50.200 alone, and in forward iteration alone. 00:17:50.200 --> 00:17:52.760 In a Rayleigh Quotient iteration, we converge to an 00:17:52.760 --> 00:17:57.440 eigen pair, but which one is not guaranteed. 00:17:57.440 --> 00:18:01.880 So this is a class of techniques, iteration methods 00:18:01.880 --> 00:18:03.560 that we are using. 00:18:03.560 --> 00:18:06.610 Another class of techniques, again iteration methods, we 00:18:06.610 --> 00:18:09.890 are stuck with iteration, as I mentioned earlier, is the 00:18:09.890 --> 00:18:11.320 polynomial iteration methods. 00:18:11.320 --> 00:18:14.630 Here, we are saying that this problem is equivalent to that, 00:18:14.630 --> 00:18:16.920 we want to solve for these 0's. 00:18:16.920 --> 00:18:20.880 Therefore, for the 0's of this polynomial, sketched out here, 00:18:20.880 --> 00:18:23.710 and one way is to use Newton iteration. 00:18:23.710 --> 00:18:25.830 Here we have p, here we have p prime. 00:18:25.830 --> 00:18:29.740 Well, of course, if we have p here, we have to also 00:18:29.740 --> 00:18:32.110 recognize we should have the coefficients here. 00:18:32.110 --> 00:18:37.330 Well, these coefficients are given here, A0 to a n. 00:18:37.330 --> 00:18:42.200 And we could rewrite this polynomial in this form, once 00:18:42.200 --> 00:18:44.260 we know the zeroes. 00:18:44.260 --> 00:18:48.420 Now, if we don't know the zeroes, of course, we would be 00:18:48.420 --> 00:18:52.840 stuck on iterating on this matrix, where we would have to 00:18:52.840 --> 00:18:55.580 calculate these coefficients. 00:18:55.580 --> 00:19:00.040 However, there are now really, two different techniques. 00:19:00.040 --> 00:19:03.680 The first technique is explicit polynomial iteration. 00:19:03.680 --> 00:19:07.260 In the explicit polynomial iteration, we expand the 00:19:07.260 --> 00:19:10.710 polynomial and iterate for the 0's, by that I mean, we 00:19:10.710 --> 00:19:14.200 actually expand this polynomial. 00:19:14.200 --> 00:19:18.090 However, this technique is not suitable for larger problems 00:19:18.090 --> 00:19:21.150 because first of all, there's much work to obtain the ai's, 00:19:21.150 --> 00:19:23.720 and then it's an unstable process. 00:19:23.720 --> 00:19:27.940 Unstable because small errors in these coefficients here, 00:19:27.940 --> 00:19:30.350 give large errors in the zeroes. 00:19:30.350 --> 00:19:34.980 So this is not a suitable technique, and we are better 00:19:34.980 --> 00:19:37.720 off using implicit polynomial iteration. 00:19:37.720 --> 00:19:40.410 In the implicit polynomial iteration, we evaluate the 00:19:40.410 --> 00:19:44.280 determinant all this matrix here. 00:19:44.280 --> 00:19:50.290 In others, p at mu i, via the L D L transpose factorization, 00:19:50.290 --> 00:19:55.770 which we discussed in the Gauss elimination solution, 00:19:55.770 --> 00:19:57.920 when we talked about Gauss elimination solution. 00:19:57.920 --> 00:20:01.530 To obtain the determinant of this matrix here then, we need 00:20:01.530 --> 00:20:04.410 simply to multiply all the diagonal elements. 00:20:04.410 --> 00:20:08.450 This technique is accurate, provided we do not encounter 00:20:08.450 --> 00:20:12.030 large multipliers L matrix. 00:20:12.030 --> 00:20:14.400 And if we do, we have to shift away, we have to 00:20:14.400 --> 00:20:16.220 change our mu i. 00:20:16.220 --> 00:20:19.570 Also we directly solve for lambda 1 and so on. 00:20:19.570 --> 00:20:22.760 And for the more we can use effectively a secant 00:20:22.760 --> 00:20:25.710 iteration, instead of calculating, in other words, 00:20:25.710 --> 00:20:28.940 the p prime at mu. 00:20:28.940 --> 00:20:33.930 This is here approximately equal to p prime, at mu i. 00:20:33.930 --> 00:20:38.010 Instead of calculating this derivative, we are putting 00:20:38.010 --> 00:20:41.330 this value here, which gives us an approximation to that 00:20:41.330 --> 00:20:44.200 derivative, this is what we call the secant iteration. 00:20:44.200 --> 00:20:47.260 Once we have calculated one eigenvalue, of course, we will 00:20:47.260 --> 00:20:51.140 have to deflate to calculate the next eigenvalue, and this 00:20:51.140 --> 00:20:53.360 is the process that he would pursue. 00:20:53.360 --> 00:20:57.480 Pictorially, we have here, our p lambda, this is the 00:20:57.480 --> 00:20:58.470 polynomial. 00:20:58.470 --> 00:20:59.860 Let me show it to you here. 00:20:59.860 --> 00:21:03.240 Here, we have the polynomial, here we have lambda 1, here we 00:21:03.240 --> 00:21:04.260 have lambda 2. 00:21:04.260 --> 00:21:07.560 In the secant iteration, we start off with the value of mu 00:21:07.560 --> 00:21:11.500 i minus 1 here, mu i, these two values are given. 00:21:11.500 --> 00:21:15.930 Usually, in a practical analysis, mu 0 would be this 00:21:15.930 --> 00:21:20.660 value here, and mu 1 would be a lower bound on lambda 1. 00:21:20.660 --> 00:21:24.990 And plugging these two values in, we get mu i plus 1, and 00:21:24.990 --> 00:21:30.620 you can see that this process would then converge to this 00:21:30.620 --> 00:21:32.210 eigenvalue. 00:21:32.210 --> 00:21:35.880 Once we have lambda 1, we can repeat the iteration. 00:21:35.880 --> 00:21:39.590 However, working now on p lambda divided by lambda 00:21:39.590 --> 00:21:41.110 minus lambda 1. 00:21:41.110 --> 00:21:44.480 And that means we are operating on this polynomial 00:21:44.480 --> 00:21:48.640 here, and polynomial then gives us the iteration on that 00:21:48.640 --> 00:21:51.190 polynomial then gives us the second eigenvalue. 00:21:51.190 --> 00:21:57.380 So we are deflating the actual p lambda by the lambda 1 value 00:21:57.380 --> 00:22:00.070 in order to obtain a convergence 00:22:00.070 --> 00:22:02.940 to the second value. 00:22:02.940 --> 00:22:06.420 Well, this is the one class of methods. 00:22:06.420 --> 00:22:09.750 Another class of methods, the third one, are the Sturm 00:22:09.750 --> 00:22:10.930 sequence methods. 00:22:10.930 --> 00:22:15.060 And here, the basic idea is the following-- 00:22:15.060 --> 00:22:17.980 if we look at this eigenproblem, and I use a 4x4 00:22:17.980 --> 00:22:24.230 matrix here, then if we take a shift on this eigenproblem, mu 00:22:24.230 --> 00:22:27.100 S, and this is the shift here. 00:22:27.100 --> 00:22:33.700 Then if we factorize K minus mu S M, equal to L D L 00:22:33.700 --> 00:22:36.540 transpose, again, we're using our Gauss elimination 00:22:36.540 --> 00:22:38.120 procedures basically. 00:22:38.120 --> 00:22:42.420 Then the number of negative elements in D is equal to the 00:22:42.420 --> 00:22:45.110 number eigenvalues smaller than mu S, a 00:22:45.110 --> 00:22:47.480 very important point. 00:22:47.480 --> 00:22:51.100 Let us look at this example here of the 4x4 matrix. 00:22:51.100 --> 00:22:56.370 Here we have plotted characteristic polynomials. 00:22:56.370 --> 00:23:00.290 The characteristic polynomial of this 4x4 matrix, with that 00:23:00.290 --> 00:23:05.310 M matrix, of course, is drawn along here, four eigenvalues. 00:23:05.310 --> 00:23:09.150 And I have put mu S here as the red line. 00:23:09.150 --> 00:23:15.430 What this theorem says is if I factorize K minus mu S, M, 00:23:15.430 --> 00:23:19.160 this being the K matrix, the total K matrix, and the 4x4 M 00:23:19.160 --> 00:23:21.020 matrix goes in here. 00:23:21.020 --> 00:23:25.740 Then there would be one negative element in D. One 00:23:25.740 --> 00:23:30.900 only because mu S is large than lambda 1. 00:23:30.900 --> 00:23:34.470 If mu S were on this side here, there would be two 00:23:34.470 --> 00:23:38.640 negative elements in D. If mu S were on this side, there 00:23:38.640 --> 00:23:41.260 would be three negative elements in D-- 00:23:41.260 --> 00:23:46.360 very important point because qw can use this procedure to 00:23:46.360 --> 00:23:50.670 directly calculate eigenvalues the way I will be discussing 00:23:50.670 --> 00:23:51.430 it just now. 00:23:51.430 --> 00:23:55.470 However, let us look a little bit more at why 00:23:55.470 --> 00:23:56.730 this hole is here. 00:23:56.730 --> 00:24:01.940 Well, it derives from the associated constraint problems 00:24:01.940 --> 00:24:04.720 with the original eigenproblem. 00:24:04.720 --> 00:24:09.440 The original eigenproblem corresponds to this beam 00:24:09.440 --> 00:24:14.030 problem, which we considered earlier already in the 00:24:14.030 --> 00:24:15.730 discussion of Gauss elimination. 00:24:15.730 --> 00:24:19.870 If I constrain this degree of freedom here, I get a new 00:24:19.870 --> 00:24:22.690 polynomial corresponding to this problem, and the 00:24:22.690 --> 00:24:26.860 important point is that the lowest root, the lowest value 00:24:26.860 --> 00:24:32.000 of this problem here lies in between these two-- lambda 1 00:24:32.000 --> 00:24:33.170 and lambda 2. 00:24:33.170 --> 00:24:37.180 The second one here, lies between lambda 2 and lambda 3, 00:24:37.180 --> 00:24:41.350 and the third one here lies between lambda 3 and lambda 4. 00:24:41.350 --> 00:24:44.340 The second associated constraint problem has two 00:24:44.340 --> 00:24:48.000 eigenvalues, and this one lies between these two, this one 00:24:48.000 --> 00:24:49.370 lies between these two. 00:24:49.370 --> 00:24:51.850 And the final eigenvalue problem here, the third 00:24:51.850 --> 00:24:54.630 associate constraint problem obtained by knocking out these 00:24:54.630 --> 00:24:57.430 degrees of freedom here, has one eigenvalue, and that one 00:24:57.430 --> 00:25:01.790 lies in between these two values here. 00:25:01.790 --> 00:25:09.970 It's this important fact that is used to 00:25:09.970 --> 00:25:13.630 derive this final result. 00:25:13.630 --> 00:25:18.820 Now, we have not time to go through the actual details, 00:25:18.820 --> 00:25:22.220 but this is the important fact here. 00:25:22.220 --> 00:25:27.470 In other words, that the [? i's ?] eigenvalue of this 00:25:27.470 --> 00:25:31.430 associated constraint problem, for example, the second 00:25:31.430 --> 00:25:35.420 eigenvalue of this associated constraint problem here, 00:25:35.420 --> 00:25:38.630 bisects or lies in between, not bisect-- 00:25:38.630 --> 00:25:44.340 lies in between lambda 2 and lambda 3 of this eigenproblem 00:25:44.340 --> 00:25:45.930 here, and so on. 00:25:45.930 --> 00:25:50.700 That result is used to derive this final result, which we 00:25:50.700 --> 00:25:57.180 use effectively then in the eigen solution. 00:25:57.180 --> 00:26:03.950 And the class of techniques then that we are referring to 00:26:03.950 --> 00:26:07.470 using this fact are the Sturm sequence methods. 00:26:07.470 --> 00:26:11.050 So basically, what we are saying is the following-- 00:26:11.050 --> 00:26:15.800 take this coefficient matrix, now I make the mu x available. 00:26:15.800 --> 00:26:20.060 We factorize that coefficient matrix into L D L transposed, 00:26:20.060 --> 00:26:24.740 and the number of negative elements in D is equal to the 00:26:24.740 --> 00:26:29.630 number of eigenvalues smaller than mu S, I. How 00:26:29.630 --> 00:26:30.560 do we use this fact? 00:26:30.560 --> 00:26:34.885 Well, if this is here, the polynomial, we could start off 00:26:34.885 --> 00:26:40.950 by putting mu S 1 into this point here, and we would find 00:26:40.950 --> 00:26:44.510 out the number of eigenvalues smaller than mu S 1. 00:26:44.510 --> 00:26:49.170 Then we are taking, mu S 2, and we find out the number of 00:26:49.170 --> 00:26:51.320 eigenvalues smaller than mu S 2. 00:26:51.320 --> 00:26:58.580 Of course, we would find that they are one more eigenvalue 00:26:58.580 --> 00:26:59.860 smaller than mu S 2. 00:26:59.860 --> 00:27:02.460 Then mu S 1, namely this one. 00:27:02.460 --> 00:27:06.150 So we know there's one eigenvalue in between mu S 1 00:27:06.150 --> 00:27:07.910 and mu S 2. 00:27:07.910 --> 00:27:11.090 Now we could bisect this interval, and we would find 00:27:11.090 --> 00:27:16.080 that, in fact, there is still one more eigenvalue smaller 00:27:16.080 --> 00:27:20.970 than mu S 3, then there was smaller than mu S 1. 00:27:20.970 --> 00:27:23.750 Because this eigenvalue still lies on the left 00:27:23.750 --> 00:27:25.780 side of mu S 3. 00:27:25.780 --> 00:27:28.650 So we bisect again, we get mu S 4. 00:27:28.650 --> 00:27:33.550 Notice this distance here, is equal to that distance here. 00:27:33.550 --> 00:27:36.970 And now, we would find that the eigenvalue still lies in 00:27:36.970 --> 00:27:40.780 between these two lines, between mu S 3 and mu S 4. 00:27:40.780 --> 00:27:43.560 And like that. we could continue bisecting the 00:27:43.560 --> 00:27:48.400 interval until we finally find a very small interval in which 00:27:48.400 --> 00:27:51.550 this eigenvalue must lie. 00:27:51.550 --> 00:27:54.490 We, of course, have to take care in the L D L transpose 00:27:54.490 --> 00:27:57.920 factorization, because we are not working on a positive 00:27:57.920 --> 00:27:58.760 definite matrix. 00:27:58.760 --> 00:28:03.800 The multipliers must be small, the L, I, J elements must 00:28:03.800 --> 00:28:06.280 remain small in the factorization. 00:28:06.280 --> 00:28:09.700 But provided that is the case, the round of errors are small, 00:28:09.700 --> 00:28:11.930 and the procedure is very stable. 00:28:11.930 --> 00:28:15.100 However, convergence can be very slow, can be extremely 00:28:15.100 --> 00:28:17.990 slow because you have to cut down the intervals further and 00:28:17.990 --> 00:28:20.990 further until you finally get a very small interval, in 00:28:20.990 --> 00:28:23.560 which you now know the eigenvalue does lie. 00:28:23.560 --> 00:28:26.340 Of course, at some point, you might want to switch to 00:28:26.340 --> 00:28:27.810 another iteration scheme. 00:28:27.810 --> 00:28:31.210 Once you know there is an eigenvalue in between say mu S 00:28:31.210 --> 00:28:34.100 3 and mu S 4, you might want to switch to a iteration 00:28:34.100 --> 00:28:38.590 scheme that converges to this particular value much faster. 00:28:38.590 --> 00:28:41.570 Finally, a whole class of methods that I would like to 00:28:41.570 --> 00:28:45.640 briefly mention to you are the methods of transformations, 00:28:45.640 --> 00:28:49.230 where we are operating on this eigenvalue problem, 00:28:49.230 --> 00:28:53.860 recognizing that phi transpose K phi is equal to lambda. 00:28:53.860 --> 00:28:59.170 And phi transpose M phi is equal to I. Where phi is 00:28:59.170 --> 00:29:02.660 defined as shown here, lambda is defined as shown here. 00:29:02.660 --> 00:29:05.960 This fact, of course, we use already earlier in the 00:29:05.960 --> 00:29:06.360 [? multiple ?] 00:29:06.360 --> 00:29:07.750 position analysis. 00:29:07.750 --> 00:29:12.020 Now how can we get phi and lambda? 00:29:12.020 --> 00:29:13.430 That is our objective. 00:29:13.430 --> 00:29:17.900 Well, if we know that this matrix diagonalizes the K 00:29:17.900 --> 00:29:21.340 matrix, and it also diagonalizes the M matrix, 00:29:21.340 --> 00:29:25.370 then why not try to construct it by iteration? 00:29:25.370 --> 00:29:29.520 And this is the basic idea in the transformation methods. 00:29:29.520 --> 00:29:33.410 What we do is the following-- we take the original K matrix, 00:29:33.410 --> 00:29:37.720 we are operating with P1 transpose K P1. 00:29:37.720 --> 00:29:45.060 On K, and on M, such as to make this resultant matrix 00:29:45.060 --> 00:29:50.020 here closer to a diagonal form. 00:29:50.020 --> 00:29:52.470 Ultimately, we want to get a diagonal matrix. 00:29:52.470 --> 00:29:56.320 Well, one step in the iteration is to make this P1 00:29:56.320 --> 00:30:00.610 transpose K P1 matrix, and P1 transpose M P1 matrix, both of 00:30:00.610 --> 00:30:03.550 these matrices a little closer to diagonal form. 00:30:03.550 --> 00:30:06.030 Then K and M were by themselves. 00:30:06.030 --> 00:30:08.260 So this is the first step. 00:30:08.260 --> 00:30:12.360 Now, we have here, let's call this a K2 matrix, and this 00:30:12.360 --> 00:30:15.870 here an M2 matrix. 00:30:18.760 --> 00:30:23.240 Now, we are operating with P2 transpose on K2, and also 00:30:23.240 --> 00:30:24.400 there's a P2 here. 00:30:24.400 --> 00:30:28.250 So let me take a different color to show that now we are 00:30:28.250 --> 00:30:34.020 operating with this part here, and that part there. 00:30:38.850 --> 00:30:41.560 And this, of course, we are doing in such a way as to have 00:30:41.560 --> 00:30:47.730 that P2 transposed, K2, P2 is a little bit closer to 00:30:47.730 --> 00:30:52.440 diagonal form than was K2. 00:30:52.440 --> 00:30:55.320 I hope you can see these lines-- let me put 00:30:55.320 --> 00:30:57.570 them in once more. 00:30:57.570 --> 00:31:01.720 What I'm saying is we have K2 after having applied P1 00:31:01.720 --> 00:31:07.200 transpose on P1, on K, and similar for M. 00:31:07.200 --> 00:31:11.010 And now we want to make K2 still a little bit closer to 00:31:11.010 --> 00:31:12.380 diagonal form. 00:31:12.380 --> 00:31:17.780 Well, we are now applying P2 transposed on K2 and the P2 on 00:31:17.780 --> 00:31:21.460 K2, this one here has been replaced by K2, this one here 00:31:21.460 --> 00:31:23.240 has been replaced by M2. 00:31:23.240 --> 00:31:26.880 And what we're saying is that this final result, in this 00:31:26.880 --> 00:31:30.760 final result should be closer to diagonal form. 00:31:30.760 --> 00:31:33.510 In fact, this should be closer to I matrix, and this should 00:31:33.510 --> 00:31:36.790 be closer to the lambda matrix, ideally. 00:31:36.790 --> 00:31:41.280 Well, this is the basic process that is used in the 00:31:41.280 --> 00:31:43.750 generalized Jacobi method. 00:31:43.750 --> 00:31:46.750 We can show that if we choose specific 00:31:46.750 --> 00:31:50.770 matrices, P I matrices. 00:31:50.770 --> 00:31:55.220 That zero always one [? of ?] diagonal element. 00:31:55.220 --> 00:31:58.970 In other words, the P I matrix, zeroes a particular 00:31:58.970 --> 00:32:03.720 off diagonal element in this current coefficient matrix, 00:32:03.720 --> 00:32:09.320 similarly for the M2, for the current M coefficient matrix. 00:32:09.320 --> 00:32:13.920 Then we know in the generalized Jacobi method that 00:32:13.920 --> 00:32:15.795 the process finally converges. 00:32:15.795 --> 00:32:21.070 It converges and we get a diagonal matrix here. 00:32:21.070 --> 00:32:23.860 With proper scaling, it will be the identity matrix. 00:32:23.860 --> 00:32:27.670 And here, we are getting also a diagonal matrix, and if this 00:32:27.670 --> 00:32:31.590 is an identity matrix, this will actually be the matrix 00:32:31.590 --> 00:32:33.330 storing the eigenvalues. 00:32:33.330 --> 00:32:36.310 Well, this is the generalized Jacobi method here. 00:32:36.310 --> 00:32:39.270 They are other techniques that can also be used, of course, 00:32:39.270 --> 00:32:44.990 this is a very common approach taken when we have M being 00:32:44.990 --> 00:32:47.460 equal to the identity matrix. 00:32:47.460 --> 00:32:50.270 However, the disadvantage of this approach is that we are 00:32:50.270 --> 00:32:53.060 calculating all eigen pairs simultaneously. 00:32:53.060 --> 00:32:55.570 In finite element analysis, we are only interested in 00:32:55.570 --> 00:33:00.240 specific eigenvalues, so lowest P, or a certain number 00:33:00.240 --> 00:33:01.880 of eigenvalues in an interval. 00:33:01.880 --> 00:33:05.810 Here, we have to calculate all of the eigenvalues stored in 00:33:05.810 --> 00:33:08.790 the lambda matrix, and all of the eigenvectors. 00:33:08.790 --> 00:33:14.676 And the eigenvectors here, would be these vectors here. 00:33:14.676 --> 00:33:21.692 These are the the matrix phi, which is P1 times P2, up to P 00:33:21.692 --> 00:33:24.390 K, after K steps of iteration. 00:33:24.390 --> 00:33:28.295 This gives us the matrix phi, storing all of the 00:33:28.295 --> 00:33:30.270 eigenvectors. 00:33:30.270 --> 00:33:33.860 So the disadvantage is that we are calculating all of that 00:33:33.860 --> 00:33:37.720 eigenvectors, which in finite element analysis is hardly 00:33:37.720 --> 00:33:38.770 ever required. 00:33:38.770 --> 00:33:43.810 All we want is the total number of eigenvalues and 00:33:43.810 --> 00:33:44.760 eigenvectors. 00:33:44.760 --> 00:33:47.740 Remember, also if K is of [? order ?] 00:33:47.740 --> 00:33:51.310 1,000, this would be a tremendous amount of work 00:33:51.310 --> 00:33:54.500 involved and, in fact, it would be beyond the current 00:33:54.500 --> 00:33:57.060 state of computing capabilities to use a 00:33:57.060 --> 00:34:02.700 generalized Jacobi method on a K matrix 1,000x1,000 and an M 00:34:02.700 --> 00:34:04.420 matrix 1,000x1,000. 00:34:04.420 --> 00:34:08.320 You will certainly stretch the current computing capabilities 00:34:08.320 --> 00:34:12.770 tremendously by using the generalized Jacobi method on 00:34:12.770 --> 00:34:15.080 large eigenproblems. 00:34:15.080 --> 00:34:20.199 For large eigenproblems, it is really best to combine the 00:34:20.199 --> 00:34:22.560 basic approaches that I just mentioned. 00:34:22.560 --> 00:34:25.070 Let's recall the basic approaches. 00:34:25.070 --> 00:34:28.920 There were iteration methods, vector iteration methods, 00:34:28.920 --> 00:34:33.040 there were Sturm sequence methods, there was the 00:34:33.040 --> 00:34:35.580 transformation method. 00:34:35.580 --> 00:34:38.050 And then, of course, there were 00:34:38.050 --> 00:34:39.590 polynomial iteration methods. 00:34:39.590 --> 00:34:41.980 So really, four classes of methods that 00:34:41.980 --> 00:34:43.500 we are talking about. 00:34:43.500 --> 00:34:46.060 Each of them, we find, have some advantages. 00:34:46.060 --> 00:34:50.889 Well, as engineer, surely we want to now combine the 00:34:50.889 --> 00:34:54.350 advantages of each of these methods to come up with 00:34:54.350 --> 00:34:57.430 optimum techniques for our specific problems. 00:34:57.430 --> 00:35:01.410 Well, there are a number of thoughts 00:35:01.410 --> 00:35:02.380 that we would go through. 00:35:02.380 --> 00:35:05.600 First of all, a determinant search, a polynomial iteration 00:35:05.600 --> 00:35:09.640 method would be effective to get near a root. 00:35:09.640 --> 00:35:12.930 Vector iteration then could be used to calculate an 00:35:12.930 --> 00:35:14.880 eigenvalue and an eigenvector. 00:35:14.880 --> 00:35:18.180 The transformation method would be good to orthogonalize 00:35:18.180 --> 00:35:22.060 the current iteration vectors if we iterate with more than 00:35:22.060 --> 00:35:23.280 just one vector. 00:35:23.280 --> 00:35:26.870 And the Sturm sequence method is extremely powerful to 00:35:26.870 --> 00:35:30.360 assure that we have actually calculated the required 00:35:30.360 --> 00:35:33.520 eigenvalues and eigenvectors. 00:35:33.520 --> 00:35:36.480 Well, the first method then that I would like to discuss 00:35:36.480 --> 00:35:38.580 with you is the determinant search method. 00:35:38.580 --> 00:35:42.370 In this method, we have combined some of these 00:35:42.370 --> 00:35:44.510 techniques that I pointed out to you earlier 00:35:44.510 --> 00:35:46.340 in an optimum way. 00:35:46.340 --> 00:35:49.760 Here, we have the polynomial P lambda, we want to, of course, 00:35:49.760 --> 00:35:51.770 solve for the roots of that polynomial-- 00:35:51.770 --> 00:35:53.660 that is our objective. 00:35:53.660 --> 00:35:55.690 Here is lambda 1, there is lambda 2. 00:35:55.690 --> 00:35:58.770 I also show a point, A, because I will 00:35:58.770 --> 00:36:00.750 refer to it just now. 00:36:00.750 --> 00:36:04.200 The basic scheme in the determinant search method is 00:36:04.200 --> 00:36:09.140 to look for the zeroes of the determinant, of the 00:36:09.140 --> 00:36:14.260 determinant K minus mu I M. Of course, these zeroes give us R 00:36:14.260 --> 00:36:15.590 equal to the eigenvalues. 00:36:15.590 --> 00:36:17.560 That's what we're interested in calculating. 00:36:17.560 --> 00:36:20.740 Well, we are calculating the determinant by the 00:36:20.740 --> 00:36:24.970 factorization of this matrix into L D L transpose, and we 00:36:24.970 --> 00:36:27.990 notice that this determinant is simply the product of the 00:36:27.990 --> 00:36:30.230 [? D I I's, ?] as I mentioned earlier. 00:36:30.230 --> 00:36:36.560 Now, having a certain starting value, mu 0, say, and mu 1. 00:36:36.560 --> 00:36:40.170 We can directly use the secant iteration, that I mentioned to 00:36:40.170 --> 00:36:43.700 you earlier to get close to an eigenvalue. 00:36:43.700 --> 00:36:47.380 Now, I also slipped in here an acceleration factor. 00:36:47.380 --> 00:36:52.010 This acceleration factor would be equal to 1 in the normal 00:36:52.010 --> 00:36:55.320 secant iteration that I showed to you earlier already. 00:36:55.320 --> 00:36:58.540 However, since we want to use this polynomial iteration 00:36:58.540 --> 00:37:04.420 technique only to get close to an eigenvalue, and we might 00:37:04.420 --> 00:37:07.260 actually obtain with eta equal to 1, sometimes slow 00:37:07.260 --> 00:37:08.340 convergence. 00:37:08.340 --> 00:37:11.810 What we do is we simply set it larger to 00:37:11.810 --> 00:37:14.120 accelerate the process. 00:37:14.120 --> 00:37:16.780 And this is what we're doing here. 00:37:16.780 --> 00:37:22.600 If we do jump beyond an eigenvalue because we have 00:37:22.600 --> 00:37:24.840 accelerated the secant iteration. 00:37:24.840 --> 00:37:29.130 In other words, if we do a jump across here, then since 00:37:29.130 --> 00:37:35.510 we are factorizing the matrix here, we can also look at the 00:37:35.510 --> 00:37:38.540 number of negative elements in the D matrix. 00:37:38.540 --> 00:37:43.140 These would tell us that we, in fact, have jumped across an 00:37:43.140 --> 00:37:46.250 unknown eigenvalue. 00:37:46.250 --> 00:37:49.510 So, in this particular case, of course, there would be one 00:37:49.510 --> 00:37:51.280 negative element in the D matrix. 00:37:51.280 --> 00:37:54.040 And we would know now we have jumped across. 00:37:54.040 --> 00:37:58.480 The jumping across an unknown eigenvalue is not possible 00:37:58.480 --> 00:38:01.370 when eta is equal to 1. 00:38:01.370 --> 00:38:05.610 However, it is possible when eta is greater than 1. 00:38:05.610 --> 00:38:11.600 And however, that jumping does not hurt us. 00:38:11.600 --> 00:38:19.540 In fact, all we will find is that the Sturm sequence count 00:38:19.540 --> 00:38:22.950 will tell us that we have jumped across an eigenvalue, 00:38:22.950 --> 00:38:27.000 and then, of course, we have to go to another strategy of 00:38:27.000 --> 00:38:29.280 solving for the eigenvalue more accurately. 00:38:29.280 --> 00:38:33.110 But the important point is that all we want to do via 00:38:33.110 --> 00:38:36.420 this procedure is to get into the vicinity of the 00:38:36.420 --> 00:38:37.190 eigenvalue. 00:38:37.190 --> 00:38:40.780 We do not want to calculate it accurately, all we want to do 00:38:40.780 --> 00:38:44.290 is get into the vicinity of an unknown eigenvalue. 00:38:44.290 --> 00:38:49.080 First, of course, lambda 1, then lambda 2, and so on. 00:38:49.080 --> 00:38:53.380 There's also a concern, of course, that we do not want to 00:38:53.380 --> 00:38:58.310 jump beyond A, because if we, in our process of looking for 00:38:58.310 --> 00:39:05.090 lambda 1, we are jumping too close to lambda 2, then via 00:39:05.090 --> 00:39:07.400 the other solution techniques that we are using later on 00:39:07.400 --> 00:39:10.050 then to calculate the eigenvalue more accurately, we 00:39:10.050 --> 00:39:11.430 would actually converge to lambda 2. 00:39:11.430 --> 00:39:13.200 And, of course, that we don't want, we want to 00:39:13.200 --> 00:39:14.450 calculate lambda 1. 00:39:16.940 --> 00:39:22.950 So, the process is then to start off with eta equal to 1, 00:39:22.950 --> 00:39:26.220 iterate, see how we are progressing towards an unknown 00:39:26.220 --> 00:39:27.160 eigenvalue. 00:39:27.160 --> 00:39:29.410 If that progressing is very slow, we are 00:39:29.410 --> 00:39:31.740 choosing eta larger. 00:39:31.740 --> 00:39:34.470 There are specific strategies that we're using to choose the 00:39:34.470 --> 00:39:37.730 right amount all eta. 00:39:37.730 --> 00:39:41.670 And until we find that we are close to an eigenvalue, either 00:39:41.670 --> 00:39:46.010 from the left or we have jumped across it already, but 00:39:46.010 --> 00:39:47.940 you should never jumped to far. 00:39:47.940 --> 00:39:50.480 And we are making sure that we never jump beyond that point 00:39:50.480 --> 00:39:54.930 A. and having jumped, we find that the Sturm sequence count 00:39:54.930 --> 00:39:57.680 tells us about it. 00:39:57.680 --> 00:40:03.310 And now, we are going on to another strategy, and that 00:40:03.310 --> 00:40:06.100 strategy is the vector iteration, vector inverse 00:40:06.100 --> 00:40:09.040 iteration method that I discussed with you earlier. 00:40:09.040 --> 00:40:12.300 What we are doing now, is we have a new coefficient matrix, 00:40:12.300 --> 00:40:13.180 that's the one. 00:40:13.180 --> 00:40:15.550 We are sitting here, say, and now I have assumed that we 00:40:15.550 --> 00:40:17.240 have, in fact, jumped. 00:40:17.240 --> 00:40:18.170 It's not necessary. 00:40:18.170 --> 00:40:20.980 We might have also approached from this side and come very 00:40:20.980 --> 00:40:27.010 close to lambda J. And now, we are using vector iteration to 00:40:27.010 --> 00:40:31.120 get this eigenvalue accurately. 00:40:31.120 --> 00:40:35.150 This coefficient here, or this value here 00:40:35.150 --> 00:40:37.830 corresponds to that length. 00:40:37.830 --> 00:40:42.700 And if you are adding this value here to this shift, we 00:40:42.700 --> 00:40:45.630 are getting an approximation to lambda J. And that, of 00:40:45.630 --> 00:40:46.990 course, is our objective. 00:40:46.990 --> 00:40:51.680 Once we have lambda J, we have calculated this eigenvalue-- 00:40:51.680 --> 00:40:56.960 the eigenvector comes out, in x k plus 1 as k increases, so 00:40:56.960 --> 00:40:59.750 we have calculated this eigen pair, and now, of course, we 00:40:59.750 --> 00:41:03.230 would use deflation method, just the way I showed it to 00:41:03.230 --> 00:41:06.200 your earlier, for the polynomial iteration method. 00:41:06.200 --> 00:41:11.010 To deflate the polynomial of this value, implicitly-- 00:41:11.010 --> 00:41:12.560 implicitly, that's important-- 00:41:12.560 --> 00:41:15.820 and we would continue our iteration process towards that 00:41:15.820 --> 00:41:18.640 eigenvalue next. 00:41:18.640 --> 00:41:20.970 In this process, we also want to make sure that the 00:41:20.970 --> 00:41:24.160 iteration vector is deflated of the previously 00:41:24.160 --> 00:41:27.110 eigenvectors, for example, Gram-Schmidt 00:41:27.110 --> 00:41:28.770 orthogonalization. 00:41:28.770 --> 00:41:32.800 If the convergence is slow in this particular iteration 00:41:32.800 --> 00:41:35.170 here, we're using also the really Rayleigh-quotient 00:41:35.170 --> 00:41:40.310 iteration to speed up iteration itself. 00:41:40.310 --> 00:41:44.620 The advantage of the method-- well, it calculates only the 00:41:44.620 --> 00:41:48.260 eigen pairs that are actually required., there's no prior 00:41:48.260 --> 00:41:51.100 transformation of the eigenprobelm necessary. 00:41:51.100 --> 00:41:53.970 The disadvantage of the method is that we have to perform 00:41:53.970 --> 00:41:56.100 many triangular factorizations. 00:41:56.100 --> 00:41:58.800 Therefore, the method is really only effective for 00:41:58.800 --> 00:42:00.760 small banded systems. 00:42:00.760 --> 00:42:04.110 If we want to calculate eigenvalues of large banded 00:42:04.110 --> 00:42:07.860 systems, we really need an algorithm with less 00:42:07.860 --> 00:42:11.320 factorizations, and more vector iterations. 00:42:11.320 --> 00:42:14.240 And that is the subspace iteration method. 00:42:14.240 --> 00:42:16.980 In fact, we have now extended the subspace iteration method 00:42:16.980 --> 00:42:20.590 to be also most effective for small banded systems. 00:42:20.590 --> 00:42:23.840 And I will just refer to it very briefly later. 00:42:23.840 --> 00:42:26.940 The basic steps of the subspace iteration method are 00:42:26.940 --> 00:42:28.380 the following. 00:42:28.380 --> 00:42:33.090 Here, we have the first equation used. 00:42:33.090 --> 00:42:34.670 We have on the right hand side, the 00:42:34.670 --> 00:42:35.950 current iteration vectors. 00:42:35.950 --> 00:42:38.290 Now, notice that in this particular case, we're 00:42:38.290 --> 00:42:42.760 iterating with q vectors, where q ip larger than p 00:42:42.760 --> 00:42:44.350 simultaneously. 00:42:44.350 --> 00:42:49.130 This gives us here, q right hand side vector. 00:42:49.130 --> 00:42:52.200 We are solving for q vectors here. 00:42:52.200 --> 00:42:54.930 This is an inverse iteration step on q vector 00:42:54.930 --> 00:42:56.360 simultaneously. 00:42:56.360 --> 00:43:00.860 Then we are calculating the projections of these k and m 00:43:00.860 --> 00:43:03.410 matrices on these vectors. 00:43:03.410 --> 00:43:07.750 Notice that this pat here, of course, is equal to 00:43:07.750 --> 00:43:08.780 this right hand side. 00:43:08.780 --> 00:43:10.850 So we would not calculate this here. 00:43:10.850 --> 00:43:13.290 We simply use the right hand side, which we have evaluated 00:43:13.290 --> 00:43:15.560 already and plug it right in there. 00:43:15.560 --> 00:43:19.390 But having calculated these two matrices, which are now of 00:43:19.390 --> 00:43:23.960 order q by q, we are solving eigenvalue problem all of 00:43:23.960 --> 00:43:26.110 these q by q matrices. 00:43:26.110 --> 00:43:27.780 And the solution of that eigenvalue 00:43:27.780 --> 00:43:29.230 problem is shown here. 00:43:29.230 --> 00:43:34.030 Notice this matrix here is a q by q matrix. 00:43:34.030 --> 00:43:36.020 We're using q vectors. 00:43:36.020 --> 00:43:40.970 This is a q by q matrix storing the eigenvectors of 00:43:40.970 --> 00:43:44.470 this eigen problem here, of the eigenproblem corresponding 00:43:44.470 --> 00:43:46.510 to these matrices. 00:43:46.510 --> 00:43:51.490 This here is the diagonal matrix listing the eigenvalues 00:43:51.490 --> 00:43:52.830 of the eigenproblem 00:43:52.830 --> 00:43:55.090 corresponding to these matrices. 00:43:55.090 --> 00:43:58.930 Having calculated these matrices here, and of course, 00:43:58.930 --> 00:44:02.070 that one too, we are performing this operation to 00:44:02.070 --> 00:44:05.390 get a new set of iteration vectors, which then are 00:44:05.390 --> 00:44:07.190 plugged back into here. 00:44:07.190 --> 00:44:12.580 Notice that the solution of this eigenproblem, q by q, 00:44:12.580 --> 00:44:17.660 might be a 50x50 matrix here, 100x100 matrix, maybe a 00:44:17.660 --> 00:44:20.720 200x200, but that already would be a large number of 00:44:20.720 --> 00:44:21.930 vectors if we're using-- 00:44:21.930 --> 00:44:23.730 a very large number of vectors. 00:44:23.730 --> 00:44:26.250 This solution of this eigenproblem can be 00:44:26.250 --> 00:44:29.290 effectively achieved via the generalized Jacobi method, 00:44:29.290 --> 00:44:32.270 which I referred to briefly earlier. 00:44:32.270 --> 00:44:35.515 Of course here, we are using our solution algorithm for 00:44:35.515 --> 00:44:36.290 static analysis. 00:44:36.290 --> 00:44:41.200 There's L D L transpose factorization involved here. 00:44:41.200 --> 00:44:45.260 And here, this is a simple matrix multiplication that has 00:44:45.260 --> 00:44:46.880 to be carried out. 00:44:46.880 --> 00:44:50.620 This part here, is the eigenproblem solution, and 00:44:50.620 --> 00:44:51.450 there's another vector 00:44:51.450 --> 00:44:53.980 multiplication carried out here. 00:44:53.980 --> 00:44:59.320 Now, under specific conditions, this eigensolution 00:44:59.320 --> 00:45:04.690 algorithm converges and the vectors stored in x k plus 1, 00:45:04.690 --> 00:45:08.940 they are q vectors now here, with appropriate ordering, 00:45:08.940 --> 00:45:11.750 converge to the phi matrix. 00:45:11.750 --> 00:45:17.800 This is now an n by q matrix, storing the lowest 00:45:17.800 --> 00:45:22.200 eigenvectorr, or I should say, the q eigenvectors 00:45:22.200 --> 00:45:25.020 corresponding to the lowest eigenvalues. 00:45:27.800 --> 00:45:30.820 And here, we have it written out, and of course, this is a 00:45:30.820 --> 00:45:31.960 diagonal matrix. 00:45:31.960 --> 00:45:35.530 The condition that we have to satisfy here is that the 00:45:35.530 --> 00:45:42.800 starting subspace spanned by these vectors here must not be 00:45:42.800 --> 00:45:52.840 orthogonal to the subspace spanned by the p eigenvectors 00:45:52.840 --> 00:45:54.680 that we are interested in. 00:45:54.680 --> 00:45:59.160 Of course, here, I'm talking about q eigenvectors. 00:45:59.160 --> 00:46:03.360 I really am only interested in p eigenvectors, in the lowest 00:46:03.360 --> 00:46:05.510 p eigenvectors. 00:46:05.510 --> 00:46:09.390 So in general, what happens is that I will only be interested 00:46:09.390 --> 00:46:13.080 in these eigenvectors here, and the lowest eigenvalues. 00:46:13.080 --> 00:46:16.850 One comes out here, I'm not too much concerned about. 00:46:16.850 --> 00:46:21.820 In fact, we are carrying these along only to accelerate the 00:46:21.820 --> 00:46:23.040 iteration scheme. 00:46:23.040 --> 00:46:28.060 And I will show you just know what the convergence rate is, 00:46:28.060 --> 00:46:29.770 and why we are carrying them along. 00:46:29.770 --> 00:46:34.800 So since we only interested in the lowest p eigenvalues and 00:46:34.800 --> 00:46:38.870 corresponding eigenvectors, we want to get convergence 00:46:38.870 --> 00:46:48.330 towards those, and what that means is that our starting 00:46:48.330 --> 00:46:57.290 vectors in here must not be m orthogonal to the eigenvectors 00:46:57.290 --> 00:47:00.780 that we want to calculate, phi 1 to phi p. 00:47:00.780 --> 00:47:02.960 If that condition is satisfied, one can show, 00:47:02.960 --> 00:47:06.300 theoretically, and of course, we see it in practice also, 00:47:06.300 --> 00:47:09.230 that this iteration scheme will always converge to the 00:47:09.230 --> 00:47:12.385 lowest eigenvectors and corresponding eigenvalues. 00:47:16.010 --> 00:47:20.760 One important point here, once we have calculated, say, a 00:47:20.760 --> 00:47:24.560 certain number of iteration vectors, and we have tested 00:47:24.560 --> 00:47:27.810 for convergence, we are saying that convergence is reached, 00:47:27.810 --> 00:47:32.670 say, when the current iterates on the eigenvalues don't 00:47:32.670 --> 00:47:34.450 change very much anymore. 00:47:34.450 --> 00:47:38.720 This tol might be 10 to the minus 6, as an example. 00:47:38.720 --> 00:47:41.820 Then once this condition is satisfied, we stop the 00:47:41.820 --> 00:47:46.440 iteration, and we really want to make then sure that we 00:47:46.440 --> 00:47:48.910 really converge to the eigenvalues of interest. 00:47:48.910 --> 00:47:50.580 And how is that accomplished? 00:47:50.580 --> 00:47:53.310 Well, if these are the p eigenvalues, we are 00:47:53.310 --> 00:47:57.060 calculating error bounds on these eigenvalues, indicated 00:47:57.060 --> 00:48:01.690 here by the blue lines here. 00:48:01.690 --> 00:48:06.670 And then just to the right of the error bound. 00:48:06.670 --> 00:48:11.300 Or we can really use that error bound, we apply a Sturm 00:48:11.300 --> 00:48:15.520 sequence shift, mu S. What this then says that k minus mu 00:48:15.520 --> 00:48:18.150 S factorize into L D L transpose. 00:48:23.090 --> 00:48:29.870 And we must now have in the D matrix, p negative elements. 00:48:29.870 --> 00:48:35.390 And only p, not more, not less, just p negative elements 00:48:35.390 --> 00:48:39.050 must be occurring in the D matrix. 00:48:39.050 --> 00:48:44.410 This check is absolutely sure. 00:48:44.410 --> 00:48:47.510 If this check a satisfied, we can be absolutely sure that we 00:48:47.510 --> 00:48:50.970 have calculated the eigenvalues of interest, the 00:48:50.970 --> 00:48:53.830 lowest p eigenvalues, in this particular case. 00:48:53.830 --> 00:48:57.910 So I repeat, after we have satisfied this convergence 00:48:57.910 --> 00:49:03.510 rate, or this conversions limit, I should say-- 00:49:03.510 --> 00:49:06.400 typically 10 to the minus 6, or 10 to the minus 8, 00:49:06.400 --> 00:49:08.700 depending on what accuracy you want. 00:49:08.700 --> 00:49:12.920 We calculating eigenvalue bounds, as shown here. 00:49:12.920 --> 00:49:18.500 We are putting a mu S Sturm sequence shift just to the 00:49:18.500 --> 00:49:24.460 right of the last eigenvalue bound, and if this is the p's 00:49:24.460 --> 00:49:28.540 eigenvalue, we must have here, exactly p 00:49:28.540 --> 00:49:30.670 negative elements indeed. 00:49:30.670 --> 00:49:35.390 The convergence rate of the iteration, we find is given by 00:49:35.390 --> 00:49:38.110 lambda i divided by lambda q plus 1. 00:49:40.750 --> 00:49:44.650 And the convergence rate to the eigenvalues of interest 00:49:44.650 --> 00:49:47.140 are just this value squared. 00:49:47.140 --> 00:49:50.280 Notice that we have here, q plus 1. 00:49:50.280 --> 00:49:53.840 So if we are interested in the p's eigenvalue,== the black 00:49:53.840 --> 00:49:56.850 pen doesn't work too well, let me take the blue one here-- 00:49:56.850 --> 00:50:01.330 if you're interested in the p's eigenvalue as the largest 00:50:01.330 --> 00:50:04.220 eigenvalue, notice that the convergence rate is that 00:50:04.220 --> 00:50:07.660 lambda p plus-- lambda p divided by lambda q plus 1. 00:50:07.660 --> 00:50:11.030 That is the worst convergence rate that we are reaching 00:50:11.030 --> 00:50:14.840 because the other eigenvalues below lambda p are smaller, 00:50:14.840 --> 00:50:16.740 giving us better convergence rate. 00:50:16.740 --> 00:50:22.110 So this is the worst that we can reach, if we are only 00:50:22.110 --> 00:50:25.280 interested in p eigenvalues and corresponding 00:50:25.280 --> 00:50:26.790 eigenvectors. 00:50:26.790 --> 00:50:29.850 Notice however, that lambda q plus 1 must be 00:50:29.850 --> 00:50:31.290 greater than lambda p. 00:50:31.290 --> 00:50:33.440 If lambda q plus 1 is equal to lambda p, 00:50:33.440 --> 00:50:34.880 we would never converge. 00:50:34.880 --> 00:50:38.930 So this is the reason why we are choosing q to be larger 00:50:38.930 --> 00:50:42.030 than p, substantially large than p. 00:50:42.030 --> 00:50:49.640 For example, we might choose q to be typically, as an 00:50:49.640 --> 00:50:57.850 example, q being the maximum of 2p and p plus 8. 00:50:57.850 --> 00:51:00.940 This is a formula that we have used with 00:51:00.940 --> 00:51:02.700 quite a lot of success. 00:51:02.700 --> 00:51:05.920 In other words, when p is equal to 4, we would use-- 00:51:11.060 --> 00:51:15.710 no, sorry, I should say the minimum here. 00:51:15.710 --> 00:51:19.100 the In other words, when p is equal to 4, we 00:51:19.100 --> 00:51:21.430 would have 8 here. 00:51:21.430 --> 00:51:24.410 When p is equal to 8, it's the same. 00:51:24.410 --> 00:51:26.390 Because these two values give us the same. 00:51:26.390 --> 00:51:30.470 And when p is equal to 100, then we just have 8 more. 00:51:30.470 --> 00:51:32.850 Q is equal to 108. 00:51:32.850 --> 00:51:37.870 This 8 here is a safeguard, so that this cannot be equal to 1 00:51:37.870 --> 00:51:39.130 in practice. 00:51:39.130 --> 00:51:43.960 It would only be equal to 1 if we would have, basically. 00:51:43.960 --> 00:51:47.390 9 time the same root, and of course, that is a case which 00:51:47.390 --> 00:51:50.290 we would hardly encounter in practice. 00:51:50.290 --> 00:51:55.060 Now, regarding the choice of the starting values. 00:51:55.060 --> 00:51:57.860 For the vectors, there are two choices. 00:51:57.860 --> 00:52:03.070 x1 is simply a vector with 1 only. 00:52:03.070 --> 00:52:10.220 xj would be a vector with zeroes everywhere, but with a 00:52:10.220 --> 00:52:11.300 1 somewhere. 00:52:11.300 --> 00:52:12.720 And that one is in the j's entry-- 00:52:15.590 --> 00:52:18.310 sorry, the k's entry. 00:52:18.310 --> 00:52:23.050 The k's entry for ek, ek has in the case, entry of 1. 00:52:23.050 --> 00:52:25.630 And otherwise, everywhere zeroes. 00:52:25.630 --> 00:52:29.390 This would be the second to the q minus first vector. 00:52:29.390 --> 00:52:31.150 And the last vector, we take a random vector. 00:52:31.150 --> 00:52:33.190 That is one choice of starting vectors. 00:52:33.190 --> 00:52:36.310 The Lanczos method can also be used to generate starting 00:52:36.310 --> 00:52:37.560 vectors very effectively. 00:52:41.260 --> 00:52:45.450 Once we have selected the starting vector, we go through 00:52:45.450 --> 00:52:48.770 the iteration schemes the way I have been displaying it. 00:52:48.770 --> 00:52:53.420 And there, it is important, once again, to perform after 00:52:53.420 --> 00:52:56.510 the convergence- the Sturm system sequence check, and 00:52:56.510 --> 00:53:00.610 also we might want to calculate error bounds on the 00:53:00.610 --> 00:53:01.970 eigenvalues and eigenvectors. 00:53:01.970 --> 00:53:05.240 These error bounds would be calculated by simply saying, 00:53:05.240 --> 00:53:08.070 well, if we want to satisfy-- 00:53:08.070 --> 00:53:11.660 or if we want to solve for lambdas and phis that 00:53:11.660 --> 00:53:17.400 satisfies this equation, then why not see how 00:53:17.400 --> 00:53:18.650 close we will satisfy. 00:53:18.650 --> 00:53:23.780 Make that a minus and put here an error, say, R, on 00:53:23.780 --> 00:53:25.050 the right hand side. 00:53:25.050 --> 00:53:28.260 Substituting then the last values for lambda and phi that 00:53:28.260 --> 00:53:32.220 we just kind calculated, the last iterative values, we get 00:53:32.220 --> 00:53:34.180 this top line here. 00:53:34.180 --> 00:53:36.500 We're taking a norm the way I've been 00:53:36.500 --> 00:53:38.310 talking about earlier. 00:53:38.310 --> 00:53:43.730 And we are also taking a norm at the bottom to get, 00:53:43.730 --> 00:53:46.690 basically, a relative value, top to bottom. 00:53:46.690 --> 00:53:49.240 And then relative value gives us epsilon i. 00:53:49.240 --> 00:53:54.070 So epsilon i should be of order, 10 to the minus 3, if 00:53:54.070 --> 00:53:56.740 we want to have three-digit accuracy in the eigenvalues 00:53:56.740 --> 00:53:57.770 and eigenvectors. 00:53:57.770 --> 00:54:00.020 In fact, the eigenvalues will generally be more accurate 00:54:00.020 --> 00:54:02.040 than the eigenvectors. 00:54:02.040 --> 00:54:04.130 That is, of course, a very low accuracy. 00:54:04.130 --> 00:54:06.260 Ideally, we have much more accuracy after 00:54:06.260 --> 00:54:07.330 the iteration schemes. 00:54:07.330 --> 00:54:11.330 We would have to tighten the tolerance, the tol value that 00:54:11.330 --> 00:54:15.010 I referred to earlier sufficiently to make epsilon i 00:54:15.010 --> 00:54:16.180 small enough. 00:54:16.180 --> 00:54:19.370 Now, I like to refer, very briefly also, to an 00:54:19.370 --> 00:54:21.740 accelerated subspace iteration method that we have been 00:54:21.740 --> 00:54:26.800 developing during the last years, which is a very 00:54:26.800 --> 00:54:30.340 substantial extension of this subspace iteration method. 00:54:30.340 --> 00:54:36.080 The idea here is, basically, that we wanted to combine the 00:54:36.080 --> 00:54:39.690 sum of the effective techniques that we are using 00:54:39.690 --> 00:54:44.115 in the polynomial iteration method, with the standard 00:54:44.115 --> 00:54:46.230 subspace iteration method that I have 00:54:46.230 --> 00:54:48.340 discussed with you here. 00:54:48.340 --> 00:54:51.710 Basically, in the accelerated subspace iteration method, we 00:54:51.710 --> 00:54:55.550 are not keeping a constant coefficient matrix, as 00:54:55.550 --> 00:54:56.870 we have done here. 00:54:56.870 --> 00:54:59.160 The k matrix never changed. 00:54:59.160 --> 00:55:02.640 Here, we are shifting the coefficient matrix, we are 00:55:02.640 --> 00:55:04.480 shifting through the spectrum. 00:55:04.480 --> 00:55:08.860 Much in the same way as we are pursuing it in the determined 00:55:08.860 --> 00:55:10.070 search method. 00:55:10.070 --> 00:55:11.450 And that is number one. 00:55:11.450 --> 00:55:15.720 And number two-- we are also not iterating with vectors 00:55:15.720 --> 00:55:18.610 that are always larger than the number of eigenvectors 00:55:18.610 --> 00:55:19.710 that we want to calculated. 00:55:19.710 --> 00:55:22.990 In fact, if we want to calculate, say, 100 00:55:22.990 --> 00:55:27.700 eigenvalues, we might just use four iteration vectors at a 00:55:27.700 --> 00:55:31.420 time, and shift through the spectrum, always with four 00:55:31.420 --> 00:55:35.020 iteration vectors, capturing the eigenvalues values and 00:55:35.020 --> 00:55:35.980 eigenvectors that we are really 00:55:35.980 --> 00:55:37.800 interested in, in bundles. 00:55:37.800 --> 00:55:42.080 So this is a way how we have very substantially accelerated 00:55:42.080 --> 00:55:45.310 the standard subspace iteration method that I have 00:55:45.310 --> 00:55:49.970 been talking to you about just now. 00:55:49.970 --> 00:55:52.315 If you are interested in seeing some solution times, 00:55:52.315 --> 00:55:56.180 please we look at this reference. 00:55:56.180 --> 00:55:59.960 Ladies and gentlemen, this then brings me to the 00:55:59.960 --> 00:56:01.640 conclusion of this lecture. 00:56:01.640 --> 00:56:05.630 And in fact, it also brings us to the conclusion of the whole 00:56:05.630 --> 00:56:09.300 set of lectures that I have been presenting to you. 00:56:09.300 --> 00:56:14.180 These set up lectures have been brief and very compact. 00:56:14.180 --> 00:56:16.900 They've been a survey introduction to the finite 00:56:16.900 --> 00:56:17.760 element method-- 00:56:17.760 --> 00:56:21.670 a very brief and compact introduction and survey. 00:56:21.670 --> 00:56:24.520 Brief because there are only 12 lectures. 00:56:24.520 --> 00:56:27.800 Compact because I wanted to present to you as much as 00:56:27.800 --> 00:56:30.570 possible in these 12 lectures. 00:56:30.570 --> 00:56:34.390 And there is, of course, a lot of information that we could 00:56:34.390 --> 00:56:35.900 have talked about. 00:56:35.900 --> 00:56:40.560 However, I hope that the information that I had given 00:56:40.560 --> 00:56:44.760 in the lectures has been of interest to you, and it will 00:56:44.760 --> 00:56:52.070 also be valuable to you in your work and general decision 00:56:52.070 --> 00:56:54.680 making process of whether you want to use a finite element 00:56:54.680 --> 00:56:58.120 method, how you want to use it, and so on. 00:56:58.120 --> 00:57:01.600 I like to take the opportunity at this point to thank the 00:57:01.600 --> 00:57:04.670 Center of Advanced Engineering Studies for the terrific 00:57:04.670 --> 00:57:07.260 effort they've made in preparing the tapes. 00:57:07.260 --> 00:57:12.150 I'm very happy that they have been cooperative with me in 00:57:12.150 --> 00:57:17.750 the way they did, in the preparation of the tapes. 00:57:17.750 --> 00:57:22.145 I also would like to thank you as listeners for your patience 00:57:22.145 --> 00:57:25.590 of listening to the tapes and, once again, I hope you found 00:57:25.590 --> 00:57:27.950 the tapes interesting and have obtained some value 00:57:27.950 --> 00:57:29.260 information on them. 00:57:29.260 --> 00:57:31.590 Best of luck in your further studies of the 00:57:31.590 --> 00:57:32.840 finite element method.