Mathematica很强大

shunfly

装了6.0，弄了个Numerical Laplace Inversion的CODE
Mathematica太强大，基于Mathematica的很强大的高精度计算贼爽，可惜的是，这种高精度计算无法在实际的
数值计算中得以很好应用。

shunfly

1. Unprotect[GWR];
   
2. 
   
3. Begin["User`NumericalLaplaceTransformInversion`Private`"]
   
4. GWR[F_, t_, M_:12, precin_:0] := Module[
   
5.     {M1, G0, Gm, Gp, best, expr, \[Tau] = Log[2]/t, Fi, broken, \
   
6. prec},
   
7.     If[precin <= 0, prec = 21 M/10, prec = precin];
   
8.     If[prec <= $MachinePrecision, prec = $MachinePrecision];
   
9.     broken = False;
   
10.     If[Precision[\[Tau]]<prec,\[Tau]=SetPrecision[\[Tau],prec]];
    
11.     Do[Fi[i] = N[F[i \[Tau]], prec], {i, 1, 2 M}];
    
12.     (* Fi[I]为F[Ln2/T*N] *)
    
13.     M1 = M;
    
14.     Do[
    
15.       G0[n - 1] = \[Tau] (2n)!/(n!(n - 1)!)Sum[
    
16.             Binomial[n, i](-1)^i  Fi[n + i], {i, 0, n}];
    
17.       If[Not[NumberQ[G0[n - 1]]], M1 = n - 1; G0[n - 1] =.; Break[]];
    
18.       , {n, 1, M}];
    
19.   Print["参数 ",G0[0],"参数 ",G0[M-1]];
    
20.     Do[Gm[n] = 0, {n, 0, M1}];
    
21.     best = G0[M1 - 1];
    
22.    
    
23.     Do[
    
24.       Do[
    
25.         expr = G0[n + 1] - G0[n];
    
26.         If[Or[Not[NumberQ[expr]], expr == 0], broken = True; \
    
27. Break[]];
    
28.         expr = Gm[n + 1] + (k + 1)/expr;
    
29.         Gp[n] = expr;
    
30.         If[OddQ[k],
    
31.           If[n == M1 - 2 - k, best = expr]
    
32.           ];
    
33.         , {n, M1 - 2 - k, 0, -1}];
    
34.       If[broken, Break[]];
    
35.       Do[Gm[n] = G0[n]; G0[n] = Gp[n], {n, 0, M1 - k}];
    
36.       , {k, 0, M1 - 2}];
    
37.     best
    
38.     ]
    
39. 
    
40. End[]  (* User`NumericalLaplaceTransformInversion` *)

复制代码

shunfly

http://library.wolfram.com/infocenter/MathSource/4738/
拉普拉斯逆变换数值计算程序

FreddyMusic

library.wolfram.com is some old archive information.

I am looking at demo project / blog, that's more creative.

http://demonstrations.wolfram.com/

shunfly

感谢FreddyMusic兄