42 | | == The Interpolation == |
| 42 | == Preparing for the Interpolation == |
| 43 | There are several steps that need to occur before the actual interpolation. When Zcooling.tab is read in, the cooling rates are assigned to a properly allocated 3D array called Zcoolingtab(i,j,k). The indices i, j, k of Zcoolingtab are integers ranging from 1:nDensities, 1:nTemps, and 1:nXs respectively. These indices can also be thought of as coordinates. Before the cooling rate can be calculated, values for ne, T, and X are calculated. The electron density ne and ionization fraction X are converted to log(ne) and log(X). Then these values for log(ne), T, and log(X) can be scaled to the integer ranges of the coordinates i, j, and k. This is where the minimum and maximum values in Zcooling.tab are used. |
| 44 | |
| 45 | In order to do a 3D interpolation, we need to find the 8 closest points to our newly scaled point of interest. You can think of these 8 points as the vertices of a cube or rectangular prism, and our point of interest is somewhere inside this box. First, we find the 8 closest coordinates. Since the coordinates i, j, k are all integers, the 8 closest coordinates are easily found by using the FLOOR and CEILING functions. Using these functions, you get a "low" and "high" value for each coordinate. Since there are 3 coordinates, 2*2*2 = 8 unique sets of coordinates. |
| 46 | |
| 47 | The next step is to match these 8 closest coordinates with their corresponding cooling rates. This is done easily with 8 simple assignments. For example, rate1 = Zcoolingtab(lognelow, templow, logxlow) would be the rate that corresponds to the lower bound coordinates. Now we have the 8 closest points. |
| 48 | |
| 49 | The final step before interpolating is to do another rescaling of the coordinates. We take the 8 closest coordinates to be the coordinates of a unit cube aligned with the origin. The coordinates of our point of interest are then scaled accordingly. This makes the interpolation very simple. Remember that in all these rescaling steps, nothing happens to the cooling rate because this is the value that we are interested in. |
| 52 | |
| 53 | == Trilinear Interpolation == |
| 54 | Now we have 8 points that are positioned on the vertices of a unit cube aligned with the origin. Each point has a unique cooling rate. Our point of interest is inside this cube, but we do not yet know what its cooling rate should be. Trilinear interpolation is the method used to find this rate. The formula used is from Paul Bourke's website: [http://paulbourke.net/miscellaneous/interpolation/]. In cooling.f90, this is done on 8 lines to avoid having one lengthy, confusing line... |
| 55 | {{{ |
| 56 | ZCoolingRate = rate1*(1d0-lognein)*(1d0-tempin)*(1d0-logxin) |
| 57 | ZCoolingRate = ZCoolingRate + rate2*lognein*(1d0-tempin)*(1d0-logxin) |
| 58 | ZCoolingRate = ZCoolingRate + rate3*(1d0-lognein)*tempin*(1d0-logxin) |
| 59 | ZCoolingRate = ZCoolingRate + rate4*(1d0-lognein)*(1d0-tempin)*logxin |
| 60 | ZCoolingRate = ZCoolingRate + rate5*lognein*(1d0-tempin)*logxin |
| 61 | ZCoolingRate = ZCoolingRate + rate6*(1d0-lognein)*tempin*logxin |
| 62 | ZCoolingRate = ZCoolingRate + rate7*lognein*tempin*(1d0-logxin) |
| 63 | ZCoolingRate = ZCoolingRate + rate8*lognein*tempin*logxin |
| 64 | }}} |
| 65 | |
| 66 | Where rate1,rate2,... are assigned as aforementioned, and the "in" values represent the coordinates of our point of interest inside the unit cube. |