Remove part of (fsurf) surface

Hi All, I'm plotting a symbolic function using fsurf, and it yields the result I want. However, roughly half of the graph is irrelevant and I would like to not display that half of the results:

I started to use fill3 in a crude attempt (below) to hide the information I'm not interested in, but then it also obscures the data I am interested in and doesn't look very clean.

Does anyone have any clever solutions to this problem? The code used is below:

clc, clear all, close all 


% DEFINE GEOMETRY CONSTs.
r = 20E-3; % absorber radius, [m]
Ah = pi*r^2; % heating surface area, [m^2]
Ac = Ah; % cooling surface area, [m^2]
l_stroke = 20E-3; % piston stroke, [m]
V1 = 250/1E6; % volume @ position1, [m^3]
% ^^ let's say V1 = 250 mL or one Red Bull can
V2 = Ah*l_stroke+V1; % volume @ position 2 (end of stroke)

V3 = V2; % for now, assuming very little expansion here
V4 = V1; 


% AIR CONSTs.
R = 8.314; % IG const., [J/mol-K]
T_amb = 20+273.15; % ambient temp., [K]
P_atm = 101325; % atm. pressure, [Pa]
n = (P_atm*V1)/(R*T_amb); % moles of air, arbitrarily chosen


% THERMAL BOUNDARY CONDITIONS
q = 1E3; % input energy, [W/m^2]
h = 10; % convection coeff., [W/m^2-K]
k = 0.8; % condunction coeff. glass cylinder, [W/m-K]
t = 3E-3; % piston wall thickness, [m]


% DEFINING WORK FUNCTIONS
syms T12 T34

W12(T12) = n*R*T12*log(V2/V1); % Work, [J]
dW12(T12) = q*Ah - (T12-T_amb)*(t/(k*Ac) + 1/(h*Ac))^-1; % Work rate, [J/s]
t12 = W12/dW12; % Time, [s]

W34(T34) = n*R*T34*log(V4/V3); % Work, [J]
dW34(T34) = -(T34-T_amb)*(t/(k*Ac) + 1/(h*Ac))^-1; % Work rate, [J/s]
t34 = W34/dW34; % Time, [s]

W_tot(T12,T34) = W12(T12) + W34(T34); 
t_stroke(T12,T34) = t12(T12) + t34(T34);

pwr = W_tot/t_stroke; % our final expression of power, [J/s]


%% 3D PLOT: Power Surface Plot
figure 
xylims = [290 1E3 290 1E3];
ph = fsurf(pwr,xylims,'ShowContours','on');
ylabel('T_1_-_2 (K)'); xlabel('T_3_-_4 (K)'); zlabel('Power (J s^-^1)')
zlim([0 q*Ah].*3)

% Image Quality Stuff
brighten(0.3)
ph.EdgeColor = [1 1 1];
ph.AmbientStrength = 0.8; 
ax = gca;
cb_ctrl = colorbar;
cb_ctrl.Limits = [0 q*Ah];
caxis([0 q*Ah]);
ax.Box = 'on';

Many thanks in advance,

Luke G.

[darova]

Is it possible to see this?

[Luke G.]

@darova, definitely! See below:

clc, clear all, close all 

% DEFINE GEOMETRY CONSTs.
r = 20E-3; % absorber radius, [m]
Ah = pi*r^2; % heating surface area, [m^2]
Ac = Ah; % cooling surface area, [m^2]
l_stroke = 20E-3; % piston stroke, [m]
V1 = 250/1E6; % volume @ position1, [m^3]
% ^^ let's say V1 = 250 mL or one Red Bull can
V2 = Ah*l_stroke+V1; % volume @ position 2 (end of stroke)
V3 = V2; % for now, assuming very little expansion here
V4 = V1; 


% AIR CONSTs.
R = 8.314; % IG const., [J/mol-K]
T_amb = 20+273.15; % ambient temp., [K]
P_atm = 101325; % atm. pressure, [Pa]
n = (P_atm*V1)/(R*T_amb); % moles of air, arbitrarily chosen


% THERMAL BOUNDARY CONDITIONS
q = 1E3; % input energy, [W/m^2]
h = 10; % convection coeff., [W/m^2-K]
k = 0.8; % condunction coeff. glass cylinder, [W/m-K]
t = 3E-3; % piston wall thickness, [m]


% DEFINING WORK FUNCTIONS
syms T12 T34

W12(T12) = n*R*T12*log(V2/V1); % Work, [J]
dW12(T12) = q*Ah - (T12-T_amb)*(t/(k*Ac) + 1/(h*Ac))^-1; % Work rate, [J/s]
t12 = W12/dW12; % Time, [s]

W34(T34) = n*R*T34*log(V4/V3); % Work, [J]
dW34(T34) = -(T34-T_amb)*(t/(k*Ac) + 1/(h*Ac))^-1; % Work rate, [J/s]
t34 = W34/dW34; % Time, [s]

W_tot(T12,T34) = W12(T12) + W34(T34); 
t_stroke(T12,T34) = t12(T12) + t34(T34);

pwr = W_tot/t_stroke; % our final expression of power, [J/s]

Answer 1

[Ameer Hamza]

The data for the surfaces plotted using fsurf is Read-only, so it is not possible to modified. However, you can use surf() like this and mask the unwanted region using nan

%% 3D PLOT: Power Surface Plot
figure 
xylims = [290 1E3 290 1E3];
[T34g, T12g] = meshgrid(linspace(xylims(1), xylims(2), 100), linspace(xylims(3), xylims(4), 100));
PHg = double(subs(pwr, [T12, T34], {T12g, T34g}));
mask = T34g < T12g;
PHg(mask) = nan;
ph = surf(T12g, T34g, PHg);
ylabel('T_1_-_2 (K)'); xlabel('T_3_-_4 (K)'); zlabel('Power (J s^-^1)')
zlim([0 q*Ah].*3)

% Image Quality Stuff
brighten(0.3)
ph.EdgeColor = [1 1 1];
ph.AmbientStrength = 0.8; 
ax = gca;
cb_ctrl = colorbar;
cb_ctrl.Limits = [0 q*Ah];
caxis([0 q*Ah]);
ax.Box = 'on';

[Luke G.]

This is brilliant! And very concise. Thanks for answering!

Answer 2

[darova]

Try numerical solution

F = matlabFunction(pwr);

[X,Y] = meshgrid(300:20:1000);
Z = F(X,Y);
Z1 = (Y-3/4*X)/10;          % line equations y = 3/4*x
Z(abs(Z) > 100) = nan;      % remove high values
Z(Z1<0) = nan;              % remove values above the line

contour(X,Y,Z1,[0 0],'b','linew',3)
surface(X,Y,Z,'facecolor','none')

view(3)
axis vis3d

[darova]

Or bettet using initmesh

F = matlabFunction(pwr);
%%
x = [300 400 1000 1000 300];        % contour
y = [300 300  750 1000 1000];
gd = [2; length(x); x(:); y(:)];    % geometry
dl = decsg(gd);                     % decomposition
[p,e,t] = initmesh(dl,'hmax',20);   % triangulate inside contour
Z = F(p(1,:),p(2,:));               % calculate Z value
ff.vertices = [p' Z(:)];
ff.faces    = t(1:3,:)';
ff.facecolor= 'y';

plot(x,y,'linew',3)
patch(ff)

view(3)
axis vis3d

[Luke G.]

@darova, thanks for these suggestions!