Formatting Double Sided Axis Issues / Compiling Formulas on One Plot Problem

Hello All,

I am working on an assignment and I am currently stuck with several issues that may be simple, however I am just not seeing it.

The problem is as follows:

We want to understand the fine details of the relationships between the various properties of a material and its orientation. This is best done by developing your own code. Using MATLAB, recreate a figure similar to Figure 2.11 in Jones using the material properties listed below. They are representative of a boron/epoxy composite. It is ok to create a separate plot for Poisson's ratio. Please download and examine the attached (the code I have written is exactly how it is supposed to be structured) for expectations and hints on how to complete MATLAB assignments for this course.

• E1 = 45.0 x 106 psi
• E2 = 5.0 x 106 psi
• G12 = 1.5 x 106 psi
• v12 = 0.3

These are the equations we have to implement:

This is what the graph should resemble:

I am having difficulties with formatting and there may possibly be an error in the functions themselves. The formulas should be correct, but my graph does not look like the book at all. I would also like to have different colors for the lines and make it neater. The given MATLAB file is structured as requested. Please do not suggest removing the for loop, the linspace feature, or how the functions are at the bottom because that is how it is supposed to be. Any help is appreciated!

Here is my code:

% The purpose of this code is to

% (1) understand the fine details of the relationships between the 

% various properties of a material and its orientation.
%start by closing all existing windows and plots that are open and clearing
%all variables
close all;
clear all;
% Material Name
E1 = 45*10^6;
E2 = 5.0*10^6;
G12 = 1.5*10^6;
v12 = .3;
v21 = v12*(E2/E1);
figure %create a new figure but hold off showing it
hold on;
theta = linspace(0,90,18);
Ex = zeros(size(theta));
Vxy = zeros(size(theta));
Gxy = zeros(size(theta));
for i = 1: length(theta)
    theta(i)
    [Ex(i)] = Ex_Longitudinal(E1, E2, G12, v12, theta(i)); 
    [Gxy] = Gxy_Longitudinal(E1, E2, G12, v12, theta(i));
    [nxy] = Nxy_Longitudinal(E1, E2, G12, v12, theta(i));
    [Vxy] = Vxy_Longitudinal(E1, E2, G12, v12, theta(i));
    yyaxis left;
    pLeft = plot(theta, Ex/E2, '-.r*', theta, Gxy/G12, '-mo', theta, -nxy, 'g+');
    yyaxis right;
    pRight = plot(theta, Vxy, 'bs');
end
xlabel ('\theta (degree)')
yyaxis left
ylabel ('E_x/E_2, Gxy/G12, -\etaxy,x')
yyaxis right
ylabel ('\nuxy');
h = [pLeft;pRight];
legend (h, 'E_x/E_2', 'Gxy/G12', '-\etaxy,x', '\nuxy');
hold off;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%Put your User Defined functions at the end of the script - that is a
% MATLAB requirement
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [Ex] = Ex_Longitudinal(E1, E2, G12, v12, theta)
    E_inverse = (1/E1)*(cosd(theta)^4) + (1/G12 - 2*v12/E1) * (sind(theta)^2) * (cosd(theta)^2) + (1/E2) * (sind(theta)^4);
    Ex = 1/E_inverse;
end

function [Vxy] = Vxy_Longitudinal(E1, E2, G12, v12, theta)
    E_inverse = (1/E1)*(cosd(theta)^4) + (1/G12 - 2*v12/E1) * (sind(theta)^2) * (cosd(theta)^2) + (1/E2) * (sind(theta)^4);
    Ex = 1/E_inverse;
    Vxy = Ex * ((v12/E1) * ((sind(theta))^4 + (cosd(theta))^4) - ((1/E1) + (1/E2) - (1/G12)) * (sind(theta))^2 * (cosd(theta))^2);
    %Vxy = Ex*((v12/E1)*((sind(theta))^4+(cosd(theta))^4)-((1/E1)+(1/E2)-(1/G12))*(sind(theta))^2*(cosd(theta))^2);
end

function [Gxy] = Gxy_Longitudinal(E1, E2, G12, v12, theta)
    Gxy_inverse = 2 * (2/E1 + 2/E2 + (4 * v12)/E1 - 1/G12) * sind(theta)^2 * cosd(theta)^2 + (1/G12) * (sind(theta)^4 + cosd(theta)^4);
    Gxy = 1/Gxy_inverse;
end

function [nxy] = Nxy_Longitudinal(E1, E2, G12, v12, theta)
    E_inverse = (1/E1)*(cosd(theta)^4) + (1/G12 - 2*v12/E1) * (sind(theta)^2) * (cosd(theta)^2) + (1/E2) * (sind(theta)^4);
    Ex = 1/E_inverse;
    nxy = Ex * ((2/E1 + (2 * v12)/E1 - 1/G12) * sind(theta) * cosd(theta)^3 - (2/E2 + (2 * v12)/E1 - 1/G12) * sind(theta)^3 + cosd(theta));
end

[dpb]

Would be a better teaching/learning experience to use MATLAB as intended and to vectorize the functions instead of looping. Is a disservice to the student otherwise imo.

However, your basic problem to begin with is in

for i = 1: length(theta)
    [Ex(i)] = Ex_Longitudinal(E1, E2, G12, v12, theta(i)); 
    [Gxy] = Gxy_Longitudinal(E1, E2, G12, v12, theta(i));
    [nxy] = Nxy_Longitudinal(E1, E2, G12, v12, theta(i));
    [Vxy] = Vxy_Longitudinal(E1, E2, G12, v12, theta(i));

NB: what is different from first expression that the rest????

BTW, the extra and superfluous [] do nothing and are expensive computationally. Write just

Ex(i)=Ex_Longitudinal(E1, E2, G12, v12, theta(i)); 

Similarly and as for style, I'd eliminate most if not all of the superfluous () in the equations--they obfuscate more than clarify when aren't needed for grouping terms--

function Ex = Ex_Longitudinal(E1,E2,G12,v12,theta)
  Ex=cosd(theta)^4/E1 + (1/G12-2*v12/E1)*sind(theta)^2*cosd(theta)^2+sind(theta)^4/E2;
  Ex=1/Ex;
end

function Vxy = Vxy_Longitudinal(E1, E2, G12, v12, theta)
  Ex_Longitudinal(E1,E2,G12,v12,theta);             % just wrote a function, use it...
  Vxy=Ex*(v12/E1*(sind(theta)^4 + cosd(theta)^4) - (1/E1+1/E2-1/G12)*sind(theta)^2*cosd(theta)^2);
end