Проблема программирования с Mathematica

Я взаимодействую с внешней программой с Mathematica. Я создаю входной файл для внешней программы. Речь идет о преобразовании данных геометрии из графики, сгенерированной Mathematica, в предопределенный формат. Вот пример геометрии.

Рисунок 1

Рисунок 1

Геометрия может быть описана многими способами в Mathematica. Один из трудоемких способов заключается в следующем.

dat={{1.,-1.,0.},{0.,-1.,0.5},{0.,-1.,-0.5},{1.,-0.3333,0.},{0.,-0.3333,0.5},
{0.,-0.3333,-0.5},{1.,0.3333,0.},{0.,0.3333,0.5},{0.,0.3333,-0.5},{1.,1.,0.},
{0.,1.,0.5},{0.,1.,-0.5},{10.,-1.,0.},{10.,-0.3333,0.},{10.,0.3333,0.},{10.,1.,0.}};

Show[ListPointPlot3D[dat,PlotStyle->{{Red,PointSize[Large]}}],Graphics3D[{Opacity[.8],
Cyan,GraphicsComplex[dat,Polygon[{{1,2,5,4},{1,3,6,4},{2,3,6,5},{4,5,8,7},{4,6,9,7},
{5,6,9,8},{7,8,11,10},{7,9,12,10},{8,9,12,11},{1,2,3},{10,12,11},{1,4,14,13},
{4,7,15,14},{7,10,16,15}}]]}],AspectRatio->GoldenRatio]

При этом создается необходимая 3D-геометрия в формате GraphicsComplex MMA. введите здесь описание изображения

Эта геометрия описана как следующий входной файл для моей внешней программы.

# GEOMETRY
# x y z [m]
NODES 16
1. -1. 0.
0. -1. 0.5
0. -1. -0.5
1. -0.3333 0.
0. -0.3333 0.50. -0.3333 -0.5
1. 0.3333 0.
0. 0.3333 0.5
0. 0.3333 -0.5
1. 1. 0.
0. 1. 0.5
0. 1. -0.5
10. -1. 0.
10. -0.3333 0.
10. 0.3333 0.
10. 1. -0.
# type node_id1 node_id2 node_id3 node_id4  elem_id1 elem_id2 elem_id3 elem_id4
PANELS 14
1 1 4 5 2 4 2 10 0
1 2 5 6 3 1 5 3 10
1 3 6 4 1 2 6 10 0
1 4 7 8 5 7 5 1 0
1 5 8 9 6 4 8 6 2
1 6 9 7 4 5 9 3 0
1 7 10 11 8 8 4 11 0
1 8 11 12 9 7 9 5 11
1 9 12 10 7 8 6 11 0
2 1 2 3 1 2 3
2 10 12 11 9 8 7
10 4 1 13 14 1 3
10 7 4 14 15 4 6
10 10 7 15 16 7 9
# end of input file

Теперь описание, которое я получил из документации этой внешней программы, довольно короткое. Я цитирую это здесь.


  1. Первое ключевое слово NODES указывает общее количество узлов. После этой строки не должно быть комментариев или пустых строк. Следующие строки состоят из трех значений координат узлов x, y и z, а количество строк должно совпадать с количеством узлов.
  2. Следующее ключевое слово — ПАНЕЛЬ. Оно указывает, сколько панелей у нас есть. После этого у нас есть линии, определяющие каждую панель. Первое целое число определяет тип панели
  3. ID 1четырехугольная панель — определяется четырьмя узлами и четырьмя соседними панелями. Соседние панели — это панели, которые имеют одни и те же стороны (пара узлов) и необходимы для расчета скорости и давления (методы 1 и 2). Отсутствующие соседи (например, для панелей рядом с задней кромкой) заполняются значением 0 (см. Рисунок 1).
  4. ID 2треугольная панель — определяется тремя узлами и тремя соседними панелями.
  5. ID 10панель следа — четырехугольная панель, определенная четырьмя узлами и двумя (соседними) панелями, расположенными на задней кромке (панели, к которым применяется панель следа Кутта). состояние).
  6. Типы панелей 1 и 2 должны быть определены до типа 10 во входном файле. Важно отметить, что поверхность нормальна; порядок узлов, определяющих панели, должен быть против часовой стрелки. По правилу правой руки, если пальцы согнуты, чтобы следовать нумерации, большой палец покажет вектор нормали, который должен указывать «внешнюю» геометрию.

