Hello, Sorry for the bother but I am wondering if anyone had any advice on my problem. The code can be found at the following: https://github.com/francispoulin/firedrakeQG/blob/master/basin_modes_sw.py There are two problems: 1) I can't seem to plot the eigenfunctions. Maybe I didn't split them correctly? I'm not sure. 2) The eigenvalues are not at all the ones I'm looking for. I know from the fenics version that there is an eigenvalue with an imaginary part of 3.14. These might be valid but I suspect that they are very noisy and not physically interesting. Any thoughts or suggestions as to what i can try would be greatly appreciated. Cheers, Francis ------------------ Francis Poulin Associate Professor Department of Applied Mathematics University of Waterloo email: fpoulin@uwaterloo.ca Web: https://uwaterloo.ca/poulin-research-group/ Telephone: +1 519 888 4567 x32637 ________________________________________ From: firedrake-bounces@imperial.ac.uk [firedrake-bounces@imperial.ac.uk] on behalf of Francis Poulin [fpoulin@uwaterloo.ca] Sent: Friday, December 02, 2016 3:41 PM To: firedrake@imperial.ac.uk Subject: Re: [firedrake] eigenvalue problems with SLEPc Hello again, I played around a bit more with my code, and after reading the documented I have restructured it so that it is more similar to what is discussed in the docs. The updated code is on the same github repo. Now I seem to be having problems plotting the function. I beileve that I have managed (very inefficiently) to store the real and imaginary parts of the free-surface heights in a proper function. See below. em_real, em_imag = Function(Z), Function(Z) u_real, u_imag = Function(Vdg), Function(Vdg) eta_real, eta_imag = Function(Vcg), Function(Vcg) eval_real, eval_imag = [], [] # Find eigenvalues an eigenvectors vr, vi = petsc_a.getVecs() em_real, em_imag = MixedFunctionSpace(Z, vr), MixedFunctionSpace(Z, vr) for i in range(nconv): lam = es.getEigenpair(i, vr, vi) print lam eval_real.append(lam.real) eval_imag.append(lam.imag) u_real, eta_real = em_real.split() #u_imag, eta_imag = em_imag.split() p = plot(eta_real) p.show() When I try plotting as above I get the following error. Sorry if this is a silly mistake but I am hoping that if i can view the eigenfunction I can see whether it looks reasonable. Given the eigenvalues I don't believe it is since it differs from my calculated values in fenics by almost 200. Cheers, Francis Traceback (most recent call last): File "basin_modes_sw.py", line 85, in <module> p = plot(eta_real) File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/plot.py", line 125, in plot raise TypeError("Unexpected type") TypeError: Unexpected type ------------------ Francis Poulin Associate Professor Department of Applied Mathematics University of Waterloo email: fpoulin@uwaterloo.ca Web: https://uwaterloo.ca/poulin-research-group/ Telephone: +1 519 888 4567 x32637 ________________________________________ From: firedrake-bounces@imperial.ac.uk [firedrake-bounces@imperial.ac.uk] on behalf of Francis Poulin [fpoulin@uwaterloo.ca] Sent: Friday, December 02, 2016 7:37 AM To: firedrake@imperial.ac.uk Subject: Re: [firedrake] eigenvalue problems with SLEPc Thank you for both suggestions! I am happy to say that I now get a solution, not a correct one mind you, but a solution. After I solve the problem I need to split the solution to get the velocity or the free-surface vr, vi = petsc_a.getVecs() for i in range(nconv): lam = es.getEigenpair(i, vr, vi) print lam eigenvaluer.append(lam.real) eigenvaluei.append(lam.imag) eigenmodes_real.vector()[:], eigenmodes_imag.vector()[:] = vr, vi The last line is not correct and yields the error Traceback (most recent call last): File "basin_modes_sw.py", line 86, in <module> eigenmodes_real.vector()[:], eigenmodes_imag.vector()[:] = vr, vi File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/vector.py", line 182, in __setitem__ self.dat.data[idx] = value TypeError: 'tuple' object does not support item assignment Can someone tell me the correct syntax for this case? Cheers, Francis ------------------ Francis Poulin Associate Professor Department of Applied Mathematics University of Waterloo email: fpoulin@uwaterloo.ca Web: https://uwaterloo.ca/poulin-research-group/ Telephone: +1 519 888 4567 x32637 ________________________________________ From: firedrake-bounces@imperial.ac.uk [firedrake-bounces@imperial.ac.uk] on behalf of Lawrence Mitchell [lawrence.mitchell@imperial.ac.uk] Sent: Friday, December 02, 2016 5:05 AM To: firedrake@imperial.ac.uk Subject: Re: [firedrake] eigenvalue problems with SLEPc
On 2 Dec 2016, at 09:56, Stephan Kramer <s.kramer@imperial.ac.uk> wrote:
On 02/12/16 04:23, Francis Poulin wrote:
Hello David,
I have an undergraduate student working with me for another few weeks. I will get her to put together some demos and we will send them along as they progress.
Sorry to ask another question but since the recent success, we were inspired to try and do something similar for the Shallow Water model. This is more complicated because it not only has the Rossby basin modes but it also has gravity wave basin modes. We have a rough version of the code below,
https://github.com/francispoulin/firedrakeQG/blob/master/basin_modes_sw.py
At the moment we have set the Coriolis parameter to zero, so we are essentially solving a 2D wave equation but with a momentum equation and a continuity equation.
I realize that I need to learn SLEPc and understand all of the options that you need for convergence but the problem that I'm having now is more fundamental.
Traceback (most recent call last): File "basin_modes_sw.py", line 64, in <module> es.solve() File "SLEPc/EPS.pyx", line 1107, in slepc4py.SLEPc.EPS.solve (src/slepc4py.SLEPc.c:28270) petsc4py.PETSc.Error: error code 56 [0] EPSSolve() line 80 in /home/fpoulin/software/firedrake/src/slepc/src/eps/interface/epssolve.c [0] EPSSetUp() line 268 in /home/fpoulin/software/firedrake/src/slepc/src/eps/interface/epssetup.c [0] STSetUp() line 310 in /home/fpoulin/software/firedrake/src/slepc/src/sys/classes/st/interface/stsolve.c [0] STSetUp_Shift() line 139 in /home/fpoulin/software/firedrake/src/slepc/src/sys/classes/st/impls/shift/shift.c [0] KSPSetUp() line 393 in /tmp/pip-SPrpSR-build/src/ksp/ksp/interface/itfunc.c [0] PCSetUp() line 968 in /tmp/pip-SPrpSR-build/src/ksp/pc/interface/precon.c [0] PCSetUp_LU() line 92 in /tmp/pip-SPrpSR-build/src/ksp/pc/impls/factor/lu/lu.c [0] MatGetOrdering() line 260 in /tmp/pip-SPrpSR-build/src/mat/order/sorder.c [0] MatGetOrdering_ND() line 19 in /tmp/pip-SPrpSR-build/src/mat/order/spnd.c [0] No support for this operation for this object type [0] Cannot get rows for matrix type nest
I think you need mat_type='aij' as an extra argument to your assemble() call
Yes. Selecting matrix types finally has its own section in the manual: http://firedrakeproject.org/solving-interface.html#specifying-solution-metho... Lawrence _______________________________________________ firedrake mailing list firedrake@imperial.ac.uk https://mailman.ic.ac.uk/mailman/listinfo/firedrake _______________________________________________ firedrake mailing list firedrake@imperial.ac.uk https://mailman.ic.ac.uk/mailman/listinfo/firedrake _______________________________________________ firedrake mailing list firedrake@imperial.ac.uk https://mailman.ic.ac.uk/mailman/listinfo/firedrake