Subscript indices must either be real positive integers or logicals.

4 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

Categories

Find more on MATLAB in Help Center and File Exchange

Products


Release

R2017a

Community Treasure Hunt

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

Start Hunting!