//example 4.12.
clc
disp("From the given information we can write,")
disp("    A = -16*10^6/j*omega  and  beta = 10^3/[2*10^3+j*omega]^2")
disp("To verify the Barkhausen condition means to verify whether |A*beta| = 1 at a frequency for which A*beta = 0 degree. Let us express, A*beta in its rectangluar form.")
disp("  A*beta = -16*10^6*10^3 / j*omega*[2*10^3+j*omega]^2 = -16*10^9 / j*omega*[4*10^6+4*10^3*j*omega+(j*omega)^2]")
disp("         = -16*10^9 / j*omega*[4*10^6+4*10^3*j*omega-omega^2]      as j*2 = -1")
disp("         = -16*10^9 / 4*10^6*j*omega+4*10^3*j^2*omega^2-j*omega^3]")
disp("         = -16*10^9 / j*omega*[4*10^6-omega^2]-[omega^2*4*10^3]")
disp("Rationalising the denominator function we get,")
disp("  A*beta = -16*10^9[-omega^2*4*10^3 - j*omega*[4*10^6-omega^2]] / [-[omega^2*4*10^3]-j*omega*[4*10^6-omega^2]]*[-omega^2*4*10^3 - j*omega*[4*10^6-omega^2]]")
disp("Using (a-b)(a+b) = a^2 - b^2 in the denominator,")
disp("  A*beta = 16*10^9[omega^2*4*10^3+j*omega*[4*10^6-omega^2]] / [-omega^2*4*10^3]^2 - [j*omega*[4*10^6-omega^2]^2")
disp("  A*beta = 16*10^9[omega^2*4*10^3+j*omega*[4*10^6-omega^2]] / 16*10^6*omega^4 + omega^2(4*10^6-omega^2)^2")
disp("Now to have A*beta = 0 degree, the imaginary part of A*beta must be zero. This is possible when,")
disp("Therefore,  omega*(4*10^6 - omega^2) = 0")
disp("Therefore,  omega = 0  or  4*10^6 - omega^2 = 0")
disp("Therefore,  omega^2 = 4*10^6          Neglecting zero value of frequency")
disp("Therefore,  omega = 2*10^3 rad/sec")
disp("At this frequency |A*beta| can be obtained as,")
disp("  |A*beta| = 16*10^9[4*10^3*omega^2] / 16*10^6*omega^4+omega^2[4*10^6-omega^2]^2       at omega = 2*10^3")
ab=(2.56*10^20)/(2.56*10^20)
disp(ab,"  |A*beta| =")
disp("Therefore,  At omega = 2*10^3 rad/sec, A*beta = 0 degree as imaginary part is zero while |A*beta| = 1. Thus Barkhausen Criterion is satisfied.")
disp("The frequency at which circuit will oscillate is the value of omega for which |A*beta|  = 1 and A*beta = 0 degree at the same time")
disp("i.e.      omega = 2*10^3  rad/sec")
disp("But       omega = 2*pi*f")
f=(2*10^3)/(2*%pi)  // in Hz
format(9)
disp(f,"Therefore,  f(in Hz) = omega / 2pi =")