Giter Club home page Giter Club logo

hippyflow's Introduction

CI Doc Status status DOI

              Inverse Problem PYthon library
 __        ______  _______   _______   __      __  __  __  __       
/  |      /      |/       \ /       \ /  \    /  |/  |/  |/  |      
$$ |____  $$$$$$/ $$$$$$$  |$$$$$$$  |$$  \  /$$/ $$ |$$/ $$ |____  
$$      \   $$ |  $$ |__$$ |$$ |__$$ | $$  \/$$/  $$ |/  |$$      \ 
$$$$$$$  |  $$ |  $$    $$/ $$    $$/   $$  $$/   $$ |$$ |$$$$$$$  |
$$ |  $$ |  $$ |  $$$$$$$/  $$$$$$$/     $$$$/    $$ |$$ |$$ |  $$ |
$$ |  $$ | _$$ |_ $$ |      $$ |          $$ |    $$ |$$ |$$ |__$$ |
$$ |  $$ |/ $$   |$$ |      $$ |          $$ |    $$ |$$ |$$    $$/ 
$$/   $$/ $$$$$$/ $$/       $$/           $$/     $$/ $$/ $$$$$$$/  
              https://hippylib.github.io

hIPPYlib implements state-of-the-art scalable algorithms for deterministic and Bayesian inverse problems governed by partial differential equations (PDEs). It builds on FEniCS (a parallel finite element element library) for the discretization of the PDE and on PETSc for scalable and efficient linear algebra operations and solvers.

For building instructions, see the file INSTALL.md. Copyright information and licensing restrictions can be found in the file COPYRIGHT.

The best starting point for new users interested in hIPPYlib's features are the interactive tutorials in the tutorial folder.

Conceptually, hIPPYlib can be viewed as a toolbox that provides the building blocks for experimenting new ideas and developing scalable algorithms for PDE-constrained deterministic and Bayesian inverse problems.

In hIPPYlib the user can express the forward PDE and the likelihood in weak form using the friendly, compact, near-mathematical notation of FEniCS, which will then automatically generate efficient code for the discretization. Linear and nonlinear, and stationary and time-dependent PDEs are supported in hIPPYlib. For stationary problems, gradient and Hessian information can be automatically generated by hIPPYlib using FEniCS symbolic differentiation of the relevant weak forms. For time-dependent problems, instead, symbolic differentiation can only be used for the spatial terms, and the contribution to gradients and Hessians arising from the time dynamics needs to be provided by the user.

Noise and prior covariance operators are modeled as inverses of elliptic differential operators allowing us to build on existing fast multigrid solvers for elliptic operators without explicitly constructing the dense covariance operator.

The key property of the algorithms underlying hIPPYlib is that solution of the deterministic and Bayesian inverse problem is computed at a cost, measured in forward PDE solves, that is independent of the parameter dimension.

hIPPYlib provides a robust implementation of the inexact Newton-conjugate gradient algorithm to compute the maximum a posterior (MAP) point. The gradient and Hessian actions are computed via their weak form specification in FEniCS by constraining the state and adjoint variables to satisfy the forward and adjoint problem. The Newton system is solved inexactly by early termination of CG iterations via Eisenstat-Walker (to prevent oversolving) and Steihaug (to avoid negative curvature) criteria. Two globalization techniques are available to the user: Armijo back-tracking line search and trust region.

In hIPPYlib, the posterior covariance is approximated by the inverse of the Hessian of the negative log posterior evaluated at the MAP point. This Gaussian approximation is exact when the parameter-to-observable map is linear; otherwise, its logarithm agrees to two derivatives with the log posterior at the MAP point, and thus it can serve as a proposal for Hessian-based Markov chain Monte Carlo (MCMC) methods. hIPPYlib makes the construction of the posterior covariance tractable by invoking a low-rank approximation of the Hessian of the log likelihood.

hIPPYlib also offers scalable methods for sample generation. To sample large scale spatially correlated Gaussian random fields from the prior distribution, hIPPYlib implements a new method that strongly relies on the structure of the covariance operator defined as the inverse of a differential operator: by exploiting the assembly procedure of finite element matrices hIPPYlib constructs a sparse Cholesky-like rectangular decomposition of the precision operator. To sample from a local Gaussian approximation to the posterior (such as at the MAP point) hIPPYlib exploits the low rank factorization of the Hessian of the log likelihood to correct samples from the prior distribution. Finally, to explore the posterior distribution, hIPPYlib implements dimension independent MCMC sampling methods enchanted by Hessian information.

