arrays - MATLAB and Forward Euler Method -

i'm writing code uses forward euler method calculate following formula:

dv(t)/dt = g - (c_drag)(p_air/p_water)(3v(t)^2 / 4*d )

i've outlined various constants in code below.

clear;  clc;  close all;  % 1. initialize dt =  input('enter time step'); % choose delta t t_start = 0; % starting time t_end = 100; g = 9.81; y(1) = 0; p_air = 1.2; p_water = 1000; c_drag = 0.517; d = 0.002; nstep = ( t_end - t_start ) / dt;% number of time steps y = 1; % y(0) = 1, initial condition all_t = ( 1 : nstep ) * dt; % time points all_fe_y = zeros( nstep , 1 ); % approximate solution - pre-allocation all_exact_y = zeros( nstep , 1 ); % exact solution - pre-allocation     % 2. forward euler time integration     n = 1 : nstep       y(n+1) = y(n) + (dt * (g - (c_drag * (p_air / p_water) * ((3*(y(dt))^2)/ (4*d))))); % compute solution @ time-step n+1  all_fe_y (n) = y(n+1); % store solution     end     % 3. compute exact solution     n = 1 : nstep      t = n * dt;      all_exact_y (n) = exp(-t);     end     % 4. write approximate solution on file     % 4.1. open file     output_id = fopen('output_forward_euler.txt', 'w');     % 4.2. write header     fprintf(output_id, '%10s %10s %10s \r\n', 'time', 'approx', 'exact' );     % 4.3. write output     fprintf(output_id, '%10.4f %10.6f %10.6f \r\n', [all_t; all_fe_y'; all_exact_y'] );     % 4.4. close file     fclose(output_id);     % 5. plot results     % 5.1 plot settings     set(0 ,'defaulttextinterpreter' , 'latex' );     set(0 ,'defaultaxesunits' , 'pixels');     myfontsize = 14;     myfontname = 'helvetica' ;     % screensize property value four-element vector: [left bottom width height]     figure('position', [100 100 500 500]);     hold on     % 5.2 plot exact , approximate solution     plot(all_t,all_exact_y,'b', 'linewidth',1)     plot(all_t,all_fe_y, '--r','linewidth',0.5)     % 5.3 more plot settings     set(gca,'box' , 'on' ,...      'color' , 'white' ,...      'linewidth' , 0.5 ,...      'tickdir' , 'in' ,...      'xminortick' , 'on' , ...      'yminortick' , 'on' , ...      'xgrid' , 'on' , ...      'ygrid' , 'on' , ...      'xminorgrid' , 'off' , ...      'yminorgrid' , 'off' , ...      'xlim' , [ 0 4] , ...      'ylim' , [ 0 1] , ...      'xtick' ,linspace(0,4,5),...      'ytick' ,linspace(0,1,6),...      'ticklength' , [.02 .02] ,...      'fontname' , myfontname ,...      'fontsize' , myfontsize );     % 5.4 background color     set(gcf, 'color', 'white');     % 5.5 labels     xlabel('time (s)', 'fontname', myfontname, 'fontsize', myfontsize)     ylabel('y (m)', 'fontname', myfontname, 'fontsize', myfontsize)     % 5.6 legend     hleg1 = legend('exact','approximation');     set(hleg1,'location','northeast') 

i'm not getting output expect code..does pertain euler method or part of function? , appreciated!