Poisson+internal Neumann: Why I am getting NaN and Inf over all cells ?


109
views
0
12 weeks ago by
rooz  
Hi there

I am trying to solve a nonlinear Poisson equation by Picard iteration,  However I am getting no errors but nonsense values from the first round, I cant see whether I have made a mistake or should I play with the solver options,

I have tried to summarize the problem into a short format, I hope it is clear enough

thank you for your help in advance

Cheers


 
from dolfin import * 
import numpy as np

#-----------------------------------------------
eps0 = 8.85E-12  
q = 1.60217657E-19
KB = 1.3806488E-23
T =300
REV = 0  
Q_interface =  2.12E-2 

#-----------------------------------------------
width_A = 2E-9 
width_B = 2E-9
x_A   = 0
x_B = width_A +x_A
x_end   = width_B + x_B 

mesh = IntervalMesh( int(abs(x_end-x_A)/.5E-10) , x_A , x_end)
V = FunctionSpace( mesh, "DG" , 0)  
u , v = TrialFunction(V) , TestFunction(V) 
#-----------------------------------------------
tol = 1E-14  ; 
class Subdomain_1(SubDomain):
     def inside(self , x , on_boundary):
         return    (x[0] >= x_A - tol )   and  (x[0] <= x_B + tol)  
 
class Subdomain_2(SubDomain):
     def inside(self , x , on_boundary):
         return     (x[0] >= x_B - tol )   and   (x[0] <= x_end + tol)  
         
subdomains = CellFunction("size_t" , mesh)  
subdomains.set_all(0)  
Subdomain_1().mark(subdomains,1) 
Subdomain_2().mark(subdomains,2) 

boundary_markers = FacetFunction("size_t" , mesh)
boundary_markers.set_all(0)  
Neuman_boundary= CompiledSubDomain(  "fabs(x[0]-"+ str(x_B) +")<DOLFIN_EPS" ) 
Neuman_boundary.mark(boundary_markers,1)
Dirich_boundary= CompiledSubDomain(  "fabs(x[0]-"+ str(x_end) +")<DOLFIN_EPS" ) 
Dirich_boundary.mark(boundary_markers,2)
bc = DirichletBC(V , Constant(REV) , boundary_markers , 2 ) 

dx =Measure('dx',domain=mesh,subdomain_data=subdomains)
ds =Measure('ds',domain=mesh,subdomain_data=boundary_markers) 
dS =Measure('dS',domain=mesh,subdomain_data=boundary_markers) 

#-----------------------------------------------
dielectric_v = eps0* np.array([10.4,10.34 ] )
eps  = Function(V)
for cell_no in range (len(subdomains.array())):
    subdomain_no = subdomains.array()[cell_no]
    eps.vector()[cell_no] =dielectric_v[subdomain_no-1]

nominal_doping_profile = np.array([1E14 ,2E14 ] )
Nd  = Function(V)
for cell_no in range (len(subdomains.array())):
    subdomain_no = subdomains.array()[cell_no]
    Nd.vector()[cell_no] =nominal_doping_profile[subdomain_no-1]
#-----------------------------------------------    
a =   inner(eps*grad(u) , grad(v))*dx(1) \
    + inner(eps*grad(u) , grad(v))*dx(2) 
         
L =  q*(Nd )*v * dx(1) \
    +q*(Nd )*v * dx(2) \
    -avg(Q_interface)*v('+')*dS(1)   

    
u= Function(V)
solve(a==L,u,bc)

 


​
Community: FEniCS Project

2 Answers


1
12 weeks ago by
You define DG0 function space, and then your energy involves a gradient. Did you try with CG1 for instance?
0
12 weeks ago by
rooz  
Hi Martin,
As you already mentioned,I had made a mistake using DG0 and using CG1 resolved that specific problem.
Thank you for your help,
Cheers  
Please login to add an answer/comment or follow this question.

Similar posts:
Search »