8
$\begingroup$

I've read here that Mathematica 10 can Obtain symbolic solutions to delay differential equations.

Would that help in numericaly evaluating Dickman's function $\rho(u)$ ?

It is a delay differential equation defined for $u\in\mathbb R^+$ by $$\begin{cases} \rho(u)=1 & \text{if }0\le u\le1\\ u\cdot\rho'(u)+\rho(u-1)=0& \text{if }u>1 \end{cases}$$

Otherwise said, can Mathematica 10 (or something) help improve/simplify this (working) code?

Dickm$pr=40 (*desired ABSOLUTE (rather than relative) accuracy*);

Dickm$ex=Ceiling[Dickm$pr*Log[3,10]]+1 (*order of expansions; empirical*);

Dickm$ma=Ceiling[Dickm$ex/Log[Dickm$ex]] (*domain limit*);

(*Dickm$dv[n] approximates DickmanRho[t+n-1/2] for -1/2\[LessEqual]t\[LessEqual]1/2*)

Dickm$dv[1]=N[1,Dickm$pr+1];

Dickm$dv[n_]:=Dickm$dv[n]=Module[{u=Dickm$dv[n-1],v},v=Integrate[(u+O[\.b2]^Ceiling[Dickm$ex-n])/(\.b2+n-1/2),\.b2];(Normal[u]/.\.b2->1/2)+(Normal[v]/.\.b2->-1/2)-v];

DickmanRho[y_?NumberQ]:=Which[y<=0,1,y>Dickm$ma,N[1,y Log[10.,y]]-1,True,Normal[Dickm$dv[Ceiling[y]]]/.\.b2->y-Ceiling[y]+1/2];

(* Example use *)
DickmanRho[15]
Improve this question
$\endgroup$
1
  • 3
    $\begingroup$ Ack. I was going to suggest this: "Might want to check this old MathGroup post." Then I realized that was your work. $\endgroup$
    – Daniel Lichtblau
    Commented Sep 11, 2014 at 20:31

1 Answer 1

Reset to default
7
$\begingroup$

Yes!

fun = NDSolve[{u ρ'[u] + ρ[u - 1] == 0, ρ[u /; u < 1] == 1}, ρ, {u, 0, 20}, 
    WorkingPrecision -> 40, PrecisionGoal -> ∞, AccuracyGoal -> 40, 
    Method -> "StiffnessSwitching"][[1, 1, 2]];

fun[15]
(* 7.5899080042980595046528227779709126741*10^-20 *)

DickmanRho[15]
(* 7.5899080042980595047*10^-20 *)

It works not only in V10 (in V8 at least).

Improve this answer
$\endgroup$
3
  • 2
    $\begingroup$ @fgrieu, I tune the options for NDSolve and now it works much better. However, it takes several seconds which is ~100 times slower than your clever method! $\endgroup$
    – ybeltukov
    Commented Sep 12, 2014 at 7:41
  • $\begingroup$ Looks like Method made a difference. $\;$ I wonder if we could make use the new DSolve capabilities. $\endgroup$
    – fgrieu
    Commented Sep 12, 2014 at 8:06
  • 1
    $\begingroup$ @fgrieu I tried DSolve, but it produces very huge formulas so we have to manually simplify and rounding them. But it is what you have already done... $\endgroup$
    – ybeltukov
    Commented Sep 12, 2014 at 8:14

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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