Jacobi And Gauss-seidel 632w6d
This document was ed by and they confirmed that they have the permission to share it. If you are author or own the copyright of this book, please report to us by using this report form. Report 3i3n4
Overview 26281t
& View Jacobi And Gauss-seidel as PDF for free.
More details 6y5l6z
- Words: 3,019
- Pages: 10
ITERATIVE METHODS FOR SOLVING LINEAR EQUATIONS There are other methods that can be used to solve a set of linear equations that are based on iteration. In these cases, an initial estimate of the parameters is estimated and then the equations are solved, yielding an updated version of the parameters. These new values are then inserted back into the equations and the process continues until the desired solution is reached. The two iterative methods discussed here are the Jacobi method and the Gauss-Seidel method.
Jacobi Iteration Method Given a set of linear equations, a 1 x1 + a 2 x 2 + L + a n x n = s1 b1 x1 + b 2 x 2 +L + bn x n = s 2 M d1 x1 + d 2 x 2 +L + d n x n = s n the problem is one of solving for x1 , x2 , …, xn . The right hand side of these equations, si, represents the solution. We begin by rearranging these equations in the form of solving for the unknown parameters one equation at a time. Thus, s 1 a 2 ( 0 ) a 3 (0 ) a − x2 − x 3 − L − n x (n0 ) a 1 a1 a1 a1 s b b b x 2 = 2 − 1 x 1(0 ) − 3 x 3(0 ) − L − n x (n0 ) b 12 b 2 b2 b2 M s d d d x n = n − 1 x 1(0 ) − 2 x (20 ) − L − n −1 x (n0−)1 an dn dn dn x1 =
The superscript (0) indicates the initial estimate of the parameters. For the first , these parameters are given the value zero. The equations are then solved which results in an updated value of the parameters. These current estimates are then inserted back into the equations and a newer set of parameters is arrived at by solving these equations. The process continues until the solution converges. As an example, take the following linear equations: 7x 1 + 3x 2 + x 3 = 18 2 x 1 − 9x 2 + 4 x 3 = 12 x 1 − 4 x 2 + 12 x 3 = 6
Rearrange these equations x1 =
18 3 1 − x2 − x3 7 7 7
⇒
x 1 = 2.571 − 0.429x 2 − 0.143x 3
12 2 4 + x1 + x3 9 9 9
⇒
x 2 = −1.333 + 0.222x 1 + 0.444 x 3
6 1 4 − x1 + x 2 12 12 12
⇒
x 3 = 0.500 − 0.083x 1 + 0.333x 2
x2 = −
x3 =
Use as the initial estimates: x1 (0) = x2 (0) = x3 (0) = 0. Insert these estimates into these equations yielding new estimates of the parameters. x 1(1) = 2.571 − 0.429 (0) − 0.143 (0) = 2.571
x (21) = −1.333 + 0.222 (0 ) + 0.444 (0 ) = −1.333 x (31 ) = 0.500 − 0.083 (0 ) + 0.333 (0) = 0.500
Insert these updated estimates back into original equation again, yielding x 1( 2 ) = 2.571 − 0.429 (− 1.333) − 0.143 (0.500 ) = 3.071
x (22 ) = −1.333 + 0.222 (2.571) + 0.444 (0.500) = −0.540 x (32 ) = 0.500 − 0.083 (2.571) + 0.333 (−1.333) = −0.159
This process is continued until the desired results are obtained. The following table shows the solutions arrived at after each iteration. These results are from the attached Fortran program. The output shows x(i) that are the parameters xi. Also output is the change in the parameters, dx(i), between each iteration. SOLVING LINEAR EQUATION USING THE JACOBI ITERATION METHOD The estimated results after each iteration are shown as: Iteration 1 2 3 4 5 6
x(1) 2.57143 3.07143 2.82540 2.87141 2.85811 2.85979
x(2) -1.33333 -.53968 -.72134 -.67695 -.68453 -.68261
x(3) .50000 -.15873 .06415 .02410 .03506 .03365
dx(1) 2.57143 .50000 -.24603 .04601 -.01330 .00168
dx(2) -1.33333 .79365 -.18166 .04439 -.00757 .00192
dx(3) .50000 -.65873 .22288 -.04005 .01096 -.00142
How good are these results? Lets take our equations and put them into an augmented matrix and solve using Gauss-Jordan elimination. Iterative Methods for Solving Linear Equations
Page 2
6 R + 31R 7 3 1 18 R1 ↔ R 3 1 − 4 12 6 R 2 − 2 R 1 1 − 4 12 2 2 − 9 4 12 2 − 9 4 12 0 −1 − 20 3 0 → R→ → −7 R −R2 1 − 4 12 6 7 3 1 18 3 1 0 31 − 83 − 24 12 6 R +4 R 1 0 92 6 R −92 R 1 0 0 2.859 1 − 4 1 2 0 1 1 3 0 1 0 − 0.683 20 0 0 1 20 0 R→ r→ − 20 R 3 0 0 − 703 − 24 3 703 0 0 1 0.034 2 0 0 1 0.034 As one can see, the values using the Jacobi iterative method are very close. Following is a Fortran program that can be used to use the Jacobi iteration to solve a set of equations. The limitation now is that it is restricted to only a 3 x 3 matrix, due to formatting procedures currently used in the program. c c c c
Program Jacobi Program to solve a linear equation using the Jacobi Iteration method IMPLICIT none REAL*8 coef(3,4), d, dx(3), x(3,4), xn(3), xnp(3) INTEGER i, iter, iterate,j, no, nv DATA iterate /0/
c c c c c c c c c
The data are entered into the program using a data file called jacobi.dat. It has the following row values number of equations number of variables x(1) x(2) x(3) solution for the first equation x(1) x(2) x(3) solution for the second equation x(1) x(2) x(3) solution for the third equation OPEN (4, file = 'jacobi.dat') OPEN (6, file = 'results')
c c c
5
c c c c
no is the number of equations and nv is the number of variables read(4,*) no do 5 i=1,no xn(i) = 0.d0 continue read(4,*) nv write(6,901) The coefficients for the variables are read in the matrix x with the solution to the equations being the last column do 10 i=1,no read(4,*)(x(i,j),j=1,no+1)
c c c c c
d is the coefficient for the variable that is being solved for it forms the denominator to compute the real number for the remaining coefficients d = x(i,i) do 7 j=1,no+1
Iterative Methods for Solving Linear Equations
Page 3
7 c c c c c c
10
13 15 c c c
coef(i,j) = x(i,j)/d end do Because the Jacobi method solves for the respect to the current estimates of the coefficient for the variable is made to use in the loop to compute the adjusted
unknown variable with other variables, the be zero for subsequent estimates
coef(i,i) = 0.d0 write(6,900)(x(i,j),j=1,nv+1) end do write(6,902) do 13 i=1,no write(6,900)(coef(i,j),j=1,nv+1) end do write(6,903) iter = 0 iterate is just a counter to keep track of the number of iterations iterate = iterate+1
c c c
18 20 c c c c
Solve for the estimate of the unknown parameters do 20 i=1,no xnp(i) = coef(i,nv+1) do 18 j=1,nv xnp(i) = xnp(i) - coef(i,j)*xn(j) end do end do dx is a vector showing the change in the estimate of the variable with respect to the estimated value used in the previous iteration do 50 i=1,no dx(i) = xnp(i) - xn(i)
c c c c c c
Test to see if the change is greater than the threshold If it is, then the variable iter is made equal to 1 At the beginning of each loop, this value is made equal to 0 If iter is 1 then this means to iterate again if (dabs(dx(i)).gt.0.01d0) iter = 1
c c c 50
900 901
902 903
904
Update the estimated parameter value xn(i) = xnp(i) continue write(6,904)iterate,(xn(i),i=1,nv),(dx(i),i=1,nv) if (iter.gt.0) go to 15 FORMAT(5x,4(f10.4,5x)) FORMAT(15x,'SOLVING LINEAR EQUATION USING THE JACOBI ITERATION MET 1HOD',//,'The coefficients to the equations with the solution at th 1e end are:',/) FORMAT(//,'Rearranging the equations to solve for the unknown vari 1ables yields',/,5x,'the following coefficients: ',/) FORMAT(//,'The estimated results after each iteration are shown as 1:'//,'Iteration',2x,'x(1)',9x,'x(2)',9x,'x(3)',7x,'dx(1)',7x,'dx(2 2)',7x,'dx(3)') FORMAT(i3,4x,6(f10.5,2x)) stop end
Iterative Methods for Solving Linear Equations
Page 4
Gauss-Seidel Iterative Method The Gauss-Seidel iterative method of solving for a set of linear equations can be thought of as just an extension of the Jacobi method. Start out using an initial value of zero for each of the parameters. Then, solve for x1 as in the Jacobi method. When solving for x2 , insert the just computed value for x1 . In other words, for each calculation, the most current estimate of the parameter value is used. To see how the Gauss-Seidel method works, lets use the values in the last example. The initial estimates are set to zero. Then, the results from the first iteration are shown as: x1(1) = 2.571 − 0.429 x (20 ) − 0.143 x (30 ) = 2.571 − 0.429 (0 ) − 0.143 (0 ) = 2.571
x (21) = −1.333 + 0.222 x (11 ) + 0.444 x (30 ) = −1.333 + 0.222 (2.571) + 0.444 (0 ) = −0.762
x (31 ) = 0.500 − 0.083 x1(1) + 0.333 x (21 ) = 0.500 − 0.083 (2.571) + 0.333 (− 0.762 ) = 0.033 The next iteration is performed in a similar fashion. It can be shown as: x (12 ) = 2.571 − 0.429x (21 ) − 0.143x (31 ) = 2.571 − 0.429 (− 0.762) − 0.143 (0.033) = 2.893
x (22 ) = −1.333 + 0.222 x (12 ) + 0.444 x (31) = −1.333 + 0.222 (2.893) + 0.444 (0.033) = −0.676 x (32 ) = 0.500 − 0.083x1( 2 ) + 0.333x (22 ) = 0.500 − 0.083 (2.893) + 0.333 (− 0.676 ) = 0.034
A Fortran program was written to solve this problem. The results are shown in the next table. SOLVING LINEAR EQUATION USING THE GAUSS-SEIDEL ITERATION METHOD The estimated results after each iteration are shown as: Iteration 1 2 3 4
x(1) 2.57143 2.89342 2.85646 2.85957
x(2) -.76190 -.67624 -.68369 -.68273
x(3) .03175 .03347 .03406 .03412
dx(1) -2.57143 -.32200 .03696 -.00311
dx(2) .76190 -.08566 .00745 -.00096
dx(3) -.03175 -.00172 -.00060 -.00006
As with the Jacobi method, the results from the Gauss-Seidel method are essentially correct. The Fortran program used to compute the Jacobi iteration method was modified to solve for the Gauss-Seidel iterative method. The program is shown as follows: c c c c
Program Gaus_sdl.for Program to solve a linear equation using the Gauss-Seidel Iteration method IMPLICIT none REAL*8 coef(3,4), d, dx(3), x(3,4), xn(3), xnp(3)
Iterative Methods for Solving Linear Equations
Page 5
INTEGER i, iter, iterate,j, no, nv DATA iterate /0/ c c c c c c c c c
The data are entered into the program using a data file called gauss.dat. It has the following row values number of equations number of variables x(1) x(2) x(3) solution for the first equation x(1) x(2) x(3) solution for the second equation x(1) x(2) x(3) solution for the third equation OPEN (4, file = 'gauss.dat') OPEN (6, file = 'results')
c c c
5
c c c c
no is the number of equations and nv is the number of variables read(4,*) no do 5 i=1,no xn(i) = 0.d0 xnp(i) = 0.d0 continue read(4,*) nv write(6,901) The coefficients for the variables are read in the matrix x with the solution to the equations being the last column do 10 i=1,no read(4,*)(x(i,j),j=1,no+1)
c c c c c
7 c c c c c c
10
13 15 c c c
d is the coefficient for the variable that is being solved for it forms the denominator to compute the real number for the remaining coefficients d = x(i,i) do 7 j=1,no+1 coef(i,j) = x(i,j)/d end do Because the Jacobi method solves for the respect to the current estimates of the coefficient for the variable is made to use in the loop to compute the adjusted
unknown variable with other variables, the be zero for subsequent estimates
coef(i,i) = 0.d0 write(6,900)(x(i,j),j=1,nv+1) end do write(6,902) do 13 i=1,no write(6,900)(coef(i,j),j=1,nv+1) end do write(6,903) iter = 0 iterate is just a counter to keep track of the number of iterations iterate = iterate+1
c c c
Solve for the estimate of the unknown parameters do 20 i=1,no xn(i) = coef(i,nv+1) do 18 j=1,nv
Iterative Methods for Solving Linear Equations
Page 6
18 20 c c c c
xn(i) = xn(i) - coef(i,j)*xn(j) end do end do dx is a vector showing the change in the estimate of the variable with respect to the estimated value used in the previous iteration do 50 i=1,no dx(i) = xnp(i) - xn(i) xnp(i) = xn(i)
c c c c c c
Test to see if the change is greater than the threshold If it is, then the variable iter is made equal to 1 At the beginning of each loop, this value is made equal to 0 If iter is 1 then this means to iterate again if (dabs(dx(i)).gt.0.01d0) iter = 1
c c c 50
900 901
902 903
904
Update the estimated parameter value xn(i) = xnp(i) continue write(6,904)iterate,(xn(i),i=1,nv),(dx(i),i=1,nv) if (iter.gt.0) go to 15 FORMAT(5x,4(f10.4,5x)) FORMAT(15x,'SOLVING LINEAR EQUATION USING THE GAUSS-SEIDEL ITERATI 1ON METHOD',//,'The coefficients to the equations with the solution 1 at the end are:',/) FORMAT(//,'Rearranging the equations to solve for the unknown vari 1ables yields',/,5x,'the following coefficients: ',/) FORMAT(//,'The estimated results after each iteration are shown as 1:'//,'Iteration',2x,'x(1)',9x,'x(2)',9x,'x(3)',7x,'dx(1)',7x,'dx(2 2)',7x,'dx(3)') FORMAT(i3,4x,6(f10.5,2x)) stop end
A more sophisticated subroutine for solving the Gauss-Seidel is shown as [source unknown]: subroutine gsitrn(a,b,x,n,ndim,niter,tol,ierr) c----------------------------------------------------------------------c c Gauss-Seidel Iterative Method c ***************************** c c This subroutine obtaines the solution to n linear equations by Gaussc Seidel iteration. An initial approximation is sent to the subroutine c in the vector x. The solution, as approximated by the subroutine is c returned in x. The iterations are continued until the maximum change c in any x component is less than tol. If this cannot be accomplished c in niter iterations, a non-zero error flag is returned. The matrix c is to be arranged so as to have the largest values on the diagonal. c (from "Applied Numerical Analysis," C.F. Gerald, p 138) c c c c INPUT/OUTPUT VARIABLES: c Iterative Methods for Solving Linear Equations
Page 7
c a(n,n) coefficient matrix with largest values on diagonal c b(n) right hand side vector c x(n) solution vector (initial guess) c n Dimension of the system you're solving c ndim Dimension of matrix a (Note: In the main program, c matrix a may have been dimensioned larger than c necessary, i.e. n, the size of the system you're c decomposing, may be smaller than ndim.) c niter Number of iterations c tol tolerance for solution c ierr Error code: ierr=0 - no errors; ierr=1 - the c solution was not obtained in maximum iterations c c----------------------------------------------------------------------dimension a(ndim,ndim),b(ndim),x(ndim) ierr = 0 c c We can save some divisions by making all the diagonal c elements equal to unity: c do 1 i=1,n temp = 1.0 / a(i,i) b(i) = b(i) * temp do 2 j=1,n a(i,j) = a(i,j) * temp 2 continue 1 continue c c Now we perform the iterations. Store max change in x values for c testing against tol. Outer loop limits iterations to niter: c do 3 iter=1,niter xmax = 0.0 do 4 i=1,n temp = x(i) x(i) = b(i) do 5 j=1,n if(i.ne.j) then x(i) = x(i) - a(i,j)*x(j) endif 5 continue if(abs(x(i)-temp).gt.xmax) xmax = abs(x(i)-temp) 4 continue if(xmax.le.tol) return 3 continue c c Normal exit from the loop means non-convergent solution. Flag the c error and control back to the calling routine: c ierr = 1 return end
Iterative Methods for Solving Linear Equations
Page 8
These iterative methods can also be very effectively programmed in a spreadsheet like Excel. For example, the Jacobi method of solving for linear equations can be shown as:
The formulas for computing x1 and x2 within the spreadsheet are shown as:
It is a simple process of copying and pasting to add more lines to solve the equations. In a similar fashion, the Gauss-Seidel method can also be programmed within Excel to arrive at the same results, as shown in the following figure.
Iterative Methods for Solving Linear Equations
Page 9
While the Gauss-Seidel method appears to be the best in this example, this is not always the case. In fact, it is very possible that the solution from either of these methods could be in error. For example, lets look at the following equations: 12x1 + 3x 2 + 4x 3 = 48 6 x1 + 15x 2 − 4 x 3 = 54.6 9x 1 − 4 x 2 + 6x 3 = 22.3 Then, using a criteria of 0.1 (this is the significant figures for the input values) to stop the iterations, the solution using the Jacobi method can be shown to be: SOLVING LINEAR EQUATION USING THE JACOBI ITERATION METHOD The estimated results after each iteration are shown as: Iteration 1 2 ! 13 14
x(1) 4.00000 1.85111
x(2) 3.64000 3.03111
x(3) 3.71667 .14333
dx(1) 4.00000 -2.14889
dx(2) 3.64000 -.60889
dx(3) 3.71667 -3.57333
2.65235 2.61850
3.05960 3.07234
1.84979 1.77788
.04225 -.03384
-.01281 .01274
.11114 -.07191
The same results using the Gauss-Seidel method, same criteria, are: SOLVING LINEAR EQUATION USING THE GAUSS-SEIDEL ITERATION METHOD The estimated results after each iteration are shown as: Iteration x(1) 1 4.00000 2 3.79778 3 3.77474
x(2) 2.04000 1.87467 1.93538
x(3) -.92333 -.73022 -.65519
dx(1) -4.00000 .20222 .02304
dx(2) -2.04000 .16533 -.06071
dx(3) .92333 -.19311 -.07503
These are markedly different results.
Iterative Methods for Solving Linear Equations
Page 10
Jacobi Iteration Method Given a set of linear equations, a 1 x1 + a 2 x 2 + L + a n x n = s1 b1 x1 + b 2 x 2 +L + bn x n = s 2 M d1 x1 + d 2 x 2 +L + d n x n = s n the problem is one of solving for x1 , x2 , …, xn . The right hand side of these equations, si, represents the solution. We begin by rearranging these equations in the form of solving for the unknown parameters one equation at a time. Thus, s 1 a 2 ( 0 ) a 3 (0 ) a − x2 − x 3 − L − n x (n0 ) a 1 a1 a1 a1 s b b b x 2 = 2 − 1 x 1(0 ) − 3 x 3(0 ) − L − n x (n0 ) b 12 b 2 b2 b2 M s d d d x n = n − 1 x 1(0 ) − 2 x (20 ) − L − n −1 x (n0−)1 an dn dn dn x1 =
The superscript (0) indicates the initial estimate of the parameters. For the first , these parameters are given the value zero. The equations are then solved which results in an updated value of the parameters. These current estimates are then inserted back into the equations and a newer set of parameters is arrived at by solving these equations. The process continues until the solution converges. As an example, take the following linear equations: 7x 1 + 3x 2 + x 3 = 18 2 x 1 − 9x 2 + 4 x 3 = 12 x 1 − 4 x 2 + 12 x 3 = 6
Rearrange these equations x1 =
18 3 1 − x2 − x3 7 7 7
⇒
x 1 = 2.571 − 0.429x 2 − 0.143x 3
12 2 4 + x1 + x3 9 9 9
⇒
x 2 = −1.333 + 0.222x 1 + 0.444 x 3
6 1 4 − x1 + x 2 12 12 12
⇒
x 3 = 0.500 − 0.083x 1 + 0.333x 2
x2 = −
x3 =
Use as the initial estimates: x1 (0) = x2 (0) = x3 (0) = 0. Insert these estimates into these equations yielding new estimates of the parameters. x 1(1) = 2.571 − 0.429 (0) − 0.143 (0) = 2.571
x (21) = −1.333 + 0.222 (0 ) + 0.444 (0 ) = −1.333 x (31 ) = 0.500 − 0.083 (0 ) + 0.333 (0) = 0.500
Insert these updated estimates back into original equation again, yielding x 1( 2 ) = 2.571 − 0.429 (− 1.333) − 0.143 (0.500 ) = 3.071
x (22 ) = −1.333 + 0.222 (2.571) + 0.444 (0.500) = −0.540 x (32 ) = 0.500 − 0.083 (2.571) + 0.333 (−1.333) = −0.159
This process is continued until the desired results are obtained. The following table shows the solutions arrived at after each iteration. These results are from the attached Fortran program. The output shows x(i) that are the parameters xi. Also output is the change in the parameters, dx(i), between each iteration. SOLVING LINEAR EQUATION USING THE JACOBI ITERATION METHOD The estimated results after each iteration are shown as: Iteration 1 2 3 4 5 6
x(1) 2.57143 3.07143 2.82540 2.87141 2.85811 2.85979
x(2) -1.33333 -.53968 -.72134 -.67695 -.68453 -.68261
x(3) .50000 -.15873 .06415 .02410 .03506 .03365
dx(1) 2.57143 .50000 -.24603 .04601 -.01330 .00168
dx(2) -1.33333 .79365 -.18166 .04439 -.00757 .00192
dx(3) .50000 -.65873 .22288 -.04005 .01096 -.00142
How good are these results? Lets take our equations and put them into an augmented matrix and solve using Gauss-Jordan elimination. Iterative Methods for Solving Linear Equations
Page 2
6 R + 31R 7 3 1 18 R1 ↔ R 3 1 − 4 12 6 R 2 − 2 R 1 1 − 4 12 2 2 − 9 4 12 2 − 9 4 12 0 −1 − 20 3 0 → R→ → −7 R −R2 1 − 4 12 6 7 3 1 18 3 1 0 31 − 83 − 24 12 6 R +4 R 1 0 92 6 R −92 R 1 0 0 2.859 1 − 4 1 2 0 1 1 3 0 1 0 − 0.683 20 0 0 1 20 0 R→ r→ − 20 R 3 0 0 − 703 − 24 3 703 0 0 1 0.034 2 0 0 1 0.034 As one can see, the values using the Jacobi iterative method are very close. Following is a Fortran program that can be used to use the Jacobi iteration to solve a set of equations. The limitation now is that it is restricted to only a 3 x 3 matrix, due to formatting procedures currently used in the program. c c c c
Program Jacobi Program to solve a linear equation using the Jacobi Iteration method IMPLICIT none REAL*8 coef(3,4), d, dx(3), x(3,4), xn(3), xnp(3) INTEGER i, iter, iterate,j, no, nv DATA iterate /0/
c c c c c c c c c
The data are entered into the program using a data file called jacobi.dat. It has the following row values number of equations number of variables x(1) x(2) x(3) solution for the first equation x(1) x(2) x(3) solution for the second equation x(1) x(2) x(3) solution for the third equation OPEN (4, file = 'jacobi.dat') OPEN (6, file = 'results')
c c c
5
c c c c
no is the number of equations and nv is the number of variables read(4,*) no do 5 i=1,no xn(i) = 0.d0 continue read(4,*) nv write(6,901) The coefficients for the variables are read in the matrix x with the solution to the equations being the last column do 10 i=1,no read(4,*)(x(i,j),j=1,no+1)
c c c c c
d is the coefficient for the variable that is being solved for it forms the denominator to compute the real number for the remaining coefficients d = x(i,i) do 7 j=1,no+1
Iterative Methods for Solving Linear Equations
Page 3
7 c c c c c c
10
13 15 c c c
coef(i,j) = x(i,j)/d end do Because the Jacobi method solves for the respect to the current estimates of the coefficient for the variable is made to use in the loop to compute the adjusted
unknown variable with other variables, the be zero for subsequent estimates
coef(i,i) = 0.d0 write(6,900)(x(i,j),j=1,nv+1) end do write(6,902) do 13 i=1,no write(6,900)(coef(i,j),j=1,nv+1) end do write(6,903) iter = 0 iterate is just a counter to keep track of the number of iterations iterate = iterate+1
c c c
18 20 c c c c
Solve for the estimate of the unknown parameters do 20 i=1,no xnp(i) = coef(i,nv+1) do 18 j=1,nv xnp(i) = xnp(i) - coef(i,j)*xn(j) end do end do dx is a vector showing the change in the estimate of the variable with respect to the estimated value used in the previous iteration do 50 i=1,no dx(i) = xnp(i) - xn(i)
c c c c c c
Test to see if the change is greater than the threshold If it is, then the variable iter is made equal to 1 At the beginning of each loop, this value is made equal to 0 If iter is 1 then this means to iterate again if (dabs(dx(i)).gt.0.01d0) iter = 1
c c c 50
900 901
902 903
904
Update the estimated parameter value xn(i) = xnp(i) continue write(6,904)iterate,(xn(i),i=1,nv),(dx(i),i=1,nv) if (iter.gt.0) go to 15 FORMAT(5x,4(f10.4,5x)) FORMAT(15x,'SOLVING LINEAR EQUATION USING THE JACOBI ITERATION MET 1HOD',//,'The coefficients to the equations with the solution at th 1e end are:',/) FORMAT(//,'Rearranging the equations to solve for the unknown vari 1ables yields',/,5x,'the following coefficients: ',/) FORMAT(//,'The estimated results after each iteration are shown as 1:'//,'Iteration',2x,'x(1)',9x,'x(2)',9x,'x(3)',7x,'dx(1)',7x,'dx(2 2)',7x,'dx(3)') FORMAT(i3,4x,6(f10.5,2x)) stop end
Iterative Methods for Solving Linear Equations
Page 4
Gauss-Seidel Iterative Method The Gauss-Seidel iterative method of solving for a set of linear equations can be thought of as just an extension of the Jacobi method. Start out using an initial value of zero for each of the parameters. Then, solve for x1 as in the Jacobi method. When solving for x2 , insert the just computed value for x1 . In other words, for each calculation, the most current estimate of the parameter value is used. To see how the Gauss-Seidel method works, lets use the values in the last example. The initial estimates are set to zero. Then, the results from the first iteration are shown as: x1(1) = 2.571 − 0.429 x (20 ) − 0.143 x (30 ) = 2.571 − 0.429 (0 ) − 0.143 (0 ) = 2.571
x (21) = −1.333 + 0.222 x (11 ) + 0.444 x (30 ) = −1.333 + 0.222 (2.571) + 0.444 (0 ) = −0.762
x (31 ) = 0.500 − 0.083 x1(1) + 0.333 x (21 ) = 0.500 − 0.083 (2.571) + 0.333 (− 0.762 ) = 0.033 The next iteration is performed in a similar fashion. It can be shown as: x (12 ) = 2.571 − 0.429x (21 ) − 0.143x (31 ) = 2.571 − 0.429 (− 0.762) − 0.143 (0.033) = 2.893
x (22 ) = −1.333 + 0.222 x (12 ) + 0.444 x (31) = −1.333 + 0.222 (2.893) + 0.444 (0.033) = −0.676 x (32 ) = 0.500 − 0.083x1( 2 ) + 0.333x (22 ) = 0.500 − 0.083 (2.893) + 0.333 (− 0.676 ) = 0.034
A Fortran program was written to solve this problem. The results are shown in the next table. SOLVING LINEAR EQUATION USING THE GAUSS-SEIDEL ITERATION METHOD The estimated results after each iteration are shown as: Iteration 1 2 3 4
x(1) 2.57143 2.89342 2.85646 2.85957
x(2) -.76190 -.67624 -.68369 -.68273
x(3) .03175 .03347 .03406 .03412
dx(1) -2.57143 -.32200 .03696 -.00311
dx(2) .76190 -.08566 .00745 -.00096
dx(3) -.03175 -.00172 -.00060 -.00006
As with the Jacobi method, the results from the Gauss-Seidel method are essentially correct. The Fortran program used to compute the Jacobi iteration method was modified to solve for the Gauss-Seidel iterative method. The program is shown as follows: c c c c
Program Gaus_sdl.for Program to solve a linear equation using the Gauss-Seidel Iteration method IMPLICIT none REAL*8 coef(3,4), d, dx(3), x(3,4), xn(3), xnp(3)
Iterative Methods for Solving Linear Equations
Page 5
INTEGER i, iter, iterate,j, no, nv DATA iterate /0/ c c c c c c c c c
The data are entered into the program using a data file called gauss.dat. It has the following row values number of equations number of variables x(1) x(2) x(3) solution for the first equation x(1) x(2) x(3) solution for the second equation x(1) x(2) x(3) solution for the third equation OPEN (4, file = 'gauss.dat') OPEN (6, file = 'results')
c c c
5
c c c c
no is the number of equations and nv is the number of variables read(4,*) no do 5 i=1,no xn(i) = 0.d0 xnp(i) = 0.d0 continue read(4,*) nv write(6,901) The coefficients for the variables are read in the matrix x with the solution to the equations being the last column do 10 i=1,no read(4,*)(x(i,j),j=1,no+1)
c c c c c
7 c c c c c c
10
13 15 c c c
d is the coefficient for the variable that is being solved for it forms the denominator to compute the real number for the remaining coefficients d = x(i,i) do 7 j=1,no+1 coef(i,j) = x(i,j)/d end do Because the Jacobi method solves for the respect to the current estimates of the coefficient for the variable is made to use in the loop to compute the adjusted
unknown variable with other variables, the be zero for subsequent estimates
coef(i,i) = 0.d0 write(6,900)(x(i,j),j=1,nv+1) end do write(6,902) do 13 i=1,no write(6,900)(coef(i,j),j=1,nv+1) end do write(6,903) iter = 0 iterate is just a counter to keep track of the number of iterations iterate = iterate+1
c c c
Solve for the estimate of the unknown parameters do 20 i=1,no xn(i) = coef(i,nv+1) do 18 j=1,nv
Iterative Methods for Solving Linear Equations
Page 6
18 20 c c c c
xn(i) = xn(i) - coef(i,j)*xn(j) end do end do dx is a vector showing the change in the estimate of the variable with respect to the estimated value used in the previous iteration do 50 i=1,no dx(i) = xnp(i) - xn(i) xnp(i) = xn(i)
c c c c c c
Test to see if the change is greater than the threshold If it is, then the variable iter is made equal to 1 At the beginning of each loop, this value is made equal to 0 If iter is 1 then this means to iterate again if (dabs(dx(i)).gt.0.01d0) iter = 1
c c c 50
900 901
902 903
904
Update the estimated parameter value xn(i) = xnp(i) continue write(6,904)iterate,(xn(i),i=1,nv),(dx(i),i=1,nv) if (iter.gt.0) go to 15 FORMAT(5x,4(f10.4,5x)) FORMAT(15x,'SOLVING LINEAR EQUATION USING THE GAUSS-SEIDEL ITERATI 1ON METHOD',//,'The coefficients to the equations with the solution 1 at the end are:',/) FORMAT(//,'Rearranging the equations to solve for the unknown vari 1ables yields',/,5x,'the following coefficients: ',/) FORMAT(//,'The estimated results after each iteration are shown as 1:'//,'Iteration',2x,'x(1)',9x,'x(2)',9x,'x(3)',7x,'dx(1)',7x,'dx(2 2)',7x,'dx(3)') FORMAT(i3,4x,6(f10.5,2x)) stop end
A more sophisticated subroutine for solving the Gauss-Seidel is shown as [source unknown]: subroutine gsitrn(a,b,x,n,ndim,niter,tol,ierr) c----------------------------------------------------------------------c c Gauss-Seidel Iterative Method c ***************************** c c This subroutine obtaines the solution to n linear equations by Gaussc Seidel iteration. An initial approximation is sent to the subroutine c in the vector x. The solution, as approximated by the subroutine is c returned in x. The iterations are continued until the maximum change c in any x component is less than tol. If this cannot be accomplished c in niter iterations, a non-zero error flag is returned. The matrix c is to be arranged so as to have the largest values on the diagonal. c (from "Applied Numerical Analysis," C.F. Gerald, p 138) c c c c INPUT/OUTPUT VARIABLES: c Iterative Methods for Solving Linear Equations
Page 7
c a(n,n) coefficient matrix with largest values on diagonal c b(n) right hand side vector c x(n) solution vector (initial guess) c n Dimension of the system you're solving c ndim Dimension of matrix a (Note: In the main program, c matrix a may have been dimensioned larger than c necessary, i.e. n, the size of the system you're c decomposing, may be smaller than ndim.) c niter Number of iterations c tol tolerance for solution c ierr Error code: ierr=0 - no errors; ierr=1 - the c solution was not obtained in maximum iterations c c----------------------------------------------------------------------dimension a(ndim,ndim),b(ndim),x(ndim) ierr = 0 c c We can save some divisions by making all the diagonal c elements equal to unity: c do 1 i=1,n temp = 1.0 / a(i,i) b(i) = b(i) * temp do 2 j=1,n a(i,j) = a(i,j) * temp 2 continue 1 continue c c Now we perform the iterations. Store max change in x values for c testing against tol. Outer loop limits iterations to niter: c do 3 iter=1,niter xmax = 0.0 do 4 i=1,n temp = x(i) x(i) = b(i) do 5 j=1,n if(i.ne.j) then x(i) = x(i) - a(i,j)*x(j) endif 5 continue if(abs(x(i)-temp).gt.xmax) xmax = abs(x(i)-temp) 4 continue if(xmax.le.tol) return 3 continue c c Normal exit from the loop means non-convergent solution. Flag the c error and control back to the calling routine: c ierr = 1 return end
Iterative Methods for Solving Linear Equations
Page 8
These iterative methods can also be very effectively programmed in a spreadsheet like Excel. For example, the Jacobi method of solving for linear equations can be shown as:
The formulas for computing x1 and x2 within the spreadsheet are shown as:
It is a simple process of copying and pasting to add more lines to solve the equations. In a similar fashion, the Gauss-Seidel method can also be programmed within Excel to arrive at the same results, as shown in the following figure.
Iterative Methods for Solving Linear Equations
Page 9
While the Gauss-Seidel method appears to be the best in this example, this is not always the case. In fact, it is very possible that the solution from either of these methods could be in error. For example, lets look at the following equations: 12x1 + 3x 2 + 4x 3 = 48 6 x1 + 15x 2 − 4 x 3 = 54.6 9x 1 − 4 x 2 + 6x 3 = 22.3 Then, using a criteria of 0.1 (this is the significant figures for the input values) to stop the iterations, the solution using the Jacobi method can be shown to be: SOLVING LINEAR EQUATION USING THE JACOBI ITERATION METHOD The estimated results after each iteration are shown as: Iteration 1 2 ! 13 14
x(1) 4.00000 1.85111
x(2) 3.64000 3.03111
x(3) 3.71667 .14333
dx(1) 4.00000 -2.14889
dx(2) 3.64000 -.60889
dx(3) 3.71667 -3.57333
2.65235 2.61850
3.05960 3.07234
1.84979 1.77788
.04225 -.03384
-.01281 .01274
.11114 -.07191
The same results using the Gauss-Seidel method, same criteria, are: SOLVING LINEAR EQUATION USING THE GAUSS-SEIDEL ITERATION METHOD The estimated results after each iteration are shown as: Iteration x(1) 1 4.00000 2 3.79778 3 3.77474
x(2) 2.04000 1.87467 1.93538
x(3) -.92333 -.73022 -.65519
dx(1) -4.00000 .20222 .02304
dx(2) -2.04000 .16533 -.06071
dx(3) .92333 -.19311 -.07503
These are markedly different results.
Iterative Methods for Solving Linear Equations
Page 10
Related Documents 3h463d
Jacobi And Gauss-seidel 632w6d
December 2019 72
Jacobi Method.pdf 145i5j
December 2019 55
Jacobi Program In Fortran 6y4u4y
December 2019 41
Gauss And Jacobi C Program 6y3r5t
December 2019 39
Metode Jacobi 1k1x19
December 2019 41