“常规”数据的ListPlot3D包含在“不规则”结构中的晶格。地区

发布于 2024-09-27 02:21:25 字数 574 浏览 1 评论 0原文

当尝试在 Mathematica 中绘制一些 3D 数据时，我遇到以下问题。数据最初是在正则格子上计算的，即 I 计算。

data = Flatten[Table[{x,y,f[x,y]},{x,x0,x1,dx},{y,y0,y1,dy}],1]

问题是 f 仅在平面上的非凸子集 U 中取实值。 U 实际上相当令人讨厌：一个“三角形”区域，所有角都是尖角的（想象一下三个相同的圆之间的区域，其中任意两个圆彼此相切）。当我尝试使用 ListPlot3D 绘制 data 时，后者在 U 的凸包上绘制一个曲面。

我想知道是否有办法告诉 Mathematica 仅在 U 内限制自身。我在想，既然我的点已经位于“常规”格子上，这应该不会太困难，但我还没有找到解决方案。

我知道 RegionFunction 选项，但在我的情况下计算起来非常昂贵，因此我试图找出一种仅使用 data 中已计算的点的方法。

原文

I have the following problem when trying to plot some 3D data in Mathematica. The data are initially computed on a regular lattice, that is I compute

data = Flatten[Table[{x,y,f[x,y]},{x,x0,x1,dx},{y,y0,y1,dy}],1]

The problem is that f takes real values only in a non-convex subset U on the plane. U is actually quite nasty: a "triangular" region where all corners are cuspy (imagine the region in between three identical circles with any two of them tangent to each other). When I try to plot data with ListPlot3D, the latter plots a surface over the convex hull of U.

I was wondering if there is a way to tell Mathematica to restrict itself only inside U. I was thinking that since my points already lie on a "regular" lattice this should not be too difficult but I have not found yet a solution.

I am aware of the RegionFunction option but it is very expensive to compute in my case so I am trying to figure out a way that uses only the already computed points in data.

分享到QQ

分享到微博

如果你对这篇内容有疑问，欢迎到本站社区发帖提问参与讨论，获取更多帮助，或者扫码二维码加入 Web 技术交流群。

发布评论

需要登录才能够评论，你可以免费注册一个本站的账号。

尴尬癌患者 2024-10-04 02:21:25

要生成非常漂亮的图片，而不使用无数的点，您可能需要使用包含所需区域边界的非均匀分布的点。这是一个有点类似于您所描述的示例。我们从三个相互相切的圆开始。

circPic = Graphics[{Circle[{0, Sqrt[3]}, 1],
  Circle[{-1, 0}, 1], Circle[{1, 0}, 1]}]

Mathematicagraphics

我们编写一个布尔函数，用于确定矩形 {-1/2,1/2} 中的点是否为{0,Sqrt[3]/2} 位于所有圆之外，并使用它来生成感兴趣区域中的一些点。

inRegionQ[p:{x_, y_}] := Norm[p - {1, 0}] > 1 &&
  Norm[p + {1, 0}] > 1 && Norm[p - {0, Sqrt[3]}] > 1;
rectPoints = N[Flatten[Table[{x, y},
  {x, -1/2, 1/2, 0.02}, {y, 0.05, Sqrt[3]/2, 0.02}], 1]];
regionPoints = Select[rectPoints, inRegionQ];

现在我们生成边界。参数 n 决定了我们在边界上放置多少个点。

n = 120;
boundary = N[Join[
  Table[{1 - Cos[t], Sin[t]}, {t, Pi/n, Pi/3, Pi/n}],
  Table[{Cos[t], Sqrt[3] - Sin[t]}, {t, Pi/3 + Pi/n, 2 Pi/3, Pi/n}],
  Table[{Cos[t] - 1, Sin[t]}, {t, Pi/3 - Pi/n, 0, -Pi/n}]]];
points = Join[boundary, regionPoints];

我们来看一下。

Show[circPic, Graphics[Point[points]],
  PlotRange -> {{-3/4, 3/4}, {-0.3, 1.3}}]

Mathematicagraphics

现在，我们定义一个函数并使用 ListPlot3D 尝试绘制它。

f[x_, y_] := -(1 - Norm[{x - 1, y}]) (1 - Norm[{x + 1, y}])*
  (1 - Norm[{x, y - Sqrt[3]}]);
