Repository: CUQI-DTU/CUQIpy-FEniCS (language: python)

Issue:
Add readme for fenics plug in



Discussion / hints:

Fix patch:
diff --git a/README.md b/README.md
index ee9b173..2aa6000 100644
--- a/README.md
+++ b/README.md
@@ -1,2 +1,47 @@
 # CUQIpy-FEniCS
-CUQIpy plug-in for FEniCS PDE modeling
+CUQIpy-FEniCS is a plugin for [CUQIpy](https://github.com/CUQI-DTU/CUQIpy) software package. It provides an interface between FEniCS PDE models and CUQIpy modules.
+
+## Installation
+First install [FEniCS](https://fenicsproject.org/download/archive/) (we recommend using Anaconda from the available installation options). Then install CUQIpy-FEniCS with pip:
+```bash
+pip install cuqipy-fenics
+```
+If CUQIpy is not installed, it will be installed automatically.
+
+## Quickstart
+```python
+import numpy as np
+import matplotlib.pyplot as plt
+import cuqi
+import cuqipy_fenics
+
+# Load a fenics forward model and data from testproblem library
+model, y_data, info = cuqipy_fenics.testproblem.FEniCSDiffusion1D.get_components(
+    dim=20,
+    endpoint=1,
+    exactSolution='smooth_step',
+    mapping='exponential',
+    SNR=10000,
+    left_bc=0,
+    right_bc=8
+)
+
+# Set up Bayesian model
+x = cuqi.distribution.GMRF(np.zeros(model.domain_dim),
+                           25, 1, 'zero', geometry=model.domain_geometry)
+# y ~ N(model(x), 0.01^2)
+y = cuqi.distribution.GaussianCov(model(x), 0.05**2)
+
+# Set up Bayesian Problem object
+BP = cuqi.problem.BayesianProblem(y, x).set_data(y=y_data)
+
+# Sample from the posterior
+samples = BP.sample_posterior(5000)
+
+# Analyze the samples
+samples.burnthin(1000).plot_ci(95, plot_par=True,
+                               exact=info.exactSolution, linestyle='-', marker='.')
+```
+
+For more examples, see the [demos](demos) folder.
+ 
\ No newline at end of file


Test patch:
diff --git a/cuqipy_fenics/testproblem.py b/cuqipy_fenics/testproblem.py
index db649e6..a9f6d7a 100644
--- a/cuqipy_fenics/testproblem.py
+++ b/cuqipy_fenics/testproblem.py
@@ -85,8 +85,11 @@ def form(m,u,p):
 
         dirichlet_bc_expression = dl.Expression("left_bc*(x[0]<eps)+right_bc*(x[0]>endpoint-eps)", eps=dl.DOLFIN_EPS, endpoint=endpoint, left_bc=left_bc, right_bc=right_bc, degree=1)
         dirichlet_bc = dl.DirichletBC(solution_function_space, dirichlet_bc_expression, u_boundary)
+        adjoint_dirichlet_bc = dl.DirichletBC(
+            solution_function_space, dl.Constant(0), u_boundary)
 
-        PDE = cuqipy_fenics.pde.SteadyStateLinearFEniCSPDE( form, mesh, solution_function_space, parameter_function_space, dirichlet_bc,observation_operator=observation_operator)
+        PDE = cuqipy_fenics.pde.SteadyStateLinearFEniCSPDE(
+            form, mesh, solution_function_space, parameter_function_space, dirichlet_bc, adjoint_dirichlet_bc, observation_operator=observation_operator)
         
         # Create PDE model
         domain_geometry = cuqipy_fenics.geometry.FEniCSContinuous(parameter_function_space)
@@ -115,18 +118,26 @@ def form(m,u,p):
             fun = lambda grid:  0.8*np.exp( -( (grid -endpoint/2.0 )**2 ) / 0.02)
             grid = np.linspace(0,endpoint,N)
             exactSolution = np.ones(N)*.8
-            exactSolution[np.where(grid>endpoint/2.0)] = fun(grid[np.where(grid>endpoint/2.0)]) 
+            exactSolution[np.where(grid > endpoint/2.0)
+                          ] = fun(grid[np.where(grid > endpoint/2.0)])
+            exactSolution = cuqi.samples.CUQIarray(
+                exactSolution, geometry=domain_geometry)
 
         # Generate exact data
         b_exact = model.forward(domain_geometry.par2fun(exactSolution),is_par=False)
 
         # Add noise to data
-        sigma = np.linalg.norm(b_exact)/SNR
-        sigma2 = sigma*sigma # variance of the observation Gaussian noise
-        data = b_exact + np.random.normal( 0, sigma, b_exact.shape )
-        
+        # Reference: Adding noise with a desired signal-to-noise ratio
+        # https://sites.ualberta.ca/~msacchi/SNR_Def.pdf
+        noise = np.random.normal(0, 1, b_exact.shape)
+        alpha = np.linalg.norm(b_exact)/(np.sqrt(SNR)*np.linalg.norm(noise))
+
+        data = cuqi.samples.CUQIarray(
+            b_exact + alpha*noise, geometry=range_geometry)
+
         # Create likelihood
-        y = cuqi.distribution.GaussianCov(model(x), sigma2*np.eye(range_geometry.par_dim))
+        y = cuqi.distribution.GaussianCov(
+            model(x), alpha*np.eye(range_geometry.par_dim))
 
         # Initialize FEniCSDiffusion1D as BayesianProblem problem
         super().__init__(y, x, y=data)