June 2016 - firedrake - imperial.ac.uk

Re: [firedrake] Saving txt files while running in parallel
by Andrew McRae 30 Aug '16

30 Aug '16

Yes, I guess that in Floriane's original code, each process is computing its own max_eta (since eta.dat.data contains the values for degrees of freedom on that MPI process). Then all process are trying to write (or append?) to the same file(?). You could get the separate max values by making the filename depend on op2.MPI.comm.rank, so that each process writes to a separate file. Otherwise, computing the max over all subdomains indeed requires a parallel operation. On 8 June 2016 at 09:53, Shipton, Jemma <j.shipton(a)imperial.ac.uk> wrote: > Hi Floriane, > > > We have defined a max method that works in parallel for our code... try > something like: > > > def max(f): fmax = op2.Global(1, [-1000], dtype=float) op2.par_loop(op2.Kernel("""void maxify(double *a, double *b) { a[0] = a[0] < fabs(b[0]) ? fabs(b[0]) : a[0]; }""", "maxify"), f.dof_dset.set, fmax(op2.MAX), f.dat(op2.READ)) return fmax.data[0] > > > But I agree, it would be nice if that weren't necessary. > > > Hope that helps, > > > Jemma > > ------------------------------ > *From:* firedrake-bounces(a)imperial.ac.uk <firedrake-bounces(a)imperial.ac.uk> > on behalf of Floriane Gidel [RPG] <mmfg(a)leeds.ac.uk> > *Sent:* 08 June 2016 09:33:46 > *To:* firedrake > *Subject:* [firedrake] Saving txt files while running in parallel > > > Dear all, > > > I am running a Firedrake code in parallel with 4 cores, in which I save > the maximal value of the amplitude of the wave at each time, with the > command: > > > max_eta = max(eta.dat.data) > > Eta_file.write('%-10s %-10s %-10s\n' % (t,'', eta_max)) > > > > But when opening the .txt file, I notice two issues: > > - the data are not saved continuously: time goes from 0 to 100 (let's say) > and then starts from 90 again; > > - the values of max_eta are not the maximal values of eta, except at the > beginning. So it looks like the function max takes the maximum value of eta > only in one subdomain, so that after the wave has crossed the subdomain, > the maximal value of eta goes back to the depth of rest value. > > > How can I force this command to be applied on the full domain, even if the > code is run in parallel? > > > Thanks, > > Floriane > > > >

5 18

Re: [firedrake] problems updating and running firedrake
by David Ham 05 Aug '16

05 Aug '16

Hi Francis, Please provide the full output, because the following is based on guesswork because I can't see what is going on. It is true that curl depends on librtmp. If your curl can't find your rtmp, then something is broken on your system. This has nothing to do with Firedrake. I *think* that the only time you use curl is to download the firedrake-install script. If you really can't unbreak curl on your system then you can obtain this script by any of the usual ways (wget, download it with a web browser etc.) Regards, David On Mon, 27 Jun 2016 at 19:35 Francis Poulin <fpoulin(a)uwaterloo.ca> wrote: > Hello David, > > Thanks for finding that. > > I tried a) and that didn't help the update. > > I then tried b) and get an curl error, "error while loading shared > libraries: librtmp.so.0". > > I have done this before but not sure why this is not working now. > > Cheers, Francis > > > > > ------------------ > Francis Poulin > Associate Professor > Department of Applied Mathematics > University of Waterloo > > email: fpoulin(a)uwaterloo.ca > Web: https://uwaterloo.ca/poulin-research-group/ > Telephone: +1 519 888 4567 x32637 > > ------------------------------ > *From:* firedrake-bounces(a)imperial.ac.uk [firedrake-bounces(a)imperial.ac.uk] > on behalf of David Ham [David.Ham(a)imperial.ac.uk] > *Sent:* Monday, June 27, 2016 10:48 AM > *To:* firedrake > *Subject:* Re: [firedrake] problems updating and running firedrake > > Hi Francis, > > It looks like you might have hit the issue they describe here: > http://askubuntu.com/questions/509283/python-no-module-named-datetime > > Basically datetime was apparently rolled into the python binary in a 2.7 > point release. I don't know why this was changed but it makes the python > binary in your virtualenv incompatible with your system python install. > > According to the thread on askubuntu, you have two options. > > a) Apparently you can re-initialise the virtualenv by running: > > virtualenv . > > in your top Firedrake directory. I'm guessing you should do this with the > virtualenv deactivated. > > or > > b) A fresh Firedrake install should work because it will create a fresh > virtualenv based on current Python binary. > > Regards, > > David > > On Mon, 27 Jun 2016 at 15:20 Francis Poulin <fpoulin(a)uwaterloo.ca> wrote: > >> Hello, >> >> I recently updated my ubuntu from 14.04 to 16.04. I tried running >> firedrake and it failed. >> >> I also tried updating firedrake and it failed, see below. >> >> Do I need to install it from scratch or should it be possible to update >> it with a new OS? >> >> Cheers, Francis >> >> >> >> Installing pip dependencies for COFFEE >> /home/fpoulin/software/firedrake/bin/pip install -r >> COFFEE/requirements.txt >> >> Traceback (most recent call last): >> File "/home/fpoulin/software/firedrake/bin/pip", line 7, in <module> >> from pip import main >> File >> "/home/fpoulin/software/firedrake/local/lib/python2.7/site-packages/pip/__init__.py", >> line 14, in <module> >> from pip.utils import get_installed_distributions, get_prog >> File >> "/home/fpoulin/software/firedrake/local/lib/python2.7/site-packages/pip/utils/__init__.py", >> line 27, in <module> >> from pip._vendor import pkg_resources >> File >> "/home/fpoulin/software/firedrake/local/lib/python2.7/site-packages/pip/_vendor/pkg_resources/__init__.py", >> line 36, in <module> >> import plistlib >> File "/usr/lib/python2.7/plistlib.py", line 62, in <module> >> import datetime >> ImportError: No module named datetime >> Traceback (most recent call last): >> File "../bin/firedrake-update", line 974, in <module> >> pip_requirements(p) >> File "../bin/firedrake-update", line 377, in pip_requirements >> run_pip_install(["-r", "%s/requirements.txt" % package]) >> File "../bin/firedrake-update", line 364, in run_pip_install >> check_call(pipinstall + pipargs) >> File "../bin/firedrake-update", line 211, in check_call >> subprocess.check_call(arguments, env=env) >> File "/usr/lib/python2.7/subprocess.py", line 541, in check_call >> raise CalledProcessError(retcode, cmd) >> subprocess.CalledProcessError: Command >> '['/home/fpoulin/software/firedrake/bin/pip', 'install', '-r', >> 'COFFEE/requirements.txt']' returned non-zero exit status 1 >> >> ------------------ >> Francis Poulin >> Associate Professor >> Department of Applied Mathematics >> University of Waterloo >> >> email: fpoulin(a)uwaterloo.ca >> Web: https://uwaterloo.ca/poulin-research-group/ >> Telephone: +1 519 888 4567 x32637 >> >>

