18 The forces on the bridge truss described in the opening t

18. The forces on the bridge truss described in the opening to this chapter satisty the equations in the following table: Joint Horizontal Component Vertical Component -10,000 = 0 This linear system can be placed in the matrix form 1 00 0 0 0 0- /a 0 10,000 0 0 o o o 0-9-1 0 0000 -1 a. Explain why the system of equations was reordered. b. Approximate the solution of the resulting linear system to within 10in thenorm using as initial approximation the vector all of whose entries are Is with (i) the Jacobi method and (ii) the Gauss-Seidel method.

Solution

function x = jacobi( M, b, N, e ) % Solve Mx = b % The diagonal entries of M and their inverses d = diag( M ); if ~all( d ) error \'at least one diagonal entry is zero\'; end invd = d.^-1; % Matrix of off-diagonal entires of N Moff = M - diag( d ); % Use d.^-1*b as the first approximation to x invdb = invd.*b; x = db; % -1 % Iterate x = D (b - M *x) % off for k = 1:N xprev = x; x = invdb - invd.*(Moff*x); if norm( x - xprev, inf ) < e return; end end error \'the method did not converge\'; end function x = gauss_seidel( M, b, N, e ) % Solve Mx = b % The diagonal entries of M and their inverses n = length( b ); d = diag( M ); if ~all( d ) error \'at least one diagonal entry is zero\'; end invd = d.^-1; % Matrix of off-diagonal entires of N Moff = M - diag( d ); % Use d.^-1*b as the first approximation to x invdb = invd.*b; x = db; % -1 % Iterate x = D (b - M *x) % off for k = 1:N xprev = x; for i = 1:n x(i) = invdb(i) - invd(i).*(Moff(i,:)*x); end if norm( x - xprev, inf ) < e return; end end error \'the method did not converge\'; end