d5a9339d1121b4da7e7f6d747f2b55e9f8fd49e8
[charm.git] / src / libs / ck-libs / ParFUM / ParFUM.h
1 /*Charm++ Finite Element Framework:
2 C interface file
3 */
4 #ifndef __PARFUM_H
5 #define __PARFUM_H
6 #include "pup_c.h"  /* for pup_er */
7 #include "idxlc.h"
8 #include "collidec.h"  // from collision framework
9 #include "charm-api.h"
10 #include "ckvector3d.h"
11 #include "charm++.h"
12 #include "tcharm.h"
13 #include "charm-api.h"
14 #include "ckvector3d.h"
15
16
17 // Forward declaration
18 class FEM_Entity;
19 class FEM_Mesh;
20 class FEM_Elem;
21 class FEM_Node;
22 class femMeshModify;
23 class FEM_Adapt_Algs;
24 class FEM_Adapt;
25 class FEM_AdaptL;
26 class IDXL_Chunk;
27 class l2g_t;
28 class FEM_ElemAdj_Layer;
29 class chunkListMsg;
30
31 /* BUG: this should not be used */
32 void _registerParFUM(void);
33
34 #ifdef __cplusplus
35 extern "C" {
36 #endif
37
38 /* datatypes: keep in sync with ParFUMf.h and idxl */
39 #define FEM_BYTE   IDXL_BYTE
40 #define FEM_INT    IDXL_INT
41 #define FEM_REAL   IDXL_REAL
42 #define FEM_FLOAT FEM_REAL /*alias*/
43 #define FEM_DOUBLE IDXL_DOUBLE
44 #define FEM_INDEX_0  IDXL_INDEX_0
45 #define FEM_INDEX_1  IDXL_INDEX_1 
46 #define FEM_VAR_INDEX (IDXL_FIRST_DATATYPE+6)
47
48 /* reduction operations: keep in sync with ParFUMf.h */
49 #define FEM_SUM IDXL_SUM
50 #define FEM_PROD IDXL_PROD
51 #define FEM_MAX IDXL_MAX
52 #define FEM_MIN IDXL_MIN
53
54 /* element types, by their number of nodes */
55 #define FEM_TRIANGULAR    3
56 #define FEM_TETRAHEDRAL   4
57 #define FEM_HEXAHEDRAL    8
58 #define FEM_QUADRILATERAL 4
59
60 /* initialization flags */
61 #define FEM_INIT_READ    2
62 #define FEM_INIT_WRITE   4
63
64 #define FEM_MESH_OUTPUT 0
65 #define FEM_MESH_UPDATE 1
66 #define FEM_MESH_FINALIZE 2
67   typedef void (*FEM_Update_mesh_fn)(int userTag);
68   typedef void (*FEM_Update_mesh_fortran_fn)(int *userTag);
69
70   typedef void (*FEM_PupFn)(pup_er, void*);
71
72   typedef void (*FEM_Mesh_alloc_fn)(void *param,int *size,int *maxSize);
73   
74   /* This should be MPI_Comm, but I want it for Fortran too */
75   typedef int FEM_Comm_t; 
76
77   /* Initialize the FEM framework (must have called MPI_Init) */
78   void FEM_Init(FEM_Comm_t defaultCommunicator);
79   void FEM_Done(void);
80
81   /*Utility*/
82   int FEM_My_partition(void);
83   int FEM_Num_partitions(void);
84   double FEM_Timer(void);
85   void FEM_Print(const char *str);
86   void FEM_Print_partition(void);
87
88 /* Mesh manipulation */
89 #define FEM_MESH_FIRST 1650000000 /*This is the first mesh ID:*/
90   /* mesh creation */
91   int FEM_Mesh_allocate(void); /* build new mesh */
92   int FEM_Mesh_copy(int fem_mesh); /* copy existing mesh */
93   void FEM_Mesh_deallocate(int fem_mesh); /* delete this mesh */
94
95   int FEM_Mesh_read(const char *prefix,int partNo,int nParts);
96   void FEM_Mesh_write(int fem_mesh,const char *prefix,int partNo,int nParts); 
97
98   int FEM_Mesh_assemble(int nParts,const int *srcMeshes);
99   void FEM_Mesh_partition(int fem_mesh,int nParts,int *destMeshes);
100   
101   int FEM_Mesh_recv(int fromRank,int tag,FEM_Comm_t comm_context);
102   void FEM_Mesh_send(int fem_mesh,int toRank,int tag,FEM_Comm_t comm_context);
103
104   int FEM_Mesh_reduce(int fem_mesh,int toRank,FEM_Comm_t comm_context);
105   int FEM_Mesh_broadcast(int fem_mesh,int fromRank,FEM_Comm_t comm_context);
106
107   void FEM_Mesh_copy_globalno(int src_mesh,int dest_mesh);
108   void FEM_Mesh_print(int fem_mesh);
109   
110 /* Mesh entity codes: (keep in sync with ParFUMf.h) */
111 #define FEM_ENTITY_FIRST 1610000000 /*This is the first entity code:*/
112 #define FEM_NODE (FEM_ENTITY_FIRST+0) /*The unique node type*/
113 #define FEM_ELEM (FEM_ENTITY_FIRST+1000) /*First element type (can add the user-defined element type) */
114 #define FEM_ELEMENT FEM_ELEM /*alias*/
115 #define FEM_SPARSE (FEM_ENTITY_FIRST+2000) /* First sparse entity (face) type */
116 #define FEM_EDGE FEM_SPARSE /* alias */
117 #define FEM_FACE FEM_SPARSE /* alias */
118 #define FEM_GHOST 10000  /* (entity add-in) Indicates we want the ghost values, not real values */
119 #define FEM_ENTITY_LAST (FEM_ENTITY_FIRST+3000+FEM_GHOST)
120
121 /* Mesh entity "attributes": per-entity data */
122 #define FEM_DATA   0  /* Backward-compatability routines' solution data: tag 0 */
123 #define FEM_ATTRIB_TAG_MAX 1000000000 /*Largest allowable user "tag" attribute*/
124 #define FEM_ATTRIB_FIRST 1620000000 /*This is the first system attribute code: one of*/
125 #define FEM_CONN   (FEM_ATTRIB_FIRST+1) /* Element-node connectivity (FEM_ELEM or FEM_SPARSE, FEM_INDEX only) */
126 #define FEM_CONNECTIVITY FEM_CONN /*alias*/
127
128   /* rarely-used external attributes */
129 #define FEM_SPARSE_ELEM (FEM_ATTRIB_FIRST+2) /* Elements each sparse data record applies to (FEM_SPARSE, 2*FEM_INDEX only) */
130 #define FEM_COOR   (FEM_ATTRIB_FIRST+3) /* Node coordinates (FEM_NODE, FEM_DOUBLE only) */
131 #define FEM_COORD FEM_COOR /*alias*/
132 #define FEM_COORDINATES FEM_COOR /*alias*/
133 #define FEM_GLOBALNO  (FEM_ATTRIB_FIRST+4) /* Global item numbers (width=1, datatype=FEM_INDEX) */
134 #define FEM_PARTITION (FEM_ATTRIB_FIRST+5) /* Destination chunk numbers (elements only; width=1, datatype=FEM_INDEX) */
135 #define FEM_SYMMETRIES (FEM_ATTRIB_FIRST+6) /* Symmetries present (width=1, datatype=FEM_BYTE) */
136 #define FEM_NODE_PRIMARY (FEM_ATTRIB_FIRST+7) /* This chunk owns this node (nodes only; width=1, datatype=FEM_BYTE) */
137 #define FEM_CHUNK (FEM_ATTRIB_FIRST+8) /* For Nodes and Elements. Used during ghost creation
138 to mark the chunk to which a ghost node or element belongs datatype=FEM_INDEX*/
139 #define FEM_BOUNDARY (FEM_ATTRIB_FIRST+9) /*provides the boundary flag for nodes, elements and sparse elements FEM_INT*/
140 #define FEM_NODE_ELEM_ADJACENCY (FEM_ATTRIB_FIRST+10) /*node to element adjacency FEM_VAR_INDEX only */
141 #define FEM_NODE_NODE_ADJACENCY (FEM_ATTRIB_FIRST+11) /*node to node adjacency FEM_VAR_INDEX only */
142 #define FEM_ELEM_ELEM_ADJACENCY (FEM_ATTRIB_FIRST+12) /*element to element adjacency FEM_VAR_INDEX only */
143 #define FEM_ELEM_ELEM_ADJ_TYPES (FEM_ATTRIB_FIRST+13) /*stores element types for those element id's listed in 
144                                                                                                                 FEM_ELEM_ELEM_ADJACENCY, needed when using 
145                                                                                                                 multiple element types*/
146 #define FEM_IS_VALID_ATTR (FEM_ATTRIB_FIRST+14) /* Stores a flag(an IDXL_BYTE) for each element or node specifying whether the entity 
147                                                                                    exists or is valid. It may be 0 whenever a mesh modification occurs that deletes the 
148                                                                                    corresponding node or element */
149
150 #define FEM_MESH_SIZING (FEM_ATTRIB_FIRST+15) /* Target edge length attr. */
151 #define FEM_ATTRIB_LAST (FEM_ATTRIB_FIRST+16) /*This is the last valid attribute code*/
152
153   /* Specialized routines: */
154   void FEM_Mesh_set_conn(int fem_mesh,int entity,
155         const int *conn, int firstItem, int length, int width);
156   void FEM_Mesh_get_conn(int fem_mesh,int entity,
157         int *conn, int firstItem, int length, int width);
158
159   void FEM_Mesh_set_data(int fem_mesh,int entity,int attr,
160         const void *data, int firstItem, int length, int datatype,int width);
161   void FEM_Mesh_get_data(int fem_mesh,int entity,int attr,
162         void *data, int firstItem, int length, int datatype,int width);
163   void FEM_Mesh_conn(int fem_mesh,int entity,
164         int *conn, int firstItem, int length, int width);
165   
166   int FEM_Mesh_get_length(int fem_mesh,int entity);
167   
168   /* General purpose routines: */
169   void FEM_Mesh_data(int fem_mesh,int entity,int attr,
170         void *data, int firstItem, int length, int datatype,int width);
171   void FEM_Mesh_data_layout(int fem_mesh,int entity,int attr,
172         void *data, int firstItem, int length, IDXL_Layout_t layout);
173   void FEM_Mesh_data_offset(int fem_mesh,int entity,int attr,
174         void *data, int firstItem, int length, 
175         int type,int width, int offsetBytes,int distanceBytes,int skewBytes);
176         
177   void FEM_Register_array(int fem_mesh,int entity,int attr,
178                           void *data, int datatype,int width);
179
180   void FEM_Register_array_layout(int fem_mesh,int entity,int attr,      
181                                  void *data, IDXL_Layout_t layout);     
182   
183   //TODO:add the most important parameter.. the function pointer to the resize function
184   void FEM_Register_entity(int fem_mesh,int entity,void *data,int len,int max,FEM_Mesh_alloc_fn fn);    
185   
186   void FEM_Mesh_set_length(int fem_mesh,int entity,int newLength);
187   int FEM_Mesh_get_width(int fem_mesh,int entity,int attr);
188   void FEM_Mesh_set_width(int fem_mesh,int entity,int attr,int newWidth);
189   int FEM_Mesh_get_datatype(int fem_mesh,int entity,int attr);
190   int FEM_Mesh_get_entities(int fem_mesh, int *entities);
191   int FEM_Mesh_get_attributes(int fem_mesh,int entity,int *attributes);
192   
193   const char *FEM_Get_entity_name(int entity,char *storage);
194   const char *FEM_Get_attr_name(int attr,char *storage);
195   const char *FEM_Get_datatype_name(int datatype,char *storage);
196
197   int FEM_Mesh_is_get(int fem_mesh); /* return 1 if this is a readable mesh */
198   int FEM_Mesh_is_set(int fem_mesh); /* return 1 if this is a writing mesh */
199   void FEM_Mesh_become_get(int fem_mesh); /* Make this a readable mesh */
200   void FEM_Mesh_become_set(int fem_mesh); /* Make this a writing mesh */
201
202   typedef void (*FEM_Userdata_fn)(pup_er p,void *data);
203   void FEM_Mesh_pup(int fem_mesh,int dataTag,FEM_Userdata_fn fn,void *data);
204
205 /* ghosts and spatial symmetries */
206 #define FEM_Is_ghost_index(idx) ((idx)<-1)
207 #define FEM_To_ghost_index(idx) (-(idx)-2)
208 #define FEM_From_ghost_index(idx) (-(idx)-2)
209
210   void FEM_Add_ghost_layer(int nodesPerTuple,int doAddNodes);
211   void FEM_Add_ghost_elem(int elType,int tuplesPerElem,const int *elem2tuple);
212
213   void FEM_Add_ghost_stencil(int nElts,int addNodes,
214         const int *ends,const int *adj);
215   void FEM_Add_ghost_stencil_type(int elType,int nElts,int addNodes,
216         const int *ends,const int *adj2);
217
218   void FEM_Add_elem2face_tuples(int fem_mesh, int elem_type, int nodesPerTuple, int tuplesPerElem,const int *elem2tuple);
219
220   void FEM_Add_linear_periodicity(int nFaces,int nPer,
221         const int *facesA,const int *facesB,
222         int nNodes,const double *nodeLocs);
223   void FEM_Sym_coordinates(int who,double *d_locs);
224   
225   void FEM_Set_sym_nodes(const int *canon,const int *sym);
226   void FEM_Get_sym(int who,int *destSym);
227   /**
228     Based on shared node communication list, compute 
229     FEM_NODE FEM_GLOBALNO and FEM_NODE_PRIMARY
230   */
231   void FEM_Make_node_globalno(int fem_mesh,FEM_Comm_t comm_context);
232
233 /* Communication: see idxlc.h */
234   IDXL_Layout_t FEM_Create_simple_field(int base_type,int vec_len);
235   IDXL_Layout_t FEM_Create_field(int base_type, int vec_len, int init_offset, 
236                        int distance);
237   
238   IDXL_t FEM_Comm_shared(int fem_mesh,int entity);
239   IDXL_t FEM_Comm_ghost(int fem_mesh,int entity);
240
241   void FEM_Get_roccom_pconn_size(int fem_mesh,int *total_len,int *ghost_len);
242   void FEM_Get_roccom_pconn(int fem_mesh,const int *paneFmChunk,int *pconn);
243   void FEM_Set_roccom_pconn(int fem_mesh,const int *paneFmChunk,const int *src,int total_len,int ghost_len);
244
245   /*Migration */
246   int FEM_Register(void *userData,FEM_PupFn _pup_ud);
247   void FEM_Migrate(void);
248   void *FEM_Get_userdata(int n);
249   
250   void FEM_Barrier(void);
251   
252   /* to be provided by the application */
253   void init(void);
254   void driver(void);
255   
256   /* Create additional mesh adjacency information */
257   void FEM_Mesh_create_node_elem_adjacency(int fem_mesh);
258   void FEM_Mesh_create_node_node_adjacency(int fem_mesh);
259   void FEM_Mesh_create_elem_elem_adjacency(int fem_mesh);
260
261   void FEM_Print_n2n(int mesh, int nodeid);
262   void FEM_Print_n2e(int mesh, int nodeid);
263   void FEM_Print_e2e(int mesh, int eid);
264   void FEM_Print_e2n(int mesh, int eid);
265
266   /* Create and modify the FEM_IS_VALID Attribute */
267   void FEM_Mesh_allocate_valid_attr(int fem_mesh, int entity_type);
268   void FEM_set_entity_valid(int mesh, int entityType, int entityIdx);
269   void FEM_set_entity_invalid(int mesh, int entityType, int entityIdx);
270   int FEM_is_valid(int mesh, int entityType, int entityIdx);
271   unsigned int FEM_count_valid(int mesh, int entityType);
272
273   /* Easy set method for coordinates, that may be helpful when creating a mesh */
274   void FEM_set_entity_coord2(int mesh, int entityType, int entityIdx, double x, double y);
275   void FEM_set_entity_coord3(int mesh, int entityType, int entityIdx, double x, double y, double z);
276   
277
278   /* Backward compatability routines: */
279   int FEM_Mesh_default_read(void);  /* return mesh used for get calls below */
280   int FEM_Mesh_default_write(void); /* return mesh used for set calls below */
281   void FEM_Mesh_set_default_read(int fem_mesh);
282   void FEM_Mesh_set_default_write(int fem_mesh);
283   
284   void FEM_Exchange_ghost_lists(int who,int nIdx,const int *localIdx);
285   int FEM_Get_ghost_list_length(void);
286   void FEM_Get_ghost_list(int *dest);
287
288   void FEM_Update_field(int fid, void *nodes);
289   void FEM_Update_ghost_field(int fid, int elTypeOrMinusOne, void *nodes);
290   void FEM_Reduce_field(int fid, const void *nodes, void *outbuf, int op);
291   void FEM_Reduce(int fid, const void *inbuf, void *outbuf, int op);
292   
293   void FEM_Read_field(int fid, void *nodes, const char *fname);
294
295   void FEM_Set_node(int nNodes,int doublePerNode);
296   void FEM_Set_node_data(const double *data);
297   void FEM_Set_elem(int elType,int nElem,int doublePerElem,int nodePerElem);
298   void FEM_Set_elem_data(int elType,const double *data);
299   void FEM_Set_elem_conn(int elType,const int *conn);
300   void FEM_Set_sparse(int uniqueIdentifier,int nRecords,
301         const int *nodes,int nodesPerRec,
302         const void *data,int dataPerRec,int dataType);
303   void FEM_Set_sparse_elem(int uniqueIdentifier,const int *rec2elem);
304
305   void FEM_Get_node(int *nNodes,int *doublePerNode);
306   void FEM_Get_node_data(double *data);
307   void FEM_Get_elem(int elType,int *nElem,int *doublePerElem,int *nodePerElem);
308   void FEM_Get_elem_data(int elType,double *data);
309   void FEM_Get_elem_conn(int elType,int *conn);
310   int  FEM_Get_sparse_length(int uniqueIdentifier); 
311   void FEM_Get_sparse(int uniqueIdentifier,int *nodes,void *data);
312   
313   void FEM_Set_mesh(int nelem, int nnodes, int nodePerElem, int* conn);
314   
315   int FEM_Get_node_ghost(void);
316   int FEM_Get_elem_ghost(int elemType);  
317
318   void FEM_Update_mesh(FEM_Update_mesh_fn callFn,int userValue,int doWhat);
319   
320   void FEM_Set_partition(int *elem2chunk);
321
322
323   /* Public functions that modify the mesh */
324   int FEM_add_node(int mesh, int* adjacent_nodes=0, int num_adjacent_nodes=0, int *chunks=0, int numChunks=0, int forceShared=0, int upcall=0);
325   int FEM_add_element(int mesh, int* conn, int conn_size, int elem_type=0, int chunkNo=-1);
326   void FEM_remove_node(int mesh,int node);
327   int FEM_remove_element(int mesh, int element, int elem_type=0, int permanent=-1);
328   int FEM_Modify_Lock(int mesh, int* affectedNodes, int numAffectedNodes, int* affectedElts=0, int numAffectedElts=0, int elemtype=0);
329   int FEM_Modify_Unlock(int mesh);
330   int FEM_Modify_LockN(int mesh, int nodeId, int readLock);
331   int FEM_Modify_UnlockN(int mesh, int nodeId, int readLock);
332   void FEM_REF_INIT(int mesh, int dim);
333   
334  
335
336   
337   // To help debugging:
338   void FEM_Print_Mesh_Summary(int mesh);
339
340 /* Routines we wish didn't exist: */
341   void FEM_Serial_split(int nchunks);
342   void FEM_Serial_begin(int chunkNo);
343   
344   void FEM_Serial_read(int chunkNo,int nChunks);
345   void FEM_Serial_assemble(void);
346   
347   int FEM_Get_comm_partners(void);
348   int FEM_Get_comm_partner(int partnerNo);
349   int FEM_Get_comm_count(int partnerNo);
350   void FEM_Get_comm_nodes(int partnerNo,int *nodeNos);
351
352
353 /* Routines that no longer exist:
354   void FEM_Composite_elem(int newElType);
355    -> Replace with IDXL_Create
356   void FEM_Combine_elem(int srcElType,int destElType,
357         int newInteriorStartIdx,int newGhostStartIdx);  
358    -> Replace with IDXL_Combine
359 */
360
361
362
363
364 /* 
365 ParFUM Collision Interface File
366
367 A few outstanding questions:
368
369 Is the use of an element based attribute to store any needed collision data a good thing?
370 Perhaps we should just use the user data attributes for the element. This may require there 
371 to be consequtive user data attributes. I.e. no FEM_DATA+5,FEM_DATA+82, without those inbetween.
372 Do we need to transmit nodal data for each element?
373 Does the user need anything beyond just some data attributes for the one remote element which is 
374 colliding locally?
375
376 THESE FUNCTIONS ARE NOT YET IMPLEMENTED!
377
378 Author: Isaac Dooley 11-09-2005
379 */
380
381   struct ParFUM_collider {
382         collide_t collide_grid;
383         double box_padding;
384         int dimension;
385
386         unsigned int *boxToElementMapping;
387         unsigned int numCollidableElements; // size of boxToElementMapping array
388   };
389
390
391
392   /* ParFUM_Collide_init() will initialize the collision library. 
393      It should be called once in driver after mesh has been loaded.
394      
395      dimension should reflect the number of coordinates associated 
396      with a node. This cannot exceed 3 with the current Collision
397      Library. The user's nodal coordinates must be registered as a 
398      particular attribute in order to determine the optimal grid sizing.
399
400      Algorithm:
401        Determine Grid Sizing
402        Call COLLIDE_Init()
403      
404   */   
405   ParFUM_collider ParFUM_Collide_Init(int dimension);
406
407
408   /* ParFUM_Collide() will create bounding boxes for each element in the local mesh chunk.
409      It will then collide these bounding boxes with those both locally and remotely.
410      It should be called at each timestep for which collisions are being tested.
411     
412      Algorithm: 
413        Create Bounding boxes for all valid elements, and priority array
414        Call COLLIDE_Boxes_prio()
415        return the number of collisions which involve a local element
416   */  
417   int ParFUM_Collide(ParFUM_collider *c, double box_padding = 0.0);
418
419   /* ParFUM_Collide_GetCollisions() is used to get the data for any remote elements which 
420      It should be called after Collide even if ParFUM_Collide returned 0
421
422      The data it returns will be double precision values associated with the
423      element attribute ParFUM_COLLISION_DATA
424
425      results should be an array allocated by the user with length equal to the number of 
426      collisions times the amount of space needed for each item in the ParFUM_COLLISION_DATA 
427      attribute
428           
429      Algorithm: 
430
431
432
433   */  
434   void ParFUM_Collide_GetCollisions(ParFUM_collider *c, void* results);
435
436   void ParFUM_Collide_Destroy(ParFUM_collider *c);
437
438
439 // End of Collision interface
440
441
442
443 // User functions for adaptivity
444
445 void FEM_ADAPT_Init(int meshID);
446 FDECL void FTN_NAME(FEM_ADAPT_INIT,fem_adapt_init)(int *meshID);
447
448
449 void FEM_ADAPT_Refine(int meshID, int qm, int method, double factor, double *sizes);
450 FDECL void FTN_NAME(FEM_ADAPT_REFINE,fem_adapt_refine)(int* meshID, 
451         int *qm, int *method, double *factor, double *sizes);
452
453
454 void FEM_ADAPT_Coarsen(int meshID, int qm, int method, double factor, 
455         double *sizes);
456 FDECL void FTN_NAME(FEM_ADAPT_COARSEN,fem_adapt_coarsen)(int* meshID, 
457         int *qm, int *method, double *factor, double *sizes);
458
459 void FEM_ADAPT_AdaptMesh(int meshID, int qm, int method, double factor, double *sizes);
460 FDECL void FTN_NAME(FEM_ADAPT_ADAPTMESH,fem_adapt_adaptmesh)(int* meshID, 
461         int *qm, int *method, double *factor, double *sizes);
462
463 void FEM_ADAPT_SetElementSizeField(int meshID, int elem, double size);
464 FDECL void FTN_NAME(FEM_ADAPT_SETELEMENTSIZEFIELD,fem_adapt_setelementsizefield)(int *meshID, int *elem, double *size);
465
466
467 void FEM_ADAPT_SetElementsSizeField(int meshID, double *sizes);
468 FDECL void FTN_NAME(FEM_ADAPT_SETELEMENTSSIZEFIELD,fem_adapt_setelementssizefield)(int *meshID, double *sizes);
469
470
471 void FEM_ADAPT_SetReferenceMesh(int meshID);
472 FDECL void FTN_NAME(FEM_ADAPT_SETREFERENCEMESH, fem_adapt_setreferencemesh)(int* meshID);
473
474
475 void FEM_ADAPT_GradateMesh(int meshID, double smoothness);
476 FDECL void FTN_NAME(FEM_ADAPT_GRADATEMESH, fem_adapt_gradatemesh)(int* meshID, double* smoothness);
477
478   // End Adaptivity interface
479
480 }
481 #endif
482
483