NGSolve
4.9
|
00001 #ifndef FILE_POSTPROC 00002 #define FILE_POSTPROC 00003 00004 /*********************************************************************/ 00005 /* File: postproc.hpp */ 00006 /* Author: Joachim Schoeberl */ 00007 /* Date: 25. Mar. 2000 */ 00008 /*********************************************************************/ 00009 00010 namespace ngcomp 00011 { 00012 00013 /* 00014 Postprocessing functions 00015 */ 00016 00017 template <class SCAL> 00018 extern NGS_DLL_HEADER void CalcFlux (const MeshAccess & ma, 00019 const S_GridFunction<SCAL> & u, 00020 S_GridFunction<SCAL> & flux, 00021 const BilinearFormIntegrator & bli, 00022 bool applyd, bool add, 00023 int domain); 00024 00025 extern NGS_DLL_HEADER 00026 void CalcFluxProject (const MeshAccess & ma, 00027 const GridFunction & u, 00028 GridFunction & flux, 00029 const BilinearFormIntegrator & bli, 00030 bool applyd, int domain, 00031 LocalHeap & lh); 00032 00033 template <class SCAL> 00034 extern NGS_DLL_HEADER void CalcFluxProject (const MeshAccess & ma, 00035 const S_GridFunction<SCAL> & u, 00036 S_GridFunction<SCAL> & flux, 00037 const BilinearFormIntegrator & bli, 00038 bool applyd, const BitArray & domains, LocalHeap & lh); 00039 00040 00041 extern NGS_DLL_HEADER 00042 void SetValues (const MeshAccess & ma, 00043 const CoefficientFunction & coef, 00044 GridFunction & u, 00045 bool bound, 00046 DifferentialOperator * diffop, // NULL is FESpace evaluator 00047 LocalHeap & clh); 00048 00049 00050 00051 template <class SCAL> 00052 extern NGS_DLL_HEADER 00053 int CalcPointFlux (const MeshAccess & ma, 00054 const GridFunction & u, 00055 const FlatVector<double> & point, 00056 FlatVector<SCAL> & flux, 00057 const BilinearFormIntegrator & bli, 00058 bool applyd, 00059 LocalHeap & lh, 00060 int component = 0); 00061 00062 00063 template <class SCAL> 00064 extern NGS_DLL_HEADER 00065 int CalcPointFlux (const MeshAccess & ma, 00066 const GridFunction & u, 00067 const FlatVector<double> & point, 00068 const Array<int> & domains, 00069 FlatVector<SCAL> & flux, 00070 const BilinearFormIntegrator & bli, 00071 bool applyd, 00072 LocalHeap & lh, 00073 int component = 0); 00074 00075 00076 extern NGS_DLL_HEADER 00077 void CalcError (const MeshAccess & ma, 00078 const GridFunction & bu, 00079 const GridFunction & bflux, 00080 const BilinearFormIntegrator & bli, 00081 FlatVector<double> & err, 00082 int domain, 00083 LocalHeap & lh); 00084 00085 template <class SCAL> 00086 extern NGS_DLL_HEADER void CalcError (const MeshAccess & ma, 00087 const S_GridFunction<SCAL> & u, 00088 const S_GridFunction<SCAL> & flux, 00089 const BilinearFormIntegrator & bli, 00090 FlatVector<double> & err, 00091 const BitArray & domains, LocalHeap & lh); 00092 00093 00094 template <class SCAL> 00095 NGS_DLL_HEADER void CalcDifference (const MeshAccess & ma, 00096 const S_GridFunction<SCAL> & u1, 00097 const S_GridFunction<SCAL> & u2, 00098 const BilinearFormIntegrator & bli1, 00099 const BilinearFormIntegrator & bli2, 00100 FlatVector<double> & diff, 00101 int domain, LocalHeap & lh); 00102 00103 NGS_DLL_HEADER void CalcDifference (const MeshAccess & ma, 00104 const GridFunction & u1, 00105 const BilinearFormIntegrator & bli1, 00106 const CoefficientFunction * coef, 00107 FlatVector<double> & diff, 00108 int domain, LocalHeap & lh); 00109 00110 00111 00112 template <class SCAL> 00113 extern NGS_DLL_HEADER void CalcGradient (const MeshAccess & ma, 00114 const FESpace & fesh1, 00115 const S_BaseVector<SCAL> & vech1, 00116 const FESpace & feshcurl, 00117 S_BaseVector<SCAL> & vechcurl); 00118 00119 template <class SCAL> 00120 extern NGS_DLL_HEADER void CalcGradientT (const MeshAccess & ma, 00121 const FESpace & feshcurl, 00122 const S_BaseVector<SCAL> & vechcurl, 00123 const FESpace & fesh1, 00124 S_BaseVector<SCAL> & vech1); 00125 00126 00127 template <class SCAL> 00128 extern NGS_DLL_HEADER void CalcErrorHierarchical (const MeshAccess & ma, 00129 const S_BilinearForm<SCAL> & bfa, 00130 const S_BilinearForm<SCAL> & bfa2, 00131 const S_LinearForm<SCAL> & lff, 00132 S_GridFunction<SCAL> & gfu, 00133 const FESpace & festest, 00134 FlatVector<double> & err, 00135 LocalHeap & lh); 00136 00137 00138 /* 00140 extern void CalcFlux (const MeshAccess & ma, const FESpace & fes, 00141 const BDBIntegrator<> & bli, 00142 const BaseVector & u, BaseVector & flux, int applyd = 1); 00143 00145 00146 00147 extern void CalcFluxNodal (const MeshAccess & ma, 00148 const FESpace & fespaceu, 00149 const FESpace & fespaceflux, 00150 const BDBIntegrator<> & bli, 00151 const BaseVector & u, BaseVector & flux, int applyd = 1, 00152 int dom = 0); 00153 00155 00156 00157 extern void CalcFluxElement (const MeshAccess & ma, 00158 const FESpace & fespaceu, 00159 const ElementFESpace & fespaceflux, 00160 const BDBIntegrator<> & bli, 00161 const BaseVector & u, BaseVector & flux, int applyd = 1, 00162 int dom = 0); 00163 00164 00166 00167 00168 extern void CalcBoundaryFlux (const MeshAccess & ma, 00169 const FESpace & fespaceu, 00170 const BDBBoundaryIntegrator<> & bli, 00171 const BaseVector & u, BaseVector & flux, int applyd = 1); 00172 00174 00175 00176 extern void CalcBoundaryFluxNodal (const MeshAccess & ma, 00177 const FESpace & fespaceu, 00178 const FESpace & fespaceflux, 00179 const BDBBoundaryIntegrator<> & bli, 00180 const BaseVector & u, BaseVector & flux, int applyd = 1, 00181 int dom = 0); 00182 00183 00185 00186 00187 extern void ZZErrorEstimator2 (const MeshAccess & ma, 00188 const FESpace & fespaceu, 00189 const FESpace & fespaceflux, 00190 const BDBIntegrator<> & bli, 00191 const BaseVector & u, const BaseVector & flux, 00192 Vector & elerr, int dom = 0); 00193 00195 00196 00197 extern void ZZBoundaryErrorEstimator2 (const MeshAccess & ma, 00198 const FESpace & fespaceu, 00199 const FESpace & fespaceflux, 00200 const BDBBoundaryIntegrator<> & bli, 00201 const BaseVector & u, const BaseVector & flux, 00202 Vector & elerr, int dom = 0); 00203 00204 00205 00206 00208 00209 extern void AverageElementData (const MeshAccess & ma, 00210 const BaseSystemVector & eldata, 00211 BaseSystemVector & nodaldata, 00212 int dom = 0); 00213 00214 00215 00217 00218 extern void ZZErrorEstimate (const MeshAccess & ma, 00219 const BaseSystemVector & eldata, 00220 const BaseSystemVector & nodaldata, 00221 Vector & elerr, 00222 int dom = 0); 00223 00225 00226 extern void Interpolate2Nodal (const MeshAccess & ma, 00227 const FESpace & fespace, 00228 const FESpace & fespacenodal, 00229 const BaseSystemVector & data, 00230 BaseSystemVector & nodaldata, 00231 int dom = 0); 00232 00234 00235 extern void ZZErrorEstimate (const MeshAccess & ma, 00236 const FESpace & fespace, 00237 const FESpace & fespacenodal, 00238 const BaseSystemVector & eldata, 00239 const BaseSystemVector & nodaldata, 00240 Vector & elerr, 00241 int dom = 0); 00242 00244 00245 extern int PointEvaluation (const MeshAccess & ma, 00246 const FESpace & fespace, 00247 const BaseVector & vec, 00248 BDBIntegrator<> & bli, 00249 double * point, 00250 Vector & result, 00251 int applyd = 0); 00252 */ 00253 00254 } 00255 00256 #endif 00257