Docstoc

UNCERTAINTY QUANTIFICATION

Document Sample
UNCERTAINTY QUANTIFICATION Powered By Docstoc
					 Marta D’Elia                                       Emory University, Math&CS Dept.
 joint work with M. Perego and A. Veneziani




A DATA ASSIMILATION TECHNIQUE FOR
    INCLUDING NOISY VELOCITY MEASURES
              INTO NAVIER-STOKES SIMULATIONS




                       UNCERTAINTY QUANTIFICATION
                                     International Centre for Mathematical Sciences
                                                          Edinburgh, May 24-28 2010
    motivation
        input:    biomedical imaging provides efficient techniques for the collection of a
                  large amount of data

          goal:   use these data to have an accurate approximation of the blood flow in vessels

ultimate goal:    estimation of reliability of our results
                  prediction of the occurrence of diseases in the patient

          how:    including the data in the numerical simulations combining measurements and
                  dynamic principles governing the system




   introduction      linear case        nonlinear case        numerical results       conclusions
    motivation
        input:       biomedical imaging provides efficient techniques for the collection of a
                     large amount of data

          goal:      use these data to have an accurate approximation of the blood flow in vessels

ultimate goal:       estimation of reliability of our results
                     prediction of the occurrence of diseases in the patient

          how:       including the data in the numerical simulations combining measurements and
                     dynamic principles governing the system


                  observed data + mathematical model = Data Assimilation (DA)




          Source: Dr. Brummer, Emory CHOA



   introduction          linear case        nonlinear case       numerical results       conclusions
     state of the art

history: firstly introduced in the 50’s in fluid geophysics

           data: sparse, irregularly distributed and noisy
           governing principles: known and described by a mathematical model


            “Analysis”            Estimation Theory        Dynamic         Control Theory (CT)     time
           Interpolation       Kalman-type techniques     relaxation      Stochastic techniques    line




   introduction            linear case       nonlinear case       numerical results       conclusions
     state of the art

history: firstly introduced in the 50’s in fluid geophysics

           data: sparse, irregularly distributed and noisy
           governing principles: known and described by a mathematical model


            “Analysis”            Estimation Theory         Dynamic         Control Theory (CT)     time
           Interpolation       Kalman-type techniques      relaxation      Stochastic techniques    line




proposed approach: availability of efficient PDE’s solvers and numerical techniques for the
                           solution of optimization problems: choice of CT approaches

                           among those, on the basis of preliminary comparison results, we propose a
                           Discretize-then-Optimize (DO) technique




   introduction            linear case       nonlinear case        numerical results       conclusions
     formulation
notations :        vessel         domain             with boundaries         ,       ,

                   variables      u,p
                   data                    , vector of measured velocites on sites




                                                                                 source: Dr. Brummer




    introduction    linear case     nonlinear case       numerical results       conclusions
     formulation
notations :              vessel         domain             with boundaries         ,       ,

                         variables      u,p
                         data                    , vector of measured velocites on sites



state equations:


                   (S)




                                                                                       source: Dr. Brummer




    introduction          linear case     nonlinear case       numerical results       conclusions
     formulation
notations :              vessel             domain             with boundaries         ,       ,

                         variables          u,p
                         data                        , vector of measured velocites on sites



state equations:


                   (S)




data assimilation:         find      s.t.
                                                                                           source: Dr. Brummer


                           under the constraint of (S)


    introduction          linear case         nonlinear case       numerical results       conclusions
      Numerical solution of the linear pb: Oseen
discretize using the Finite Element (FE) method the state eq.s and the functional

                               discretized u and p

                               restriction matrix

                               selection matrix

                               boundary mass matrix

                               Oseen matrix




    introduction       linear case         nonlinear case     numerical results     conclusions
      Numerical solution of the linear pb: Oseen
discretize using the Finite Element (FE) method the state eq.s and the functional

                               discretized u and p

                               restriction matrix

                               selection matrix

                               boundary mass matrix

                               Oseen matrix




optimize solving the KKT system induced by the Lagrangian with the Reduced Hessian method:


                                                       Reduced Hessian Hr



                                                        Sensitivity matrix


    introduction       linear case         nonlinear case         numerical results   conclusions
 DA procedure for Navier-Stokes eq.s
nonlinear constraint:




the nonlinearity introduced by the advection
term makes the forward problem challenging




introduction         linear case         nonlinear case   numerical results   conclusions
   DA procedure for Navier-Stokes eq.s
 nonlinear constraint:




 the nonlinearity introduced by the advection
 term makes the forward problem challenging




idea:       exploit Picard/Newton method for treating the nonlinearity
            use the DA procedure proposed for linear problems


                   Newton-like iterative procedure




 introduction         linear case         nonlinear case    numerical results   conclusions
    well-posedness
note:     the regularized formulation is always well-posed
          sufficient condition for an equality PDE constrained opt. pb: Hr positive definite
          goal: find conditions on the selection matrix sufficient for the well-posedness




  introduction        linear case        nonlinear case         numerical results        conclusions
     well-posedness
