Can someone help me display the area bounded by the curve and the linear line. Pls see attached photo.

[John D'Errico]

I already showed you how to solve that problem in your last question. Thre is no need to keep on asking it over and over again.

[Mark Cuanan]

I am so sorry but I actually use different approach. The method is different from the one you posted that is why I am finding difficulty in finding the intersect areas. Thank you so much for your help.

Answer 1

[Star Strider]

One approach —

clear all;
close all

format shortG

filename = 'GT+10-ubc-X-Axis - Static.csv'; 


T1 = readtable(filename, 'VariableNamingRule','preserve')

T1 = 31x2 table
    displacement->    Base_Shear
    ______________    __________

            0              0.1  
         0.03           1256.2  
         0.05           2512.4  
         0.08           3768.5  
         0.11             4998  
         0.13           5772.7  
         0.16           6167.4  
         0.19           6478.3  
         0.21           6687.1  
         0.24           6783.4  
         0.27           6873.1  
         0.29           6962.3  
         0.32           6994.1  
         0.35             7014  
         0.37           7005.1  
          0.4           7084.8  

% T1 = readtable('GT+10-ubc-X-Axis - Static.csv', 'VariableNamingRule','preserve')

% clear all 
% filename = 'GT+10-ubc-X-Axis - Static.csv';
data1 = readmatrix(filename);
x=data1(:,1);
y=data1(:,2);
% start point
xs=0;
ys=0;
% end point
maxy=max(y)

maxy = 
       7238.3

ip= find(y ==maxy);
%ip= find(x =maxy=displacement);
dist=1

dist = 
     1

figure(4)
plot(x(1:ip+dist),y(1:ip+dist),'LineWidth',2, 'DisplayName','Non-Elastic Response')

xe=x(ip);
ye=y(ip);

yin = @(t) interp1(x(1:ip), y(1:ip), t);
ynn=yin

ynn = function_handle with value:
    @(t)interp1(x(1:ip),y(1:ip),t)

% Generate t values for plotting
t_values = linspace(x(1), x(ip), 200);
% Evaluate yin for each t value
y_values = yin(t_values);

% Plot the graph

%plot(t_values, y_values);

x=t_values

x = 1x200

0 0.0024121 0.0048241 0.0072362 0.0096482 0.01206 0.014472 0.016884 0.019296 0.021709 0.024121 0.026533 0.028945 0.031357 0.033769 0.036181 0.038593 0.041005 0.043417 0.045829 0.048241 0.050653 0.053065 0.055477 0.057889 0.060302 0.062714 0.065126 0.067538 0.06995

y=y_values

y = 1x200

1.0e+00 * 0.1 101.09 202.08 303.08 404.07 505.06 606.05 707.05 808.04 909.03 1010 1111 1212 1341.4 1492.9 1644.4 1795.9 1947.4 2098.9 2250.4 2401.9 2539.7 2640.7 2741.7 2842.7 2943.7 3044.7 3145.7 3246.7 3347.7

d=length(x)

d = 
   200

ym=0

ym = 
     0

i=1

i = 
     1

while ym<y(end)-400
%intersectingpoint
xi=x(i);
yi=y(i);
%midpoint 
xm=xi/(0.6);
ym=yi/(0.6);
% Calculate slope
% Define the equation of the line
syms xx % Define x as a symbolic variable

m1 = (ym - yi) / (xm - xi);
m2=(ye-ym)/(xe-xm);
eqn2 = ym + m2*(xx - xm);
dxx=0.01;
eqn = yi + m1*(xx - xi);

domx=linspace(xi,xm,100);%
dxx=domx(2)-domx(1);
domx2=xm:dxx:xe-0.1;
curvev = yin(domx);
linev = subs(eqn,xx,domx);
dx1=domx(2)-domx(1);
curvev2 = yin(domx2);
linev2 = subs(eqn,xx,domx2);
%cld=min(abs(curvev2-linev2));
%intpo= find(abs(curvev2-linev2) ==cld);
%x3=domx2(intpo);
%y3=curvev2(intpo);
domx3=linspace(xm,xe,100);%xm:dxx:x3
curvev2 = yin(domx3);
linev2 = subs(eqn2,xx,domx3);
dx2=domx3(2)-domx3(1);
area1(i)=sum(abs(curvev-linev))*dx1;
area2(i)=sum(curvev2-linev2)*dx2;

tot_diff(i)=abs((area1(i)-area2(i)));


i=i+1;
end

cld=min(tot_diff);
minindex= find( tot_diff==cld);
figure(4)
xi=x(minindex);
yi=y(minindex);
xm=xi/(0.6);
ym=yi/(0.6);
yn = ynn(xm);
yn2=(yn/xm)*xe;

%hold on
%plot([xm,xe],[ym,yn2],'red','DisplayName','Linear Elastic Response')

xy =  [x; y];           % Non-Linear Elastic Response Curve

% [b1, b2] = bounds(xy,2)
% s = [xs ys]
% m = [xm ym]
% e = [xe ye]

B_s_m = [xs 1; xm 1] \ [ys;ym]                      % Linear Regression

B_s_m = 2x1

1.0e+00 * 46957 0

s_m_x = x(x >= xs & x <= xm);
s_m_y = [s_m_x(:) ones(size(s_m_x(:)))] * B_s_m;

B_m_e = [xm 1; xe 1] \ [ym;ye]                      % Linear Regression

B_m_e = 2x1

1.0e+00 * 2388.9 6091.7

m_e_x = x(x >= xm & x <= xe);
m_e_y = [m_e_x(:) ones(size(m_e_x(:)))] * B_m_e;

xint = x;
yc = [s_m_y; m_e_y];
dy = yc - y(:);

Lv1 = x >= xm;
y2 = y(Lv1);
Lv1 = y2(:) >= m_e_y;                           % Region Of 'A_2'
A_2 = trapz(x(Lv1), dy(Lv1))                    % 'A_2'

A_2 = 
       7.1642

cumulative_area = cumtrapz(xint, abs(dy));
Area = cumulative_area(end);

A_1 = Area - A_2

A_1 = 
       47.147

hold on
plot([xs,xm,xe],[ys,ym,ye],'red','DisplayName','Billinear Curve Fit')
hold on
plot(xs, ys, '.', 'DisplayName','ys')%, 'Marker','p', 'MarkerSize',20)
plot(xm, ym, '.', 'DisplayName','ym')%, 'Marker','p', 'MarkerSize',20)
plot(xe, ye, '.', 'DisplayName','ye')%, 'Marker','p', 'MarkerSize',20)
% hold off
%Plotting
title('Billinear Curve Fitting of Pushover Curve')
xlabel('Displacement (m)')
ylabel('Baseshear Force (kN)')
ax=gca;
ax.YAxis.Exponent = 0;

