constant FMT = " %7.3f"

procedure ivp_euler(integer f, atom y, integer step, integer end_t)
integer t = 0;

    printf(1, " Step %2d: ", step);
    while t<=end_t do
        if remainder(t,10)==0 then printf(1, FMT, y) end if
        y += step * call_func(f,{t, y});
        t += step
    end while
    printf(1, "\n");
end procedure

procedure analytic()
    printf(1, "    Time: ");
    for t = 0 to 100 by 10 do printf(1," %7g", t) end for
    printf(1, "\nAnalytic: ");
    for t = 0 to 100 by 10 do
        printf(1, FMT, 20 + 80 * exp(-0.07 * t))
    end for
    printf(1,"\n");
end procedure

function cooling(atom /*t*/, atom temp)
    return -0.07 * (temp - 20);
end function
constant r_cooling = routine_id("cooling")

analytic();
ivp_euler(r_cooling, 100, 2, 100);
ivp_euler(r_cooling, 100, 5, 100);
ivp_euler(r_cooling, 100, 10, 100);