32

Okay, this is a bit of an embarassment. Here is a very small modification of the original code. I simply made explicit option settings, made a denominator to Sin explicitly real, that kind of thing. My tests show the same timing as the original, give or take an iota. ie = 200; ez = ConstantArray[0., {ie + 1}]; hy = ConstantArray[0., {ie}]; fdtd1d = ...


28

Perhaps you are looking to build a bifurcation diagram. There are a few approaches in Mathematica mentioned in Documentation, which I give below. Also please take a look at apps of similar nature at the Wolfram Demonstration Project. I do not have time to dive into your specific problem, and give classic examples of logistic map which also a quadratic ...


20

This is a "visual proof" of the Archimedean limiting regular polygons. You could implement the recursion and it would progressively approach $\pi$. The proof lies in the "squeezing" argument. $\pi$ is transcendental (no solution to this recurrence in a closed algebraic expression). Whatever implementation of recursion with approach $\pi$ that is reassuring ...


17

Since a fully NDSolve-based solution is acceptable for you, let me give you one. You simply need the magic of "Pseudospectral" or a dense enough 4th order spatial discretization: mol[n_, o_:"Pseudospectral"] := {"MethodOfLines", "SpatialDiscretization" -> {"TensorProductGrid", "MaxPoints" -> n, "MinPoints" -> n, "DifferenceOrder" -> o}} ...


15

Using Annuity: pmt /.Solve[TimeValue[Annuity[pmt, 52, 1], .02, 0] == 5000, pmt] 155.545 Exercise left for the reader...


13

As Bob Hanlon's answer points out, RecurrenceTable does not hold its arguments, but most especially, it does not hold its iterator arguments. This must surely be viewed as a bug as it departs completely from the iterator idiom established by Table, etc. On the other hand, there is precedent for evaluating the iterator specifications. For example, ...


13

