close all
clear all

batchfilename = 'new_batch6.m';
dir = 'Batch6Results\'

ratio = [.25 .5 1 2 4 6 8];
names = {'BuoyA2','BuoyA10','BuoyA20'};
areas = [2,10,20];
wave_heights = [1 3 5];

%ratio = 4;
%names = {'BuoyA10'}
%areas = 10;
%wave_heights = [5];

b = 00:1000:12000;
Lc = 0.4;
tf = 500;
dt = 0.05;

for Hs = wave_heights
    for i = 1:length(names)
      name = names{i}; A = areas(i);
      u_mean  = -ones(length(ratio),length(b));
      u_force_mean  = -ones(length(ratio),length(b));
      pto_mean  = -ones(length(ratio),length(b));
      pow_disp_mean = -ones(length(ratio),length(b));
      prop_power_mean = -ones(length(ratio),length(b));

      for j = 1:length(ratio)
          Ls = ratio(j)*A/Lc
          for k = 1:length(b)
            dtt = dt;
            while (pto_mean(j,k) < 0)
               [u_mean(j,k) u_force_mean(j,k) pto_mean(j,k) pow_disp_mean(j,k) prop_power_mean(j,k)] = Compute_FlappingFoil_Performance(dir,name,Hs,Ls,Lc,dtt,tf,b(k));

               if(pto_mean(j,k) == -1)  %Results filenm not found

                   fid = fopen(batchfilename,'a');
                   fprintf(fid, 'clear all;close all;Hs=%f;Ls=%f;Lc=%f;dt=%f;tf=%f;b=%f;TD_FlappingFoil_%s\n',Hs,Ls,Lc,dtt,tf,b(k),name);
                   fclose(fid);
                   u_mean(j,k) = NaN;
                   pto_mean(j,k) = NaN;
                   pow_disp_mean(j,k) = NaN;
                   prop_power_mean(j,k) = NaN;
                   u_force_mean(j,k) = NaN;
                   break;
               end
               if(isnan(pto_mean(j,k)))  %Results not converged, look for finer timestep results
                   dtt = dtt/2;
                   pto_mean(j,k) = -1;
               end
            end
          end
      end
      
    close all;
    figure('Position', [55 100 1175 844],'PaperOrientation','landscape','PaperPosition',[0 0 11 8.5])
    subplot(4,1,1)
    set(gca,'NextPlot','replacechildren');
    set(gca,'LineStyleOrder',{'-x',':x'});
    plot(ratio,pto_mean)
    title(['FlappingFoil\_' name ' : Hs = ' num2str(Hs)]);
    ylabel('Mean PTO (watts)');
    legend(['b = ' num2str(b(1)) ' N/m/s'],num2str(b(2)),num2str(b(3)),num2str(b(4)),num2str(b(5)),num2str(b(6)),num2str(b(7)),num2str(b(8)),num2str(b(9)),num2str(b(10)),num2str(b(11)),num2str(b(12)),num2str(b(13)));
    subplot(4,1,2)
    set(gca,'NextPlot','replacechildren');
    set(gca,'LineStyleOrder',{'-x',':x'});
    plot(ratio,u_force_mean/4.45)
    ylabel('mean propulsion force (lbs)');
    subplot(4,1,3)
    set(gca,'NextPlot','replacechildren');
    set(gca,'LineStyleOrder',{'-x',':x'});
    plot(ratio,prop_power_mean)
    ylabel('mean propulsion power (watts)');
    subplot(4,1,4)
    set(gca,'NextPlot','replacechildren');
    set(gca,'LineStyleOrder',{'-x',':x'});
    plot(ratio,u_mean)
    ylabel('Mean Forward Speed (m/s)');
    xlabel('Foil to Buoy Area Ratio'); 
    
    print('-dpdf',[dir 'FlappingFoil_' name '_Hs' num2str(Hs) '_Vr'])
    print('-dill',[dir 'FlappingFoil_' name '_Hs' num2str(Hs) '_Vr'])
     
    
    figure('Position', [55 100 1175 844],'PaperOrientation','landscape','PaperPosition',[0 0 11 8.5])
    subplot(4,1,1)
    set(gca,'NextPlot','replacechildren');
    set(gca,'LineStyleOrder',{'-x',':x'});
    plot(b,pto_mean')
    title(['FlappingFoil\_' name ' : Hs = ' num2str(Hs)]);
    ylabel('Mean PTO (watts)');
    legend(['ratio = ' num2str(ratio(1))],['ratio = ' num2str(ratio(2))],['ratio = ' num2str(ratio(3))],['ratio = ' num2str(ratio(4))],['ratio = ' num2str(ratio(5))],['ratio = ' num2str(ratio(6))],['ratio = ' num2str(ratio(7))]);
    subplot(4,1,2)
    set(gca,'NextPlot','replacechildren');
    set(gca,'LineStyleOrder',{'-x',':x'});
    plot(b,u_force_mean'/4.45)
    ylabel('mean propulsion force (lbs)');
    subplot(4,1,3)
    set(gca,'NextPlot','replacechildren');
    set(gca,'LineStyleOrder',{'-x',':x'});
    plot(b,prop_power_mean')
    ylabel('mean propulsion power (watts)');
    subplot(4,1,4)
    set(gca,'NextPlot','replacechildren');
    set(gca,'LineStyleOrder',{'-x',':x'});
    plot(b,u_mean')
    ylabel('Mean Forward Speed (m/s)');
    xlabel('damping coeff N/(m/s)'); 
    print('-dpdf',[dir 'FlappingFoil_' name '_Hs' num2str(Hs) '_Vb'])
    print('-dill',[dir 'FlappingFoil_' name '_Hs' num2str(Hs) '_Vb'])
    end
end