3 3

optimisation of the running time
by Floriane Gidel [RPG] 12 Jul '16

12 Jul '16

Dear all, I would like some advice to decrease the running time in my code 3D_NL.py . Details are given in the attached pdf and below is a quick summary: In the attached files, I run a 3D nonlinear code for potential flow water wave equations. The vertical and horizontal space discretisation are separated, so that Firedrake solves the equations in the horizontal plane only, and the vertical dependency is involved through matrices computed by evaluating the z-dependent integrals with the functions defined in the file 'discr_sigma.py'. Basically, the potential flow phi(x,y,z,t) is defined as phi(x,y,t)phi(z), where phi(z) is discretised as a n^th order Lagrange expansion in one element: phi(z_i) = product_{k=0}^n (z-z_k)/(z_i-z_k) for k!= i There are therefore 3 unknowns, which are : - h (the total depth of water) defined in the horizontal plane - phi_s (the potential flow at the surface) defined in the horizontal plane - hat_psi a vector function containing all the phi_i for i = 1:n , which are the values of phi in the vertical nodes (excluding the surface ones), where n is the order of the Lagrange expansion. Each phi_i is defined in the horizontal plane. The time scheme is defined like this : step 1 : update phi_s and hat_psi simultaneously (that is, sum of n+1 weak formulations) step 2 : update h and hat_psi simultaneously (that is, sum of n+1 weak formulations) step 3 : update phi_s (one weak formulation). (those steps are written in the pdf file as well) At the moment, the equations are solved in a small domain with resolution 0.1, but it's running very slowly. Running in parallel does not improve it since the number of nodes is too small. Is there another way to increase the simulation speed? Maybe by changing the iteration criteria in the solver ? I hope this is clear enough, otherwise I can send you more details. Best regards, Floriane

