Subscript indices must either be real positive integers or logicals.

2 views (last 30 days)
clear;
clc;
% For Graphite-Eboxy
S_1m = 204000; % maximum Sigma 1 (psi)
S_2m = 9050; % maximum Sigma 2 (psi)
S_12m = 11500; % maximum Tau 12 (psi)
St_1m = 0.00888; % maximum normal strain 1
St_2m = 0.00725; % maximum normal strain 2
St_12m = 0.00725; % maximum shear strain 12
% Initial Modulus :
E_1o = 23E6; % (psi)
E_2o = 1.25E6; % (psi)
G_12o = 0.975E6; % (psi)
% Initial Vlaue for Secant Modulus (0.9 * Initial Modulus):
E_1s = E_1o; % (psi)
E_2s = 0.9*E_2o; % (psi)
G_12s = 0.9*G_12o; % (psi)
% Material constant:
B1 = 0; C1 = 1; D1 = 0;
B2 = 0; C2 = 1; D2 = 0;
B3 = 0.05139; C3 = 0.067339; D3 = -0.003119;
%Fiber Orientation:
Theta = 0;
C = cosd(Theta);
S = sind(Theta);
%Total Strain Energy (psi):
U_1m = 0.5 * S_1m * St_1m ;
U_2m = 0.5 * S_2m * St_2m ;
U_12m = 0.5 * S_12m * St_12m ;
% Tau = Tau _xy First trial:
t_xy = 0 ;
S_x = 0;
s_y = 0;
for i = 0:0.01:200 % Ry range
a1 = C^2 + i*S^2;
a2 = S^2 + i*C^2;
a12 = S*C*(i - 1);
t1 = 2 * t_xy * S * C;
t2 = -t1;
t12 = t_xy *( C^2 - S^2);
Error = 1;
U_p_total_o = 0;
S_x =0;
S_y = 0;
while Error > 0.00000001 %To calculate the plastic strain energy
C_1 = 2 * E_1s * U_1m ;
C_2 = 2 * E_2s * U_2m ;
C_12 = 2 * G_12s * U_12m ;
b = (2*a1*t1/C_1) + (2*a2*t2/C_2)+(2*a12*t12/C_12);
a = (a1^2/C_1) + (a2^2/C_2)+(a12^2/C_12);
c = -1 + ((t1^2/C_1)+(t2^2/C_2)+(t12^2/C_12));
d = sqrt (b^2 - 4*a*c) ;
% For Sigma X, Y, and txy
S_x = (-b+d)/(2*a);
S_y = i * S_x ;
% For Sigma 1, 2, and t12
S_1 = S_x * C^2 + S_y * S^2 + t_xy *(2 *C * S);
S_2 = S_x * S^2 + S_y * C^2 - t_xy *(2 * C * S);
t_12 = (-S_x * C * S) + (S_y * C * S) + t_xy * (C^2 - S^2);
% For plastic strain energy:
U_p1 = ((S_1)^2 / 2 )*((E_1o-E_1s)/(E_1o*E_1s));
U_p2 = ((S_2)^2 / 2 )*((E_2o-E_2s)/(E_2o*E_2s));
U_p12 = ((t_12)^2 / 2 )*((G_12o-G_12s)/(G_12o*G_12s));
U_p_total_n = U_p1 + U_p2 + U_p12 ;
%For New Secant modulus:
E_1s = E_1o*(1- B1*(U_p_total_n)^C1 + D1*(U_p_total_n) );
E_2s = E_2o*(1- B2*(U_p_total_n)^C2 + D2*(U_p_total_n) );
G_12s = G_12o*(1- B3*(U_p_total_n)^C3 + D3*(U_p_total_n) );
%Error Calculation:
Error = (( U_p_total_n - U_p_total_o )/ U_p_total_n);
U_p_total_o = U_p_total_n ;
end
S_x(i) = S_x ;
S_y(i) = i * S_x ;
end
Array indices must be positive integers or logical values.
plot (S_x, S_y);
The error came from the last two line ( S_x(i) = S_x; and S_y(i) = S_x * i; ) I'm trying to save all the output result from For Loop but the error: "Subscript indices must either be real positive integers or logicals" and plot it, any help please ??

Accepted Answer

