Finding a value and marking it in a graph

Question

Shimon Katzman 2019-12-5

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/495069-finding-a-value-and-marking-it-in-a-graph

评论： VBBV 2022-12-25

Hi everyone,

I am trying to find the value of M when epss=epssh=0.009 for every graph (there are 4 of them) and to mark it on the graph unsucssesfully.

Can anyone help me please with this please?

b=400;  %mm
d=500;  %mm
Asv = [3000 5000 3000 5000]; %mm^2
fckv = [30  30  90  90]; %Mpa
num_fckv = numel(fckv);
epscmv = linspace(0.1, 100, 5000)*1E-3;
num_epscvm = numel(epscmv);
fsolve_opts = optimoptions('fsolve','Display','off');
M = zeros(num_epscvm, num_fckv);
phi = zeros(num_epscvm, num_fckv);
idx = zeros(1, num_fckv);
Mmax = zeros(1, num_fckv);
phiAtMmax = zeros(1, num_fckv);
epsAtMmax = zeros(1, num_fckv);
tic
for k = 1:num_fckv
    fck = fckv(k);
    Ecshah=57000/145*(fck*145)^0.5;  %Mpa
    Es=200000;  %Mpa
    Esh=8500; %Mpa
    fy=500;  %Mpa
    fsu=750; %Mpa
    epssh=0.009;
    epssu=0.075;
    eps0=1.027*10^-7*fck*145+0.00195;
    kshah=0.025*fck*10^3;
    A=Ecshah*eps0/fck;
    P=Esh*((epssu-epssh)/(fsu-fy));
    epsy=fy/Es;
    As = Asv(k);
    c0 = 1000;
    c = zeros(1, num_epscvm);
    for i=1:num_epscvm
        epscm = epscmv(i);
        epss=@(c) (d-c)/c*epscm;
        funCshah=@(epsc) (1-(1-epsc./eps0).^A) .* (epsc<=eps0) + exp(-kshah*(epsc-eps0).^1.15) .* (epsc>eps0);
        compression=@(c) b*fck*c/epscm*integral(funCshah,0,epscm)/1000;
        sigmaSteel=@(c) Es*epss(c) .* (epss(c)<=epsy) + fy .* (epss(c)>epsy & epss(c)<=epssh) + (fsu+(fy-fsu)*abs((epssu-epss(c))./(epssu-epssh)).^(1/P)) .* (epss(c)>epssh & epss(c)<=epssu) + 0 .* (epss(c)>epssu);
        tension=@(c) sigmaSteel(c).*As/1000;
        c(i)=fsolve(@(c) compression(c)-tension(c), c0, fsolve_opts);
        c0 = c(i);
        funM=@(epsc) (1-(1-epsc./eps0).^A).*(d-c(i)+(c(i)./epscm).*epsc) .* (epsc<=eps0) + exp(-kshah*(epsc-eps0).^1.15).*(d-c(i)+(c(i)./epscm).*epsc) .* (epsc>eps0);
        M(i,k)=b*fck*c(i)/epscm*integral(funM,0,epscm)/1000000;
        phi(i,k)=epscm/c(i);
        epss1(k,i)= (d-c(i))/c(i)*epscm;
        sigmaSteel1(k,i)= Es*epss1(k,i) .* (epss1(k,i)<=epsy) + fy .* (epss1(k,i)>epsy & epss1(k,i)<=epssh) + (fsu+(fy-fsu)*abs((epssu-epss1(k,i))./(epssu-epssh)).^(1/P)) .* (epss1(k,i)>epssh & epss1(k,i)<=epssu) + 0 .* (epss1(k,i)>epssu);
        c_mtx(i,k) = c(i);
    end
    idx(k) = find(diff(M(:,k)) < 0, 1, 'first')
    [Mmax(k),idx(k)]=max(M(:,k)) %[kNm]
    phiAtMmax(k)=phi(idx(k),k)  %[1/mm]
    epsAtMmax(k)=epscmv(idx(k))*1000 %[promil]
end
toc
% epss_c = fsolve(@(c) epss(c)-epssh, 100, fsolve_opts)
% c_idx = find(c >= epss_c, 1, 'first')
% epscm_exact = interp1(c([-5 5]+c_idx), epscmv([-5 5]+c_idx), epss_c)
figure
hold all
hp = zeros(1, num_fckv);
lgdstr = cell(1, num_fckv);
for k = 1:num_fckv
    hp(k) = plot(phi(1:idx(k)+50,k), M(1:idx(k)+50,k));