3 5

Assembling Jacobian within TAO
by Justin Chang 30 Jun '16

30 Jun '16

Hi all, So I am attempting to use TAO's mixed complementarity solvers and it basically requires a formGradient/Constraints routine and a formJacobian routine. I am solving a variational inequality via TAO and want the linear operator A to be my jacobian. However, I don't think it was correctly done because my semi smooth linear solver does not converge (but it should!). Attached is the code. Line 108 is what TaoSetJacobianRoutine() calls. Basically I want to assemble the linear operator into the matrix petsc_J. Is what I have the correct implementation? run the code as: > python MCP_poisson.py 100 100 2 Thanks, Justin

2 2

Re: [firedrake] FEniCS implementation works but Firedrake does not, why?
by Fabio Luporini 29 Jun '16

29 Jun '16

Hi Justin, I should be able to come up with a non-verbose mode in the next few days, thanks -- Fabio 2016-06-01 19:12 GMT+02:00 Justin Chang <jychang48(a)gmail.com>: > A couple more related questions: > > 1) Every time I run this (or any firedrake code) I get these notifications: > > COFFEE finished in XXX seconds (flops: X -> X) > > Is there a way to disable these? > > 2) Also, if I use quadrilateral elements, i also get this notification (I > believe it happens during the errornorm() call) > > Discontinuous Lagrange element requested on quadrilateral, creating DQ > element. > > Is there a way to also disable this? > > Thanks, > Justin > > On Wed, Jun 1, 2016 at 11:56 AM, Justin Chang <jychang48(a)gmail.com> wrote: > >> Ah this works wonderfully, thanks Lawrence! >> >> On Wed, Jun 1, 2016 at 4:23 AM, Lawrence Mitchell < >> lawrence.mitchell(a)imperial.ac.uk> wrote: >> >>> Hi Justin, >>> >>> On 01/06/16 08:44, Justin Chang wrote: >>> >>> ... >>> >>> > #===================== >>> > # Medium properties >>> > #===================== >>> > D = Function(Q) >>> > alpha = Constant(0.0) >>> > d1 = 1.0 >>> > d2 = 0.0001 >>> > theta = pi/6.0 >>> > co = cos(theta) >>> > si = sin(theta) >>> > >>> D.interpolate(Expression(("co*co*d1+si*si*d2","-co*si*(d1-d2)","-co*si*(d1-d2)","si*si*d1+co*co*d2"),co=co,si=si,d1=d1,d2=d2)) >>> >>> >>> FWIW, I recommend moving to using interpolation of UFL expressions for >>> these kind of initialisations (you get, amongst other things) better >>> error messages when mistyping something. See >>> http://firedrakeproject.org/interpolation.html#ufl-expressions >>> >>> > #===================== >>> > # Volumetric source >>> > #===================== >>> > def f(u): >>> > return Constant(10)*u*u*u >>> > >>> > def df(u): >>> > return Constant(30)*u*u >>> > >>> > #===================== >>> > # Variational form >>> > #===================== >>> > a = alpha * inner(u,v) * dx + inner(D*nabla_grad(u), nabla_grad(v))*dx >>> + inner(df(u_k) * u,v)*dx >>> >>> So a is a bilinear form, but it depends on the value of the previous >>> solution. >>> >>> > A = assemble(a,bcs=bcs) >>> > solver = LinearSolver(A,options_prefix="solver_") >>> >>> Here you assemble A and build a LinearSolver. >>> >>> > #======================== >>> > # Optimization >>> > #======================== >>> > lb = Function(V) >>> > ub = Function(V) >>> > ub.assign(1000) >>> > taoSolver = PETSc.TAO().create(PETSc.COMM_WORLD) >>> > with lb.dat.vec_ro as lb_vec, ub.dat.vec as ub_vec: >>> > taoSolver.setVariableBounds(lb_vec,ub_vec) >>> > >>> > def ObjGrad(tao, petsc_x, petsc_g): >>> > with b.dat.vec as b_vec: >>> > A.M.handle.mult(petsc_x,petsc_g) >>> > xtHx = petsc_x.dot(petsc_g) >>> > xtf = petsc_x.dot(b_vec) >>> > petsc_g.axpy(-1.0,b_vec) >>> > return 0.5*xtHx - xtf >>> > >>> > taoSolver.setObjectiveGradient(ObjGrad) >>> > taoSolver.setType(PETSc.TAO.Type.BLMVM) >>> > >>> > #======================== >>> > # Picard iteration >>> > #======================== >>> > while it < maxit and eps > tol: >>> > A = assemble(a,bcs=bcs) >>> >>> Here you reassemble a and assign it to a new variable "A". >>> >>> > # Standard solver >>> > if opt == 0: >>> > b = assemble(L,bcs=bcs) >>> > solver.solve(u_k1,b) >>> >>> This is never seen inside the solver, so this solver always runs with >>> the originally assembled A (with the initial value for u_k). >>> >>> > # Optimization solver >>> > else: >>> > b = assemble(L) >>> > b -= rhs_opt >>> > for bc in bcs: >>> > bc.apply(b) >>> > with u_k1.dat.vec as u_k1_vec: >>> > taoSolver.solve(u_k1_vec) >>> g >>> I'm actually slightly surprised that this one works, but that appears >>> to be because I don't understand how python's scoping rules work. >>> >>> The LinearSolver object isn't designed to work if your jacobian >>> changes hence the problem. >>> >>> I recommend you use a LinearVariationalSolver, but you need to be >>> careful, since by default it doesn't rebuild the jacobian every solve. >>> >>> If I do: >>> >>> problem = LinearVariationalProblem(a, L, u_k1, bcs=bcs, >>> constant_jacobian=False) >>> solver = LinearVariationalSolver(problem, options_prefix="solver_") >>> >>> if opt == 0: >>> solver.solve() >>> >>> >>> Then I get: >>> >>> python foo.py 50 50 0 >>> >>> Error norm: 2.193e-01 >>> Error norm: 1.705e-02 >>> Error norm: 5.105e-04 >>> Error norm: 4.980e-07 >>> Error norm: 5.332e-13 >>> >>> We could add the option to rebuild the jacobian inside a LinearSolver >>> (maybe if you've called something like solver.reset_operator() first), >>> but it's not available right now. >>> >>> I've attached a modified version that gets quadratic convergence in >>> both cases. I've also modified the taoSolver setup slightly so that >>> you reuse the same storage for the operator and RHS (rather than >>> building a new one every time). >>> >>> Cheers, >>> >>> Lawrence >>> >>> _______________________________________________ >>> firedrake mailing list >>> firedrake(a)imperial.ac.uk >>> https://mailman.ic.ac.uk/mailman/listinfo/firedrake >>> >>> >> >

