Issues with plotting/best method

1 vue (au cours des 30 derniers jours)
Eddy Ramirez le 17 Avr 2021
Modifié(e) : per isakson le 17 Avr 2021
Greetings,
I would like to plot the following, but I get a blank graph and I am not sure if it is the way I am inputting the values
any recommendations? or a different approach?
close all
clear all
clc
%%%FIBER
Ef=220e9;%[N/m] GPA to Newton/square meter
Vf=.63; %fiber volume fraction
vf=.33;%fiber poissions ratio
%%MATRIX
Em=10e9;%[N/m] GPA to Newton/square meter
Vm=1-Vf; %%matrix volume fraction
vm=.33;%%matrix poissons ratio
Em2=Em/(1-vm^2); %equation 3.31
%%Lamina's Thickness
h=1e-3;%[N/m] mm to Newton/square meter
%%Shear modulus
Gf=Ef/(2*(1+vf));
Gm=Em/(2*(1+vm));
%%Lamina Properties
E1=(Vf*Ef)+(Vm*Em);
E2=(Ef*Em)/(Vf*Gm+Vm*Gf);
E_2=Ef*Em2/(Vf*Em2+Vm*Ef); %equation 3.32
v12=(Vf*vf)+(Vm*vm);
v21=(E2/E1)*v12;
G12=(Gf*Gm)/((Vf*Gm)+(Vm*Gf));
% Reduced local in plane stiffness Q
Q11=E1/(1-v12*v21);
Q12=(v21*E1)/(1-v12*v21);
Q22=E2/(1-v12*v21);
Q66=G12;
Q=[Q11,Q12,0;
Q12,Q22,0;
0,0,Q66];
%%Laminate rotations in degrees
theta=[0,45,-45,90,30,60,-30,-60,-60,-30,60,30,90,-45,45,0];
% theta=[theta;theta(end:-1:1)];
z=(-length(theta)/2+(0:length(theta)))'*h;
% Global Q matrix:
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for i=1:length(theta)
m=cosd(theta(i));
n=sind(theta(i));
T=[m.^2,n.^2,2.*m.*n;
n.^2,m.^2,-2.*m.*n;
-m.*n,m.*n,m.^2-n.^2];
% Global transformed reduced stiffness coefficients Q_bar
Q_bar{i}=inv(T).*Q.*T;
end
Aij=0;
Bij=0;
Dij=0;
for j=1:length(theta)
A_j=Q_bar{j}*(z(j+1)-z(j));
Aj{j}=A_j;
Aij=Aij+Aj{j};
B_j=Q_bar{j}*((z(j+1))^2-(z(j))^2);
Bj{j}=B_j;
Bij=Bij+Bj{j};
D_j=Q_bar{j}*((z(j+1))^3-(z(j))^3);
Dj{j}=D_j;
Dij=Dij+Dj{j};
end
A=Aij;
B=(1/2)*Bij;
D=(1/3)*Dij;
%%Forces
Nx=2e6;
Ny=4.6e6;
Ns=0;
N=[Nx;Ny;Ns];
%%Moments
Mx=3e3;
My=0;
Ms=-1e-3;
M=[Mx;My;Ms];
%%Shear
e_x=sym('Epsilonx');
e_y=sym('Epsilony');
gamma_xy=sym('Gammaxy');
shear=[e_x;e_y; gamma_xy];
%%Bending Twist
k_x=sym('Kx');
k_y=sym('Ky');
k_xy=sym('Kxy');
bending=[k_x;k_y;k_xy];
%%Shear Extension Coupling
SEC=N==A*shear;
SEC_F=solve(SEC);
e_xf=vpa(SEC_F.Epsilonx);
e_yf=vpa(SEC_F.Epsilony);
gamma_xyf=(SEC_F.Gammaxy);
shearf=[e_xf;e_yf; gamma_xyf];
%%Bending
BEC=M==D*bending;
BEC_F=solve(BEC);
k_xf=vpa(BEC_F.Kx);
k_yf=vpa(BEC_F.Ky);
k_xyf=vpa(BEC_F.Kxy);
bendingf=[k_xf;k_yf;k_xyf];
test=Q_bar{1}*shearf+z(1,1)*Q_bar{1}*bendingf;
test1=Q_bar{1}.*shearf+z(2,1)*Q_bar{1}.*bendingf;
% plot(test(1,1), z(1,1));
for k=1:length(theta)
stress1=Q_bar{k}*shearf+z(k)*Q_bar{k}*bendingf;%%BOTTOM LAYER
stress2=Q_bar{k}*shearf+z(k+1)*Q_bar{k}*bendingf;%%TOP LAYER
stress1f{k}=stress1; %Bottom
stress2f{k}=stress2; %Top
end
b0=stress1f{1};
b1=stress1f{2};
b2=stress1f{3};
b3=stress1f{4};
b4=stress1f{5};
b5=stress1f{6};
b6=stress1f{1};
b7=stress1f{7};
b8=stress1f{8};
b9=stress1f{9};
b10=stress1f{10};
b11=stress1f{11};
b12=stress1f{12};
b13=stress1f{13};
b14=stress1f{14};
b15=stress1f{15};
t0=stress2f{1};
t1=stress2f{2};
t2=stress2f{3};
t3=stress2f{4};
t4=stress2f{5};
t5=stress2f{6};
t6=stress2f{1};
t7=stress2f{7};
t8=stress2f{8};
t9=stress2f{9};
t10=stress2f{10};
t11=stress2f{11};
t12=stress2f{12};
t13=stress2f{13};
t14=stress2f{14};
t15=stress2f{15};
plot(b0(1,1), z(1,1))
keyboard
hold on
plot(t0(1,1), z(2,1))
plot(b1(1,1), z(2,1))
plot(t1(1,1), z(3,1))
plot(b2(1,1), z(3,1))
plot(t2(1,1), z(4,1))
plot(b3(1,1), z(4,1))
plot(t3(1,1), z(5,1))
plot(b4(1,1), z(5,1))
plot(t4(1,1), z(6,1))
plot(b5(1,1), z(6,1))
plot(t5(1,1), z(7,1))
plot(b6(1,1), z(7,1))
plot(t6(1,1), z(8,1))
plot(b7(1,1), z(8,1))
plot(t7(1,1), z(9,1))
plot(b8(1,1), z(9,1))
plot(t8(1,1), z(10,1))
plot(b9(1,1), z(10,1))
plot(t9(1,1), z(11,1))
plot(b10(1,1), z(11,1))
plot(t10(1,1), z(12,1))
plot(b11(1,1), z(12,1))
plot(t11(1,1), z(13,1))
plot(b12(1,1), z(13,1))
plot(t12(1,1), z(14,1))
plot(b13(1,1), z(14,1))
plot(t13(1,1), z(15,1))
plot(b14(1,1), z(15,1))
plot(t14(1,1), z(16,1))
plot(b15(1,1), z(16,1))
plot(t15(1,1), z(17,1))
hold off
keyboard
8 commentairesAfficher 6 commentaires plus anciensMasquer 6 commentaires plus anciens
Eddy Ramirez le 17 Avr 2021
should I run a stairs instead of a plot?
Clayton Gotberg le 17 Avr 2021
To connect the dots, just plot all of the points at once in the order you want them to be connected, using a line.
If I have three points (x_a, y_a), (x_b,y_b), and (x_c,y_c), I can plot them separately or together:
%Separately
plot(x_a, y_a,'k.')
plot(x_b, y_b,'k.')
plot(x_c, y_c,'k.')
%Together
plot([x_a x_b x_c],[y_a y_b y_c],'k-')
% Or
X = [x_a x_b x_c];
Y = [y_a y_b y_c];
plot(X,Y,'k-')

