MATLAB隐函数三维绘图：可视化技巧

前言：

        小白路过，只是想在学习高等数学三重积分时让图像更直观一些，在网上粗略学习了一下MATLAB，本图来自张宇高数18讲例18.11

代码：

        首先创建一个函数文件isocurve3，用以绘制交线：

1.     function [nX,nY,nZ]=isocurve3(X,Y,Z,f1,f2)
2.     % 获取f1隐函数的三角面和三角顶点
3.     V1=f1(X,Y,Z);
4.     hel=isosurface(X,Y,Z,V1,0);
5.     % 将f1获取的三角顶点带入f2求得数值
6.     V2=f2(hel.vertices(:,1),hel.vertices(:,2),hel.vertices(:,3));
7.     % 检查三个顶点的数值是否同时含有有大于0数值及小于0数值
8.     mask=V2>0;
9.     outcount=sum(mask(hel.faces),2);
10.     cross=(outcount==2)|(outcount==1);
11.     crossing_tris=hel.faces(cross,:);
12.     % 通过旋转交换三个顶点次序，将小于0的点放在第一列
13.     out_vert=mask(crossing_tris);
14.     flip=sum(out_vert,2)==1;
15.     out_vert(flip,:)=1-out_vert(flip,:);
16.     ntri=size(out_vert,1);
17.     overt=zeros(ntri,3);
18.     for i=1:ntri
19.         v1i=find(~out_vert(i,:));
20.         v2i=1+mod(v1i,3);
21.         v3i=1+mod(v1i+1,3);
22.         overt(i,:)=crossing_tris(i,[v1i v2i v3i]);
23.     end
24.     % 类似于求重心
25.     u=(-V2(overt(:,1)))./(V2(overt(:,2))-V2(overt(:,1)));
26.     v=(-V2(overt(:,1)))./(V2(overt(:,3))-V2(overt(:,1)));
27.     uverts=repmat((1-u),[1 3]).*hel.vertices(overt(:,1),:)+repmat(u,[1 3]).*hel.vertices(overt(:,2),:);
28.     vverts=repmat((1-v),[1 3]).*hel.vertices(overt(:,1),:)+repmat(v,[1 3]).*hel.vertices(overt(:,3),:);
29.     % 因为可能含有多条曲线，因此逐段连线
30.     nX=nan(3,ntri);
31.     nX(1,:)=uverts(:,1)';
32.     nX(2,:)=vverts(:,1)';
33.     nY=nan(3,ntri);
34.     nY(1,:)=uverts(:,2)'; 
35.     nY(2,:)=vverts(:,2)';
36.     nZ=nan(3,ntri);
37.     nZ(1,:)=uverts(:,3)';
38.     nZ(2,:)=vverts(:,3)';
39.     end

        然后，创建主文件，绘制曲面及交线：

1.     clear all
2.     clc
3.     
4.     u = 0;                                          %控制坐标区间
5.     lmax = 2; hmax = 4;
6.     if u == 0
7.         lx = lmax; tx = -lx;                        %tx:x负半轴长；lx：x正半轴长
8.         ly = lx; ty = -ly; % 全轴
9.         lz = hmax; tz = 0;
10.     elseif u == 1 
11.         lx = lmax; tx = 0;
12.         ly = lx; ty = 0;  % 正半轴
13.         lz = hmax; tz = 0;   
14.     else 
15.         lx = 0; tx = -lmax;
16.         ly = 0; ty = tx;  % 负半轴
17.         lz = 0; tz = -hmax;  
18.     end
19.       
20.     for i = 1:4
21.         subplot(2,2,i)                                           %绘制四个视图
22.     
23.         f1 = @(x,y,z) 2*x - z;                                   %绘制隐函数曲面
24.         f2 = @(x,y,z) x.^2 + y.^2 - z;
25.         fimplicit3(f1,[tx lx ty ly tz lz],'EdgeColor','none',...
26.             'FaceColor','cyan')
27.         hold on
28.         fimplicit3(f2,[tx lx ty ly tz lz],'EdgeColor','none',...
29.             'FaceColor','y')
30.         alpha(0.3)                                                %曲面透明度
31.         axis equal
32.     
33.         [X,Y,Z]=meshgrid(-5:.1:5);                                 %绘制交线
34.         [nX,nY,nZ]=isocurve3(X,Y,Z,f1,f2);
35.         line(nX(:),nY(:),nZ(:),'LineWidth',1.5,'Color',[0.5,0,1]); 
36.     
37.         
38.     hold on                                                         %绘制XYZ轴
39.         plot3([tx,0,lx],[0,0,0],[0,0,0],'r',[0,0,0],[ty,0,ly],...
40.             [0,0,0],'g',[0,0,0],[0,0,0],[tz,0,lz],'b')
41.         xlabel('x轴');ylabel('y轴');zlabel('z轴');
42.     
43.         if i == 1
44.             title('三维图')
45.         elseif i == 2
46.             view([1,0,0])
47.             title('yOz平面')
48.             
49.         elseif i == 3
50.             view([0,-1,0])
51.             title('xOz平面')
52.         else
53.             view([0,0,1])
54.             title('xOy平面')
55.         end
56.         axis equal
57.         hold off
58.     end

展示效果：

PS：交线部分的函数来自这位大佬https://blog.csdn.net/slandarer/article/details/125710021

我也是小白，很多错误不是很了解