Error during implementation of KrylovSolver

Question

Hi All,

The following use of Krylov Solver gives me this error:

Solving linear system of size 41473 x 41473 (PETSc Krylov solver).
PETSc Krylov solver starting to solve 41473 x 41473 system.
Traceback (most recent call last):
File "wave_newmark.py", line 124, in
solver.solve(e.vector(), b)
RuntimeError:
Error: Unable to solve linear system using PETSc Krylov solver.
* Reason: Solution failed to converge in 0 iterations (PETSc reason DIVERGED_NANORINF, residual norm ||r|| = nan).
* Where: This error was encountered inside PETScKrylovSolver.cpp.
* Process: 0

Can anyone please help in getting rid of this? Is it being implemented wrong?

a = some expression ; L = some expression
S  = SystemAssembler(a, L, bc)
A  = PETScMatrix()
b  = PETScVector()
S.assemble(A)

e = Function(V)

solver = PETScKrylovSolver("cg")
solver.set_operator(A) 

while t < T:
  S.assemble(b)  
  solver.solve(e.vector(), b) 
  e0.assign(e1) 
  e1.assign(e)
  t += dt

Answer 1

It seems that both e and b are zero in your case. This makes the
residual zero and any kind of relative convergence criteria will break.
Try to print the norm of e and b. Try also to assemble something into b.

Comments

Hi Kent-Andre,

Thank you very much for your reply. The above code is preceded by the following:

e = TrialFunction(V) 
v = TestFunction(V)  

a =   epsilon/dt/dt*inner(v, e)*dx \
    + nbeta*nu*inner(curl(v), curl(e))*dx

L =   2.0*epsilon/dt/dt*inner(v,e1)*dx \
    - epsilon/dt/dt*inner(v, e0)*dx \
    - (1.0 - 2.0*nbeta)*nu*inner(curl(v), curl(e1))*dx \
    - nbeta*nu*inner(curl(v), curl(e0))*dx 

Is it because I am defining 'e' twice? I am trying to assemble the vector from Systemassembler 'S' into 'b' at the beginning of the while loop (refer code). Should that not work?

---

Hey, do reply whenever you find time. Thanks!

---

If you define 'e' twice you should be careful. The first forms
are in terms of an 'e' that is later 'lost' in the namespace. As long
as you know what you are doing, you are fine.

Anyway, I would have assembled
S.assemble(A, b) before the time loop and printed
out the norm of b. If it is zero then the solution will be zero
which is not what you wanted. Also print out the norm of e.
In general, make sure all functions are initialized to something.
It is not clear to from your code what e1 and e0 are in the first solve. If they are zero
then there will be no fun, right.

---

Hey,

Thanks a lot for your reply!

No, my e0 and e1 are not zero. I tried substituting the Krylov solver with LUsolver and it works. I don't know what's going on, but with the exact same syntax, it worked!

Do you have any experience with modelling wave equations by any chance?

Thanks,
Chaitanya