Finally, randomized and probing algorithms are available to compute the pointwise variance of the prior/posterior distribution and the trace of the covariance operator.

hippyflow's People

Contributors

dc-luo avatar tomoleary avatar

Stargazers

 avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar

Watchers

 avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar  avatar

hippyflow's Issues

Mesh parallelism correctness for decompositions

  • What needs to be done to ensure the decompositions and numpy representations of multivectors make sense for mesh partitioned problems? So far instance parallelism works, but there could be issues for mesh parallelism.

More efficient Jacobian data generation

In the case that the full Jacobian is requested, and it is either short-fat or tall-skinny randomized SVD is inefficient, and its better to compute action on identity. This case can be detected and implemented in ActiveSubspaceProjector object.

Correct projections / bases

In the hf.DataGenerator class the input basis should be used and the output projector should be used (including Riesz representer). The naming conventions in POD and KLE and AS should be updated to reflect this. Also the key word arguments should be updated in the data generator to reflect this.

Memory footprints of ActiveSubspace and POD for large problems.

The way that these decompositions are currently computed (as well as data generation) are memory intensive operations, since they require storing many vectors simultaneously.

Reducing the memory footprint for data generation is easy, and already mostly fixed on tom_dev branch.

Reducing the memory footprint for the active subspace decomposition requires refactoring all the way from the randomized eigensolver down. The strategy is to iterate over samples one by one, solve the forward problem, set linearization points and then perform the entire Gaussian random matrix for that sample, updating the target matrix in place, and scaling accordingly for the averaging operation. The issue is that MatMVMult in hippylib iteratively applies a given operator to different vectors to perform a matrix-matrix product, but we want to put this MatMVMult itself inside of a for loop, so we need to create a custom averagedMatMVMult function that then can apply the entire matrix-matrix action for one sample, and then iterate over samples until the entire action is complete.

Additionally because the doublePass algorithm calls this operation twice, the samples need to be exactly the same for successive calls to this operator, but we want the option to not have to save the samples directly, so we need a deterministic sampling strategy (this requires custom instance of Random with deterministic burn-ins). This raises a question about the parRandom when it is on a communicator mesh grid.

pip installation / setup.py

@tomoleary, it'd be nice to be able to pip install hippyflow instead of having to create an environment variable to find the library.

  1. How much can we copy from hippylib's setup.py and INSTALL.md?
  2. What other dependencies & optionals would I need to add?

Mass and stiffness matrices

This is a minor addition, but adding methods that get dense mass and stiffness matrices from the POD object would be very useful.

Recommend Projects

  • React photo React

    A declarative, efficient, and flexible JavaScript library for building user interfaces.

  • Vue.js photo Vue.js

    ๐Ÿ–– Vue.js is a progressive, incrementally-adoptable JavaScript framework for building UI on the web.

  • Typescript photo Typescript

    TypeScript is a superset of JavaScript that compiles to clean JavaScript output.

  • TensorFlow photo TensorFlow

    An Open Source Machine Learning Framework for Everyone

  • Django photo Django

    The Web framework for perfectionists with deadlines.

  • D3 photo D3

    Bring data to life with SVG, Canvas and HTML. ๐Ÿ“Š๐Ÿ“ˆ๐ŸŽ‰

Recommend Topics

  • javascript

    JavaScript (JS) is a lightweight interpreted programming language with first-class functions.

  • web

    Some thing interesting about web. New door for the world.

  • server

    A server is a program made to process requests and deliver data to clients.

  • Machine learning

    Machine learning is a way of modeling and interpreting data that allows a piece of software to respond intelligently.

  • Game

    Some thing interesting about game, make everyone happy.

Recommend Org

  • Facebook photo Facebook

    We are working to build community through open source technology. NB: members must have two-factor auth.

  • Microsoft photo Microsoft

    Open source projects and samples from Microsoft.

  • Google photo Google

    Google โค๏ธ Open Source for everyone.

  • D3 photo D3

    Data-Driven Documents codes.