Hello,
I've never had this problem before, but the last thing preventing me from using a program I wrote is that the term 1 that I've isolated gives the wrong answer for some unknown reason. It's line 56 where I just made the variable "term 1" used in line 54 in the real formula. I put in the correct number at the end of the Function_M into wolframalpha to find the correct plot and it just seems to not give the correct answer when function = 0.
clear;clc;%Input Values
%Shear_max=560*(10^6); %yield stress from data sheet
g_max_A=172.369*(10^6); %max cyclic shear stress austentite MPa (25ksi)
g_max_M=68.9476*(10^6); %max cyclic shear stress martensite MPa (10ksi)
G_M=10.8*(10^9); %Shear Modulus of the cold Material
G_A=28.8*(10^9); %Shear Modulus of the hot Material
v=0.33; %Poisson's ratio
%Wire and Coil Diameter
d_max=0.058444; %Stroke Displacement (increased by 1/(%usable strain))
d=150*(10^-6); %The diameter of the Wire (microns converted to m) cools in 1 second with mesh geometry
F_M=0.9807/2; %100 grams in newtons / 2 for 50 grams martensite pulling force
F_A=0.9807/2*(g_max_A/g_max_M); %Force adjusted for austenite max stress at same D
%w=[1 2 4 6 8 10 12 14 16 18 20 22 24 26 28 30];
w=10;D=g_max_M*pi*d^3.*w./(8*F_M); %Spring diameter using martensite numbers
%D_A=g_max_A*pi*d^3.*w./(8*F_A);
D_Mandrel=D-d;g_eight=0.71;g_twelve=1;d_net=0;c=1;length_f=0.12; angle_A_f=sym('angle_A_f'); angle_M_f=sym('angle_M_f'); Function_A=sym('Function_A'); Function_M=sym('Function_M');%Iterating final length to find a diameter for force
%while (length_f >= 0.120)
% length_f=length_f-0.001;
angle_i=40;while length_f >= 0.100 angle_i=angle_i-1; %iterated initial angle
angle_A_f=sym('angle_A_f'); angle_M_f=sym('angle_M_f');% Function_A=0;
Function_A=G_A*d^4*pi/(8*D^2)*... cosd(angle_i)^2*... (sind(angle_A_f)-sind(angle_i))/... (cosd(angle_A_f)^2*... (cosd(angle_A_f)^2+sind(angle_A_f)^2/... (1+v)))-... (F_A/w); %angle_A_f=double(angle_A_f);
% if angle_i==1
angle_A_f=double(vpasolve(Function_A==0, angle_A_f, [0 90]));% else
%angle_A_f=double(vpasolve(Function_A==0, angle_A_f, [0 90]));
% end
Function_M=G_M*d^4*pi/(8*D^2).*cosd(angle_i)^2*(sind(angle_M_f)-sind(angle_i))/... (cosd(angle_M_f)^2*... (cosd(angle_M_f)^2+(sind(angle_M_f)^2/... (1+v))))-(pi*d^3/(8*D)*G_M*0.06*g_eight) -... (F_M/w); term1=-(G_M.*g_eight); term1_2=(pi*d^3)/(8*D); term2=-(F_M/w); angle_M_f=double(vpasolve(Function_M==0, angle_M_f, [0 90])); n=cosd(angle_i)*d_max/(pi*D.*(sind(angle_M_f)-sind(angle_A_f))); length_f=n*pi*D*(tand(angle_M_f)); length_i=n*pi*D*(tand(angle_i)); %pitch_f=tand(angle_f)*pi*D; %Coil pitch
%Length of unsprung spring set to forearm
%length_i=length_f-dl; %Final length = initial + displacement
%n=length_f/pitch_f; %Number of loops in coil
%angle_i=atand(length_i/(pi*n*D));
%Shear Strain and Detwinned Martensite Fraction
% angle_A_f0 = 0;
% for c=1:length(w)
% angle_A_f = fzero(@(angle_A_f) findangle_A(angle_A_f,G_A,d,D_A,angle_i,v,w,c)...
% ,angle_A_f0);
%
% angle_M_f0 = 0;
% angle_M_f = fzero(@(angle_M_f) findangle_M(angle_M_f,G_M,d,D_M,angle_i,v,g_eight,w,c)...
% ,angle_M_f0);
% for c=1:length(w)% end if angle_i <= 3 break endendHere are some pictures of the formula used on line 51 and proof from wolframalpha that the numbers are different.
Formula at angle_i=5:
0=10.8^9*(0.15/1000)^4/(8*0.00186^2)*(cos(5 degrees)^2)*(sin(x degrees)-(sin(5 degrees)))/(sin(x degrees)^2*(cos(x degrees)^2+sin(x degrees)^2/(1+0.33)))-0.3762
Formula that I'm trying to reproduce:
everything besides yL and ESt is referenced clearly in the code, and yL = 0.06, and ESt = 0.71.
Best Answer