Генерация меня sh из неявной функции - PullRequest
Новая
       40

Генерация меня sh из неявной функции

0 голосов
/ 08 февраля 2020

Я хочу использовать Python (или Python 3) для создания тома (3D) me sh из неявной функции: 

def func(x,y,z):
    q = 0.25
    mu = q/(1.+q)
    return -(1-mu)*pow(x*x + y*y + z*z,-1./2.)  - mu*pow(pow(x-1,2) + y*y + z*z,-1./2.) - 0.5*(pow(x-mu,2) + y*y) + 1.9023266381531847

Эта функция имеет сложную изоповерхность, но Я хочу ограничить поверхность между x = -0.615 и x = 1.4, а y = -0.6, y = 0.6, ограничений по направлению z нет, но интересная часть находится между z = +/- 1 .

Я пытался pygalme sh, но я не мог заставить их пример адаптироваться к моей функции. Вылетает мое ядро ​​Python без вывода. Возможно ли получить pygalme sh для этого? Если нет, то что было бы лучше?

1 Ответ

1 голос
/ 10 февраля 2020

Только для справки, для меня не получается sh без вывода: 

import numpy
import pygalmesh


class GrimReaper(pygalmesh.DomainBase):
    def __init__(self):
        super().__init__()

    def eval(self, x):
        q = 0.25
        mu = q / (1.0 + q)
        x, y, z = x
        return (
            -(1 - mu) / numpy.sqrt(x ** 2 + y ** 2 + z ** 2)
            - mu / numpy.sqrt((x - 1) ** 2 + y ** 2 + z ** 2)
            - 0.5 * ((x - mu) ** 2 + y ** 2)
            + 1.9023266381531847
        )

    def get_bounding_sphere_squared_radius(self):
        return 2.0


d = GrimReaper()
mesh = pygalmesh.generate_mesh(d, cell_size=0.1)

mesh.write("out.xmf")

Inserting protection balls...
  refine_balls = true
  min_balls_radius = 0
  min_balls_weight = 0
insert_corners() done. Nb of points in triangulation: 0
insert_balls_on_edges() done. Nb of points in triangulation: 0
refine_balls() done. Nb of points in triangulation: 0

construct initial points (nb_points: 12)
s.py:17: RuntimeWarning: divide by zero encountered in double_scalars
  + 1.9023266381531847
12/12 initial point(s) found...
Start surface scan...Scanning triangulation for bad facets (sequential) - number of finite facets = 50...
Number of bad facets: 0
scanning edges (lazy)
scanning vertices (lazy)
end scan. [Bad facets:0]

Refining Surface...
Legend of the following line: (#vertices,#steps,#facets to refine,#tets to refine)
(12,0,0,0)
Total refining surface time: 1.90735e-05s

Start volume scan...Scanning triangulation for bad cells (sequential)... 20 cells scanned, done.
Number of bad cells: 1
end scan. [Bad tets:1]

Refining...
Legend of the following line: (#vertices,#steps,#facets to refine,#tets to refine)
(23,11,0,70) (18839.3 vertices/s)Segmentation fault (core dumped)

Он отлично работает без термина - 0.5 * ((x - mu) ** 2 + y ** 2).

. Segfault указывает на проблему в CGAL. Возможно, полезно подать ошибку там.

...