Re: ValueError: math domain error

[Daniel Wheeler]

On Nov 15, 2006, at 4:42 PM, Omer Yassur wrote:
Hi Daniel,

I tried to work with the debugger in order to trace the error,  
unfortunately with no success.
You can find below, the original script and the Traceback.
Short glance at the division terms show that order of math operators
is rational and seems o.k.

Can you please advice me on how to debug this problem?

File "C:\Python23\Lib\site-packages\fipy\tools\numerix.py", line  
767, in exp
    return NUMERIC.exp(arr)


It is clear from the traceback below that the Numeric.exp function is  
exploding.

File "C:\Python23\Lib\site-packages\fipy\tools\numerix.py", line  
767, in exp
    return NUMERIC.exp(arr)

You can check this by going to line 767 of numerix.py and sticking in  
a print statement
before the return. One of the values in the array is becoming large  
and positive.
Only print the maximum value. You'll get loads of values, but the  
last one will be the
one that fails. This will be coming from your coefficients in the  
places where you use the exp function.


For example you have exp(-B*psy) in one of your coeffiecents. B or  
psy is probably becoming negative
you can print the maximum of these arrays before you go into sweep()  
at each iteration. Should be clear if
you do that.

Also you can plot during the sweep to get a feel for how the solution  
is behaving during the sweep.
psy is probably becoming negative.

Many thanks,
Omer



#!/usr/bin/env python


# define a 1D mesh
# total length 3um
dx = 5e-8
nx = 6000

# define mesh
from fipy.meshes.grid1D import Grid1D
mesh = Grid1D(dx=dx, nx=nx)

# declare the problem variable (potential across the device)
from fipy.variables.cellVariable import CellVariable
psy = CellVariable(name="psy", mesh=mesh, hasOld=1)

from fipy.tools.numerix import exp
Kb = 1.38e-23
T = 300
q = 1.6e-19
B = q/(Kb*T)
ni = 1.5e10
Na = 5e15
np0 = (ni**2)/Na
eps_0 = 8.854e-14
eps_si = 11.8*eps_0
Vref=0.2
Ho=6.02e17
eps_aqa=81*eps_0
dbulk=1e-4
dh=1.001e-4

# calculate the charge density in the semiconductor
x = mesh.getCellCenters()[...,0]
s1_semi =( -((q*B)/(eps_si))*(Na*exp(-B*psy)+np0*exp(B*psy))) * (x  
<= dbulk)
s0_semi= ((q/eps_si)*(Na*(exp(-B*psy)-1)-np0*(exp(B*psy)-1))- 
s1_semi*psy) * (x <= dbulk)
from fipy.terms.implicitSourceTerm import ImplicitSourceTerm
implicitSource_semi = ImplicitSourceTerm(coeff = s1_semi)


# calculate the charge density in the diffuse layer
psyD=psy-Vref
s1_electrolyte=-(((q*Ho*B)/(eps_aqa))*(exp(-B*psyD)+exp(B*psyD))) *  
(x >= dh)
s0_electrolyte=(((q*Ho)/(eps_aqa))*(exp(-B*psyD)-exp(B*psyD))- 
s1_electrolyte*psyD) * (x >= dh)
from fipy.terms.implicitSourceTerm import ImplicitSourceTerm
implicitSource_electrolyte = ImplicitSourceTerm(coeff =  
s1_electrolyte)


# boundary condition
from fipy.boundaryConditions.fixedValue import FixedValue
from fipy.boundaryConditions.fixedFlux import FixedFlux
bcs = (FixedValue(faces=mesh.getFacesRight(), value=Vref),
       FixedFlux(faces=mesh.getFacesLeft(), value=0))

from fipy.terms.implicitDiffusionTerm import ImplicitDiffusionTerm
diffCoeff = 1.
diffTerm = ImplicitDiffusionTerm(coeff=diffCoeff)

# poisson equation
samfet_eq = diffTerm + s0_semi + implicitSource_semi +  
s0_electrolyte + implicitSource_electrolyte == 0

# solve the equation
res = 1e+10
while res > 1e-6:
        res = samfet_eq.sweep(var=psy,boundaryConditions=bcs)

if __name__ == '__main__':
   from fipy.viewers import make
   viewer = make(vars=(psy))
   viewer.plot()

from fipy.viewers.tsvViewer import TSVViewer
TSVViewer(vars=(psy)).plot(filename=".\out\samfet.tsv")



cg not converged
Traceback (most recent call last):
  File "C:\Documents and Settings\oo\My Documents\omer\Thesis\Fipy 
\work\samfet.py", line 65, in ?
    res = samfet_eq.sweep(var=psy,boundaryConditions=bcs)
  File "C:\Python23\Lib\site-packages\fipy\terms\term.py", line  
149, in sweep
    matrix, RHSvector = self.__buildMatrix(var, boundaryConditions,  
dt)
  File "C:\Python23\Lib\site-packages\fipy\terms\term.py", line 95,  
in __buildMatrix
    return self._buildMatrix(var, boundaryConditions, dt)
  File "C:\Python23\Lib\site-packages\fipy\terms\binaryTerm.py",  
line 64, in _buildMatrix
    matrix, RHSvector = self.term1._buildMatrix(var,  
boundaryConditions, dt = dt)
  File "C:\Python23\Lib\site-packages\fipy\terms\binaryTerm.py",  
line 64, in _buildMatrix
    matrix, RHSvector = self.term1._buildMatrix(var,  
boundaryConditions, dt = dt)
  File "C:\Python23\Lib\site-packages\fipy\terms\binaryTerm.py",  
line 64, in _buildMatrix
    matrix, RHSvector = self.term1._buildMatrix(var,  
boundaryConditions, dt = dt)
  File "C:\Python23\Lib\site-packages\fipy\terms\binaryTerm.py",  
