adds pyramids to l2 projection

nasa · Apr 22, 2021 · e469c46 · e469c46
1 parent d017e02
commit e469c46
Showing 1 changed file with 35 additions and 0 deletions.
diff --git a/src/ref_recon.c b/src/ref_recon.c
@@ -145,6 +145,41 @@ REF_STATUS ref_recon_l2_projection_grad(REF_GRID ref_grid, REF_DBL *scalar,
           case REF_CELL_TR3:
           case REF_CELL_QUA:
           case REF_CELL_PYR:
+            tet_nodes[0] = nodes[0];
+            tet_nodes[1] = nodes[4];
+            tet_nodes[2] = nodes[1];
+            tet_nodes[3] = nodes[2];
+            vol_status = ref_node_tet_vol(ref_node, tet_nodes, &cell_vol);
+            grad_status =
+                ref_node_tet_grad_nodes(ref_node, tet_nodes, scalar, cell_grad);
+            if (REF_SUCCESS == vol_status && REF_SUCCESS == grad_status) {
+              for (cell_node = 0; cell_node < 4; cell_node++)
+                for (i = 0; i < 3; i++)
+                  grad[i + 3 * tet_nodes[cell_node]] += cell_vol * cell_grad[i];
+              for (cell_node = 0; cell_node < 4; cell_node++)
+                vol[tet_nodes[cell_node]] += cell_vol;
+            } else {
+              printf("%s: %d: %s: vol status %d grad status %d\n", __FILE__,
+                     __LINE__, __func__, vol_status, grad_status);
+            }
+            tet_nodes[0] = nodes[0];
+            tet_nodes[1] = nodes[3];
+            tet_nodes[2] = nodes[4];
+            tet_nodes[3] = nodes[2];
+            vol_status = ref_node_tet_vol(ref_node, tet_nodes, &cell_vol);
+            grad_status =
+                ref_node_tet_grad_nodes(ref_node, tet_nodes, scalar, cell_grad);
+            if (REF_SUCCESS == vol_status && REF_SUCCESS == grad_status) {
+              for (cell_node = 0; cell_node < 4; cell_node++)
+                for (i = 0; i < 3; i++)
+                  grad[i + 3 * tet_nodes[cell_node]] += cell_vol * cell_grad[i];
+              for (cell_node = 0; cell_node < 4; cell_node++)
+                vol[tet_nodes[cell_node]] += cell_vol;
+            } else {
+              printf("%s: %d: %s: vol status %d grad status %d\n", __FILE__,
+                     __LINE__, __func__, vol_status, grad_status);
+            }
+            break;
           case REF_CELL_HEX:
             RSS(REF_IMPLEMENT, "implement cell type");
             break;