proteus  1.8.1
C/C++/Fortran libraries
postprocessing.h
Go to the documentation of this file.
1 #ifndef POSTPROCESSING_H
2 #define POSTPROCESSING_H
3 #include PROTEUS_LAPACK_H
4 
16 /***********************************************************************
17  data structure for node star solves, same as ASM smoother basically
18 ***********************************************************************/
19 typedef struct
20 {
21  int N;
23  double **subdomain_L,
26  PROTEUS_LAPACK_INTEGER** subdomain_pivots;
27  PROTEUS_LAPACK_INTEGER** subdomain_column_pivots;
29 
30 extern void invertLocal(
31  int nSpace,
32  double (*A)[3],
33  double (*AI)[3]
34  );
35 
36 extern void updateSelectedExteriorElementBoundaryFlux(int nExteriorElementBoundaries_global,
37  int nElementBoundaries_element,
38  int nQuadraturePoints_elementBoundary,
39  int nDOF_test_element,
40  int* exteriorElementBoundaries,
41  int* elementBoundaryElements,
42  int* elementBoundaryLocalElementBoundaries,
43  int* skipflag_elementBoundaries,
44  double* flux,
45  double* w_dS,
46  double* residual);
47 extern void updateRT0velocityWithAveragedPotentialP1nc(int nElements_global,
48  int nQuadraturePoints_element,
49  int nSpace,
50  double * detJ,
51  double * quad_a,
52  double * phi,
53  double * gradphi,
54  double * a,
55  double * rt0vdofs);
57  int nElements_global,
58  int nQuadraturePoints_element,
59  int nDOF_test_element,
60  int nElementBoundaries_element,
61  int nQuadraturePoints_elementBoundary,
62  int nSpace,
63  int * nFreeDOF_element,
64  int * freeLocal_element,
65  double * detJ,
66  double * sqrt_det_g,
67  double * n,
68  double * elementBarycenters,
69  double * quad_a,
70  double * quad_f,
71  double * w_dV_r,
72  double * w_dV_m,
73  double * u,
74  double * gradu,
75  double * a,
76  double * f,
77  double * r,
78  double * mt,
79  double * rt0vdofs
80  );
82  int nElements_global,
83  int nQuadraturePoints_element,
84  int nDOF_test_element,
85  int nElementBoundaries_element,
86  int nQuadraturePoints_elementBoundary,
87  int nSpace,
88  int * nFreeDOF_element,
89  int * freeLocal_element,
90  double * detJ,
91  double * sqrt_det_g,
92  double * n,
93  double * elementBarycenters,
94  double * quad_a,
95  double * quad_f,
96  double * w_dV_r,
97  double * u,
98  double * gradu,
99  double * a,
100  double * f,
101  double * r,
102  double * rt0vdofs
103  );
105  int nElements_global,
106  int nQuadraturePoints_element,
107  int nDOF_test_element,
108  int nElementBoundaries_element,
109  int nQuadraturePoints_elementBoundary,
110  int nSpace,
111  int *nFreeDOF_element,
112  int *freeLocal_element,
113  double *detJ,
114  double *sqrt_det_g,
115  double *n,
116  double *elementBarycenters,
117  double *quad_a,
118  double *quad_f,
119  double *w_dV_r,
120  double *w_dV_m,
121  double *u,
122  double *gradu,
123  double *a,
124  double *f,
125  double *r,
126  double *mt,
127  double *rt0vdofs
128  );
130  int nElements_global,
131  int nQuadraturePoints_element,
132  int nDOF_test_element,
133  int nElementBoundaries_element,
134  int nQuadraturePoints_elementBoundary,
135  int nSpace,
136  int *nFreeDOF_element,
137  int *freeLocal_element,
138  double *detJ,
139  double *sqrt_det_g,
140  double *n,
141  double *elementBarycenters,
142  double *quad_a,
143  double *quad_f,
144  double *w_dV_r,
145  double *u,
146  double *gradu,
147  double *a,
148  double *f,
149  double *r,
150  double *rt0vdofs
151  );
152 extern void getElementRT0velocityValues(
153  int nElements_global,
154  int nPoints_element,
155  int nSpace,
156  double *x_element,
157  double *rt0vdofs_element,
158  double *v_element
159  );
161  int nElements_global,
162  int nElementBoundaries_element,
163  int nPoints_elementBoundary,
164  int nSpace,
165  double *x_elementBoundary,
166  double *rt0vdofs_element,
167  double *v_elementBoundary
168  );
170  int nElementBoundaries_global,
171  int nPoints_elementBoundary,
172  int nSpace,
173  int *elementBoundaryElementsArray,
174  double *x_elementBoundary_global,
175  double *rt0vdofs_element,
176  double *v_elementBoundary_global
177  );
179  int nElements_global,
180  int nQuadraturePoints_element,
181  int nElementBoundaries_element,
182  int nQuadraturePoints_elementBoundary,
183  int nSpace,
184  double *uQuadratureWeights_element,
185  double *elementBarycenters,
186  double *aElementQuadratureWeights,
187  double *detJ,
188  double *uQuadratureWeights_elementBoundary,
189  double *x,
190  double *u,
191  double *gradu,
192  double *x_elementBoundary,
193  double *u_elementBoundary,
194  double *n,
195  double *a,
196  double *f,
197  double *r,
198  double *rt0vdofs,
199  double *rt0potential
200  );
202  int nElements_global,
203  int nElementBoundaries_element,
204  int nQuadraturePoints_elementBoundary,
205  int nSpace,
206  double *elementBoundaryQuadratureWeights,
207  double *n,
208  double *v_elementBoundary,
209  double *rt0vdofs_element
210  );
212  int nElements_global,
213  int nElementBoundaries_element,
214  int nQuadraturePoints_elementBoundary,
215  int nDOF_RT0V_element,
216  int *elementBoundaryElementsArray,
217  int *elementBoundariesArray,
218  double *elementBoundaryQuadratureWeights,
219  double *flux_elementBoundary,
220  double *rt0vdofs_element
221  );
223  int nElements_global,
224  int nElementBoundaries_element,
225  int nPoints_element,
226  int nSpace,
227  int nDetVals_element,
228  double *nodeArray,
229  int *elementNodesArray,
230  double *abs_det_J,
231  double *x_element,
232  double *rt0vdofs_element,
233  double *v_element
234  );
236  int nElements_global,
237  int nElementBoundaries_element,
238  int nPoints_elementBoundary,
239  int nSpace,
240  int nDetVals_element,
241  double *nodeArray,
242  int *elementNodesArray,
243  double *abs_det_J,
244  double *x_elementBoundary,
245  double *rt0vdofs_element,
246  double *v_elementBoundary
247  );
249  int nElementBoundaries_global,
250  int nPoints_elementBoundary_global,
251  int nSpace,
252  int nDetVals_element,
253  double *nodeArray,
254  int *elementNodesArray,
255  int *elementBoundaryElementsArray,
256  double *abs_det_J,
257  double *x_elementBoundary_global,
258  double *rt0vdofs_element,
259  double *v_elementBoundary_global
260  );
262  int nElements_global,
263  int nElementBoundaries_element,
264  int nQuadraturePoints_elementBoundary,
265  int nSpace,
266  int nDOFs_test_element,
267  int nDOFs_trial_element,
268  int nVDOFs_element,
269  double *w_dS_f,
270  double *ebq_n,
271  double *ebq_v,
272  double *BDMprojectionMat_element
273  );
275  int nElements_global,
276  int nVDOFs_element,
277  double *BDMprojectionMat_element,
278  int *BDMprojectionMatPivots_element
279  );
280 extern void solveLocalBDM1projection(
281  int nElements_global,
282  int nElementBoundaries_element,
283  int nQuadraturePoints_elementBoundary,
284  int nSpace,
285  int nDOFs_test_element,
286  int nVDOFs_element,
287  double *BDMprojectionMatFact_element,
288  int *BDMprojectionMatPivots_element,
289  double *w_dS_f,
290  double *ebq_n,
291  double *ebq_velocity,
292  double *p1_velocity_dofs
293  );
295  int nElements_global,
296  int nElementBoundaries_element,
297  int nQuadraturePoints_elementBoundary,
298  int nDOFs_test_element,
299  int nVDOFs_element,
300  double *BDMprojectionMatFact_element,
301  int *BDMprojectionMatPivots_element,
302  int *elementBoundaryElementsArray,
303  int *elementBoundariesArray,
304  double *w_dS_f,
305  double *ebq_global_flux,
306  double *p1_velocity_dofs
307  );
309  int nElements_global,
310  int nQuadraturePoints_element,
311  int nSpace,
312  int nDOF_trial_element,
313  int nVDOF_element,
314  double *q_v,
315  double *p1_velocity_dofs,
316  double *q_velocity
317  );
319  int nElements_global,
320  int nInteriorElementBoundaries_global,
321  int nExteriorElementBoundaries_global,
322  int nElementBoundaries_element,
323  int nQuadraturePoints_elementBoundary,
324  int nNodes_element,
325  int nSpace,
326  int *interiorElementBoundaries,
327  int *exteriorElementBoundaries,
328  int *elementBoundaryElements,
329  int *elementBoundaryLocalElementBoundaries,
330  int *exteriorElementBoundariesToSkip,
331  double *dS,
332  double *normal,
333  double *elementResidual,
334  double *velocity,
335  double *conservationResidual
336  );
337 extern void sunWheelerGSsweep(
338  int nElements_global,
339  int nElementBoundaries_global,
340  int nInteriorElementBoundaries_global,
341  int nExteriorElementBoundaries_global,
342  int nElementBoundaries_element,
343  int nQuadraturePoints_elementBoundary,
344  int nSpace,
345  int *interiorElementBoundaries,
346  int *exteriorElementBoundaries,
347  int *elementBoundaryElements,
348  int *elementBoundaryLocalElementBoundaries,
349  double *dS,
350  double *normal,
351  double *sqrt_det_g,
352  double *alpha,
353  double *fluxCorrection,
354  double *conservationResidual
355  );
356 extern void fluxCorrectionVelocityUpdate(
357  int nElements_global,
358  int nElementBoundaries_global,
359  int nInteriorElementBoundaries_global,
360  int nExteriorElementBoundaries_global,
361  int nElementBoundaries_element,
362  int nQuadraturePoints_elementBoundary,
363  int nSpace,
364  int *interiorElementBoundaries,
365  int *exteriorElementBoundaries,
366  int *elementBoundaryElements,
367  int *elementBoundaryLocalElementBoundaries,
368  double *dS,
369  double *normal,
370  double *fluxCorrection,
371  double *vConservative,
372  double *vConservative_element
373  );
374 extern void computeFluxCorrectionPWC(
375  int nElementBoundaries_global,
376  int nInteriorElementBoundaries_global,
377  int nExteriorElementBoundaries_global,
378  int *interiorElementBoundaries,
379  int *exteriorElementBoundaries,
380  int *elementBoundaryElements,
381  double *pwcW,
382  double *pwcV,
383  double *fluxCorrection
384  );
385 extern int nodeStar_init(
386  int nElements_global,
387  int nNodes_element,
388  int nNodes_global,
389  int *nElements_node,
390  int *nodeStarElementsArray,
391  int *nodeStarElementNeighborsArray,
392  int *N_p,
393  int **subdomain_dim_p,
394  double ***subdomain_L_p,
395  double ***subdomain_R_p,
396  double ***subdomain_U_p,
397  PROTEUS_LAPACK_INTEGER *** subdomain_pivots_p,
398  PROTEUS_LAPACK_INTEGER *** subdomain_column_pivots_p
399  );
400 extern int nodeStar_free(
401  int N,
402  int *subdomain_dim,
403  double **subdomain_L,
404  double **subdomain_R,
405  double **subdomain_U,
406  PROTEUS_LAPACK_INTEGER ** subdomain_pivots,
407  PROTEUS_LAPACK_INTEGER ** subdomain_column_pivots_p
408  );
409 extern int nodeStar_setU(
410  NodeStarFactorStruct * nodeStarFactor,
411  double val
412  );
413 extern int nodeStar_copy(
414  int other_N,
415  int *other_subdomain_dim,
416  double **other_subdomain_L,
417  double **other_subdomain_R,
418  double **other_subdomain_U,
419  PROTEUS_LAPACK_INTEGER ** other_subdomain_pivots,
420  PROTEUS_LAPACK_INTEGER ** other_subdomain_column_pivots,
421  int *N_p,
422  int **subdomain_dim_p,
423  double ***subdomain_L_p,
424  double ***subdomain_R_p,
425  double ***subdomain_U_p,
426  PROTEUS_LAPACK_INTEGER *** subdomain_pivots_p,
427  PROTEUS_LAPACK_INTEGER *** subdomain_column_pivots_p
428  );
429 extern
430 void calculateConservationResidualPWL(int nElements_global,
431  int nInteriorElementBoundaries_global,
432  int nExteriorElementBoundaries_global,
433  int nElementBoundaries_element,
434  int nQuadraturePoints_elementBoundary,
435  int nDOF_element,
436  int nNodes_element,
437  int nSpace,
438  int* interiorElementBoundaries,
439  int* exteriorElementBoundaries,
440  int* elementBoundaryElements,
441  int* elementBoundaryLocalElementBoundaries,
442  int* elementNodes,
443  int* dofMapl2g,
444  int* nodeStarElements,
445  int* nodeStarElementNeighbors,
446  int* nElements_node,
447  int* fluxElementBoundaries,
448  double* elementResidual,
449  double* vAverage,
450  double* dX,
451  double* w,
452  double* normal,
453  NodeStarFactorStruct* nodeStarFactor,
454  double* conservationResidual,
455  double* vConservative,
456  double* vConservative_element);
457 
458 
459 extern
460 void calculateConservationJacobianPWL(int nNodes_global,
461  int nNodes_internal,
462  int nElements_global,
463  int nInteriorElementBoundaries_global,
464  int nExteriorElementBoundaries_global,
465  int nElementBoundaries_element,
466  int nQuadraturePoints_elementBoundary,
467  int nNodes_element,
468  int nDOF_element,
469  int nSpace,
470  int* interiorElementBoundaries,
471  int* exteriorElementBoundaries,
472  int* elementBoundaryElements,
473  int* elementBoundaryLocalElementBoundaries,
474  int* elementNodes,
475  int* dofMapl2g,
476  int* nodeStarElements,
477  int* nodeStarElementNeighbors,
478  int* nElements_node,
479  int* internalNodes,
480  int* fluxElementBoundaries,
481  int* fluxBoundaryNodes,
482  double* w,
483  double* normal,
484  NodeStarFactorStruct* nodeStarFactor);
485 
486 extern
487 void calculateConservationFluxPWL(int nNodes_global,
488  int nNodes_internal,
489  int* nElements_node,
490  int* internalNodes,
491  int* fluxBoundaryNodes,
492  NodeStarFactorStruct* nodeStarFactor);
493 extern
494 void calculateConservationResidualPWL_opt(int nNodes_owned,
495  int nElements_global,
496  int nInteriorElementBoundaries_global,
497  int nExteriorElementBoundaries_global,
498  int nElementBoundaries_element,
499  int nQuadraturePoints_elementBoundary,
500  int nNodes_element,
501  int nSpace,
502  int* interiorElementBoundaries,
503  int* exteriorElementBoundaries,
504  int* elementBoundaryElements,
505  int* elementBoundaryLocalElementBoundaries,
506  int* elementNodes,
507  int* nodeStarElements,
508  int* nodeStarElementNeighbors,
509  int* nElements_node,
510  int* fluxElementBoundaries,
511  double* elementResidual,
512  double* vAverage,
513  double* dX,
514  double* w,
515  double* normal,
516  NodeStarFactorStruct* nodeStarFactor,
517  double* conservationResidual,
518  double* vConservative,
519  double* vConservative_element);
520 
521 
522 extern
523 void calculateConservationJacobianPWL_opt(int nNodes_owned,
524  int nNodes_global,
525  int nNodes_internal,
526  int nElements_global,
527  int nInteriorElementBoundaries_global,
528  int nExteriorElementBoundaries_global,
529  int nElementBoundaries_element,
530  int nQuadraturePoints_elementBoundary,
531  int nNodes_element,
532  int nSpace,
533  int* interiorElementBoundaries,
534  int* exteriorElementBoundaries,
535  int* elementBoundaryElements,
536  int* elementBoundaryLocalElementBoundaries,
537  int* elementNodes,
538  int* nodeStarElements,
539  int* nodeStarElementNeighbors,
540  int* nElements_node,
541  int* internalNodes,
542  int* fluxElementBoundaries,
543  int* fluxBoundaryNodes,
544  double* w,
545  double* normal,
546  NodeStarFactorStruct* nodeStarFactor);
547 
548 extern
549 void calculateConservationFluxPWL_opt(int nNodes_owned,
550  int nNodes_global,
551  int nNodes_internal,
552  int* nElements_node,
553  int* internalNodes,
554  int* fluxBoundaryNodes,
555  NodeStarFactorStruct* nodeStarFactor);
556 
558  int nElements_global,
559  int nInteriorElementBoundaries_global,
560  int nExteriorElementBoundaries_global,
561  int nElementBoundaries_element,
562  int nQuadraturePoints_elementBoundary,
563  int nNodes_element,
564  int nSpace,
565  int *interiorElementBoundaries,
566  int *exteriorElementBoundaries,
567  int *elementBoundaryElements,
568  int *elementBoundaryLocalElementBoundaries,
569  int *elementNodes,
570  int *nodeStarElements,
571  int *nodeStarElementNeighbors,
572  int *nElements_node,
573  int *fluxElementBoundaries,
574  double *elementResidual,
575  double *vAverage,
576  double *dX,
577  double *w,
578  double *normal,
579  NodeStarFactorStruct * nodeStarFactor,
580  double *conservationResidual,
581  double *vConservative,
582  double *vConservative_element
583  );
585  int nNodes_global,
586  int nNodes_internal,
587  int nElements_global,
588  int nInteriorElementBoundaries_global,
589  int nExteriorElementBoundaries_global,
590  int nElementBoundaries_element,
591  int nQuadraturePoints_elementBoundary,
592  int nNodes_element,
593  int nSpace,
594  int *interiorElementBoundaries,
595  int *exteriorElementBoundaries,
596  int *elementBoundaryElements,
597  int *elementBoundaryLocalElementBoundaries,
598  int *elementNodes,
599  int *nodeStarElements,
600  int *nodeStarElementNeighbors,
601  int *nElements_node,
602  int *internalNodes,
603  int *fluxElementBoundaries,
604  int *fluxBoundaryNodes,
605  double *w,
606  double *normal,
607  NodeStarFactorStruct * nodeStarFactor
608  );
609 extern void calculateConservationFluxPWLv3(int nNodes_global,
610  int nNodes_internal,
611  int* nElements_node,
612  int* internalNodes,
613  int* fluxBoundaryNodes,
614  NodeStarFactorStruct* nodeStarFactor);
615 
616 extern void postprocessAdvectiveVelocityPointEval(int nPoints,
617  int nSpace,
618  double updateCoef,
619  const double* f,
620  double * velocity);
621 
622 extern void postprocessDiffusiveVelocityPointEval(int nPoints,
623  int nSpace,
624  double updateCoef,
625  const double* a,
626  const double* grad_phi,
627  double * velocity);
628 extern void calculateConservationResidualPWL_primative(int nElements_global,
629  int nInteriorElementBoundaries_global,
630  int nExteriorElementBoundaries_global,
631  int nElementBoundaries_element,
632  int nQuadraturePoints_elementBoundary,
633  int nNodes_element,
634  int nSpace,
635  int* interiorElementBoundaries,
636  int* exteriorElementBoundaries,
637  int* elementBoundaryElements,
638  int* elementBoundaryLocalElementBoundaries,
639  int* skipflag_elementBoundaries,
640  double* elementResidual,
641  double* dX,
642  double* normal,
643  double* conservationResidual,
644  double* vConservative);
645 extern void postProcessRT0velocityFromP1nc_sd(int nElements_global,
646  int nQuadraturePoints_element,
647  int nDOF_test_element,
648  int nElementBoundaries_element,
649  int nQuadraturePoints_elementBoundary,
650  int nSpace,
651  int* rowptr,
652  int* colind,
653  int * nFreeDOF_element,
654  int * freeLocal_element,
655  double * detJ,
656  double * sqrt_det_g,
657  double * n,
658  double * elementBarycenters,
659  double * quad_a,
660  double * quad_f,
661  double * w_dV_r,
662  double * w_dV_m,
663  double * u,
664  double * gradu,
665  double * a,
666  double * f,
667  double * r,
668  double * mt,
669  double * rt0vdofs);
670 extern void postProcessRT0velocityFromP1ncNoMass_sd(int nElements_global,
671  int nQuadraturePoints_element,
672  int nDOF_test_element,
673  int nElementBoundaries_element,
674  int nQuadraturePoints_elementBoundary,
675  int nSpace,
676  int* rowptr,
677  int* colind,
678  int * nFreeDOF_element,
679  int * freeLocal_element,
680  double * detJ,
681  double * sqrt_det_g,
682  double * n,
683  double * elementBarycenters,
684  double * quad_a,
685  double * quad_f,
686  double * w_dV_r,
687  double * u,
688  double * gradu,
689  double * a,
690  double * f,
691  double * r,
692  double * rt0vdofs);
693 extern void updateRT0velocityWithAveragedPotentialP1nc_sd(int nElements_global,
694  int nQuadraturePoints_element,
695  int nSpace,
696  int* rowptr,
697  int* colind,
698  double * detJ,
699  double * quad_a,
700  double * phi,
701  double * gradphi,
702  double * a,
703  double * rt0vdofs);
704 extern void postProcessRT0potentialFromP1nc_sd(int nElements_global,
705  int nQuadraturePoints_element,
706  int nElementBoundaries_element,
707  int nQuadraturePoints_elementBoundary,
708  int nSpace,
709  int* rowptr,
710  int* colind,
711  double * uQuadratureWeights_element,
712  double * elementBarycenters,
713  double * aElementQuadratureWeights,
714  double * detJ,
715  double * uQuadratureWeights_elementBoundary,
716  double * x,
717  double * u,
718  double * gradu,
719  double * x_elementBoundary,
720  double * u_elementBoundary,
721  double * n,
722  double * a,
723  double * f,
724  double * r,
725  double * rt0vdofs,
726  double * rt0potential);
727 extern void postprocessDiffusiveVelocityPointEval_sd(int nPoints,
728  int nSpace,
729  double updateCoef,
730  int* rowptr,
731  int* colind,
732  const double* a,
733  const double* grad_phi,
734  double * velocity);
735 extern void getGlobalExteriorElementBoundaryRT0velocityValues(int nExteriorElementBoundaries_global,
736  int nPoints_elementBoundary,
737  int nSpace,
738  int * elementBoundaryElementsArray,
739  int * exteriorElementBoundariesArray,
740  double * x_elementBoundary_global,
741  double * rt0vdofs_element,
742  double * v_elementBoundary_global);
743 extern void buildLocalBDM2projectionMatrices(int degree,
744  int nElements_global,
745  int nElementBoundaries_element,
746  int nQuadraturePoints_elementBoundary,
747  int nQuadraturePoints_elementInterior,
748  int nSpace,
749  int nDOFs_test_element,
750  int nDOFs_trial_boundary_element,
751  int nDOFs_trial_interior_element,
752  int nVDOFs_element,
753  int *edgeFlags,
754  double * w_dS_f,
755  double * ebq_n,
756  double * ebq_v,
757  double * BDMprojectionMat_element,
758  double * q_basis_vals,
759  double * w_int_test_grads,
760  double * w_int_div_free,
761  double * piola_trial_fun);
762 extern void factorLocalBDM2projectionMatrices(int nElements_global,
763  int nVDOFs_element,
764  double *BDMprojectionMat_element,
765  int *BDMprojectionMatPivots_element);
766 extern void solveLocalBDM2projection(int nElements_global,
767  int nElementBoundaries_element,
768  int nQuadraturePoints_elementBoundary,
769  int nSpace,
770  int nDOFs_test_element,
771  int nVDOFs_element,
772  double * BDMprojectionMatFact_element,
773  int* BDMprojectionMatPivots_element,
774  double * w_dS_f,
775  double * ebq_n,
776  double * w_interior_gradients,
777  double * q_velocity,
778  double * ebq_velocity,
779  double * p1_velocity_dofs);
780 extern void buildBDM2rhs(int nElements_global,
781  int nElementBoundaries_element,
782  int nQuadraturePoints_elementBoundary,
783  int nQuadraturePoints_elementInterior,
784  int nSpace,
785  int nDOFs_test_element,
786  int nVDOFs_element,
787  int nDOFs_trial_interior_element,
788  double * BDMprojectionMatFact_element,
789  int* BDMprojectionMatPivots_element,
790  int *edgeFlags,
791  double * w_dS_f,
792  double * ebq_n,
793  double * w_interior_grads,
794  double * w_interior_divfree,
795  double * ebq_velocity,
796  double * q_velocity,
797  double * p1_velocity_dofs);
798 extern void getElementBDM2velocityValuesLagrangeRep(int nElements_global,
799  int nQuadraturePoints_element,
800  int nSpace,
801  int nDOF_trial_element,
802  int nVDOF_element,
803  double * q_v, /*scalar P^1 shape fncts*/
804  double * p1_velocity_dofs,
805  double * q_velocity);
806 extern void getElementLDGvelocityValuesLagrangeRep(int nElements_global,
807  int nQuadraturePoints_element,
808  int nSpace,
809  int nDOF_trial_element,
810  int nVDOF_element,
811  double * q_v, /*scalar shape fncts*/
812  double * velocity_dofs,
813  double * q_velocity);
814 extern void getGlobalExteriorElementBoundaryBDM1velocityValuesLagrangeRep(int nExteriorElementBoundaries_global,
815  int nQuadraturePoints_elementBoundary,
816  int nSpace,
817  int nDOF_trial_element,
818  int nVDOF_element,
819  int *elementBoundaryElementsArray,
820  int *exteriorElementBoundariesArray,
821  double * ebqe_v, /*scalar P^1 shape fncts*/
822  double * p1_velocity_dofs,
823  double * ebqe_velocity);
824 
825 extern void getElementBoundaryBDM1velocityValuesLagrangeRep(int nElements_global,
826  int nBoundaries_Element,
827  int nQuadraturePoints_elementBoundary,
828  int nSpace,
829  int nDOF_trial_element,
830  int nVDOF_element,
831  int *elementBoundaryElementsArray, /* not used */
832  int *exteriorElementBoundariesArray, /*not used */
833  double * ebq_v, /*scalar P^1 shape fncts*/
834  double * p1_velocity_dofs,
835  double * ebq_velocity);
836 extern void getGlobalExteriorElementBoundaryRT0velocityValuesFluxRep(int nExteriorElementBoundaries_global,
837  int nPoints_elementBoundary_global,
838  int nSpace,
839  int nDetVals_element,
840  double * nodeArray,
841  int *elementNodesArray,
842  int *elementBoundaryElementsArray,
843  int* exteriorElementBoundariesArray,
844  double * abs_det_J,
845  double * x_ebqe,
846  double * rt0vdofs_element,
847  double * v_ebqe);
848 extern void getRT0velocityValuesFluxRep_arbitraryElementMembership(int nElements_global,
849  int nElementBoundaries_element,
850  int nPoints,
851  int nSpace,
852  int nDetVals_element,
853  const double * nodeArray,
854  const int * elementNodesArray,
855  const double * abs_det_J,
856  const double * x,
857  const int * element_locations,
858  const double * rt0vdofs_element,
859  double * v_element);
860 extern void calculateConservationResidualPWL_interiorBoundaries(int nElements_global,
861  int nInteriorElementBoundaries_global,
862  int nExteriorElementBoundaries_global,
863  int nElementBoundaries_element,
864  int nQuadraturePoints_elementBoundary,
865  int nNodes_element,
866  int nSpace,
867  int* interiorElementBoundaries,
868  int* exteriorElementBoundaries,
869  int* elementBoundaryElements,
870  int* elementBoundaryLocalElementBoundaries,
871  int* elementNodes,
872  int* nodeStarElements,
873  int* nodeStarElementNeighbors,
874  int* nElements_node,
875  int* fluxElementBoundaries,
876  double* elementResidual,
877  double* vAverage,
878  double* dX,
879  double* w,
880  double* normal,
881  NodeStarFactorStruct* nodeStarFactor,
882  double* conservationResidual,
883  double* vConservative,
884  double* vConservative_element);
885 extern void calculateConservationJacobianPWL_interiorBoundaries(int nNodes_global,
886  int nElements_global,
887  int nInteriorElementBoundaries_global,
888  int nExteriorElementBoundaries_global,
889  int nElementBoundaries_element,
890  int nQuadraturePoints_elementBoundary,
891  int nNodes_element,
892  int nSpace,
893  int* interiorElementBoundaries,
894  int* exteriorElementBoundaries,
895  int* elementBoundaryElements,
896  int* elementBoundaryLocalElementBoundaries,
897  int* elementNodes,
898  int* nodeStarElements,
899  int* nodeStarElementNeighbors,
900  int* nElements_node,
901  int* fluxElementBoundaries,
902  double* w,
903  double* normal,
904  NodeStarFactorStruct* nodeStarFactor);
905 extern void subdomain_U_copy_global2local(int max_nN_owned,
906  int nElements_global,
907  int nNodes_element,
908  int* elementNodes,
909  int* nodeStarElements,
910  NodeStarFactorStruct* nodeStarFactor,
911  double* subdomain_U);
912 extern void subdomain_U_copy_local2global(int max_nN_owned,
913  int nElements_global,
914  int nNodes_element,
915  int* elementNodes,
916  int* nodeStarElements,
917  NodeStarFactorStruct* nodeStarFactor,
918  double* subdomain_U);
919 
920 extern void calculateElementResidualPWL(int nElements,
921  int nDOF_element_res,
922  int nDOF_element_resPWL,
923  double* alpha,
924  double* elementResidual,
925  double* elementResidualPWL);
927 #endif
sunWheelerGSsweep
void sunWheelerGSsweep(int nElements_global, int nElementBoundaries_global, int nInteriorElementBoundaries_global, int nExteriorElementBoundaries_global, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nSpace, int *interiorElementBoundaries, int *exteriorElementBoundaries, int *elementBoundaryElements, int *elementBoundaryLocalElementBoundaries, double *dS, double *normal, double *sqrt_det_g, double *alpha, double *fluxCorrection, double *conservationResidual)
Definition: postprocessing.c:3746
getGlobalExteriorElementBoundaryRT0velocityValues
void getGlobalExteriorElementBoundaryRT0velocityValues(int nExteriorElementBoundaries_global, int nPoints_elementBoundary, int nSpace, int *elementBoundaryElementsArray, int *exteriorElementBoundariesArray, double *x_elementBoundary_global, double *rt0vdofs_element, double *v_elementBoundary_global)
Definition: postprocessing.c:1636
nodeStar_setU
int nodeStar_setU(NodeStarFactorStruct *nodeStarFactor, double val)
Definition: postprocessing.c:4028
solveLocalBDM1projection
void solveLocalBDM1projection(int nElements_global, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nSpace, int nDOFs_test_element, int nVDOFs_element, double *BDMprojectionMatFact_element, int *BDMprojectionMatPivots_element, double *w_dS_f, double *ebq_n, double *ebq_velocity, double *p1_velocity_dofs)
Definition: postprocessing.c:2964
postProcessRT0potentialFromP1nc
void postProcessRT0potentialFromP1nc(int nElements_global, int nQuadraturePoints_element, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nSpace, double *uQuadratureWeights_element, double *elementBarycenters, double *aElementQuadratureWeights, double *detJ, double *uQuadratureWeights_elementBoundary, double *x, double *u, double *gradu, double *x_elementBoundary, double *u_elementBoundary, double *n, double *a, double *f, double *r, double *rt0vdofs, double *rt0potential)
Definition: postprocessing.c:1682
postProcessRT0velocityFromP1ncV2
void postProcessRT0velocityFromP1ncV2(int nElements_global, int nQuadraturePoints_element, int nDOF_test_element, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nSpace, int *nFreeDOF_element, int *freeLocal_element, double *detJ, double *sqrt_det_g, double *n, double *elementBarycenters, double *quad_a, double *quad_f, double *w_dV_r, double *w_dV_m, double *u, double *gradu, double *a, double *f, double *r, double *mt, double *rt0vdofs)
getElementBoundaryBDM1velocityValuesLagrangeRep
void getElementBoundaryBDM1velocityValuesLagrangeRep(int nElements_global, int nBoundaries_Element, int nQuadraturePoints_elementBoundary, int nSpace, int nDOF_trial_element, int nVDOF_element, int *elementBoundaryElementsArray, int *exteriorElementBoundariesArray, double *ebq_v, double *p1_velocity_dofs, double *ebq_velocity)
Definition: postprocessing.c:3597
w
#define w(x)
Definition: jf.h:22
calculateConservationJacobianPWL_opt
void calculateConservationJacobianPWL_opt(int nNodes_owned, int nNodes_global, int nNodes_internal, int nElements_global, int nInteriorElementBoundaries_global, int nExteriorElementBoundaries_global, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nNodes_element, int nSpace, int *interiorElementBoundaries, int *exteriorElementBoundaries, int *elementBoundaryElements, int *elementBoundaryLocalElementBoundaries, int *elementNodes, int *nodeStarElements, int *nodeStarElementNeighbors, int *nElements_node, int *internalNodes, int *fluxElementBoundaries, int *fluxBoundaryNodes, double *w, double *normal, NodeStarFactorStruct *nodeStarFactor)
Definition: postprocessing.c:5170
updateRT0velocityWithAveragedPotentialP1nc_sd
void updateRT0velocityWithAveragedPotentialP1nc_sd(int nElements_global, int nQuadraturePoints_element, int nSpace, int *rowptr, int *colind, double *detJ, double *quad_a, double *phi, double *gradphi, double *a, double *rt0vdofs)
Definition: postprocessing.c:1409
invertLocal
void invertLocal(int nSpace, double(*A)[3], double(*AI)[3])
postprocessAdvectiveVelocityPointEval
void postprocessAdvectiveVelocityPointEval(int nPoints, int nSpace, double updateCoef, const double *f, double *velocity)
Definition: postprocessing.c:5529
getGlobalElementBoundaryRT0velocityValues
void getGlobalElementBoundaryRT0velocityValues(int nElementBoundaries_global, int nPoints_elementBoundary, int nSpace, int *elementBoundaryElementsArray, double *x_elementBoundary_global, double *rt0vdofs_element, double *v_elementBoundary_global)
Definition: postprocessing.c:1593
getElementBoundaryRT0velocityValuesFluxRep
void getElementBoundaryRT0velocityValuesFluxRep(int nElements_global, int nElementBoundaries_element, int nPoints_elementBoundary, int nSpace, int nDetVals_element, double *nodeArray, int *elementNodesArray, double *abs_det_J, double *x_elementBoundary, double *rt0vdofs_element, double *v_elementBoundary)
Definition: postprocessing.c:2209
factorLocalBDM1projectionMatrices
void factorLocalBDM1projectionMatrices(int nElements_global, int nVDOFs_element, double *BDMprojectionMat_element, int *BDMprojectionMatPivots_element)
Definition: postprocessing.c:2906
subdomain_U_copy_global2local
void subdomain_U_copy_global2local(int max_nN_owned, int nElements_global, int nNodes_element, int *elementNodes, int *nodeStarElements, NodeStarFactorStruct *nodeStarFactor, double *subdomain_U)
Definition: postprocessing.c:5425
f
Double f
Definition: Headers.h:64
nodeStar_copy
int nodeStar_copy(int other_N, int *other_subdomain_dim, double **other_subdomain_L, double **other_subdomain_R, double **other_subdomain_U, PROTEUS_LAPACK_INTEGER **other_subdomain_pivots, PROTEUS_LAPACK_INTEGER **other_subdomain_column_pivots, int *N_p, int **subdomain_dim_p, double ***subdomain_L_p, double ***subdomain_R_p, double ***subdomain_U_p, PROTEUS_LAPACK_INTEGER ***subdomain_pivots_p, PROTEUS_LAPACK_INTEGER ***subdomain_column_pivots_p)
Definition: postprocessing.c:4038
calculateElementResidualPWL
void calculateElementResidualPWL(int nElements, int nDOF_element_res, int nDOF_element_resPWL, double *alpha, double *elementResidual, double *elementResidualPWL)
Definition: postprocessing.c:5616
NodeStarFactorStruct::N
int N
Definition: postprocessing.h:21
NodeStarFactorStruct::subdomain_column_pivots
PROTEUS_LAPACK_INTEGER ** subdomain_column_pivots
Definition: postprocessing.h:27
NodeStarFactorStruct::subdomain_L
double ** subdomain_L
Definition: postprocessing.h:23
subdomain_U_copy_local2global
void subdomain_U_copy_local2global(int max_nN_owned, int nElements_global, int nNodes_element, int *elementNodes, int *nodeStarElements, NodeStarFactorStruct *nodeStarFactor, double *subdomain_U)
Definition: postprocessing.c:5449
postprocessDiffusiveVelocityPointEval
void postprocessDiffusiveVelocityPointEval(int nPoints, int nSpace, double updateCoef, const double *a, const double *grad_phi, double *velocity)
Definition: postprocessing.c:5552
n
Int n
Definition: Headers.h:28
phi
Double phi
Definition: Headers.h:76
postProcessRT0velocityFromP1nc_sd
void postProcessRT0velocityFromP1nc_sd(int nElements_global, int nQuadraturePoints_element, int nDOF_test_element, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nSpace, int *rowptr, int *colind, int *nFreeDOF_element, int *freeLocal_element, double *detJ, double *sqrt_det_g, double *n, double *elementBarycenters, double *quad_a, double *quad_f, double *w_dV_r, double *w_dV_m, double *u, double *gradu, double *a, double *f, double *r, double *mt, double *rt0vdofs)
Definition: postprocessing.c:414
getGlobalExteriorElementBoundaryRT0velocityValuesFluxRep
void getGlobalExteriorElementBoundaryRT0velocityValuesFluxRep(int nExteriorElementBoundaries_global, int nPoints_elementBoundary_global, int nSpace, int nDetVals_element, double *nodeArray, int *elementNodesArray, int *elementBoundaryElementsArray, int *exteriorElementBoundariesArray, double *abs_det_J, double *x_ebqe, double *rt0vdofs_element, double *v_ebqe)
Definition: postprocessing.c:2340
getElementRT0velocityValues
void getElementRT0velocityValues(int nElements_global, int nPoints_element, int nSpace, double *x_element, double *rt0vdofs_element, double *v_element)
Definition: postprocessing.c:1504
NodeStarFactorStruct::subdomain_U
double ** subdomain_U
Definition: postprocessing.h:25
calculateConservationJacobianPWLv3
void calculateConservationJacobianPWLv3(int nNodes_global, int nNodes_internal, int nElements_global, int nInteriorElementBoundaries_global, int nExteriorElementBoundaries_global, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nNodes_element, int nSpace, int *interiorElementBoundaries, int *exteriorElementBoundaries, int *elementBoundaryElements, int *elementBoundaryLocalElementBoundaries, int *elementNodes, int *nodeStarElements, int *nodeStarElementNeighbors, int *nElements_node, int *internalNodes, int *fluxElementBoundaries, int *fluxBoundaryNodes, double *w, double *normal, NodeStarFactorStruct *nodeStarFactor)
postProcessRT0velocityFromP1ncNoMass_sd
void postProcessRT0velocityFromP1ncNoMass_sd(int nElements_global, int nQuadraturePoints_element, int nDOF_test_element, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nSpace, int *rowptr, int *colind, int *nFreeDOF_element, int *freeLocal_element, double *detJ, double *sqrt_det_g, double *n, double *elementBarycenters, double *quad_a, double *quad_f, double *w_dV_r, double *u, double *gradu, double *a, double *f, double *r, double *rt0vdofs)
Definition: postprocessing.c:1046
getElementBDM2velocityValuesLagrangeRep
void getElementBDM2velocityValuesLagrangeRep(int nElements_global, int nQuadraturePoints_element, int nSpace, int nDOF_trial_element, int nVDOF_element, double *q_v, double *p1_velocity_dofs, double *q_velocity)
Definition: postprocessing.c:3427
postProcessRT0velocityFromP1nc
void postProcessRT0velocityFromP1nc(int nElements_global, int nQuadraturePoints_element, int nDOF_test_element, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nSpace, int *nFreeDOF_element, int *freeLocal_element, double *detJ, double *sqrt_det_g, double *n, double *elementBarycenters, double *quad_a, double *quad_f, double *w_dV_r, double *w_dV_m, double *u, double *gradu, double *a, double *f, double *r, double *mt, double *rt0vdofs)
Definition: postprocessing.c:51
getElementBoundaryRT0velocityValues
void getElementBoundaryRT0velocityValues(int nElements_global, int nElementBoundaries_element, int nPoints_elementBoundary, int nSpace, double *x_elementBoundary, double *rt0vdofs_element, double *v_elementBoundary)
Definition: postprocessing.c:1543
calculateConservationResidualPWL_opt
void calculateConservationResidualPWL_opt(int nNodes_owned, int nElements_global, int nInteriorElementBoundaries_global, int nExteriorElementBoundaries_global, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nNodes_element, int nSpace, int *interiorElementBoundaries, int *exteriorElementBoundaries, int *elementBoundaryElements, int *elementBoundaryLocalElementBoundaries, int *elementNodes, int *nodeStarElements, int *nodeStarElementNeighbors, int *nElements_node, int *fluxElementBoundaries, double *elementResidual, double *vAverage, double *dX, double *w, double *normal, NodeStarFactorStruct *nodeStarFactor, double *conservationResidual, double *vConservative, double *vConservative_element)
Definition: postprocessing.c:4834
getGlobalElementBoundaryRT0velocityValuesFluxRep
void getGlobalElementBoundaryRT0velocityValuesFluxRep(int nElementBoundaries_global, int nPoints_elementBoundary_global, int nSpace, int nDetVals_element, double *nodeArray, int *elementNodesArray, int *elementBoundaryElementsArray, double *abs_det_J, double *x_elementBoundary_global, double *rt0vdofs_element, double *v_elementBoundary_global)
Definition: postprocessing.c:2277
updateSelectedExteriorElementBoundaryFlux
void updateSelectedExteriorElementBoundaryFlux(int nExteriorElementBoundaries_global, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nDOF_test_element, int *exteriorElementBoundaries, int *elementBoundaryElements, int *elementBoundaryLocalElementBoundaries, int *skipflag_elementBoundaries, double *flux, double *w_dS, double *residual)
Update the element boundary flux on exterior element boundaries.
Definition: postprocessing.c:5482
calculateConservationResidualPWL_interiorBoundaries
void calculateConservationResidualPWL_interiorBoundaries(int nElements_global, int nInteriorElementBoundaries_global, int nExteriorElementBoundaries_global, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nNodes_element, int nSpace, int *interiorElementBoundaries, int *exteriorElementBoundaries, int *elementBoundaryElements, int *elementBoundaryLocalElementBoundaries, int *elementNodes, int *nodeStarElements, int *nodeStarElementNeighbors, int *nElements_node, int *fluxElementBoundaries, double *elementResidual, double *vAverage, double *dX, double *w, double *normal, NodeStarFactorStruct *nodeStarFactor, double *conservationResidual, double *vConservative, double *vConservative_element)
Definition: postprocessing.c:5649
calculateConservationResidualPWL
void calculateConservationResidualPWL(int nElements_global, int nInteriorElementBoundaries_global, int nExteriorElementBoundaries_global, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nDOF_element, int nNodes_element, int nSpace, int *interiorElementBoundaries, int *exteriorElementBoundaries, int *elementBoundaryElements, int *elementBoundaryLocalElementBoundaries, int *elementNodes, int *dofMapl2g, int *nodeStarElements, int *nodeStarElementNeighbors, int *nElements_node, int *fluxElementBoundaries, double *elementResidual, double *vAverage, double *dX, double *w, double *normal, NodeStarFactorStruct *nodeStarFactor, double *conservationResidual, double *vConservative, double *vConservative_element)
Definition: postprocessing.c:4185
calculateConservationFluxPWL
void calculateConservationFluxPWL(int nNodes_global, int nNodes_internal, int *nElements_node, int *internalNodes, int *fluxBoundaryNodes, NodeStarFactorStruct *nodeStarFactor)
Definition: postprocessing.c:4738
u
Double u
Definition: Headers.h:89
updateRT0velocityWithAveragedPotentialP1nc
void updateRT0velocityWithAveragedPotentialP1nc(int nElements_global, int nQuadraturePoints_element, int nSpace, double *detJ, double *quad_a, double *phi, double *gradphi, double *a, double *rt0vdofs)
Definition: postprocessing.c:1315
getElementBDM1velocityValuesLagrangeRep
void getElementBDM1velocityValuesLagrangeRep(int nElements_global, int nQuadraturePoints_element, int nSpace, int nDOF_trial_element, int nVDOF_element, double *q_v, double *p1_velocity_dofs, double *q_velocity)
Definition: postprocessing.c:3386
calculateConservationResidualPWLv3
void calculateConservationResidualPWLv3(int nElements_global, int nInteriorElementBoundaries_global, int nExteriorElementBoundaries_global, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nNodes_element, int nSpace, int *interiorElementBoundaries, int *exteriorElementBoundaries, int *elementBoundaryElements, int *elementBoundaryLocalElementBoundaries, int *elementNodes, int *nodeStarElements, int *nodeStarElementNeighbors, int *nElements_node, int *fluxElementBoundaries, double *elementResidual, double *vAverage, double *dX, double *w, double *normal, NodeStarFactorStruct *nodeStarFactor, double *conservationResidual, double *vConservative, double *vConservative_element)
calculateConservationResidualPWL_primative
void calculateConservationResidualPWL_primative(int nElements_global, int nInteriorElementBoundaries_global, int nExteriorElementBoundaries_global, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nNodes_element, int nSpace, int *interiorElementBoundaries, int *exteriorElementBoundaries, int *elementBoundaryElements, int *elementBoundaryLocalElementBoundaries, int *skipflag_elementBoundaries, double *elementResidual, double *dX, double *normal, double *conservationResidual, double *vConservative)
Definition: postprocessing.c:5079
postProcessRT0potentialFromP1nc_sd
void postProcessRT0potentialFromP1nc_sd(int nElements_global, int nQuadraturePoints_element, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nSpace, int *rowptr, int *colind, double *uQuadratureWeights_element, double *elementBarycenters, double *aElementQuadratureWeights, double *detJ, double *uQuadratureWeights_elementBoundary, double *x, double *u, double *gradu, double *x_elementBoundary, double *u_elementBoundary, double *n, double *a, double *f, double *r, double *rt0vdofs, double *rt0potential)
Definition: postprocessing.c:1842
postProcessRT0velocityFromP1ncNoMass
void postProcessRT0velocityFromP1ncNoMass(int nElements_global, int nQuadraturePoints_element, int nDOF_test_element, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nSpace, int *nFreeDOF_element, int *freeLocal_element, double *detJ, double *sqrt_det_g, double *n, double *elementBarycenters, double *quad_a, double *quad_f, double *w_dV_r, double *u, double *gradu, double *a, double *f, double *r, double *rt0vdofs)
Definition: postprocessing.c:778
NodeStarFactorStruct::subdomain_pivots
PROTEUS_LAPACK_INTEGER ** subdomain_pivots
Definition: postprocessing.h:26
NodeStarFactorStruct
Definition: postprocessing.h:20
getRT0velocityValuesFluxRep_arbitraryElementMembership
void getRT0velocityValuesFluxRep_arbitraryElementMembership(int nElements_global, int nElementBoundaries_element, int nPoints, int nSpace, int nDetVals_element, const double *nodeArray, const int *elementNodesArray, const double *abs_det_J, const double *x, const int *element_locations, const double *rt0vdofs_element, double *v_element)
Definition: postprocessing.c:2405
calculateConservationFluxPWL_opt
void calculateConservationFluxPWL_opt(int nNodes_owned, int nNodes_global, int nNodes_internal, int *nElements_node, int *internalNodes, int *fluxBoundaryNodes, NodeStarFactorStruct *nodeStarFactor)
Definition: postprocessing.c:5361
getGlobalExteriorElementBoundaryBDM1velocityValuesLagrangeRep
void getGlobalExteriorElementBoundaryBDM1velocityValuesLagrangeRep(int nExteriorElementBoundaries_global, int nQuadraturePoints_elementBoundary, int nSpace, int nDOF_trial_element, int nVDOF_element, int *elementBoundaryElementsArray, int *exteriorElementBoundariesArray, double *ebqe_v, double *p1_velocity_dofs, double *ebqe_velocity)
Definition: postprocessing.c:3506
NodeStarFactorStruct::subdomain_R
double ** subdomain_R
Definition: postprocessing.h:24
postprocessDiffusiveVelocityPointEval_sd
void postprocessDiffusiveVelocityPointEval_sd(int nPoints, int nSpace, double updateCoef, int *rowptr, int *colind, const double *a, const double *grad_phi, double *velocity)
Definition: postprocessing.c:5583
buildLocalBDM2projectionMatrices
void buildLocalBDM2projectionMatrices(int degree, int nElements_global, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nQuadraturePoints_elementInterior, int nSpace, int nDOFs_test_element, int nDOFs_trial_boundary_element, int nDOFs_trial_interior_element, int nVDOFs_element, int *edgeFlags, double *w_dS_f, double *ebq_n, double *ebq_v, double *BDMprojectionMat_element, double *q_basis_vals, double *w_int_test_grads, double *w_int_div_free, double *piola_trial_fun)
Definition: postprocessing.c:2709
solveLocalBDM1projectionFromFlux
void solveLocalBDM1projectionFromFlux(int nElements_global, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nDOFs_test_element, int nVDOFs_element, double *BDMprojectionMatFact_element, int *BDMprojectionMatPivots_element, int *elementBoundaryElementsArray, int *elementBoundariesArray, double *w_dS_f, double *ebq_global_flux, double *p1_velocity_dofs)
Definition: postprocessing.c:3261
calculateConservationResidualGlobalBoundaries
void calculateConservationResidualGlobalBoundaries(int nElements_global, int nInteriorElementBoundaries_global, int nExteriorElementBoundaries_global, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nNodes_element, int nSpace, int *interiorElementBoundaries, int *exteriorElementBoundaries, int *elementBoundaryElements, int *elementBoundaryLocalElementBoundaries, int *exteriorElementBoundariesToSkip, double *dS, double *normal, double *elementResidual, double *velocity, double *conservationResidual)
Definition: postprocessing.c:3651
calculateConservationFluxPWLv3
void calculateConservationFluxPWLv3(int nNodes_global, int nNodes_internal, int *nElements_node, int *internalNodes, int *fluxBoundaryNodes, NodeStarFactorStruct *nodeStarFactor)
computeFluxCorrectionPWC
void computeFluxCorrectionPWC(int nElementBoundaries_global, int nInteriorElementBoundaries_global, int nExteriorElementBoundaries_global, int *interiorElementBoundaries, int *exteriorElementBoundaries, int *elementBoundaryElements, double *pwcW, double *pwcV, double *fluxCorrection)
Definition: postprocessing.c:3877
projectElementBoundaryVelocityToRT0fluxRep
void projectElementBoundaryVelocityToRT0fluxRep(int nElements_global, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nSpace, double *elementBoundaryQuadratureWeights, double *n, double *v_elementBoundary, double *rt0vdofs_element)
Definition: postprocessing.c:2003
r
Double r
Definition: Headers.h:83
calculateConservationJacobianPWL
void calculateConservationJacobianPWL(int nNodes_global, int nNodes_internal, int nElements_global, int nInteriorElementBoundaries_global, int nExteriorElementBoundaries_global, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nNodes_element, int nDOF_element, int nSpace, int *interiorElementBoundaries, int *exteriorElementBoundaries, int *elementBoundaryElements, int *elementBoundaryLocalElementBoundaries, int *elementNodes, int *dofMapl2g, int *nodeStarElements, int *nodeStarElementNeighbors, int *nElements_node, int *internalNodes, int *fluxElementBoundaries, int *fluxBoundaryNodes, double *w, double *normal, NodeStarFactorStruct *nodeStarFactor)
Definition: postprocessing.c:4547
NodeStarFactorStruct::subdomain_dim
int * subdomain_dim
Definition: postprocessing.h:22
buildLocalBDM1projectionMatrices
void buildLocalBDM1projectionMatrices(int nElements_global, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nSpace, int nDOFs_test_element, int nDOFs_trial_element, int nVDOFs_element, double *w_dS_f, double *ebq_n, double *ebq_v, double *BDMprojectionMat_element)
Definition: postprocessing.c:2571
nodeStar_free
int nodeStar_free(int N, int *subdomain_dim, double **subdomain_L, double **subdomain_R, double **subdomain_U, PROTEUS_LAPACK_INTEGER **subdomain_pivots, PROTEUS_LAPACK_INTEGER **subdomain_column_pivots_p)
Definition: postprocessing.c:3996
factorLocalBDM2projectionMatrices
void factorLocalBDM2projectionMatrices(int nElements_global, int nVDOFs_element, double *BDMprojectionMat_element, int *BDMprojectionMatPivots_element)
Definition: postprocessing.c:2935
calculateConservationJacobianPWL_interiorBoundaries
void calculateConservationJacobianPWL_interiorBoundaries(int nNodes_global, int nElements_global, int nInteriorElementBoundaries_global, int nExteriorElementBoundaries_global, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nNodes_element, int nSpace, int *interiorElementBoundaries, int *exteriorElementBoundaries, int *elementBoundaryElements, int *elementBoundaryLocalElementBoundaries, int *elementNodes, int *nodeStarElements, int *nodeStarElementNeighbors, int *nElements_node, int *fluxElementBoundaries, double *w, double *normal, NodeStarFactorStruct *nodeStarFactor)
Definition: postprocessing.c:6005
postProcessRT0velocityFromP1ncV2noMass
void postProcessRT0velocityFromP1ncV2noMass(int nElements_global, int nQuadraturePoints_element, int nDOF_test_element, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nSpace, int *nFreeDOF_element, int *freeLocal_element, double *detJ, double *sqrt_det_g, double *n, double *elementBarycenters, double *quad_a, double *quad_f, double *w_dV_r, double *u, double *gradu, double *a, double *f, double *r, double *rt0vdofs)
buildBDM2rhs
void buildBDM2rhs(int nElements_global, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nQuadraturePoints_elementInterior, int nSpace, int nDOFs_test_element, int nVDOFs_element, int nDOFs_trial_interior_element, double *BDMprojectionMatFact_element, int *BDMprojectionMatPivots_element, int *edgeFlags, double *w_dS_f, double *ebq_n, double *w_interior_grads, double *w_interior_divfree, double *ebq_velocity, double *q_velocity, double *p1_velocity_dofs)
Definition: postprocessing.c:3054
nodeStar_init
int nodeStar_init(int nElements_global, int nNodes_element, int nNodes_global, int *nElements_node, int *nodeStarElementsArray, int *nodeStarElementNeighborsArray, int *N_p, int **subdomain_dim_p, double ***subdomain_L_p, double ***subdomain_R_p, double ***subdomain_U_p, PROTEUS_LAPACK_INTEGER ***subdomain_pivots_p, PROTEUS_LAPACK_INTEGER ***subdomain_column_pivots_p)
Definition: postprocessing.c:3921
solveLocalBDM2projection
void solveLocalBDM2projection(int nElements_global, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nSpace, int nDOFs_test_element, int nVDOFs_element, double *BDMprojectionMatFact_element, int *BDMprojectionMatPivots_element, double *w_dS_f, double *ebq_n, double *w_interior_gradients, double *q_velocity, double *ebq_velocity, double *p1_velocity_dofs)
Definition: postprocessing.c:3213
fluxCorrectionVelocityUpdate
void fluxCorrectionVelocityUpdate(int nElements_global, int nElementBoundaries_global, int nInteriorElementBoundaries_global, int nExteriorElementBoundaries_global, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nSpace, int *interiorElementBoundaries, int *exteriorElementBoundaries, int *elementBoundaryElements, int *elementBoundaryLocalElementBoundaries, double *dS, double *normal, double *fluxCorrection, double *vConservative, double *vConservative_element)
Definition: postprocessing.c:3795
getElementRT0velocityValuesFluxRep
void getElementRT0velocityValuesFluxRep(int nElements_global, int nElementBoundaries_element, int nPoints_element, int nSpace, int nDetVals_element, double *nodeArray, int *elementNodesArray, double *abs_det_J, double *x_element, double *rt0vdofs_element, double *v_element)
Definition: postprocessing.c:2143
getElementLDGvelocityValuesLagrangeRep
void getElementLDGvelocityValuesLagrangeRep(int nElements_global, int nQuadraturePoints_element, int nSpace, int nDOF_trial_element, int nVDOF_element, double *q_v, double *velocity_dofs, double *q_velocity)
Definition: postprocessing.c:3468
projectElementBoundaryFluxToRT0fluxRep
void projectElementBoundaryFluxToRT0fluxRep(int nElements_global, int nElementBoundaries_element, int nQuadraturePoints_elementBoundary, int nDOF_RT0V_element, int *elementBoundaryElementsArray, int *elementBoundariesArray, double *elementBoundaryQuadratureWeights, double *flux_elementBoundary, double *rt0vdofs_element)
Definition: postprocessing.c:2084