clc
clear

m1=10;
k1=500000;
ksi1=0.005;
c1=ksi1*2*m1*sqrt(k1/m1);
mu=0.1;
m2=mu*m1;
ome1=2*pi*10;
ome2=2*pi*100;
C11=(ome1+ome2)/2;
omega=0:0.5:2*ome2;

c20=76.2917;
k20=3.9386e+004;
om=236;


C2 =c20-50:c20+50;
K2 =k20-5e4:0.5e4:k20+5e4;                
[c2,k2] = meshgrid(C2,K2); 

for ii=1:length(k2)
    for jj=1:length(c2)
        for iomeg=1:length(omega)
            H11(iomeg)=(-omega(iomeg)^2*m2+k2(ii,jj)+i*omega(iomeg)*c2(ii,jj))./((-omega(iomeg)^2*m2+k2(ii,jj)).*(-omega(iomeg)^2*m1+k1)- ...
                  omega(iomeg)^2*m2.*k2(ii,jj)+i*omega(iomeg).*c2(ii,jj).*(-omega(iomeg)^2*m1+k1-omega(iomeg)^2*m2));
        end
        
        H12(ii,jj)=max(abs(H11));
    end    
end

plot3(c2,k2,H12)

hold on
H110=(-om^2*m2+k20+i*om*c20)./((-om^2*m2+k20).*(-om^2*m1+k1)-om^2*m2.*k20+i*om.*c20.*(-om^2*m1+k1-om^2*m2));
plot3(c20,k20,abs(H110),'go')

xlabel('c2')
ylabel('k2')
zlabel('abs(H11)')
rotate3d
