Separated fortran interface into header and implementation. Now it should link
authorAaron Becker <abecker3@illinois.edu>
Thu, 15 Dec 2005 17:21:49 +0000 (17:21 +0000)
committerAaron Becker <abecker3@illinois.edu>
Thu, 15 Dec 2005 17:21:49 +0000 (17:21 +0000)
properly without breaking the fortran interface.

src/libs/ck-libs/fem/Makefile
src/libs/ck-libs/fem/fem_adapt_if.C [new file with mode: 0644]
src/libs/ck-libs/fem/fem_adapt_if.h

index 902d9957c157bdbcbb2cbb509910173be0b5b141..92c03ab171a4d706b39b7c8106bd8743c8eb38c6 100644 (file)
@@ -8,7 +8,8 @@ OBJS=fem.o fem_mesh.o symmetries.o \
        partition.o map.o fem_compat.o call_init.o \
        parallel_part.o fem_mesh_modify.o fem_mesh_adjacency.o \
        fem_adapt_new.o fem_adapt_algs.o fem_interpolate.o \
-       fem_lock.o fem_util.o fem_lock_node.o fem_adapt_lock.o
+       fem_lock.o fem_util.o fem_lock_node.o fem_adapt_lock.o \
+       fem_adapt_if.o
 COMPAT=compat_init.o compat_finit.o compat_driver.o compat_fdriver.o 
 LIB=libmodulefem
 
@@ -97,6 +98,9 @@ fem_lock_node.o: fem_lock_node.C $(HEADDEP)
 fem_adapt_lock.o: fem_adapt_lock.C $(HEADDEP)
        $(CHARMC) -c fem_adapt_lock.C
 
+fem_adapt_if.o: fem_adapt_if.C $(HEADDEP)
+       $(CHARMC) -c fem_adapt_if.C
+
 fem_comm.o: fem_comm.C $(HEADDEP)
        $(CHARMC) -c fem_comm.C
 
diff --git a/src/libs/ck-libs/fem/fem_adapt_if.C b/src/libs/ck-libs/fem/fem_adapt_if.C
new file mode 100644 (file)
index 0000000..28e520e
--- /dev/null
@@ -0,0 +1,89 @@
+#include "fem_adapt_if.h"
+
+
+void FEM_ADAPT_Init(int meshID) {
+  const int triangleFaces[6] = {0,1,1,2,2,0};
+  FEM_Add_elem2face_tuples(meshID, 0, 2, 3, triangleFaces);
+  FEM_Mesh_create_elem_elem_adjacency(meshID);
+  FEM_Mesh_create_node_elem_adjacency(meshID);
+  FEM_Mesh_create_node_node_adjacency(meshID);
+  _registerFEMMeshModify();
+  FEM_REF_INIT(meshID, 2);  // dim=2
+  FEM_Mesh *meshP = FEM_Mesh_lookup(meshID, "FEM_ADAPT_Init");
+  CtvInitialize(FEM_Adapt_Algs *, _adaptAlgs);
+  CtvAccess(_adaptAlgs) = meshP->getfmMM()->getfmAdaptAlgs();
+  CtvAccess(_adaptAlgs)->FEM_Adapt_Algs_Init(FEM_DATA+0, FEM_BOUNDARY);
+}
+FDECL void FTN_NAME(FEM_ADAPT_INIT,fem_adapt_init)(int *meshID)
+{
+  FEM_ADAPT_Init(*meshID);
+}
+
+
+void FEM_ADAPT_Refine(int meshID, int qm, int method, double factor,
+        double *sizes) {
+    FEM_Mesh* mesh = FEM_Mesh_lookup(meshID, "FEM_ADAPT_Init");
+    mesh->getfmMM()->getfmAdaptAlgs()->FEM_Refine(qm, method, factor, sizes);
+}
+FDECL void FTN_NAME(FEM_ADAPT_REFINE,fem_adapt_refine)(int* meshID, 
+        int *qm, int *method, double *factor, double *sizes)
+{
+  FEM_ADAPT_Refine(*meshID, *qm, *method, *factor, sizes);
+}
+
+
+ void FEM_ADAPT_Coarsen(int meshID, int qm, int method, double factor, 
+        double *sizes) {
+    FEM_Mesh* mesh = FEM_Mesh_lookup(meshID, "FEM_ADAPT_Init");
+    mesh->getfmMM()->getfmAdaptAlgs()->FEM_Coarsen(qm, method, factor, sizes);
+}
+FDECL  void FTN_NAME(FEM_ADAPT_COARSEN,fem_adapt_coarsen)(int* meshID, 
+        int *qm, int *method, double *factor, double *sizes)
+{
+  FEM_ADAPT_Coarsen(*meshID, *qm, *method, *factor, sizes);
+}
+
+
+ void FEM_ADAPT_SetElementSizeField(int meshID, int elem, double size) {
+  FEM_Mesh *meshP = FEM_Mesh_lookup(meshID, "FEM_ADAPT_Init");
+  meshP->elem[0].setMeshSizing(elem, size);
+}
+FDECL  void FTN_NAME(FEM_ADAPT_SETELEMENTSIZEFIELD,fem_adapt_setelementsizefield)(int *meshID, int *elem, double *size)
+{
+  FEM_ADAPT_SetElementSizeField(*meshID, *elem, *size);
+}
+
+
+ void FEM_ADAPT_SetElementsSizeField(int meshID, double *sizes) {
+  FEM_Mesh *meshP = FEM_Mesh_lookup(meshID, "FEM_ADAPT_Init");
+  int numElements = meshP->elem[0].size();
+  for (int i=0; i<numElements; i++) {
+    meshP->elem[0].setMeshSizing(i, sizes[i]);
+  }
+}
+FDECL  void FTN_NAME(FEM_ADAPT_SETELEMENTSSIZEFIELD,fem_adapt_setelementssizefield)(int *meshID, double *sizes)
+{
+  FEM_ADAPT_SetElementsSizeField(*meshID, sizes);
+}
+
+
+void FEM_ADAPT_SetReferenceMesh(int meshID) {
+    FEM_Mesh* mesh = FEM_Mesh_lookup(meshID, "FEM_ADAPT_Init");
+    mesh->getfmMM()->getfmAdaptAlgs()->SetReferenceMesh();
+}
+FDECL void FTN_NAME(FEM_ADAPT_SETREFERENCEMESH, fem_adapt_setreferencemesh)(int* meshID)
+{
+    FEM_ADAPT_SetReferenceMesh(*meshID);
+}
+
+
+void FEM_ADAPT_GradateMesh(int meshID, double smoothness)
+{
+    FEM_Mesh* mesh = FEM_Mesh_lookup(meshID, "FEM_ADAPT_Init");
+    mesh->getfmMM()->getfmAdaptAlgs()->GradateMesh(smoothness);
+}
+FDECL void FTN_NAME(FEM_ADAPT_GRADATEMESH, fem_adapt_gradatemesh)(int* meshID, double* smoothness)
+{
+    FEM_ADAPT_GradateMesh(*meshID, *smoothness);
+}
+
index 86a70554caa3941fd00e5b4ccd92e094f7e8c2e9..29e46aa2229643751736dda477a0b3991f660fd1 100644 (file)
@@ -2,95 +2,39 @@
 #define __CHARM_FEM_ADAPT_IF_H
 
 #include "charm-api.h"
