Took another wack at some semi-Lagrangian Navier-Stokes code thats I abaondoned some time back. Basically, it worked, but I never bothered to make the LHS symmetric. Instead, being lazy, I just solved with BiCGSTAB_l.
So, after a few really stupid errors, it now works again, solving with conjugate gradients. Unfortunatly, the bloody thing now takes forever to converge. Arg! Guess this means that I have to write preconditioning code. Oof.
It is some consolation that now, at least, each iteration is about 3x as fast