使用FG方法将光盘映射到正方形

Question

上下文

我正在处理计算流体力学中的网格。我想在一个圆圈周围生成一个结构化的网格。我的计划是生成一个极坐标网格(下图中的左网格)，然后使用FG公式得到最终的网格(下图中的右网格)。我使用FG方法从第4页的这篇文章地图光盘到正方形。不幸的是，这篇文章没有提到如何处理公式中的奇点。以下是几种表达方式：

x = sgn(uv)/(v*sqrt 2)*sqrt(u**2+v**2-sqrt((u**2+v**2)(u**2+v**2-4u**2v**2)))

y = sgn(uv)/(u*sqrt 2)*sqrt(u**2+v**2-sqrt((u**2+v**2)(u**2+v**2-4u**2v**2)))

​

​

在编写程序之前，我正与这些格式的一些问题作斗争。

问题

• 为什么这些公式将下列点映射为：(1,0)、(0,1)、(-1,0)、(-1,-1)和(0,0)到(0,0)点？
• 我应该如何得到一个圆和一个正方形之间的中间形状，如下图所示。
• 是否有可能提供一个算法来从左边的网格中得到正确的网格？

​

​

以下是我的尝试：

"""Map a circular computational domain with structured mesh around a circle (circular cylinder in 3D) to
Rectangular domain"""

import numpy as np
from numpy import sqrt, sign, pi, cos, sin
import matplotlib.pyplot as plt

def FGsquircle(u, v):
    SMALL = 1e-15
    t0 = u**2+v**2
    t1 = (u**2+v**2)*(u**2+v**2-4*u**2*v**2)
    t2 = u**2+v**2
    t3 = (u**2+v**2)*(u**2+v**2-4*u**2*v**2)

    x = sign(u*v)/(v*sqrt(2.0)+SMALL)*sqrt(t0-sqrt(t1))
    y = sign(u*v)/(u*sqrt(2.0)+SMALL)*sqrt(t2-sqrt(t3))
    return x, y

R0 = 1.0 # radius of the disc
RMAX = 5.0 # the radius of the outer circle in the domain
NT = 360 # num of division in the theta direction
NR = 10 # num of radial divisions
r = [R0+(RMAX-R0)/NR*k for k in range(NR)] # the radii of circles
theta = np.array([2*pi/NT*k for k in range(NT+1)])

u = [r[k]*cos(theta) for k in range(NR)]
v = [r[k]*sin(theta) for k in range(NR)]

u = np.array(u)
v = np.array(v)

x, y = FGsquircle(u, v)

我得到了以下错误：

utils.py:21: RuntimeWarning: invalid value encountered in sqrt
  x = sign(u*v)/(v*sqrt(2.0)+SMALL)*sqrt(t0-sqrt(t1))
utils.py:22: RuntimeWarning: invalid value encountered in sqrt
  y = sign(u*v)/(u*sqrt(2.0)+SMALL)*sqrt(t2-sqrt(t3))

我很感谢你的帮助。

Answer 1

为什么这些公式把(1,0)，(0,1)，(-1,0)，(-1，-1)和(0,0)映射到点(0,0)？

在这两个表达式中，都用u或v除以，因此当两个表达式中的一个为0时，表达式将变为未定义(而不是零)。

我应该如何得到一个圆和一个正方形之间的中间形状，如下图所示。

只需映射一个半径小于1的圆。

是否有可能提供一个算法来从左边的网格中得到正确的网格？

你只需要转换网格点。细胞可以保持不变。

示例代码：

from numpy import sqrt, sign
import numpy
import matplotlib.pyplot as plt


def f(x):
    u, v = x
    alpha = sqrt(
        u ** 2
        + v ** 2
        - sqrt((u ** 2 + v ** 2) * (u ** 2 + v ** 2 - 4 * u ** 2 * v ** 2))
    )
    return numpy.array(
        [sign(u * v) / (v * sqrt(2)) * alpha, sign(u * v) / (u * sqrt(2)) * alpha]
    )


for r in numpy.linspace(0.1, 1.0, 10):
    theta = numpy.linspace(0.0, 2 * numpy.pi, 1000, endpoint=True)
    uv = r * numpy.array([numpy.cos(theta), numpy.sin(theta)])
    xy = f(uv)

    plt.plot(xy[0], xy[1], "-")

plt.gca().set_aspect("equal")
plt.show()

​

​