Hi tcfisher,

Here's one less-tedious way, using Python's list comprehension, string formatting, and an alternate format for the `matrix`

command:

```
sage: plist = ['p%s%s'%(i,j) for i in range(4) for j in range(4)]
sage: plist
['p00', 'p01', 'p02', 'p03', 'p10', 'p11', 'p12', 'p13', 'p20', 'p21', 'p22', 'p23', 'p30', 'p31', 'p32', 'p33']
sage: m = matrix(SR,4,4,plist); m
[p00 p01 p02 p03]
[p10 p11 p12 p13]
[p20 p21 p22 p23]
[p30 p31 p32 p33]
```

Note that this does not inject the variables '`pij`

' into the global namespace, so trying to write them directly will result in an error:

```
sage: p00
Traceback (most recent call last):
...
NameError: name 'p00' is not defined
```

But if all you want to do is manipulate these matrices, that will work fine:

```
sage: m^2
[ p00^2 + p01*p10 + p02*p20 + p03*p30 p00*p01 + p01*p11 + p02*p21 + p03*p31 p00*p02 + p01*p12 + p02*p22 + p03*p32 p00*p03 + p01*p13 + p02*p23 + p03*p33]
[p00*p10 + p10*p11 + p12*p20 + p13*p30 p01*p10 + p11^2 + p12*p21 + p13*p31 p02*p10 + p11*p12 + p12*p22 + p13*p32 p03*p10 + p11*p13 + p12*p23 + p13*p33]
[p00*p20 + p10*p21 + p20*p22 + p23*p30 p01*p20 + p11*p21 + p21*p22 + p23*p31 p02*p20 + p12*p21 + p22^2 + p23*p32 p03*p20 + p13*p21 + p22*p23 + p23*p33]
[p00*p30 + p10*p31 + p20*p32 + p30*p33 p01*p30 + p11*p31 + p21*p32 + p31*p33 p02*p30 + p12*p31 + p22*p32 + p32*p33 p03*p30 + p13*p31 + p23*p32 + p33^2]
```

And of course if you need the '`pij`

' in the global namespace, you can always declare the variables by looping over `plist`

(or maybe something more clever, if you give it some thought, but I haven't really :)

```
sage: for p in plist:
....: var(p)
....:
p00
p01
p02
p03
p10
p11
p12
p13
p20
p21
p22
p23
p30
p31
p32
p33
sage: p00
p00
sage: p00 + 2*p13
p00 + 2*p13
```