hold on
plot([0,xm],[ym,ym],'--','HandleVisibility','off')
plot([xm,xm],[0,ym],'--','HandleVisibility','off')
%plot([xe,xe],[0,yn2],'--','HandleVisibility','off')
%plot([0,xe],[yn2,yn2],'--','HandleVisibility','off')
%plot([0,xe],[ye,ye],'--','HandleVisibility','off')

text(xm+0.01, ym-200, sprintf('\\leftarrow A_1 = %.3f', A_1))
text(0.3, 6900, sprintf('\\leftarrow A_2 = %.3f', A_2), 'Rotation',-10)

% hold on
legend('Location','best')

% lgd = legend;
% lgd.NumColumns = 1;

% Display the results
disp('Values:');

Values:

disp(['y coordinate of yield point umax = ', num2str(xe)]);

y coordinate of yield point umax = 0.48

disp(['x cordinate of yield point uy = ', num2str(xm)]);

x cordinate of yield point uy = 0.13668

disp(['y coordinate of yield point Vy = ', num2str(ym)]);

y coordinate of yield point Vy = 6418.1787

disp(['y coordinate of yield point Ve = ', num2str(yn2)]);

y coordinate of yield point Ve = 20581.2369

disp(['y coordinate of yield point Vmax = ', num2str(ye)]);

y coordinate of yield point Vmax = 7238.34

figure
yyaxis left
plot(x, y, 'DisplayName','Non-Elastic Response')
hold on
plot(s_m_x, s_m_y,'DisplayName','s - m')
plot(m_e_x, m_e_y,'DisplayName','m - e')
plot(x, dy, 'DisplayName','Curve Difference')
plot(x, abs(dy), '.', 'DisplayName','Curve Difference - Absolute Value')
xline(xm, '--k', 'xm', 'DisplayName','xm')
hold off
yyaxis right
plot(x, cumulative_area, 'DisplayName','Cumulative Area')
hold off
grid
legend('Location','southoutside', 'NumColumns',2)
title('Plots Demonstraating ''Cumulative Area'' Calculation')

To use trapz or cumtrapz you need to interpolate the lines defined by the three points and then subtract them from the curve. The linear regressions do that here. Since I am not certain how you want to calculate the ares, I provided both the difference and the absolute value of the difference between the olines and the curve, and did the integration on the absolute value. Change that if necessary.

It took a bit of time to figure out how your code works. My apologies for the delay.

EDIT — (10 Mar 2024 at 04:28)

Added specific and calculations.

EDIT — (10 Mar 2024 at 04:45)

Added text objects for and .

.

[Mark Cuanan]

Hello sir. Can you somehow display the final area A1 and A2 please? corrsponding to the final plot? I hope this could be done. Thank you

[Star Strider]

That looks like different data.

A couple text calls are all that would be necessary for that. I will add them to my original code. Check it in about five minutes.