Вызов!!

Нам предоставляется 3D-модель CAD в файле с именем One.obj, и она отлично экспортируется в ММА.

cd = Import["One.obj"]

Результатом является объект MMA Graphics3D введите здесь описание изображения

Теперь я могу легко получить доступ к данным геометрии, так как MMA считывает их изнутри.

{ver1, pol1} = cd[[1]][[2]] /. GraphicsComplex -> List;
MyPol = pol1 // First // First;
Graphics3D[GraphicsComplex[ver1,MyPol],Axes-> True]

введите здесь описание изображения

  1. Как мы можем использовать информацию о вершинах и полигонах, содержащуюся в ver1 и pol1, и записать их в текстовый файл, как описано в примере входного файла выше. В этом случае у нас будут только панели типа ID2 (треугольные).
  2. Используя триангуляцию Mathematica, как найти площадь поверхности этого трехмерного объекта. Есть ли встроенная функция, которая может вычислять площадь поверхности в ММА?
  3. Сейчас нет необходимости создавать панель пробуждения или элементы типа ID10. Входной файл только с треугольными элементами будет в порядке.

Извините за такой длинный пост, но это головоломка, которую я пытаюсь решить в течение длительного времени. Надеюсь, что у некоторых из вас, экспертов, может быть правильное понимание, чтобы взломать его.

BR


person PlatoManiac    schedule 14.09.2011    source источник
comment
Обычно я стараюсь получать деньги только за то, что читаю такие задачи, как эта.   -  person Dr. belisarius    schedule 15.09.2011
comment
@belisarius К сожалению, я из гораздо менее привилегированного общества, где даже тщательно продуманное решение этой проблемы наверняка не принесет никакой оплаты за решатель. Но еще раз приношу свои извинения за столь пространный вопрос.   -  person PlatoManiac    schedule 15.09.2011
comment
Знаете ли вы, что Mathematica 7 (и, возможно, более ранние версии) поддерживает как импорт, так и экспорт графических файлов формата .obj?   -  person Codie CodeMonkey    schedule 21.09.2011
comment
@DeepYellow Я был в этом уверен. Но экспорт-импорт для файлов .obj не решает мою проблему.   -  person PlatoManiac    schedule 22.09.2011


Ответы (1)


Q1 и Q2 достаточно просты, чтобы вы могли отказаться от ярлыков «вызов» в своем вопросе. Вопрос 3 может потребовать некоторых разъяснений.

Q1

edges = cd[[1, 2, 1]];

polygons = cd[[1, 2, 2, 1, 1, 1]];

Обновление Q1

Основная проблема состоит в том, чтобы найти соседа каждого полигона. Это делает следующее:

