THMMY.gr

function dx=odefun(t,x,flag)

if x<= 21
    dx=30-x;
    
elseif x>=27
    dx=-x;
   
elseif x>= 21 && x<=27 && flag==0
    dx=-x;
elseif x>= 21 && x<=27 && flag==1
    dx=30-x;
end

end

[a b]=ode45(@(t,x)odefun(t,x,flag),[0 10],16);

function heater
    global my_flag
    
    my_flag = 1;
    options = odeset('MaxStep',0.1);
    [t, degrees] = ode45(@odefun,[0 10],16,options);
    
    plot(t,degrees);
    grid
    xlabel('t (sec)')
    ylabel('degrees')
end

function dx = odefun(~,x)
    global my_flag
    
    if x >= 27
        my_flag = 0;
    elseif x < 21
        my_flag = 1;
    end
    dx = 30*my_flag - x;
end

function heater_2
    
    t_start = 0;
    t_final = 10;
    x_start = 16;
    options = odeset('Events',@events,'MaxStep',0.1);
    t = [];
    degrees = [];
    odefun = @odefun_on;

    while true
        [t_temp, degrees_temp, t_ev, degree_ev] = ode45(odefun,[t_start t_final],x_start,options);
        if degree_ev >= 27
            odefun = @odefun_off;
        else
            odefun = @odefun_on;
        end
        t_start = t_ev;
        x_start = degree_ev;
        t = [t ; t_temp];
        degrees = [degrees ; degrees_temp];
        if t(end) >= t_final
            break;
        end
    end
    
    plot(t,degrees);
    grid
    xlabel('t (sec)')
    ylabel('degrees')
end

function dx = odefun_on(~,x)
    dx = 30 - x;
end

function dx = odefun_off(~,x)
    dx = -x;
end

function [ value, isterminal, direction ] = events( ~, x )
    value(1) = x - 21;
    value(2) = 27 - x;
    isterminal(1) = 1;
    isterminal(2) = 1;
    direction(1) = -1;
    direction(2) = -1;
end