Last night I was thinking again about my Sine Approximation routine and holy crud, I must have been really tired not to recognize this:

In the approximation, I was taking the upper 8 bits, truncating propagation. This was so that I had a fixed-point routine. But if I had taken the lower 8 bits, I would have had a useful 8-bit integer routine for -x²-x.

But wait.

If I start the accumulator with x (which is actually 3 cycles faster to do), then I end up with -x².

Well then, anybody wanna optimize this? ^^ :
L_sqrd:
;Input: L
;Output: L*L->A
;159 t-states
;39 bytes
ld h,l
ld c,l
rr h
sbc a,a
xor l
add a,c
ld c,a
rl l

rr h
sbc a,a
xor l
and %11111000
add a,c
ld c,a
rl l

rr h
sbc a,a
xor l
and %11100000
add a,c
ld c,a
rl l

ld a,h
rrca
xor l
and 128
xor c
neg
ret

It might be more easy if I explain that it is computing:
=(hhhhhhhh^abcdefgh)+(gggggg00^bcdefg00)+(ffff0000^cdef0000)+(ee000000^de000000)
=(hhhhhhhh^abcdefgh)+(ggggg000^bcdef000)+(fff00000^cde00000)+(e0000000^d0000000)
=((hhhhhhhh^abcdefgh)+(ggggg000^bcdef000)+(fff00000^cde00000))^e0000000^d0000000

^ is XOR
+ is addition modulo 256
EDIT: Cross-posted:
I thought of a way to optimize the first iteration. Saved 2 bytes, 8 t-states.
Basically, at the first iteration, C=-1 or C=2L, then I shift L up one bit for the next iteration. I got rid of the initial ld c,a and use the first iteration to compute c. To do this, I just shift L at the beginning, then after the "sbc a,a" I just OR that with L. If a was $FF, the result is FF (which is -1), else it is 2*L:
L_sqrd:
;Input: L
;Output: L*L->A
;151 t-states
;37 bytes
ld h,l
;First iteration, get the lowest 3 bits
sla l
rrh
sbc a,a
or l
;second iteration, get the next 2 bits
ld c,a
rr h
sbc a,a
xor l
and $F8
add a,c
;third iteration, get the next 2 bits
ld c,a
sla l
rr h
sbc a,a
xor l
and $E0
add a,c
;fourth iteration, get the last bit
ld c,a
ld a,l
add a,a
rrc h
xor h
and $80
xor c
neg
ret
Also, as a note, if you stop at any iteration and perform NEG, you can mask to get the lower bits of the square. So for example, the first iteration gives you -L*L mod 8->A, the second returns -L*L mod 32->A, the third gives it mod 128.

Yep, that is what I am using. By the way, I am just running through the program 3 times and timing it with the given rand seeds. (I am doing it all in one program, though)

[discussion in IRC]
[time passes]
Here is my program, it just does one at a time for competition and it spits out a new flake at each iteration with probability 1/2.

It generates 2 each cycle by default.
There are 128 positions.
This will never happen

And in case you wanted to do a number that doesn't evenly divide into 128, this piece of the code allows exiting as soon as 128 flakes are on screen:
E+1→E
If E=128
S→B
That exits the For() loop to prevent generating any more snow flakes, then if you look at the main Repeat loop conditions:
Repeat E>=128 or getKey=45
This exits if [Clear] is pressed or E is 128 or larger. E counts how many flakes have been drawn.

Uh, yup, I forgot to add one line:
Ans→L2(Ans
This should be inserted before (or after, doesn't really matter) the line redrawing the new location of the snowflake.

From Source Coder:
{16,8->dim([A]
DelVar EE[A]→[A]
ClrHome
2→S
Repeat E>=128 or getKey=45
For(B,1,S
Repeat not([A](Ans,1
randInt(1,16
End
Ans→[A](Ans,1
Output(1,Ans,"*
E+1→E
If E=128
S→B
End
Matr>list([A],8,L2
For(A,7,1,-1
Matr>list([A],A,L1
L1*(L1 and not(L2→L3
L1-Ans→L1
While max(L3
max(L3→B
0→L3(Ans
0→[A](B,A
Output(A,B,"
min(16,max(1,B+1-2int(2rand
If L2(Ans
B
Ans→[A](Ans,A+1
Ans→L2(Ans
Output(A+1,Ans,"*
End
L1→L2
End
End

EDIT: Forgot a line.

Mine draws the flakes when they are newly added at the top of the screen, then it searches the matrix for snowflakes that can move down. It erases those ones and draws their new position.

Also, DJ_O, you didn't need the "1=" part of those If statements.
The way I test for flakes that can move is I store the matrix rotated so that column 1 is the first row of snow flakes. Then each snow flake is represented by the column number it is in (so if it is homescreen column 6, it has a 6 in the matrix). My code is then:
Matr>list([A],8,L2 ;get the last row
For(A,7,1,-1
Matr>list([A],A,L1
L1*(L1 and not(L2→L3 ;This checks if anything in L1 can move down into L2
;"not(L2" leaves a 0 where there is already a flake, else 1 if empty
;Then "L1 and not(L2" leaves a 1 if there is a snowflake above an empty space
L1-Ans→L1 ;remove the snowflakes that can move down from L1
While max(L3
max(L3→B ;location of the snowflake furthest to the right
0→L3(Ans ;remove the snowflake from L3, the list of moveable flakes
0→[A](B,A ;remove it from the matrix
Output(A,B," ;Erase it
min(16,max(1,B+1-2int(2rand ;randomly move left/right
If L2(Ans ;check if the space is occupied
B
Ans→[A](Ans,A+1 ;write the snowflake to the new coordinate
Output(A+1,Ans,"*
End
L1→L2 ;now This row becomes the new lower row
End


Judging by your screenshot, your program would actually run a lot faster with a matrix, since it seems like you are only keeping track of one snowflake at a time. I didn't realize that until just now as I was about to post my attempt Mine handles multiple flakes at a time

EDIT: Oops, now I see that it does handle multiple flakes. I am too sleepy.