Connectez-vous pour commenter.

Réponses (2)

LO le 17 Avr 2021
short answer: you are plotting points singularly so you can't see them
short solution: pool them together in arrays then use the command to plot
arrays can be made using the line
[a,b,c,d,.....,z]
where letters are your values for bottom or top
I am not familiar with the sym data format but I think the issue might be just in the way you plot there. Also comment the keyboard otherwise it gets stuck in debugging.
see your code with scatter instead of plot (scatter can plot single points, plot can too but if you plot a single points one by one you will not see them, they have to be organized in arrays)
close all
clear all
clc
%%%FIBER
Ef=220e9;%[N/m] GPA to Newton/square meter
Vf=.63; %fiber volume fraction
vf=.33;%fiber poissions ratio
%%MATRIX
Em=10e9;%[N/m] GPA to Newton/square meter
Vm=1-Vf; %%matrix volume fraction
vm=.33;%%matrix poissons ratio
Em2=Em/(1-vm^2); %equation 3.31
%%Lamina's Thickness
h=1e-3;%[N/m] mm to Newton/square meter
%%Shear modulus
Gf=Ef/(2*(1+vf));
Gm=Em/(2*(1+vm));
%%Lamina Properties
E1=(Vf*Ef)+(Vm*Em);
E2=(Ef*Em)/(Vf*Gm+Vm*Gf);
E_2=Ef*Em2/(Vf*Em2+Vm*Ef); %equation 3.32
v12=(Vf*vf)+(Vm*vm);
v21=(E2/E1)*v12;
G12=(Gf*Gm)/((Vf*Gm)+(Vm*Gf));
% Reduced local in plane stiffness Q
Q11=E1/(1-v12*v21);
Q12=(v21*E1)/(1-v12*v21);
Q22=E2/(1-v12*v21);
Q66=G12;
Q=[Q11,Q12,0;
Q12,Q22,0;
0,0,Q66];
%%Laminate rotations in degrees
theta=[0,45,-45,90,30,60,-30,-60,-60,-30,60,30,90,-45,45,0];
% theta=[theta;theta(end:-1:1)];
z=(-length(theta)/2+(0:length(theta)))'*h;
% Global Q matrix:
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for i=1:length(theta)
m=cosd(theta(i));
n=sind(theta(i));
T=[m.^2,n.^2,2.*m.*n;
n.^2,m.^2,-2.*m.*n;
-m.*n,m.*n,m.^2-n.^2];
% Global transformed reduced stiffness coefficients Q_bar
Q_bar{i}=inv(T).*Q.*T;
end
Aij=0;
Bij=0;
Dij=0;
for j=1:length(theta)
A_j=Q_bar{j}*(z(j+1)-z(j));
Aj{j}=A_j;
Aij=Aij+Aj{j};
B_j=Q_bar{j}*((z(j+1))^2-(z(j))^2);
Bj{j}=B_j;
Bij=Bij+Bj{j};
D_j=Q_bar{j}*((z(j+1))^3-(z(j))^3);
Dj{j}=D_j;
Dij=Dij+Dj{j};
end
A=Aij;
B=(1/2)*Bij;
D=(1/3)*Dij;
%%Forces
Nx=2e6;
Ny=4.6e6;
Ns=0;
N=[Nx;Ny;Ns];
%%Moments
Mx=3e3;
My=0;
Ms=-1e-3;
M=[Mx;My;Ms];
%%Shear
e_x=sym('Epsilonx');
e_y=sym('Epsilony');
gamma_xy=sym('Gammaxy');
shear=[e_x;e_y; gamma_xy];
%%Bending Twist
k_x=sym('Kx');
k_y=sym('Ky');
k_xy=sym('Kxy');
bending=[k_x;k_y;k_xy];
%%Shear Extension Coupling
SEC=N==A*shear;
SEC_F=solve(SEC);
e_xf=vpa(SEC_F.Epsilonx);
e_yf=vpa(SEC_F.Epsilony);
gamma_xyf=(SEC_F.Gammaxy);
shearf=[e_xf;e_yf; gamma_xyf];
%%Bending
BEC=M==D*bending;
BEC_F=solve(BEC);
k_xf=vpa(BEC_F.Kx);
k_yf=vpa(BEC_F.Ky);
k_xyf=vpa(BEC_F.Kxy);
bendingf=[k_xf;k_yf;k_xyf];
test=Q_bar{1}*shearf+z(1,1)*Q_bar{1}*bendingf;
test1=Q_bar{1}.*shearf+z(2,1)*Q_bar{1}.*bendingf;
% scatter(test(1,1), z(1,1));
for k=1:length(theta)
stress1=Q_bar{k}*shearf+z(k)*Q_bar{k}*bendingf; %%BOTTOM LAYER
stress2=Q_bar{k}*shearf+z(k+1)*Q_bar{k}*bendingf; %%TOP LAYER
stress1f{k}=stress1; %Bottom
stress2f{k}=stress2; %Top
end
b0=stress1f{1};
b1=stress1f{2};
b2=stress1f{3};
b3=stress1f{4};
b4=stress1f{5};
b5=stress1f{6};
b6=stress1f{1};
b7=stress1f{7};
b8=stress1f{8};
b9=stress1f{9};
b10=stress1f{10};
b11=stress1f{11};
b12=stress1f{12};
b13=stress1f{13};
b14=stress1f{14};
b15=stress1f{15};
t0=stress2f{1};
t1=stress2f{2};
t2=stress2f{3};
t3=stress2f{4};
t4=stress2f{5};
t5=stress2f{6};
t6=stress2f{1};
t7=stress2f{7};
t8=stress2f{8};
t9=stress2f{9};
t10=stress2f{10};
t11=stress2f{11};
t12=stress2f{12};
t13=stress2f{13};
t14=stress2f{14};
t15=stress2f{15};
scatter(b0(1,1), z(1,1))
% keyboard
hold on
scatter(t0(1,1), z(2,1))
scatter(b1(1,1), z(2,1))
scatter(t1(1,1), z(3,1))
scatter(b2(1,1), z(3,1))
scatter(t2(1,1), z(4,1))
scatter(b3(1,1), z(4,1))
scatter(t3(1,1), z(5,1))
scatter(b4(1,1), z(5,1))
scatter(t4(1,1), z(6,1))
scatter(b5(1,1), z(6,1))
scatter(t5(1,1), z(7,1))
scatter(b6(1,1), z(7,1))
scatter(t6(1,1), z(8,1))
scatter(b7(1,1), z(8,1))
scatter(t7(1,1), z(9,1))
scatter(b8(1,1), z(9,1))
scatter(t8(1,1), z(10,1))
scatter(b9(1,1), z(10,1))
scatter(t9(1,1), z(11,1))
scatter(b10(1,1), z(11,1))
scatter(t10(1,1), z(12,1))
scatter(b11(1,1), z(12,1))
scatter(t11(1,1), z(13,1))
scatter(b12(1,1), z(13,1))
scatter(t12(1,1), z(14,1))
scatter(b13(1,1), z(14,1))
scatter(t13(1,1), z(15,1))
scatter(b14(1,1), z(15,1))
scatter(t14(1,1), z(16,1))
scatter(b15(1,1), z(16,1))
scatter(t15(1,1), z(17,1))
hold off
% keyboard
0 commentairesAfficher -2 commentaires plus anciensMasquer -2 commentaires plus anciens

