4
$\begingroup$

Bug persisting through 11.2


I'am facing a trouble with HypergeometricPFQRegularized. Mathematica hangs if imaginary part of one of the arguments has zero imaginary part with decimal point as the following example shows:

\[Nu]j = 2 + I*0
HypergeometricPFQRegularized[{1/2, 1}, {1 - \[Nu]j, 
  1 + \[Nu]j}, -14.28]
\[Nu]j = 2 + I*0.
HypergeometricPFQRegularized[{1/2, 1}, {1 - \[Nu]j, 
  1 + \[Nu]j}, -14.28]

(*
Out[404]= 2

Out[405]= 0.112647

Out[406]= 2. + 0. I

During evaluation of In[404]:= Divide::infy: Infinite expression 3/0 encountered.

During evaluation of In[404]:= Divide::infy: Infinite expression 3/0 encountered.

During evaluation of In[404]:= Divide::infy: Infinite expression 3/0 encountered.

During evaluation of In[404]:= General::stop: Further output of Divide::infy will be suppressed during this calculation.
*)

This code is used in heavy serial calculations which hangs sometimes and my investigation indicated on the above trouble with HypergeometricPFQRegularized. Could someone propose simple workaround? I need effective code that would round [Nu]j when it is close to an integer.

$\endgroup$
  • $\begingroup$ Perhaps Rationalize or Chop will solve your problem. $\endgroup$ – Bob Hanlon Feb 8 '18 at 13:30
8
$\begingroup$

Diagnosis

Upon inspecting the source code for HypergeometricPFQRegularized with

<< GeneralUtilities`
PrintDefinitions[HypergeometricPFQRegularized]

I have tracked that error in the OP's example down to a poorly constructed pattern for the 10th DownValue:

HypergeometricPFQRegularized[a_List,{b1___,Except[_Complex, c_],b2___},z_] := ...

Clearly the programmer intended for the rule to be applied in case c was real, but it fails to match if c is an approximate real number with head Complex.

Furthermore, the RHS of the rule, additionally imposes a Condition:

... := Condition[(*body*), Internal`RealValuedNumericQ[c] && 
                           Block[{rc = Round[c]}, rc <= 0 && c == rc]]

which among other things again checks for the real-valuedness of c, but failed to pass if c has head Complex.

Resolution

To fix this code, I have taken the liberty of modifying that rule (along with parts of the body) so that it functions as it was (likely) intended. From a fresh kernel, run the following:

HypergeometricPFQRegularized[{}, {}, 0]; (*loads definitions*)
Unprotect[HypergeometricPFQRegularized];

(*Modifies the 10th rule*)
DownValues[HypergeometricPFQRegularized] = 
  ReplacePart[DownValues[HypergeometricPFQRegularized], 10 ->
    HoldPattern[HypergeometricPFQRegularized[a_List, {b1___, c_, b2___}, z_] /; PossibleZeroQ[Im[c]]] :> 
      Module[{temp, v},
         temp = Select[{b1, c, b2}, (Last[Internal`TestIntegerQ[#]] && # <= 0) &];
         v = 1 - Round[Min[temp]]; 
         z^v (Times @@ Pochhammer[a, v] HypergeometricPFQRegularized[a + v,
           Append[System`HypergeometricPFQDump`reduce[{b1, c, b2}, {1 - v}] + v, 1 + v], z])
      ] /; Block[{rc = Round[c]}, rc <= 0 && c == rc]];

Protect[HypergeometricPFQRegularized];

And now, the numerical evaluation proceeds without problems:

νj = 2 + I*0.;
HypergeometricPFQRegularized[{1/2, 1}, {1 - νj, 1 + νj}, -14.28]

(*  0.112647  *)
|improve this answer|||||
$\endgroup$
  • $\begingroup$ PrintDefinitions seems like a really useful function! $\endgroup$ – chris Feb 9 '18 at 7:48

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.