3 4

Re: [firedrake] Stochastic differential equations
by Colin Cotter 28 Jun '16

28 Jun '16

Hi Francis, The model is described in this paper: https://arxiv.org/abs/1410.8311 all the best --cjc On 27 June 2016 at 14:45, Francis Poulin <fpoulin(a)uwaterloo.ca> wrote: > Hello Colin, > > I would be curious to hear more about the stochastic QG model. Should I > ask Ali about this? > > Cheers, Francis > > > > ------------------ > Francis Poulin > Associate Professor > Department of Applied Mathematics > University of Waterloo > > email: fpoulin(a)uwaterloo.ca > Web: https://uwaterloo.ca/poulin-research-group/ > Telephone: +1 519 888 4567 x32637 > > ------------------------------ > *From:* firedrake-bounces(a)imperial.ac.uk [firedrake-bounces(a)imperial.ac.uk] > on behalf of Colin Cotter [colin.cotter(a)imperial.ac.uk] > *Sent:* Monday, June 27, 2016 2:12 AM > *To:* firedrake > *Subject:* Re: [firedrake] Stochastic differential equations > > Hi! > > Yes, you can just use the numpy random number generator to write to the > basis coefficients for this. > > from numpy import random > > F = Function(V) > F.dat.data[:] = random.randn(N) > > etc. > > The main thing you have to be careful of is MPI parallelism. > > We have been using this to implement a stochastic QG equation for example. > > It might be nice to integrate random number generators into TSFL generated > code so that we can use them directly in UFL expressions - that would tidy > up a lot of our code (and help make it more parallel safe). > > cheers > --cjc > > On 26 June 2016 at 21:27, Angwenyi David <kipkoej(a)gmail.com> wrote: > >> Thanks, Francis. >> >> True, the DE is driven by white noise, that is, I assume it is gaussian >> with zero mean vector and identity covariance matrix premultiplied by the >> square root of the time-step. >> >> On 26 Jun 2016, at 22:16, Francis Poulin <fpoulin(a)uwaterloo.ca> wrote: >> >> Hello, >> >> I can't say whether there is anything specific to firedrake but I can say >> that in general it depends on what kind of noise you are considering in the >> SDE. >> >> If you are using coloured noise, then you treat your forcing term as you >> would any other term. However, given the fact that you are using W I >> imagine you are using White noise. >> >> In that case the popular thing to do is to use a splitting method. For >> the "regular" deterministic terms, you use whatever time stepping scheme >> you normally would use. However, for the white noise you probably want to >> use an Euler step. Basically, that means you have square root of dt in >> front of the random variable. The stochastic part is only first order >> accurate but it's much easier than using other techniques, but they do >> exist. >> >> Anyhow, that's my two cents worth. >> >> Francis >> >> ------------------ >> Francis Poulin >> Associate Professor >> Department of Applied Mathematics >> University of Waterloo >> >> email: fpoulin(a)uwaterloo.ca >> Web: https://uwaterloo.ca/poulin-research-group/ >> Telephone: +1 519 888 4567 x32637 >> >> ------------------------------ >> *From:* firedrake-bounces(a)imperial.ac.uk [ >> firedrake-bounces(a)imperial.ac.uk] on behalf of Angwenyi David [ >> kipkoej(a)gmail.com] >> *Sent:* Sunday, June 26, 2016 4:11 PM >> *To:* firedrake(a)imperial.ac.uk >> *Subject:* [firedrake] Stochastic differential equations >> >> Dear all, >> >> Is it possible to solve pdes driven by some noise using Firedrake? For >> example, see equation 3.2 in the attachment whose discretisation yields a >> system 3.3 and 3.4, where dW(t):=sqrt(dt)N(0,Identity). My other related >> concern is whether It is possible to implement, say, A Kalman filter with >> equation 3.2 as model equation on Firedrake. >> >> Thanks in advance. >> >> Angwenyi David >> <Screen Shot 2016-06-26 at 21.51.52.png> >> _______________________________________________ >> firedrake mailing list >> firedrake(a)imperial.ac.uk >> https://mailman.ic.ac.uk/mailman/listinfo/firedrake >> >> >> > > > -- > http://www.imperial.ac.uk/people/colin.cotter > > www.cambridge.org/9781107663916 > > > -- http://www.imperial.ac.uk/people/colin.cotter www.cambridge.org/9781107663916

