I have a complicated function called f2[re]
f2[re_] :=
1/(2 cep^(3/2) (cep - 2 cp) cp^(3/2) re)
E^(-cp inf^2 - cep inf re - ce re^2 - (cep re \[Alpha]2)/cp -
2 (inf + re) \[Alpha]2)
norm^2 \[Pi] (2 cep^(3/2) Sqrt[cp] E^((cep re \[Alpha]2)/cp) -
4 Sqrt[cep] cp^(3/2) E^((cep re \[Alpha]2)/cp) -
2 cep^(3/2) Sqrt[cp] E^(
re (2 cep inf + 4 \[Alpha]2 + (cep \[Alpha]2)/cp)) +
4 Sqrt[cep] cp^(3/2) E^(
re (2 cep inf + 4 \[Alpha]2 + (cep \[Alpha]2)/cp)) +
4 cp^(3/2) E^(
2 (inf + re) \[Alpha]2 + (2 \[Alpha]2^2)/cep + (
cep (cp (inf + re)^2 + 2 re \[Alpha]2))/(2 cp)) Sqrt[
2 \[Pi]] \[Alpha]2 Erf[(cep (inf - re) + 2 \[Alpha]2)/(
Sqrt[2] Sqrt[cep])] +
cep^(3/2) E^(
cp inf^2 + cep inf re + 2 inf \[Alpha]2 +
4 re \[Alpha]2 + ((cep^2 re^2)/4 + \[Alpha]2^2)/cp)
Sqrt[\[Pi]] (cep re - 2 (cp re + \[Alpha]2)) Erf[(
2 cp inf - cep re + 2 \[Alpha]2)/(2 Sqrt[cp])] +
cep^(5/2) E^(
cp inf^2 + cep inf re + (cep^2 re^2)/(4 cp) +
2 inf \[Alpha]2 + (2 cep re \[Alpha]2)/cp + \[Alpha]2^2/cp)
Sqrt[\[Pi]]
re Erf[(2 cp inf + cep re + 2 \[Alpha]2)/(2 Sqrt[cp])] -
2 cep^(3/2) cp E^(
cp inf^2 + cep inf re + (cep^2 re^2)/(4 cp) +
2 inf \[Alpha]2 + (2 cep re \[Alpha]2)/cp + \[Alpha]2^2/cp)
Sqrt[\[Pi]]
re Erf[(2 cp inf + cep re + 2 \[Alpha]2)/(2 Sqrt[cp])] +
2 cep^(3/2) E^(
cp inf^2 + cep inf re + (cep^2 re^2)/(4 cp) +
2 inf \[Alpha]2 + (2 cep re \[Alpha]2)/cp + \[Alpha]2^2/cp)
Sqrt[\[Pi]] \[Alpha]2 Erf[(2 cp inf + cep re + 2 \[Alpha]2)/(
2 Sqrt[cp])] -
cep^(5/2) E^(
cp inf^2 + cep inf re + (cep^2 re^2)/(4 cp) +
2 inf \[Alpha]2 + (2 cep re \[Alpha]2)/cp + \[Alpha]2^2/cp)
Sqrt[\[Pi]]
re Erf[(cep re - 2 cp re + 2 \[Alpha]2)/(2 Sqrt[cp])] +
2 cep^(3/2) cp E^(
cp inf^2 + cep inf re + (cep^2 re^2)/(4 cp) +
2 inf \[Alpha]2 + (2 cep re \[Alpha]2)/cp + \[Alpha]2^2/cp)
Sqrt[\[Pi]]
re Erf[(cep re - 2 cp re + 2 \[Alpha]2)/(2 Sqrt[cp])] -
2 cep^(3/2) E^(
cp inf^2 + cep inf re + (cep^2 re^2)/(4 cp) +
2 inf \[Alpha]2 + (2 cep re \[Alpha]2)/cp + \[Alpha]2^2/cp)
Sqrt[\[Pi]] \[Alpha]2 Erf[(cep re - 2 cp re + 2 \[Alpha]2)/(
2 Sqrt[cp])] -
4 cp^(3/2) E^(
2 (inf + re) \[Alpha]2 + (2 \[Alpha]2^2)/cep + (
cep (cp (inf + re)^2 + 2 re \[Alpha]2))/(2 cp)) Sqrt[
2 \[Pi]] \[Alpha]2 Erf[(cep (inf + re) + 2 \[Alpha]2)/(
Sqrt[2] Sqrt[cep])] +
cep^(5/2) E^(
cp inf^2 + cep inf re + 2 inf \[Alpha]2 +
4 re \[Alpha]2 + ((cep^2 re^2)/4 + \[Alpha]2^2)/cp) Sqrt[\[Pi]]
re Erf[(cep re - 2 (cp re + \[Alpha]2))/(2 Sqrt[cp])] -
2 cep^(3/2) cp E^(
cp inf^2 + cep inf re + 2 inf \[Alpha]2 +
4 re \[Alpha]2 + ((cep^2 re^2)/4 + \[Alpha]2^2)/cp) Sqrt[\[Pi]]
re Erf[(cep re - 2 (cp re + \[Alpha]2))/(2 Sqrt[cp])] -
2 cep^(3/2) E^(
cp inf^2 + cep inf re + 2 inf \[Alpha]2 +
4 re \[Alpha]2 + ((cep^2 re^2)/4 + \[Alpha]2^2)/cp)
Sqrt[\[Pi]] \[Alpha]2 Erf[(cep re - 2 (cp re + \[Alpha]2))/(
2 Sqrt[cp])]);
The values of function f2
must be equal to function f1
in different re
.
function f1
includes a numerical integration as follows
f1[re_?NumericQ] := (2*\[Pi]*norm^2)/
re*(Exp[-ce re^2]) (1/cep^(3/2)) NIntegrate[
rp *Exp[-cp rp^2]*(E^(-cep re rp) (Sqrt[
cep] (-E^(-2 (re + rp) \[Alpha]2) + E^(
2 cep re rp - 2 \[Alpha]2 Abs[re - rp])) -
E^(1/2 cep (re + rp)^2 + (2 \[Alpha]2^2)/cep) Sqrt[
2 \[Pi]] \[Alpha]2 Erf[(cep (re + rp) + 2 \[Alpha]2)/(
Sqrt[2] Sqrt[cep])] +
E^(1/2 cep (re + rp)^2 + (2 \[Alpha]2^2)/cep) Sqrt[
2 \[Pi]] \[Alpha]2 Erf[(2 \[Alpha]2 + cep Abs[re - rp])/(
Sqrt[2] Sqrt[cep])])), {rp, 0, inf}];
So I specify the needed constants as follows
ce = 0.0006533997022026266`;
cp = 2.8846533997022026`;
cep = -0.026693200595594747`;
inf = 15;
norm = 0.5293704593281444`;
{\[Alpha]2, \[Beta]2} = {0.9948987042561277`, 0.0002930460049474671`};
when I try
f1[0.1]
f2[0.1]
I get
0.0968872575537 + 0. I
0.0968872575541 + 0. I
which is OK as I expect. But as soon as I set inf=16
it returns
General::munfl: Exp[-770.463385382] is too small to represent as a normalized machine number; precision may be lost.
0. + 0. I
I googled this error and it said it relates to MachinePrecision
, It's obvious if inf
>15 some parts of f2[re]
become so small so that Mathematica makes them zero, so the whole function becomes zero. However I don't any idea to address this problem and increase the precision of calculations. Any idea?
Note: inf
is finite number (its maximum is 20) which I used as upper limit of integration to obtain f2[re]
NIntegrate
to useWorkingPrecision->25
(or anything modestly lower than the precision of the constants), and using a bignum instead of .1 in the inputs. $\endgroup$