%     plot([1 1]*epss_c, ylim, ':k')
    lgdstr{k} = sprintf('fck = %2d [Mpa], As = %4d [mm^2]',fckv(k), Asv(k));
    plot(phi(idx(k),k),M(idx(k),k),'r*')  % marking the Mmax
end
hold off
grid on
xlabel('phi  [1/mm]')
ylabel('Moment  [kNm]')
hl = legend(hp,lgdstr);
hl.FontSize = 7;
% text(epss_c, 400, sprintf('\\leftarrow epss_c = %.1f', epss_c), 'HorizontalAlignment','left')

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Star Strider 2019-12-5

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/495069-finding-a-value-and-marking-it-in-a-graph#answer_405019

在 MATLAB Online 中打开

The only way I can see to calculate ‘M’ where ‘epss’ is equal to 0.009 is to duplicate the first nested loop, defining:

epss = @(c) epssh;

If you do that:

for k = 1:numel(fckv)
    fck = fckv(k);
    Ecshah=57000/145*(fck*145)^0.5;  %Mpa
    Es=200000;  %Mpa
    Esh=8500; %Mpa
    fy=500;  %Mpa
    fsu=750; %Mpa
    epssh=0.009;
    epssu=0.075;
    eps0=1.027*10^-7*fck*145+0.00195;
    kshah=0.025*fck*10^3;
    A=Ecshah*eps0/fck;
    P=Esh*((epssu-epssh)/(fsu-fy));
    epsy=fy/Es;
    epscmv = linspace(0.1, 100, 5000)*1E-3;
    As = Asv(k);
    for i=1:numel(epscmv);
        epscm = epscmv(i);
        epss = @(c) epssh;
%         epss=@(c) (d-c)/c*epscm;
        funCshah=@(epsc) (1-(1-epsc./eps0).^A) .* (epsc<=eps0) + exp(-kshah*(epsc-eps0).^1.15) .* (epsc>eps0);
        compression=@(c) b*fck*c/epscm*integral(funCshah,0,epscm)/1000;
        sigmaSteel=@(c) Es*epss(c) .* (epss(c)<=epsy) + fy .* (epss(c)>epsy & epss(c)<=epssh) + (fsu+(fy-fsu)*abs((epssu-epss(c))./(epssu-epssh)).^(1/P)) .* (epss(c)>epssh & epss(c)<=epssu) + 0 .* (epss(c)>epssu);
        tension=@(c) sigmaSteel(c).*As/1000;
        c(i)=fsolve(@(c) compression(c)-tension(c),1000);
        funM=@(epsc) (1-(1-epsc./eps0).^A).*(d-c(i)+(c(i)./epscm).*epsc) .* (epsc<=eps0) + exp(-kshah*(epsc-eps0).^1.15).*(d-c(i)+(c(i)./epscm).*epsc) .* (epsc>eps0);
        M009(i,k)=b*fck*c(i)/epscm*integral(funM,0,epscm)/1000000;
        phi(i,k)=epscm/c(i);
        c_mtx(i,k) = c(i);
    end
end
toc
figure
mesh(M009)
grid on

and also do the surface plot, you can see the result.

I have no idea what to do with this result, so I leave that to you.

11 个评论
显示 9更早的评论隐藏 9更早的评论

Shimon Katzman 2019-12-7

在 MATLAB Online 中打开

Hi Star,

I am trying this but it errors me. Can you please help me with the error?

Thank you very much.

