3512c690e3039005a9cdf1ea3bb26a204699e9e6
[charm.git] / src / libs / ck-libs / fem / fem_adapt_if.h
1 #ifndef __CHARM_FEM_ADAPT_IF_H
2 #define __CHARM_FEM_ADAPT_IF_H
3
4 #include "charm-api.h"
5 #include "fem_mesh_modify.h"
6 #include "fem_adapt_algs.h"
7
8 extern void _registerFEMMeshModify(void);
9
10 void FEM_ADAPT_Init(int meshID) {
11   _registerFEMMeshModify();
12   FEM_REF_INIT(meshID, 2);  // dim=2
13   FEM_Mesh *meshP = FEM_Mesh_lookup(meshID, "FEM_ADAPT_Init");
14   CtvInitialize(FEM_Adapt_Algs *, _adaptAlgs);
15   CtvAccess(_adaptAlgs) = meshP->getfmMM()->getfmAdaptAlgs();
16   CtvAccess(_adaptAlgs)->FEM_Adapt_Algs_Init(FEM_DATA+0, FEM_BOUNDARY);
17 }
18 FDECL void FTN_NAME(FEM_ADAPT_INIT,fem_adapt_init)(int *meshID)
19 {
20   FEM_ADAPT_Init(*meshID);
21 }
22
23
24 void FEM_ADAPT_Refine(int qm, int method, double factor, double *sizes) {
25   CtvAccess(_adaptAlgs)->FEM_Refine(qm, method, factor, sizes);
26 }
27 FDECL void FTN_NAME(FEM_ADAPT_REFINE,fem_adapt_refine)(int *qm, int *method, double *factor, double *sizes)
28 {
29   FEM_ADAPT_Refine(*qm, *method, *factor, sizes);
30 }
31
32
33 void FEM_ADAPT_Coarsen(int qm, int method, double factor, double *sizes) {
34   CtvAccess(_adaptAlgs)->FEM_Coarsen(qm, method, factor, sizes);
35 }
36 FDECL void FTN_NAME(FEM_ADAPT_COARSEN,fem_adapt_coarsen)(int *qm, int *method, double *factor, double *sizes)
37 {
38   FEM_ADAPT_Coarsen(*qm, *method, *factor, sizes);
39 }
40
41
42 void FEM_ADAPT_SetElementSizeField(int meshID, int elem, double size) {
43   FEM_Mesh *meshP = FEM_Mesh_lookup(meshID, "FEM_ADAPT_Init");
44   meshP->elem[0].setMeshSizing(elem, size);
45 }
46 FDECL void FTN_NAME(FEM_ADAPT_SETELEMENTSIZEFIELD,fem_adapt_setelementsizefield)(int *meshID, int *elem, double *size)
47 {
48   FEM_ADAPT_SetElementSizeField(*meshID, *elem, *size);
49 }
50
51
52 void FEM_ADAPT_SetElementsSizeField(int meshID, double *sizes) {
53   FEM_Mesh *meshP = FEM_Mesh_lookup(meshID, "FEM_ADAPT_Init");
54   int numElements = meshP->elem[0].size();
55   for (int i=0; i<numElements; i++) {
56     meshP->elem[0].setMeshSizing(i, sizes[i]);
57   }
58 }
59 FDECL void FTN_NAME(FEM_ADAPT_SETELEMENTSSIZEFIELD,fem_adapt_setelementssizefield)(int *meshID, double *sizes)
60 {
61   FEM_ADAPT_SetElementsSizeField(*meshID, sizes);
62 }
63
64
65 void FEM_ADAPT_SetReferenceMesh(int meshID) {
66     FEM_Mesh* mesh = FEM_Mesh_lookup(meshID, "FEM_ADAPT_Init");
67
68     // for each element, set its size to its average edge length
69     // TODO: do we need to run this loop for element types other than 0?
70     
71     double avgLength;
72     int width = mesh->elem[0].getConn().width();
73     int* eConn = (int*)malloc(width*sizeof(int));
74     int numElements = mesh->elem[0].size();
75     
76     for (int i=0; i<numElements; ++i, avgLength=0) {
77         mesh->e2n_getAll(i, eConn);
78         
79         for (int j=0; j<width-1; ++j) {
80             avgLength += CtvAccess(_adaptAlgs)->length(eConn[j], eConn[j+1]);
81         }
82         avgLength += CtvAccess(_adaptAlgs)->length(eConn[0], eConn[width-1]);
83         avgLength /= width;
84         mesh->elem[0].setMeshSizing(i, avgLength);      
85     }
86     free(eConn);
87 }
88 FDECL void FTN_NAME(FEM_ADAPT_SETREFERENCEMESH, fem_adapt_setreferencemesh)(int* meshID)
89 {
90     FEM_ADAPT_SetReferenceMesh(*meshID);
91 }
92
93 #endif