🍍C++(Python)肥皂泡沫普拉托边界膜曲面模型算法

Python | C++ | MATLAB | 数学 | 算法 | 物理 | 流体力学 | 泡沫 | 皂膜 | 曲面模型 | 厚度

🎯要点

🎯肥皂泡二维流体模拟 | 🎯泡沫普拉托边界膜曲面模型算法演化厚度变化 | 🎯螺旋曲面三周期最小结构生成

📜皂膜用例：Python计算物理粒子及拉格朗日和哈密顿动力学

🍇Python(MATLAB)三维曲面图

普拉托定律描述了肥皂膜的结构。这些定律是由比利时物理学家普拉托在 19 世纪根据他的实验观察制定的。自然界中的许多图案都是基于遵守这些定律的泡沫。此定律描述了皂膜的形状和构造如下：

肥皂膜由整个（完整的）光滑表面制成
肥皂膜的一部分的平均曲率在同一块肥皂膜上的任何点处处处恒定
肥皂膜总是沿着称为普拉托边界的边缘成三块相交，并且它们以 arccos(− 1 / 2 ) = 120°‘
这些普拉托边界在一个顶点处以四面体角 arccos(− 1 / 3 ) ≈ 109.47°。

除普拉托定律之外的结构是不稳定的，并且薄膜将很快倾向于重新排列自身以符合这些定律。美国数学家吉恩·艾伦·泰勒使用几何测量理论在数学上证明了这些定律适用于最小曲面。

肥皂膜是被空气包围的薄层液体（通常为水基）。例如，如果两个肥皂泡接触，它们会合并并在其间形成一层薄膜。因此，泡沫由通过普拉托边界连接的薄膜网络组成。肥皂膜可用作极小曲面的模型系统，极小曲面在数学中被广泛使用。

从数学角度来看，肥皂膜是最小表面。表面张力是单位面积产生表面所需的能量。薄膜——与任何物体或结构一样——倾向于以最小势能状态存在。为了最小化其能量，自由空间中的液滴自然呈现球形，对于给定的体积，其表面积最小。水坑和薄膜可以在其他力的存在下存在，例如重力和对基质原子的分子间吸引力。后一种现象称为润湿：基质原子和薄膜原子之间的结合力会导致总能量降低。在这种情况下，物体的最低能量配置是尽可能多的薄膜原子尽可能靠近基质。这将导致无限薄的薄膜，无限广泛地分布在基质上。实际上，粘附润湿效应（导致表面最大化）和表面张力效应（导致表面最小化）会相互平衡：稳定的结构可以是液滴、水坑或薄膜，具体取决于作用于身体的力。

💦Python绘制三维曲面：

Matplotlib 的 mpl_toolkits.mplot3d 工具包中的 axis3d 提供了用于创建三维曲面图的必要函数。曲面图是通过使用 ax.plot_surface() 函数创建的。

语法：

 ax.plot_surface(X, Y, Z)

其中 X 和 Y 是 x 和 y 点的二维数组，而 Z 是高度的二维数组。

示例：

 # Import libraries
 from mpl_toolkits import mplot3d
 import numpy as np
 import matplotlib.pyplot as plt
 
 x = np.outer(np.linspace(-3, 3, 32), np.ones(32))
 y = x.copy().T # transpose
 z = (np.sin(x **2) + np.cos(y **2) )
 
 fig = plt.figure(figsize =(14, 9))
 ax = plt.axes(projection ='3d')
 
 ax.plot_surface(x, y, z)
 plt.show()

梯度曲面图是三维曲面图与二维轮廓图的组合。在此图中，三维曲面的颜色与二维轮廓图相同。曲面较高的部分与曲面较低的部分颜色不同。

语法：

 surf = ax.plot_surface(X, Y, Z, cmap=, linewidth=0, antialiased=False)

属性 cmap= 设置表面的颜色。还可以通过调用Fig.colorbar来添加颜色条。下面的代码创建一个梯度曲面图：

 from mpl_toolkits import mplot3d
 import numpy as np
 import matplotlib.pyplot as plt
 
 x = np.outer(np.linspace(-3, 3, 32), np.ones(32))
 y = x.copy().T # transpose
 z = (np.sin(x **2) + np.cos(y **2) )
 
 fig = plt.figure(figsize =(14, 9))
 ax = plt.axes(projection ='3d')
 my_cmap = plt.get_cmap('hot')
 
 surf = ax.plot_surface(x, y, z,
                     cmap = my_cmap,
                     edgecolor ='none')
 
 fig.colorbar(surf, ax = ax,
             shrink = 0.5, aspect = 5)
 
 ax.set_title('Surface plot')
 plt.show()

