MATLAB Answers

How to solve this problem ?

2 views (last 30 days)
I have added another program to the first program but I hve obtained this error
Assignment has more non-singleton rhs dimensions than non-singleton subscripts
Error in Try_var (line 237)
numerator(c1,km)= ( j*omega(c1)*fi_r(:,km)*sigma_r(:,km) )/( wr(km).^2
-omega(c1).^2+j^2*zetak(km)*wr(km)*omega(c1));
Please help me
The first program is rectified
L = 0.035;
hs = 0.0005;
hp = 0.0001;
rhop = 7800;
rhos = 8700 ;
b0 = 0.01;
x = linspace(0,L,100);
ratio = 1;
nlambda=11;
b = b0*ratio + (1/L)*b0*(1-ratio)*(L-x);
m_const = (rhos * hs + rhop* hp);
m = b * m_const;
f= @(lambda) 1+(cos(lambda)*cosh(lambda)) ;
lamda=zeros(1,nlambda); % preallocate
for i=1:nlambda
lambda(i) = fzero(f, (pi/2)+pi*(i-1));
end
x = linspace(0,L,100);
psi=zeros(1,nlambda);
sigma_r=zeros(100,nlambda);
for km =1:nlambda
psi(km)= ( -cos (lambda(km))- cosh(lambda(km)) )/ ( -sin(lambda(km))+ sinh(lambda(km)) ) ;
Mode_shape = sin ( (lambda(km)/L)*x)- sinh( (lambda(km)/L)*x)+ psi(km) *( cos( (lambda(km)/L)*x) - cosh( (lambda(km)/L)*x));
sigma_r(:,km) = trapz(x,m.* Mode_shape.');
end
The second program need to be rectified :
n_mode=10;
R_load= 1e2;
w_r = 1678,50389473644 10519,0043516527 29453,5167737505 57717,1375770109 95410,5639715215 142526,861582819 199066,444575967 265029,286979108 340415,398708447 425220,754783481 519536,710722902];
e_S_33 = 25.55e-9;
Ep = 6.7e10 ;
d31 = -247.76e-12;
hpc= 3e-04;
omega = linspace (0,5024, 500)';
sz = zeros(length(omega),n_mode);
tot_V = sz; numerator=sz; denomenator=sz;
V_freq = sz; voltage_freq = sz; voltage_phase_freq = sz;
% position of maximal deflection
[~,pos] = max(abs(Mode_shape(:,1)));
%The slope or the derived function of the mode shape
for km=1:nlambda
slope(:,km) = (lambda(km)/L) * ( cos((lambda(km)/L)*x) - cosh ((lambda(km)/L)*x) + psi(km)* ( -sin((lambda(km)/L)*x) - sinh((lambda(km)/L)*x) ) );
end
for c = 1:length(R_matrix)
% Set current resistance value
R_l = R_matrix(c);
% Time constant
to_c = R_l*e_S_33*b*L / hp;
%
teta = - Ep*b*d31* hpc;
%
for km =1:n_mode
xi_r(:,km) = teta* slope(pos,km);
end
%
zetak = 0,00820668626018796 0,0163274909149453 0,0434769520932623 0,0847227579442571 0,139877591855777 0,208869541213683 0,291680805021625 0,388304307061774 0,498736899482713 0,622971093973019];
%
for km=1:n_mode
fi_r (:,km)= (-d31*Ep*hpc*hp / e_S_33*L) * slope(pos,km);
end
for km=1:n_mode
numerator(c1,km)= ( j*omega(c1)*fi_r(:,km)*sigma_r(:,km) )/( wr(km).^2 -omega(c1).^2+j^2*zetak(km)*wr(km)*omega(c1));
denominator(c1,km) = (j*omega(c1)*xi_r(:,km)*fi_r(:,km)^2)/( wr(km)^2-omega(c1)^2+j^2*zetak(km)*wr(km)*omega(c1));
end
V_freq (c1,:) = sum(numerator(c1,:))/( sum(denominator(c1,:)) + ( (1+ j *omega(c1)* to_c) / to_c) );
voltage_freq(c1,:) = abs(V_freq(c1,:));
voltage_phase_freq(c1,:) = angle(V_freq(c1,:));
end
  9 Comments
