Changes between Version 7 and Version 8 of ThinPlateShell


Ignore:
Timestamp:
Dec 9, 2014, 5:30:04 PM (6 years ago)
Author:
sommeria
Comment:

--

Legend:

Unmodified
Added
Removed
Modified
  • ThinPlateShell

    v7 v8  
    99We consider a set of measurement points ${\bf r_i}, i = 1,2, \ldots,N$, and the corresponding set of measurements values for a quantity $f$ (a velocity component for PIV data).
    1010
    11 $ f({\bf r_i})=f_i, i = 1,2, \ldots,N} \; [equ. 1] $
     11 $ f({\bf r_i})=f_i, i = 1,2, \ldots,N} \; [equ. 1] $
    1212
    1313A pure interpolation function $f({\bf r})$ must exactly reach these values $f_i $ at the measurement points, while a  smoothing function should approach these values within a range smaller than the estimated error bar.
     
    1515Considering first the case of a pure interpolation, the function $f({\bf r})$ must minimise the 'elastic energy' or mean squared curvature
    1616
    17                             $ E = \int\left[\left(\frac{\partial^2 f}{\partial x^2}\right)^2 + 2\left(\frac{\partial^2 f}{\partial xy}\right)^2 + \left(\frac{\partial^2 f}{\partial y^2}\right)^2 \right] \textrm{d} x \, \textrm{d}y  \;\textrm{[equ. 2]} $
     17 $ E = \int\left[\left(\frac{\partial^2 f}{\partial x^2}\right)^2 + 2\left(\frac{\partial^2 f}{\partial xy}\right)^2 + \left(\frac{\partial^2 f}{\partial y^2}\right)^2 \right] \textrm{d} x \, \textrm{d}y  \;\textrm{[equ. 2]} $
    1818
    19   with the constraints [->equ. 1]. These constraints can be also written as a domain integral with Dirac functions at measurement points:
     19with the constraints [->equ. 1]. These constraints can be also written as a domain integral with Dirac functions at measurement points:
    2020
    21 $ {\int f \delta({\bf r-r_i}) \textrm{d} x \, \textrm{d}y=f_i , i = 1,2, \ldots,N} \;\textrm{[equ. 3]}  $
     21 $ {\int f \delta({\bf r-r_i}) \textrm{d} x \, \textrm{d}y=f_i , i = 1,2, \ldots,N} \;\textrm{[equ. 3]}  $
    2222
    2323This variational problem is solved by introducing a Lagrange parameter $8\pi S_i$ for each constraint (the factor $8\pi$ is put for further convenience). The variation of the energy can be rewritten by a double integration by part, leading to $\int \Delta^2 f \delta f dx dy$. The optimum function $f$ is thus solution of the  double Laplacian equation with pointwise source terms at locations $r_i$,
     
    2727The solution $\phi$ for a single source point is such that $\Delta^2 \phi=0$ everywhere outside the source. Using the axisymmetric form of the Laplacian $\Delta \phi=\frac{1}{r}\frac{d}{dr}(r\frac{d\phi}{dr})$, we can show that the general axisymmetric solution is
    2828
    29 $\phi=C_0 \log(r)+C_1 r^2 + C_2 r^2 \log (r)+C_3 \;\textrm{[equ. 5]}$$ whose Laplacian writes $$\Delta\phi=\frac{1}{r}\frac{d}{dr}(r\frac{d\phi}{dr})=4 C_1+ 4 C_2 (1+\log(r))\;\textrm{[equ. 6]}$
     29 $\phi=C_0 \log(r)+C_1 r^2 + C_2 r^2 \log (r)+C_3 \;\textrm{[equ. 5]}$
     30whose Laplacian writes
     31
     32 $\Delta\phi=\frac{1}{r}\frac{d}{dr}(r\frac{d\phi}{dr})=4 C_1+ 4 C_2 (1+\log(r))\;\textrm{[equ. 6]}$
    3033
    3134The coefficient $C_0$ must be zero to avoid divergence of $\phi$ at the origin. The constant $C_3$ can be included in a global constant $a_0$. The constant $C_1$ can be viewed as a change of unit for $r$ (replacing $r$ by $r/r_0$ with $\log(r_0)=-C_1/C_2)$), and can be set to 0 without loss of generality. The integral of $\Delta\phi$ over a small neighborhood of the source must be equal to $8\pi S_i$. By the divergence theorem, this is also equal to the flux of $\nabla \Delta \phi=4 C_2/r$ around the source, equal to $8\pi C_2$. This must be equal to $8\pi S_i$, which sets $C_2=S_i$. The elementary function is thus $\phi(r)=S_i r^2\log (r)$.^
     
    3336The general form of $f$ is thus
    3437
    35 $\label{sol_gene} f({\bf r})=\sum S_i \phi({\bf|r-r_i}|)+a_0+a_1x+a_2y\;\textrm{[equ. 7]} $
     38 $\label{sol_gene} f({\bf r})=\sum S_i \phi({\bf|r-r_i}|)+a_0+a_1x+a_2y\;\textrm{[equ. 7]} $
    3639
    3740with $\phi(r)=r^2\log (r)$. The  values at the measurement points ${\bf r_j}$ are
    3841
    39 $ f({\bf r_j})=\sum S_i \phi({\bf|r_j-r_i}|)+a_0+a_1x_j+a_2y_j;\;\textrm{[equ. 8]} $
     42 $ f({\bf r_j})=\sum S_i \phi({\bf|r_j-r_i}|)+a_0+a_1x_j+a_2y_j;\;\textrm{[equ. 8]} $
    4043
    4144In other words, defining the vector $F=f({\bf r_j}) \{ j=1,...N \}, a_0, a_1,a_2$ and the matrix
    4245
    43 $ M=\left[  \begin{array}{cc}^
     46 $ M=\left[  \begin{array}{cc}^
    4447  \phi(|r_j-r_i|) , 1,\;  x_j,\;  y_j \\ ... \;\;\;\;\;\;\;\;\;\;\;\;\; ... \;\;\; ... \;\;\; ...\\
    4548
     
    5558The variational problem then leads to the equation
    5659
    57 $ \rho\Delta^2 f= \sum  (f_i-f) \delta ({\bf r-r_i}) $
     60 $ \rho\Delta^2 f= \sum  (f_i-f) \delta ({\bf r-r_i}) $
    5861
    5962 The solution is still obtained as a sum of radial basis functions $\phi({\bf|r-r_i}|)$, whose weight $S_i$ now satisfies $M_\rho*S=F$, with the matrix
    6063
    61 $ M_\rho=M+\rho I_{000} $
     64 $ M_\rho=M+\rho I_{000} $
    6265
    6366where $I_{000}$ is the NxN identity matrix extended by three columns of 0.^
     
    6669Spatial derivatives of the interpolated quantity $f$ can be obtained by direct differentiation of the result.  For any function $\phi(r)$, with radial distance $r=|{\bf r-r_i}|$, $r^2=X^2+Y^2$, $\partial_X \phi=(d\phi/dr) \partial_X r$, and $\partial_X r=X/r$. This yields $\partial_X \phi=X (2 \log(r) +1)$, so that,
    6770
    68 $\partial_x f({\bf r})=\sum S_i (2 \log|{\bf r-r_i}|+1)(x-x_i)+a_1;$
    69 $\partial_y f({\bf r})=\sum S_i (2 \log|{\bf r-r_i}|+1)(y-y_i)+a_2;$
     71 $\partial_x f({\bf r})=\sum S_i (2 \log|{\bf r-r_i}|+1)(x-x_i)+a_1;$
     72 $\partial_y f({\bf r})=\sum S_i (2 \log|{\bf r-r_i}|+1)(y-y_i)+a_2;$
    7073
    7174= Sub-domains =