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);
        }
}