b=400;  %mm
d=500;  %mm
Asv = [3000 5000 3000 5000]; %mm^2
fckv = [30  30  90  90]; %Mpa
num_fckv = numel(fckv);
fsolve_opts = optimoptions('fsolve','Display','off');
cAtMy = zeros(1, num_epscvm);
c0 = 1000;
tic
for k = 1:num_fckv
    fck = fckv(k);
    Ecshah=57000/145*(fck*145)^0.5;  %Mpa
    Es=200000;  %Mpa
    Esh=8500; %Mpa
    fy=500;  %Mpa
    fsu=750; %Mpa
    epssh=0.009;
    epssu=0.075;
    eps0=1.027*10^-7*fck*145+0.00195;
    kshah=0.025*fck*10^3;
    A=Ecshah*eps0/fck;
    P=Esh*((epssu-epssh)/(fsu-fy));
    epsy=fy/Es;
    As = Asv(k);
        %Values at yield
    epssAtMy= epsy;
    epscmAtMy=@(cAtMy) epssAtMy*cAtMy./(d-cAtMy);
    funCshahAtMy=@(epsc) (1-(1-epsc./eps0).^A) .* (epsc<=eps0) + exp(-kshah*(epsc-eps0).^1.15) .* (epsc>eps0);
    compressionAtMy=@(cAtMy) b*fck*cAtMy/epscmAtMy*integral(funCshahAtMy,0,epscmAtMy)/1000;
    sigmaSteelAtMy= fy;
    tensionAtMy= sigmaSteelAtMy*As/1000;
    cAtMy(k)=fsolve(@(cAtMy) compressionAtMy(cAtMy)-tensionAtMy, c0, fsolve_opts)
    c0=cAtMy(k)
    %funM=@(epsc) (1-(1-epsc./eps0).^A).*(d-cAtMy+(cAtMy./epscmAtMy).*epsc) .* (epsc<=eps0) + exp(-kshah*(epsc-eps0).^1.15).*(d-cAtMy+(cAtMy./epscmAtMy).*epsc) .* (epsc>eps0);
    %M(k)=b*fck*cAtMy/epscmAtMy*integral(funM,0,epscmAtMy)/1000000;
        
     
end 
toc

Shimon Katzman 2019-12-7

在 MATLAB Online 中打开

You can delete it.. still doesnt work.

b=400;  %mm
d=500;  %mm
Asv = [3000 5000 3000 5000]; %mm^2
fckv = [30  30  90  90]; %Mpa
num_fckv = numel(fckv);
fsolve_opts = optimoptions('fsolve','Display','off');
c0 = 1000;
tic
for k = 1:num_fckv
    fck = fckv(k);
    Ecshah=57000/145*(fck*145)^0.5;  %Mpa
    Es=200000;  %Mpa
    Esh=8500; %Mpa
    fy=500;  %Mpa
    fsu=750; %Mpa
    epssh=0.009;
    epssu=0.075;
    eps0=1.027*10^-7*fck*145+0.00195;
    kshah=0.025*fck*10^3;
    A=Ecshah*eps0/fck;
    P=Esh*((epssu-epssh)/(fsu-fy));
    epsy=fy/Es;
    As = Asv(k);
        %Values at yield
    epssAtMy= epsy;
    epscmAtMy=@(cAtMy) epssAtMy*cAtMy./(d-cAtMy);
    funCshahAtMy=@(epsc) (1-(1-epsc./eps0).^A) .* (epsc<=eps0) + exp(-kshah*(epsc-eps0).^1.15) .* (epsc>eps0);
    compressionAtMy=@(cAtMy) b*fck*cAtMy/epscmAtMy*integral(funCshahAtMy,0,epscmAtMy)/1000;
    sigmaSteelAtMy= fy;
    tensionAtMy= sigmaSteelAtMy*As/1000;
    cAtMy(k)=fsolve(@(cAtMy) compressionAtMy(cAtMy)-tensionAtMy, c0, fsolve_opts)
    c0=cAtMy(k)
    %funM=@(epsc) (1-(1-epsc./eps0).^A).*(d-cAtMy+(cAtMy./epscmAtMy).*epsc) .* (epsc<=eps0) + exp(-kshah*(epsc-eps0).^1.15).*(d-cAtMy+(cAtMy./epscmAtMy).*epsc) .* (epsc>eps0);
    %M(k)=b*fck*cAtMy/epscmAtMy*integral(funM,0,epscmAtMy)/1000000;
        
     
end 
toc

Shimon Katzman 2019-12-8

在 MATLAB Online 中打开

Hi Star.

1) It gives me back only one value of cAtMy (instead of 4) what is missing?

