Define the γ vector (‘gammav’ here), use a loop, and see the documentation section on Passing Extra Parameters.
tspan = [0 50];
y0 = [2 0]';
gammav = logspace(-1, 1, 5)
for k = 1:numel(gammav)
[t{k},y{k}] = ode45(@(t,y)hw1(t,y,gammav(k)),tspan,y0);
end
figure
hold on
for k = 1:numel(gammav)
plot(t{k},y{k}(:,1))
end
grid
legend(compose('$\\gamma = %.2f$',gammav), 'Location','best', 'Interpreter','latex')
function dy = hw1(t,y,gamma)
dy(1)= y(2);
dy(2)= -gamma*y(2)-y(1);
dy = [dy(1);dy(2)];
end
.