I am trying to implement an algorithm described for example in this paper, Section III, or referenced here, from where numerical parameters have been taken. I added the snippet I use to compute the matrix from Eq. 3.1 in the first reference at the end, if it could be of any benefit. The need arises to compute very small numbers. Concretely, set
Km = 0.025;
Dm = 0.063;
am = 4.2;
kb = 8.617*^-5;
V[delta_] := Dm*(1 - Exp[-delta])^2;
numbers such as
Exp[-V[-40]]
must be calculated. I read these answers How to check for underflow, and find a constant to correct it? and here New General::munfl error and loss of precision but I am not getting anywhere, probably because I do not have yet any understanding of the underlying theory.
I tried
Activate@SetPrecision[Inactivate[Exp[-V[-40]]], 15]
as well as
Exp[SetPrecision[-V[-40], 15]]
or
Quiet[Exp[-V[-40]], General::unfl]
but I still get underflows errors. The question whether these figures are to be computed at all is also interesting, but I am trying to replicate the results of the paper.
How should these overflows be handled?
As a layman, I would not in principle mind something simple as Julia does
julia> dm = 0.063
0.063
julia> morse(x)=dm*(1-exp(-x))^2
morse (generic function with 1 method)
julia> morse(-40)
3.4905921021679106e33
julia> exp(-morse(-40))
0.0
although I understand this could be insufficient for more complex applications. I assume underflow messages a warning about the achieved precision, and automatic truncation to 0 might be undesirable, but I do not know.
** Matrix to be diagonalised**
lmin = -40;
lmax = 5000;
ndiv = 20;
divis = (lmax - lmin)/ndiv;
W[delta1_, delta2_] := (Km/(2*am^2))*(delta1 - delta2)^2;
mater[tempe_] :=
Table[(divis)*
Exp[ -(1/(kb*tempe))*(0.5*V[lmin + divis*i] +
0.5*V[lmin + divis*j] +
W[lmin + divis*i, lmin + divis*j])] , {i, 0, ndiv}, {j, 0,
ndiv}];
Quiet[Exp[-V[-40]], General::munfl]
. OrOff[General::munfl]
. $\endgroup$ – Michael E2 Apr 12 '20 at 12:48Exp[-V[-40]]
by itself givesGeneral::munfl
for me, notGeneral::unfl
. I'm using V12.1.0, in case it's a version issue. -- BTW, theSetPrecision
does not work because it does not inspect the code forV[]
, which contains a machine-precision number. $\endgroup$ – Michael E2 Apr 12 '20 at 13:39