我有一个参数曲面,其定义如下:朗伯W函数我想要用 Asymptote 来显示。为此,Lambert W 函数是使用牛顿法实现的,并且必须将原始(封闭)表面分成两个开放表面以避免除以零的问题。以下是 MWE:
settings.render=8;
settings.prc=false;
settings.outformat="pdf";
import graph3;
size(200);
currentprojection=orthographic(40,10,10);
// contour value
real c = 0.006;
// parameter ranges
real umin = asin(1.5*exp(1)*sqrt(c*(sqrt(2*pi))));
real umax = pi-asin(1.5*exp(1)*sqrt(c*(sqrt(2*pi))));
real vmin = 0;
real vmax = 1;
// Lambert W function
real w1(real w, real z, int i){return z<-1/exp(1) - 0.00001 ? 1/0 : z<-1/exp(1) ? -1 : i>0 && abs((w*exp(w)-z)/(exp(w)+w*exp(w))) > 1e-7 ? w1(w-(w*exp(w)-z)/(exp(w)+w*exp(w)),z,i-1) : w-(w*exp(w)-z)/(exp(w)+w*exp(w));};
// auxiliary functions
real y5(real h, real p){return (1/4.)*sqrt(15./pi) * sin(2*p) * sin(h)**2;};
real r1(real y){return -6*w1(-2,-sqrt(c*9*sqrt(30)/abs(y))/4,200);};
real r2(real y){return -6*w1(1,-sqrt(c*9*sqrt(30)/abs(y))/4,200);};
// x, y, and z coordinates of the surfaces
real x11(real u, real v){return r1(y5(u,v*pi/2+(0.5-v)*asin(exp(1)**2*9.*c*sqrt(2.*pi)/(4.*sin(u)**2))))*sin(u)*cos(v*pi/2+(0.5-v)*asin(exp(1)**2*9.*c*sqrt(2.*pi)/(4.*sin(u)**2)));};
real y11(real u, real v){return r1(y5(u,v*pi/2+(0.5-v)*asin(exp(1)**2*9.*c*sqrt(2.*pi)/(4.*sin(u)**2))))*sin(u)*sin(v*pi/2+(0.5-v)*asin(exp(1)**2*9.*c*sqrt(2.*pi)/(4.*sin(u)**2)));};
real z11(real u, real v){return r1(y5(u,v*pi/2+(0.5-v)*asin(exp(1)**2*9.*c*sqrt(2.*pi)/(4.*sin(u)**2)) ))*cos(u);};
real x12(real u, real v){return r2(y5(u,v*pi/2+(0.5-v)*asin(exp(1)**2*9.*c*sqrt(2.*pi)/(4.*sin(u)**2))))*sin(u)*cos(v*pi/2+(0.5-v)*asin(exp(1)**2*9.*c*sqrt(2.*pi)/(4.*sin(u)**2)));};
real y12(real u, real v){return r2(y5(u,v*pi/2+(0.5-v)*asin(exp(1)**2*9.*c*sqrt(2.*pi)/(4.*sin(u)**2))))*sin(u)*sin(v*pi/2+(0.5-v)*asin(exp(1)**2*9.*c*sqrt(2.*pi)/(4.*sin(u)**2)));};
real z12(real u, real v){return r2(y5(u,v*pi/2+(0.5-v)*asin(exp(1)**2*9.*c*sqrt(2.*pi)/(4.*sin(u)**2)) ))*cos(u);};
triple f11(pair p){return (x11(p.x,p.y),y11(p.x,p.y),z11(p.x,p.y));};
triple f12(pair p){return (x12(p.x,p.y),y12(p.x,p.y),z12(p.x,p.y));};
surface s11 = surface(f=f11,a=(umin, vmin),b=(umax,vmax)); // this works
surface s12 = surface(f=f12,a=(umin, vmin),b=(umax,vmax)); // this works
// surface s11 = surface(f=f11,a=(umin, vmin),b=(umax,vmax),Spline); // this doesn't work
// surface s12 = surface(f=f12,a=(umin, vmin),b=(umax,vmax),Spline); // this doesn't work
draw(s11, red+opacity(0.5));
draw(s12, red+opacity(0.5));
如果我添加Spline
指令使表面看起来平滑,Asymptote 会因错误而崩溃some path/graph_splinetype.asy: 89.10: function values are not periodic
。我试图了解出了什么问题graph_splinetype.asy
,graph3.asy
但不幸的是,我不够熟练,无法成功。所以我的问题是:是否有机会使用Spline
这个参数曲面或者用另一种方法使它看起来光滑?
更令人费解的是,Spline
对于类似的参数曲面来说,它可以很好地工作(尽管这个曲面绕 z 轴具有旋转对称性,这可能很重要),即这个:
settings.render=8;
settings.prc=false;
settings.outformat="pdf";
import graph3;
size(200);
currentprojection=orthographic(40,10,10);
// contour value
real c = 0.006;
// parameter ranges
real umin = 0;
real umax = acos(2*exp(1)*c*sqrt(2*pi))-0.0000001;
real vmin = 0;
real vmax = 2*pi;
// Lambert W function
real w1(real w, real z, int i){return z<-1/exp(1) - 0.00001 ? 1/0 : z<-1/exp(1) ? -1 : i>0 && abs((w*exp(w)-z)/(exp(w)+w*exp(w))) > 1e-7 ? w1(w-(w*exp(w)-z)/(exp(w)+w*exp(w)),z,i-1) : w-(w*exp(w)-z)/(exp(w)+w*exp(w));};
// auxiliary functions
real y1(real h){return sqrt(3./pi)/2.*cos(h);};
real r1(real y){return -2*w1(-2,-c*sqrt(6)/abs(y),200);};
real r2(real y){return -2*w1(1,-c*sqrt(6)/abs(y),200);};
// x, y, and z coordinates of the surfaces
real x11(real u, real v){return r1(y1(u))*sin(u)*cos(v);};
real y11(real u, real v){return r1(y1(u))*sin(u)*sin(v);};
real z11(real u, real v){return r1(y1(u))*cos(u);};
real x12(real u, real v){return r2(y1(u))*sin(u)*cos(v);};
real y12(real u, real v){return r2(y1(u))*sin(u)*sin(v);};
real z12(real u, real v){return r2(y1(u))*cos(u);};
triple f11(pair p){return (x11(p.x,p.y),y11(p.x,p.y),z11(p.x,p.y));};
triple f12(pair p){return (x12(p.x,p.y),y12(p.x,p.y),z12(p.x,p.y));};
surface s11 = surface(f11,(umin,vmin),(umax,vmax),50,Spline);
surface s12 = surface(f12,(umin,vmin),(umax,vmax),50,Spline);
draw(s11, red+opacity(0.5));
draw(s12, red+opacity(0.5));
答案1
您能详细说明一下您的系统和 Asymptote 版本吗?在我的计算机(Linux 64 位、Sid、Asymptote svn)上我没有遇到任何问题,您的代码产生了一个平滑的表面。
由于在表面例程中实现了周期性样条线或非结样条线之间的选择,因此它似乎比最终的参数样条线例程限制更少(我很久以前编写了此平滑表面例程的一部分)。可以使用以下代码强制选择“非结”
splinetype[] Notaknot={notaknot,notaknot,notaknot};
surface s11 = surface(f=f11,a=(umin,vmin),b=(umax,vmax),8,16,Notaknot,Notaknot);
surface s12 = surface(f=f12,a=(umin,vmin),b=(umax,vmax),8,16,Notaknot,Notaknot);
我希望它能在你的系统上运行。
奥格