458   PetscInt         Nq, Nc, dimC, qNc, N;  in DMFieldEvaluateFV_DS()  local
468   PetscCall(DMGetCoordinateDim(field->dm, &dimC));  in DMFieldEvaluateFV_DS()
492   if (D) PetscCall(DMGetWorkArray(field->dm, N * dimC, mpitype, &qD));  in DMFieldEvaluateFV_DS()
493   if (H) PetscCall(DMGetWorkArray(field->dm, N * dimC * dimC, mpitype, &qH));  in DMFieldEvaluateFV_DS()
538         for (j = 0; j < Nc * dimC; j++) sD[i * Nc * dimC + j] = 0.;  in DMFieldEvaluateFV_DS()
544             for (l = 0; l < dimC; l++) {  in DMFieldEvaluateFV_DS()
545 …for (m = 0; m < dim; m++) pD[l] += geom->invJ[((i * Nq + k) * dimC + m) * dimC + l] * sqD[((i * Nq…  in DMFieldEvaluateFV_DS()
547 …for (l = 0; l < dimC; l++) sD[(i * Nc + j) * dimC + l] += geom->detJ[i * Nq + k] * weights[k] * pD…  in DMFieldEvaluateFV_DS()
550         for (k = 0; k < Nc * dimC; k++) sD[i * Nc * dimC + k] /= vol;  in DMFieldEvaluateFV_DS()
559         for (j = 0; j < Nc * dimC; j++) rD[i * Nc * dimC + j] = 0.;  in DMFieldEvaluateFV_DS()
565             for (l = 0; l < dimC; l++) {  in DMFieldEvaluateFV_DS()
566 …for (m = 0; m < dim; m++) pD[l] += geom->invJ[((i * Nq + k) * dimC + m) * dimC + l] * rqD[((i * Nq…  in DMFieldEvaluateFV_DS()
568 …for (l = 0; l < dimC; l++) rD[(i * Nc + j) * dimC + l] += geom->detJ[i * Nq + k] * weights[k] * pD…  in DMFieldEvaluateFV_DS()
571         for (k = 0; k < Nc * dimC; k++) rD[i * Nc * dimC + k] /= vol;  in DMFieldEvaluateFV_DS()
585         for (j = 0; j < Nc * dimC * dimC; j++) sH[i * Nc * dimC * dimC + j] = 0.;  in DMFieldEvaluateFV_DS()
587           const PetscReal *invJ = &geom->invJ[(i * Nq + k) * dimC * dimC];  in DMFieldEvaluateFV_DS()
596             const PetscScalar *spH = &sqH[((i * Nq + k) * Nc + j) * dimC * dimC];  in DMFieldEvaluateFV_DS()
598             for (l = 0; l < dimC; l++) {  in DMFieldEvaluateFV_DS()
599               for (m = 0; m < dimC; m++) {  in DMFieldEvaluateFV_DS()
601 …  for (r = 0; r < dim; r++) pH[l][m] += invJ[q * dimC + l] * invJ[r * dimC + m] * spH[q * dim + r];  in DMFieldEvaluateFV_DS()
605             for (l = 0; l < dimC; l++) {  in DMFieldEvaluateFV_DS()
606 …for (m = 0; m < dimC; m++) sH[(i * Nc + j) * dimC * dimC + l * dimC + m] += geom->detJ[i * Nq + k]…  in DMFieldEvaluateFV_DS()
610         for (k = 0; k < Nc * dimC * dimC; k++) sH[i * Nc * dimC * dimC + k] /= vol;  in DMFieldEvaluateFV_DS()
619         for (j = 0; j < Nc * dimC * dimC; j++) rH[i * Nc * dimC * dimC + j] = 0.;  in DMFieldEvaluateFV_DS()
621           const PetscReal *invJ = &geom->invJ[(i * Nq + k) * dimC * dimC];  in DMFieldEvaluateFV_DS()
630             const PetscReal *rpH = &rqH[((i * Nq + k) * Nc + j) * dimC * dimC];  in DMFieldEvaluateFV_DS()
632             for (l = 0; l < dimC; l++) {  in DMFieldEvaluateFV_DS()
633               for (m = 0; m < dimC; m++) {  in DMFieldEvaluateFV_DS()
635 …  for (r = 0; r < dim; r++) pH[l][m] += invJ[q * dimC + l] * invJ[r * dimC + m] * rpH[q * dim + r];  in DMFieldEvaluateFV_DS()
639             for (l = 0; l < dimC; l++) {  in DMFieldEvaluateFV_DS()
640 …for (m = 0; m < dimC; m++) rH[(i * Nc + j) * dimC * dimC + l * dimC + m] += geom->detJ[i * Nq + k]…  in DMFieldEvaluateFV_DS()
644         for (k = 0; k < Nc * dimC * dimC; k++) rH[i * Nc * dimC * dimC + k] /= vol;  in DMFieldEvaluateFV_DS()
649   if (D) PetscCall(DMRestoreWorkArray(field->dm, N * dimC, mpitype, &qD));  in DMFieldEvaluateFV_DS()
650   if (H) PetscCall(DMRestoreWorkArray(field->dm, N * dimC * dimC, mpitype, &qH));  in DMFieldEvaluateFV_DS()