note:      the regularized formulation is always well-posed
           sufficient condition for an equality PDE constrained opt. pb: Hr positive definite
           goal: find conditions on the selection matrix sufficient for the well-posedness


                                                 T
result:   for a=0, if Null(D) ∩ Range(S-1Rin Min)={0} then Hr is positive definite and the
          the problem is well-posed.




   introduction        linear case        nonlinear case         numerical results        conclusions
     well-posedness
note:      the regularized formulation is always well-posed
           sufficient condition for an equality PDE constrained opt. pb: Hr positive definite
           goal: find conditions on the selection matrix sufficient for the well-posedness


                                                 T
result:   for a=0, if Null(D) ∩ Range(S-1Rin Min)={0} then Hr is positive definite and the
          the problem is well-posed.
          This condition is satisfied by choosing Q s.t. its restriction to rows corresponding to
          sites on Γin has rank Nin (degrees of freedom, DOF, of U on Γin)




   introduction        linear case        nonlinear case         numerical results        conclusions
     well-posedness
note:      the regularized formulation is always well-posed
           sufficient condition for an equality PDE constrained opt. pb: Hr positive definite
           goal: find conditions on the selection matrix sufficient for the well-posedness


                                                   T
result:   for a=0, if Null(D) ∩ Range(S-1Rin Min)={0} then Hr is positive definite and the
          the problem is well-posed.
          This condition is satisfied by choosing Q s.t. its restriction to rows corresponding to
          sites on Γin has rank Nin (degrees of freedom, DOF, of U on Γin)




practically:           sites on grid nodes                          sparse sites on
                       (P1bubble-P1)                                inflow boundary
                                                                    (P1bubble-P1)




   introduction        linear case           nonlinear case      numerical results        conclusions
    well-posedness
idea:      given sparse measurement on inflow (not satisfying assumptions for well-posedness)
           recover the well-posedness with approximated data on grid nodes on inflow




   introduction       linear case      nonlinear case        numerical results      conclusions
    well-posedness
idea:      given sparse measurement on inflow (not satisfying assumptions for well-posedness)
           recover the well-posedness with approximated data on grid nodes on inflow




approximation:       interpolation of each velocity component of the given data
                     dj = Πkd(xj)   where Πkd is the interpolator recovered from k data.




   introduction       linear case        nonlinear case        numerical results       conclusions
     well-posedness
idea:      given sparse measurement on inflow (not satisfying assumptions for well-posedness)
           recover the well-posedness with approximated data on grid nodes on inflow




approximation:         interpolation of each velocity component of the given data
                       dj = Πkd(xj)       where Πkd is the interpolator recovered from k data.




result:    if the original data is s.t.
            di = uex(yi) + εi where εi ≈ N (0, σ2)
           then, the new data on inflow grid nodes is s.t
           dj = uex(yj) + ηj where η j ≈ N (μj, ν2j)




   introduction         linear case            nonlinear case        numerical results       conclusions
     well-posedness
idea:       given sparse measurement on inflow (not satisfying assumptions for well-posedness)
            recover the well-posedness with approximated data on grid nodes on inflow




approximation:          interpolation of each velocity component of the given data
                        dj = Πkd(xj)       where Πkd is the interpolator recovered from k data.




result:     if the original data is s.t.
             di = uex(yi) + εi where εi ≈ N (0, σ2)
            then, the new data on inflow grid nodes is s.t
            dj = uex(yj) + ηj where η j ≈ N (μj, ν2j)


mean:       affected by discretization error of the interpolation process
variance:   affected by the error caused by interpolation of noisy data



   introduction          linear case            nonlinear case              numerical results   conclusions
    analytic test case - setting
test case: the advection vector corresponds to uex




  introduction         linear case          nonlinear case   numerical results   conclusions
    analytic test case - setting
test case: the advection vector corresponds to uex           mesh generation:
                                                             software - FreeFem++

                                                                                    (1.5, 2)




                                                             (-.5, 0)




  introduction         linear case          nonlinear case      numerical results       conclusions
     analytic test case - setting
test case: the advection vector corresponds to uex              mesh generation:
                                                                software - FreeFem++

                                                                                       (1.5, 2)




                                                                (-.5, 0)




data generation:

di = uex(yi) + εi

simplified case: εi ≈ σ U(-0.5, 0.5)

variance determined s.t. the signal-noise ratio (SNR)
is fixed a priori according to experimental results

estimated ratio SNR = 8.3


   introduction          linear case           nonlinear case      numerical results       conclusions
     analytic test case - setting
test case: the advection vector corresponds to uex               mesh generation:
                                                                 software - FreeFem++

                                                                                           (1.5, 2)




                                                                 (-.5, 0)




data generation:                                                implementation details:

di = uex(yi) + εi                                               DA solver – lifeV
                                                                external library – AZTEC
simplified case: εi ≈ σ U(-0.5, 0.5)                            postprocessing – PARAVIEW
                                                                FE spaces – P1bubble-P1
variance determined s.t. the signal-noise ratio (SNR)
                                                                Oseen solver – Monolithic PCD