1 0

Re: [firedrake] problems updating and running firedrake
by David Ham 27 Jun '16

27 Jun '16

Hi Francis, It looks like you might have hit the issue they describe here: http://askubuntu.com/questions/509283/python-no-module-named-datetime Basically datetime was apparently rolled into the python binary in a 2.7 point release. I don't know why this was changed but it makes the python binary in your virtualenv incompatible with your system python install. According to the thread on askubuntu, you have two options. a) Apparently you can re-initialise the virtualenv by running: virtualenv . in your top Firedrake directory. I'm guessing you should do this with the virtualenv deactivated. or b) A fresh Firedrake install should work because it will create a fresh virtualenv based on current Python binary. Regards, David On Mon, 27 Jun 2016 at 15:20 Francis Poulin <fpoulin(a)uwaterloo.ca> wrote: > Hello, > > I recently updated my ubuntu from 14.04 to 16.04. I tried running > firedrake and it failed. > > I also tried updating firedrake and it failed, see below. > > Do I need to install it from scratch or should it be possible to update it > with a new OS? > > Cheers, Francis > > > > Installing pip dependencies for COFFEE > /home/fpoulin/software/firedrake/bin/pip install -r COFFEE/requirements.txt > > Traceback (most recent call last): > File "/home/fpoulin/software/firedrake/bin/pip", line 7, in <module> > from pip import main > File > "/home/fpoulin/software/firedrake/local/lib/python2.7/site-packages/pip/__init__.py", > line 14, in <module> > from pip.utils import get_installed_distributions, get_prog > File > "/home/fpoulin/software/firedrake/local/lib/python2.7/site-packages/pip/utils/__init__.py", > line 27, in <module> > from pip._vendor import pkg_resources > File > "/home/fpoulin/software/firedrake/local/lib/python2.7/site-packages/pip/_vendor/pkg_resources/__init__.py", > line 36, in <module> > import plistlib > File "/usr/lib/python2.7/plistlib.py", line 62, in <module> > import datetime > ImportError: No module named datetime > Traceback (most recent call last): > File "../bin/firedrake-update", line 974, in <module> > pip_requirements(p) > File "../bin/firedrake-update", line 377, in pip_requirements > run_pip_install(["-r", "%s/requirements.txt" % package]) > File "../bin/firedrake-update", line 364, in run_pip_install > check_call(pipinstall + pipargs) > File "../bin/firedrake-update", line 211, in check_call > subprocess.check_call(arguments, env=env) > File "/usr/lib/python2.7/subprocess.py", line 541, in check_call > raise CalledProcessError(retcode, cmd) > subprocess.CalledProcessError: Command > '['/home/fpoulin/software/firedrake/bin/pip', 'install', '-r', > 'COFFEE/requirements.txt']' returned non-zero exit status 1 > > ------------------ > Francis Poulin > Associate Professor > Department of Applied Mathematics > University of Waterloo > > email: fpoulin(a)uwaterloo.ca > Web: https://uwaterloo.ca/poulin-research-group/ > Telephone: +1 519 888 4567 x32637 > >

