PETSs Error Code 63 when trying to do Nonlinear Problem

1

Hello,

I am new to fenics and am trying to do some 3d nonlinear finite element analysis on an equation for cancer modeling. I am running into an error that I can’t seem to find a fix for on the internet. I was hoping someone can help - it’s probably something easy that I just don’t understand. If a fix is found I would love some explanation on why it solves the problem (and what this error even means :smiley:).

Here is the error verbatim:
*** Error: Unable to successfully call PETSc function ‘MatSetValuesLocal’.
*** Reason: PETSc error code is: 63 (Argument out of range).
*** Where: This error was encountered inside /tmp/dolfin/dolfin/la/PETScMatrix.cpp.
*** Process: 0

And here is my code:

from __future__ import print_function
from fenics import *
import matplotlib.pyplot as plt
import time

T = 2.0            # final time
num_steps = 50     # number of time steps
dt = T / num_steps # time step size

# Create mesh and define function space
# mesh = RectangleMesh(Point(-2, -2), Point(2, 2), nx, ny)
mesh = UnitCubeMesh(16, 16, 16)
V = FunctionSpace(mesh, 'P', 1)

# Define boundary condition
# def boundary(x, on_boundary):
#     return on_boundary
#
# bc = DirichletBC(V, Constant(0), boundary)

# Define initial value
u_0 = Expression('exp(1 / (1 - pow(pow(x[0], 2) + pow(x[1], 2) + pow(x[2], 2), 2) ) )', degree=2)
u_n = interpolate(u_0, V)

# Define variational problem
u = Function(V)
v = TestFunction(V)

D = Constant((1/1000)/.18)
a = Constant(1/.18)

beta = Constant(0.9)
one  = Constant(1)

F = u*v*dx + dt*D*dot(grad(u), grad(v))*dx - (u_n + dt*a*(one - beta*u)*u)*v*dx

# Create VTK file for saving solution
vtkfile = File('CarTCellModel/solution.pvd')

# Time-stepping
u = Function(V)
t = 0
for n in range(num_steps):

    # Update current time
    t += dt

    # Compute solution
    solve(F == 0, u)

    # Save to file and plot solution
    vtkfile << (u, t)
    plot(u)

    # Update previous solution
    u_n.assign(u)

# Hold plot
#interactive()
# hold plot
#plt.show()

Let me know if you need any more info. If it helps, most of this code was adapted from the Gaussian Heat equation example here: fenics-tutorial/ft04_heat_gaussian.py at master · hplgit/fenics-tutorial · GitHub. I then tried to modify it to work on 3D and for the nonlinear weak form.

2

You’ve overwritten the variable u such that when you compute the Frechet derivative necessary for Newton’s method in solve(F == 0, u) you get a zero Jacobian since F does not contain your new u.

If you delete the line

...
# Time-stepping
u = Function(V)  # Delete this
...

you will at least have a system which can be assembled.

Futhermore, for some hints and tips regarding nonlinear problems, see here.

2 Likes
3

Thank you! Really appreciate your quick response. That makes a lot of sense. I was confused about that line, but kept it in since the heat gaussian example I linked had 2 definitions of u like that (one was a TrialFunction though).

Thanks a ton for the link as well.