• API
• FAQ
• Tools
• Archive
SHARE
TWEET # Untitled a guest Jul 5th, 2017 65 Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
1. rscalc[rm2m_] := rm2m + 2 m*(1 + Log[(rm2m)/(2 m)])
2. xmax = 700; xmin = -700; delx = 0.03; imax = Floor[(xmax - xmin)/delx];
3. rstar = xmax; m = 1.0; l = 2.0; rcour = 0.99; ntot = 3;
4. If[rstar > 4 m, r = rstar, r = 2 m*E^((rstar/(2 m)) - 1)];
5. rm2m = r - 2 m;
6. Do[rstest = rscalc[rm2m];
7.   drsdr = (rm2m + 2 m)/rm2m;
8.   rtempm2m = rm2m + (rstar - rstest)/drsdr;
9.   rm2m = rtempm2m;, {i, 7}];
10. mesh = NDSolve[{rm2mc'[rs] == rm2mc[rs]/(rm2mc[rs] + 2 m),
11.    rm2mc[rstar] == rm2m}, rm2mc, {rs, xmax, xmin}, AccuracyGoal -> 25,
12.    PrecisionGoal -> 25, WorkingPrecision -> 50]
13. Plot[Evaluate[rm2mc[lex] /. mesh], {lex, -20, 20}, PlotRange -> All]
14. xtest = -250;
15. b = Evaluate[rm2mc[xtest] /. mesh]
16. xchck = rscalc[b];
17. err = xchck - xtest
18. V[rm2m_] := (rm2m/(rm2m +
19.        2 m))*((l (l + 1))/(rm2m + 2 m)^2 + (2 m/(rm2m + 2 m)^3));
20. Plot[V[nax - 2 m], {nax, 0, 10}]
21. Vsx[x_] := V[Evaluate[rm2mc[x] /. mesh]];
22. Plot[Vsx[x], {x, -20, 20}]
23. XTab = Table[xmin + i*delx, {i, imax + 1}];
RAW Paste Data
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy.
Top