2 1

problems updating and running firedrake
by Francis Poulin 27 Jun '16

27 Jun '16

Hello, I recently updated my ubuntu from 14.04 to 16.04. I tried running firedrake and it failed. I also tried updating firedrake and it failed, see below. Do I need to install it from scratch or should it be possible to update it with a new OS? Cheers, Francis Installing pip dependencies for COFFEE /home/fpoulin/software/firedrake/bin/pip install -r COFFEE/requirements.txt Traceback (most recent call last): File "/home/fpoulin/software/firedrake/bin/pip", line 7, in <module> from pip import main File "/home/fpoulin/software/firedrake/local/lib/python2.7/site-packages/pip/__init__.py", line 14, in <module> from pip.utils import get_installed_distributions, get_prog File "/home/fpoulin/software/firedrake/local/lib/python2.7/site-packages/pip/utils/__init__.py", line 27, in <module> from pip._vendor import pkg_resources File "/home/fpoulin/software/firedrake/local/lib/python2.7/site-packages/pip/_vendor/pkg_resources/__init__.py", line 36, in <module> import plistlib File "/usr/lib/python2.7/plistlib.py", line 62, in <module> import datetime ImportError: No module named datetime Traceback (most recent call last): File "../bin/firedrake-update", line 974, in <module> pip_requirements(p) File "../bin/firedrake-update", line 377, in pip_requirements run_pip_install(["-r", "%s/requirements.txt" % package]) File "../bin/firedrake-update", line 364, in run_pip_install check_call(pipinstall + pipargs) File "../bin/firedrake-update", line 211, in check_call subprocess.check_call(arguments, env=env) File "/usr/lib/python2.7/subprocess.py", line 541, in check_call raise CalledProcessError(retcode, cmd) subprocess.CalledProcessError: Command '['/home/fpoulin/software/firedrake/bin/pip', 'install', '-r', 'COFFEE/requirements.txt']' returned non-zero exit status 1 ------------------ Francis Poulin Associate Professor Department of Applied Mathematics University of Waterloo email: fpoulin(a)uwaterloo.ca Web: https://uwaterloo.ca/poulin-research-group/ Telephone: +1 519 888 4567 x32637

