MFEM v4.7.0
Finite element discretization library
Loading...
Searching...
No Matches
fmsconvert.cpp
Go to the documentation of this file.
1// Copyright (c) 2010-2024, Lawrence Livermore National Security, LLC. Produced
2// at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3// LICENSE and NOTICE for details. LLNL-CODE-806117.
4//
5// This file is part of the MFEM library. For more information and source code
6// availability visit https://mfem.org.
7//
8// MFEM is free software; you can redistribute it and/or modify it under the
9// terms of the BSD-3 license. We welcome feedback and contributions, see file
10// CONTRIBUTING.md for details.
11
12#include "../config/config.hpp"
13
14#ifdef MFEM_USE_FMS
15
16#include "fmsconvert.hpp"
17#include <climits>
18
19using std::endl;
20
21// #define DEBUG_FMS_MFEM 1
22// #define DEBUG_MFEM_FMS 1
23
24namespace mfem
25{
26
27static inline int
28FmsBasisTypeToMfemBasis(FmsBasisType b)
29{
30 int retval = -1;
31 switch (b)
32 {
33 case FMS_NODAL_GAUSS_OPEN:
34 retval = mfem::BasisType::GaussLegendre;
35 break;
36 case FMS_NODAL_GAUSS_CLOSED:
37 retval = mfem::BasisType::GaussLobatto;
38 break;
39 case FMS_POSITIVE:
40 retval = mfem::BasisType::Positive;
41 break;
42 case FMS_NODAL_UNIFORM_OPEN:
43 retval = mfem::BasisType::OpenUniform;
44 break;
45 case FMS_NODAL_UNIFORM_CLOSED:
46 retval = mfem::BasisType::ClosedUniform;
47 break;
48 case FMS_NODAL_CHEBYSHEV_OPEN:
49 case FMS_NODAL_CHEBYSHEV_CLOSED:
50 mfem::out <<
51 "FMS_NODAL_CHEBYSHEV_OPEN, FMS_NODAL_CHEBYSHEV_CLOSED need conversion to MFEM types."
52 << endl;
53 break;
54 }
55 return retval;
56}
57
58// The following function is unused (for now), so it is commented out to
59// suppress compilation warning.
60#if 0
61/**
62@brief Get the order and layout of the field.
63*/
64static int
65FmsFieldGetOrderAndLayout(FmsField f, FmsInt *f_order, FmsLayoutType *f_layout)
66{
67 int err = 0;
68 FmsFieldDescriptor fd;
69 FmsLayoutType layout;
70 FmsScalarType data_type;
71 const void *data = nullptr;
72 FmsInt order = 0;
73
74 FmsFieldGet(f, &fd, NULL, &layout, &data_type,
75 &data);
76
77 FmsFieldDescriptorType f_fd_type;
78 FmsFieldDescriptorGetType(fd, &f_fd_type);
79 if (f_fd_type != FMS_FIXED_ORDER)
80 {
81 err = 1;
82 }
83 else
84 {
85 FmsFieldType field_type;
86 FmsBasisType basis_type;
87 FmsFieldDescriptorGetFixedOrder(fd, &field_type,
88 &basis_type, &order);
89 }
90
91 *f_order = order;
92 *f_layout = layout;
93
94 return err;
95}
96#endif
97
98/**
99@brief This function converts an FmsField to an MFEM GridFunction.
100@note I took some of the Pumi example code from the mesh conversion function
101 that converted coordinates and am trying to make it more general.
102 Coordinates are just another field so it seems like a good starting
103 point. We still have to support a bunch of other function types, etc.
104*/
105template <typename DataType>
106int
107FmsFieldToGridFunction(FmsMesh fms_mesh, FmsField f, Mesh *mesh,
108 GridFunction &func, bool setFE)
109{
110 int err = 0;
111
112 // NOTE: transplanted from the FmsMeshToMesh function
113 // We should do this work once and save it.
114 //--------------------------------------------------
115 FmsInt dim, n_elem, space_dim;
116
117 // Find the first component that has coordinates - that will be the new mfem
118 // mesh.
119 FmsInt num_comp;
120 FmsMeshGetNumComponents(fms_mesh, &num_comp);
121 FmsComponent main_comp = NULL;
122 FmsField coords = NULL;
123 for (FmsInt comp_id = 0; comp_id < num_comp; comp_id++)
124 {
125 FmsComponent comp;
126 FmsMeshGetComponent(fms_mesh, comp_id, &comp);
127 FmsComponentGetCoordinates(comp, &coords);
128 if (coords) { main_comp = comp; break; }
129 }
130 if (!main_comp) { return 1; }
131 FmsComponentGetDimension(main_comp, &dim);
132 FmsComponentGetNumEntities(main_comp, &n_elem);
133 FmsInt n_ents[FMS_NUM_ENTITY_TYPES];
134 FmsInt n_main_parts;
135 FmsComponentGetNumParts(main_comp, &n_main_parts);
136 for (FmsInt et = FMS_VERTEX; et < FMS_NUM_ENTITY_TYPES; et++)
137 {
138 n_ents[et] = 0;
139 for (FmsInt part_id = 0; part_id < n_main_parts; part_id++)
140 {
141 FmsInt num_ents;
142 FmsComponentGetPart(main_comp, part_id, (FmsEntityType)et, NULL, NULL,
143 NULL, NULL, &num_ents);
144 n_ents[et] += num_ents;
145 }
146 }
147 //--------------------------------------------------
148
149 // Interrogate the field.
150 FmsFieldDescriptor f_fd;
151 FmsLayoutType f_layout;
152 FmsScalarType f_data_type;
153 const void *f_data;
154 FmsFieldGet(f, &f_fd, &space_dim, &f_layout, &f_data_type,
155 &f_data);
156 // FmsFieldGet(coords, NULL, &space_dim, NULL, NULL, NULL);
157
158 FmsInt f_num_dofs;
159 FmsFieldDescriptorGetNumDofs(f_fd, &f_num_dofs);
160
161 // Access FMS data through this typed pointer.
162 auto src_data = reinterpret_cast<const DataType *>(f_data);
163
164 FmsFieldDescriptorType f_fd_type;
165 FmsFieldDescriptorGetType(f_fd, &f_fd_type);
166 if (f_fd_type != FMS_FIXED_ORDER)
167 {
168 return 9;
169 }
170 FmsFieldType f_field_type;
171 FmsBasisType f_basis_type;
172 FmsInt f_order;
173 FmsFieldDescriptorGetFixedOrder(f_fd, &f_field_type,
174 &f_basis_type, &f_order);
175
176 if (f_field_type != FMS_CONTINUOUS && f_field_type != FMS_DISCONTINUOUS &&
177 f_field_type != FMS_HDIV)
178 {
179 return 10;
180 }
181
182 int btype = FmsBasisTypeToMfemBasis(f_basis_type);
183 if (btype < 0)
184 {
185 mfem::out << "\tInvalid BasisType: " << BasisType::Name(btype) << std::endl;
186 return 11;
187 }
188
189 //------------------------------------------------------------------
190 if (setFE)
191 {
192 // We could assemble a name based on fe_coll.hpp rules and pass to
193 // FiniteElementCollection::New()
194
195 mfem::FiniteElementCollection *fec = nullptr;
196 switch (f_field_type)
197 {
198 case FMS_DISCONTINUOUS:
199 fec = new L2_FECollection(f_order, dim, btype);
200 break;
201 case FMS_CONTINUOUS:
202 fec = new H1_FECollection(f_order, dim, btype);
203 break;
204 case FMS_HDIV:
205 fec = new RT_FECollection(f_order, dim);
206 break;
207 case FMS_HCURL:
208 case FMS_DISCONTINUOUS_WEIGHTED:
209 MFEM_ABORT("Field types FMS_HCURL and FMS_DISCONTINUOUS_WEIGHTED"
210 " are not supported yet.");
211 break;
212 }
213
214 int ordering = (f_layout == FMS_BY_VDIM) ? Ordering::byVDIM : Ordering::byNODES;
215 auto fes = new FiniteElementSpace(mesh, fec, space_dim, ordering);
216 func.SetSpace(fes);
217 }
218 //------------------------------------------------------------------
219 const FmsInt nstride = (f_layout == FMS_BY_VDIM) ? space_dim : 1;
220 const FmsInt vstride = (f_layout == FMS_BY_VDIM) ? 1 : f_num_dofs;
221
222 // Data reordering to store the data into func.
223 if ((FmsInt)(func.Size()) != f_num_dofs*space_dim)
224 {
225 return 12;
226 }
227
228 mfem::FiniteElementSpace *fes = func.FESpace();
229 const int vdim = fes->GetVDim();
230 const mfem::FiniteElementCollection *fec = fes->FEColl();
231 const int vert_dofs = fec->DofForGeometry(mfem::Geometry::POINT);
232 const int edge_dofs = fec->DofForGeometry(mfem::Geometry::SEGMENT);
233 const int tri_dofs = fec->DofForGeometry(mfem::Geometry::TRIANGLE);
234 const int quad_dofs = fec->DofForGeometry(mfem::Geometry::SQUARE);
235 const int tet_dofs = fec->DofForGeometry(mfem::Geometry::TETRAHEDRON);
236 const int hex_dofs = fec->DofForGeometry(mfem::Geometry::CUBE);
237 int ent_dofs[FMS_NUM_ENTITY_TYPES];
238 ent_dofs[FMS_VERTEX] = vert_dofs;
239 ent_dofs[FMS_EDGE] = edge_dofs;
240 ent_dofs[FMS_TRIANGLE] = tri_dofs;
241 ent_dofs[FMS_QUADRILATERAL] = quad_dofs;
242 ent_dofs[FMS_TETRAHEDRON] = tet_dofs;
243 ent_dofs[FMS_HEXAHEDRON] = hex_dofs;
244 FmsInt fms_dof_offset = 0;
245 int mfem_ent_cnt[4] = {0,0,0,0}; // mfem entity counters, by dimension
246 int mfem_last_vert_cnt = 0;
247 mfem::HashTable<mfem::Hashed2> mfem_edge;
248 mfem::HashTable<mfem::Hashed4> mfem_face;
249 if (dim >= 2 && edge_dofs > 0)
250 {
251 mfem::Array<int> ev;
252 for (int i = 0; i < mesh->GetNEdges(); i++)
253 {
254 mesh->GetEdgeVertices(i, ev);
255 int id = mfem_edge.GetId(ev[0], ev[1]);
256 if (id != i) { return 13; }
257 }
258 }
259 if (dim >= 3 &&
260 ((n_ents[FMS_TRIANGLE] > 0 && tri_dofs > 0) ||
261 (n_ents[FMS_QUADRILATERAL] > 0 && quad_dofs > 0)))
262 {
263 mfem::Array<int> fv;
264 for (int i = 0; i < mesh->GetNFaces(); i++)
265 {
266 mesh->GetFaceVertices(i, fv);
267 if (fv.Size() == 3) { fv.Append(INT_MAX); }
268 // HashTable uses the smallest 3 of the 4 indices to hash Hashed4
269 int id = mfem_face.GetId(fv[0], fv[1], fv[2], fv[3]);
270 if (id != i) { return 14; }
271 }
272 }
273
274 // Loop over all parts of the main component
275 for (FmsInt part_id = 0; part_id < n_main_parts; part_id++)
276 {
277 // Loop over all entity types in the part
278 for (FmsInt et = FMS_VERTEX; et < FMS_NUM_ENTITY_TYPES; et++)
279 {
280 FmsDomain domain;
281 FmsIntType ent_id_type;
282 const void *ents;
283 const FmsOrientation *ents_ori;
284 FmsInt num_ents;
285 FmsComponentGetPart(main_comp, part_id, (FmsEntityType)et, &domain,
286 &ent_id_type, &ents, &ents_ori, &num_ents);
287 if (num_ents == 0) { continue; }
288 if (ent_dofs[et] == 0)
289 {
290 if (et == FMS_VERTEX) { mfem_last_vert_cnt = mfem_ent_cnt[et]; }
291 mfem_ent_cnt[FmsEntityDim[et]] += num_ents;
292 continue;
293 }
294
295 if (ents != NULL &&
296 (ent_id_type != FMS_INT32 && ent_id_type != FMS_UINT32))
297 {
298 return 15;
299 }
300 if (ents_ori != NULL)
301 {
302 return 16;
303 }
304
305 if (et == FMS_VERTEX)
306 {
307 const int mfem_dof_offset = mfem_ent_cnt[0]*vert_dofs;
308 for (FmsInt i = 0; i < num_ents*vert_dofs; i++)
309 {
310 for (int j = 0; j < vdim; j++)
311 {
312 const int idx = i*nstride+j*vstride;
313 func(mfem_dof_offset*nstride+idx) =
314 static_cast<double>(src_data[fms_dof_offset*nstride+idx]);
315 }
316 }
317 fms_dof_offset += num_ents*vert_dofs;
318 mfem_last_vert_cnt = mfem_ent_cnt[et];
319 mfem_ent_cnt[0] += num_ents;
320 continue;
321 }
322 mfem::Array<int> dofs;
323 if (FmsEntityDim[et] == dim)
324 {
325 for (FmsInt e = 0; e < num_ents; e++)
326 {
327 fes->GetElementInteriorDofs(mfem_ent_cnt[dim]+e, dofs);
328 for (int i = 0; i < ent_dofs[et]; i++, fms_dof_offset++)
329 {
330 for (int j = 0; j < vdim; j++)
331 {
332 func(fes->DofToVDof(dofs[i],j)) =
333 static_cast<double>(src_data[fms_dof_offset*nstride+j*vstride]);
334 }
335 }
336 }
337 mfem_ent_cnt[dim] += num_ents;
338 continue;
339 }
340 const FmsInt nv = FmsEntityNumVerts[et];
341 mfem::Array<int> ents_verts(num_ents*nv), m_ev;
342 const int *ei = (const int *)ents;
343 if (ents == NULL)
344 {
345 FmsDomainGetEntitiesVerts(domain, (FmsEntityType)et, NULL, FMS_INT32,
346 0, ents_verts.GetData(), num_ents);
347 }
348 else
349 {
350 for (FmsInt i = 0; i < num_ents; i++)
351 {
352 FmsDomainGetEntitiesVerts(domain, (FmsEntityType)et, NULL,
353 FMS_INT32, ei[i], &ents_verts[i*nv], 1);
354 }
355 }
356 for (int i = 0; i < ents_verts.Size(); i++)
357 {
358 ents_verts[i] += mfem_last_vert_cnt;
359 }
360 const int *perm;
361 switch ((FmsEntityType)et)
362 {
363 case FMS_EDGE:
364 {
365 for (FmsInt part_ent_id = 0; part_ent_id < num_ents; part_ent_id++)
366 {
367 const int *ev = &ents_verts[2*part_ent_id];
368 int mfem_edge_id = mfem_edge.FindId(ev[0], ev[1]);
369 if (mfem_edge_id < 0)
370 {
371 return 17;
372 }
373 mesh->GetEdgeVertices(mfem_edge_id, m_ev);
374 int ori = (ev[0] == m_ev[0]) ? 0 : 1;
375 perm = fec->DofOrderForOrientation(mfem::Geometry::SEGMENT, ori);
376 fes->GetEdgeInteriorDofs(mfem_edge_id, dofs);
377 for (int i = 0; i < edge_dofs; i++)
378 {
379 for (int j = 0; j < vdim; j++)
380 {
381 func(fes->DofToVDof(dofs[i],j)) =
382 static_cast<double>(src_data[(fms_dof_offset+perm[i])*nstride+j*vstride]);
383 }
384 }
385 fms_dof_offset += edge_dofs;
386 }
387 break;
388 }
389 case FMS_TRIANGLE:
390 {
391 for (FmsInt part_ent_id = 0; part_ent_id < num_ents; part_ent_id++)
392 {
393 const int *tv = &ents_verts[3*part_ent_id];
394 int mfem_face_id = mfem_face.FindId(tv[0], tv[1], tv[2], INT_MAX);
395 if (mfem_face_id < 0)
396 {
397 return 18;
398 }
399 mesh->GetFaceVertices(mfem_face_id, m_ev);
400 int ori = 0;
401 while (tv[ori] != m_ev[0]) { ori++; }
402 ori = (tv[(ori+1)%3] == m_ev[1]) ? 2*ori : 2*ori+1;
403 perm = fec->DofOrderForOrientation(mfem::Geometry::TRIANGLE, ori);
404 fes->GetFaceInteriorDofs(mfem_face_id, dofs);
405 for (int i = 0; i < tri_dofs; i++)
406 {
407 for (int j = 0; j < vdim; j++)
408 {
409 func(fes->DofToVDof(dofs[i],j)) =
410 static_cast<double>(src_data[(fms_dof_offset+perm[i])*nstride+j*vstride]);
411 }
412 }
413 fms_dof_offset += tri_dofs;
414 }
415 break;
416 }
417 case FMS_QUADRILATERAL:
418 {
419 for (FmsInt part_ent_id = 0; part_ent_id < num_ents; part_ent_id++)
420 {
421 const int *qv = &ents_verts[4*part_ent_id];
422 int mfem_face_id = mfem_face.FindId(qv[0], qv[1], qv[2], qv[3]);
423 if (mfem_face_id < 0) { return 19; }
424 mesh->GetFaceVertices(mfem_face_id, m_ev);
425 int ori = 0;
426 while (qv[ori] != m_ev[0]) { ori++; }
427 ori = (qv[(ori+1)%4] == m_ev[1]) ? 2*ori : 2*ori+1;
428 perm = fec->DofOrderForOrientation(mfem::Geometry::SQUARE, ori);
429 fes->GetFaceInteriorDofs(mfem_face_id, dofs);
430 for (int i = 0; i < quad_dofs; i++)
431 {
432 for (int j = 0; j < vdim; j++)
433 {
434 func(fes->DofToVDof(dofs[i],j)) =
435 static_cast<double>(src_data[(fms_dof_offset+perm[i])*nstride+j*vstride]);
436 }
437 }
438 fms_dof_offset += quad_dofs;
439 }
440 break;
441 }
442 default: break;
443 }
444 mfem_ent_cnt[FmsEntityDim[et]] += num_ents;
445 }
446 }
447
448 return err;
449}
450
451int
452FmsMeshToMesh(FmsMesh fms_mesh, Mesh **mfem_mesh)
453{
454 FmsInt dim, n_vert, n_elem, n_bdr_elem, space_dim;
455
456 // Find the first component that has coordinates - that will be the new mfem
457 // mesh.
458 FmsInt num_comp;
459 FmsMeshGetNumComponents(fms_mesh, &num_comp);
460 FmsComponent main_comp = NULL;
461 FmsField coords = NULL;
462 for (FmsInt comp_id = 0; comp_id < num_comp; comp_id++)
463 {
464 FmsComponent comp;
465 FmsMeshGetComponent(fms_mesh, comp_id, &comp);
466 FmsComponentGetCoordinates(comp, &coords);
467 if (coords) { main_comp = comp; break; }
468 }
469 if (!main_comp) { return 1; }
470 FmsComponentGetDimension(main_comp, &dim);
471 FmsComponentGetNumEntities(main_comp, &n_elem);
472 FmsInt n_ents[FMS_NUM_ENTITY_TYPES];
473 FmsInt n_main_parts;
474 FmsComponentGetNumParts(main_comp, &n_main_parts);
475
476#define RENUMBER_ENTITIES
477#ifdef RENUMBER_ENTITIES
478 // I noticed that to get domains working right, since they appear to be
479 // defined in a local vertex numbering scheme, we have to offset the vertex
480 // ids that MFEM makes for shapes to move them to the coordinates in the
481 // current domain.
482
483 // However, parts would just be a set of element ids in the current domain
484 // and it does not seem appropriate to offset the points in that case.
485 // Should domains be treated specially?
486 int *verts_per_part = new int[n_main_parts];
487#endif
488
489 // Sum the counts for each entity type across parts.
490 for (FmsInt et = FMS_VERTEX; et < FMS_NUM_ENTITY_TYPES; et++)
491 {
492 n_ents[et] = 0;
493 for (FmsInt part_id = 0; part_id < n_main_parts; part_id++)
494 {
495 FmsInt num_ents;
496 FmsComponentGetPart(main_comp, part_id, (FmsEntityType)et, NULL, NULL,
497 NULL, NULL, &num_ents);
498 n_ents[et] += num_ents;
499#ifdef RENUMBER_ENTITIES
500 if (et == FMS_VERTEX)
501 {
502 verts_per_part[part_id] = num_ents;
503 }
504#endif
505 }
506 }
507 n_vert = n_ents[FMS_VERTEX];
508
509#ifdef RENUMBER_ENTITIES
510 int *verts_start = new int[n_main_parts];
511 verts_start[0] = 0;
512 for (FmsInt i = 1; i < n_main_parts; ++i)
513 {
514 verts_start[i] = verts_start[i-1] + verts_per_part[i-1];
515 }
516#endif
517
518 // The first related component of dimension dim-1 will be the boundary of the
519 // new mfem mesh.
520 FmsComponent bdr_comp = NULL;
521 FmsInt num_rel_comps;
522 const FmsInt *rel_comp_ids;
523 FmsComponentGetRelations(main_comp, &rel_comp_ids, &num_rel_comps);
524 for (FmsInt i = 0; i < num_rel_comps; i++)
525 {
526 FmsComponent comp;
527 FmsMeshGetComponent(fms_mesh, rel_comp_ids[i], &comp);
528 FmsInt comp_dim;
529 FmsComponentGetDimension(comp, &comp_dim);
530 if (comp_dim == dim-1) { bdr_comp = comp; break; }
531 }
532 if (bdr_comp)
533 {
534 FmsComponentGetNumEntities(bdr_comp, &n_bdr_elem);
535 }
536 else
537 {
538 n_bdr_elem = 0;
539 }
540
541 FmsFieldGet(coords, NULL, &space_dim, NULL, NULL, NULL);
542 int err = 0;
543 Mesh *mesh = nullptr;
544 mesh = new Mesh(dim, n_vert, n_elem, n_bdr_elem, space_dim);
545
546 // Element tags
547 FmsInt num_tags;
548 FmsMeshGetNumTags(fms_mesh, &num_tags);
549 FmsTag elem_tag = NULL, bdr_tag = NULL;
550 for (FmsInt tag_id = 0; tag_id < num_tags; tag_id++)
551 {
552 FmsTag tag;
553 FmsMeshGetTag(fms_mesh, tag_id, &tag);
554 FmsComponent comp;
555 FmsTagGetComponent(tag, &comp);
556 if (!elem_tag && comp == main_comp)
557 {
558#if DEBUG_FMS_MFEM
559 const char *tn = NULL;
560 FmsTagGetName(tag, &tn);
561 mfem::out << "Found element tag " << tn << std::endl;
562#endif
563 elem_tag = tag;
564 }
565 else if (!bdr_tag && comp == bdr_comp)
566 {
567#if DEBUG_FMS_MFEM
568 const char *tn = NULL;
569 FmsTagGetName(tag, &tn);
570 mfem::out << "Found boundary tag " << tn << std::endl;
571#endif
572 bdr_tag = tag;
573 }
574 }
575 FmsIntType attr_type;
576 const void *v_attr, *v_bdr_attr;
577 mfem::Array<int> attr, bdr_attr;
578 FmsInt num_attr;
579 // Element attributes
580 if (elem_tag)
581 {
582 FmsTagGet(elem_tag, &attr_type, &v_attr, &num_attr);
583 if (attr_type == FMS_UINT8)
584 {
585 mfem::Array<uint8_t> at((uint8_t*)v_attr, num_attr);
586 attr = at;
587 }
588 else if (attr_type == FMS_INT32 || attr_type == FMS_UINT32)
589 {
590 attr.MakeRef((int*)v_attr, num_attr);
591 }
592 else
593 {
594 err = 1; // "attribute type not supported!"
595 goto func_exit;
596 }
597 }
598 // Boundary attributes
599 if (bdr_tag)
600 {
601 FmsTagGet(bdr_tag, &attr_type, &v_bdr_attr, &num_attr);
602 if (attr_type == FMS_UINT8)
603 {
604 mfem::Array<uint8_t> at((uint8_t*)v_bdr_attr, num_attr);
605 bdr_attr = at;
606 }
607 else if (attr_type == FMS_INT32 || attr_type == FMS_UINT32)
608 {
609 bdr_attr.MakeRef((int*)v_bdr_attr, num_attr);
610 }
611 else
612 {
613 err = 2; // "bdr attribute type not supported!"
614 goto func_exit;
615 }
616 }
617
618 // Add elements
619 for (FmsInt part_id = 0; part_id < n_main_parts; part_id++)
620 {
621 for (int et = FMS_VERTEX; et < FMS_NUM_ENTITY_TYPES; et++)
622 {
623 if (FmsEntityDim[et] != dim) { continue; }
624
625 FmsDomain domain;
626 FmsIntType elem_id_type;
627 const void *elem_ids;
628 const FmsOrientation *elem_ori;
629 FmsInt num_elems;
630 FmsComponentGetPart(main_comp, part_id, (FmsEntityType)et, &domain,
631 &elem_id_type, &elem_ids, &elem_ori, &num_elems);
632
633 if (num_elems == 0) { continue; }
634
635 if (elem_ids != NULL &&
636 (elem_id_type != FMS_INT32 && elem_id_type != FMS_UINT32))
637 {
638 err = 3; goto func_exit;
639 }
640 if (elem_ori != NULL)
641 {
642 err = 4; goto func_exit;
643 }
644
645 const FmsInt nv = FmsEntityNumVerts[et];
646 mfem::Array<int> ents_verts(num_elems*nv);
647 if (elem_ids == NULL)
648 {
649 FmsDomainGetEntitiesVerts(domain, (FmsEntityType)et, NULL, FMS_INT32,
650 0, ents_verts.GetData(), num_elems);
651 }
652 else
653 {
654 const int *ei = (const int *)elem_ids;
655 for (FmsInt i = 0; i < num_elems; i++)
656 {
657 FmsDomainGetEntitiesVerts(domain, (FmsEntityType)et, NULL, FMS_INT32,
658 ei[i], &ents_verts[i*nv], 1);
659 }
660 }
661 const int elem_offset = mesh->GetNE();
662 switch ((FmsEntityType)et)
663 {
664 case FMS_EDGE:
665 err = 5;
666 goto func_exit;
667 break;
668 case FMS_TRIANGLE:
669#ifdef RENUMBER_ENTITIES
670 // The domain vertices/edges were defined in local ordering. We
671 // now have a set of triangle vertices defined in terms of local
672 // vertex numbers. Renumber them to a global numbering.
673 for (FmsInt i = 0; i < num_elems*3; i++)
674 {
675 ents_verts[i] += verts_start[part_id];
676 }
677#endif
678
679 for (FmsInt i = 0; i < num_elems; i++)
680 {
681 mesh->AddTriangle(
682 &ents_verts[3*i], elem_tag ? attr[elem_offset+i] : 1);
683 }
684 break;
685 case FMS_QUADRILATERAL:
686#ifdef RENUMBER_ENTITIES
687 for (FmsInt i = 0; i < num_elems*4; i++)
688 {
689 ents_verts[i] += verts_start[part_id];
690 }
691#endif
692 for (FmsInt i = 0; i < num_elems; i++)
693 {
694 mesh->AddQuad(&ents_verts[4*i], elem_tag ? attr[elem_offset+i] : 1);
695 }
696 break;
697 case FMS_TETRAHEDRON:
698#ifdef RENUMBER_ENTITIES
699 for (FmsInt i = 0; i < num_elems*4; i++)
700 {
701 ents_verts[i] += verts_start[part_id];
702 }
703#endif
704 for (FmsInt i = 0; i < num_elems; i++)
705 {
706 mesh->AddTet(&ents_verts[4*i], elem_tag ? attr[elem_offset+i] : 1);
707 }
708 break;
709
710 // TODO: What about wedges and pyramids?
711
712
713 case FMS_HEXAHEDRON:
714#ifdef RENUMBER_ENTITIES
715 for (FmsInt i = 0; i < num_elems*8; i++)
716 {
717 ents_verts[i] += verts_start[part_id];
718 }
719
720#endif
721 for (FmsInt i = 0; i < num_elems; i++)
722 {
723 const int *hex_verts = &ents_verts[8*i];
724#if 0
725 const int reorder[8] = {0, 1, 2, 3, 5, 4, 6, 7};
726 const int new_verts[8] = {hex_verts[reorder[0]],
727 hex_verts[reorder[1]],
728 hex_verts[reorder[2]],
729 hex_verts[reorder[3]],
730 hex_verts[reorder[4]],
731 hex_verts[reorder[5]],
732 hex_verts[reorder[6]],
733 hex_verts[reorder[7]]
734 };
735 hex_verts = new_verts;
736#endif
737 mesh->AddHex(hex_verts, elem_tag ? attr[elem_offset+i] : 1);
738 }
739 break;
740 default:
741 break;
742 }
743 }
744 }
745
746 // Add boundary elements
747 if (bdr_comp && n_bdr_elem > 0)
748 {
749 FmsInt n_bdr_parts;
750 FmsComponentGetNumParts(bdr_comp, &n_bdr_parts);
751
752 for (FmsInt part_id = 0; part_id < n_bdr_parts; part_id++)
753 {
754 for (int et = FMS_VERTEX; et < FMS_NUM_ENTITY_TYPES; et++)
755 {
756 if (FmsEntityDim[et] != dim-1) { continue; }
757
758 FmsDomain domain;
759 FmsIntType elem_id_type;
760 const void *elem_ids;
761 const FmsOrientation *elem_ori;
762 FmsInt num_elems;
763 FmsComponentGetPart(bdr_comp, part_id, (FmsEntityType)et, &domain,
764 &elem_id_type, &elem_ids, &elem_ori, &num_elems);
765 if (num_elems == 0) { continue; }
766
767 if (elem_ids != NULL &&
768 (elem_id_type != FMS_INT32 && elem_id_type != FMS_UINT32))
769 {
770 err = 6; goto func_exit;
771 }
772 if (elem_ori != NULL)
773 {
774 err = 7; goto func_exit;
775 }
776
777 const FmsInt nv = FmsEntityNumVerts[et];
778 mfem::Array<int> ents_verts(num_elems*nv);
779 if (elem_ids == NULL)
780 {
781 FmsDomainGetEntitiesVerts(domain, (FmsEntityType)et, NULL, FMS_INT32,
782 0, ents_verts.GetData(), num_elems);
783 }
784 else
785 {
786 const int *ei = (const int *)elem_ids;
787 for (FmsInt i = 0; i < num_elems; i++)
788 {
789 FmsDomainGetEntitiesVerts(domain, (FmsEntityType)et, NULL,
790 FMS_INT32, ei[i], &ents_verts[i*nv], 1);
791 }
792 }
793 const int elem_offset = mesh->GetNBE();
794 switch ((FmsEntityType)et)
795 {
796 case FMS_EDGE:
797 for (FmsInt i = 0; i < num_elems; i++)
798 {
799 mesh->AddBdrSegment(
800 &ents_verts[2*i], bdr_tag ? bdr_attr[elem_offset+i] : 1);
801 }
802 break;
803 case FMS_TRIANGLE:
804 for (FmsInt i = 0; i < num_elems; i++)
805 {
806 mesh->AddBdrTriangle(
807 &ents_verts[3*i], bdr_tag ? bdr_attr[elem_offset+i] : 1);
808 }
809 break;
810 case FMS_QUADRILATERAL:
811 for (FmsInt i = 0; i < num_elems; i++)
812 {
813 mesh->AddBdrQuad(
814 &ents_verts[4*i], bdr_tag ? bdr_attr[elem_offset+i] : 1);
815 }
816 break;
817 default:
818 break;
819 }
820 }
821 }
822 }
823
824#ifdef RENUMBER_ENTITIES
825 delete [] verts_per_part;
826 delete [] verts_start;
827#endif
828
829 // Transfer coordinates
830 {
831 // Set the vertex coordinates to zero
832 const double origin[3] = {0.,0.,0.};
833 for (FmsInt vi = 0; vi < n_vert; vi++)
834 {
835 mesh->AddVertex(origin);
836 }
837
838 // Finalize the mesh topology
839 mesh->FinalizeTopology();
840
841 FmsFieldDescriptor coords_fd = NULL;
842 FmsLayoutType coords_layout;
843 FmsFieldGet(coords, &coords_fd, NULL, &coords_layout, NULL, NULL);
844 if (!coords_fd)
845 {
846 mfem::err << "Error reading the FMS mesh coords' FieldDescriptor." << std::endl;
847 err = 8;
848 goto func_exit;
849 }
850 FmsInt coords_order = 0;
851 FmsBasisType coords_btype = FMS_NODAL_GAUSS_CLOSED;
852 FmsFieldType coords_ftype = FMS_CONTINUOUS;
853 FmsFieldDescriptorGetFixedOrder(coords_fd, &coords_ftype, &coords_btype,
854 &coords_order);
855 // Maybe this is extra but it seems mesh->SetCurvature assumes
856 // btype=1. Maybe protects us against corrupt data.
857 if (coords_btype != FMS_NODAL_GAUSS_CLOSED)
858 {
859 mfem::err << "Error reading FMS mesh coords." << std::endl;
860 err = 9;
861 goto func_exit;
862 }
863
864 // Switch to mfem::Mesh with nodes (interpolates the linear coordinates)
865 const bool discont = (coords_ftype == FMS_DISCONTINUOUS);
866 mesh->SetCurvature(coords_order, discont, space_dim,
867 (coords_layout == FMS_BY_VDIM) ?
868 mfem::Ordering::byVDIM : mfem::Ordering::byNODES);
869
870 // Finalize mesh construction
871 mesh->Finalize();
872
873 // Set the high-order mesh nodes
874 mfem::GridFunction &nodes = *mesh->GetNodes();
875 FmsFieldToGridFunction<double>(fms_mesh, coords, mesh, nodes, false);
876 }
877
878func_exit:
879
880 if (err)
881 {
882 delete mesh;
883 }
884 else
885 {
886 *mfem_mesh = mesh;
887 }
888 return err;
889}
890
891bool
892BasisTypeToFmsBasisType(int bt, FmsBasisType &btype)
893{
894 bool retval = false;
895 switch (bt)
896 {
897 case mfem::BasisType::GaussLegendre:
898 // mfem::out << "mfem::BasisType::GaussLegendre -> FMS_NODAL_GAUSS_OPEN" << endl;
899 btype = FMS_NODAL_GAUSS_OPEN;
900 retval = true;
901 break;
902 case mfem::BasisType::GaussLobatto:
903 // mfem::out << "mfem::BasisType::GaussLobato -> FMS_NODAL_GAUSS_CLOSED" << endl;
904 btype = FMS_NODAL_GAUSS_CLOSED;
905 retval = true;
906 break;
907 case mfem::BasisType::Positive:
908 // mfem::out << "mfem::BasisType::Positive -> FMS_POSITIVE" << endl;
909 btype = FMS_POSITIVE;
910 retval = true;
911 break;
912 case mfem::BasisType::OpenUniform:
913 // mfem::out << "mfem::BasisType::OpenUniform -> ?" << endl;
914 btype = FMS_NODAL_UNIFORM_OPEN;
915 retval = true;
916 break;
917 case mfem::BasisType::ClosedUniform:
918 // mfem::out << "mfem::BasisType::ClosedUniform -> ?" << endl;
919 btype = FMS_NODAL_UNIFORM_CLOSED;
920 retval = true;
921 break;
922 case mfem::BasisType::OpenHalfUniform:
923 // mfem::out << "mfem::BasisType::OpenHalfUniform -> ?" << endl;
924 break;
925 case mfem::BasisType::Serendipity:
926 // mfem::out << "mfem::BasisType::Serendipity -> ?" << endl;
927 break;
928 case mfem::BasisType::ClosedGL:
929 // mfem::out << "mfem::BasisType::ClosedGL -> ?" << endl;
930 break;
931
932 }
933 /*
934 Which MFEM types map to:?
935 FMS_NODAL_CHEBYSHEV_OPEN,
936 FMS_NODAL_CHEBYSHEV_CLOSED,
937 */
938
939 return retval;
940}
941
942/**
943@note We add the FMS field descriptor and field in here so we can only do it
944 after successfully validating the inputs (handling certain grid function
945 types, etc.)
946*/
947int
948GridFunctionToFmsField(FmsDataCollection dc,
949 FmsComponent comp,
950 const std::string &fd_name,
951 const std::string &field_name,
952 const Mesh *mesh,
953 const GridFunction *gf,
954 FmsField *outfield)
955{
956 if (!dc) { return 1; }
957 if (!comp) { return 2; }
958 if (!mesh) { return 3; }
959 if (!gf) { return 4; }
960 if (!outfield) { return 5; }
961
962 double *c = gf->GetData();
963
964 const mfem::FiniteElementSpace *fespace = gf->FESpace();
965 const mfem::FiniteElementCollection *fecoll = fespace->FEColl();
966
967#ifdef DEBUG_MFEM_FMS
968 mfem::out << "Adding FMS field for " << field_name << "..." << endl;
969#endif
970
971 /* Q: No getter for the basis, do different kinds of FECollection have
972 implied basis? There are two subclasses that actually have the getter,
973 maybe those aren't implied? */
974 FmsInt order = 1;
975 int vdim = 1;
976 FmsFieldType ftype = FMS_CONTINUOUS;
977 FmsBasisType btype = FMS_NODAL_GAUSS_CLOSED;
978 switch (fecoll->GetContType())
979 {
980 case mfem::FiniteElementCollection::CONTINUOUS:
981 {
982 ftype = FMS_CONTINUOUS;
983 order = static_cast<FmsInt>(fespace->GetOrder(0));
984 vdim = gf->VectorDim();
985 auto fec = dynamic_cast<const mfem::H1_FECollection *>(fecoll);
986 if (fec != nullptr)
987 {
988 if (!BasisTypeToFmsBasisType(fec->GetBasisType(), btype))
989 {
990 mfem::err << "Error converting MFEM basis type to FMS for"
991 " FMS_CONTINUOUS." << std::endl;
992 return 6;
993 }
994 }
995 break;
996 }
997 case mfem::FiniteElementCollection::DISCONTINUOUS:
998 {
999 ftype = FMS_DISCONTINUOUS;
1000 order = static_cast<FmsInt>(fespace->GetOrder(0));
1001 vdim = gf->VectorDim();
1002 auto fec = dynamic_cast<const mfem::L2_FECollection *>(fecoll);
1003 if (fec != nullptr)
1004 {
1005 if (!BasisTypeToFmsBasisType(fec->GetBasisType(), btype))
1006 {
1007 mfem::err << "Error converting MFEM basis type to FMS for"
1008 " FMS_DISCONTINUOUS." << std::endl;
1009 return 7;
1010 }
1011 }
1012 break;
1013 }
1014 case mfem::FiniteElementCollection::TANGENTIAL:
1015 {
1016 mfem::out << "Warning, unsupported ContType (TANGENTIAL) for "
1017 << field_name << ". Using FMS_CONTINUOUS." << std::endl;
1018 break;
1019 }
1020 case mfem::FiniteElementCollection::NORMAL:
1021 {
1022 ftype = FMS_HDIV;
1023 // This RT_FECollection type seems to arise from "RT" fields such as "RT_3D_P1".
1024 // Checking fe_coll.hpp, this contains verbiage about H_DIV so we assign it the
1025 // FMS type FMS_HDIV.
1026
1027 // I've seen RT_3D_P1 return the wrong order so get it from the name.
1028 int idim, iorder;
1029 if (sscanf(fecoll->Name(), "RT_%dD_P%d", &idim, &iorder) == 2)
1030 {
1031 order = (FmsInt)iorder;
1032 }
1033 else
1034 {
1035 order = static_cast<FmsInt>(fespace->GetOrder(0));
1036 }
1037
1038 // Get the vdim from the fespace since the grid function is returning
1039 // 3 but we need it to be what was read from the file so we can pass the
1040 // right vdim to the FMS field descriptor to compute the expected number of dofs.
1041 vdim = fespace->GetVDim();
1042 break;
1043 }
1044 default:
1045 mfem::out << "Warning, unsupported ContType for field " << field_name
1046 << ". Using FMS_CONTINUOUS." << std::endl;
1047 ftype = FMS_CONTINUOUS;
1048 break;
1049 }
1050
1051 // Now that we're not failing, create the fd and field.
1052 FmsFieldDescriptor fd = NULL;
1053 FmsField f = NULL;
1054 FmsDataCollectionAddFieldDescriptor(dc, fd_name.c_str(), &fd);
1055 FmsDataCollectionAddField(dc, field_name.c_str(), &f);
1056 *outfield = f;
1057
1058 /* Q: Why is order defined on a per element basis? */
1059 FmsFieldDescriptorSetComponent(fd, comp);
1060 FmsFieldDescriptorSetFixedOrder(fd, ftype, btype, order);
1061
1062 FmsInt ndofs;
1063 FmsFieldDescriptorGetNumDofs(fd, &ndofs);
1064
1065 const char *name = NULL;
1066 FmsFieldGetName(f, &name);
1067 FmsLayoutType layout = fespace->GetOrdering() == mfem::Ordering::byVDIM ?
1068 FMS_BY_VDIM : FMS_BY_NODES;
1069
1070#ifdef DEBUG_MFEM_FMS
1071 switch (ftype)
1072 {
1073 case FMS_CONTINUOUS:
1074 mfem::out << "\tFMS_CONTINUOUS" << std::endl;
1075 break;
1076 case FMS_DISCONTINUOUS:
1077 mfem::out << "\tFMS_DISCONTINUOUS" << std::endl;
1078 break;
1079 case FMS_HDIV:
1080 mfem::out << "\tFMS_HDIV" << std::endl;
1081 break;
1082 }
1083 mfem::out << "\tField is order " << order << " with vdim " << vdim <<
1084 " and nDoFs " << ndofs << std::endl;
1085 mfem::out << "\tgf->size() " << gf->Size() << " ndofs * vdim " << ndofs * vdim
1086 << std::endl;
1087 mfem::out << "\tlayout " << layout << " (0 = BY_NODES, 1 = BY_VDIM)" <<
1088 std::endl;
1089#endif
1090
1091 if (FmsFieldSet(f, fd, vdim, layout, FMS_DOUBLE, c))
1092 {
1093 mfem::err << "Error setting field " << field_name << " in FMS." << std::endl;
1094 return 8;
1095 }
1096 return 0;
1097}
1098
1099bool
1100MfemMetaDataToFmsMetaData(DataCollection *mdc, FmsDataCollection fdc)
1101{
1102 if (!mdc) { return false; }
1103 if (!fdc) { return false; }
1104
1105 int *cycle = NULL;
1106 double *time = NULL, *timestep = NULL;
1107 FmsMetaData top_level = NULL;
1108 FmsMetaData *cycle_time_timestep = NULL;
1109 int mdata_err = 0;
1110 mdata_err = FmsDataCollectionAttachMetaData(fdc, &top_level);
1111 if (!top_level || mdata_err)
1112 {
1113 mfem::err << "Failed to attach metadata to the FmsDataCollection" << std::endl;
1114 return 40;
1115 }
1116
1117 mdata_err = FmsMetaDataSetMetaData(top_level, "MetaData", 3,
1118 &cycle_time_timestep);
1119 if (!cycle_time_timestep || mdata_err)
1120 {
1121 mfem::err << "Failed to acquire FmsMetaData array" << std::endl;
1122 return false;
1123 }
1124
1125 if (!cycle_time_timestep[0])
1126 {
1127 mfem::err << "The MetaData pointer for cycle is NULL" << std::endl;
1128 return false;
1129 }
1130 mdata_err = FmsMetaDataSetIntegers(cycle_time_timestep[0], "cycle",
1131 FMS_INT32, 1, (void**)&cycle);
1132 if (!cycle || mdata_err)
1133 {
1134 mfem::err << "The data pointer for cycle is NULL" << std::endl;
1135 return false;
1136 }
1137 *cycle = mdc->GetCycle();
1138
1139 if (!cycle_time_timestep[1])
1140 {
1141 mfem::err << "The FmsMetaData pointer for time is NULL" << std::endl;
1142 return false;
1143 }
1144 mdata_err = FmsMetaDataSetScalars(cycle_time_timestep[1], "time", FMS_DOUBLE,
1145 1, (void**)&time);
1146 if (!time || mdata_err)
1147 {
1148 mfem::err << "The data pointer for time is NULL." << std::endl;
1149 return false;
1150 }
1151 *time = mdc->GetTime();
1152
1153 if (!cycle_time_timestep[2])
1154 {
1155 mfem::err << "The FmsMetData pointer for timestep is NULL" << std::endl;
1156 return false;
1157 }
1158 mdata_err = FmsMetaDataSetScalars(cycle_time_timestep[2], "timestep",
1159 FMS_DOUBLE, 1, (void**)&timestep);
1160 if (!timestep || mdata_err)
1161 {
1162 mfem::err << "The data pointer for timestep is NULL" << std::endl;
1163 return false;
1164 }
1165 *timestep = mdc->GetTimeStep();
1166
1167 return true;
1168}
1169
1170//---------------------------------------------------------------------------
1171bool
1172FmsMetaDataGetInteger(FmsMetaData mdata, const std::string &key,
1173 std::vector<int> &values)
1174{
1175 if (!mdata) { return false; }
1176
1177 bool retval = false;
1178 FmsMetaDataType type;
1179 FmsIntType int_type;
1180 FmsInt i, size;
1181 FmsMetaData *children = nullptr;
1182 const void *data = nullptr;
1183 const char *mdata_name = nullptr;
1184 if (FmsMetaDataGetType(mdata, &type) == 0)
1185 {
1186 switch (type)
1187 {
1188 case FMS_INTEGER:
1189 if (FmsMetaDataGetIntegers(mdata, &mdata_name, &int_type, &size, &data) == 0)
1190 {
1191 if (strcasecmp(key.c_str(), mdata_name) == 0)
1192 {
1193 retval = true;
1194
1195 // Interpret the integers and store them in the std::vector<int>
1196 switch (int_type)
1197 {
1198 case FMS_INT8:
1199 for (i = 0; i < size; i++)
1200 {
1201 values.push_back(static_cast<int>(reinterpret_cast<const int8_t*>(data)[i]));
1202 }
1203 break;
1204 case FMS_INT16:
1205 for (i = 0; i < size; i++)
1206 {
1207 values.push_back(static_cast<int>(reinterpret_cast<const int16_t*>(data)[i]));
1208 }
1209 break;
1210 case FMS_INT32:
1211 for (i = 0; i < size; i++)
1212 {
1213 values.push_back(static_cast<int>(reinterpret_cast<const int32_t*>(data)[i]));
1214 }
1215 break;
1216 case FMS_INT64:
1217 for (i = 0; i < size; i++)
1218 {
1219 values.push_back(static_cast<int>(reinterpret_cast<const int64_t*>(data)[i]));
1220 }
1221 break;
1222 case FMS_UINT8:
1223 for (i = 0; i < size; i++)
1224 {
1225 values.push_back(static_cast<int>(reinterpret_cast<const uint8_t*>(data)[i]));
1226 }
1227 break;
1228 case FMS_UINT16:
1229 for (i = 0; i < size; i++)
1230 {
1231 values.push_back(static_cast<int>(reinterpret_cast<const uint16_t*>(data)[i]));
1232 }
1233 break;
1234 case FMS_UINT32:
1235 for (i = 0; i < size; i++)
1236 {
1237 values.push_back(static_cast<int>(reinterpret_cast<const uint32_t*>(data)[i]));
1238 }
1239 break;
1240 case FMS_UINT64:
1241 for (i = 0; i < size; i++)
1242 {
1243 values.push_back(static_cast<int>(reinterpret_cast<const uint64_t*>(data)[i]));
1244 }
1245 break;
1246 default:
1247 retval = false;
1248 break;
1249 }
1250 }
1251 }
1252 break;
1253 case FMS_META_DATA:
1254 if (FmsMetaDataGetMetaData(mdata, &mdata_name, &size, &children) == 0)
1255 {
1256 // Recurse to look for the key we want.
1257 for (i = 0; i < size && !retval; i++)
1258 {
1259 retval = FmsMetaDataGetInteger(children[i], key, values);
1260 }
1261 }
1262 break;
1263 default:
1264 break;
1265 }
1266 }
1267
1268 return retval;
1269}
1270
1271//---------------------------------------------------------------------------
1272bool
1273FmsMetaDataGetScalar(FmsMetaData mdata, const std::string &key,
1274 std::vector<double> &values)
1275{
1276 if (!mdata) { return false; }
1277
1278 bool retval = false;
1279 FmsMetaDataType type;
1280 FmsScalarType scal_type;
1281 FmsInt i, size;
1282 FmsMetaData *children = nullptr;
1283 const void *data = nullptr;
1284 const char *mdata_name = nullptr;
1285 if (FmsMetaDataGetType(mdata, &type) == 0)
1286 {
1287 switch (type)
1288 {
1289 case FMS_SCALAR:
1290 if (FmsMetaDataGetScalars(mdata, &mdata_name, &scal_type, &size, &data) == 0)
1291 {
1292 if (strcasecmp(key.c_str(), mdata_name) == 0)
1293 {
1294 retval = true;
1295
1296 // Interpret the integers and store them in the std::vector<int>
1297 switch (scal_type)
1298 {
1299 case FMS_FLOAT:
1300 for (i = 0; i < size; i++)
1301 {
1302 values.push_back(static_cast<double>(reinterpret_cast<const float*>(data)[i]));
1303 }
1304 break;
1305 case FMS_DOUBLE:
1306 for (i = 0; i < size; i++)
1307 {
1308 values.push_back(reinterpret_cast<const double*>(data)[i]);
1309 }
1310 break;
1311 default:
1312 retval = false;
1313 break;
1314 }
1315 }
1316 }
1317 break;
1318 case FMS_META_DATA:
1319 if (FmsMetaDataGetMetaData(mdata, &mdata_name, &size, &children) == 0)
1320 {
1321 // Recurse to look for the key we want.
1322 for (i = 0; i < size && !retval; i++)
1323 {
1324 retval = FmsMetaDataGetScalar(children[i], key, values);
1325 }
1326 }
1327 break;
1328 default:
1329 break;
1330 }
1331 }
1332
1333 return retval;
1334}
1335
1336//---------------------------------------------------------------------------
1337bool
1338FmsMetaDataGetString(FmsMetaData mdata, const std::string &key,
1339 std::string &value)
1340{
1341 if (!mdata) { return false; }
1342
1343 bool retval = false;
1344 FmsMetaDataType type;
1345 FmsInt i, size;
1346 FmsMetaData *children = nullptr;
1347 const char *mdata_name = nullptr;
1348 const char *str_value = nullptr;
1349
1350 if (FmsMetaDataGetType(mdata, &type) == 0)
1351 {
1352 switch (type)
1353 {
1354 case FMS_STRING:
1355 if (FmsMetaDataGetString(mdata, &mdata_name, &str_value) == 0)
1356 {
1357 if (strcasecmp(key.c_str(), mdata_name) == 0)
1358 {
1359 retval = true;
1360 value = str_value;
1361 }
1362 }
1363 break;
1364 case FMS_META_DATA:
1365 if (FmsMetaDataGetMetaData(mdata, &mdata_name, &size, &children) == 0)
1366 {
1367 // Recurse to look for the key we want.
1368 for (i = 0; i < size && !retval; i++)
1369 {
1370 retval = FmsMetaDataGetString(children[i], key, value);
1371 }
1372 }
1373 break;
1374 default:
1375 break;
1376 }
1377 }
1378
1379 return retval;
1380}
1381
1382/* -------------------------------------------------------------------------- */
1383/* FMS to MFEM conversion function */
1384/* -------------------------------------------------------------------------- */
1385
1386int FmsDataCollectionToDataCollection(FmsDataCollection dc,
1387 DataCollection **mfem_dc)
1388{
1389 int retval = 0;
1390 FmsMesh fms_mesh;
1391 FmsDataCollectionGetMesh(dc, &fms_mesh);
1392
1393 // NOTE: The MFEM data collection has a single Mesh. Mesh has a constructor
1394 // to take multiple Mesh objects but it appears to glue them together.
1395 Mesh *mesh = nullptr;
1396 int err = FmsMeshToMesh(fms_mesh, &mesh);
1397 if (err == 0)
1398 {
1399 std::string collection_name("collection");
1400 const char *cn = nullptr;
1401 FmsDataCollectionGetName(dc, &cn);
1402 if (cn != nullptr)
1403 {
1404 collection_name = cn;
1405 }
1406
1407 // Make a data collection that contains the mesh.
1408 DataCollection *mdc = new DataCollection(collection_name, mesh);
1409 mdc->SetOwnData(true);
1410
1411 // Now do fields, etc. and add them to mdc.
1412 FmsField *fields = nullptr;
1413 FmsInt num_fields = 0;
1414 if (FmsDataCollectionGetFields(dc, &fields, &num_fields) == 0)
1415 {
1416 for (FmsInt i = 0; i < num_fields; ++i)
1417 {
1418 const char *name = nullptr;
1419 FmsFieldGetName(fields[i], &name);
1420
1421 GridFunction *gf = new GridFunction;
1422
1423 // Get the data type.
1424 FmsFieldDescriptor f_fd;
1425 FmsLayoutType f_layout;
1426 FmsScalarType f_data_type;
1427 const void *f_data;
1428 FmsFieldGet(fields[i], &f_fd, NULL, &f_layout, &f_data_type,
1429 &f_data);
1430
1431 // Interpret the field according to its data type.
1432 int err = 1;
1433 switch (f_data_type)
1434 {
1435 case FMS_FLOAT:
1436 err = FmsFieldToGridFunction<float>(fms_mesh, fields[i], mesh, *gf, true);
1437 break;
1438 case FMS_DOUBLE:
1439 err = FmsFieldToGridFunction<double>(fms_mesh, fields[i], mesh, *gf, true);
1440 break;
1441 case FMS_COMPLEX_FLOAT:
1442 case FMS_COMPLEX_DOUBLE:
1443 // Does MFEM support complex?
1444 break;
1445 default:
1446 break;
1447 }
1448
1449 if (err == 0)
1450 {
1451 mdc->RegisterField(name, gf);
1452 }
1453 else
1454 {
1455 mfem::out << "There was an error converting " << name << " code: " << err <<
1456 std::endl;
1457 delete gf;
1458 }
1459
1460 const char *fname = NULL;
1461 FmsFieldGetName(fields[i], &fname);
1462 }
1463 }
1464
1465 // If we have metadata in FMS, pass what we can through to MFEM.
1466 FmsMetaData mdata = NULL;
1467 FmsDataCollectionGetMetaData(dc, &mdata);
1468 if (mdata)
1469 {
1470 std::vector<int> ivalues;
1471 std::vector<double> dvalues;
1472 std::string svalue;
1473 if (FmsMetaDataGetInteger(mdata, "cycle", ivalues))
1474 {
1475 if (!ivalues.empty())
1476 {
1477 mdc->SetCycle(ivalues[0]);
1478 }
1479 }
1480 if (FmsMetaDataGetScalar(mdata, "time", dvalues))
1481 {
1482 if (!dvalues.empty())
1483 {
1484 mdc->SetTime(dvalues[0]);
1485 }
1486 }
1487 if (FmsMetaDataGetScalar(mdata, "timestep", dvalues))
1488 {
1489 if (!dvalues.empty())
1490 {
1491 mdc->SetTimeStep(dvalues[0]);
1492 }
1493 }
1494 }
1495
1496 *mfem_dc = mdc;
1497 }
1498 else
1499 {
1500 mfem::out << "FmsDataCollectionToDataCollection: mesh failed to convert. err="
1501 << err
1502 << endl;
1503
1504 retval = 1;
1505 }
1506
1507#if DEBUG_FMS_MFEM
1508 if (*mfem_dc)
1509 {
1510 VisItDataCollection visit_dc(std::string("DEBUG_DC"), mesh);
1511 visit_dc.SetOwnData(false);
1512 const auto &fields = (*mfem_dc)->GetFieldMap();
1513 for (const auto &field : fields)
1514 {
1515 visit_dc.RegisterField(field.first, field.second);
1516 }
1517 visit_dc.Save();
1518 }
1519#endif
1520
1521 return retval;
1522}
1523
1524
1525
1526/* -------------------------------------------------------------------------- */
1527/* MFEM to FMS conversion function */
1528/* -------------------------------------------------------------------------- */
1529
1530int
1531MeshToFmsMesh(const Mesh *mmesh, FmsMesh *fmesh, FmsComponent *volume)
1532{
1533 if (!mmesh) { return 1; }
1534 if (!fmesh) { return 2; }
1535 if (!volume) { return 3; }
1536
1537
1538 int err = 0;
1539 const int num_verticies = mmesh->GetNV();
1540#ifdef DEBUG_MFEM_FMS
1541 const int num_edges = mmesh->GetNEdges();
1542#endif
1543 const int num_faces = mmesh->GetNFaces();
1544 const int num_elements = mmesh->GetNE();
1545
1546#ifdef DEBUG_MFEM_FMS
1547 mfem::out << "nverts: " << num_verticies << std::endl;
1548 mfem::out << "nedges: " << num_edges << std::endl;
1549 mfem::out << "nfaces: " << num_faces << std::endl;
1550 mfem::out << "nele: " << num_elements << std::endl;
1551#endif
1552
1553 FmsMeshConstruct(fmesh);
1554 FmsMeshSetPartitionId(*fmesh, 0, 1);
1555
1556 FmsDomain *domains = NULL;
1557 FmsMeshAddDomains(*fmesh, "Domain", 1, &domains);
1558 FmsDomainSetNumVertices(domains[0], num_verticies);
1559
1560 const int edge_reorder[2] = {1, 0};
1561 const int quad_reorder[4] = {0,1,2,3};
1562 const int tet_reorder[4] = {3,2,1,0};
1563 const int hex_reorder[6] = {0,5,1,3,4,2};
1564 const int *reorder[8] = {NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL};
1565 reorder[FMS_EDGE] = edge_reorder;
1566 reorder[FMS_QUADRILATERAL] = quad_reorder;
1567 reorder[FMS_TETRAHEDRON] = tet_reorder;
1568 reorder[FMS_HEXAHEDRON] = hex_reorder;
1569
1570 const mfem::Table *edges = mmesh->GetEdgeVertexTable();
1571 if (!edges)
1572 {
1573 mfem::err << "Error, mesh has no edges." << std::endl;
1574 return 1;
1575 }
1576 mfem::Table *faces = mmesh->GetFaceEdgeTable();
1577 if (!faces && num_faces > 0)
1578 {
1579 mfem::err <<
1580 "Error, mesh contains faces but the \"GetFaceEdgeTable\" returned NULL." <<
1581 std::endl;
1582 return 2;
1583 }
1584
1585 // Build edges
1586 std::vector<int> edge_verts(edges->Size() * 2);
1587 for (int i = 0; i < edges->Size(); i++)
1588 {
1589 mfem::Array<int> nids;
1590 edges->GetRow(i, nids);
1591 for (int j = 0; j < 2; j++)
1592 {
1593 edge_verts[i*2 + j] = nids[j];
1594 }
1595 }
1596
1597 // TODO: Move this code to after the for loop so edges can be added at top level entities
1598 FmsDomainSetNumEntities(domains[0], FMS_EDGE, FMS_INT32, edge_verts.size() / 2);
1599 FmsDomainAddEntities(domains[0], FMS_EDGE, reorder, FMS_INT32,
1600 edge_verts.data(), edge_verts.size() / 2);
1601#ifdef DEBUG_MFEM_FMS
1602 mfem::out << "EDGES: ";
1603 for (int i = 0; i < edge_verts.size(); i++)
1604 {
1605 if (i % 2 == 0) { mfem::out << std::endl << "\t" << i/2 << ": "; }
1606 mfem::out << edge_verts[i] << " ";
1607 }
1608 mfem::out << std::endl;
1609#endif
1610
1611 // Build faces
1612 if (faces)
1613 {
1614 // TODO: Support Triangles and Quads, and move this code after the for
1615 // loop so these can be added as top level entities
1616 int rowsize = faces->RowSize(0);
1617 std::vector<int> face_edges(faces->Size() * rowsize);
1618 for (int i = 0; i < faces->Size(); i++)
1619 {
1620 mfem::Array<int> eids;
1621 faces->GetRow(i, eids);
1622 for (int j = 0; j < rowsize; j++)
1623 {
1624 face_edges[i*rowsize + j] = eids[j];
1625 }
1626 }
1627 FmsEntityType ent_type = (rowsize == 3) ? FMS_TRIANGLE : FMS_QUADRILATERAL;
1628 FmsDomainSetNumEntities(domains[0], ent_type, FMS_INT32,
1629 face_edges.size() / rowsize);
1630 FmsDomainAddEntities(domains[0], ent_type, NULL, FMS_INT32, face_edges.data(),
1631 face_edges.size() / rowsize);
1632#ifdef DEBUG_MFEM_FMS
1633 mfem::out << "FACES: ";
1634 for (int i = 0; i < face_edges.size(); i++)
1635 {
1636 if (i % rowsize == 0) { mfem::out << std::endl << "\t" << i/rowsize << ": "; }
1637 mfem::out << "(" << edge_verts[face_edges[i]*2] << ", " <<
1638 edge_verts[face_edges[i]*2+1] << ") ";
1639 }
1640 mfem::out << std::endl;
1641#endif
1642 }
1643
1644 // Add top level elements
1645 std::vector<int> tags;
1646 std::vector<int> tris;
1647 std::vector<int> quads;
1648 std::vector<int> tets;
1649 std::vector<int> hexes;
1650 for (int i = 0; i < num_elements; i++)
1651 {
1652 auto etype = mmesh->GetElementType(i);
1653 tags.push_back(mmesh->GetAttribute(i));
1654 switch (etype)
1655 {
1656 case mfem::Element::POINT:
1657 {
1658 // TODO: ?
1659 break;
1660 }
1661 case mfem::Element::SEGMENT:
1662 {
1663 // TODO: ?
1664 break;
1665 }
1666 case mfem::Element::TRIANGLE:
1667 {
1668 mfem::Array<int> eids, oris;
1669 mmesh->GetElementEdges(i, eids, oris);
1670 for (int e = 0; e < 3; e++)
1671 {
1672 tris.push_back(eids[e]);
1673 }
1674 break;
1675 }
1676 case mfem::Element::QUADRILATERAL:
1677 {
1678 mfem::Array<int> eids, oris;
1679 mmesh->GetElementEdges(i, eids, oris);
1680 for (int e = 0; e < 4; e++)
1681 {
1682 quads.push_back(eids[e]);
1683 }
1684 break;
1685 }
1686 case mfem::Element::TETRAHEDRON:
1687 {
1688 mfem::Array<int> fids, oris;
1689 mmesh->GetElementFaces(i, fids, oris);
1690 for (int f = 0; f < 4; f++)
1691 {
1692 tets.push_back(fids[f]);
1693 }
1694 break;
1695 }
1696 case mfem::Element::HEXAHEDRON:
1697 {
1698 mfem::Array<int> fids, oris;
1699 mmesh->GetElementFaces(i, fids, oris);
1700 for (int f = 0; f < 6; f++)
1701 {
1702 hexes.push_back(fids[f]);
1703 }
1704 break;
1705 }
1706 default:
1707 mfem::err << "Error, element not implemented." << std::endl;
1708 return 3;
1709 }
1710 }
1711
1712 if (tris.size())
1713 {
1714 FmsDomainSetNumEntities(domains[0], FMS_TRIANGLE, FMS_INT32, tris.size() / 3);
1715 FmsDomainAddEntities(domains[0], FMS_TRIANGLE, reorder, FMS_INT32, tris.data(),
1716 tris.size() / 3);
1717#ifdef DEBUG_MFEM_FMS
1718 mfem::out << "TRIS: ";
1719 for (int i = 0; i < tris.size(); i++)
1720 {
1721 if (i % 3 == 0) { mfem::out << std::endl << "\t" << i/3 << ": "; }
1722 mfem::out << tris[i] << " ";
1723 }
1724 mfem::out << std::endl;
1725#endif
1726 }
1727
1728 if (quads.size())
1729 {
1730 // TODO: Not quite right either, if there are hexes and quads then this
1731 // will overwrite the faces that made up the hexes
1732 FmsDomainSetNumEntities(domains[0], FMS_QUADRILATERAL, FMS_INT32,
1733 quads.size() / 4);
1734 FmsDomainAddEntities(domains[0], FMS_QUADRILATERAL, reorder, FMS_INT32,
1735 quads.data(), quads.size() / 4);
1736#ifdef DEBUG_MFEM_FMS
1737 mfem::out << "QUADS: ";
1738 for (int i = 0; i < quads.size(); i++)
1739 {
1740 if (i % 4 == 0) { mfem::out << std::endl << "\t" << i/4 << ": "; }
1741 mfem::out << quads[i] << " ";
1742 }
1743 mfem::out << std::endl;
1744#endif
1745 }
1746
1747 if (tets.size())
1748 {
1749 FmsDomainSetNumEntities(domains[0], FMS_TETRAHEDRON, FMS_INT32,
1750 tets.size() / 4);
1751 FmsDomainAddEntities(domains[0], FMS_TETRAHEDRON, reorder, FMS_INT32,
1752 tets.data(), tets.size() / 4);
1753#ifdef DEBUG_MFEM_FMS
1754 mfem::out << "TETS: ";
1755 for (int i = 0; i < tets.size(); i++)
1756 {
1757 if (i % 4 == 0) { mfem::out << std::endl << "\t" << i/4 << ": "; }
1758 mfem::out << tets[i] << " ";
1759 }
1760 mfem::out << std::endl;
1761#endif
1762 }
1763
1764 if (hexes.size())
1765 {
1766 FmsDomainSetNumEntities(domains[0], FMS_HEXAHEDRON, FMS_INT32,
1767 hexes.size() / 6);
1768 FmsDomainAddEntities(domains[0], FMS_HEXAHEDRON, reorder, FMS_INT32,
1769 hexes.data(), hexes.size() / 6);
1770#ifdef DEBUG_MFEM_FMS
1771 mfem::out << "HEXES: ";
1772 for (int i = 0; i < hexes.size(); i++)
1773 {
1774 if (i % 6 == 0) { mfem::out << std::endl << "\t" << i/6 << ": "; }
1775 mfem::out << hexes[i] << " ";
1776 }
1777 mfem::out << std::endl;
1778#endif
1779 }
1780
1781 err = FmsMeshFinalize(*fmesh);
1782 if (err)
1783 {
1784 mfem::err << "FmsMeshFinalize returned error code " << err << std::endl;
1785 return 4;
1786 }
1787
1788 err = FmsMeshValidate(*fmesh);
1789 if (err)
1790 {
1791 mfem::err << "FmsMeshValidate returned error code " << err << std::endl;
1792 return 5;
1793 }
1794
1795 FmsComponent v = NULL;
1796 FmsMeshAddComponent(*fmesh, "volume", &v);
1797 FmsComponentAddDomain(v, domains[0]);
1798
1799 FmsTag tag;
1800 FmsMeshAddTag(*fmesh, "element_attribute", &tag);
1801 FmsTagSetComponent(tag, v);
1802 FmsTagSet(tag, FMS_INT32, FMS_INT32, tags.data(), tags.size());
1803
1804 // Add boundary component
1805 std::vector<int> bdr_eles[FMS_NUM_ENTITY_TYPES];
1806 std::vector<int> bdr_attributes;
1807 const int NBE = mmesh->GetNBE();
1808 for (int i = 0; i < NBE; i++)
1809 {
1810 const Element::Type betype = mmesh->GetBdrElementType(i);
1811 bdr_attributes.push_back(mmesh->GetBdrAttribute(i));
1812 switch (betype)
1813 {
1814 case Element::POINT:
1815 bdr_eles[FMS_VERTEX].push_back(mmesh->GetBdrElementFaceIndex(i));
1816 break;
1817 case Element::SEGMENT:
1818 bdr_eles[FMS_EDGE].push_back(mmesh->GetBdrElementFaceIndex(i));
1819 break;
1820 case Element::TRIANGLE:
1821 bdr_eles[FMS_TRIANGLE].push_back(mmesh->GetBdrElementFaceIndex(i));
1822 break;
1823 case Element::QUADRILATERAL:
1824 bdr_eles[FMS_QUADRILATERAL].push_back(mmesh->GetBdrElementFaceIndex(i));
1825 break;
1826 case Element::TETRAHEDRON:
1827 bdr_eles[FMS_TETRAHEDRON].push_back(mmesh->GetBdrElementFaceIndex(i));
1828 break;
1829 case Element::HEXAHEDRON:
1830 bdr_eles[FMS_HEXAHEDRON].push_back(mmesh->GetBdrElementFaceIndex(i));
1831 break;
1832 default:
1833 MFEM_WARNING("Unsupported boundary element " << betype << " at boundary index "
1834 << i);
1835 break;
1836 }
1837 }
1838
1839 if (NBE)
1840 {
1841 FmsComponent boundary = NULL;
1842 FmsMeshAddComponent(*fmesh, "boundary", &boundary);
1843 FmsInt part_id;
1844 FmsComponentAddPart(boundary, domains[0], &part_id);
1845 for (int i = FMS_NUM_ENTITY_TYPES - 1; i > 0; i--)
1846 {
1847 if (bdr_eles[i].size())
1848 {
1849 FmsComponentAddPartEntities(boundary, part_id, (FmsEntityType)i,
1850 FMS_INT32, FMS_INT32, FMS_INT32, NULL,
1851 bdr_eles[i].data(),
1852 NULL, bdr_eles[i].size());
1853 break;
1854 }
1855 }
1856 FmsComponentAddRelation(v, 1);
1857 FmsTag boundary_tag = NULL;
1858 FmsMeshAddTag(*fmesh, "boundary_attribute", &boundary_tag);
1859 FmsTagSetComponent(boundary_tag, boundary);
1860 FmsTagSet(boundary_tag, FMS_INT32, FMS_INT32, bdr_attributes.data(),
1861 bdr_attributes.size());
1862 }
1863 *volume = v;
1864 return 0;
1865}
1866
1867int
1868DataCollectionToFmsDataCollection(DataCollection *mfem_dc,
1869 FmsDataCollection *dc)
1870{
1871 // TODO: Write me
1872 int err = 0;
1873 const Mesh *mmesh = mfem_dc->GetMesh();
1874
1875 FmsMesh fmesh = NULL;
1876 FmsComponent volume = NULL;
1877 err = MeshToFmsMesh(mmesh, &fmesh, &volume);
1878 if (!fmesh || !volume || err)
1879 {
1880 mfem::err << "Error converting mesh topology from MFEM to FMS" << std::endl;
1881 if (fmesh)
1882 {
1883 FmsMeshDestroy(&fmesh);
1884 }
1885 return 1;
1886 }
1887
1888 err = FmsDataCollectionCreate(fmesh, mfem_dc->GetCollectionName().c_str(), dc);
1889 if (!*dc || err)
1890 {
1891 mfem::err << "There was an error creating the FMS data collection." <<
1892 std::endl;
1893 FmsMeshDestroy(&fmesh);
1894 if (*dc)
1895 {
1896 FmsDataCollectionDestroy(dc);
1897 }
1898 return 3;
1899 }
1900
1901 // Add the coordinates field to the data collection
1902 const mfem::GridFunction *mcoords = mmesh->GetNodes();
1903 if (mcoords)
1904 {
1905 FmsField fcoords = NULL;
1906 err = GridFunctionToFmsField(*dc, volume, "CoordsDescriptor", "Coords", mmesh,
1907 mcoords, &fcoords);
1908 err |= FmsComponentSetCoordinates(volume, fcoords);
1909 }
1910 else
1911 {
1912 // Sometimes the nodes are stored as just a vector of vertex coordinates
1913 mfem::Vector mverts;
1914 mmesh->GetVertices(mverts);
1915 FmsFieldDescriptor fdcoords = NULL;
1916 FmsField fcoords = NULL;
1917 err = FmsDataCollectionAddFieldDescriptor(*dc, "CoordsDescriptor", &fdcoords);
1918 err |= FmsFieldDescriptorSetComponent(fdcoords, volume);
1919 err |= FmsFieldDescriptorSetFixedOrder(fdcoords, FMS_CONTINUOUS,
1920 FMS_NODAL_GAUSS_CLOSED, 1);
1921 err |= FmsDataCollectionAddField(*dc, "Coords", &fcoords);
1922 err |= FmsFieldSet(fcoords, fdcoords, mmesh->SpaceDimension(), FMS_BY_NODES,
1923 FMS_DOUBLE, mverts.HostRead());
1924 err |= FmsComponentSetCoordinates(volume, fcoords);
1925 }
1926
1927 if (err)
1928 {
1929 mfem::err << "There was an error setting the mesh coordinates." << std::endl;
1930 FmsMeshDestroy(&fmesh);
1931 FmsDataCollectionDestroy(dc);
1932 return 4;
1933 }
1934
1935 const auto &fields = mfem_dc->GetFieldMap();
1936 for (const auto &pair : fields)
1937 {
1938 std::string fd_name(pair.first + "Descriptor");
1939 FmsField field;
1940 err = GridFunctionToFmsField(*dc, volume, fd_name, pair.first.c_str(), mmesh,
1941 pair.second, &field);
1942 if (err)
1943 {
1944 mfem::err << "WARNING: There was an error adding the " << pair.first <<
1945 " field. Continuing..." << std::endl;
1946 }
1947 }
1948
1949 // /* TODO:
1950 // const auto &qfields = mfem_dc->GetQFieldMap();
1951 // for(const auto &pair : qfields) {
1952 // FmsFieldDescriptor fd = NULL;
1953 // FmsField f = NULL;
1954 // std::string fd_name(pair.first + "Collection");
1955 // FmsDataCollectionAddFieldDescriptor(*dc, fd_name.c_str(), &fd);
1956 // FmsDataCollectionAddField(*dc, pair.first.c_str(), &f);
1957 // GridFunctionToFmsField(*dc, fd, f, volume, pair.second); // TODO: Volume isn't always going to be correct
1958 // } */
1959
1960 MfemMetaDataToFmsMetaData(mfem_dc, *dc);
1961
1962 return 0;
1963}
1964
1965} // namespace mfem
1966#endif
mfem::Array::Size
int Size() const
Return the logical size of the array.
Definition array.hpp:144
mfem::Array::MakeRef
void MakeRef(T *data_, int size_, bool own_data=false)
Make this Array a reference to a pointer.
Definition array.hpp:882
mfem::Array::Append
int Append(const T &el)
Append element 'el' to array, resize if necessary.
Definition array.hpp:769
mfem::Array::GetData
T * GetData()
Returns the data.
Definition array.hpp:118
mfem::BasisType::ClosedGL
@ ClosedGL
Closed GaussLegendre.
Definition fe_base.hpp:38
mfem::BasisType::OpenHalfUniform
@ OpenHalfUniform
Nodes: x_i = (i+1/2)/n, i=0,...,n-1.
Definition fe_base.hpp:36
mfem::BasisType::Serendipity
@ Serendipity
Serendipity basis (squares / cubes)
Definition fe_base.hpp:37
mfem::BasisType::GaussLobatto
@ GaussLobatto
Closed type.
Definition fe_base.hpp:32
mfem::BasisType::GaussLegendre
@ GaussLegendre
Open type.
Definition fe_base.hpp:31
mfem::BasisType::Positive
@ Positive
Bernstein polynomials.
Definition fe_base.hpp:33
mfem::BasisType::OpenUniform
@ OpenUniform
Nodes: x_i = (i+1)/(n+1), i=0,...,n-1.
Definition fe_base.hpp:34
mfem::BasisType::ClosedUniform
@ ClosedUniform
Nodes: x_i = i/(n-1), i=0,...,n-1.
Definition fe_base.hpp:35
mfem::BasisType::Name
static const char * Name(int b_type)
Check and convert a BasisType constant to a string identifier.
Definition fe_base.hpp:92
mfem::DataCollection::GetCollectionName
const std::string & GetCollectionName() const
Get the name of the collection.
Definition datacollection.hpp:334
mfem::DataCollection::GetCycle
int GetCycle() const
Get time cycle (for time-dependent simulations)
Definition datacollection.hpp:327
mfem::DataCollection::SetTimeStep
void SetTimeStep(real_t ts)
Set the simulation time step (for time-dependent simulations)
Definition datacollection.hpp:324
mfem::DataCollection::RegisterField
virtual void RegisterField(const std::string &field_name, GridFunction *gf)
Add a grid function to the collection.
Definition datacollection.hpp:244
mfem::DataCollection::SetCycle
void SetCycle(int c)
Set time cycle (for time-dependent simulations)
Definition datacollection.hpp:319
mfem::DataCollection::SetTime
void SetTime(real_t t)
Set physical time (for time-dependent simulations)
Definition datacollection.hpp:321
mfem::DataCollection::GetTimeStep
real_t GetTimeStep() const
Get the simulation time step (for time-dependent simulations)
Definition datacollection.hpp:331
mfem::DataCollection::GetFieldMap
const FieldMapType & GetFieldMap() const
Get a const reference to the internal field map.
Definition datacollection.hpp:294
mfem::DataCollection::GetMesh
Mesh * GetMesh()
Get a pointer to the mesh in the collection.
Definition datacollection.hpp:304
mfem::DataCollection::SetOwnData
void SetOwnData(bool o)
Set the ownership of collection data.
Definition datacollection.hpp:336
mfem::DataCollection::GetTime
real_t GetTime() const
Get physical time (for time-dependent simulations)
Definition datacollection.hpp:329
mfem::Element::Type
Type
Constants for the classes derived from Element.
Definition element.hpp:41
mfem::FiniteElementCollection
Collection of finite elements from the same family in multiple dimensions. This class is used to matc...
Definition fe_coll.hpp:27
mfem::FiniteElementCollection::DofOrderForOrientation
virtual const int * DofOrderForOrientation(Geometry::Type GeomType, int Or) const =0
Returns an array, say p, that maps a local permuted index i to a local base index: base_i = p[i].
mfem::FiniteElementCollection::DISCONTINUOUS
@ DISCONTINUOUS
Field is discontinuous across element interfaces.
Definition fe_coll.hpp:48
mfem::FiniteElementCollection::CONTINUOUS
@ CONTINUOUS
Field is continuous across element interfaces.
Definition fe_coll.hpp:45
mfem::FiniteElementCollection::NORMAL
@ NORMAL
Normal component of vector field.
Definition fe_coll.hpp:47
mfem::FiniteElementCollection::TANGENTIAL
@ TANGENTIAL
Tangential components of vector field.
Definition fe_coll.hpp:46
mfem::FiniteElementCollection::GetContType
virtual int GetContType() const =0
mfem::FiniteElementCollection::DofForGeometry
virtual int DofForGeometry(Geometry::Type GeomType) const =0
mfem::FiniteElementCollection::Name
virtual const char * Name() const
Definition fe_coll.hpp:79
mfem::FiniteElementSpace
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Definition fespace.hpp:220
mfem::FiniteElementSpace::GetEdgeInteriorDofs
void GetEdgeInteriorDofs(int i, Array< int > &dofs) const
Returns the indices of the degrees of freedom for the interior of the specified edge.
Definition fespace.cpp:3149
mfem::FiniteElementSpace::GetOrdering
Ordering::Type GetOrdering() const
Return the ordering method.
Definition fespace.hpp:725
mfem::FiniteElementSpace::GetOrder
int GetOrder(int i) const
Returns the polynomial degree of the i'th finite element.
Definition fespace.hpp:696
mfem::FiniteElementSpace::FEColl
const FiniteElementCollection * FEColl() const
Definition fespace.hpp:727
mfem::FiniteElementSpace::GetVDim
int GetVDim() const
Returns vector dimension.
Definition fespace.hpp:706
mfem::FiniteElementSpace::GetFaceInteriorDofs
void GetFaceInteriorDofs(int i, Array< int > &dofs) const
Returns the indices of the degrees of freedom for the interior of the specified face.
Definition fespace.cpp:3125
mfem::FiniteElementSpace::GetElementInteriorDofs
void GetElementInteriorDofs(int i, Array< int > &dofs) const
Returns the indices of the degrees of freedom for the interior of the specified element.
Definition fespace.cpp:3104
mfem::FiniteElementSpace::DofToVDof
int DofToVDof(int dof, int vd, int ndofs=-1) const
Compute a single vdof corresponding to the index dof and the vector index vd.
Definition fespace.cpp:249
mfem::GridFunction
Class for grid function - Vector with associated FE space.
Definition gridfunc.hpp:31
mfem::GridFunction::FESpace
FiniteElementSpace * FESpace()
Definition gridfunc.hpp:696
mfem::GridFunction::VectorDim
int VectorDim() const
Definition gridfunc.cpp:323
mfem::GridFunction::SetSpace
virtual void SetSpace(FiniteElementSpace *f)
Associate a new FiniteElementSpace with the GridFunction.
Definition gridfunc.cpp:195
mfem::H1_FECollection
Arbitrary order H1-conforming (continuous) finite elements.
Definition fe_coll.hpp:260
mfem::HashTable::GetId
int GetId(int p1, int p2)
Get id of item whose parents are p1, p2... Create it if it doesn't exist.
Definition hash.hpp:611
mfem::HashTable::FindId
int FindId(int p1, int p2) const
Find id of item whose parents are p1, p2... Return -1 if it doesn't exist.
Definition hash.hpp:702
mfem::L2_FECollection
Arbitrary order "L2-conforming" discontinuous finite elements.
Definition fe_coll.hpp:330
mfem::Mesh
Mesh data type.
Definition mesh.hpp:56
mfem::Mesh::GetNEdges
int GetNEdges() const
Return the number of edges.
Definition mesh.hpp:1232
mfem::Mesh::AddBdrQuad
int AddBdrQuad(int v1, int v2, int v3, int v4, int attr=1)
Definition mesh.cpp:2063
mfem::Mesh::GetElementType
Element::Type GetElementType(int i) const
Returns the type of element i.
Definition mesh.cpp:7316
mfem::Mesh::GetBdrElementType
Element::Type GetBdrElementType(int i) const
Returns the type of boundary element i.
Definition mesh.cpp:7321
mfem::Mesh::GetAttribute
int GetAttribute(int i) const
Return the attribute of element i.
Definition mesh.hpp:1333
mfem::Mesh::AddQuad
int AddQuad(int v1, int v2, int v3, int v4, int attr=1)
Adds a quadrilateral to the mesh given by 4 vertices v1 through v4.
Definition mesh.cpp:1743
mfem::Mesh::GetBdrAttribute
int GetBdrAttribute(int i) const
Return the attribute of boundary element i.
Definition mesh.hpp:1339
mfem::Mesh::AddTriangle
int AddTriangle(int v1, int v2, int v3, int attr=1)
Adds a triangle to the mesh given by 3 vertices v1 through v3.
Definition mesh.cpp:1729
mfem::Mesh::GetBdrElementFaceIndex
int GetBdrElementFaceIndex(int be_idx) const
Return the local face (codimension-1) index for the given boundary element index.
Definition mesh.hpp:1518
mfem::Mesh::GetVertices
void GetVertices(Vector &vert_coord) const
Definition mesh.cpp:8902
mfem::Mesh::GetNFaces
int GetNFaces() const
Return the number of faces in a 3D mesh.
Definition mesh.hpp:1235
mfem::Mesh::FinalizeTopology
void FinalizeTopology(bool generate_bdr=true)
Finalize the construction of the secondary topology (connectivity) data of a Mesh.
Definition mesh.cpp:3135
mfem::Mesh::AddVertex
int AddVertex(real_t x, real_t y=0.0, real_t z=0.0)
Definition mesh.cpp:1658
mfem::Mesh::GetNE
int GetNE() const
Returns number of elements.
Definition mesh.hpp:1226
mfem::Mesh::AddTet
int AddTet(int v1, int v2, int v3, int v4, int attr=1)
Adds a tetrahedron to the mesh given by 4 vertices v1 through v4.
Definition mesh.cpp:1757
mfem::Mesh::AddBdrSegment
int AddBdrSegment(int v1, int v2, int attr=1)
Definition mesh.cpp:2035
mfem::Mesh::GetElementFaces
void GetElementFaces(int i, Array< int > &faces, Array< int > &ori) const
Return the indices and the orientations of all faces of element i.
Definition mesh.cpp:7201
mfem::Mesh::SpaceDimension
int SpaceDimension() const
Dimension of the physical space containing the mesh.
Definition mesh.hpp:1163
mfem::Mesh::GetNodes
void GetNodes(Vector &node_coord) const
Definition mesh.cpp:8973
mfem::Mesh::AddBdrTriangle
int AddBdrTriangle(int v1, int v2, int v3, int attr=1)
Definition mesh.cpp:2049
mfem::Mesh::GetNV
int GetNV() const
Returns number of vertices. Vertices are only at the corners of elements, where you would expect them...
Definition mesh.hpp:1223
mfem::Mesh::GetEdgeVertices
void GetEdgeVertices(int i, Array< int > &vert) const
Returns the indices of the vertices of edge i.
Definition mesh.cpp:7069
mfem::Mesh::GetFaceVertices
void GetFaceVertices(int i, Array< int > &vert) const
Returns the indices of the vertices of face i.
Definition mesh.hpp:1456
mfem::Mesh::AddHex
int AddHex(int v1, int v2, int v3, int v4, int v5, int v6, int v7, int v8, int attr=1)
Adds a hexahedron to the mesh given by 8 vertices v1 through v8.
Definition mesh.cpp:1806
mfem::Mesh::GetNBE
int GetNBE() const
Returns number of boundary elements.
Definition mesh.hpp:1229
mfem::Mesh::Finalize
virtual void Finalize(bool refine=false, bool fix_orientation=false)
Finalize the construction of a general Mesh.
Definition mesh.cpp:3241
mfem::Mesh::GetFaceEdgeTable
Table * GetFaceEdgeTable() const
Definition mesh.cpp:7078
mfem::Mesh::GetElementEdges
void GetElementEdges(int i, Array< int > &edges, Array< int > &cor) const
Return the indices and the orientations of all edges of element i.
Definition mesh.cpp:6985
mfem::Mesh::SetCurvature
virtual void SetCurvature(int order, bool discont=false, int space_dim=-1, int ordering=1)
Set the curvature of the mesh nodes using the given polynomial degree.
Definition mesh.cpp:6211
mfem::Mesh::GetEdgeVertexTable
Table * GetEdgeVertexTable() const
Definition mesh.cpp:7106
mfem::RT_FECollection
Arbitrary order H(div)-conforming Raviart-Thomas finite elements.
Definition fe_coll.hpp:386
mfem::Table::RowSize
int RowSize(int i) const
Definition table.hpp:108
mfem::Table::GetRow
void GetRow(int i, Array< int > &row) const
Return row i in array row (the Table must be finalized)
Definition table.cpp:187
mfem::Table::Size
int Size() const
Returns the number of TYPE I elements.
Definition table.hpp:92
mfem::Vector
Vector data type.
Definition vector.hpp:80
mfem::Vector::HostRead
virtual const real_t * HostRead() const
Shortcut for mfem::Read(vec.GetMemory(), vec.Size(), false).
Definition vector.hpp:478
mfem::Vector::Size
int Size() const
Returns the size of the vector.
Definition vector.hpp:218
mfem::Vector::GetData
real_t * GetData() const
Return a pointer to the beginning of the Vector data.
Definition vector.hpp:227
mfem::VisItDataCollection
Data collection with VisIt I/O routines.
Definition datacollection.hpp:420
mfem::VisItDataCollection::Save
virtual void Save() override
Save the collection and a VisIt root file.
Definition datacollection.cpp:471
mfem::VisItDataCollection::RegisterField
virtual void RegisterField(const std::string &field_name, GridFunction *gf) override
Add a grid function to the collection and update the root file.
Definition datacollection.cpp:415
dim
int dim
Definition ex24.cpp:53
b
real_t b
Definition lissajous.cpp:42
mfem::FmsMetaDataGetScalar
bool FmsMetaDataGetScalar(FmsMetaData mdata, const std::string &key, std::vector< double > &values)
Definition fmsconvert.cpp:1273
mfem::FmsMetaDataGetString
bool FmsMetaDataGetString(FmsMetaData mdata, const std::string &key, std::string &value)
Definition fmsconvert.cpp:1338
mfem::GridFunctionToFmsField
int GridFunctionToFmsField(FmsDataCollection dc, FmsComponent comp, const std::string &fd_name, const std::string &field_name, const Mesh *mesh, const GridFunction *gf, FmsField *outfield)
Definition fmsconvert.cpp:948
mfem::out
OutStream out(std::cout)
Global stream used by the library for standard output. Initially it uses the same std::streambuf as s...
Definition globals.hpp:66
mfem::FmsFieldToGridFunction
int FmsFieldToGridFunction(FmsMesh fms_mesh, FmsField f, Mesh *mesh, GridFunction &func, bool setFE)
This function converts an FmsField to an MFEM GridFunction.
Definition fmsconvert.cpp:107
mfem::BasisTypeToFmsBasisType
bool BasisTypeToFmsBasisType(int bt, FmsBasisType &btype)
Definition fmsconvert.cpp:892
mfem::FmsDataCollectionToDataCollection
int FmsDataCollectionToDataCollection(FmsDataCollection dc, DataCollection **mfem_dc)
Definition fmsconvert.cpp:1386
mfem::MeshToFmsMesh
int MeshToFmsMesh(const Mesh *mmesh, FmsMesh *fmesh, FmsComponent *volume)
Definition fmsconvert.cpp:1531
mfem::err
OutStream err(std::cerr)
Global stream used by the library for standard error output. Initially it uses the same std::streambu...
Definition globals.hpp:71
mfem::DataCollectionToFmsDataCollection
int DataCollectionToFmsDataCollection(DataCollection *mfem_dc, FmsDataCollection *dc)
Definition fmsconvert.cpp:1868
mfem::MfemMetaDataToFmsMetaData
bool MfemMetaDataToFmsMetaData(DataCollection *mdc, FmsDataCollection fdc)
Definition fmsconvert.cpp:1100
mfem::f
std::function< real_t(const Vector &)> f(real_t mass_coeff)
Definition lor_mms.hpp:30
mfem::FmsMetaDataGetInteger
bool FmsMetaDataGetInteger(FmsMetaData mdata, const std::string &key, std::vector< int > &values)
Definition fmsconvert.cpp:1172
mfem::FmsMeshToMesh
int FmsMeshToMesh(FmsMesh fms_mesh, Mesh **mfem_mesh)
Definition fmsconvert.cpp:452
doxygen 1.11.0