U0 MyDerivative(CMathODE *ode, F64, COrder2D3 *, COrder2D3 *)
{
	//The forces due to springs and drag are
	//automatically handled by the
	//ode code.  We can add new forces
	//here.
	CTask	*task = ode->win_task;
	F64		 d, dd;
	CD3		 p, p2;
	MyMass	*tmpm1, *tmpm2;

	//Collisions
	tmpm1 = ode->next_mass;
	while (tmpm1 != &ode->next_mass)
	{
		tmpm2 = tmpm1->next;
		while (tmpm2 != &ode->next_mass)
		{
			D3Sub(&p, &tmpm2->state->x, &tmpm1->state->x);
			dd = D3NormSqr(&p);
			if (dd <= Sqr(tmpm1->radius + tmpm2->radius)) {
				d = Sqrt(dd) + 0.0001;
				dd = 10.0 * Sqr(Sqr(Sqr(tmpm1->radius + tmpm2->radius) - dd));
				D3MulEqu(&p, dd / d);
				D3AddEqu(&tmpm2->DstateDt->DxDt, &p);
				D3SubEqu(&tmpm1->DstateDt->DxDt, &p);
			}
			tmpm2 = tmpm2->next;
		}
		tmpm1 = tmpm1->next;
	}

	tmpm1 = ode->next_mass;
	while (tmpm1 != &ode->next_mass)
	{
		if (!(tmpm1->flags & MSF_FIXED))
			tmpm1->DstateDt->DyDt += 10.0 * tmpm1->mass; //Gravity
		tmpm1 = tmpm1->next;
	}

	if (cursor_mass)
	{
		p2.x = mouse.pos.x - task->pix_left - task->scroll_x;
		p2.y = mouse.pos.y - task->pix_top  - task->scroll_y;
		p2.z = 0;
		D3Sub(&p, &p2, &cursor_mass->state->x);
		d = 10.0 * D3NormSqr(&p);
		D3MulEqu(&p, d);
		D3AddEqu(&cursor_mass->DstateDt->DxDt, &p);
	}
}
