Matrix Op Storage and Speed :)
Warning do not make the matrices too big!

Matrices and Sparse Matrices and Speed (Oh My!)
Differencing and Sparse Matrices

I am going to implement a few difference operators and give them semi-meaningful names.  First I am going to assume that I have a table of values that are naturally arranged in a grid of x,y values.  I am going to call the grid of values U[[i,j]] and we are going to assume that U[[i,j]]≈u[(Subscript[x, i]),(Subscript[y, j])] for some function u.

AδPlusX[InMat_]:= Module[
{Blah,i},
Table[
InMat[[i+1]]-InMat[[i]],
{i,1, Length[InMat]-1}
]
]
AδMinusX[InMat_]:= Module[
{Blah,i},
Table[
InMat[[i]]-InMat[[i-1]],
{i,2, Length[InMat]}
]
]

Or (probably better and measurably faster)

δPlusX[InMat_]:= Module[
{Blah},
RotateLeft[InMat]-InMat
]
δMinusX[InMat_]:= Module[
{Blah,i},
InMat-RotateRight[InMat]
]

What about the y-derivatives

δPlusY[InMat_]:= Module[
{Temp},
Temp=Transpose[InMat];
Transpose[δPlusX[Temp]]
]
δMinusY[InMat_]:= Module[
{Temp},
Temp=Transpose[InMat];
Transpose[δMinusX[Temp]]
]

TδPlusY[InMat_]:= Module[
{Blah},
Transpose[δPlusX[Transpose[InMat]]]
]
TδMinusY[InMat_]:= Module[
{Blah},
Transpose[δMinusX[Transpose[InMat]]]
]

U=Table[u[x,y],{x,0, 10,0.01},{y,2,3,0.02}];
Dimensions[U]
Timing[AδPlusX[U];]
Timing[δPlusX[U];]
Timing[δPlusY[U];]
Timing[TδPlusY[U];]

Of course we can use matrices. (first the version that does not wrap around!)

δPlusXMat=SparseArray[
Flatten[
Table[ {
{i,i}->1,
{i,i+1}->-1
},
{i, Length[U]-1}]
]
]
Dimensions[U]
Timing[Dimensions[δPlusXMat.U]]
NonSparseδ=Normal[δPlusXMat];
Timing[Dimensions[NonSparseδ.U]]

The conclusion is that you should not use full matrices for sparse things!    Now a version that does wrap around but without the full matrix computation!

δMinusXMat=SparseArray[
Flatten[
{
{{1,1}->-1,{1,Length[U]}->1},
Table[ {
{i,i-1}->1,
{i,i}->-1
},
{i, 2,Length[U]}]
}
]
]
Dimensions[U]
Timing[Dimensions[δMinusXMat.U]]

You can do the same thing with the y derivatives but there is not much point since it is slower than the "rotate" version.

Sparse vs Full Matrices

The question came up why not just use sparse matrices.

n=1000;
A= Table[RandomReal[],{n},{n}];
AA=SparseArray[A]
B=Table[RandomReal[],{n},{n}];
BB=SparseArray[B]
Timing[A.B;]
Timing[A.BB;]
Timing[AA.B;]
Timing[AA.BB;]

The memory storage is also worse.

ByteCount[A]
ByteCount[AA]

But for genuinely sparse matrices the sparse matrix stuff is a bit faster.

n=7*(10^3);
i=1;
SparseSample=SparseArray[
Flatten[
{
Table[{k,k}->RandomReal[{10, 100}],
{k,1,n}],
Table[
i=RandomInteger[{i,n-1}];
{i, RandomInteger[{i+1,n}]}->RandomReal[{-1,1}],
{100}
]
}
],
{n,n}
]
b= Table[RandomReal[{0,1}],{n}];
Timing[LinearSolve[SparseSample,b];]
Timing[LinearSolve[Normal[SparseSample],b];]