Mallouli Marwa
Mallouli Marwa on 1 Jan 2019
The first program is rectified :
L = 0.035;
hs = 0.0005;
hp = 0.0001;
rhop = 7800;
rhos = 8700 ;
b0 = 0.01;
x = linspace(0,L,100);
ratio = 1;
nlambda=11;
b = b0*ratio + (1/L)*b0*(1-ratio)*(L-x);
m_const = (rhos * hs + rhop* hp);
m = b * m_const;
f= @(lambda) 1+(cos(lambda)*cosh(lambda)) ;
lamda=zeros(1,nlambda); % preallocate
for i=1:nlambda
lambda(i) = fzero(f, (pi/2)+pi*(i-1));
end
x = linspace(0,L,100);
psi=zeros(1,nlambda);
sigma_r=zeros(100,nlambda);
for km =1:nlambda
psi(km)= ( -cos (lambda(km))- cosh(lambda(km)) )/ ( -sin(lambda(km))+ sinh(lambda(km)) ) ;
Mode_shape(:,km) = sin ( (lambda(km)/L)*x)- sinh( (lambda(km)/L)*x)+ psi(km) *( cos( (lambda(km)/L)*x) - cosh( (lambda(km)/L)*x));
sigma_r(:,km) = trapz(x,m.* Mode_shape(:,km).');
end
The second program need to be rectified :
R_matrix = R_load;
n_mode=10;
R_load= 1e2;
w_r = 1678,50389473644 10519,0043516527 29453,5167737505 57717,1375770109 95410,5639715215 142526,861582819 199066,444575967 265029,286979108 340415,398708447 425220,754783481 519536,710722902];
e_S_33 = 25.55e-9;
Ep = 6.7e10 ;
d31 = -247.76e-12;
hpc= 3e-04;
omega = linspace (0,5024, 500)';
sz = zeros(length(omega),n_mode);
tot_V = sz; numerator=sz; denomenator=sz;
V_freq = sz; voltage_freq = sz; voltage_phase_freq = sz;
% position of maximal deflection
[~,pos] = max(abs(Mode_shape(:,1)));
%The slope or the derived function of the mode shape
for km=1:nlambda
slope(:,km) = (lambda(km)/L) * ( cos((lambda(km)/L)*x) - cosh ((lambda(km)/L)*x) + psi(km)* ( -sin((lambda(km)/L)*x) - sinh((lambda(km)/L)*x) ) );
end
for c = 1:length(R_matrix)
% Set current resistance value
R_l = R_matrix(c);
% Time constant
to_c = R_l*e_S_33*b*L / hp;
%
teta = - Ep*b*d31* hpc;
%
for km =1:n_mode
xi_r(:,km) = teta* slope(pos,km);
end
%
zetak = 0,00820668626018796 0,0163274909149453 0,0434769520932623 0,0847227579442571 0,139877591855777 0,208869541213683 0,291680805021625 0,388304307061774 0,498736899482713 0,622971093973019];
%
for km=1:n_mode
fi_r (:,km)= (-d31*Ep*hpc*hp / e_S_33*L) * slope(pos,km);
end
for c1 = 1:length(omega)
for km=1:n_mode
numerator(c1,km)= ( j*omega(c1)*fi_r(:,km)*sigma_r(:,km) )/( wr(km).^2 -omega(c1).^2+j^2*zetak(km)*wr(km)*omega(c1));
denominator(c1,km) = (j*omega(c1)*xi_r(:,km)*fi_r(:,km)^2)/( wr(km)^2-omega(c1)^2+j^2*zetak(km)*wr(km)*omega(c1));
end
V_freq (c1,:) = sum(numerator(c1,:))/( sum(denominator(c1,:)) + ( (1+ j *omega(c1)* to_c) / to_c) );
voltage_freq(c1,:) = abs(V_freq(c1,:));
voltage_phase_freq(c1,:) = angle(V_freq(c1,:));
end
% Total voltage
tot_V = voltage_freq * 9.81;
end
end

Sign in to comment.

Accepted Answer

Image Analyst
Image Analyst on 1 Jan 2019
Edited: Image Analyst on 1 Jan 2019
First, your vectors for w_r and zetak are missing left brackets.
Next, R_load is not defined. If it came from loading a mat file, attach the mat file, and all your m-files, with the paper clip icon like madhan already asked you.
Next, is anything in
numerator(c1,km)= ( j*omega(c1)*fi_r(:,km)*sigma_r(:,km) )/( wr(km).^2
-omega(c1).^2+j^2*zetak(km)*wr(km)*omega(c1));
a vector instead of a scalar? Use the debugger to find out. If so, then you might need ./ instead of / or .^ instead of ^ or .* instead of *.
  8 Comments
Stephan
Stephan on 2 Jan 2019
49990 +2

Sign in to comment.

More Answers (0)

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by