특정 3D 볼록 세트를 둘러싸고 새기는 Graph / Construct (John) 타원체

Nov 13 2020

Fritz John 의 유명한 정리 JohnEllipsoids 는 볼록한 몸체와 관련된 최소 및 최대 볼륨의 외접 및 내접 타원체를 알려줍니다.

이제, SpectraConvexity에 대한 Nathaniel Johnston의 대답에서 볼록이라고 주장되는 몸체 는 절대적으로 분리 가능한 2 큐 비트 상태의 정렬 된 스펙트럼 집합입니다. 이 세트는 제약 조건에 의해 정의됩니다.

1 > x && x >= y && y >= z && z >= 1 - x - y - z >= 0 && 
 x - z < 2 Sqrt[y (1 - x - y - z)

관련 John 타원체의 명시적인 구성을 제외하고 (그리고 검색을 지원할 수 있음) 그래픽 탐색 (RegionPlot3D, Ellipsoid 및 RegionMeasure를 포함한 Mathematica의 많은 도구 사용)을 찾아서 근사화하는 것이 흥미로울 것입니다.

두 개의 다른 관련 관심 세트가 있으며, 또한 검사중인 볼록 체 내부에 포함되어 있습니다. 이것들은 제약에 의해 주어집니다

 1 > x && x >= y && y >= z && z >= 1 - x - y - z >= 0 && 
  x^2 + y^2 + (1 - x - y - z)^2 + z^2 < 3/8]

 1 > x && x >= y && y >= z && z >= 1 - x - y - z >= 0 && 
  x^2 + y^2 + (1 - x - y - z)^2 + z^2 < 1/3]

이것이 문제의 타원체 일 수 있으며, 그렇지 않다면 기하학적 모양은 무엇일까요?

다음은 위에 주어진 세 가지 제약 조건과 관련된 세 세트의 RegionPlot3D를 사용하는 플롯입니다. Ellipsoid 명령을 이러한 그래픽에 통합하고 볼륨을 찾는 데 RegionMeasure를 사용하려고합니다.

ThreeSetPlot

다음은 Mathematica를 사용하여 제기 된 질문을 탐색하기위한 매우 거칠고 예비적인 노력입니다. 플롯은 절대적으로 분리 가능한 2- 큐 비트 상태의 정렬 된 스펙트럼의 볼록 세트를 둘러싸는 타원체 "가까운"을 보여줍니다. 그러나 최소한의 부피의 외접 타원체를 구성하는 것은 매우 어려운 일입니다. 그 존재는 프리츠 존 정리에 의해 주어집니다. ( "John 타원체는 계산하기 어렵다" M- 타원체 .) 최소화 할 목적 함수는 무엇입니까? 또한 물론 최대 볼륨 문제의 "이중"으로 새겨진 타원체가 있습니다.

EllipsoidPlot

이 퀘스트에서 Ellipsoid 및 RegionMeasure 명령 (다른 명령 중)을 활용할 수 있는지 확실하지 않습니다.

절대적으로 분리 가능한 2- 큐 비트 상태의 정렬 된 스펙트럼의 볼록 세트 (여기서 주요 관심 사항)의 (유클리드) 부피는 다음과 같습니다. $\approx 0.00227243$ (정확한 값을 얻을 수 있어야 함) 마지막 플롯에 표시된 타원체의 부피는 $\frac{\pi }{150 \sqrt{15}} \approx 0.0054077$.

답변

3 GeorgeVarnavides Nov 16 2020 at 21:25

더 많은 확장 주석이지만 BoundingRegion기능을 인식하지 못한 경우 :

rm=RegionMember[ImplicitRegion[conditionABS,{x,y,z}]];
pts=RandomVariate[UniformDistribution[N[{{1/4,1/8 (2+Sqrt[6])},{1/24 (6-Sqrt[6]),1/8 (2+Sqrt[2])},{1/8 (2-Sqrt[2]),1/3}}]],10^5];
insidePts=Select[pts,rm];
fastEllipsoid=BoundingRegion[insidePts,"FastEllipsoid"]
RegionMeasure[fastEllipsoid]
Graphics3D[{{Opacity[0.5],fastEllipsoid},Point[insidePts]}]