line 66, in _buildMatrix
    termMatrix, termRHSvector = self.term2._buildMatrix(var,  
boundaryConditions, dt = dt)
  File "C:\Python23\Lib\site-packages\fipy\terms\cellTerm.py", line  
121, in _buildMatrix
    inline._optionalInline(self._buildMatrixIn,  
self._buildMatrixPy, L, var.getOld(), b, dt, coeffVectors)
  File "C:\Python23\Lib\site-packages\fipy\tools\inline\inline.py",  
line 10, in _optionalInline
    return pythonFn(*args)
  File "C:\Python23\Lib\site-packages\fipy\terms\cellTerm.py", line  
89, in _buildMatrixPy
    b += Numeric.ones([N]) * Numeric.array(coeffVectors['b vector'])
  File "C:\Python23\Lib\site-packages\fipy\variables\variable.py",  
line 162, in __array__
    return Numeric.array(self.getValue(), t)
  File "C:\Python23\Lib\site-packages\fipy\variables 
\cellVariable.py", line 108, in getValue
    return Variable.getValue(self)
  File "C:\Python23\Lib\site-packages\fipy\variables\variable.py",  
line 342, in getValue
    value = self._calcValue()
  File "C:\Python23\Lib\site-packages\fipy\variables\variable.py",  
line 1708, in _calcValue
    return self.op(self.var[0].getValue(), val1)
  File "C:\Python23\Lib\site-packages\fipy\variables 
\cellVariable.py", line 108, in getValue
    return Variable.getValue(self)
  File "C:\Python23\Lib\site-packages\fipy\variables\variable.py",  
line 342, in getValue
    value = self._calcValue()
  File "C:\Python23\Lib\site-packages\fipy\variables\variable.py",  
line 1708, in _calcValue
    return self.op(self.var[0].getValue(), val1)
  File "C:\Python23\Lib\site-packages\fipy\variables 
\cellVariable.py", line 108, in getValue
    return Variable.getValue(self)
  File "C:\Python23\Lib\site-packages\fipy\variables\variable.py",  
line 342, in getValue
    value = self._calcValue()
  File "C:\Python23\Lib\site-packages\fipy\variables\variable.py",  
line 1708, in _calcValue
    return self.op(self.var[0].getValue(), val1)
  File "C:\Python23\Lib\site-packages\fipy\variables 
\cellVariable.py", line 108, in getValue
    return Variable.getValue(self)
  File "C:\Python23\Lib\site-packages\fipy\variables\variable.py",  
line 342, in getValue
    value = self._calcValue()
  File "C:\Python23\Lib\site-packages\fipy\variables\variable.py",  
line 1702, in _calcValue
    val1 = self.var[1].getValue()
  File "C:\Python23\Lib\site-packages\fipy\variables 
\cellVariable.py", line 108, in getValue
    return Variable.getValue(self)
  File "C:\Python23\Lib\site-packages\fipy\variables\variable.py",  
line 342, in getValue
    value = self._calcValue()
  File "C:\Python23\Lib\site-packages\fipy\variables\variable.py",  
line 1708, in _calcValue
    return self.op(self.var[0].getValue(), val1)
  File "C:\Python23\Lib\site-packages\fipy\variables 
\cellVariable.py", line 108, in getValue
    return Variable.getValue(self)
  File "C:\Python23\Lib\site-packages\fipy\variables\variable.py",  
line 342, in getValue
    value = self._calcValue()
  File "C:\Python23\Lib\site-packages\fipy\variables\variable.py",  
line 1708, in _calcValue
    return self.op(self.var[0].getValue(), val1)
  File "C:\Python23\Lib\site-packages\fipy\variables 
\cellVariable.py", line 108, in getValue
    return Variable.getValue(self)
  File "C:\Python23\Lib\site-packages\fipy\variables\variable.py",  
line 342, in getValue
    value = self._calcValue()
  File "C:\Python23\Lib\site-packages\fipy\variables\variable.py",  
line 1708, in _calcValue
    return self.op(self.var[0].getValue(), val1)
  File "C:\Python23\Lib\site-packages\fipy\variables 
\cellVariable.py", line 108, in getValue
    return Variable.getValue(self)
  File "C:\Python23\Lib\site-packages\fipy\variables\variable.py",  
line 342, in getValue
    value = self._calcValue()
  File "C:\Python23\Lib\site-packages\fipy\variables\variable.py",  
line 1702, in _calcValue
    val1 = self.var[1].getValue()
  File "C:\Python23\Lib\site-packages\fipy\variables 
\cellVariable.py", line 108, in getValue
    return Variable.getValue(self)
  File "C:\Python23\Lib\site-packages\fipy\variables\variable.py",  
line 342, in getValue
    value = self._calcValue()
  File "C:\Python23\Lib\site-packages\fipy\variables\variable.py",  
line 1708, in _calcValue
    return self.op(self.var[0].getValue(), val1)
  File "C:\Python23\Lib\site-packages\fipy\variables 
\cellVariable.py", line 108, in getValue
    return Variable.getValue(self)
  File "C:\Python23\Lib\site-packages\fipy\variables\variable.py",  
line 342, in getValue
    value = self._calcValue()
  File "C:\Python23\Lib\site-packages\fipy\variables\variable.py",  
line 745, in _calcValue
    return self.op(self.var[0].getValue())
  File "C:\Python23\Lib\site-packages\fipy\variables\variable.py",  
line 1951, in <lambda>
    return self._getUnaryOperatorVariable(lambda a: numerix.exp(a))
  File "C:\Python23\Lib\site-packages\fipy\tools\numerix.py", line  
767, in exp
    return NUMERIC.exp(arr)
ValueError: math domain error






--
Daniel Wheeler