is fixed a priori according to experimental results
                                                                 ν – 0.035
estimated ratio SNR = 8.3


   introduction          linear case           nonlinear case       numerical results           conclusions
     analytic test case – error behavior

sanity check: with noise-free data the FE convergence rate is recovered




    introduction      linear case       nonlinear case       numerical results   conclusions
     analytic test case – error behavior
noisy data: error behavior for the Oseen test case
           1. error of the sample mean over Nr realizations, UmNr, vs number of realizations, Nr :

              ||UmNr - Uex||2  E∆x     as Nr      , O(Nr-0.5)




                                                8
              since Uex = Um+ E∆x        where E∆x is the noise-free discretization error
                                         and Um is the mean value for the variable UmNr




   introduction      linear case         nonlinear case           numerical results     conclusions
     analytic test case – error behavior
noisy data: error behavior for the Oseen test case
           1. error of the sample mean over Nr realizations, UmNr, vs number of realizations, Nr :

              ||UmNr - Uex||2  E∆x     as Nr      , O(Nr-0.5)




                                                8
              since Uex = Um+ E∆x        where E∆x is the noise-free discretization error
                                         and Um is the mean value for the variable UmNr

              ∆x = 0.16




   introduction           linear case    nonlinear case           numerical results     conclusions
     analytic test case – error behavior
noisy data: error behavior for the Oseen test case
           1. error of the sample mean over Nr realizations, UmNr, vs number of realizations, Nr :

              ||UmNr - Uex||2  E∆x     as Nr      , O(Nr-0.5)




                                                8
              since Uex = Um+ E∆x        where E∆x is the noise-free discretization error
                                         and Um is the mean value for the variable UmNr

              ∆x = 0.072




   introduction        linear case       nonlinear case           numerical results     conclusions
     analytic test case – error behavior
noisy data: error behavior for the Oseen test case
           1. error of the sample mean over Nr realizations, UmNr, vs number of realizations, Nr :

              ||UmNr - Uex||2  E∆x     as Nr      , O(Nr-0.5)




                                                8
              since Uex = Um+ E∆x        where E∆x is the noise-free discretization error
                                         and Um is the mean value for the variable UmNr

              ∆x = 0.05




   introduction           linear case    nonlinear case           numerical results     conclusions
     analytic test case – error behavior
noisy data: error behavior for the Oseen test case
           2. average error vs Ns , O(Ns-0.5) with ∆x = 0.16, SNR = 10




   introduction      linear case      nonlinear case       numerical results   conclusions
     analytic test case – interpolation validation
comparison: regularization vs interpolation


parameters: SNR = 20; ∆x = 0.072




   introduction     linear case      nonlinear case   numerical results   conclusions
   towards real geometries – velocity & vorticity
 curved domain: approximation of a section of the carotid

data generation: reference solution (FE solution on very fine grid) with addition of white noise




  introduction       linear case       nonlinear case        numerical results      conclusions
   towards real geometries – velocity & vorticity
 curved domain: approximation of a section of the carotid

data generation: reference solution (FE solution on very fine grid) with addition of white noise




  introduction       linear case       nonlinear case        numerical results      conclusions
     final remarks

conclusion

1.    the DA procedure proved to be consisitent in case of non-noisy data

2.    in presence of noise it is a tool for filtering the unavoidable uncertainties in the
      measured data: the assimilation of boundary and internal measures enhances the
      accuracy and reliability of the recovered variables

3.    the combination of the DA procedure with the interpolation technique is a
      competitive tool with respect to common regularization approaches




state of the art       approaches            linear case         nonlinear case         conclusions
     final remarks

conclusion

1.    the DA procedure proved to be consisitent in case of non-noisy data

2.    in presence of noise it is a tool for filtering the unavoidable uncertainties in the
      measured data: the assimilation of boundary and internal measures enhances the
      accuracy and reliability of the recovered variables

3.    the combination of the DA procedure with the interpolation technique is a
      competitive tool with respect to common regularization approaches



future guidelines

1.    precise estimation of statistics related to assimilated variables via more
      sophisticated stochastic tools
2.    assimilation of velocities from real geometries and data




state of the art       approaches            linear case         nonlinear case         conclusions
         thank you for your attention



         questions?                 …

                                    ...




         special thanks to          Alessandro Veneziani

                                    Mauro Perego

                                    Michele Benzi

                                    Max Gunzburger




state of the art      approaches          linear case      nonlinear case   conclusions
     more...
comparison: (accuracy and conditioning) vs (data location)




   introduction      linear case     nonlinear case     numerical results   conclusions
 more...

                                    outer iterations: GMRES




                                    system in S to be solved twice (transpose)




                                    decomposed using exact factorization




                                    1st inner iteration: system in C: GMRES
                                    2nd inner iteration: system in Σ (*) : GMRES
                                    preconditioner needed (Elman, Silverster, Wathen)
          (*)



state of the art   approaches   linear case          nonlinear case              conclusions
   more...
results with Stream line Diffusion stabilization for noise/free Oseen




 state of the art          approaches               linear case         nonlinear case   conclusions