문서 자체는 다음과 같이 경고합니다.

"FastEllipsoid" gives a bounding Ellipsoid, not necessarily with the minimal volume
1 DanielHuber Nov 18 2020 at 23:17

다음은 4 개의 점이 주어지면 타원체를 외 접하는 방법입니다.

조건 형성 :

cond = 1 > x && x >= y && y >= z && z >= 1 - x - y - z >= 0 && 
   x - z <= 2 Sqrt[y (1 - x - y - z)];

, 여기서 <를 <=로 변경 한 경우 먼저 Minimizeand를 사용하여 4 개의 극한 점을 결정합니다 Maximize. 예 Maximize[{y, cond}, {x, y, z}]. 이것은 4 가지 포인트를줍니다 :

pts={{1/3, 1/3, 1/3}, {1/4, 1/4, 1/4}, {1/2, 1/6, 1/
  6}, {1/8 (2 + Sqrt[2]), 1/8 (2 + Sqrt[2]), 
  1/2 (1 + 1/4 (-2 - Sqrt[2]))}}//N;

다음으로 우리는 가장 멀리 떨어진 두 지점을 결정합니다. 우리의 경우 이것은 pts[[2]]pts[[4]]입니다. 우리는 선의 중간 점을 선택 pts[[2]]하는 pts[[4]]COM (질량 중심) : 우리의 타원체의 중심. 그리고 거리의 절반이 가장 큰 절반 축이 될 것입니다 : 타원체의 a3 :

com = (pts[[2]] + pts[[4]])/2 // N;
a3 = Norm[pts[[2]] - pts[[4]]]/2 // N;

다음 계산을 더 쉽게하기 위해 점이 원점에 있도록 변환합니다. 그런 다음 반축 a3가 z 방향을 가리 키도록 좌표계를 회전합니다.

