matlab等值面

clear; clc;
x=-100:0.5:100;
y=0:0.5:100;
z=0:0.5:100;
[X,Y,Z]=meshgrid(x,y,z);
r1=sqrt(X.^2+Y.^2);
r2=sqrt((X-Z).^2+Y.^2);
r=r1+r2+Z;
k1x=Z./r2-X./r1-X./r2;
k1y=-Y./r1-Y./r2;
k2x=X./r2-Z./r2-1;
k2y=Y./r2;
k1=sqrt(k1x.^2+k1y.^2);
k2=sqrt(k2x.^2+k2y.^2);
eta=(sqrt(k1)+sqrt(k2))/sqrt(2);
figure; hold on;
isosurface(X,Y,Z,eta,1.5);
isosurface(X,Y,Z,r,200);