Remove part of (fsurf) surface

23 ビュー (過去 30 日間)
Luke G.
Luke G. 2020 年 5 月 20 日
編集済み: Ameer Hamza 2020 年 5 月 20 日
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.
  2 件のコメント
darova
darova 2020 年 5 月 20 日
Is it possible to see this?
Luke G.
Luke G. 2020 年 5 月 20 日
@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]

サインインしてコメントする。

採用された回答

Ameer Hamza
Ameer Hamza 2020 年 5 月 20 日
編集済み: Ameer Hamza 2020 年 5 月 20 日
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';
  1 件のコメント
Luke G.
Luke G. 2020 年 5 月 20 日
This is brilliant! And very concise. Thanks for answering!

サインインしてコメントする。

その他の回答 (1 件)

darova
darova 2020 年 5 月 20 日
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
  2 件のコメント
darova
darova 2020 年 5 月 20 日
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
75
Luke G.
Luke G. 2020 年 5 月 20 日
@darova, thanks for these suggestions!

サインインしてコメントする。

カテゴリ

Help Center および File ExchangeSurface and Mesh Plots についてさらに検索

タグ

製品


リリース

R2020a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by