doc/doxygen/interptheory.dox

   1 /*!
   2 \page InterpKerRemapGlobal Linear remapping
   3
   4 For fields with polynomial representation on each cell, the components of the discretized field  \f$ \phi_s \f$ on the source side can be expressed as linear combinations of the components of the discretized field \f$ \phi_t \f$ on the target side, in terms of a matrix-vector product:
   5
   6 \f[
   7  \phi_t=W.\phi_s.
   8 \f]
   9
  10 \f$W\f$ is called the \anchor interpolationmatrix interpolation matrix.
  11 The objective of interpolators is to compute the matrix W depending on their physical
  12 properties (\ref IntExtFields) and their mesh discretisation (on cells P0, on nodes P1,...).
  13
  14 \section ConsInterp Conservative interpolation
  15
  16 At the basis of many CFD numerical schemes is the fact that physical
  17 quantities such as density, momentum per unit volume or energy per
  18 unit volume obey some balance laws that should be preserved at the
  19 discrete level on every cell.
  20
  21 It is therefore often desired that the process interpolation preserve the
  22 integral of \f$ \phi \f$ on any domain. At the discrete level, for any
  23 target cell \f$ T_i \f$, the following \b general \b interpolation \b
  24 formula \anchor InterpKerGenralEq has to
  25 be satisfied :
  26
  27 \f[
  28 \int_{T_i} \phi_t = \sum_{S_j\cap T_i \neq \emptyset} \int_{T_i\cap S_j} \phi_s.
  29 \f]
  30
  31 This equation is used to compute \f$ W_{ij} \f$, based on the fields representation ( P0, P1, P1d etc..) and the
  32 geometry of source and target mesh cells.
  33
  34 \section MeshOverlap Mesh overlapping
  35
  36 Another important property of the interpolation process is the maximum principle: the field values resulting from the interpolation should remain between the upper and lower bounds of the original field.
  37 When interpolation is performed between a source mesh S and a target
  38 mesh T the aspect of overlapping is important. In fact if any cell of
  39 of S is fully overlapped by cells of T and inversely any cell of T is
  40 fully overlapped by cells of S that is
  41 \f[
  42 \sum_{S_j} Vol(T_i\cap S_j) = Vol(T_i),\hspace{1cm} and \hspace{1cm} \sum_{T_i} Vol(S_j\cap T_i) = Vol(S_j)
  43 \f]
  44 then the meshes S and T are said to be \b
  45 overlapping. In this case the two formulas in a given column in the table below give the same
  46 result. All intensive formulas result in the same output, and all the extensive formulas give also the same output.
  47
  48 The ideal interpolation algorithm should be conservative and respect the maximum principle. However such an algorithm can be impossible to design if the two meshes do not overlap. When the meshes do not overlap, using either \f$Vol(T_i)\f$ or \f$\sum_{S_j} Vol(T_i\cap S_j)\f$ one obtains an algorithm that respects either the conservativity or the maximum principle (see the nature of field \ref TableNatureOfField "summary table").
  49
  50
  51 \section InterpKerRemapInt Linear conservative remapping of P0 (cell based) fields
  52
  53 We assume that the field is represented by a vector with a discrete value on each cell.
  54 This value can represent either
  55 - an average value of the field in the cell (average density, velocity or temperature in the cell) in which case the representation is said to be \b intensive,
  56 - an integrated value over the cell (total mass, power of the cell) in which case the representation is said to be \b extensive
  57
  58 \section InterpKerP0P0Int cell-cell (P0->P0) conservative remapping of intensive fields
  59
  60 For intensive fields such as mass density or power density, the
  61 left hand side in the \ref InterpKerGenralEq "general interpolation equation" becomes :
  62
  63 \f[
  64 \int_{T_i} \phi = Vol(T_i).\phi_{T_i}.
  65 \f]
  66
  67 Here Vol represents the volume when the mesh dimension is equal to 3, the
  68 area when mesh dimension is equal to 2, and length when mesh dimension is equal to 1.
  69
  70 In the \ref InterpKerGenralEq "general interpolation equation" the
  71 right hand side becomes :
  72
  73 \f[
  74 \sum_{S_j\cap T_i \neq \emptyset} \int_{T_i\cap S_j} \phi = \sum_{S_j\cap T_i \neq \emptyset} {Vol(T_i\cap S_j)}.\phi_{S_j}.
  75 \f]
  76
  77 As the field values are constant on each
  78 cell, the coefficients of the linear remapping matrix \f$ W \f$ are
  79 given by the formula :
  80
  81 \f[
  82  W_{ij}=\frac{Vol(T_i\cap S_j)}{ Vol(T_i) }.
  83 \f]
  84
  85
  86 \section InterpKerP0P0Ext cell-cell (P0->P0) conservative remapping of extensive fields
  87
  88 In code coupling from neutronics to hydraulics, \b extensive field
  89 of power is exchanged and the total power should remain the same.
  90 The discrete values of the field represent the total power contained in the cell.
  91 Hence in the \ref InterpKerGenralEq "general interpolation equation" the
  92 left hand side becomes :
  93
  94 \f[
  95 \int_{T_i} \phi = P_{T_i},
  96 \f]
  97
  98 while the right hand side is now :
  99
 100 \f[
 101 \sum_{S_j\cap T_i \neq \emptyset} \int_{T_i\cap S_j} \phi =
 102 \sum_{S_j\cap T_i \neq \emptyset} \frac{Vol(T_i\cap S_j)}{ Vol(S_j)}.P_{S_j}.
 103 \f]
 104
 105 The coefficients of the linear remapping matrix \f$ W \f$ are then
 106 given by the formula :
 107
 108 \f[
 109  W_{ij}=\frac{Vol(T_i\cap S_j)}{  Vol(S_j) }.
 110 \f]
 111
 112 \section TableNatureOfField Summary
 113 In the case of fields with a P0 representation (cell based) and when the meshes do not overlap, the scheme is either conservative or maximum preserving (not both). Depending on the \ref NatureOfField the interpolation coefficients take the following value:
 114
 115  * <TABLE BORDER=1 >
 116  * <TR><TD> </TD><TD>Intensive</TD><TD> extensive </TD></TR>
 117  * <TR><TD> Conservation</TD><TD> \f[\frac{Vol(T_i\cap S_j)}{ Vol(T_i)}\f] <br /> \ref TableNatureOfFieldExampleRevIntegral "RevIntegral" </TD><TD> \f[ \frac{Vol(T_i\cap S_j)}{ \sum_{T_i} Vol(S_j\cap T_i) }\f] <br /> \ref TableNatureOfFieldExampleIntegralGlobConstraint "IntegralGlobConstraint" </TD></TR>
 118  * <TR><TD> Maximum principle </TD><TD> \f[\frac{Vol(T_i\cap S_j)}{ \sum_{S_j} Vol(T_i\cap S_j)}\f] <br /> \ref TableNatureOfFieldExampleConservVol "ConservativeVolumic" </TD><TD>  \f[\frac{Vol(T_i\cap S_j)}{  Vol(S_j) }\f] <br /> \ref TableNatureOfFieldExampleIntegral "Integral"</TD></TR>
 119  *</TABLE>
 120
 121 \section TableNatureOfFieldExample Illustration of a non overlapping P0P0 interpolation
 122
 123 Let's consider the following case with a source mesh containing two cells and a target mesh containing one cell.
 124 Let's consider a field FS on cells on the source mesh that we want to interpolate on the target mesh.
 125
 126 The value of FS on the cell#0 is 4 and the value on the cell#1 is 100.
 127
 128 The aim here is to compute the interpolated field FT on the target mesh of field FS depending on the \ref NatureOfField "nature of the field".
 129
 130 \anchor TableNatureOfFieldEx1
 131 \image html NonOverlapping.png "An example of non overlapping intersection of two meshes."
 132
 133 The first step of the interpolation leads to the following M1 matrix :
 134
 135 \f[
 136     M1=\left[\begin{tabular}{cc}
 137     0.125 & 0.75 \\
 138     \end{tabular}\right]
 139     \f]
 140
 141 \subsection TableNatureOfFieldExampleConservVol Conservative volumic case
 142
 143 If we apply the formula \ref TableNatureOfField "above" it leads to the following \f$ M_{Conservative Volumic} \f$ matrix :
 144
 145 \f[
 146     M_{Conservative Volumic}=\left[\begin{tabular}{cc}
 147     $\displaystyle{\frac{0.125}{0.125+0.75}}$ &
 148     $\displaystyle{\frac{0.75}{0.125+0.75}}$ \\
 149     \end{tabular}\right]=\left[\begin{tabular}{cc}
 150     0.14286 & 0.85714 \\
 151     \end{tabular}\right]
 152 \f]
 153 \f[
 154     FT=\left[\begin{tabular}{cc}
 155     $\displaystyle\frac{0.125}{0.875}$ & $\displaystyle\frac{0.75}{0.875}$ \\
 156     \end{tabular}\right].\left[\begin{tabular}{c}
 157     4 \\
 158     100 \\
 159     \end{tabular}\right]
 160     =\left[\begin{tabular}{c}
 161     86.28571\\
 162     \end{tabular}\right]
 163 \f]
 164
 165 As we can see here the maximum principle is respected.This nature of field is particulary recommended to interpolate an intensive
 166 field such as \b temperature or \b pression.
 167
 168 \subsection TableNatureOfFieldExampleIntegral Integral case
 169
 170 If we apply the formula \ref TableNatureOfField "above" it leads to the following \f$ M_{Integral} \f$ matrix :
 171
 172 \f[
 173     M_{Integral}=\left[\begin{tabular}{cc}
 174     $\displaystyle{\frac{0.125}{9}}$ & $\displaystyle{\frac{0.75}{3}}$ \\
 175     \end{tabular}\right]=\left[\begin{tabular}{cc}
 176     0.013888 & 0.25 \\
 177     \end{tabular}\right]
 178 \f]
 179 \f[
 180     FT=\left[\begin{tabular}{cc}
 181     $\displaystyle{\frac{0.125}{9}}$ & $\displaystyle{\frac{0.75}{3}}$ \\
 182     \end{tabular}\right].\left[\begin{tabular}{c}
 183     4 \\
 184     100 \\
 185     \end{tabular}\right]
 186     =\left[\begin{tabular}{c}
 187     25.055\\
 188     \end{tabular}\right]
 189 \f]
 190
 191 This type of interpolation is typically recommended for the interpolation of \b power (\b NOT \b power \b density !) for
 192 a user who wants to conserve the quantity \b only on the intersecting part of the source mesh (the green part on the \ref TableNatureOfFieldEx1 "example")
 193
 194 This type of interpolation is equivalent to the computation of \f$ FS_{vol} \f$ followed by a multiplication by \f$ M1 \f$ where \f$ FS_{vol} \f$ is given by :
 195
 196 \f[
 197    FS_{vol}=\left[\begin{tabular}{c}
 198     $\displaystyle{\frac{4}{9}}$ \\
 199     $\displaystyle{\frac{100}{3}}$ \\
 200     \end{tabular}\right]
 201 \f]
 202
 203 In the particular case treated \ref TableNatureOfFieldEx1 "here", it means that only a power of 25.055 W is intercepted by the target cell !
 204
 205 So from the 104 W of the source field \f$ FS \f$, only 25.055 W are transmited in the target field using this nature of field.
 206 In order to treat differently a power field, another policy, \ref TableNatureOfFieldExampleIntegralGlobConstraint "integral global constraint nature" is available.
 207
 208 \subsection TableNatureOfFieldExampleIntegralGlobConstraint Integral with global constraints case
 209
 210 If we apply the formula \ref TableNatureOfField "above" it leads to the following \f$ M_{IntegralGlobConstraint} \f$ matrix :
 211
 212 \f[
 213     M_{IntegralGlobConstraint}=\left[\begin{tabular}{cc}
 214     $\displaystyle{\frac{0.125}{0.125}}$ & ${\displaystyle\frac{0.75}{0.75}}$ \\
 215     \end{tabular}\right]=\left[\begin{tabular}{cc}
 216     1 & 1 \\
 217     \end{tabular}\right]
 218 \f]
 219 \f[
 220     FT=\left[\begin{tabular}{cc}
 221     1 & 1 \\
 222     \end{tabular}\right].\left[\begin{tabular}{c}
 223     4 \\
 224     100 \\
 225     \end{tabular}\right]
 226     =\left[\begin{tabular}{c}
 227     104\\
 228     \end{tabular}\right]
 229 \f]
 230
 231 This type of interpolation is typically recommended for the interpolation of \b power (\b NOT \b power \b density !) for
 232 a user who wants to \b conserve \b all \b the \b power in its source field. Here we have 104 W in source field, we have 104 W too,
 233 in the output target interpolated field.
 234
 235 \b BUT, As we can see here, the maximum principle is \b not respected here, because the target cell #0 has a value higher than the two
 236 intercepted source cells.
 237
 238 \subsection TableNatureOfFieldExampleRevIntegral Reverse integral case
 239
 240 If we apply the formula \ref TableNatureOfField "above" it leads to the following \f$ M_{RevIntegral} \f$ matrix :
 241
 242 \f[
 243     M_{RevIntegral}=\left[\begin{tabular}{cc}
 244     $\displaystyle{\frac{0.125}{1.5}}$ & $\displaystyle{\frac{0.75}{1.5}}$ \\
 245     \end{tabular}\right]=\left[\begin{tabular}{cc}
 246     0.083333 & 0.5 \\
 247     \end{tabular}\right]
 248 \f]
 249 \f[
 250     FT=\left[\begin{tabular}{cc}
 251     $\displaystyle{\frac{0.125}{1.5}}$ & $\displaystyle{\frac{0.75}{1.5}}$ \\
 252     \end{tabular}\right].\left[\begin{tabular}{c}
 253     4 \\
 254     100 \\
 255     \end{tabular}\right]
 256     =\left[\begin{tabular}{c}
 257     50.333\\
 258     \end{tabular}\right]
 259 \f]
 260
 261 This type of nature is particulary recommended to interpolate an intensive \b density
 262 field (moderator density, power density).
 263 The difference with \ref TableNatureOfFieldExampleConservVol "conservative volumic" seen above is that here the
 264 target field is homogeneized to the \b whole target cell. It explains why this nature of field does not follow the maximum principle.
 265
 266 To illustrate the case, let's consider that \f$ FS \f$ is a power density field in \f$ W/m^2 \f$.
 267 With this nature of field the target cell #0 cumulates 0.125*4=0.5 W of power from the source cell #0 and 0.75*100=75 W of power from the source cell #1.
 268 It leads to 75.5 W of power on the \b whole target cell #0. So, the final power density is equal to 75.5/1.5=50.333 W/m^2.
 269
 270 */
 271