← Back to team overview

dolfin team mailing list archive

[Question #149194]: convergence problem with coupled poisson equations

 

New question #149194 on DOLFIN:
https://answers.launchpad.net/dolfin/+question/149194

Dear all,

I need to solve the drift-diffusion equations for a degenerate semiconductor. The following code works well when the doping is not too high +/- 1000. However at doping of +/- 10000 the problem does not converge and I get the error below. I think this has something to do with the numbers becoming too big but I'd like to know exactly where the error happens. 

Also any ideas on how to scale the problem for this not to happen? I am using dolfin 0.9.9. 

Newton iteration 1000: r (abs) = -nan (tol = 1.000e-10) r (rel) = -nan (tol = 1.000e-09)
---------------------------------------------------------------------------
RuntimeError                              Traceback (most recent call last)
/home/chaffra/Projects/Fenics/Semifem/test/ddtest.py in <module>()
     59 newton_solver = problem.newton_solver()
     60 newton_solver.parameters['maximum_iterations'] = 1000
---> 61 problem.solve(u_func)
/usr/local/lib/python2.6/dist-packages/dolfin/fem/variationalproblem.py in solve(self, u)
     64 
     65         # Solve
---> 66         cpp.VariationalProblem.solve(self, u)
     67 
     68         # Return solution
RuntimeError: *** Error: Newton solver did not converge.


------------------------My example
from dolfin import *
#from semifem.utils.fermi_dirac import fermi_dirac_integral, form_compiler_parameters


m = UnitInterval(100)
Q = FunctionSpace(m,'CG',1)


W = MixedFunctionSpace([Q,Q,Q])
u_func = Function(W)
u, un, up = split(u_func)
v, vn, vp = TestFunctions(W)

du_func = TrialFunction(W)
du, dun, dup = split(du_func)

class Doping(Expression):
    
    def eval(self, v, x):
        if x[0] < 0.5:
            v[0] = 10000
        else:
            v[0] = -10000

n = exp(u-un)
p = exp(-(u-up))
d = Doping()
rho = p+d-n
mun = Constant(10000)
mup = Constant(1000)

L1 = inner(grad(u),grad(v))*dx -rho*v*dx
a1 = derivative(L1,u_func,du_func)

L2 = inner(mun*n*grad(un), grad(vn))*dx
a2 = derivative(L2, u_func, du_func)

L3 = inner(mup*p*grad(up), grad(vp))*dx
a3 = derivative(L3,u_func,du_func)

bc_left0 = DirichletBC(W.sub(0),0.0,"x[0]<DOLFIN_EPS")
bc_left1 = DirichletBC(W.sub(1),1e-1,"x[0]<DOLFIN_EPS")
bc_left2 = DirichletBC(W.sub(2),1e-1, "x[0]<DOLFIN_EPS")

bc_right0 = DirichletBC(W.sub(0),0.0,"x[0] == 1")
bc_right1 = DirichletBC(W.sub(1),0.0,"x[0] == 1")
bc_right2 = DirichletBC(W.sub(2),0.0,"x[0] == 1")

problem = VariationalProblem(a1+a2+a3,L1+L2+L3,bcs=[bc_left0, bc_left1, bc_left2, bc_right0, bc_right1, bc_right2],nonlinear=True, 
                             form_compiler_parameters={})

newton_solver = problem.newton_solver()
newton_solver.parameters['maximum_iterations'] = 1000
problem.solve(u_func)



-- 
You received this question notification because you are a member of
DOLFIN Team, which is an answer contact for DOLFIN.