l = 1;
B = [0.0250; 0.0500; 0.0750; 0.1000; 0.1250; 0.1500; 0.1750; 0.2000; 0.2250; 0.2500; 0.2750; 0.3000; 0.3250; 0.3500; 0.3750; 0.4000];
B = 0.025:0.025:0.4;
kV = nan(1,1,numel(B));
CM = nan(2,4,numel(B));
rig = 0.01;
for j = 1:numel(B)
k = B(j)/(2*rig);
TMP = [cos(k*l)^2, (sin(k*l)*cos(k*l))/k, sin(k*l)*cos(k*l), (sin(k*l)^2)/k;
-k*sin(k*l)*cos(k*l), cos(k*l)^2, -k*sin(k*l)^2, sin(k*l)*cos(k*l);
-sin(k*l)*cos(k*l), -(sin(k*l)^2)/k, cos(k*l)^2, (sin(k*l)*cos(k*l))/k;
k*sin(k*l)^2, -sin(k*l)*cos(k*l), -k*sin(k*l)*cos(k*l), cos(k*l)^2];
TMP = TMP*TMP;
CM(:,:,j) = [TMP(1,:);TMP(3,:)];
kV(j) = k;
end
Best Answer