(* Split every triangle in 3 edges, with nodes in each edge sorted *)
triangleEdges = (Sort /@ Subsets[#, {2}]) & /@ polygons;

(* Generate a list of edges *)
singleEdges = Union[Flatten[triangleEdges, 1]];

(* Define a function which, given an edge (node number list), returns the bordering  *)
(* triangle numbers. It's done by working through each of the triangles' edges       *)
ClearAll[edgesNeighbors]
edgesNeighbors[_] = {};
MapIndexed[(
   edgesNeighbors[#1[[1]]] = Flatten[{edgesNeighbors[#1[[1]]], #2[[1]]}];
   edgesNeighbors[#1[[2]]] = Flatten[{edgesNeighbors[#1[[2]]], #2[[1]]}];
   edgesNeighbors[#1[[3]]] = Flatten[{edgesNeighbors[#1[[3]]], #2[[1]]}];
   ) &, triangleEdges
];

(* Build a triangle relation table. Each '1' indicates a triangle relation *)
relations = ConstantArray[0, {triangleEdges // Length, triangleEdges // Length}];
Scan[
  (n = edgesNeighbors[##]; 
     If[Length[n] == 2, 
        {n1, n2} = n; 
        relations[[n1, n2]] = 1;  relations[[n2, n1]] = 1];
   ) &, singleEdges
]

MatrixPlot[relations]

отношения треугольника

(* Build a neighborhood list *)
triangleNeigbours = 
    Table[Flatten[Position[relations[[i]], 1]], {i,triangleEdges // Length}];

(* Test: Which triangles border on triangle number 1? *)
triangleNeigbours[[1]]

(* ==> {32, 61, 83} *)

(* Check this *)
polygons[[{1, 32, 61, 83}]]

(* ==> {{1, 2, 3}, {3, 2, 52}, {1, 3, 50}, {19, 2, 1}} *)
(* Indeed, they all share an edge with #1 *)

Для их вывода можно использовать функции вывода низкого уровня, описанные здесь. Я оставлю детали вам (это мой вызов вам).

Q2
Площадь крыла – это сумма площадей отдельных полигонов. Отдельные площади можно рассчитать следующим образом:

ClearAll[polygonArea];
polygonArea[pts_List] :=
 Module[{dtpts = Append[pts, pts[[1]]]},
   If[Length[pts] < 3, 
      0, 
      1/2 Sum[Det[{dtpts[[i]], dtpts[[i + 1]]}], {i, 1, Length[dtpts] - 1}]
   ]
 ]

на основе этой страницы Mathworld.

Область подписана BTW, поэтому вы можете использовать Abs.

ИСПРАВЛЕНИЕ
Приведенная выше функция площади может использоваться только для общих многоугольников в 2D. Для области треугольника в 3D можно использовать следующее:

ClearAll[polygonArea];
polygonArea[pts_List?(Length[#] == 3 &)] := 
    Norm[Cross[pts[[2]] - pts[[1]], pts[[3]] - pts[[1]]]]/2
person Sjoerd C. de Vries    schedule 14.09.2011
comment
На вопрос 1 вообще нельзя ответить, просто получив ребра и полигоны. Он уже извлечен, когда мы формируем GraphicsComplex. Основная проблема заключается в том, что нужно также предоставить информацию о том, как ребра в каждом полигоне связаны с другими полигонами. Посмотрите на столбцы элементов примера входного файла. Итак, нам нужна маркировка для каждого полигона. Проверьте отмеченные красным цветом многоугольники на рисунке 1. - person PlatoManiac; 15.09.2011
comment
Большое спасибо за материал для полигона. - person PlatoManiac; 15.09.2011
comment
@PlatoManiac Приношу свои извинения за слишком поверхностное прочтение вашего вопроса. Я нырну глубже, но это займет некоторое время. А пока вопрос: не должен ли 0 в строке «1 1 4 5 2 4 2 10 0», непосредственно следующей за PANEL, быть 12 (или даже 3)? Обе поверхности кажутся связанными с узлами 1 и 4 и, следовательно, являются соседними панелями. - person Sjoerd C. de Vries; 15.09.2011
comment
Извините, что отвечаю после столь долгого перерыва. Получил травму во время игры и действительно плохо провел время. Кстати, ваше решение все еще не работает для внешнего решателя. Я думаю, что заказ панелей происходит неправильно. Вот пример входного файла, который работает. Это сфера с 24 панелями. mediafire.com/?0hu5ql817h2ghiy . Здесь используются четырехугольные панели, что означает, что вместо треугольников у них четырехугольник. Надеюсь, вы мне поможете... - person PlatoManiac; 23.09.2011
comment
@PlatoManiac Без подробностей о том, как вы использовали мои функции для создания списка панелей, я не могу ничем помочь. Я предполагаю, что узлы панели должны быть упорядочены против часовой стрелки, и я предполагаю, что соседние плоскости должны быть упорядочены в том же порядке, что и ребра. Есть еще вопрос выше, на который вы не ответили. о ребрах, которые имеют более 1 соединительной панели, например, панели 1, 3 и 12? Вы действительно должны понять это сами. Детали — это не проблема Mathematica, и у меня, к сожалению, не так много свободного времени. - person Sjoerd C. de Vries; 24.09.2011
comment
если вы заинтересованы, проверьте stackoverflow. com/questions/7676032/, если вы можете дать больше информации. - person PlatoManiac; 10.10.2011