diff --git a/src/ref_recon.c b/src/ref_recon.c index 6dee97345..09c8f91fc 100644 --- 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;