pts1 = (# - com) & /@ pts // N;
pts2 = (r2 = RotationMatrix[{pts1[[2]] - pts1[[4]], {0, 0, 1}}]).# & /@
    pts1;

이제 점 1 또는 3 (이 경우 점 3)이 원점에서 더 먼 지점을 결정하고이 점이 yz 평면에 놓 이도록 z 축을 중심으로 회전합니다.

pts3 = (r3 = 
       RotationMatrix[
        ArcTan[pts2[[3, 1]], pts2[[3, 2]]], {0, 0, 1}]).# & /@ pts2;

다음으로 점 3이 타원 위에 놓 이도록 y 축을 따라 절반 축을 결정하고 yz 평면이 타원체를 잘라냅니다.

a2 = Sqrt[pts3[[3, 2]]^2/(1 - (pts3[[3, 3]]/a3)^2)]

이제 우리는 마지막 점 1이 타원체 위에 놓 이도록 x 좌표 방향으로 절반 축 a1을 결정합니다.

a1 = Sqrt[
  pts3[[1, 1]]^2/(1 - (pts3[[1, 2]]/a2)^2 - (pts3[[1, 3]]/a3)^2)]

이제 새 좌표계에서 타원체와 변환 된 점을 그리는 데 필요한 모든 데이터가 있습니다.

마지막으로 타원체 공식을 이전 좌표에 작성하여 타원체를 원래 좌표로 다시 변환해야합니다.

fun[{x_, y_, z_}] = Total[((r3.r2.({x, y, z} - com))/{a1, a2, a3})^2];

이제 원래 좌표에 타원체를 그릴 수 있습니다.

Show[
 ContourPlot3D[
  fun[{x, y, z}] == 1, {x, .1, .6}, {y, .1, .55}, {z, -.1, .4}, 
  AxesLabel -> {"x", "y", "z"}, ContourStyle -> Opacity[0.5], 
  Mesh -> None]
 , Graphics3D[{PointSize[0.03], Point[pts]}, Axes -> True]
 , reg
 ]

편의를 위해 모든 코드를 한 조각으로 :

cond = 1 > x && x >= y && y >= z && z >= 1 - x - y - z >= 0 && 
  x - z <= 2 Sqrt[y (1 - x - y - z)]; pts = {{1/3, 1/3, 1/3}, {1/4, 
    1/4, 1/4}, {1/2, 1/6, 1/6}, {1/8 (2 + Sqrt[2]), 1/8 (2 + Sqrt[2]),
     1/2 (1 + 1/4 (-2 - Sqrt[2]))}} // N;
com = (pts[[2]] + pts[[4]])/2 // N;
a3 = Norm[pts[[2]] - pts[[4]]]/2 // N;
pts1 = (# - com) & /@ pts // N;
pts2 = (r2 = RotationMatrix[{pts1[[2]] - pts1[[4]], {0, 0, 1}}]).# & /@
    pts1;
pts3 = (r3 = 
       RotationMatrix[
        ArcTan[pts2[[3, 1]], pts2[[3, 2]]], {0, 0, 1}]).# & /@ pts2;
a2 = Sqrt[pts3[[3, 2]]^2/(1 - (pts3[[3, 3]]/a3)^2)]
a1 = Sqrt[
  pts3[[1, 1]]^2/(1 - (pts3[[1, 2]]/a2)^2 - (pts3[[1, 3]]/a3)^2)]
Show[ContourPlot3D[
  Total[({x, y, z}/{a1, a2, a3})^2] == 
   1, {x, -.2, .2}, {y, -.25, .25}, {z, -.2, .21}, 
  AxesLabel -> {"x", "y", "z"}], 
 Graphics3D[{PointSize[0.03], Point[pts3], 
   Line[{pts3[[2]], pts3[[4]]}], Line[{{0, 0, 0}, pts3[[3]]}]}, 
  Axes -> True]]
fun[{x_, y_, z_}] = Total[((r3.r2.({x, y, z} - com))/{a1, a2, a3})^2];
reg = RegionPlot3D[
   cond, {x, 1/4, 1/8 (2 + Sqrt[6])}, {y, 1/24 (6 - Sqrt[6]), 
    1/8 (2 + Sqrt[2])}, {z, 1/3, 1/8 (2 - Sqrt[2])}, 
   PlotPoints -> 100];
Show[ContourPlot3D[
  fun[{x, y, z}] == 1, {x, .1, .6}, {y, .1, .55}, {z, -.1, .4}, 
  AxesLabel -> {"x", "y", "z"}, ContourStyle -> Opacity[0.5], 
  Mesh -> None], 
 Graphics3D[{PointSize[0.03], Point[pts]}, Axes -> True], reg]
1 Dominic Nov 20 2020 at 01:16

내접 타원체를 찾는 방법은 다음과 같습니다.

  1. 먼저 경계 영역에서 100 만 개의 포인트를 생성하고 conditionABS 내부의 포인트를 선택합니다.
  2. ConvexHull of points 생성,
  3. LinearOptimization을 사용하여 폴리 토프를 구성하고,
  4. ConicOptimization을 실행하여 타원체를 찾습니다.

그러나 번역 벡터 d의 부호를 변경해야했습니다. 타원체의 부피 (기계 정밀도)는 0.001442입니다. 자세한 내용 은 다각형의 가장 큰 타원에 대한 PF링크를 참조하십시오 .

conditionABS = 
  1 > x && x >= y && y >= z && z >= 1 - x - y - z >= 0 && 
   x - z < 2 Sqrt[y (1 - x - y - z)];
(*
  generate one million points in bounding region and select points \
inside conditionABS
*)
rm = RegionMember[ImplicitRegion[conditionABS, {x, y, z}]];
pts = RandomVariate[
   UniformDistribution[
    N[{{1/4, 1/8 (2 + Sqrt[6])}, {1/24 (6 - Sqrt[6]), 
       1/8 (2 + Sqrt[2])}, {1/8 (2 - Sqrt[2]), 1/3}}]], 10^6];
insidePts = Select[pts, rm];
(*
 generate a convex hull for the points
*)
mesh = ConvexHullMesh[insidePts];
meshP = Show[Graphics3D@{Opacity[0.02, Blue], mesh}, Axes -> True]
(*
 Obtain polytope inequalities to represent the region 
*)
{A, b} = LinearOptimization[0, {}, x \[Element] mesh, 
   "LinearInequalityConstraints"];
Length@A
(* 
use ConicOptimization to find max ellipsoid
*)
polyA = A;
polyB = b; constraints = 
 Table[Norm[polyA[[i]].c] + polyA[[i]].d <= polyB[[i]], {i, 
   Length[polyA]}]; {cEllipse, dEllipse} = {c, d} /. 
  ConicOptimization[-Tr[c], 
   constraints, {c \[Element] Matrices[{3, 3}], d}]
(*
 compute volume
*)
eVolume = 
 4 Pi/3 (Norm[cEllipse[[All, 1]]] Norm[cEllipse[[All, 2]]] 
    Norm[cEllipse[[All, 3]]])
(*
 construct affine paramaterization for ellipsoid
*)
aFine[d_, m_, \[Theta]_, \[Phi]_] := 
  d + m[[All, 1]] Cos[\[Theta]] Cos[\[Phi]] + 
   m[[All, 2]] Cos[\[Theta]] Sin[\[Phi]] + m[[All, 3]] Sin[\[Theta]];
(*
 generate plots
*)
pp1 = ParametricPlot3D[
  aFine[-dEllipse, cEllipse, t, p], {t, -Pi/2, Pi/2}, {p, 0, 2 Pi}]
Show[{meshP, pp1}, Axes -> True, BoxRatios -> {1, 1, 1}]

Dominic Nov 15 2020 at 01:40

예를 들어, 아래 지역에 대한 최소 및 최대 볼륨의 외접 및 내접 엘 립소 이드를 찾으려고합니까 (코드를 약간 수정했습니다)?

    conditionABS = 
      1 > x && x >= y && y >= z && z >= 1 - x - y - z >= 0 && 
       x - z < 2 Sqrt[y (1 - x - y - z)];
    RegionPlot3D[conditionABS, {x, 1/4, 1/8 (2 + Sqrt[6])}, {y, 
  1/24 (6 - Sqrt[6]), 1/8 (2 + Sqrt[2])}, {z, 1/3, 
  1/8 (2 - Sqrt[2])}, 
 AxesLabel -> {Style["x", 16, Bold, Black], 
   Style[ "y", 16, Bold, Black], Style[ "z", 16, Bold, Black]}, 
 PlotPoints -> 100]

PaulB.Slater Nov 16 2020 at 23:15

그 자체로 완전한 답은 아니지만 단순히 두 가지 문제의 증폭입니다.

첫째로, 문제의 볼록 세트의 면적 / 체적 비율이 6 인 AreaVolumeRatio 가 나타납니다 . 만약 그렇다면, 알려진 볼록 세트 패밀리 중에서 세트의 특성을 식별하는 데 도움이 될 수 있습니다.

둘째, 중심 ( "문제가되는") 불평등 제약

 x - z < 2 Sqrt[y (1 - x - y - z)

의 양의 반 정확성과 동일합니다. $2 \times 2$ 매트릭스,

P = {{2 (1-x-y-z), -x + z}, {-x + z, 2 y}},

Nathaniel Johnston (R. Hildebrand의 인용 작업)이 PositiveSemidefiniteness 에 대한 답변 끝에 지적한대로 .

자,이 행렬 P가 파이썬 코드 "내부 및 외부 Löwner-John Ellipsoids"에서 필요한 것 (P로도 표시됨) 일까요? 사용자 Dominic이 언급 한 PythonCode 는이 질문에 대한 그의 의견 중 하나입니까?

그렇다면 (이 시점에서 P가 어떤 식 으로든 폴리 토프를 나타내는 것이 다소 회의적입니다), 파이썬 코드 구현을 시도했습니다 (ConfigurePythonForExternalEvaluate를 활용하면 보일 것입니다).

아마도 양의 유사 유한 조건이 제약 조건을 정의하는 행렬을 구성 할 수 있습니다.

1 > x && x >= y && y >= z && z >= 1 - x - y - z >= 0 &&  x - z < 2 Sqrt[y (1 - x - y - z)

질문의 시작 부분에 주어집니다. 아마도 그러한 행렬은 파이썬 코드에 입력하기에 적합한 행렬 일 것입니다.

그러한 행렬을 얻는 다소 사소한 방법은 표시된 $2 \times 2$ 행렬 (주된 부등식 제약 조건 생성)

 P=  {{2 (1 - x - y - z), -x + z}, {-x + z, 2 y}}

원래 null의 상단 모서리 $6 \times 6$ 행렬을 만들고 1-x, xy, yz 및 z- (1-xyz) 항목을 나머지 4 개의 대각선 위치에 삽입합니다.

PaulB.Slater Nov 18 2020 at 20:16

다시 말하지만, 문제가되는 두 개의 타원체의 요청 된 구성이 아니라이 문제에 대해 관심이있는 두 가지 개발 사항에주의를 기울이려는 노력입니다.

첫째로, 사용자 Dominic은 여기에있는 주석에서 "Inner and outer Löwner-John Ellipsoids" Mosekpythoncode 라는 제목의 정교한 (Mosek- 소프트웨어 패키지) 파이썬 코드를 언급했습니다 . 파이썬 사용자도 아니고 기본 최적화 절차의 전문가도 아니고 pythonQuestion이라는 질문을 올렸습니다 .

Mosek과 분명히 제휴 한 사용자 Michal Adamaszek은 다음과 같이 언급했습니다.

"Mosek 코드는 폴리 토프 P에 새겨진 타원체를위한 것입니다. P가 볼록하지만 폴리 토프가 아니라면"for all u "부분을 더 관리하기 쉬운 것으로 다시 작성할 수 있는지 여부에 따라 가능할 수도 있고 아닐 수도 있습니다. 세트에 SDP 표현이있는 것 같습니다. 따라서 최소한 충분히 많은 u를 샘플링하고 해당 Cu + d를 P에 속하도록 제한하여 근사치를 얻을 수 있습니다. " ( "SDP 표현"은$6 \times 6$ 매트릭스

{{2 (1 - x - y - z), -x + z, 0, 0, 0, 0}, {-x + z, 2 y, 0, 0, 0, 0}, {0, 0, 1 - x, 0, 0, 0}, {0, 0, 0, x - y, 0, 0}, {0, 0, 0, 0, y - z, 0}, {0, 0, 0, 0, 0, -1 + x + y + 2 z}}

내 이전 "답변"에서 구성되었습니다.

나는 대답했다 :

"Michal Adamaszek에게 정말 감사합니다. 질문을 제기하면서 얻고 자하는 전문 지식입니다. 저는 Python 사용자가 아니기 때문에 제안 된 접근 방식을 구현하는 데 더 많은 노력을 기울여야 할 수도 있습니다.이 시점에서 저는 P가 폴리 토프인지 아닌지에 대한 확고한 지식이 없습니다. "너무 좋았 기 때문에 사실 일 수 없습니다."라고 생각합니다. 제한된 이해 내에서 P가 폴리 토프인지 여부는 그 자체로 어려운 질문입니다. "

내가 여기서 강조하고 싶은 두 번째 개발에 관해서는, 그것은 지금 설정이가 획득 (사용자가 JimB)의 즉각적인 결과이다 6입니다 ( "스펙트럼을 주문") 볼록의 면적 / 부피 비율 것으로 알려져있다 AreaVolumeRatio 의 세트의 볼륨

1/576 (8 - 6 Sqrt[2] - 9 Sqrt[2] π + 24 Sqrt[2] ArcCos[1/3])  ,

이 표현의 6 배에 해당하는 지역의 이전 발견과 결합됩니다.