使用 Matplotlib 绘制的三维曲面图可以投影到二维曲面上。下面的代码创建三维图并可视化其在二维等高线图上的投影：

 from mpl_toolkits import mplot3d
 import numpy as np
 import matplotlib.pyplot as plt
 
 x = np.outer(np.linspace(-3, 3, 32), np.ones(32))
 y = x.copy().T # transpose
 z = (np.sin(x **2) + np.cos(y **2) )
 
 fig = plt.figure(figsize =(14, 9))
 ax = plt.axes(projection ='3d')
 
 my_cmap = plt.get_cmap('hot')
 
 surf = ax.plot_surface(x, y, z, 
                     rstride = 8,
                     cstride = 8,
                     alpha = 0.8,
                     cmap = my_cmap)
 cset = ax.contourf(x, y, z,
                 zdir ='z',
                 offset = np.min(z),
                 cmap = my_cmap)
 cset = ax.contourf(x, y, z,
                 zdir ='x',
                 offset =-5,
                 cmap = my_cmap)
 cset = ax.contourf(x, y, z, 
                 zdir ='y',
                 offset = 5,
                 cmap = my_cmap)
 fig.colorbar(surf, ax = ax, 
             shrink = 0.5,
             aspect = 5)
 
 
 ax.set_xlabel('X-axis')
 ax.set_xlim(-5, 5)
 ax.set_ylabel('Y-axis')
 ax.set_ylim(-5, 5)
 ax.set_zlabel('Z-axis')
 ax.set_zlim(np.min(z), np.max(z))
 ax.set_title('3D surface having 2D contour plot projections')
 
 plt.show()

💦MATLAB最小曲面

构造最小曲面

 n=35;
 [X,Y,Z]=meshgrid(linspace(-pi,pi,n));
 [F,V] = isosurface(X,Y,Z,S,0);

曲面可视化

 cFigure;
 hold on;
 title('Schwarz P-surface','FontSize',fontSize);
 gpatch(F,V,'kw','k',1);
 axisGeom;
 camlight headlight;
 drawnow;

可视化曲面所有变体：

 cFigure;
 
 subplot(2,3,1);
 title('Schwarz P-surface','FontSize',fontSize);
 hold on;
 gpatch(F,V,pColors(1,:),'none',1);
 axisGeom;
 camlight headlight; lighting gouraud;
 view(-50,30);
 
 [X,Y,Z]=meshgrid(linspace(-pi,pi,n));
 [F,V] = isosurface(X,Y,Z,S,0.1);
 
 subplot(2,3,2);
 title('d','FontSize',fontSize);
 hold on;
 gpatch(F,V,pColors(2,:),'none',1);
 axisGeom;
 camlight headlight; lighting gouraud;
 view(-50,30);
 
 [X,Y,Z]=meshgrid(linspace(-2*pi,2*pi,n));
 [F,V] = isosurface(X,Y,Z,S,0.6);
 
 subplot(2,3,3);
 title('Gyroid','FontSize',fontSize);
 hold on;
 gpatch(F,V,pColors(3,:),'none',1);
 axisGeom;
 camlight headlight; lighting gouraud;
 view(-50,30);
 
 [X,Y,Z]=meshgrid(linspace(-pi,pi,n));
 [F,V] = isosurface(X,Y,Z,S,0);
 
 subplot(2,3,4);
 title('Neovius','FontSize',fontSize);
 hold on;
 gpatch(F,V,pColors(4,:),'none',1);
 axisGeom;
 camlight headlight; lighting gouraud;
 view(-50,30);
 
 [X,Y,Z]=meshgrid(linspace(-pi,pi,n));
 [F,V] = isosurface(X,Y,Z,S,0);
 
 subplot(2,3,5);
 title('w','FontSize',fontSize);
 hold on;
 gpatch(F,V,pColors(5,:),'none',1);
 axisGeom;
 camlight headlight; lighting gouraud;
 view(-50,30);
 
 [X,Y,Z]=meshgrid(linspace(-pi,pi,n));
 [F,V] = isosurface(X,Y,Z,S,0.5);
 
 subplot(2,3,6);
 title('pw','FontSize',fontSize);
 hold on;
 gpatch(F,V,pColors(6,:),'none',1);
 axisGeom;
 camlight headlight; lighting gouraud;
 view(-50,30);
 drawnow;