1 0

Re: [firedrake] Stochastic differential equations
by Colin Cotter 27 Jun '16

27 Jun '16

Hi! Yes, you can just use the numpy random number generator to write to the basis coefficients for this. from numpy import random F = Function(V) F.dat.data[:] = random.randn(N) etc. The main thing you have to be careful of is MPI parallelism. We have been using this to implement a stochastic QG equation for example. It might be nice to integrate random number generators into TSFL generated code so that we can use them directly in UFL expressions - that would tidy up a lot of our code (and help make it more parallel safe). cheers --cjc On 26 June 2016 at 21:27, Angwenyi David <kipkoej(a)gmail.com> wrote: > Thanks, Francis. > > True, the DE is driven by white noise, that is, I assume it is gaussian > with zero mean vector and identity covariance matrix premultiplied by the > square root of the time-step. > > On 26 Jun 2016, at 22:16, Francis Poulin <fpoulin(a)uwaterloo.ca> wrote: > > Hello, > > I can't say whether there is anything specific to firedrake but I can say > that in general it depends on what kind of noise you are considering in the > SDE. > > If you are using coloured noise, then you treat your forcing term as you > would any other term. However, given the fact that you are using W I > imagine you are using White noise. > > In that case the popular thing to do is to use a splitting method. For > the "regular" deterministic terms, you use whatever time stepping scheme > you normally would use. However, for the white noise you probably want to > use an Euler step. Basically, that means you have square root of dt in > front of the random variable. The stochastic part is only first order > accurate but it's much easier than using other techniques, but they do > exist. > > Anyhow, that's my two cents worth. > > Francis > > ------------------ > Francis Poulin > Associate Professor > Department of Applied Mathematics > University of Waterloo > > email: fpoulin(a)uwaterloo.ca > Web: https://uwaterloo.ca/poulin-research-group/ > Telephone: +1 519 888 4567 x32637 > > ------------------------------ > *From:* firedrake-bounces(a)imperial.ac.uk [firedrake-bounces(a)imperial.ac.uk] > on behalf of Angwenyi David [kipkoej(a)gmail.com] > *Sent:* Sunday, June 26, 2016 4:11 PM > *To:* firedrake(a)imperial.ac.uk > *Subject:* [firedrake] Stochastic differential equations > > Dear all, > > Is it possible to solve pdes driven by some noise using Firedrake? For > example, see equation 3.2 in the attachment whose discretisation yields a > system 3.3 and 3.4, where dW(t):=sqrt(dt)N(0,Identity). My other related > concern is whether It is possible to implement, say, A Kalman filter with > equation 3.2 as model equation on Firedrake. > > Thanks in advance. > > Angwenyi David > <Screen Shot 2016-06-26 at 21.51.52.png> > _______________________________________________ > firedrake mailing list > firedrake(a)imperial.ac.uk > https://mailman.ic.ac.uk/mailman/listinfo/firedrake > > > -- http://www.imperial.ac.uk/people/colin.cotter www.cambridge.org/9781107663916

3 3

Re: [firedrake] Stochastic differential equations
by Colin Cotter 27 Jun '16

27 Jun '16