points3D = {#[[1]], #[[2]], f[#[[1]], #[[2]]]} & /@ points;
pic = ListPlot3D[points3D, Mesh -> All]

Mathematicagraphics

无论如何，我们必须删除该区域之外的内容。在这个特定的例子中，我们可以利用函数在边界上为零的事实。

DeleteCases[Normal[pic], Polygon[{
  {x1_, y1_, z1_?(Abs[#] < 1/10.0^6 &)},
  {x2_, y2_, z2_?(Abs[#] < 1/10.0^6 &)},
  {x3_, y3_, z3_?(Abs[#] < 1/10.0^6 &)}}, ___],
  Infinity]

Mathematicagraphics

相当不错，但在尖点附近有一些问题，而且它绝对不是很通用，因为它使用了特定的函数的属性。如果您检查 pic 的结构，您会发现它包含一个 GraphicsComplex，并且该 GraphicsComplex 第一个参数中的前 n 个点正是边界。证据如下：

Most /@ pic[[1, 1, 1 ;; n]] == boundary

现在边界由三个分量组成，我们想要删除由仅从这些分量之一中选择的点形成的任何三角形。以下代码执行此操作。请注意，boundaryComponents 包含形成边界的那些点的索引，并且如果 A 是任一 B 的子集，则 someSubsetQ[A,Bs] 返回 true。我们想要删除多多边形中作为边界组件之一的子集的那些三角形索引。这是通过以下代码中的 DeleteCases 命令完成的。

哦，让我们也添加一些装饰。

subsetQ[A_, B_] := Complement[A, B] == {};
someSubsetQ[A_, Bs_] := Or @@ Map[subsetQ[A, #] &, Bs];
boundaryComponents = Partition[Prepend[Range[n], n], 1 + n/3, n/3];
Show[pic /. Polygon[triangles_] :> {EdgeForm[Opacity[0.3]],
  Polygon[DeleteCases[triangles, _?(someSubsetQ[#, boundaryComponents] &)]]},
Graphics3D[{Thick, Line[Table[Append[pt, 0],
  {pt, Prepend[boundary, Last[boundary]]}]]}]]

To generate a really nice picture, without using zillions of points, you might want to use a non-uniform distribution of points that includes the boundary of the region you want. Here's an example somewhat like what you describe. We start with three mutually tangent circles.

circPic = Graphics[{Circle[{0, Sqrt[3]}, 1],
  Circle[{-1, 0}, 1], Circle[{1, 0}, 1]}]

Mathematica graphics

We write a Boolean function that determines whether a point in the rectangle {-1/2,1/2} by {0,Sqrt[3]/2} lies outside of all the circles and use this to generate some points in the region of interest.

inRegionQ[p:{x_, y_}] := Norm[p - {1, 0}] > 1 &&
  Norm[p + {1, 0}] > 1 && Norm[p - {0, Sqrt[3]}] > 1;
rectPoints = N[Flatten[Table[{x, y},
  {x, -1/2, 1/2, 0.02}, {y, 0.05, Sqrt[3]/2, 0.02}], 1]];
regionPoints = Select[rectPoints, inRegionQ];

Now we generate the boundary. The parameter n determines how many points we place on the boundary.

n = 120;
boundary = N[Join[
  Table[{1 - Cos[t], Sin[t]}, {t, Pi/n, Pi/3, Pi/n}],
  Table[{Cos[t], Sqrt[3] - Sin[t]}, {t, Pi/3 + Pi/n, 2 Pi/3, Pi/n}],
  Table[{Cos[t] - 1, Sin[t]}, {t, Pi/3 - Pi/n, 0, -Pi/n}]]];
points = Join[boundary, regionPoints];

Let's take a look.

Show[circPic, Graphics[Point[points]],
  PlotRange -> {{-3/4, 3/4}, {-0.3, 1.3}}]

Mathematica graphics

Now, we define a function and use ListPlot3D to try to plot it.

f[x_, y_] := -(1 - Norm[{x - 1, y}]) (1 - Norm[{x + 1, y}])*
  (1 - Norm[{x, y - Sqrt[3]}]);
points3D = {#[[1]], #[[2]], f[#[[1]], #[[2]]]} & /@ points;
pic = ListPlot3D[points3D, Mesh -> All]

Mathematica graphics

Somehow, we've got to delete that stuff that lies outside the region. In this particular example, we can use the fact that the function is zero on the boundary.

DeleteCases[Normal[pic], Polygon[{
  {x1_, y1_, z1_?(Abs[#] < 1/10.0^6 &)},
  {x2_, y2_, z2_?(Abs[#] < 1/10.0^6 &)},
  {x3_, y3_, z3_?(Abs[#] < 1/10.0^6 &)}}, ___],
  Infinity]

Mathematica graphics

Pretty good, but there are a couple of problems near the cusps and it's definitely not very general since it used a specific property of the function. If you examine the structure of pic, you'll find that it contains a GraphicsComplex and the first n points in the first argument of that GraphicsComplex is exactly the boundary. Here's proof:

Most /@ pic[[1, 1, 1 ;; n]] == boundary

Now the boundary comes in three components and we want to delete any triangle that is formed by points chosen from just one of those components. The following code does this. Note that boundaryComponents contains the indices of those points that form the boundary and someSubsetQ[A,Bs] returns true if A is a subset of any one of the Bs. We want to delete those triangle indices in the multi-Polygon that are subsets of one of the boundaryComponents. That's accomplished in the following code by the DeleteCases command.

Oh, and let's add some decoration, too.

subsetQ[A_, B_] := Complement[A, B] == {};
someSubsetQ[A_, Bs_] := Or @@ Map[subsetQ[A, #] &, Bs];
boundaryComponents = Partition[Prepend[Range[n], n], 1 + n/3, n/3];
Show[pic /. Polygon[triangles_] :> {EdgeForm[Opacity[0.3]],
  Polygon[DeleteCases[triangles, _?(someSubsetQ[#, boundaryComponents] &)]]},
Graphics3D[{Thick, Line[Table[Append[pt, 0],
  {pt, Prepend[boundary, Last[boundary]]}]]}]]

回复收藏 0 原文

伪心 2024-10-04 02:21:25

这可能不是最佳解决方案，所以我将保留这个问题，以防有人有更好的想法。

这是我解决问题中描述的问题的方法。首先，我替换了 data 中的点 {x,y,f[x,y]}，其中 f[x,y] 是复杂的通过{x,y,None}。然后以下函数根据我的数据点创建 3D 表面。请注意，data 是

data = Table[{x,y,f[x,y]},{x,x0,x1,dx},{y,y0,y1,dy}]

没有扁平化的结果（这对于以下函数效果更好）。功能是：

makeSurface[data_] := Module[{len1, len2, polys, a, b, c, d, p},
  len1 = Length[data];
  len2 = Length[data[[1]]];
  polys = Table[
    a = data[[i, j]];
    b = data[[i + 1, j]];
    c = data[[i + 1, j + 1]];
    d = data[[i, j + 1]];
    p = Select[{a, b, c, d}, #[[3]] =!= None &];
    If[Length[p] >= 2, Polygon[p], None],
    {i, 1, len1 - 1}, {j, 1, len2 - 1}];
  Graphics3D[Join[{EdgeForm[None],FaceForm[Directive[GrayLevel[0.5]]]},
    Select[Flatten[polys, 1], # =!= None &]]]]

上面的代码可能可以优化，但对我来说已经足够好了。

This is probably not the optimal solution, so I will keep the question open in case someone has a better idea.

Here is how I solved the problem I described in my question. First, I replaced points {x,y,f[x,y]} in data for which f[x,y] was complex by {x,y,None}. Then the following function creates a 3D surface out of my data points. Note here that data is the result of

data = Table[{x,y,f[x,y]},{x,x0,x1,dx},{y,y0,y1,dy}]

that is, no flattening (which works better for the following function). The function is:

makeSurface[data_] := Module[{len1, len2, polys, a, b, c, d, p},
  len1 = Length[data];
  len2 = Length[data[[1]]];
  polys = Table[
    a = data[[i, j]];
    b = data[[i + 1, j]];
    c = data[[i + 1, j + 1]];
    d = data[[i, j + 1]];
    p = Select[{a, b, c, d}, #[[3]] =!= None &];
    If[Length[p] >= 2, Polygon[p], None],
    {i, 1, len1 - 1}, {j, 1, len2 - 1}];
  Graphics3D[Join[{EdgeForm[None],FaceForm[Directive[GrayLevel[0.5]]]},
    Select[Flatten[polys, 1], # =!= None &]]]]

The above code can be probably optimized but it worked well enough for me.

回复收藏 0 原文

~没有更多了~