2) i added the M(k) and again it errors :(((

b=400;  %mm
d=500;  %mm
Asv = [3000 5000 3000 5000]; %mm^2
fckv = [30  30  90  90]; %Mpa
num_fckv = numel(fckv);
fsolve_opts = optimoptions('fsolve','Display','off');
c0 = 100;
tic
for k = 1:num_fckv
    fck = fckv(k);
    Ecshah=57000/145*(fck*145)^0.5;  %Mpa
    Es=200000;  %Mpa
    Esh=8500; %Mpa
    fy=500;  %Mpa
    fsu=750; %Mpa
    epssh=0.009;
    epssu=0.075;
    eps0=1.027*10^-7*fck*145+0.00195;
    kshah=0.025*fck*10^3;
    A=Ecshah*eps0/fck;
    P=Esh*((epssu-epssh)/(fsu-fy));
    epsy=fy/Es;
    As = Asv(k);
        %Values at yield
    epssAtMy= epsy;
    epscmAtMy=@(cAtMy) epssAtMy*cAtMy./(d-cAtMy);
    funCshahAtMy=@(epsc) (1-(1-epsc./eps0).^A) .* (epsc<=eps0) + exp(-kshah*(epsc-eps0).^1.15) .* (epsc>eps0);
    compressionAtMy=@(cAtMy) b*fck*cAtMy/epscmAtMy(cAtMy)*integral(funCshahAtMy,0,epscmAtMy(cAtMy))/1000;
    %compressionAtMy=@(cAtMy) b*fck*cAtMy/epscmAtMy*integral(funCshahAtMy,0,epscmAtMy)/1000;
    sigmaSteelAtMy= fy;
    tensionAtMy= sigmaSteelAtMy*As/1000;
    cAtMy(k)=fsolve(@(cAtMy) compressionAtMy(cAtMy)-tensionAtMy, c0, fsolve_opts);
    c0=cAtMy(k);
    funM=@(epsc) (1-(1-epsc./eps0).^A).*(d-cAtMy+(cAtMy./epscmAtMy).*epsc) .* (epsc<=eps0) + exp(-kshah*(epsc-eps0).^1.15).*(d-cAtMy+(cAtMy./epscmAtMy).*epsc) .* (epsc>eps0);
    M(k)=b*fck*cAtMy(k)/epscmAtMy*integral(funM,0,epscmAtMy)/1000000
          
end 
toc

Shimon Katzman 2019-12-8

在 MATLAB Online 中打开

b=400;  %mm
d=500;  %mm
Asv = [3000 5000 3000 5000]; %mm^2
fckv = [30  30  90  90]; %Mpa
num_fckv = numel(fckv);
fsolve_opts = optimoptions('fsolve','Display','off');
c0 = 100;
tic
for k = 1:num_fckv
    fck = fckv(k);
    Ecshah=57000/145*(fck*145)^0.5;  %Mpa
    Es=200000;  %Mpa
    Esh=8500; %Mpa
    fy=500;  %Mpa
    fsu=750; %Mpa
    epssh=0.009;
    epssu=0.075;
    eps0=1.027*10^-7*fck*145+0.00195;
    kshah=0.025*fck*10^3;
    A=Ecshah*eps0/fck;
    P=Esh*((epssu-epssh)/(fsu-fy));
    epsy=fy/Es;
    As = Asv(k);
        %Values at yield
    epssAtMy= epsy;
    epscmAtMy=@(cAtMy) epssAtMy*cAtMy./(d-cAtMy);
    funCshahAtMy=@(epsc) (1-(1-epsc./eps0).^A) .* (epsc<=eps0) + exp(-kshah*(epsc-eps0).^1.15) .* (epsc>eps0);
    compressionAtMy=@(cAtMy) b*fck*cAtMy/epscmAtMy(cAtMy)*integral(funCshahAtMy,0,epscmAtMy(cAtMy))/1000;
    %compressionAtMy=@(cAtMy) b*fck*cAtMy/epscmAtMy*integral(funCshahAtMy,0,epscmAtMy)/1000;
    sigmaSteelAtMy= fy;
    tensionAtMy= sigmaSteelAtMy*As/1000;
    cAtMy(k)=fsolve(@(cAtMy) compressionAtMy(cAtMy)-tensionAtMy, c0, fsolve_opts)
    c0=cAtMy(k);
    funM=@(epsc) (1-(1-epsc./eps0).^A).*(d-cAtMy(k)+(cAtMy(k)./epscmAtMy).*epsc) .* (epsc<=eps0) + exp(-kshah*(epsc-eps0).^1.15).*(d-cAtMy(k)+(cAtMy(k)./epscmAtMy).*epsc) .* (epsc>eps0);
    M(k)=b*fck*cAtMy(k)/epscmAtMy*integral(funM,0,epscmAtMy)/1000000
          
end 
toc

Whatever i do it errors me :

Undefined operator '/' for input arguments of type 'function_handle'.
Error in Untitled (line 39)
    M(k)=b*fck*cAtMy(k)/epscmAtMy*integral(funM,0,epscmAtMy)/1000000

Shimon Katzman 2019-12-8

I know but I still cant find what's wrong :(

VBBV 2022-12-25

在 MATLAB Online 中打开

b=400;  %mm
d=500;  %mm
Asv = [3000 5000 3000 5000]; %mm^2
fckv = [30  30  90  90]; %Mpa
num_fckv = numel(fckv);
fsolve_opts = optimoptions('fsolve','Display','off');
c0 = 100;
tic
for k = 1:num_fckv
    fck = fckv(k);
    Ecshah=57000/145*(fck*145)^0.5;  %Mpa
    Es=200000;  %Mpa
    Esh=8500; %Mpa
    fy=500;  %Mpa
    fsu=750; %Mpa
    epssh=0.009;
    epssu=0.075;
    eps0=1.027*10^-7*fck*145+0.00195;
    kshah=0.025*fck*10^3;
    A=Ecshah*eps0/fck;
    P=Esh*((epssu-epssh)/(fsu-fy));
    epsy=fy/Es;
    As = Asv(k);
        %Values at yield
    epssAtMy= epsy;
    epscmAtMy=@(cAtMy) epssAtMy*cAtMy./(d-cAtMy);
    funCshahAtMy=@(epsc) (1-(1-epsc./eps0).^A) .* (epsc<=eps0) + exp(-kshah*(epsc-eps0).^1.15) .* (epsc>eps0);
    compressionAtMy=@(cAtMy) b*fck*cAtMy/epscmAtMy(cAtMy)*integral(funCshahAtMy,0,epscmAtMy(cAtMy))/1000;
    %compressionAtMy=@(cAtMy) b*fck*cAtMy/epscmAtMy*integral(funCshahAtMy,0,epscmAtMy)/1000;
    sigmaSteelAtMy= fy;
    tensionAtMy= sigmaSteelAtMy*As/1000;
    cAtMy(k)=fsolve(@(cAtMy) compressionAtMy(cAtMy)-tensionAtMy, c0, fsolve_opts);
    c0=cAtMy(k);
    funM=@(epsc) (1-(1-epsc./eps0).^A).*(d-c0+(c0./epscmAtMy(c0)).*epsc) .* (epsc<=eps0) + exp(-kshah*(epsc-eps0).^1.15).*(d-c0+(c0./epscmAtMy(c0)).*epsc) .* (epsc>eps0);
    K = integral(funM,0,epscmAtMy(c0));
    M(k)=b*fck*cAtMy(k)/epscmAtMy(c0)*K/1000000;
          
end 
M
M = 1×4
1.0e+03 *

    0.6333    0.9588    0.6696    1.0822
toc
Elapsed time is 0.099555 seconds.

you need to provide input value to variable to the function epscmAtMy defined using function handle when calling it in the lines

funM=@(epsc) (1-(1-epsc./eps0).^A).*(d-c0+(c0./epscmAtMy(c0)).*epsc) .* (epsc<=eps0) + exp(-kshah*(epsc-eps0).^1.15).*(d-c0+(c0./epscmAtMy(c0)).*epsc) .* (epsc>eps0);
M(k)=b*fck*cAtMy(k)/epscmAtMy(c0)*K/1000000;

请先登录，再进行评论。

Finding a value and marking it in a graph

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

采纳的回答

11 个评论
显示 9更早的评论隐藏 9更早的评论

更多回答（0 个）

另请参阅

类别

标签

Community Treasure Hunt

Finding a value and marking it in a graph

0 个评论 显示 -2更早的评论隐藏 -2更早的评论

采纳的回答

11 个评论 显示 9更早的评论隐藏 9更早的评论

更多回答（0 个）

另请参阅

类别

标签

Community Treasure Hunt

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

11 个评论
显示 9更早的评论隐藏 9更早的评论