# 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.

## 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...');
bigFDxx = kron(FDxx,Iy);
clear bigFDxx
```
```Constructing actual preconditioner...
```
```bigFDyy=kron(Ix,FDyy);
```
```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)
title('preconditioned var diff')
subplot(2,2,2)
title('const diff')
subplot(2,2,3)
title('var - const')
subplot(2,2,4)
title('precon - no precon')
figure(2)
clf
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
```