Root finding and plotting

How do I find the three roots of " f "

I got the plot of " f." And because the " f " and X asis has three intersections, it would have three roots.

I would like to know how to find the three roots at the same time and plot the three roots on my original plot.

Below is my coding.

Thank you very much!!

clc
clear
format long
v=0.3; 
E=209e+3; 
G=E/(2*(1+v)); 
q=-1;
h=15;
D=(E*h^3)/(12*(1-v^2)); 
I=(h^3)/12; 
a=600;b=2400; 

n =3;
[T1, T2] = meshgrid(1:2:n);
mn = [T1(:), T2(:)]

syms x y

x_value=301;
y_value=0:1:2400;

len=length(mn);
amn=zeros(1,len);

for  i=1:len
    m=mn(i,1);
    n=mn(i,2);
    amn(i)=(16*q/(m*n*D*pi^6))*(1/((m/a)^2+(n/b)^2)^2); 
end

wmn=sym(zeros(1,len));
wxx=sym(zeros(1,len));
wyy=sym(zeros(1,len));
for i=1:len
    m=mn(i,1);
    n=mn(i,2);
    wmn(i)=amn(i).*((sin(m.*pi.*x./a)).*(sin(n.*pi.*y./b))); 
    wxx(i)=diff(wmn(i),x,2);
    wyy(i)=diff(wmn(i),y,2);
end

combine_wxx=sum(wxx);
combine_wyy=sum(wyy);

My=-D*(combine_wyy+v*combine_wxx);
Differentiation_My=diff(My,y,1)
f=subs(Differentiation_My,x,301)
%fsolve(f,500)

Differentiation_My_value=double(subs(Differentiation_My,{x,y},{x_value,y_value}));
plot(Differentiation_My_value)
grid on

Answer 1

[Star Strider]

One way is to search for the indices nearest the zero-crossings, then interpolate (if necessary, to get more exact values) —

format long g
v=0.3; 
E=209e+3; 
G=E/(2*(1+v)); 
q=-1;
h=15;
D=(E*h^3)/(12*(1-v^2)); 
I=(h^3)/12; 
a=600;b=2400; 

n =3;
[T1, T2] = meshgrid(1:2:n);
mn = [T1(:), T2(:)]

mn = 4×2

1 1 1 3 3 1 3 3

syms x y

x_value=301;
y_value=0:1:2400;

len=length(mn);
amn=zeros(1,len);

for  i=1:len
    m=mn(i,1);
    n=mn(i,2);
    amn(i)=(16*q/(m*n*D*pi^6))*(1/((m/a)^2+(n/b)^2)^2); 
end

wmn=sym(zeros(1,len));
wxx=sym(zeros(1,len));
wyy=sym(zeros(1,len));
for i=1:len
    m=mn(i,1);
    n=mn(i,2);
    wmn(i)=amn(i).*((sin(m.*pi.*x./a)).*(sin(n.*pi.*y./b))); 
    wxx(i)=diff(wmn(i),x,2);
    wyy(i)=diff(wmn(i),y,2);
end

combine_wxx=sum(wxx);
combine_wyy=sum(wyy);

My=-D*(combine_wyy+v*combine_wxx);
Differentiation_My=diff(My,y,1)

Differentiation_My = 

f=subs(Differentiation_My,x,301)

f = 

%fsolve(f,500)

Differentiation_My_value=double(subs(Differentiation_My,{x,y},{x_value,y_value}));

zxi = find(diff(sign(Differentiation_My_value)))                            % Zero-Crossing Indices

zxi = 1×4

583 1200 1201 1818

for k = 1:numel(zxi)
    idxrng = max(zxi(k)-2,1):min(zxi(k)+2,numel(Differentiation_My_value)); % Index Range For Interpolation
    zc(k) = interp1(Differentiation_My_value(idxrng),idxrng,0);             % Interpolate
end
Zero_Crossings = zc

Zero_Crossings = 1×4

583.376183230044 1201 1201 1818.62381676996

plot(Differentiation_My_value)
hold on
plot(zc, zeros(size(zc)), 'rs')
hold off
grid on

To plot it against an ‘x’ vector instead of the indices, the ‘zc’ calculation becomes:

% zc(k) = interp1(Differentiation_My_value(idxrng), x(idxrng), 0);

I commented it here because it would otherwise execute in the online Run feature, and throw an error.

.

[Mark]

Hi Star Strider, Thanks for your explanation! One more thing I would like to ask is that why there are 4 Indices at the calculation process of "Zero-Crossing Indices" because it should only 3 Indices for this case.

[Star Strider]

As always, my pleasure!

One of them is a near-duplicate (1200, 1201). I have no idea how that arises unless there is ‘noise’ or a duplicated value in that region that causes a second sign reversal. (I did not look closely at the data.) Using the unique (or uniquetol) function would remove it.

[Mark]

Great! I took your suggestion using the "unique" and it remove the repeat data.

Could you breifly explain in the for loop % Index Range For Interpolation

idxrng = max(zxi(k)-2,1):min(zxi(k)+2,numel(Differentiation_My_value));

I cannot really understand the meaning of the line coding above.

Also, about :To plot it against an ‘x’ vector instead of the indices, the ‘zc’ calculation becomes:

I tried to replace "idxrng" with "x(idxrng)" as your comment did but it indicated an error:Index exceeds the number of array elements (1).

How should I fix it if I Want to plot it against an ‘x’ vector instead of the indices.

Thanks you very much!!

[Star Strider]

Sure!

This assignment:

idxrng = max(zxi(k)-2,1):min(zxi(k)+2,numel(Differentiation_My_value));

sets the beginning of the index range to the greater of 1 or the index value - 2, so that it does not address indices less than 1. The second part does the same at the opposite end, taking the lesser of the index value + 2 or the number of elements in the vector, so that it does not address indices greater than the number of elements in the vector. It is simply ‘crash-proofing’ the indexing range. (It is not strictly necessasry in this instance, however experience has taught me to always include it.)

‘How should I fix it if I Want to plot it against an ‘x’ vector instead of the indices.’

Since ‘x’ would have to be the same length as ‘Differentiation_My_value’, one way would be:

xstart = 0;
xend = 1;
x = linspace(xstart, xend, numel(Differentiation_My_value));

with ‘xstart’ and ‘xend’ being any values you want (within the ability of MATLAB to represent them).

As always, my pleasure!

.

[Mark]

Thanks, Star Strider! I got it! Really appreciate your explanation!!

[Star Strider]

As always, my pleasure!

If you have any further questions or concerns, please post them. I’ll do my best ro respond.

[Mark]

Sounds Great! It so nice of you!! I am a matlab beginner user and sometimes I may meet problems when I do the simulation work for my graduate school study. I will try my best to deal with it by myshelf at first and if I still cannot deal with it I will let you know my questions. Really thanks a lot : )

[Star Strider]

I’ve definitely been there (although by now a few decades removed).

As always, my pleasure!

Answer 2

[Sulaymon Eshkabilov]

There are a couple of ways to find roots of this eqn:

(1) ginput() --> graphical method, e.g.:

[Roots, y] = ginput(3)  % click on three crossing points

(2) fsolve()

[Mark]

Hi Sulaymon Eshkabilov,

Really thanks for your information!