So I am having trouble with a very basic physics system that I am implementing.
There is a sphere that essentially behaves like a single point in 3D space, and there is a heightmap terrain. I want the sphere to be affected by gravity, so that it moves down hills and resists movement up hills, and also be able to bounce when it hits the ground, moving off with a damped velocity.

Right now I just implemented a standard acceleration/velocity/position model for the sphere with a constant downwards acceleration for gravity, and during each update I reflect the current velocity vector across the normal vector of the point on the ground (unless the sphere is in the air). This results in pretty decent gravity effects, but the bouncing isn't the best. Is there any better way to do this stuff?

Also, what's the best way to dampen the velocity after bouncing by a certain elasticity? Simply multiplying the velocity by the elasticity doesn't produce a good result. I've had some success with reflecting the velocity and then velocity -= normal * velocity.Length() * elasticity.

Well the problem with that, for example, if it it moving horizontally and it bounces on a parallel surface, then it will slow down horizontally as well as vertically, which doesn't seem right when you're bouncing down a hill to gain velocity. Subtracting by the normal is working well for me right now though.

I don't fully understand the reflection part, maybe you're doing it right, dunno. Anyway, what should happen on an impact between a ball and a surface with some normal is:

1. Split the ball's velocity vector into the normal and tangencial subvectors. The normal velocity v[sub]n[/sub] is parallel to the surface normal (you get its magnitude using a dot product between the velocity and the normalized surface normal and its direction is equal to the surface normal) and tangencial velocity v[sub]t[/sub] is given as v[sub]t[/sub] = v - v[sub]n[/sub] (a vector operation!).
2. You do nothing with v[sub]t[/sub] (it stays as it was).
3. You replace v[sub]n[/sub] with an inverse of v[sub]n[/sub] (opposite direction). When the impact in not fully elastic, you also have to modify the magnitude of v[sub]n[/sub]. Notice that we are applying damping (elasticity) only to v[sub]n[/sub], not the the whole velocity v. This IMHO could be the reason of your strange behaviour. Damping can be something like 0.8 and then this whole operation is for example v[sub]n[/sub] *= -0.8 (change direction and magnitude at the same time).
4. You make the new total velocity by simply adding the subvectors: v = v[sub]n[/sub] + v[sub]t[/sub].
5. If you want to implement also rolling, you'll use also the v[sub]t[/sub] component, but this is a bit more complicated.

Yes, Tom, I started using that system and it works very well, thanks for the info.
As for the rolling, the tangential velocity adter the reflection handles that well enough already.

The only issue I am getting is that at elasticities greater than 0.5 the bouncing never completely dampens, but that can be dealt with using some extra processing.

So how are you doing it by reflecting the velocity vector by the normal and then applying the damping, or the way I showed? As I said, there is a difference, because in my solution you dampen only the normal component of the velocity and thus the angle of impact is different than the angle of rebound. I'm a bit confused because you mentioned reflection again.

About never completely dampening - that's expected because dampening is multiplication and not substitution, you can keep mutliplying a number by 0.5 forever and you still won't have perfect zero
In practise, you'll get to zero either because of float variable resolution or (preferably) you use some terminator which checks for a limiting velocity or energy value.

Ah I see what you are describing in your method now.
Yeah what I do is first reflect the entire velocity vector across the normal, and then dampen the normal component of this new velocity. It results in both pretty good bouncing and rolling/sliding due to gravity.