Connectez-vous pour commenter.

per isakson le 17 Avr 2021
plot(b0(1,1), z(1,1), 'd' )
keyboard
hold on
plot(t0(1,1), z(2,1), 'd' )
plot(b1(1,1), z(2,1), 'd' )
plot(t1(1,1), z(3,1), 'd' )
plot(b2(1,1), z(3,1), 'd' )
etc.
Now the script outputs
3 commentairesAfficher 1 commentaire plus ancienMasquer 1 commentaire plus ancien
Eddy Ramirez le 17 Avr 2021
thank you for the response, I was able to get the dots but I am looking for a stairs plot scenario (pic attached)
I tried to run something like this, but it is not even close
figure(1)
stairs(b0(1,1), z(1,1), 'LineWidth',2,'Marker','d','MarkerFaceColor','c')
hold on
stairs(t0(1,1), z(2,1),'LineWidth',2,'Marker','d','MarkerFaceColor','c')
stairs(b1(1,1), z(2,1), 'LineWidth',2,'Marker','d','MarkerFaceColor','c')
stairs(t1(1,1), z(3,1), 'LineWidth',2,'Marker','d','MarkerFaceColor','c')
hold off
per isakson le 17 Avr 2021
Modifié(e) : per isakson le 17 Avr 2021
To draw a line you need a vector of at least two elements. That's true for both plot and stairs. I joined only a subset of your points.
x = [t0(1,1),b1(1,1),t1(1,1),b2(1,1),t2(1,1),b3(1,1),t3(1,1),b4(1,1),t4(1,1),b5(1,1),t5(1,1),b6(1,1),t6(1,1),b7(1,1),t7(1,1)];
y = [z(2,1),z(2,1),z(3,1),z(3,1),z(4,1),z(4,1),z(5,1),z(5,1),z(6,1),z(6,1),z(7,1),z(7,1),z(8,1),z(8,1),z(9,1)];
stairs( x, y, '-dk' )
outputs
and replace '-dk' by '-k' to remove the markers
Not all of the lines of unnamed.png are vertical or horizontal. plot( x, y, '-dk' ) produces

Connectez-vous pour commenter.

Catégories

En savoir plus sur Loops and Conditional Statements dans Help Center et File Exchange

Community Treasure Hunt

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

Start Hunting!

Translated by