%=============================================================
% Uma introdução à modelagem quase-estática de automóveis
% Capítulo 6 - Diagramas de desempenho
% Exemplo página 137 da edição de Fevereiro de 2012
%=============================================================
clear all
close all
clc
%-------------------------------------------------------------------------
% dados do veículo
%-------------------------------------------------------------------------
% relações das marchas
ic=[2.909; 1.9776; 1.471; 1.0];
% relação do diferencial
id=3.091;
% rendimento da transmissão
etam=0.90;
% dados dos pneus
rd=0.32;        % raio dinâmico (m)
escorr=0.02;    % escorregamento
f=0.015;        % coeficiente de atrito de rolamento
% carroceria
A=2.0;          % área frontal (m^2)
Cx=0.42;
% peso do veículo (N)
G=16503;
% pressão atmosférica (mmHg)
patm=760;    
% temperatura ambiente (K)
Tamb=15+273.15;
% densidade do ar
ro=0.4647*patm/Tamb;
% aclive
alfa=0*pi/180;
% número de pontos nos gráficos
NP=20;
%% curva de potência do motor
% rotação (rpm)
%nrpm=[1500;  2000;  2500;  3000;  3500;   4000;   4500;   5000;   5500;   6000];
nrpm=linspace(1500,6000,NP);
% potência (W)
%PeWtab= [5560; 33100; 57500; 78500; 96000; 110100; 120800; 128000; 131800; 132100];
PeW=10^(-11)*nrpm.^4-1.7*10^(-7)*nrpm.^3-0.0058*(nrpm.^2)+77*nrpm-96000;
%% cálculo do torque (N.m)
T=PeW./(2*pi*nrpm/60);
%%
% aceleração da gravidade (m/s^2)
grav=9.81;
% número de marchas
NM=size(ic,1);
%-------------------------------------------------------------------------
%% cálculo das curvas de desempenho
%-------------------------------------------------------------------------
for j=1:NM
    for i=1:NP
        % velocidade tangencial
        vt(i,j)=2*pi/60*nrpm(i)/(ic(j)*id)*rd;
        % velocidade real
        vr(i,j)=vt(i,j)*(1-escorr);
        % resistência aerodinâmica
        Qa(i,j)=1/2*ro*vr(i,j)^2*A*Cx;
        % resistência ao rolamento
        Qr(i,j)=f*G*cos(alfa);
        % resistência ao aclive
        Qs(i,j)=G*sin(alfa);
        % potência no cubo
        Pc(i,j)=PeW(i)*etam;
        % potência consumida
        Pcons(i,j)=(Qa(i,j)+Qr(i,j)+Qs(i,j))*vt(i,j);
        % potência líquida
        PL(i,j)=Pc(i,j)-Pcons(i,j);
    end
end
%-------------------------------------------------------------------------
%% captura dos parâmetros
%-------------------------------------------------------------------------
%% aclives máximos
maxPLvr=max(PL./vr);
alfamax = asin(maxPLvr*(1-escorr)/G);
aclivemax = tan(alfamax)*100;
%% aceleração máxima
% inércia equivalente
delta=0.004+0.05*ic.^2;
% aceleração máxima em cada marcha
acelmax=maxPLvr*(1-escorr)./(1+delta)'/G*grav;
acelmaxG=acelmax/grav;
%% velocidades limites
vrmin=vr(1,:);
vrmax=vr(end,:);

%-------------------------------------------------------------------------
%% gráfico de torque e potência do motor
%-------------------------------------------------------------------------
figure('name','Torque e potência do motor')
[AX,H1,H2] = plotyy(nrpm,T,nrpm,PeW);
set(get(AX(1),'Ylabel'),'String','T (Nm)') 
set(get(AX(2),'Ylabel'),'String','Pe (W)') 
set(H1,'LineStyle',':')
set(H2,'LineStyle','--')
set(H1,'LineWidth',2)
set(H2,'LineWidth',2)
grid on
xlabel('n (rpm)')
%-------------------------------------------------------------------------
%% gráfico das potências no cubo e perdida
%-------------------------------------------------------------------------
figure('Name','Potências no cubo e consumida')    % título na moldura superior
plot(vr*3.6,Pc,'LineWidth',2)      
hold on
plot(vr*3.6,Pcons,'k','LineWidth',2)   
grid on                           % linhas de grade
xlabel('v (km/h)')                  % rótulo do eixo x
ylabel('P (W)')               % rótulo do eixo y
title('Potência no cubo e potência consumida')         % título do gráfico
xlim([0 Inf])                           % limites mínimo e máximo do eixo x
ylim([0 Inf])                       % limites mínimo e máximo do eixo y
legend('P_c 1^a marcha',...                 % legenda identificadora de cada marcha
       'P_c 2^a marcha',...
       'P_c 3^a marcha',...
       'P_c 4^a marcha',...
       'P_{cons}',4)
%-------------------------------------------------------------------------
%% gráfico da potência líquida
%-------------------------------------------------------------------------
figure('Name','Potência líquida')    % título na moldura superior
plot(vr*3.6,PL,'LineWidth',2)        
grid on                           % linhas de grade
xlabel('v (km/h)')                  % rótulo do eixo x
ylabel('P_L (W)')               % rótulo do eixo y
title('Potência líquida')         % título do gráfico
legend('1^a marcha',...                 % legenda identificadora de cada marcha
       '2^a marcha',...
       '3^a marcha',...
       '4^a marcha',1)
hold on
plot(vr*3.6,PL,'ok','LineWidth',1)
%-------------------------------------------------------------------------
%% gráfico das relações da caixa
%-------------------------------------------------------------------------
figure('Name','Relações')
plot(vr(:,1)*3.6,nrpm,vr(:,2)*3.6,nrpm,vr(:,3)*3.6,nrpm,vr(:,4)*3.6,nrpm,...
    'LineWidth',2)
axis([0 inf 0 inf])
xlabel('v (km/h')
ylabel('n (rpm)')
legend('1^a','2^a','3^a','4^a',4)
grid on
title('Relações da caixa')
% ======================================================================== 
%% === F I M ===
% ========================================================================


