# Generating noise in particular way

**URL:** <https://fenicsproject.discourse.group/t/generating-noise-in-particular-way/1505>\
**Category:** Uncategorized\
**Created:** [September 9, 2019, 12:32pm UTC](https://fenicsproject.discourse.group/t/generating-noise-in-particular-way/1505 "2019-09-09T12:32:09Z")\
**Posts on this page:** 8\
**Page:** 1

<div class="post-metadata">

**Author:** ![erja](https://yyz2.discourse-cdn.com/free1/user_avatar/fenicsproject.discourse.group/erja/32/469_2.png) [@erja](https://fenicsproject.discourse.group/u/erja)\
**Post date:** [September 9, 2019, 12:32pm UTC](https://fenicsproject.discourse.group/t/generating-noise-in-particular-way/1505/1 "2019-09-09T12:32:09Z")

</div>

Hi there,

I would like to generate noise in a certain way. Assume I have my FEM-space V\_h.  
Then my noise is to be generated via

\sum\_{k=1}^{N\_h} \xi\_k e\_k

where e\_k is a basis for the FEM space and \xi\_i are normally distributed variables. The issue is to obtain a basis for the FEM space. This method has to be generalizable, so I wonder if there is some “generate basis”-funtion for the FEM-space? There seems to be something like it,  
[https://fenicsproject.org/docs/dolfin/1.6.0/python/programmers-reference/cpp/fem/BasisFunction.html](https://fenicsproject.org/docs/dolfin/1.6.0/python/programmers-reference/cpp/fem/BasisFunction.html)  
but I don’t understand the documenation, and more importantly, how do I get this to work with fenics? I need some sort of “integrable” expression, right? Furthermore, I want the covariance matrix R of the noise to be given by

(R^T R)^{-1}=M^{-1}

where M is the mass matrix of the FEM space Vh.

---

<div class="post-metadata">

**Author:** ![Yannik](https://avatars.discourse-cdn.com/v4/letter/y/5fc32e/32.png) [@Yannik](https://fenicsproject.discourse.group/u/Yannik)\
**Post date:** [September 10, 2019, 9:03am UTC](https://fenicsproject.discourse.group/t/generating-noise-in-particular-way/1505/2 "2019-09-10T09:03:46Z")

</div>

Perhaps this works for you:

```
from fenics import *
import numpy as np

# mesh, FunctionSpace
mesh = UnitIntervalMesh(5)
Vh = FunctionSpace(mesh, 'P', 1)

# allocate function u in Vh
u = Function(Vh)
Nh = u.vector()[:].size
xi = np.random.rand(Nh, ) # or any other random array function

# set the entries on the nodes of u to (the random array) xi
u.vector()[:] = xi

# generating mass matrix
z = TestFunction(Vh)
y = TrialFunction(Vh)
m = y*z*dx
Mass = assemble(m)
# Mass.array() gives you the desired matrix
M = Mass.array()
L = np.linalg.cholesky(M)
R = np.transpose(L)

```

I do not know if the last two lines gives you the covariance matrix R. I supposed that R should be upper triangular? The Cholesky-decomposition gives a lower triangular matrix L such that L L^T = M. Define R := L^T, which is upper triangular, then we have R^T R = L^{T^T} L^T = L L^T = M.

---

<div class="post-metadata">

**Author:** ![erja](https://yyz2.discourse-cdn.com/free1/user_avatar/fenicsproject.discourse.group/erja/32/469_2.png) [@erja](https://fenicsproject.discourse.group/u/erja)\
**Post date:** [September 10, 2019, 9:26am UTC](https://fenicsproject.discourse.group/t/generating-noise-in-particular-way/1505/3 "2019-09-10T09:26:05Z")

</div>

Thank you for your reply, I appreciate it a lot!

I am sorry, I am new to fenics (not helped by the complete lack of documentation), but I don’t see how that random generation would generate according to my formula above… it is quite important that it takes a sum of normal random coefficients times an eigenbasis function, that is to say

\sum\_{k=1}^{N\_h} \xi\_k e\_k

where then e\_k is the eigenbasis.  
Furhtermore, I want to use this as a RHS in an elliptic equation. I can’t get this to a form that fenics accepts or understand.

---

<div class="post-metadata">

**Author:** ![Yannik](https://avatars.discourse-cdn.com/v4/letter/y/5fc32e/32.png) [@Yannik](https://fenicsproject.discourse.group/u/Yannik)\
**Post date:** [September 10, 2019, 9:55am UTC](https://fenicsproject.discourse.group/t/generating-noise-in-particular-way/1505/4 "2019-09-10T09:55:17Z")

</div>

Ok. So let u = \sum\_{k = 0}^{N\_h} u\_k e\_k, where e\_k is your basis of V\_h. With `u = Function(V)` I get exactly such function in FEniCS. Here `u.vector()[0]` corresponds to u\_0, `u.vector()[5]` to u\_5 and so on. With `u.vector()[:]` I get an array U = [u\_0, u\_1, \ldots, u\_{N\_h}]. What I want to do now is to set u\_k = \xi\_k for k = 0, \ldots, N\_h. Thats why I am generating a random-array `xi`. Now I have to assign the nodal values of u with the of \xi. Translating this is `u.vector()[k] = xi[k]` for all `k`. Or, in short notation `u.vector()[:] = xi`. Here, I assign an array with an array.  
Hope that helps. Otherwise I have not understood your question.

---

<div class="post-metadata">

**Author:** ![erja](https://yyz2.discourse-cdn.com/free1/user_avatar/fenicsproject.discourse.group/erja/32/469_2.png) [@erja](https://fenicsproject.discourse.group/u/erja)\
**Post date:** [September 10, 2019, 10:18am UTC](https://fenicsproject.discourse.group/t/generating-noise-in-particular-way/1505/5 "2019-09-10T10:18:38Z")

</div>

Okey! Great! I think I get it now! Just one question, is u=Function(V) any basis element or is it an eigenbasis?

---

<div class="post-metadata">

**Author:** ![Yannik](https://avatars.discourse-cdn.com/v4/letter/y/5fc32e/32.png) [@Yannik](https://fenicsproject.discourse.group/u/Yannik)\
**Post date:** [September 10, 2019, 10:24am UTC](https://fenicsproject.discourse.group/t/generating-noise-in-particular-way/1505/6 "2019-09-10T10:24:21Z")

</div>

None of both. V is your functionspace of choice (P1, P2, DG, …), `u = Function(V)` is a Function from the functionspace V, which is 0 everywhere.

---

<div class="post-metadata">

**Author:** ![erja](https://yyz2.discourse-cdn.com/free1/user_avatar/fenicsproject.discourse.group/erja/32/469_2.png) [@erja](https://fenicsproject.discourse.group/u/erja)\
**Post date:** [September 10, 2019, 10:26am UTC](https://fenicsproject.discourse.group/t/generating-noise-in-particular-way/1505/7 "2019-09-10T10:26:02Z")

</div>

So what is u.vector()[0] then?

---

<div class="post-metadata">

**Author:** ![Yannik](https://avatars.discourse-cdn.com/v4/letter/y/5fc32e/32.png) [@Yannik](https://fenicsproject.discourse.group/u/Yannik)\
**Post date:** [September 10, 2019, 12:04pm UTC](https://fenicsproject.discourse.group/t/generating-noise-in-particular-way/1505/8 "2019-09-10T12:04:06Z")

</div>

The value of u corresponding to the 0-th basisfunction of the functionspace V.