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

Question

上下文：

我正在研究计算流体动力学中的网格。我想围绕一个圆生成一个结构化的网格。我的计划是生成一个极性网格（下图为左网格），然后使用FG公式获得最终网格（下图为右网格）。我正在使用this article第4页中的FG方法将光盘映射到square。不幸的是，本文没有提到如何处理公式中的奇点。这是表达式：

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时，表达式变为undefined（不为零）。

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

仅绘制半径小于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()