Dear Angwenyi, For an ensemble filtering approach, please can you talk to my student Ali Gregory about how he did it? http://www.imperial.ac.uk/people/a.gregory14 If you are working with an ensemble, then one approach is to use a list of Firedrake fields for the different states. Then, the Kalman update needs to be implemented through a PyOP2 parloop. Ali can help you with this (he has done it for the ETPF). all the best --cjc On 27 June 2016 at 12:16, Angwenyi David <kipkoej(a)gmail.com> wrote: > The numpy random number generator works! Thanks Colin. > > Angwenyi > > On 27 Jun 2016, at 09:16, Angwenyi David <kipkoej(a)gmail.com> wrote: > > Thanks, Collins. > > Could you please, if it won’t matter to you, send me the code for the > stochastic quasi geostrophic equation you implemented. It shall help me a > great deal. > > How about implementing a Kalman filter with the system of equations 3.3 > and 3.4 as a model equations with given data? I presume it is possible; for > the system can be represented as psi^(n+1)=A^(n)*psi^(n)+noise. The > complication I see is having the spatial discretisation of phi fit into the > matrix A. Links to the articles with similar problems are appreciated. > > Thanks. > > Angwenyi > > On 27 Jun 2016, at 08:12, Colin Cotter <colin.cotter(a)imperial.ac.uk> > wrote: > > Hi! > > Yes, you can just use the numpy random number generator to write to the > basis coefficients for this. > > from numpy import random > > F = Function(V) > F.dat.data[:] = random.randn(N) > > etc. > > The main thing you have to be careful of is MPI parallelism. > > We have been using this to implement a stochastic QG equation for example. > > It might be nice to integrate random number generators into TSFL generated > code so that we can use them directly in UFL expressions - that would tidy > up a lot of our code (and help make it more parallel safe). > > cheers > --cjc > > On 26 June 2016 at 21:27, Angwenyi David <kipkoej(a)gmail.com> wrote: > >> Thanks, Francis. >> >> True, the DE is driven by white noise, that is, I assume it is gaussian >> with zero mean vector and identity covariance matrix premultiplied by the >> square root of the time-step. >> >> On 26 Jun 2016, at 22:16, Francis Poulin <fpoulin(a)uwaterloo.ca> wrote: >> >> Hello, >> >> I can't say whether there is anything specific to firedrake but I can say >> that in general it depends on what kind of noise you are considering in the >> SDE. >> >> If you are using coloured noise, then you treat your forcing term as you >> would any other term. However, given the fact that you are using W I >> imagine you are using White noise. >> >> In that case the popular thing to do is to use a splitting method. For >> the "regular" deterministic terms, you use whatever time stepping scheme >> you normally would use. However, for the white noise you probably want to >> use an Euler step. Basically, that means you have square root of dt in >> front of the random variable. The stochastic part is only first order >> accurate but it's much easier than using other techniques, but they do >> exist. >> >> Anyhow, that's my two cents worth. >> >> Francis >> >> ------------------ >> Francis Poulin >> Associate Professor >> Department of Applied Mathematics >> University of Waterloo >> >> email: fpoulin(a)uwaterloo.ca >> Web: https://uwaterloo.ca/poulin-research-group/ >> Telephone: +1 519 888 4567 x32637 >> >> ------------------------------ >> *From:* firedrake-bounces(a)imperial.ac.uk [ >> firedrake-bounces(a)imperial.ac.uk] on behalf of Angwenyi David [ >> kipkoej(a)gmail.com] >> *Sent:* Sunday, June 26, 2016 4:11 PM >> *To:* firedrake(a)imperial.ac.uk >> *Subject:* [firedrake] Stochastic differential equations >> >> Dear all, >> >> Is it possible to solve pdes driven by some noise using Firedrake? For >> example, see equation 3.2 in the attachment whose discretisation yields a >> system 3.3 and 3.4, where dW(t):=sqrt(dt)N(0,Identity). My other related >> concern is whether It is possible to implement, say, A Kalman filter with >> equation 3.2 as model equation on Firedrake. >> >> Thanks in advance. >> >> Angwenyi David >> <Screen Shot 2016-06-26 at 21.51.52.png> >> _______________________________________________ >> firedrake mailing list >> firedrake(a)imperial.ac.uk >> https://mailman.ic.ac.uk/mailman/listinfo/firedrake >> >> >> > > > -- > http://www.imperial.ac.uk/people/colin.cotter > > www.cambridge.org/9781107663916 > > > _______________________________________________ > firedrake mailing list > firedrake(a)imperial.ac.uk > https://mailman.ic.ac.uk/mailman/listinfo/firedrake > > > > -- http://www.imperial.ac.uk/people/colin.cotter www.cambridge.org/9781107663916

2 1