// symmetrical class A amplifier with BJTs
// Vcc=10V, Icq=0.5A, m=1, n=1/2, R=2.5ohm

// initialization
clear
npoint=720; // reduce if the simulation is too slow
i=1:npoint;
wt=2*%pi*(i-0.5)/npoint;
deg=wt*180/%pi;
s=sin(wt);

np=200; // number of frames, reduce if the simulation is too slow

scf(1);

for k=1:np
    
    // computation of voltages and currents
    Vinm(k)=k/np*0.025; // amplitude of the input voltage
    vin=Vinm(k)*s; // amplitude * pattern
    vout=200*vin; // A*vin, A=200 for this amplifier
    ic1=0.5+20*vin; // Icq+gm*m*vin
    ic2=0.5-20*vin; // Icq-gm*m*vin
    icc=ic1+ic2;
    vce1=10-2*vout; // Vcc-vout/n
    vce2=10+2*vout; // Vcc+vout/n
    pd1=ic1.*vce1;
    pd2=ic2.*vce2;
    pcc=10*icc; // Vcc*icc
    pout=vout.*vout/2.5;
    Pcc(k)=mean(pcc);
    Pd1(k)=mean(pd1);
    Pd2(k)=mean(pd2);
    Pout(k)=mean(pout);
    efficiency(k)=Pout(k)/Pcc(k)*100;
    
    
    drawlater()
    
    clf(1);
    
    subplot(4,3,1)
    plot2d(deg,ic1,rect=[0 0 360 1],nax=[0 7 0 5])
    xtitle('','wt [deg]','ic1 [A]')

    subplot(4,3,2)
    plot2d(deg,vce1,rect=[0 0 360 20],nax=[0 7 0 5])
    xtitle('','wt [deg]','vce1 [V]')
    
    subplot(4,3,3)
    plot2d(vce1,ic1,rect=[0 0 20 1],nax=[0 5 0 5],style=0)
    xtitle('','vce1 [V]','ic1 [A]')
    
    subplot(4,3,4)
    plot2d(deg,ic2,rect=[0 0 360 1],nax=[0 7 0 5])
    xtitle('','wt [deg]','ic2 [A]')
    
    subplot(4,3,5)
    plot2d(deg,vce2,rect=[0 0 360 20],nax=[0 7 0 5])
    xtitle('','wt [deg]','vce2 [V]')
    
    subplot(4,3,6)
    plot2d(vce2,ic2,rect=[0 0 20 1],nax=[0 5 0 5],style=0)
    xtitle('','vce2 [V]','ic2 [A]')
        
    subplot(4,3,7)
    plot2d(deg,vin,rect=[0 -0.025 360 0.025],nax=[0 7 0 5])
    xtitle('','wt [deg]','vin [V]')

    subplot(4,3,8)
    plot2d(deg,vout,rect=[0 -5 360 5],nax=[0 7 0 5])
    xtitle('','wt [deg]','vout [V]')
 
    subplot(4,3,9)
    plot2d(deg,pout,rect=[0 0 360 12],nax=[0 7 0 5])
    xtitle('','wt [deg]','pout [W]')
        
    subplot(4,3,10)
    plot2d(deg,pd1,rect=[0 0 360 12],nax=[0 7 0 5])
    xtitle('','wt [deg]','pd1 [W]')

    subplot(4,3,11)
    plot2d(deg,pd2,rect=[0 0 360 12],nax=[0 7 0 5])
    xtitle('','wt [deg]','pd2 [W]')

    subplot(4,3,12)
    plot2d(deg,pcc,rect=[0 0 360 12],nax=[0 7 0 5])
    xtitle('','wt [deg]','pcc [W]')
        
    drawnow()
    
end

scf(2);
clf(2);
plot2d(Vinm,Pcc,rect=[0 0 0.025 12],style=color("green"))
plot2d(Vinm,Pd1,rect=[0 0 0.025 12],style=color("blue"))
plot2d(Vinm,Pd2,rect=[0 0 0.025 12],style=color("magenta"))
plot2d(Vinm,Pout,rect=[0 0 0.025 12],style=color("red"))
plot2d(Vinm,Pd1+Pd2,rect=[0 0 0.025 12],style=color("yellow"))
xtitle('','Vinmm [V]','Pcc [green], Pd1 [blue], Pd2 [magenta], Pout [red], Pd1+Pd2 [yellow], all in [W]')

scf(3);
clf(3);
plot2d(Vinm,efficiency)
xtitle('','Vinm [V]','efficiency [%]')