You were correct. NestList is exactly the function you want to use. NestList[Dot[A, #]&, x0, 5] (* {{2, 0}, {1., 1.5}, {-0.4, 2.4}, {-1.64, 2.34}, {-2.224, 1.344}, {-1.9184, -0.1896}} *) Note that the first argument of NestList must be a function.


12

When a difference equation does not have a closed form solution, you cannot realistically expect Mathematica to find one. The exercises on the page you link to ask for the steady state, not for a closed form solution! We can find the steady state solution by letting r=r[n] and f=f[n] and then rewriting the relationships that hold at the fixed point. ss = {...


12

As mentioned in the comment above, handling irregular region with FDM is cumbersome and frustrating in my view, actually that's where I stopped my self-learning of FDM and turned to finite element method (FEM), which is more suitable for this task, and finally write this. (Have a look at this post for more information. ) Nevertheless, there seems to be no ...


12

Introduction A lack of time to write up an answer ironically provided time to reflect on the problem, and some nagging uncertainties about some issues contributed to the delay. The slowness of the OP's code can be seen by a simple analysis, which reveals that some expensive calculations are repeated multiple times for each step in the time integration. A ...


11

First, you can try to apply the FunctionExpand command to the DifferenceRoot object. If it is able to find a closed form of the sequence, then the Limit might be able to find an exact symbolic limit. To find a numerical approximation, you can use the SequenceLimit command. In general, it does not guarantee to give the correct result, but if your sequence ...


11

You should use RecurrenceTable: N@RecurrenceTable[{a[n + 1] == 2 a[n] b[n]/(a[n] + b[n]), b[n + 1] == Sqrt[a[n + 1] b[n]], a[0] == 2 Sqrt[3], b[0] == 3}, {a[n], b[n]}, {n, 1, 10}] N@Pi with result: {{3.21539, 3.10583}, {3.15966, 3.13263}, {3.14609, 3.13935}, {3.14271, 3.14103}, {3.14187, 3.14145}, {3.14166, 3.14156}, {3.14161, 3.14158},...


11

Your recurrence is an instance of the logistic map $$x_{n+1} = r x_n(1 - x_n).$$ The first sentence of the Solution in some cases section in the link above says: The special case of r = 4 can in fact be solved exactly, as can the case with r = 2; however the general case can only be predicted statistically. Mathematica is indeed able to solve the $r=2$ ...


11

Unlike many other functions that use an iterator, RecurrenceTable does not have the attribute HoldAll. Presumably, this indicates a different evaluation sequence (Evaluation in Iteration Functions). Grid[ Prepend[ {#, Attributes[#]} & /@ {Do, Product, Sum, Table, RecurrenceTable}, Style[#, Bold] & /@ {"Function", "Attributes"}], Frame -> ...


11

Proper memoization is the key: Q[1] = 1; Q[2] = 1; Q[n_] := Q[n] = Q[n - Q[n - 1]] + Q[n - Q[n - 2]] data = Table[Q[n]/n, {n, 10000}]; // AbsoluteTiming {0.051379, Null} ListPlot[data]


10

Edit If the corresponding model is to be used for a lifespan L, then RSolve works: RSolve[{y[0] == y0, y[1 + t] == (1 + 1/L (-1 - t)) y[t]}, y[t], t] (* {{y[t] -> -(-1 + L) (-L)^-t y0 Pochhammer[2 - L, -1 + t]}} *) I suspect that RSolve fails because 1 - (t + 1)/100 becomes zero. Although mathematically that shouldn't matter, algorithmically it may. ...


10

Alternatively, and perhaps more directly, use RSolve[{a[n] == Sqrt[3] + 1/2 a[n - 1], a[0] == 1}, a[n], n] (* {{a[n] -> 2^-n (1 - 2 Sqrt[3] + 2^(1 + n) Sqrt[3])}} *) Limit[a[n] /. %[[1]], n -> Infinity] (* 2 Sqrt[3] *)


10

You can use Mathematica to prove the induction step in a proof that $a_n, b_n \rightarrow L$ for some limit $L$ by showing that $a_n$ is decreasing and $b_n$ increasing toward each other and that the difference goes to zero. To show $L = \pi$ you would need a definition of $\pi$ you can relate to this sequence. For instance, if you define $\pi$ to be the ...


10

Your can use MatrixPower for this example: f[n_] := MatrixPower[{{.5, -.6}, {.75, 1.1}}, n].{2, 0} f /@ Range[0, 5] yields: {{2., 0.}, {1., 1.5}, {-0.4, 2.4}, {-1.64, 2.34}, {-2.224, 1.344}, {-1.9184, -0.1896}}


9

There is a hidden option to change the variables used by DifferenceRoot[]: SetSystemOptions["HolonomicOptions" -> {"IteratorSymbol" -> k, "SequenceSymbol" -> ϕ}]; after which DifferenceRootReduce[n!, n] (* DifferenceRoot[Function[{ϕ, k}, {(-1 - k) ϕ[k] + ϕ[1 + k] == 0, ϕ[0] == 1}]][n] *) For the analogous case of DifferentialRoot[], the ...


9

This is how you define the function $g$ in Mathematica: g[a_, x_] = a*Sin[π*x]; Making a color-plot of the function $g(g(g(x)))-x$ in the $a$-$x$ plane and indicating the zero-contours, we see the lobes of 3-cycles: DensityPlot[g[a, g[a, g[a, x]]] - x, {a, 0, 1}, {x, 0, 1}, MeshFunctions -> {#3 &}, Mesh -> {{0}}, PlotPoints -> 100] Now we ...


8

This is definitely a bug. It seems that the problem is caused because you increment the argument by two (your recurrence has c[k] and c[k+2]). This results in two cases (even and odd), as can be seen by your Maple output, and Mathematica apparently does not know how to deal with this properly. One way around this is to substitute the variables so that the ...


8

Here is a modified version of your code. On my PC it completes your example run in about 5 seconds. I won't try to describe every change but will point out the major features. Some of the changes are stylistic rather than performance-based. This is not a criticism of your style but a reflection of the way I broke the original code down in order to ...


8

The general formula can be derived as follows. First@RSolve[{prin[n] == (1 + int) prin[n - 1] - pay, prin[0] == loan}, prin[n], n]; First@Solve[(prin[n] /. %) == 0, pay] (* {pay -> (int (1 + int)^n loan)/(-1 + (1 + int)^n)} *) where pay is the payment per period, int is the interest per period, and loan is the original principal. For the example given ...


7

Started as a comment to the OP. Turning this into an answer was requested. The OP asked for Mathematica help solving the recurrence. I suggested recognizing that the recurrence was the sum of angles formula for tangent, then going from there. Consequently, this answer uses Mathematica to solve this problem using only one "external" ingredient right at ...


7

When building a simulation like yours you should test the performance of the individual components before incorporating them into the simulation. That is, you should know the cost of the components as well as the values they return. Here is a simple example based on your code. You use Clip in a couple of places to limit values from below. That is suspect, ...


7

Several changes are required to obtain the desired results. First, the syntax error mu[n] == mu must be replaced by mu[n] = mu. Next, initial conditions must be provided for the recurrence in n: mu[0] := mu /. FindRoot[f[1/2, mu, 2^0] == 1/2, {mu, 0.9}] mu[1] := mu /. FindRoot[f[1/2, mu, 2^1] == 1/2, {mu, 0.9}] mu[2] := mu /. FindRoot[f[1/2, mu, 2^2] == 1/...


7

This is a bug in RSolve. The problem appears to be caused by a missing check while applying a transformation to solve this example, which should have returned unevaluated with the available methods. Sorry for the confusion caused by the incorrect answer.


7

Eliminating b from the first two expressions in the question and redefining n yields c[n] == 2 c[n - 1] + c[n - 2] + 2 a[n - 1] Although RSolve can handle this expression with a as given in the question, it seems better to Rationalize it. 794367 10^-7 (460811 10^-6)^n+255972 10^-6 (-675131 10^-6)^n+664591 10^-6 (321432 10^-5)^n Then, Simplify[RSolve[{c[...


7

This is not a complete answer. It is a long comment about what I found when I looked at the inspectable source code of this function using various spelunking tools (e.g.. Spelunking package, GeneralUtilities`PrinteDefinitions, etc.). Method parsing looks unusual First, I wonder if parsing for the Method option is buggy. First it does this internally: $...


Only top voted, non community-wiki answers of a minimum length are eligible