PreviousPython自行车六自由度飞行器多连接件非线性运动方程模型 NextPython和MATLAB粘性力接触力动态模型半隐式欧拉算法

Last updated 4 months ago

# Import libraries from mpl_toolkits import mplot3d import numpy as np import matplotlib.pyplot as plt x = np.outer(np.linspace(-3, 3, 32), np.ones(32)) y = x.copy().T # transpose z = (np.sin(x **2) + np.cos(y **2) ) fig = plt.figure(figsize =(14, 9)) ax = plt.axes(projection ='3d') ax.plot_surface(x, y, z) plt.show()

from mpl_toolkits import mplot3d import numpy as np import matplotlib.pyplot as plt x = np.outer(np.linspace(-3, 3, 32), np.ones(32)) y = x.copy().T # transpose z = (np.sin(x **2) + np.cos(y **2) ) fig = plt.figure(figsize =(14, 9)) ax = plt.axes(projection ='3d') my_cmap = plt.get_cmap('hot') surf = ax.plot_surface(x, y, z, rstride = 8, cstride = 8, alpha = 0.8, cmap = my_cmap) cset = ax.contourf(x, y, z, zdir ='z', offset = np.min(z), cmap = my_cmap) cset = ax.contourf(x, y, z, zdir ='x', offset =-5, cmap = my_cmap) cset = ax.contourf(x, y, z, zdir ='y', offset = 5, cmap = my_cmap) fig.colorbar(surf, ax = ax, shrink = 0.5, aspect = 5) ax.set_xlabel('X-axis') ax.set_xlim(-5, 5) ax.set_ylabel('Y-axis') ax.set_ylim(-5, 5) ax.set_zlabel('Z-axis') ax.set_zlim(np.min(z), np.max(z)) ax.set_title('3D surface having 2D contour plot projections') plt.show()

cFigure; subplot(2,3,1); title('Schwarz P-surface','FontSize',fontSize); hold on; gpatch(F,V,pColors(1,:),'none',1); axisGeom; camlight headlight; lighting gouraud; view(-50,30); [X,Y,Z]=meshgrid(linspace(-pi,pi,n)); [F,V] = isosurface(X,Y,Z,S,0.1); subplot(2,3,2); title('d','FontSize',fontSize); hold on; gpatch(F,V,pColors(2,:),'none',1); axisGeom; camlight headlight; lighting gouraud; view(-50,30); [X,Y,Z]=meshgrid(linspace(-2*pi,2*pi,n)); [F,V] = isosurface(X,Y,Z,S,0.6); subplot(2,3,3); title('Gyroid','FontSize',fontSize); hold on; gpatch(F,V,pColors(3,:),'none',1); axisGeom; camlight headlight; lighting gouraud; view(-50,30); [X,Y,Z]=meshgrid(linspace(-pi,pi,n)); [F,V] = isosurface(X,Y,Z,S,0); subplot(2,3,4); title('Neovius','FontSize',fontSize); hold on; gpatch(F,V,pColors(4,:),'none',1); axisGeom; camlight headlight; lighting gouraud; view(-50,30); [X,Y,Z]=meshgrid(linspace(-pi,pi,n)); [F,V] = isosurface(X,Y,Z,S,0); subplot(2,3,5); title('w','FontSize',fontSize); hold on; gpatch(F,V,pColors(5,:),'none',1); axisGeom; camlight headlight; lighting gouraud; view(-50,30); [X,Y,Z]=meshgrid(linspace(-pi,pi,n)); [F,V] = isosurface(X,Y,Z,S,0.5); subplot(2,3,6); title('pw','FontSize',fontSize); hold on; gpatch(F,V,pColors(6,:),'none',1); axisGeom; camlight headlight; lighting gouraud; view(-50,30); drawnow;