Using Finite Difference preconditioners to solve the variable diffusivity Helmholtz problem
This code is a tutorial code for solving a variable coefficient Helmholtz (Poisson) problem in a doubly periodic domain. It uses GMRES preconditioned with 2nd order finite differences. Sometimes if you have a very dense matrix, it can take an extremely long time to converge using iterative methods. One way to improve this is to "precondition" your matrix by multiplying it by an approximation to the solution so that the condition number is smaller and thus will converge faster. Instead of solving Ax=b, you solve AP^{-1}y=b for y, and then solve Px=y. Please note that this script defines functions at the end, which is only supported by MATLAB 2016b or later.
Contents
- Construct grid
- GMRES parameters
- FFT related material
- Diffusivity
- Beginning the Pre-conditioner material
- Construct the Preconditioner
- Compute the LU factorization
- Solve the Problem using GMRES
- Plot the Results
- Helmholtz Operator
- Compute the derivaitves
- Compute the operator
- Return as a vector
- Precondition using LU
Construct grid
Define the various grid parameters.
Nx = 256; % uonal resolution Ny = 256; % meridional resolution Lx=3e3; Ly=3e3; dx=Lx/Nx; dy=Ly/Ny; x=dx*(1:Nx)'; y=dy*(1:Ny)'; dx2o2 = dx^2/2; dy2o2 = dy^2/2; [xx,yy] = meshgrid(x,y);
GMRES parameters
These parameters are used by the GMRES function to tell it to stop iterating.
TOL=1e-9; % GMRES convergence tolerance MAXIT = min(250,Nx); % Upper-bound for # of iterations GMRES can perform.
FFT related material
For this problem we will be solving the PDE using spectral methods. For the constant diffusivity case see adv_spectralhelmholtz.m.
Build the wavenumbers and create the variable diffusivity.
dk = 2*pi/Lx; % Nyquist wavenumbers dl = 2*pi/Ly; k=[0:Nx/2-1 Nx/2 -Nx/2+1:-1]'*dk; l=[0:Ny/2-1 Ny/2 -Ny/2+1:-1]'*dl; [kk,ll]=meshgrid(k,l); l2 = l.^2; % 1D meridional wavenumber arrays for splitting up the 1D problems il = 1i*l; k2 = k.^2; % 1D uonal wavenumber arrays for differentiating within a sub problem ik = 1i*k; ll2 = ll.^2; % 2D meridional wavenumber arrays for differentiation purposes ill = 1i*ll; kk2 = kk.^2; % 2D uonal wavenumber arrays for general differentiation purposes ikk = 1i*kk;
Diffusivity
The key difference for this case (and why we need to use GMRES) is that we have a variable diffusivity and so it is no longer trivial to solve the matrix problem in Fourier space.
myr=sqrt(4*(xx-0.5*Lx).^2 + (yy-0.5*Ly).^2); kappa = 100*(1-0.95*exp(-(myr/(0.3*Lx)).^4)); kappax= real(ifft(ikk.*fft(kappa,[],2),[],2)); kappay= real(ifft(ill.*fft(kappa,[],1),[],1));
Beginning the Pre-conditioner material
Create the finite difference derivative matrices and make sure the boundary conditions are correct.
disp('Building preconditioner...');
tic;
Building preconditioner...
Build 1st and 2nd derivative matrices from FD (for preconditioner)
disp('Calculating FD matrices...'); FDx=spalloc(Nx,Nx,3*Nx); FDxx=spalloc(Nx,Nx,3*Nx); FDy=spalloc(Ny,Ny,3*Ny); FDyy=spalloc(Ny,Ny,3*Ny); matx=[1 1 1;-dx 0 dx;dx2o2 0 dx2o2]; maty=[1 1 1;-dy 0 dy;dy2o2 0 dy2o2]; coefs1x=matx\[0;1;0]; coefs1y=maty\[0;1;0]; coefs2x=matx\[0;0;1]; coefs2y=maty\[0;0;1]; for ii=2:Nx-1 FDx(ii,ii-1:ii+1)=coefs1x; FDxx(ii,ii-1:ii+1)=coefs2x; end for ii=2:Ny-1 FDy(ii,ii-1:ii+1)=coefs1y; FDyy(ii,ii-1:ii+1)=coefs2y; end
Calculating FD matrices...
%Employ periodicity at the ends
FDx(1, [Nx 1 2]) = coefs1x;
FDx(Nx, [Nx-1 Nx 1]) = coefs1x;
FDxx(1, [Nx 1 2]) = coefs2x;
FDxx(Nx, [Nx-1 Nx 1]) = coefs2x;
FDy(1, [Ny 1 2]) = coefs1y;
FDy(Ny, [Ny-1 Ny 1]) = coefs1y;
FDyy(1, [Ny 1 2]) = coefs2y;
FDyy(Ny, [Ny-1 Ny 1]) = coefs2y;
Now build big matrices
Ix = speye(Nx,Nx); Iy = speye(Ny,Ny);
Construct the Preconditioner
Do this in pieces so we don't run out of memory.
disp('Constructing actual preconditioner...'); kappadiag = spdiags(kappa(:),0,Nx*Ny,Nx*Ny); bigFDxx = kron(FDxx,Iy); precon = kappadiag*bigFDxx; clear bigFDxx
Constructing actual preconditioner...
bigFDyy=kron(Ix,FDyy); precon = precon + kappadiag*bigFDyy; clear bigFDyy kappadiag
kappaxdiag = spdiags(kappax(:),0,Nx*Ny,Nx*Ny); bigFDx=kron(FDx,Iy); precon = precon + kappaxdiag*bigFDx; clear bigFDx kappaxdiag
kappaydiag = spdiags(kappay(:),0,Nx*Ny,Nx*Ny); bigFDy=kron(Ix,FDy); precon = precon + kappaydiag*bigFDy; clear bigFDy kappaydiag clear Ix Iy FDx FDxx FDy FDyy coefs1x coefs2x coefs1y coefs2y disp('done');
done
Compute the LU factorization
tic;
[ll,uu,pp,qq]=lu(precon);
toc;
disp('done');
UMFPACK V5.4.0 (May 20, 2009), Control: Matrix entry defined as: double Int (generic integer) defined as: UF_long 0: print level: 2 1: dense row parameter: 0.2 "dense" rows have > max (16, (0.2)*16*sqrt(n_col) entries) 2: dense column parameter: 0.2 "dense" columns have > max (16, (0.2)*16*sqrt(n_row) entries) 3: pivot tolerance: 0.1 4: block size for dense matrix kernels: 32 5: strategy: 0 (auto) 6: initial allocation ratio: 0.7 7: max iterative refinement steps: 2 12: 2-by-2 pivot tolerance: 0.01 13: Q fixed during numerical factorization: 0 (auto) 14: AMD dense row/col parameter: 10 "dense" rows/columns have > max (16, (10)*sqrt(n)) entries Only used if the AMD ordering is used. 15: diagonal pivot tolerance: 0.001 Only used if diagonal pivoting is attempted. 16: scaling: 0 (no) 17: frontal matrix allocation ratio: 0.5 18: drop tolerance: 0 19: AMD and COLAMD aggressive absorption: 1 (yes) The following options can only be changed at compile-time: 8: BLAS library used: Fortran BLAS. size of BLAS integer: 8 9: compiled for MATLAB 10: CPU timer is POSIX times ( ) routine. 11: compiled for normal operation (debugging disabled) computer/operating system: unknown size of int: 4 UF_long: 8 Int: 8 pointer: 8 double: 8 Entry: 8 (in bytes) UMFPACK V5.4.0 (May 20, 2009), Info: matrix entry defined as: double Int (generic integer) defined as: UF_long BLAS library used: Fortran BLAS. size of BLAS integer: 8 MATLAB: yes. CPU timer: POSIX times ( ) routine. number of rows in matrix A: 65536 number of columns in matrix A: 65536 entries in matrix A: 327680 memory usage reported in: 16-byte Units size of int: 4 bytes size of UF_long: 8 bytes size of pointer: 8 bytes size of numerical entry: 8 bytes strategy used: symmetric ordering used: amd on A+A' modify Q during factorization: no prefer diagonal pivoting: yes pivots with zero Markowitz cost: 0 submatrix S after removing zero-cost pivots: number of "dense" rows: 0 number of "dense" columns: 0 number of empty rows: 0 number of empty columns 0 submatrix S square and diagonal preserved pattern of square submatrix S: number rows and columns 65536 symmetry of nonzero pattern: 1.000000 nz in S+S' (excl. diagonal): 262144 nz on diagonal of matrix S: 65536 fraction of nz on diagonal: 1.000000 AMD statistics, for strict diagonal pivoting: est. flops for LU factorization: 1.38146e+09 est. nz in L+U (incl. diagonal): 5253060 est. largest front (# entries): 580644 est. max nz in any column of L: 762 number of "dense" rows/columns in S+S': 0 symbolic factorization defragmentations: 0 symbolic memory usage (Units): 1605819 symbolic memory usage (MBytes): 24.5 Symbolic size (Units): 229554 Symbolic size (MBytes): 4 symbolic factorization CPU time (sec): 0.12 symbolic factorization wallclock time(sec): 0.12 matrix scaled: no symbolic/numeric factorization: upper bound actual % variable-sized part of Numeric object: initial size (Units) 1572871 1507332 96% peak size (Units) 117714908 3567256 3% final size (Units) 114847127 2853864 2% Numeric final size (Units) 115207611 3148812 3% Numeric final size (MBytes) 1757.9 48.0 3% peak memory usage (Units) 118655709 4475289 4% peak memory usage (MBytes) 1810.5 68.3 4% numeric factorization flops 1.85846e+11 1.38146e+09 1% nz in L (incl diagonal) 60451329 2659298 4% nz in U (incl diagonal) 92422791 2659298 3% nz in L+U (incl diagonal) 152808584 5253060 3% largest front (# entries) 2608154 580644 22% largest # rows in front 1357 762 56% largest # columns in front 1922 762 40% initial allocation ratio used: 0.0438 # of forced updates due to frontal growth: 0 number of off-diagonal pivots: 0 nz in L (incl diagonal), if none dropped 2659298 nz in U (incl diagonal), if none dropped 2659298 number of small entries dropped 0 nonzeros on diagonal of U: 65536 min abs. value on diagonal of U: 2.62e-12 max abs. value on diagonal of U: 2.91e+00 estimate of reciprocal of condition number: 9.00e-13 indices in compressed pattern: 487433 numerical values stored in Numeric object: 5253060 numeric factorization defragmentations: 1 numeric factorization reallocations: 1 costly numeric factorization reallocations: 1 numeric factorization CPU time (sec): 0.55 numeric factorization wallclock time (sec): 0.33 numeric factorization mflops (CPU time): 2511.75 numeric factorization mflops (wallclock): 4186.25 symbolic + numeric CPU time (sec): 0.67 symbolic + numeric mflops (CPU time): 2061.88 symbolic + numeric wall clock time (sec): 0.45 symbolic + numeric mflops (wall clock): 3069.91 Elapsed time is 0.629444 seconds. done
Solve the Problem using GMRES
Define the RHS and solve the Helmholtz problem using GMRES. GMRES takes the inputs: the matrix A, the vector b, restart number [], tolerance of the method TOL, maximum number of iterations MAXIT, and the preconditioner matrix M. A does not have to be a matrix explicitly, here we have used an anonymous function that calls the karhelm function so that its output is the matrix. This allows more parameters to be passed to the karnhelm function than is traditionally allowed. A similar procedure is used for the preconditioner matrix.
RHShelm = sin(6*pi*xx/Lx).*sin(4*pi*yy/Ly); disp('Helmholtz inversion without preconditioner'); tic; [u00]=gmres(@(foo)func_helmholtz(foo,kappax,kappay,kappa,ill,ikk,Nx,Ny),RHShelm(:),[],TOL,MAXIT); toc; disp('Helmholtz inversion with incomplete LU preconditioner...'); tic; [u0]=gmres(@(foo)func_helmholtz(foo,kappax,kappay,kappa,ill,ikk,Nx,Ny),RHShelm(:),[],TOL,MAXIT,@(bar)func_preconLU(bar,ll,uu,pp,qq)); toc; disp('Zero diffusivity helmholtz inversion without preconditioner'); tic; [u0const]=gmres(@(foo)func_helmholtz(foo,zeros(size(kappa)),zeros(size(kappa)),mean(kappa(:))*ones(size(kappa)),ill,ikk,Nx,Ny),RHShelm(:),[],TOL,MAXIT); toc;
Helmholtz inversion without preconditioner
Plot the Results
u0=reshape(u0,Ny,Nx); % the preconditioned solution u00=reshape(u00,Ny,Nx); % the not preconditioned solution u0const=reshape(u0const,Ny,Nx); % the constant diffusivity solution ud=u0-u0const; % the difference between constant and non-constant diffusivity udo = u0-u00; % the difference between preconditioned and not maxdiff=max(abs(ud(:)))/max(abs(u0(:))) maxdiff0 = max(abs(udo(:)))/max(abs(u0(:))) figure(1) clf subplot(2,2,1) pcolor(xx,yy,u0),shading flat; colorbar; title('preconditioned var diff') subplot(2,2,2) pcolor(xx,yy,u0const),shading flat; colorbar; title('const diff') subplot(2,2,3) pcolor(xx,yy,u0-u0const),shading flat; colorbar; title('var - const') subplot(2,2,4) pcolor(xx,yy,u0-u00),shading flat; colorbar; title('precon - no precon') figure(2) clf pcolor(xx,yy,kappa),shading flat; colorbar; title('diffusivity')
maxdiff = 0.0219 maxdiff0 = 0.0219
Helmholtz Operator
This function returns the results of the Helmholtz operator taking in the variable u, the diffusivity kappa, its derivatives kappax and kappay, i times the wavenumber vector squared, ikk and ill, and the grid sizes Nx and Ny.
function result = func_helmholtz(u,kappax,kappay,kappa,ill,ikk,Nx,Ny)
u = reshape(u,Ny,Nx);
Compute the derivaitves
ux = real(ifft(ikk.*fft(u,[],2),[],2)); uxx = real(ifft(ikk.*fft(ux,[],2),[],2)); uy = real(ifft(ill.*fft(u,[],1),[],1)); uyy = real(ifft(ill.*fft(uy,[],1),[],1));
Compute the operator
result = kappax.*ux + kappa.*uxx + kappay.*uy + kappa.*uyy - u;
Return as a vector
result = result(:);
end
gmres converged at iteration 11 to a solution with relative residual 8.6e-10. Elapsed time is 0.599190 seconds. Helmholtz inversion with incomplete LU preconditioner...
gmres converged at iteration 1 to a solution with relative residual 3.2e-15. Elapsed time is 0.151477 seconds.
Precondition using LU
Pre-condition using stored sparse LU factors. (with permutation matrices P and Q)
function approxans = func_preconLU(rhs,L,U,P,Q) approxans = Q*(U\(L\(P*rhs))); end
sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. Warning: Input tol may not be achievable by GMRES. Try to use a bigger tolerance. sp\: bandwidth = 65534+1+0. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. sp\: bandwidth = 0+1+65534. sp\: is A diagonal? no. sp\: is A triangular? yes. sp\: do a triangular solve. gmres stopped at iteration 151 without converging to the desired tolerance 1e-09 because the method stagnated. The iterate returned (number 146) has relative residual 0.0061. Elapsed time is 12.284280 seconds. Zero diffusivity helmholtz inversion without preconditioner