Scott MacKenzie
Scott MacKenzie on 5 Mar 2022
Edited: Scott MacKenzie on 5 Mar 2022
You've got some funny things going on with the S_x and S_y variables. I made a small change: building up the S_x and S_y values in each pass in vectors xResult and yResult. This, I think, is what you are after:
% For Graphite-Eboxy
S_1m = 204000; % maximum Sigma 1 (psi)
S_2m = 9050; % maximum Sigma 2 (psi)
S_12m = 11500; % maximum Tau 12 (psi)
St_1m = 0.00888; % maximum normal strain 1
St_2m = 0.00725; % maximum normal strain 2
St_12m = 0.00725; % maximum shear strain 12
% Initial Modulus :
E_1o = 23E6; % (psi)
E_2o = 1.25E6; % (psi)
G_12o = 0.975E6; % (psi)
% Initial Vlaue for Secant Modulus (0.9 * Initial Modulus):
E_1s = E_1o; % (psi)
E_2s = 0.9*E_2o; % (psi)
G_12s = 0.9*G_12o; % (psi)
% Material constant:
B1 = 0; C1 = 1; D1 = 0;
B2 = 0; C2 = 1; D2 = 0;
B3 = 0.05139; C3 = 0.067339; D3 = -0.003119;
%Fiber Orientation:
Theta = 0;
C = cosd(Theta);
S = sind(Theta);
%Total Strain Energy (psi):
U_1m = 0.5 * S_1m * St_1m ;
U_2m = 0.5 * S_2m * St_2m ;
U_12m = 0.5 * S_12m * St_12m ;
% Tau = Tau _xy First trial:
t_xy = 0 ;
% S_x = 0;
% s_y = 0;
xResult = [];
yResult = [];
for i = 0:0.01:200 % Ry range
a1 = C^2 + i*S^2;
a2 = S^2 + i*C^2;
a12 = S*C*(i - 1);
t1 = 2 * t_xy * S * C;
t2 = -t1;
t12 = t_xy *( C^2 - S^2);
Error = 1;
U_p_total_o = 0;
S_x =0;
S_y = 0;
while Error > 0.00000001 %To calculate the plastic strain energy
C_1 = 2 * E_1s * U_1m ;
C_2 = 2 * E_2s * U_2m ;
C_12 = 2 * G_12s * U_12m ;
b = (2*a1*t1/C_1) + (2*a2*t2/C_2)+(2*a12*t12/C_12);
a = (a1^2/C_1) + (a2^2/C_2)+(a12^2/C_12);
c = -1 + ((t1^2/C_1)+(t2^2/C_2)+(t12^2/C_12));
d = sqrt (b^2 - 4*a*c) ;
% For Sigma X, Y, and txy
S_x = (-b+d)/(2*a);
S_y = i * S_x ;
% For Sigma 1, 2, and t12
S_1 = S_x * C^2 + S_y * S^2 + t_xy *(2 *C * S);
S_2 = S_x * S^2 + S_y * C^2 - t_xy *(2 * C * S);
t_12 = (-S_x * C * S) + (S_y * C * S) + t_xy * (C^2 - S^2);
% For plastic strain energy:
U_p1 = ((S_1)^2 / 2 )*((E_1o-E_1s)/(E_1o*E_1s));
U_p2 = ((S_2)^2 / 2 )*((E_2o-E_2s)/(E_2o*E_2s));
U_p12 = ((t_12)^2 / 2 )*((G_12o-G_12s)/(G_12o*G_12s));
U_p_total_n = U_p1 + U_p2 + U_p12 ;
%For New Secant modulus:
E_1s = E_1o*(1- B1*(U_p_total_n)^C1 + D1*(U_p_total_n) );
E_2s = E_2o*(1- B2*(U_p_total_n)^C2 + D2*(U_p_total_n) );
G_12s = G_12o*(1- B3*(U_p_total_n)^C3 + D3*(U_p_total_n) );
%Error Calculation:
Error = (( U_p_total_n - U_p_total_o )/ U_p_total_n);
U_p_total_o = U_p_total_n ;
end
xResult = [xResult, S_x];
yResult = [yResult, i*S_y];
end
plot (xResult, yResult);

More Answers (2)

AndresVar
AndresVar on 5 Mar 2022
In matlab arrays are indexed starting at 1 so S_x(0) is an error.
You can use another array to store the 'range' values
range = 0:0.1:20;
for idx = 1:length(range)
i = range(idx);
..
S_x(idx) = ..
end

Image Analyst
Image Analyst on 6 Mar 2022

Products


Release

R2017a

Community Treasure Hunt

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

Start Hunting!

Translated by