-#include "fem_mesh_modify.h"
 #include "fem_adapt_algs.h"
+#include "fem_mesh_modify.h"
 
 extern void _registerFEMMeshModify(void);
 
-inline void FEM_ADAPT_Init(int meshID) {
-  const int triangleFaces[6] = {0,1,1,2,2,0};
-  FEM_Add_elem2face_tuples(meshID, 0, 2, 3, triangleFaces);
-  FEM_Mesh_create_elem_elem_adjacency(meshID);
-  FEM_Mesh_create_node_elem_adjacency(meshID);
-  FEM_Mesh_create_node_node_adjacency(meshID);
-  _registerFEMMeshModify();
-  FEM_REF_INIT(meshID, 2);  // dim=2
-  FEM_Mesh *meshP = FEM_Mesh_lookup(meshID, "FEM_ADAPT_Init");
-  CtvInitialize(FEM_Adapt_Algs *, _adaptAlgs);
-  CtvAccess(_adaptAlgs) = meshP->getfmMM()->getfmAdaptAlgs();
-  CtvAccess(_adaptAlgs)->FEM_Adapt_Algs_Init(FEM_DATA+0, FEM_BOUNDARY);
-}
-FDECL inline void FTN_NAME(FEM_ADAPT_INIT,fem_adapt_init)(int *meshID)
-{
-  FEM_ADAPT_Init(*meshID);
-}
-
-
-inline void FEM_ADAPT_Refine(int meshID, int qm, int method, double factor,
-        double *sizes) {
-    FEM_Mesh* mesh = FEM_Mesh_lookup(meshID, "FEM_ADAPT_Init");
-    mesh->getfmMM()->getfmAdaptAlgs()->FEM_Refine(qm, method, factor, sizes);
-}
-FDECL inline void FTN_NAME(FEM_ADAPT_REFINE,fem_adapt_refine)(int* meshID, 
-        int *qm, int *method, double *factor, double *sizes)
-{
-  FEM_ADAPT_Refine(*meshID, *qm, *method, *factor, sizes);
-}
-
-
-inline void FEM_ADAPT_Coarsen(int meshID, int qm, int method, double factor, 
-        double *sizes) {
-    FEM_Mesh* mesh = FEM_Mesh_lookup(meshID, "FEM_ADAPT_Init");
-    mesh->getfmMM()->getfmAdaptAlgs()->FEM_Coarsen(qm, method, factor, sizes);
-}
-FDECL inline void FTN_NAME(FEM_ADAPT_COARSEN,fem_adapt_coarsen)(int* meshID, 
-        int *qm, int *method, double *factor, double *sizes)
-{
-  FEM_ADAPT_Coarsen(*meshID, *qm, *method, *factor, sizes);
-}
-
-
-inline void FEM_ADAPT_SetElementSizeField(int meshID, int elem, double size) {
-  FEM_Mesh *meshP = FEM_Mesh_lookup(meshID, "FEM_ADAPT_Init");
-  meshP->elem[0].setMeshSizing(elem, size);
-}
-FDECL inline void FTN_NAME(FEM_ADAPT_SETELEMENTSIZEFIELD,fem_adapt_setelementsizefield)(int *meshID, int *elem, double *size)
-{
-  FEM_ADAPT_SetElementSizeField(*meshID, *elem, *size);
-}
-
-
-inline void FEM_ADAPT_SetElementsSizeField(int meshID, double *sizes) {
-  FEM_Mesh *meshP = FEM_Mesh_lookup(meshID, "FEM_ADAPT_Init");
-  int numElements = meshP->elem[0].size();
-  for (int i=0; i<numElements; i++) {
-    meshP->elem[0].setMeshSizing(i, sizes[i]);
-  }
-}
-FDECL inline void FTN_NAME(FEM_ADAPT_SETELEMENTSSIZEFIELD,fem_adapt_setelementssizefield)(int *meshID, double *sizes)
-{
-  FEM_ADAPT_SetElementsSizeField(*meshID, sizes);
-}
-
-
-inline void FEM_ADAPT_SetReferenceMesh(int meshID) {
-    FEM_Mesh* mesh = FEM_Mesh_lookup(meshID, "FEM_ADAPT_Init");
-    mesh->getfmMM()->getfmAdaptAlgs()->SetReferenceMesh();
-}
-FDECL inline void FTN_NAME(FEM_ADAPT_SETREFERENCEMESH, fem_adapt_setreferencemesh)(int* meshID)
-{
-    FEM_ADAPT_SetReferenceMesh(*meshID);
-}
-
-
-void inline FEM_ADAPT_GradateMesh(int meshID, double smoothness)
-{
-    FEM_Mesh* mesh = FEM_Mesh_lookup(meshID, "FEM_ADAPT_Init");
-    mesh->getfmMM()->getfmAdaptAlgs()->GradateMesh(smoothness);
-}
-FDECL inline void FTN_NAME(FEM_ADAPT_GRADATEMESH, fem_adapt_gradatemesh)(int* meshID, double* smoothness)
-{
-    FEM_ADAPT_GradateMesh(*meshID, *smoothness);
-}
+void FEM_ADAPT_Init(int meshID);
+FDECL void FTN_NAME(FEM_ADAPT_INIT,fem_adapt_init)(int *meshID);
+
+
+void FEM_ADAPT_Refine(int meshID, int qm, int method, double factor, double *sizes);
+FDECL void FTN_NAME(FEM_ADAPT_REFINE,fem_adapt_refine)(int* meshID, 
+        int *qm, int *method, double *factor, double *sizes);
+
+
+void FEM_ADAPT_Coarsen(int meshID, int qm, int method, double factor, 
+        double *sizes);
+FDECL void FTN_NAME(FEM_ADAPT_COARSEN,fem_adapt_coarsen)(int* meshID, 
+        int *qm, int *method, double *factor, double *sizes);
+
+
+void FEM_ADAPT_SetElementSizeField(int meshID, int elem, double size);
+FDECL void FTN_NAME(FEM_ADAPT_SETELEMENTSIZEFIELD,fem_adapt_setelementsizefield)(int *meshID, int *elem, double *size);
+
+
+void FEM_ADAPT_SetElementsSizeField(int meshID, double *sizes);
+FDECL void FTN_NAME(FEM_ADAPT_SETELEMENTSSIZEFIELD,fem_adapt_setelementssizefield)(int *meshID, double *sizes);
+
+
+void FEM_ADAPT_SetReferenceMesh(int meshID);
+FDECL void FTN_NAME(FEM_ADAPT_SETREFERENCEMESH, fem_adapt_setreferencemesh)(int* meshID);
+
+
+void FEM_ADAPT_GradateMesh(int meshID, double smoothness);
+FDECL void FTN_NAME(FEM_ADAPT_GRADATEMESH, fem_adapt_gradatemesh)(int* meshID, double* smoothness);
 
 #endif