# An Interesting, and Perhaps New, Matrix

I do not recall having seen this matrix before, but I will not be surprised to learn that somebody else already knows all about it, especially if that person's name is Nick.

### Contents

`Q`

I've been investigating the matrices generated by this elegant one-liner.

Q = @(n) (-n:n).^2 + (-n:n)'.^2;

The `Q` is for "quadratic".

The middle column contains the squares of the integers from `-n` to `n`. So does the middle row. The apostrophe summons singleton expansion. The resulting matrix has order `2*n+1`. Here is `Q(5)`.

Q5 = Q(5)

Q5 = 50 41 34 29 26 25 26 29 34 41 50 41 32 25 20 17 16 17 20 25 32 41 34 25 18 13 10 9 10 13 18 25 34 29 20 13 8 5 4 5 8 13 20 29 26 17 10 5 2 1 2 5 10 17 26 25 16 9 4 1 0 1 4 9 16 25 26 17 10 5 2 1 2 5 10 17 26 29 20 13 8 5 4 5 8 13 20 29 34 25 18 13 10 9 10 13 18 25 34 41 32 25 20 17 16 17 20 25 32 41 50 41 34 29 26 25 26 29 34 41 50

I like the contour plot.

contourf(Q(100)) axis square colorbar title('Q(100)')

`D`

For another blog post under development, I need a logical mask that carves a circular region out of graphic. This disc does the job.

D = @(n) Q(n) <= n^2;

Here is my carver.

```
spy(D(100))
title('D(100)')
```

Did you notice the digits in the count of nonzeros in `D(100)`? It happens whenever `n` is a power of 10.

fprintf('%15s %12s\n','n','nnz(D(n))') for n = 10.^(0:4) fprintf('%15d %12d\n',n, nnz(D(n))) end

n nnz(D(n)) 1 5 10 317 100 31417 1000 3141549 10000 314159053

#### O.E.I.S.

A classic problem, described in the Online Encyclopedia of Integer Sequences, asks how many points with integer coordinates lie within a disc of increasing radius. Our nonzero count provides the answer.

fprintf('%15s %8s\n','n','a(n)') for n = [1:15 99:101 499:501 999:1001] if mod(n,100) == 99 fprintf('%15s %8s\n','-','-') end a(n) = nnz(D(n)); fprintf('%15d %8d\n',n,a(n)) end

n a(n) 1 5 2 13 3 29 4 49 5 81 6 113 7 149 8 197 9 253 10 317 11 377 12 441 13 529 14 613 15 709 - - 99 30757 100 31417 101 32017 - - 499 782197 500 785349 501 788509 - - 999 3135157 1000 3141549 1001 3147833

`R`

Taking the reciprocals of the matrix entries, and reducing the range of the anonymous index, produces a matrix that behaves a bit like the Hilbert matrix, `hilb(n)`.

R = @(n) 1./((1:n).^2 + (1:n)'.^2);

Here are the 5-by-5's.

```
format rat
R5 = R(5)
H5 = hilb(5)
```

R5 = 1/2 1/5 1/10 1/17 1/26 1/5 1/8 1/13 1/20 1/29 1/10 1/13 1/18 1/25 1/34 1/17 1/20 1/25 1/32 1/41 1/26 1/29 1/34 1/41 1/50 H5 = 1 1/2 1/3 1/4 1/5 1/2 1/3 1/4 1/5 1/6 1/3 1/4 1/5 1/6 1/7 1/4 1/5 1/6 1/7 1/8 1/5 1/6 1/7 1/8 1/9

#### Condition

Going away from the diagonal, the elements of `R(n)` decay more rapidly than those of `hilb(n)`, so `R(n)` is better conditioned than `hilb(n)`.

format short e fprintf('%15s %12s %12s\n','n','cond R','cond H') for n = 1:12 fprintf('%15d %12.2e %12.2e\n',n,cond(R(n)),cond(hilb(n))) end

n cond R cond H 1 1.00e+00 1.00e+00 2 1.53e+01 1.93e+01 3 2.04e+02 5.24e+02 4 2.59e+03 1.55e+04 5 3.21e+04 4.77e+05 6 3.89e+05 1.50e+07 7 4.67e+06 4.75e+08 8 5.54e+07 1.53e+10 9 6.53e+08 4.93e+11 10 7.65e+09 1.60e+13 11 8.92e+10 5.22e+14 12 1.04e+12 1.62e+16

#### Extra Credit

What is the rank of `Q(n)`? Why? See a paper in SIAM Review by Strang and Moler.

Why is the table of values for `nnz(D(10^k))` so short? How might you extend this table?

Investigate `R(n)`. Is it positive definite? What are its eigenvalues? What is its inverse? What is the sign pattern of the elements of its inverse? For what values of `n` can you compute the inverse reliably using floating point arithmetic? How does all this compare with `hilb(n)` and `invhilb(n)` ?

Comments in the Comments, or in email to me, are welcome.

Get
the MATLAB code

Published with MATLAB® R2022a