Gaussian processes by stochastic differential equations on manifolds

 
FEniCS code:

from dolfin import *

import os

import numpy as np

from mshr import *

mesh = BoundaryMesh(generate_mesh(Sphere(Point(0,0,0),1),50),"exterior")

V = FunctionSpace(mesh, "CG", 1)

u = TrialFunction(V)

v = TestFunction(V)

k=inner(grad(u),grad(v))*dx

K=assemble(k)

K=K.array()

m=u*v*dx

M = assemble(m)

M=M.array()

alpha = 1.0

beta = 0.2

np.random.seed(12)

L = np.linalg.cholesky(M)

file_realisation = File('results/realisation_sphere.pvd')

func = Function(V)

func.rename('realisation','realisation')

for i in range(10):

    W=np.random.randn(M.shape[0],1)

    y = np.matmul(L,W)

    e = np.linalg.solve(M+(beta**2)*K,alpha*y)

    func.vector().set_local(e)

    file_realisation << func  

 

Generation on Cardiff Dinosaur, Eddy