title: Incorrect evaluation for Thue-Morse signed harmonic series

text: I would like to evaluate $$s = 1 - \frac{1}{2} - \frac{1}{3} + \frac{1}{4} - \frac{1}{5} + \frac{1}{6}+\frac{1}{7}-\frac{1}{8} - ... + \frac{(-1)^{\textrm{binary digit sum}(n-1)}}{n} + ... $$ where the signs follow the Thue-Morse sequence. When I try to evaluate `Sum[ (-1)^Total[ IntegerDigits[ n - 1, 2]]/n, {n, 1, ∞}]` using _Mathematica 7_ I get > >     Log[2] >   This is incorrect. $$\log( 2) = \sum_{n=1}^\infty \frac{(-1)^{n-1}}{n} \approx 0.693147$$ but $$s \approx 0.39876108810841881241$$ I checked this using _C#_ and finite sums within _Mathematica_. This value is not in the inverse symbolic calculator, but since that is incomplete I still have some hope for an analytic summation. The difference between these series can be expressed as a sum of positive terms $$\frac{2}{3\times 4} + \frac{2}{5\times 6} + \frac{2}{9\times 10} + \frac{2}{15\times 16}+\ldots+\frac{2}{n \times (n+1)}+ \ldots$$ where the sum is over all odd $n$ with even binary digit sum, so it is definitely the case that $s \lt \log( 2).$ I also tried `Sum[ (-1)^Total[ IntegerDigits[ n - 1, 2]]/n, { n, 1, ∞, 2}]` which equals $\frac{3}{2}s \approx 0.59814163216262821861$, but which instead evaluates to > >     Sum::div: Sum does not converge >   >     Sum[ (-1)^(1 + n)/n, { n, 1, ∞, 2}] >   So, is this a known bug?   What sorts of related sums are affected?   Is there a way to get the right answer from _Mathematica_?


title: CDF of TransformedDistribution Not Returned And No Error Message Produced

text: I am unable to get Mathematica to return the CDF of a distribution produced by TransformedDistribution. It "executes" for 6 hours on my machine and then returns a blank plot. Could someone execute the code below and explain why this is happening and offer a solution ? Is there something wrong with my code ?               beg = AbsoluteTime[];     dist1 = ExponentialDistribution[0.1];     pdf1 = PDF[dist1, x];     cdf1 = CDF[dist1, x];     plot1 = Plot[pdf1, {x, 0, 100}, Frame -> True,      FrameLabel -> {"X", "PDF"}, PlotLabel -> "PDF of X (dist1)",      PlotStyle -> Red, PlotRange -> All]     plot2 = Plot[cdf1, {x, 0, 100}, Frame -> True,      FrameLabel -> {"X", "CDF"}, PlotLabel -> "CDF of X (dist1)",      PlotStyle -> Red, PlotRange -> All]          dist2 = NormalDistribution[53, 2.4];     pdf2 = PDF[dist2, x];     cdf2 = CDF[dist2, x];     plot3 = Plot[pdf2, {x, 0, 100}, Frame -> True,      FrameLabel -> {"X", "PDF"}, PlotLabel -> "PDF of X (dist2)",      PlotStyle -> Red, PlotRange -> All]     plot4 = Plot[cdf2, {x, 0, 100}, Frame -> True,      FrameLabel -> {"X", "CDF"}, PlotLabel -> "CDF of X (dist2)",      PlotStyle -> Red, PlotRange -> All]     end = AbsoluteTime[];     executiontime = (end - beg)/60  (* minutes *)      Here is the output, so far so good. ![Output1](http://i.stack.imgur.com/y9KHN.jpg) Now here is where the problem is.               beg = AbsoluteTime[];     dist3 = TransformedDistribution[     a + b, {a \[Distributed] dist1, b \[Distributed] dist2}]     pdf3 = PDF[dist3, x];     cdf3 = CDF[dist3, x];     plot5 = Plot[pdf3, {x, 0, 100}, Frame -> True,      FrameLabel -> {"X", "PDF"}, PlotLabel -> "PDF of X (dist3)",      PlotStyle -> Red, PlotRange -> All]     plot6 = Plot[cdf3, {x, 0, 100}, Frame -> True,      FrameLabel -> {"X", "CDF"}, PlotLabel -> "CDF of X (dist3)",      PlotStyle -> Red, PlotRange -> All]     end = AbsoluteTime[];     executiontime = (end - beg)/60  (* minutes *)      And here is the output. ![Ooutput2](http://i.stack.imgur.com/YwOkm.jpg) Thanks for any help.


title: How to plot this function?

text: I'd like to plot the function $f$ in this answer http://math.stackexchange.com/a/788818/66096 > Let $h(x) = \begin{cases} e^{ -\frac{1}{1 - (1-2x)^2}} & \mbox{ for } 0<x < > 1\\ 0 & \mbox{ otherwise} \end{cases}$ Then let $$g(x)=\sum_{n=1}^\infty > h(n^2(x-n))=\begin{cases}h(n^2(x-n))&\text{if }n\le x<n+1, n\in\mathbb > N\\0&\text{otherwise}\end{cases}$$ And let $f(x)=\int_0^xg(t)\,\mathrm dt$. I only know the basics of _Mathematica_ , so I've been trying to plot piecewise functions without success. Here is what I got               h[x_] = Piecewise[{{E^(-(1 - (1 - 2 x)^2)^(-1)), 0 < x < 1}, {0,          x > 1}, {0, x < 0}}]      But how do I deal with $g$ ? I could input partial sums with the series expression, but at the same time, I could save computations by just defining it piecewise.


title: Which functions automatically generate new functions?

text: In his video http://www.wolfram.com/broadcast/video.php?channel=362&video=1643, Wolfram gives the following example: the output from               dict = Nearest[WordData[]]      is a `NearestFunction` object -- a function that may then be evaluated at an input such as:               dict["eleephannt"]      > {elephant} How can one determine which other functions `func` similarly generate a `funcFunction`? (Not sure if the Title correctly captures the meaning of my question!)


title: Deciding whether a semi-algebraic set is empty. Which command is best?

text: [Update] I was able to prove that the "larger" set $ 1\ge x \ge y,w \ge 0$ and $p(x,y,w)>0 $ is empty with another computer algebra system in my slow laptop in few minutes but in Mathematica even if I use a Linux cluster with several processor and lots of memory, I can not get the answer with in one hour (will try with a day later). The same problem also happens with the other system if I use Cylindrical Algebraic Decomposition. My conclusion is that it is not Mathematica's fault but rather Cylindrical Algebraic Decomposition is very inefficient. It seems that Regular Chains (triangular decomposition) methods may perform better. [end of update] Question was: We have system of polynomial inequalities in three variables: $$ 1> x > y > w > 0\,\text{ and } p(x,y,w)>0 $$ where $p$ is a polynomial with integer coefficients with order 9. Mathematica has several different commands I could use to try to check if this set is empty: Reduce, FindInstances, CylindricalDecomposition,SemialgebraicComponentInstances, etc... I looked at the documentation but it does not say which algorithms each command uses. I need to choose the most efficient one. In Maple, running 8 processors with 8G, I have a code already running for 8 hours and no answer. Edit [due to Artes' comment] the expression for $p$ is the following:                f[x_,y_]:=(x^2+4*x*y+y^2)/((x+y)*(2*x+y)*(x+2*y))      p:= -1 * Numerator[Together[Simplify[(f[1,x]+f[y,z]-f[1,y]-f[x,z])/(2*(x-y)*(1-z))]]]


title: Decomposing a large text file, with alternating headers and text, into an array of headers and an array of text segments

text: I have a very large text file, `test.txt`, in the following format:               (Same specific character or word like '^' or '*') Title 1 (end-of-line)          Text 1 (with paragraphs markers, end-of-line markers spaces, all kinds of stuff)          (Same specific character or word like '^' or '*') Title 2 (end-of-line)          Text 2 (with paragraphs markers, spaces, all kinds of stuff)          (Same specific character or word like '^' or '*') Title 3 (end-of-line)          Text 3 (with paragraphs markers, spaces, all kinds of stuff)          ...      I'd like to make two arrays, one corresponding to the strings `{ "Title 1", "Title 2", "Title 3", ...}`, and the other corresponding to the strings `{ "Text 1", "Text 2", "Text 3", ...}`.   Is there a simple one-liner to do this? Here's a specific test example. Here each line starting with `^` should be a `"Title"` (of which there are three), and the material between the `^` corresponds to `"Text"` (of which there are three sections). Notice that `^` only appears as the first character in the `"Title"` string, and that each title is finished with an end-of-line, and finally that each `"Text"` section can consist of multiple lines of strings.


title: Obtaining the envelope of an oscillating function and perform fitting

text: There is an elegant way to obtain the envelope of an oscillation function explained here http://mathematica.stackexchange.com/a/27854/15966 My question is how to modify the code below so I can fit data to the obtained envelope with a parameter. How do I introduce a parameter into the module so it can be used later for fitting when the result of the module is an interpolation function. Lets' generate some data               a = RandomReal[{20, 30}];       b = RandomReal[{5, 10}];     data = Table[{x, Cos[a*x]^2*Exp[-b*x]}, {x, 0, 1, 0.001}];     ListPlot[data, PlotRange -> All]      The code generates something like this: ![data](http://i.stack.imgur.com/qmrLE.jpg) Now the module below produces envelope for the given oscillating function               FunctionEnvelope[f_, {t_, a_, b_}, n_: 40] :=       Module[{seeds, x, y, points, progress = 0, tempf, union},        seeds = Rescale[Range[0, 1, 1/n] + 1/(2 n), {0, 1}, {a, b}];       points =         Last[Last[          Reap[Monitor[            Do[Quiet@              Check[progress++; {y, x} = FindMaximum[Abs[f], {t, x0}];               x = t /. x; If[a <= x <= b, Sow[{x, y}]], Null], {x0,               seeds}],             ProgressIndicator[progress, {0, Length[seeds]}]]]]];       union[] :=         points =          Union[points,           SameTest -> (Abs[First[#1] - First[#2]]/               Replace[Max[Abs[First[#1]], Abs[First[#2]]],                 u_ /; u == 0 :> 1] < 10^-6 &)];       union[];       tempf = Interpolation[points];       points = Quiet[Join[{{a, tempf[a]}}, points, {{b, tempf[b]}}]];       union[];       Interpolation[points]]      It can be used as follows:               q = 50;     f[t_] := E^(- 10*Abs[t]) * Cos[1000*t]/2 +        E^(-10*Abs[t])*Cos[(1000 + q)*t ]/2     g = FunctionEnvelope[f[x], {x, 0, 1}, 200];     Plot[{f[x], g[x]}, {x, 0, 1}, PlotRange -> All]      ![enter image description here](http://i.stack.imgur.com/C4AIu.jpg) My question is how can I modify the module above to make it a function of the parameter q, so that the resulting interpolation function can be used to fit the data to get the optimal parameter q. (The function f[t] is in principle can be moved to the definition of the module, I'm not going to change it.) Thanks!


title: Importing adjacency matrix from file of 0/1 arrays

text: Is it possible to import adjacency matrix for graph representation, from simple text file with NxN (rows x columns) of 0s/1s:               0 1 0 0     1 0 1 1     0 1 0 0     0 1 0 0


title: Numerical Instability?

text: **Current Status** Please skip ahead to Update 2: the key questions are now: 1/ Why does Mma generally fail to find the minimum for a well behaved function when the function does not have infinite precision but the specified precision seems more than sufficient to the task 2/ Why does Mma fail to find the minimum when the precision of calculations is increased beyond the point at which it _is_ first able to minimise the function The original question need only be consulted for background, if desired. **End Current Status** Further to this initial question on controlling precision and with the benefit of rasher & Mr Wizard's (previous) input I now have control over my calculation precision (I think, but I am still very much a beginner so I beg your collective indulgence...). However, the results are wildly unstable as I increase precision, though to my eyes the equations look well-behaved in the regions of interest. I only took an interest in the precision when the machine precision results did not seem physically reasonable. I am seeking maxima and minima for four variables, t, x, y & z and have - I think - good reason (see separate section below) to believe that for x, y, and z the minima should be equal, and so should the maxima. What I find for ~machine precision is (focussing on the minima only):               x = -0.222222222222222     y = -1.242253153494817     z = -1.222222222222222      These are clearly unequal with differences of the same magnitude as the values. (NB. When 3D plotting the resulst of the tranformations, these erroneous values lead to the transformed sphere being truncated when the plot range is determined in terms of these extrema, so when Mma is plotting the equations it is generating values significantly different to those it comes up with during minimization etc.] With precision set to 50 the results are very different:               x = -1.2619796324000608024205195519983014313647470787550     y = -1.2619796324000608024205195519983014313647470787550     z = -1.2619796324000608024205195519983014313647470787550      Clearly these three are identical to the specified precision, but increase it further to 100 and one obtains               x = -1.261979632400060802420519551998301431364747078755025069028958918198601669661636873543419739259161731     y = -1.261979632400060802420519551998301431364747078755025069028958918198601669661636873543419739259161731     z = -1.222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222222      Now x = y still, and the results agree with the previous output to 50 digits, but z has changed in the 2nd decimal place. There are of course similar problems with the maxima - and it should be noted that the maxima probably ought to have the same magnitudes as the minima, which they don't. Can anyone explain what is going on and how can I reliably estimate the desired minima and maxima? **Code** Note also the output of the MinValue function... as an incidental question, why does it not achieve the same prevision as Minimize? (Thanks to belisarius for prettifying the code) To the best of my knowledge the equations are correct for the physics.               $PreRead = (# /. s_String /; StringMatchQ[s, NumberString] &&              Precision@ToExpression@s == MachinePrecision :>   s <> "`100." &);     Block[{γ, vx, vy, vz, t0, x0, y0, z0, t1},      t0 = γ (t1 - vx Sin[θ] Sin[ϕ] - vy Cos[θ] Sin[ϕ] - vz Cos[ϕ]);      x0 = (-γ t1 vx - γ^2 t1 vx +    Sin[θ] Sin[ϕ] + γ Sin[θ] Sin[ϕ] +            γ^2 vx^2 Sin[θ] Sin[ϕ] + γ^2 vx vy Cos[θ] Sin[ϕ] + γ^2 vx vz Cos[ϕ])/(1 + γ);      y0 = ((1 + γ + γ^2 vy^2) Cos[θ] Sin[ϕ] - γ vy (t1 + γ t1 - γ vx Sin[θ] Sin[ϕ]) +            γ^2 vy vz Cos[ϕ])/(1 + γ);      z0 = (γ vz (-(1 + γ) t1 + γ vx Sin[θ] Sin[ϕ] + γ vy Cos[θ] Sin[ϕ]) +            (1 + γ + γ^2 vz^2) Cos[ϕ])/(1 + γ);      vx = vy = vz = 0.8/Sqrt[3];       γ = 1/Sqrt[1 - (vx^2 + vx^2 + vx^2)];       Print["γ = ", γ];      t1 = 0;      Print["Precision of t0 = ", Precision[t0], " of x0 = ",              Precision[x0], " of y0 = ", Precision[t0], " of t0 = ", Precision[z0]];      Print["Precision of Cos[ϕ] = ", Precision[Cos[ϕ]]];      Print["MinValue of t0 = ", MinValue[{t0, 0. <= θ < 2 Pi, -Pi <= ϕ <= Pi}, {θ, ϕ}, Reals]] ;      Print[{Minimize[{t0, 0 <= θ < 2 Pi, -Pi <= ϕ < Pi}, {θ, ϕ}, Reals],             Maximize[{t0, 0 <= θ < 2 Pi, -Pi <= ϕ < Pi}, {θ, ϕ},  Reals]}];      Print[{Minimize[{x0, 0 <= θ < 2 Pi, -Pi <= ϕ < Pi}, {θ, ϕ}, Reals],             Maximize[{x0, 0 <= θ < 2 Pi, -Pi <= ϕ < Pi}, {θ, ϕ},  Reals]}];      Print[{Minimize[{y0, 0 <= θ < 2 Pi, -Pi <= ϕ < Pi}, {θ, ϕ}, Reals],             Maximize[{y0, 0 <= θ < 2 Pi, -Pi <= ϕ < Pi}, {θ, ϕ},  Reals]}];      Print[{Minimize[{z0, 0 <= θ < 2 Pi, -Pi <= ϕ < Pi}, {θ, ϕ}, Reals],              Maximize[{z0, 0 <= θ < 2 Pi, -Pi <= ϕ < Pi}, {θ, ϕ}, Reals]}];      ]      **Reasons to believe extrema should be equal** The equations given in the code sample describe the Lorentz transformation of a sphere from one frame of reference to another. For simplicity, albeit in the general case, the sphere is moving with equal velocity in all three spatial directions. As is well known, according to Special Relativity the sphere should be contracted _along the direction of motion_ - in this case in the direction {1,1,1}. Consequently by symmetry of the sphere and equality of motion with respect to x, y, z, the x, y, z minima (the coordinate distances to the "back" of the sphere, assuming it is moving away) should all be the same. Likewise for the maxima. The fact that the sphere is defined in terms of parametric equations in polar coordinates (SinSin, CosSin and Cos) obscures this, but again by symmetry the orientation of the axes in the sphere's rest frame shouldn't make any difference either. The fact that the x, y, z extrema are _bigger_ and not smaller than 1 (sphere radius) is due to the fact that they are the transformed values in the other frame and are not simultaneous, so do not represent the Lorentz contraction directly: what is actually being calculated here is the spatial envelope over time in the other frame. **Update** I have tried to apply MichaelE2's "uncertainty" function given in answer to another question here and obtain the following results (focusing on z0 only for brevity). Adding the following code to that already given (declare extra block vars as necessary - as noted in chat, I have been unable to get the M.SE toolbar to work for me, so apologies for the [Gamma], etc. again)               minz = Minimize[{z0,         0 <= \[Theta] < 2 Pi, -Pi <= \[Phi] < Pi}, {\[Theta], \[Phi]},        Reals]; maxz =       Maximize[{z0,         0 <= \[Theta] < 2 Pi, -Pi <= \[Phi] < Pi}, {\[Theta], \[Phi]},        Reals];     (* unset all variables again so that uncertainty can take the \     derivative... *)     vx =. ; vx =.; vz =.; \[Gamma] =.; t1 =.;     Print[minz[[2, 1]], maxz[[2, 2]]];     Print["Uncertainty in min z0 = ",        uncertainty[        z0  /.  \[Gamma] -> 1/Sqrt[1 - (vx^2 + vy^2 + vz^2)], {minz[[2,           1]], minz[[2, 2]], vx -> 0.8/Sqrt[3], vy -> 0.8/Sqrt[3],          vz -> 0.8/Sqrt[3], t1 -> 0}]];      I obtain               for precision =  20, Uncertainty in min z0 = 4.61791*10^-20                   =  50, Uncertainty in min z0 = 8.63608*10^-50                   = 100, Uncertainty in min z0 = 4.55654*10^-100      This suggests to me, naively perhaps!, that the answers previously given should have agreed to the precision specified, and since they did not I remained baffled as to what is going on. **Update 2 (21 Mar) - the plot thickens** Summary: having had my attention drawn to the form of z0 by @MichelE2's comment an error was detected: phi, which is a spherical polar coordinate was running from -Pi to Pi, instead of 0 to Pi, which created double minima. Corrected, there is now only one minimum, as this plot of z0 shows nicely ![Good plot](http://i.stack.imgur.com/ovNOn.jpg) We now find that if z0 has infinite precision, Minimize solves analytically and N[] of the solutions can be obtained to specified precision - and the answers agree with the application of NMinimize [why are the specs for those functions so different???] _**However..._** we also find that subjecting the variables to careful precision control, both Minimize and NMinimize fail to find the correct answer unless prec = 45 to 47. It can be seen that they are getting stuck in the corner at phi = Pi, which for the well behaved function z0 does not seem reasonable. Automatic plotting of the region of interest fails (plot range min and max become the same) but a manual plot range shows this ![Bad plot](http://i.stack.imgur.com/H8hSg.jpg) which doesn't look nice at all. [NB note that I have tweaked WorkingPrecision to keep it at least consistent with prec] Given the reasonable behaviour of z0, what is the reason that Mma is generally unable to (N)Minimize[z0] and when precision has "increased enough" to do it, how can it fail again at _higher_ precision. This does not seem right at all. Code for this update follows...               Block[{\[Gamma], vx, vy, vz, t0, prec, mint, maxt, minx, maxx, miny,        maxy, minz, maxz, K, plotMinTheta, plotMinPhi, },      z0 = (\[Gamma] vz (-(1 + \[Gamma]) t1 + \[Gamma] vx Sin[\[Theta]] \     Sin[\[Phi]] + \[Gamma] vy Cos[\[Theta]] Sin[\[Phi]]) + (1 + \[Gamma] \     + \[Gamma]^2 vz^2) Cos[\[Phi]])/(1 + \[Gamma]);      prec = 44; K = 8/10;       t1 = 0; (* prec = 45 to 47 is 1st range for which part 2 gives \     correct answer, then from 103 to..? *)      vx = vy = vz = (8/10)/Sqrt[3]; \[Gamma] =        1/Sqrt[1 - (vx^2 + vx^2 + vx^2)]; t1 = 0;      Print["For precision of t1 = ", Precision[t1], ", of K = ",        Precision[K], ", of \[Gamma] = ", Precision[\[Gamma]], ", of z0 = ",        Precision[z0]];      minz = Minimize[{z0, 0 <= \[Theta] < 2 Pi,          0 <= \[Phi] < Pi}, {\[Theta], \[Phi]}, Reals];      plotMinTheta = N[2 (  \[Theta] - Pi) /. minz[[2]][[1]], prec];       plotMinPhi =        N[(2 \[Phi] - Pi) /.  minz[[2]][[2]],         prec]; (* to centre the minimum found *)      Print["N[Min,$MachinePrecision,] = ", N[minz, $MachinePrecision]];      Print["N[Min,", prec, "] = ", N[minz, prec]];      Print["NMinimize[z0] = ",        NMinimize[{z0, 0 <= \[Theta] < 2 Pi,          0 <= \[Phi] < Pi, \[Theta] \[Element] Reals, \[Phi] \[Element]           Reals}, {\[Theta] , \[Phi]}, WorkingPrecision -> (89/100)*prec,         MaxIterations -> 200]];       Print[Plot3D[        z0, {\[Theta], plotMinTheta, 2 Pi}, {\[Phi], plotMinPhi, Pi},         PlotRange -> Automatic,  AxesLabel -> {"theta", "phi"},         WorkingPrecision -> (89/100)*prec, Mesh -> 10]];           vx = vy = vz = N[K, prec]/Sqrt[3]; \[Gamma] =        1/Sqrt[1 - (vx^2 + vx^2 + vx^2)]; t1 = 0;      Print["For precision of t1 = ", Precision[t1], ", of K = ",        Precision[K], ", of \[Gamma] = ", Precision[\[Gamma]], ", of z0 = ",        Precision[z0]];      Print["NMinimize[z0] = ",        NMinimize[{z0, 0 <= \[Theta] < 2 Pi,          0 <= \[Phi] < Pi, \[Theta] \[Element] Reals, \[Phi] \[Element]           Reals}, {\[Theta] , \[Phi]}, WorkingPrecision -> (89/100)*prec]];       minz = Minimize[{z0, 0 <= \[Theta] < 2 Pi,          0 <= \[Phi] < Pi}, {\[Theta], \[Phi]}, Reals];       plotMinTheta = N[2 (  \[Theta] - Pi) /. minz[[2]][[1]], prec];       plotMinPhi =        N[(2 \[Phi] - Pi) /.  minz[[2]][[2]],         prec]; (* depending on prec, will not work when phi of minimum \     \[Equal] Pi ! *)      Print["N[Minimize[z0],", prec, "] = ", N[minz, prec]];      Print[Plot3D[        z0, {\[Theta], plotMinTheta, 2 Pi}, {\[Phi], plotMinPhi, Pi},         PlotRange -> Automatic,  AxesLabel -> {"theta", "phi"},         WorkingPrecision -> (89/100)*prec, Mesh -> 10]];      Print[Plot3D[        z0, {\[Theta], N[6.28317, prec], 2 Pi}, {\[Phi], 0.999999 Pi, Pi},         PlotRange -> Automatic, AxesLabel -> {"theta", "phi"},         WorkingPrecision -> (89/100)*prec, Mesh -> 10]];      ]


title: Apply Collatz Funtion to solve and print the answer

text: Good day to all! The objective of my task is apply the collatz function for a number `x1` and print the exact number of iterations until `s1`. I want to create a function using `Module`:               CollatzFunction[x1_, s1_]:= Module[{.......};      Print iterations to `s1`-> "number of iterations" I barely have an idea using `Nest[]`, but I don't know how to apply `Nest` until the `s1` number is reached. Please, help in this task, it's just a step towards creating a more complicated function. Thanks.


title: Do not understand expression returned by FullSimplify

text: Why _Mathematica_ returns (sec x)^2?               (Tan[x]^2 + 1 - Sec[x]^2)/(Sin[x]^2 + Cos[x]^2 - 1) // FullSimplify      > >     Sec[x]^2 >   Because this? Then differentiate？               Integrate[(Tan[x]^2 + 1 - Sec[x]^2)/(Sin[x]^2 + Cos[x]^2 - 1), x]      > >     Tan[x] >   And `WolframAlpha` returns several answers ![enter image description here](http://i.stack.imgur.com/93AdK.png) Why these happened？It's strange。


title: Frontend cursor movement speed depends on the position in large notebook file

text: I run Mathematica v9 on my Atom laptop (don't ask) and I got used to the fact that everything is slower. Except for one annoying thing: cursor movement (in text) is slower than the keyboard key repeat rate. In other editors I am used to hold down the up/down/left/right keys and I expect the cursor to stop when I release the key. In Mathematica Frontend this is not the case, keyboard buffer gets filled with keyboard events which get executed by the Frontend way later than the release of the key. This annoying behavior trained me to use only single presses of the movement keys (or mouse clicks at the right position). You can imagine, scrolling through a lot of code in this way can be quite painful. What is even more interesting is that the speed with which Frontend executes the cursor movements depends on the overall size of the notebook file and the current relative position of the cursor the file. The same approximately holds also for edits. Having a large file with a lot of graphics, editing text or moving the cursor at the top of such notebook becomes extremely slow on such weak netbooks. On the contrary, adding code or moving cursor at the end of such file is quick enough not to be noticed. I have been following the evolution of Mathematica "speed" since v5 and I have to say that there was huge drop of performance when Wolfram switched the backend language to java and ever since the performance of the kernel or the Frontend has been getting worse (the fact substantially masked by the CPU speed increases during the "frequency wars"; now this Moore's law is over...). It also seems that all the incredible performance gains (with JIT technology and loop/execution tracing) in other interpreted languages never really reached Mathematica. Maybe Wolfram should dedicate one major version bump only to address these performance issues otherwise (this great) language may put itself into oblivion. Nevertheless, the behavior described above is quite odd, i.e. it seems edits at the top are $\mathcal{O}($file size$)$ while at the bottom they are more like $\mathcal{O}(1)$. Are there any preferences/options that can be tweaked in order to reduce the load Frontend uses for rendering?


title: Ad hoc graphics primitives-like objects

text: `Sphere` is one of the three-dimensional graphics primitives available in Mathematica and can be easily used to created very useful images. For instance, in the figure below I created three images of an ellipsoidal object that is lit from a fixed direction and it is then just rotated around the z-axis.               GraphicsRow[         Table[             Graphics3D[             Rotate[Scale[{GrayLevel[.7], Sphere[]}, {1, 1.5, 1}, {0, 0, 0}],              i Degree, {0, 0, 1}], PlotRange -> {{-2, 2}, {-2, 2}, {-2, 2}},              AspectRatio -> 1, Background -> Black, Boxed -> False,              ViewPoint -> Front, SphericalRegion -> True,              Lighting -> {{"Directional", White, ImageScaled[{0, 0, 1}]}}         ], {i, {0, 45, 90}}]]      ![ellipsoidalObject](http://i.stack.imgur.com/qSTdb.png) Suppose that instead of using the ellipsoidal object I would like to use an object like this bump:               Plot3D[Exp[-(x^2 + y^2)], {x, -2, 2}, {y, -2, 2}]      ![bump](http://i.stack.imgur.com/Hjfn5.png) Is there a way to create a "graphics primitive" object of the bump, so as to just use it instead of `Sphere` in the above example? I suppose there should be a way to extract all the polygons specifying the shape of the bump, but I couldn't find a way to do it. Any other approach is also more than welcome.


title: Evaluating Polynomials at Grid Points

text: I am continuing my **quest on B-splines**. The code below builds a 5x5 matrix out of B-splines, using the BSplineBasis[] routine. I now want to evaluate the polynomials that are stored in each matrix element. The problem that I am facing is, however, that each matrix element corresponds to a different position on my xy-grid. For instance, B[[2,2]] corresponds to the grid-pint x=-1, y = -1 . B[[2,3]] corresponds to x=-1, y = 0 and so forth. At the moment I do x = smth and then evaluate the polynomial manually. This works for a small system, but it is not optimal for larger problems. **So my question is** , is there some notation that does something like **B[[2,3]](put out the value of this matrix element at (x = smth,y = smth))**? I tried it with Evaluate[], but didn't get it to work.               knots = {-1, -1, -1 , -1, 0, 1, 1, 1, 1};     B1 = Table[        D[BSplineBasis[{3, knots}, i, x], {x, 2}] BSplineBasis[{3, knots},           j, y], {i, 0, 4}, {j, 0, 4}];     B2 = Table[        BSplineBasis[{3, knots}, i, x] D[          BSplineBasis[{3, knots}, j, y], {y, 2}], {i, 0, 4}, {j, 0, 4}];     B = B1 + B2;     (*set boundary conditions*)     B[[1]] = {0, 0, 0, 0, 0};     B[[5]] = {0, 0, 0, 0, 0};     B[[All, 5]] = {0, 0, 0, 0, 0};     B[[2, 1]] = BSplineBasis[{3, knots}, 2, x];     B[[3, 1]] = BSplineBasis[{3, knots}, 3, x];     B[[4, 1]] = BSplineBasis[{3, knots}, 4, x];      Here is a plot of what I'm doing, just for nurturing the eye.               Plot3D[Evaluate[B], {x, -1, 1}, {y, -1, 1}, Mesh -> None,       PlotStyle -> Directive[Opacity[0.8]], PlotRange -> {-2, 2},       ColorFunction -> Function[{x, y, z}, Hue[z]]]      ![enter image description here](http://i.stack.imgur.com/jcico.png)


title: ContourStyle for a particular contour line in ContourPlot

text: I have the following command:               ContourPlot[(-1 + x)^2 (-1 + y) + x^2 y,     {x, 0, 1}, {y, 0, 1},     ContourLabels -> True]      which gives me ![enter image description here](http://i.stack.imgur.com/FD2ug.png) I want to highlight the 0-contour line, e.g. make it thicker, dashed, red, etc., but _Mathematica_ 's documentation on `ContourPlot` does not seem to provide a way to specify the style for particular contour lines. How should I specify the style for particular contour lines?


title: CleanSlate[] not working?

text: In[1]:= Needs["Utilities`CleanSlate`"]          In[2]:= x := 4          In[3]:= x          Out[3]= 4          In[4]:= CleanSlate[]          During evaluation of In[4]:=   (CleanSlate) Contexts purged: {Global`}          During evaluation of In[4]:=   (CleanSlate) Approximate kernel memory recovered: 0 Kb          Out[0]= {"Utilities`CleanSlate`", "PacletManager`", "QuantityUnits`", \ "WebServices`", "System`", "Global`"}          In[1]:= x          Out[1]= 4      When I use CleanSlate[] the variable x should be cleared from memory right? Why does it still remain afterwards?


title: How do I solve $1 - (1 - (A x)^2)^\frac{3}{2} - B(1 - \cos(x))= 0$?

text: Consider the equation $$ 1 - (1 - (A x)^2)^\frac{3}{2} - B(1 - \cos(x)) = 0 $$ where $A,B \in \mathbb{R}$ are constants. What is the analytic solution?


title: webMathematica on Linux

text: I have been using _Mathematica_ and _webMathematica_ (with application server) in Windows and it was working fine. Recently we changed to _Linux_ , so I have to install both _Mathematica_ and _webMathematica_ (with application server) in Ubuntu-12.04. I manage to install _Mathematica_ in Ubuntu-12.04 and get it work. Now I have to install _webMathematica_ in Ubuntu-12.04 and deploy it in Application server. This is what I tried base on the link. I tried to install Xvnc, I downloaded vnc for Linux from site. Am not root of the Ubuntu OS, So I have to login as a user and did               cp vncviewer vncserverui vncpasswd Xvnc /usr/local/bin      and I tried to launch Xvnc as               niren@ubuntu-1204-precise-64-minimal:/$ vncserver      I got output is > New 'ubuntu-1204-precise-64-minimal:7 (niren)' desktop is > ubuntu-1204-precise-64-minimal:7 Starting applications specified in > /home/niren/.vnc/xstartup > > Log file is /home/niren/.vnc/ubuntu-1204-precise-64-minimal:7.log Instead of > New 'ubuntu-1204-precise-64-minimal:7 (niren)' desktop is > ubuntu-1204-precise-64-minimal:7 Starting applications specified in > /home/niren/.vnc/xstartup > > Log file is > /home/niren/.vnc/ubuntu-1204-precise-64-minimal.wolfram.com:7.log Since I run the `vncserver` under `niren@ubuntu-1204-precise-64-minimal`. Then I tried to test Xvnc as               mathematica -display :7      but it is not working and it is showing > mathematica: command not found. any suggestion. how can I run the X server to make _webMathemaica_ under Linux work? I am using               OS(Linux): Ubuntu-12.04.     Application server: Geronimo 3.0.1     webMathematica: 3.2     Mathematica(linux version): 9.0.1


title: Numeric calculation of Hessian

text: I want to calculate the Hessian matrix for a function that can only be evaluated numerically. So far, I have the following (where f is just for testing):               Needs["NumericalCalculus`"]          np = 3;  (* number of parameters *)     f0[{x_, y_, z_}] := x^2*y*z^3 + z + 1     (* numeric version; c to count evaluations *)     f[{pars__?NumericQ}] := (c += 1; f0[{pars}])          (* construct symmetric hessian *)     hessian[{pars__?NumericQ}] := Module[{diag, urt, temp1, temp2},        (* diagonal *)        diag = Table[          ND[f[ReplacePart[{pars}, i -> temp1]], {temp1,             2}, {pars}[[i]]], {i, 1, np}          ];        (* upper right triangle *)        urt = Table[          If[j > i,           ND[ND[f[ReplacePart[{pars}, {i -> temp1, j -> temp2}]],              temp1, {pars}[[i]]], temp2, {pars}[[j]]],           0           ],          {i, 1, np}, {j, 1, np}          ];        (* result *)        Table[         If[i == j, diag[[i]],          If[j > 1,           urt[[i, j]],           urt[[j, i]]           ]          ],         {i, 1, np}, {j, 1, np}         ]        ];      This seems to work:                c = 0;      hessian[{3, 4, 5}] // MatrixForm      yields the expected result. But the number of function evaluations appears to be rather large (I got c = 3099). Is this normal for ND, or can the above calculation be improved? Thank you for any answers/comments.


title: How can I set specific ImageSize to each cell of the Grid?

text: I want to make a `Grid` with each cell Imagesize was `{10,10}`.If I place any Gui-Element with more than that cell Imagesize,I don't want to scale that cell. I want to set each cell ImageSize is `{10,10}`. For this purpose,I tried with both `Grid` and `GraphicsGrid`. If you can use `Grpahicsrid`,it's working fine.That corresponding code is in below.               GraphicsGrid[ConstantArray["234", {3, 4}], Frame -> All,ImageSize -> {40, 30}]      If you use `Grid`,it's not working.That corresponding code is in below.                Grid[ConstantArray["234", {3, 4}], Frame -> All, ItemSize -> {10, 10}]      I think,`Itemsize` units totally different compared with `ImageSize`.That's why it's not working . How can I resolve this ?


title: Manipulate Plot and control variables using multiple controls

text: Here's a simple Manipulate example:               f[x_, a_, b_] := a x + b     Manipulate[      Plot[f[x, u[[1]], u[[2]]], {x, -1, 1},       PlotRange -> {-1, 1}], {{u, {1, 0}, "a/b"}, {-1, -1}, {1, 1}}]      ![enter image description here](http://i.stack.imgur.com/DoAo8.png) I would like to be able to control `a` and `b` with individual 1D sliders as well, such that all controls are dynamically updated. How can I achieve this?


title: Parameters for the CIE XYZ Color Space in Mathematica

text: In Mathematica, producing a color based on its RGB parameters is easy; one can use (for example)               Graphics[{RGBColor[0.8, 0.6, 0.5], Rectangle[]}, ImageSize -> 100]      or               Image[{{{0.8, 0.6, 0.5}}}, Automatic, ColorSpace -> "RGB", ImageSize -> 100]      This produces a square as follows: ![RGB Color Sample](http://i.stack.imgur.com/F3OJ9.png) The range of values for each of the `r, g, b` parameters are from `0 to 1`, as listed in the documentation. Any value of out of this range will be clipped. However, if I would like to produce a color based on its CIE XYZ color, I can still do as follows               Image[{{{0.8, 0.6, 0.5}}}, Automatic, ColorSpace -> "XYZ", ImageSize -> 100]      ![XYZ Color Sample](http://i.stack.imgur.com/6Smmo.png) Note that the color produced will be visually different since the `XYZ` color space is different from that of the `RGB` color space. My question is, do you know what are the range of values for each of the `X, Y, Z` parameters? I can't find this in the existing documentation. I was thinking of checking the range of values by trial and error (changing one parameter and seeing whether the values are clipped by comparing the results visually), but I hope there's a better way to do this or that there was some documentation that I missed. **Note** : The online EasyRGB calculator limits the `X, Y, Z` parameters from `0 to 95.047, 0 to 100, 0 to 108.883` respectively, and I'm suspecting that the limits for Mathematica would be simply 1/100 of these limits. Would like some verification on this if possible, though.


title: Animated gif form parametric plot

text: I need to create an animated gif for the following function.               rx1=((-T \[Omega] Cos[T \[Omega]] + Sin[T \[Omega]]))/(2\[Omega]^2)     ry1=(-2 + 2 Cos[T \[Omega]] + T \[Omega] Sin[T \[Omega]])/(2\[Omega]^2)     pulsefrequency=200000*\[pi]; \[tau]pulse1= 20*10^-6;      When I use                  Animate[       Show[ParametricPlot[{rx1, ry1}/. \[Omega] -> pulsefrequency,        {T,0,u},PlotStyle -> {Thick, Red}],        Graphics[        {Blue, PointSize[0.02],          Point[{rx1, ry1}  /. \[Omega] ->                pulsefrequency /. T -> u]}]],         {u, 0, \[Tau]pulse1}, AnimationRepetitions -> 1]      It works well but when for getting an animated gif, I replace Animate command with Table and then it ceases to work               Table[ Show[     ParametricPlot[{rx1, ry1}  /. \[Omega] -> pulsefrequency, {T, 0, u /. u -> s}, AspectRatio -> 1, PlotRange -> {{-0.010, 0.010}, {-0.010, 0.010}}, Background -> RGBColor[100.0, 100.0, 100.0]],Graphics[{PointSize[.025], Hue[0], Point[{rx1, ry1}  /. \[Omega] -> pulsefrequency  /. T -> s]}]] , {s, 0, \[Tau]pulse1, \[Tau]pulse1/10}]      The output that I recieve is               {Show[ParametricPlot[{rx1, ry1} /. E0 -> pulseamp /. \[Omega] ->         pulsefrequency /. m -> mass /. q -> charge, {T, 0,      u /. u -> s}, AspectRatio -> 1, PlotRange -> {{-0.01, 0.01}, {-0.01, 0.01}}, Background -> RGBColor[100., 100., 100.]], \!\(\*GraphicsBox[{Hue[0], PointSize[0.025], PointBox[{0., 0.}]}]\)], **followed by the list of the images**       I have spend hours to figure out how to just get the animated image rather than the command. I am not sure whats going wrong. Any help would be appreciated. Regards


title: How can I draw the transition diagram of a Markov chain?

text: I apologize for my heaviness as I am a _Mathematica_ novice. Is there a kind soul who can tell me how to draw diagrams like those shown on page 158, 160 and 161 of this document.Chaine de Markov I oversee several plotting examples unsuccessfully.


title: Exporting German Umlauts to PDF and importing them doesn't work - character encoding problem

text: I am using _Mathematica_ 9.0.0.0 under German Windows 32bit. I am on trying to export and import German text to PDF and re-import it. It seems that there is a confusion about the character encoding. It obviously differs between `Export` and `Import`. However, nothing related to the encoding is described in the PDF documentation, nor does specifying `CharacterEncoding` in PDF import and export seem to make any difference.               ImportString[ExportString[       "This is a string. Aber was ist mit äöü, ß und ÄÖÜ? Also x\.b2 und \y\.b3 measured in µm",       "PDF"]]      I wonder if there is a work-around which doesn't require a complete re- mapping. PDF export is slow, and doing it character by character will probably take hours.


title: Extract connected components from PNG image file and write to individual PNG files

text: Given the PNG image below with the three leaves in it, I'd like to extract each individual leaf and write it to its own PNG file, using _Mathematica_ 9. Note that in the original image, the box bounding a component may possibly intersect with another component, so that needs to be taken care of (i.e. each output PNG should have one complete leaf and no part of another leaf). The segmentation itself looks quite easy - just use the alpha channel, I'd say. Hopefully this is a simple problem for an experienced _Mathematica_ user (which I am not) - _Mathematica_ just happens to be within reach and I figured it might be a good tool to use. And if you could have the code go through a directory of images similar to this one, perform the desired operation on each PNG image in it and save the output images to another directory, that would be great! (the components of image.png should be saved as image1.png, image2.png, etc.) ![](http://i.stack.imgur.com/8enYZ.png) **Edit** : While it would obviously be unrealistic to expect a segmentation algorithm to work universally based on a single sample image that I showed, unfortunately some of the ones that have been suggested are failing for some cases they clearly shouldn't. Two examples below: ![](http://i.stack.imgur.com/OZEvk.png)![](http://i.stack.imgur.com/Tl1Pm.png) I'm hoping the authors will be able to modify their answers so that they work with these images (and hopefully most images "like" these ones), particularly if the problem is with the cropping/trimming/masking step (as opposed to the image preprocessing step, which might need to be "tuned" until it works for most examples). I should mention that I am interested in the special case of PNG images where the **AlphaChannel[]** might be used to advantage, but I'm happy that the solutions offered so far are more generally applicable. I do need to be able to deal with cases in which the components are quite close together (but not touching).


title: Save plot by pressing on button

text: I want to save a certain graphics by pressing a button. Thinks work fine for something like this               Button["Save",      export = Histogram[RandomVariate[NormalDistribution[0, 1], 10^3]];      Export["test.pdf", export]      ]      By increasing the number of data points the button press doesn`t generate test.pdf file anymore               Button["Save",      export = Histogram[RandomVariate[NormalDistribution[0, 1], 10^7]];      Export["test.pdf", export]      ]      Simply executing               Export["test.pdf", export]      after pressing the button, generats the pdf file. What is the cause of this behaviour? Thanks in advance for your help


title: Stochastic problem

text: I have to organize a small sports league and I am puzzled on how to create the game plan. We are 8 persons playing table soccer with 2 vs. 2 matches. The idea is that each person plays once with every other person. So e.g. person 1 will play once with person 2, 3, 4, 5, 6, 7 and 8. So there are 7 matchdays and 2 games on each. If a team wins both players receive a point. At the end of the season a player can have earned a max of 7 points. _E.g. matchday 1: p1 &p2; p3&p4; p5&p6; p7&p8 E.g. matchday 2: p1&p3; p2&p4; p5&p7; p6&p8_ How can I calculate the teams for each of the seven game days given the conditions above? Since the strongest two teams would play against each other depending on points collected in the earlier games I only need a set of 4 teams for each matchday. I would be grateful for some smart help on calculating a game plan. I guess I am stuck here: players = Range[1, 8]; teams = Select[Subsets[players , 2], Length[#] == 2 &]; Thanks a lot Patrick


title: OpenCLLink specify grid size when calling function

text: The Mathematica documentation (section "Details": https://reference.wolfram.com/mathematica/OpenCLLink/ref/OpenCLFunction.html) states that: > On launch, if the number of threads is not specified (as an extra argument > to the ), then the dimension of the element with largest rank and dimension > is chosen. For images, the rank is set to 2. Great - but how do I do that when I call an OpenCL Function? There is so documentation for that.. Thanks


title: Plot and Histogram together in Mathematica 7

text: I want to draw a histogram and some line plots together. I have one observed variable and several predicted ones. Observed data should be shown as a histogram and the predicted variables should be shown as line plots or smooth histograms. I use _Mathematica_ 7. Please help me to show which fitted distribution is best.


title: Using BarChart and pastel colors is there a method to get the labels to be the same color as the bar?

text: I am using version 8 and I have the following code:               totals4 = {{1093, "AUSTRALIA"}, {235, "BELGIUM"}, {474, "BERMUDA"}, {187, "BRAZIL"}};     BarChart[#, ChartStyle -> "Pastel",               ChartLabels -> (Placed[Rotate[Style[#, Blue], Pi/3] & /@ #2, {{.5, 0}, {1, 1}}])              ] & @@ Transpose[totals4]      Is there a way to get the labels to be the same color as their respective bar? Thank you in advance. ![bar chart](http://i.stack.imgur.com/3X0OS.jpg)


title: How specify PlotRange and display bitmap over plot in DateListPlot?

text: The problem came up in a plot of RGB data obtained by imaging the variable star d-Cephei with a DSLR camera on several successive nights. I tried to plot the R, G, B intensity points for each day in a column and place a bitmap of the star (taken that day) over the column. This worked fine with `ListPlot` during September, but with October coming up, it was clear I would need something like `DateListPlot`. But when I substituted `DateListPlot`, after rewriting the data structures appropriately, I found that the points would plot, but the bitmaps that I had inserted with `Epilog` disappeared. Furthermore, by adding a `PlotRange` statement, I could get the bitmaps, but the R, G, B intensity points would not display. In what follows, I have tried to boil it down using just one suitable list of points for one color and just one bitmap created by rasterizing a Disk. Somewhat to my surprise, the problem is reproducible. As written, the image plotted with `ListPlot` shows it working properly, but it lacks the formatting of the date information on the x axis. The image plotted with `DateListPlot` formats the dates nicely, but it doesn't show the bitmap. My questions are:   1. What will make the bitmap display with `DateListPlot`?   2. How can I move the bitmap display above the range of the points on the Y axis, and below the title?   3. Is it possible to specify a `PlotRange` in DateListPlot? How would one specify the range for the x axis?   4. Why does this happen?   5. When will a book be published listing all the do's and don't's of plots and graphics in Mathematica? My code follows:               Clear[disk, datePts, pts, dateListPlot, listPlot]          disk = Rasterize[Graphics[Disk[]]];          datePts = {{{2012, 9, 24}, .1}, {{2012, 9, 25}, .35}, {{2012, 9,           26}, .75}, {{2012, 9, 27}, .90}, {{2012, 9, 28}, .60}, {{2012, 9,           29}, .30}, {{2012, 9, 30}, .10}};          pts = {#[[1, 3]], #[[2]]} & /@ datePts;          ListPlot[       pts,       PlotStyle -> Directive[{Red, PointSize[Large]}],       PlotLabel -> "ListPlot: Epilog disk visible",       Epilog -> Inset[disk, {23.75, .8}, {0, 0}, .4],       PlotRange -> {{23.5, 30.5}, {0, 1}},       Background -> LightGray       ]          DateListPlot[      datePts,      PlotStyle -> Directive[{Red, PointSize[Large]}],      PlotLabel -> "DateListPlot: Epilog disk not visible",      Epilog -> Inset[disk, {23.75, .8}, {0, 0}, .4],      (* Activate PlotRange makes Epilog visible, but points will not be visible. *)      (*PlotRange->{{23,30},{0,1}},*)      Background -> LightGray      ]


title: How do I operate on this list elegantly?

text: calibData = { {200, 10}, {100, 20} }          MapAt[N, MapAt[Log10, calibData, {;; , 1}], {;; , 1}]      I want to apply a base 10 Log to the first element of each tuple. Then the output should be displayed as a floating point number. This bit of code does the job, but it is ugly and its hard to understand what it is doing. Is there a better way of completing this task?


title: Better solution than returning a list of 3 values?

text: I have a function (using `SetDelayed`) that currently returns 3 values in a list. Later on I use the result of this list along with `[[1]]`, `[[2]]`, and `[[3]]` to use the values. Is there a way I can give each value a "name" of some sort, and return only one value in such a way that all these values can be accessed by name? (Coming from an object-oriented programming perspective, I just want to return a single object with a few fields/accessors.)


title: How to measure cross sectional area of Tube

text: I have modified my question to a simple version. I have vertex clusters, polygons and a list of center points. Attached pts and polygons of 3d tube are:                pts =Import["https://dl.dropboxusercontent.com/u/68983831/npts.txt","Table"];       polys= Import["https://dl.dropboxusercontent.com/u/68983831/npolys.txt","Table"];      These are approximate center points:                cpts = Import["https://dl.dropboxusercontent.com/u/68983831/cpts.txt","Table"];           Show[Graphics3D[{Opacity[.1], EdgeForm[None], GraphicsComplex[pts, Polygon@polys]},       Boxed -> False], ListPointPlot3D[cpts]]      ![enter image description here](http://i.stack.imgur.com/4cuX0.png) How to measure cross-sectional area, major & minor axis length as move along the center points? or may be without center points?


title: Commutative algebra by generators and quadratic relations

text: My apologies if this has been asked before. Please also be gentle, as this is my first question here. I am trying to define a commutative algebra by generators and relations in the following way. For a simple example, we can consider the Z[b] algebra generated by a, modulo a^2=b. I would like to introduce a rule (or a list of rules in more complicated examples):               In[1]: TheAlgebra = a*a -> b      It would be cool to evaluate products like this:               In[2]: a^3 //. TheAlgebra      However, this input simply returns the symbol `a^3`. I would have expected the result `a*b`. The following input works though:               In[2]': (a*a//.TheAlgebra)*a //. TheAlgebra      My question is this: how can I get Mathematica to simplify an expression such as a^kb^l in a smart way? There are similar questions on the site, but maybe the answer for this one is simpler, since I'm only interested in simplifying products of the various symbols together with some integers. Maybe there is a way to force Mathematica to write `a^3` in a full-form expression that I like, such as: `Times[a,Times[a,a]]`?


title: Combining heads of lists so that you can create a nested list from two sublists

text: I have a list such as `list = {1,4,3,4,6,3,1,2,3,4 ... }` and a range, say, `range = Range[20000]`. I want to obtain the following output:               {{0,1}, {1,4}, {2,3}, {3,4}, {4,6}, ... {nth element of list, nth element of range}}      where `n` is basically `Range[Length[list]]`. The general question here would be how can I create a list from the $a_i$th index of list `a`, the $b_k$th index of list `b` and so on. I would have tried `MapThread` but I'm trying to use the `Histogram3D` function so I needed this "nice" list.


title: Manipulating infinite series

text: If I type               Sum[f[x],{x,m,Infinity}]-Sum[f[x],{x,m+3,Infinity}]      I would like Mathematica to return something like               Sum[f[x],{x,m,m+2}]      or maybe               f[m]+f[m+1]+f[m+2]      But it doesn't seem to want to do this (maybe because it's worried that the difference makes no sense if the sums don't converge?). Is there a simple way to let Mathematica know I would like it to make this kind of simplification?


title: Masked Image Registration/Alignment/Correlation

text: I'm trying to register (align) two circular grayscale images using a boolean mask to prevent them simply being aligned on the fixed circular border. There are also small fixed circular regions in the image that cause further problems (image artifacts, black dot marks). In the Mathematica documentation (https://reference.wolfram.com/mathematica/ref/ImageCorrelate.html) there is a simple example for image correlation with a circular mask. Here is an example image: ![example](http://i.stack.imgur.com/0c5E2.jpg) Here is the code I put together using `ImageCorrelate`, I want to use this, except for larger images then are in the example:               (*note this mask isn't complete*)     mask = ColorConvert[        Image@Rasterize[          Graphics[{{White, Disk[{0, 0}, 700]}, {Black,              Disk[{70, -20}, 60]}},            PlotRange -> {{-800, 800}, {-800, 800}}, Background -> Black],           ImageSize -> 500, RasterSize -> 500], "GrayScale"];          pos = Position[Flatten[ImageData@mask], 1.];     f = NormalizedSquaredEuclideanDistance[Extract[#1, pos], Extract[#2, pos]]&;          image1 = Import["https://dl.dropboxusercontent.com/u/3730003/map4__00001.jpg"];     image2 = Import["https://dl.dropboxusercontent.com/u/3730003/map5__00011.jpg"];          (*this takes forever...*)     corr = ImageCorrelate[image1, image2, f];     minOffset = PixelValuePositions[corr, Min[ImageData[corr]]][[1]] - ImageDimensions[image1]/2;      How can I make this more efficient such that it can be used on larger images? I've tried to make my own compiled version of `NormalizedSquaredEuclideanDistance`, however I haven't been able to figure out how to deal with arbitrary length vectors... Some code borrowed from: Image correlation To follow up on @bills comment, I did try `ImageCorrespondingPoints` and `ImageAlign`, however, they often find matching points on the static boundary as corresponding points which would prevent accurate alignment. I couldn't figure out how to implement some sort of boolean mask. Cropping a square (code below) from the center of the image helped the situation but the black spots still cause problems.               cropper[image_] := ImageAdjust[ImageCrop[image, 1000]];          alignerKey[image1_, image2_] := Module[{cropped1, cropped2, transform},       transform = FindGeometricTransform[cropper@image1, cropper@image2, "Transformation" -> "Translation"];       ImagePerspectiveTransformation[image2, Last@transform, DataRange -> Full]       ]


title: Numerically solving PDE with high precision

text: I want to numerically solve the PDE $\partial_t u(t,x)=c\partial_x u(t,x)+(mx-l)u(t,x)$ with some initial and boundary conditions and given parameters $c$, $m$ and $l$. Consider the code               T=55;     xMax=200;     l=.6;     m=0.0072;     c=2;          system={       D[u[t,x],t]==c*D[u[t,x],x]+(m*x-l)*u[t,x],       u[0,x]==500*Evaluate@PDF[TruncatedDistribution[{0,xMax-1},NormalDistribution[100,10]],x],       u[t,0]==u[t,xMax]==0     };          (*prepare animation of progress update*)     "Computation started: "<>DateString[]     currentTimeStep=0;     "Time step: "ProgressIndicator[Dynamic[currentTimeStep],{0,T}]          sol=First@NDSolve[system,u[t,x],{t,0,T},{x,0,xMax},StepMonitor:>(currentTimeStep=t;),Method->{"MethodOfLines","TemporalVariable"->t,"SpatialDiscretization"->{"TensorProductGrid","MinPoints"->10000}}];          "Computation stopped: "<>DateString[]          Manipulate[         Plot[Evaluate[u[t,x]/.sol/.{t->time}],{x,0,xMax},PlotRange->{{0,xMax},{-50,50}},AxesOrigin->{0,0}]     ,{time,0,T}]      As you can see, for $T=55$, the solution shows numerical artifacts. ![Numerical solution of u\(55,x\)](http://i.stack.imgur.com/xlKun.png) It seems that the cause for this is a not sufficent accuracy of the computation, since               u[t,x]/.sol/.{t->10,x->#}&/@Range[xMax-5,xMax]      yields               {2.53114*10^-18, -1.23827*10^-17, -7.49676*10^-15, -2.81533*10^-15, 2.97761*10^-14, 0.}      but it follows from the PDE that the solution stays non-negative for non- negative initial data. Using a higher `AccuracyGoal`, such as `AccuracyGoal->6` only results in huge memory consumption, even with a broader grid-spacing. Any ideas how to eliminate the produced artifacts?


title: 3D spectrogram plots

text: > **Possible Duplicate:**   >  Computing and plotting a spectrogram in Mathematica I have a notebook that reads in some time-series data and calculates a spectrogram of it using Mathematica 9.0's spectrogram function. This makes a 2D color contour plot of frequency vs. time (colors give the amplitude function). I would like to extend this to a 3D spectrogram plot [for example, see: http://en.wikipedia.org/wiki/File:Spectrogram.png]. The necessary data is already in the existing spectrogram function, but how do I pass this to something like Plot3D to get a 3D plot of amplitude vs. frequency and time?


title: Converting expression to C

text: I have an expression of the form                   -0.0000289725287527177708 - 2.52403420408155732 x +       138.677105376831122 x^2 - 3402.37981527828424 x^3 +       34440.8443628217428 x^4 + 158064.877964911022 x^5 -       8.04498826077845134*10^6 x^6      that I would like to convert over to an expression that can be used in C++. In C++ I have defined some functions for taking powers and the syntax is that $x^3=third()$, $x^4=fourth()$ and so on. Besides this, I to convert $*10^6$ to $e6$. Is there an easy way to do these things?


title: Getting Mathematica.SE reputation data

text: I'm asking the following on behalf of Mr.Wizard. On StackOverflow I used the old v1.0 SO API and _Mathematica_ to obtain user reputation data. Now we are no longer a tag on StackOverflow and the v2.1 API is the current top dog, can we still use similar techniques to get the reputation data here?


title: NDSolve Issue with initial conditions

text: I am trying to implement a Poincare section for a gravitational movement on the plane $(x,y)$. Here is the code I wrote               PoincareSection[cix_, civx_, ciy_, civy_]:=         q = 0.9;         Rc = 0.14;         v0 = 1;         Φ[x_, y_] := 1/2 v0^2 Log[Rc^2 + x^2 + (y/q)^2];         h1 = x'[t] == vx[t];         h2 = vx'[t] == -D[Φ[x, y], x] /. x -> x[t] /. y -> y[t];         h3 = y'[t] == vy[t];         h4 = vy'[t] == - D[Φ[x, y], y] /. x -> x[t] /. y -> y[t];         Return[Reap[NDSolve[{h1, h2, h3, h4,         x[0] == cix, vx[0] == civx, y[0] == ciy, vy[0] == civy,         WhenEvent[{y[t] == 0, vy[t] > 0}, Sow[{x[t], vx[t]}]]},         {x[t], vx[t], y[t], vy[t]},         {t, 0, 100}]][[2]]];      However, when I try to run this code, I always get the following error: > NDSolve::ndinnt: Initial condition civx is not a number or a rectangular > array of numbers. I have tried using `N` and pattern `_?NumericQ` for the parameters of the function, but that didn't work out... Any ideas what may be causing the problem ? (If you see any other improvements I could make to my code I would be happy to hear them...)


title: Built-ins for discrete ArgMax or ArgMin?

text: Wondering whether there are discrete versions of `ArgMax` / `ArgMin`, that is, something that will find the part of an expression at which some function `f` of the elements is maximized (or minimized). `Last[Ordering[f/@list]]` and `First[Ordering[f/@list]]` will do, but the following does not work:               ArgMax[Function[i, (f/@list)[[i]]], i]


title: Mathematica does not respect tensor order?

text: I'm trying to verify an identity involving symmetric traceless tensors over the reals. What I tried was:               In[11]:= $Assumptions = d \[Element] Matrices[{3, 3}, Reals, Symmetric[{1, 2}]]     Out[11]= d \[Element] Matrices[{3, 3}, Reals, Symmetric[{1, 2}]]     In[24]:= dd = d - Tr[d] IdentityMatrix[3]/3     Out[24]= {{d - Tr[d]/3, d, d}, {d, d - Tr[d]/3, d}, {d, d,d - Tr[d]/3}}      I think the output in 24 is assuming that d is now a scalar. That gets me into trouble later when I tried to compute outer products of dd and use TensorReduce. For example               lh = TensorReduce[TensorContract[dd\[TensorProduct]dd, {{1, 3}, {2, 4}}]]     TensorDimensions::fscl: Nonscalar expression d encountered as an argument of numeric function Power. >>     TensorDimensions::inhom: Inhomogeneous dimensions in sum d-Tr[d]/3. >>      Why is this not working, and what is the right way to do computations on traceless matrices? Sune


title: FindMinimum slower in MMA10

text: I've recently updated one of my PCs (both have the exact same build) to Mathematica 10. After running the exact same notebook (which has two `FindMinimum` calls) on both, the one with Mathematica 10 takes 831 and 19 065 seconds to converge, whereas the old version only takes 27 and 1510 seconds. This is quite a difference! The question is : has anything changed in the implementation of `FindMinimum`, or in the default settings, that would cause this? I checked `WorkingPrecision` and it hasn't changed in MMA10. Any other clues? As a minimal example I checked with the code in this question, and got the same issue.


title: GeoDirection and GeoDistance Memory Leaks: How to Recover the Memory?

text: The functions GeoDirection and GeoDistance have memory leaks. To see this, run the example functions below and watch as the Physical Memory Available in the Windows Task Manager, Performance tab, monotonically decreases with function execution.               memoryLeakGeoDirection[x_] :=      Block[{},      Do[GeoDirection[{RandomReal[{5, 15}],       RandomReal[{105, 135}]}, {RandomReal[{20, 25}],       RandomReal[{140, 155}]}], {x}]]          memoryLeakGeoDirection[10000]      and               memoryLeakGeoDistance[x_] :=      Block[{},      Do[GeoDistance[{RandomReal[{5, 15}],       RandomReal[{105, 135}]}, {RandomReal[{20, 25}],       RandomReal[{140, 155}]}], {x}]]          memoryLeakGeoDistance[10000]      If one swaps GeoDirection or GeoDistance with any of GeoDestination, GeoPosition or Log, Sin, etc. then no memory loss occurs. I have coded a large Monte Carlo simulation using GeoDirection and GeoDistance. To run the actual simulation it is necessary to call these functions on the order of 108 to 109 times. I have upgraded the RAM from 16 GB to 32 GB, and upgraded from Windows 7 Home Premium to Windows 7 Professional, all in a failed attempt to stay ahead of the leaking memory. I am using Mathematica 8.0.4.0 64 bit. I have tried `$HistoryLength = 0` and `ClearSystemCache[]`, without success in recovering the leaking memory. Do you have any suggestions as to how to recover this lost memory during a running simulation?


title: Find Roots in Do loop

text: Task: Finding roots in loop               t = List[1, 2, 3, 4, 5]     fx[x_] := a*x^2 - 5     List[Do[Print[FindRoot[fx[k] == 1, {a, 1}]], {k, 0, 5}]]      Output: Currently the output is as follows:                {a-> 6}      {a-> 1.5}      and so forth But I want to have the values of a in a List in order to perfom further tasks (Plotting root of f_a against a)


title: Named string patterns in Alternatives, why does it depend on order?

text: When using the same name for part of patterns inside `Alternatives` how come it behaves differently for strings than in normal patterns?               {StringMatchQ["ab", (x_ ~~ "c") | ("a" ~~ x_)],      StringMatchQ["ab", ("a" ~~ x_) | (x_ ~~ "c")]}     (* False, True *)          {MatchQ[{1, 2}, {x_, 3} | {1, x_}],      MatchQ[{1, 2}, {1, x_} | {x_, 3}]}     (* True, True *)          StringMatchQ["ab", (x_ ~~ "c") | (x_ ~~ "b")]     (* False *)


title: NSolve for high degree univariate polynomials

text: I am trying to solve a high degree univariate polynomial using Mathematica's NSolve command. But when I plug the solutions generated by Mathematica back to the equations, the equations are massive numbers and not even close to zero.               Clear["Global`*"];     r := Rationalize[RandomReal[NormalDistribution[0, 1]]]     randompol[n_, i_] := Sum[r Sqrt[Binomial[n, j]] x^j, {j, 0, n}]          e= randompol[100,1];          SetSystemOptions["NSolveOptions" -> {"Tolerance" -> 10^(-10)}];     sol = Solve[e == 0, x];          e/.sol      which gives some of the entries as big as `-1.8513*10^56` I don't know what other method can be used to solve high degree polynomials. Any help would be much appreciated. Regards, -John


title: Why does MapIndexed return the index in a list?

text: `MapIndex` returns the index in a list instead of as shown in my revised example. I can work with it (or work around it), but what's the design concept? It seems like an unnecessary complication, but I'm sure there is a reason I don't see. ![enter image description here](http://i.stack.imgur.com/7hqey.png) I tagged this with a design-pattern tag. I'm not sure if that's appropriate. Editors please remove the tag if it's not.


title: Why do Mathematica and MathKernel need "incoming connections"?

text: My firewall (the "Application Firewall" in OS X 10.8.2) detects that _Mathematica_ (9.0.1.0) and "MathKernel" both request incoming connections. My configuration expects all incoming connections to be blocked, so I wonder what the consequences of denying these connections is. ![enter image description here](http://i.stack.imgur.com/0FBzT.png) Why do _Mathematica_ and "MathKernel" request incoming connections? What are these connections used for and what, specifically, will be the consequences of denying them? When I allow connections, and take a quick look at them, I see that both have opened several incoming/outgoing paris of connections to "localhost", which suggests that the connections are being used for local communication between each other, as suggested in Szabolcs's answer ![enter image description here](http://i.stack.imgur.com/7TePT.png) ![enter image description here](http://i.stack.imgur.com/BY64W.png) but closer examination shows an additional connection for each that I'm not sure how to interpret ![enter image description here](http://i.stack.imgur.com/ikoyn.png) Blocking incoming connections seems to have no effect, and I wonder of that's really true (e.g., am I paying a performance price). Allowing would seem to be innocent, of all the connections are local; but I know of no way to guarantee that.


title: Is it possible to specify which range that slot takes from?

text: consider this example:               list = {a, b, c};     Table[{Sin[i], Cos[j]}, {i, 1, 3}, {j, list[[i ;;]]}]          (*{{{Sin[1], Cos[a]}, {Sin[1], Cos[b]}, {Sin[1], Cos[c]}}, {{Sin[2],         Cos[b]}, {Sin[2], Cos[c]}}, {{Sin[3], Cos[c]}}}*)      for the following code:               {Sin[#], Cos[#]} & /@ list[[# ;;]] & /@ {1, 2, 3}      we can get same result if the slot inside `Sin` `(Sin[#])` takes the values from the outer range `({1,2,3})` so in general, is it possible to specify which range that slot takes from?


title: How do I use the result of Solve in a function definition?

text: I'm trying to adapt the following piece to the part where `PowerF` is a function of `amax`.               Z = Max[Z /. Solve[FullSimplify[PowerF, Z ∈ Reals] == P0, Z]]      Now I have modified my code because I want to also make a plot that depends on the parameter `amax`. However the following does not work:               Z[amax_] := Max[Z[amax] /. Solve[FullSimplify[PowerF[amax], Z ∈ Reals] == P0, Z]]      How to fix this?


title: Replacement rule to create graphs

text: {{{x11, 31}, {1, {{x10, 19}}}, {0, {{x9, 12}}}}} /. {___,         s : {{a_, b_}, {c_, {{d_, e_}}}, {f_, {{g_, h_}}}}, ___} ->        Graph[{a \[UndirectedEdge] d, a \[UndirectedEdge] g},         VertexLabels -> {a -> a, d -> d, g -> g}, EdgeWeight -> {c, f},         ImagePadding -> 10, EdgeLabels -> {{a \[UndirectedEdge] d ->             c}, {a \[UndirectedEdge] g -> f}}]      It generates a graph but it doesn't replace the variables as is expected,variables stay same without replacement. I am also facing problem with labeling vertices using patterns. Does someone know how to achieve this replacement ?


title: Antialiased text in manipulate constantly re-evaluated, ignores TrackedSymbols

text: I am trying to display antialiased text (using the ImportString[ExportString[]] method shown here) but inside a manipulate so that the text can be changed interactively. I have found a way to do this, but problem is that the function is constantly being updated, even when the tracked symbol "i" is unchanging. A similar problem is outlined here. Below are simplified versions of my code to:               Manipulate[     Style[i,86,Bold,GrayLevel[0.4],FontFamily->"Verdana"],     {{i,36890},0,99999,1},TrackedSymbols:>{i}     ]          toCurves[a_]:=ImportString[ExportString[a,"PDF"],"PDF"][[1,1]];     toGraphics[a_,b_]:=Graphics[{a},Method->{"ShrinkWrap"->True},ImageSize->b];          Manipulate[     toGraphics[toCurves[Style[i,Bold,GrayLevel[0.4],FontFamily->"Verdana"]],{{Automatic},{64}}],     {{i,36890},0,99999,1},TrackedSymbols:>{i}     ]      The first manipulate displays the text with jaggies and works as expected. The second manipulate uses the custom functions toCurves and toGraphics and displays the text without jaggies but is constantly being updated. The cause of the problem appears to be the ExportString function, but it could also be the ImportString function. Edit 1: I found a way to stop the manipulate from constantly updating using the If function:               style[a_] := Style[a, Bold, GrayLevel[0.4], FontFamily -> "Verdana"]     toCurves[a_] := ImportString[ExportString[a, "PDF"], "PDF"][[1, 1]];     toGraphics[a_] :=        Graphics[{a}, Method -> {"ShrinkWrap" -> True},         ImageSize -> {{Automatic}, {64}}];          Manipulate[      If[       i != j,       j = i; k = toGraphics[toCurves[style[i]]],       k       ],      {{i, 36890}, 0, 99999, 1},       Initialization :> {j = ""},       TrackedSymbols :> {i}      ]      I modified the toGraphics function and added a custom style function to reduce some of the clutter in the code. It isn't as elegant as I'd like it to be, but it works.


title: Accumulate second elements in sublist

text: I want to accumulate the second elements in each sublist identified by its first element.               {       {{1, 1}}, {{2, 1}}, {{3, 1}}, {{2, 2}}, {{5, 1}},        {{2, 1}, {3, 1}},       {{7, 1}}, {{2, 3}}, {{3, 2}},       {{2, 1}, {5, 1}},       {{11, 1}},       {{2, 2}, {3, 1}},       {{13, 1}},       {{2, 1}, {7, 1}},       {{3, 1}, {5, 1}},       {{2, 4}}, {{17, 1}},       {{2, 1}, {3, 2}}      }      For the first eight elements in my list               {{{1, 1}}, {{2, 1}}, {{3, 1}}, {{2, 2}}, {{5, 1}}, {{2, 1}, {3, 1}}, {{7, 1}}, {{2, 3}}}, ...      the desired result is               {{{1, 1}}, {{2, 1}}, {{3, 1}}, {{2, 3}}, {{5, 1}}, {{2, 4}, {3, 2}}, {{7, 1}}, {2, 7}}}, ...      In the result, the first element of each sublist will be the identification element, and the second element is just how much of the first element there are upto any given postition of the list So for the sixth position in the list, `{{2, 1}, {3, 1}}`, the result would be `{{2, 4}, {3, 2}}`, because by that time in the list there will be an accumilation of 4 'twos' and 2 'threes'.


title: Manipulate that remembers its own input code and can show it by request

text: This question is a further question to my previous one. The solution for that question works well as long as the original input `Manipulate` code is still provided in the notebook. But in general, there are cases that only the output cell is included and the input code omitted. In such case, I think we could possily make the `Manipulate` remembers its own input code, and by request (e.g. a button) we can show the code. A first attempt starts with either `DynamicModule` or `Pattern`. For instance, I tried               DynamicModule[{mc = HoldComplete@Manipulate[      Column[{        Plot[Sin[x - n], {x, 0, 12.4}],        Button["show code",         NotebookPut[Notebook[{Cell[BoxData[ToBoxes[mc]], "Intput"]}],           EvaluationNotebook[]]         ]}], {n, .3, 11.8}]},      ReleaseHold[mc]]      where mc is the short for 'manipulatecode'. But the code will not be shown if we click the button, only the output instead. I also tried               mem : Manipulate[     Column[{     Plot[Sin[x - n], {x, 0, 12.4}],     Button["show code",      NotebookPut[Notebook[{Cell[BoxData[ToBoxes[mc]], "Intput"]}],        EvaluationNotebook[]]      ]}],     {n, .3, 11.8},     Initialization :> (mc := HoldComplete[mem];)]      But this will not work also. Is there a way to let the button creates a new cell and show the `Manipulate` code in this cell (better in expression format)?


title: How to update Print-out "in place"?

text: I want to run an exhaustive search algorithm that will take several hours to finish, so I'd like to have it execute `Print` statement periodically, showing the fraction of the search space covered and the best result found so far. E.g.               For[i = 1, i <= n, i++,       ...       If[Mod[i, stride] == 0, Print[{N[i/n], bestsofar}]     ];      This works OK if the number of updates is not too big. Otherwise, the voluminous print out begins to overwhelm the notebook, making it increasingly difficult to navigate through it. Closing the cell does not help here: it is opened automatically with each new print-out. (Also, when the cell remains open, the window will not necessarily auto-scroll to keep the latest print-out in view, so one needs to manually scroll to see it. No big deal; just a bit annoying.) I'm looking for a way to have each new print-out overwrite the previous one "in place". With other programming environments, at least in Unix-like systems, one achieves this effect by prepending a carriage return character (ASCII 13, aka `\r`) to the printed string, and omitting the newline character (ASCII 10, aka `\n`) from the end of the string1. A single newline character may be eventually printed, before execution terminates. The `Print` command implicitly adds the newline at the end of the output, and I have not found a way to suppress this behavior. Is there some other way to achieve the effect of updating a running cell's print-out in-place? UPDATE: based on molekyla777's answer, I modified the code above to:               Monitor[For[i = 1, i <= n, i++,       ...       Refresh[{N[i/n], bestsofar}, UpdateInterval -> 1, TrackedSymbols -> {}]     ];      I added the `Refresh` call to reduce the performance drain that would otherwise result from updating the `N[i/n]` term, which changes with each iteration. It basically does what `Mod[i, stride] == 0` did in the previous version. * * * 1If the length of the output string is not monotonically increasing, it may also be necessary to right-pad the output string with enough spaces to "erase" the previous output.)


title: Pass model parameters to NDSolve'Reinitialize

text: I am wondering if there is a way to pass parameter values during the NDSolve'Reinitialize step of the NDSolve process? I've read that one way to possibly speed up repeated calls to NDSolve is to Process the equations once, and then use NDSolve'Reinitialize to update the initial conditions. In the ODE example below, the initial conditions will never change, instead I am searching for the rate parameters ka and kb that minimize a scoring function. How could I use NDSolve'Reinitialize to pass different values of those parameters, and more efficiently search for the optimal parameter sets? k               a = .1;     kb = .01;          decay = First[          NDSolve`ProcessEquations[        {          A'[t] == -ka*A[t],         B'[t] == ka*A[t] - kb*B[t],         c'[t] == kb*B[t],         A[0] == 35.,         B[0] == 0.,         c[0] == 0.         },        {A, B, c},        t,        Method -> "Automatic"]]     NDSolve`Iterate[decay, 200]          sol = NDSolve`ProcessSolutions[decay]          Plot[      Evaluate[       {A[t], B[t], c[t]} /. sol,       {t, 0, 200}       ],      PlotStyle -> {Red, Blue, Green}]


title: Matrix constructed as a function inside a for-loop

text: I am trying to construct a matrix `A` which depends on variables `x` and `y` through iteration. My code is:               Nmax = 2;     For[iw = 1, iw <= 5, iw++,               A[x_, y_] :=        Normal[SparseArray[{{i_, j_} /; j == i :> iw + 2 (Cos[x] + Cos[y]),          {i_, j_} /; j - i == 1 :> 1, {i_, j_} /; i - j == 1 :> 1}, {Nmax,           Nmax}] ]              ;       Print["Testing A-matrix at (x,y)=(\[Pi]/2,\[Pi]/2)"];            MatrixForm[A[Pi/2, Pi/2]] // N      ]             MatrixForm[A[Pi/2, Pi/2]] // N      You can see that I am trying to construct a matrix where the diagonal elements are decided by the function `2 (Cos[x] + Cos[y])` and loop index `iw` and all offdiagonal elements are assigned to be 1 . However, I do not see the matrix when I try to see the matrix for each loop! Instead, the code shows the matrix outside the loop. Why is it so? [ PS: A request, I am probably paying the cost of being novice in Mathematica and got flagged for posting silly questions. So if I my question again bears stupidity, please let me know. I'd prefer to remove questions rather than getting negatively tagged. Also senior members, please assist me to follow the community guidelines. ] Thanks.


title: ListPlot with redundant labels for tick marks

text: Is it possible to have two sets of labels for the $x$-axis tick marks using `ListPlot`? For example, could we have tick marks labels: {1,2,3,4,...}, and then right below these labels, the values for the powers of two: {2,4,8,16,...}?


title: How to limit the argument to binary type in mathematica function

text: How do I write a function that only takes binary input?


title: Splitting list of sums into a sum of lists

text: I want to split a list of sums into a sum of lists For example, I'd like to turn               {{x + y + z}, {x}, {y}, {z}}       into               x{{1}, {1}, {0}, {0}} +y{{1}, {0}, {1}, {0}} + z{{1}, {0}, {0}, {1}}      as you would when turning a vector with a sum inside it into a sum of vectors. I've looked all over for documentation for this kind of thing, but I couldn't find anything :(. It would be even further helpful if someone could tell me how to put these three vectors into a single list, I know I can append them, I'm just not sure how to do it when they're being multiplied by something.


title: How to select the first "lowest" (first valley) number in a list?

text: I would to select the first valley of the list of coordinates below. It should be `{3, 1}`.               list = {{1, 7}, {2, 4}, {3, 1}, {4, 4}, {5, 9}, {6, 8}, {7, 3}, {8, 1}, {9, 4}, {10, 9}};      I thought it could be done like this               Select[list, #1[[2]] > #2[[2]] &, 1]      but that does not work. How can it be done then?


title: Find Numerical Derivatives (With Errors) given some data

text: I have a question. Suposse I have some data with their respective errors:                Data = {{0.067, 0.423}, {0.30, 0.408}, {0.60, 0.433}, {0.25,0.3512}, {0.37, 0.4602}, {0.44, 0.413}, {0.60, 0.390}, {0.73,0.437}, {0.8, 0.47}};                Errors = {0.055, 0.0552, 0.0662, 0.0583, 0.0378, 0.080, 0.063, 0.072, 0.08};      I want to construct somehow numerical derivatives of these data at som specific point. My idea was first to find a fit of those data (I'm not including the errors so far)               model = a*x + b*x^2 + c*x^3 + d*x^4 + f*x^5 + g*x^6 + h*x^7 + j*x^8 + k*x^9;          fit = FindFit[Dataf, model, {a, b, c, d, f, g, h, j, k}, z];          modelf = Function[{z}, Evaluate[model /. fitfs8]];          fits = Plot[modelf[x], {x, 0, 1.0}, Epilog -> Map[Point, Datafs8],       Frame -> True, ImageSize -> 700, PlotStyle -> Directive[Red],      FrameLabel -> {Style["x", Italic],       Style["y", Italic]}, Axes -> False,       LabelStyle -> Directive[16]];      ![enter image description here](http://i.stack.imgur.com/X8Lr0.jpg) The reason why I choose those model was matter of convenience, since I have some theoretical model (Black line, I don't add the expression because is quite "complicated"). After that I make a Numerical Derivative of the Model and compare with the analytical derivative of the theoretical prediction ![enter image description here](http://i.stack.imgur.com/AiEpV.jpg) The questions are:   1. How can I include the errors in the method?   2. How can I get then expressions (numbers) for y' at some point with errors as well? I appreciate your help :)


title: How to speed up the plot of NIntegrate?

text: Here is a toy example:               f[t_] := NIntegrate[Sin[x], {x, 0, t}];     Plot[f[t], {t, 0, 10}] // Timing      Even such a simple example will take 2.8 seconds on my computer. Since many of the plot family functions have the attribute `HoldAll`, _Mathematica_ evaluates `f[t]` only after assigning specific numerical values to `t`, thus causing a lot of repeated evaluations of `NIntegrate`. The integration to a smaller upper-limit is **forgotten** when integrating to a bigger upper-limit. On the other hand, I can’t benefit from wrapping the function in `Evaluate`, because of the numerical nature of `NIntegrate`. Something stuck here. So the questions are:   * Can NIntegrate remember or make full use of the result of a smaller upper-limit integral? Or generally, how to speed up the plot involving `NIntegrate` or is there any principle to do it?   * Is it possible to realize **ParametricNIntegrate** just as `ParametricNDSolve` in version 9? **Edit** Following the way of changing the integral to `NDSolve` (Thanks to Mark McClure) ,I’ve come up with one possible way to realize ParametricNIntegrate. Here is the code for a more general way compared to the toy model:               sol=ParametricNDSolveValue[{f'[x] == Sin[a x], f[0] == 0},f,{x,0,10},{a}];     Manipulate[Plot[{sol[a][t], 1/a (1 - Cos[a t])}, {t, 0, 10}], {a, 1, 5}]      Now I can change the parameter and plot the integral in real-time.Any better ideas?


title: Country borders and grids with CountryData plot?

text: I am coloring several countries in the manner of the following:               myc[name_] := If[name == "UnitedStates", Red, Lighter[Gray]]      Graphics[{myc[#], CountryData[#, {"FullPolygon", "WinkelTripel"}]} & /@ CountryData[]]       How can I add country borders and longitude-latitude grid, as we could do with WorldPlot?


title: Caching of Graph Graphics

text: In this answer, Simon Woods concludes that when you display an expression with head graph more than once, something is getting stored somewhere. I looked into it too and found out this only happens if you tie the expression with head graph to a symbol. Evaluating the symbol a second time will no longer print messages and results of random numbers are the same. Consider this example               g=Graph[{1 -> 2, 2 -> 3, 3 -> 1},         VertexShapeFunction -> ({Random[] // Hue, Disk[Print["VSF"]; #, 0.05]} &)]      --prints--> messages --displays-> a graph with pretty colors               g      --displays--> a graph with the same pretty colors (and does not print messages) even though               g // FullForm // InputForm      returns               (*output*)     FullForm[Graph[{1, 2, 3}, {DirectedEdge[1, 2], DirectedEdge[2, 3],            DirectedEdge[3, 1]}, {VertexShapeFunction ->              {{Hue[Random[]], Disk[Print["VSF"]; #1, 0.05]} & }}]]      where we can see that the colors are supposed to be generated randomly by the option `VertexShapeFunction`. However, after evaluating               ClearSystemCache[];      evaluating g will (almost almost surely ;) ) yield a graph with a new colors. this Q&A seems related. **My question is** : How can we turn off such caching? We can find systemoptions having to do with Cache by evaluating               SystemOptions["Cache*"]      but setting all the rules of the form               "Cache"-> True      to               "Cache" -> False      does not change the way the graphics behave. They do however influence the value of               Timing[N[Pi, 1000000]][[1]]      The SystemOptions can be set by evaluating               SetSystemOptions["CacheOptions" -> "Developer" -> "Cache" -> False];     SetSystemOptions["CacheOptions" -> "Numeric" -> "Cache" -> False];     SetSystemOptions[       "CacheOptions" -> "ParametricFunction" -> "Cache" -> False];     SetSystemOptions["CacheOptions" -> "Quantity" -> "Cache" -> False];     SetSystemOptions["CacheOptions" -> "Symbolic" -> "Cache" -> False];


title: Align plots with a shift

text: I've got three plots combined through Overlay:               Overlay[Table[Plot[PDF[NormalDistribution[0, ss], x], {x, -2, 2},           Axes -> {True, False},          Ticks -> None,          PlotRange -> {All, {0, 2}},           ImageSize -> 200],           {ss, {.2, .5, 1}}]]      Now I'd like to obtain something like this: ![overlaid plots with a shift](http://i.stack.imgur.com/IuuAc.png) I tried to achieve that with the Alignment option but to no avail. Thanks in advance for any ideas.


title: Algebraic expressions on pure functions

text: I would like to know how to perform algebraic operations with pure functions. **Simple version:** Here's a silly toy model: I want to transform the algebraic expression Sin*Cos into the function evaluated at r, that is, Sin[r]*Cos[r]. But I could not find any way to use Apply or Evaluate to do this. In other words, I would like to make Mathematica understand that something like this is true:               (Sin*Cos)[r_] := Sin[r_]*Cos[r_]      and analogously for all other possible pure functions. I realize that this means I need to convert               (Sin[#] &) Cos[#] &      into the version without the first '&', so that the whole expression defines only one pure function:               (Sin[#]) Cos[#] &[r]      * * * **Full version:** I have more elaborate algebraic expressions with two symbols 'Fm' and 'Fp', which I think as representing pure functions that will be defined later. These expressions involve derivatives, squares, etc., of these functions, for instance:               4 Fm Derivative[1][Fm]      or               1/4 Fp (2 Fp - Fp^3/Fm^2 + 4 Fm Derivative[1][Fm] Derivative[1][Fp])      I then want to simplify these expressions when I choose explicit pure functions, such as:               Fm = Sin[#] &;     Fp = Sin[#] Cos[#] &;      But when I declare these pure functions as above and then re-evaluate the expressions, I get:               4 (Cos[#1] &) (Sin[#1] &)     1/4 (Sin[#1] Cos[#1] &) (2 (Sin[#1] Cos[#1] &) - (Sin[#1] Cos[#1] &)^3/(Sin[#1] &)^2 + 4 (Cos[#1] &) (Sin[#1] &) (Cos[#1]^2 - Sin[#1]^2 &))      and instead I would like to get               4 Cos[#1] Sin[#1] &     1/4 (Sin[#1] Cos[#1]) (2 (Sin[#1] Cos[#1]) - (Sin[#1] Cos[#1])^3/(Sin[#1])^2 + 4 (Cos[#1]) (Sin[#1]) (Cos[#1]^2 - Sin[#1]^2)) &      Again, the solution to this is essentially erasing the '&'s systematically (except the very last one); but I don't know how to make it work. Any help would be greatly appreciated!!


title: How does one specify Neumann conditions for NDSolve?

text: I have a series of functions defined in my notebook, and then want to use this to solve a diffusion-reaction type equation. At the moment, something like this works:               eqn = D[Ef[r, t], t] - Def5*( D[Ef[r, t], r, r] + (2/(r )) D[Ef[r, t], r]) + q[r ] Ef[r, t];          ic = {Ef[r, 0] == kd, Ef[rn, t] == kd, Ef[ro, t] == kd};     s = NDSolve[{eqn == 0}~Join~ic, Ef, {r, rn, ro}, {t, 0, 60}];     Plot3D[Evaluate[Ef[r, t] /. s], {r, rn, ro}, {t, 0, 60},      AxesLabel -> Automatic, MaxRecursion -> 8, PlotPoints -> 32]       `q[r]` is a predefined function, and `rn` and `ro` are boundaries. Now, the above code works, BUT the boundary conditions are wrong and inconsistent; I want to change them to Neumann-type conditions so that $d(Ef)/dr = 0$ at both `ro` and `rn`, but I am having trouble implementing this; any ideas? All help appreciated! This is probably very simple but I'm new to _Mathematica_ and sometimes screw up syntax... DRG


title: Tolerance on image detection

text: I am trying to pick straight lines out of a photographic image to track deflections of a structure under loading. The idea is to mark certain structural elements with straight lines, as they deflect I can then extract the vectors from the image, find the centre point and measure how much the vector has displaced and rotated since the previous image. I have used the second answer to this question to get the bones of a process up and running, code below. Using a generic vector image from web as a test image until my rig is up and running.               i = ColorNegate@Import["http://goo.gl/5R4MAl"]     p = Closing[Binarize@ImageTake[i, 480], 1];     Show[p, Graphics[{Thick, Orange, Line /@ ImageLines[p, Segmented -> True]}]]      This seems to pick out the main vectors quite nicely. ![enter image description here](http://i.stack.imgur.com/aAAuk.jpg) Whilst the length of results returned is 9 which is what I expected, on closer inspection there are additional points that are very close together. I'm guessing that these are rogue parts of the arrowheads being picked up.               foo = ImageLines[p, Segmented -> True]     foo // Length      Gives the length to be 9, which is what I would have expected given there are nine primary lines, but if displayed in matrixform I'm getting additional rogue points.               foo // MatrixForm      ![enter image description here](http://i.stack.imgur.com/euLlA.png) I've tried to tinker about with limiting the MaxFeatures->9 and using DeleteSmallComponents with mixed success. I'm struggling to find any helpful guidance too on what t (threshold) and d (distinctness) actually do with regards the image recognition. What methods exist to extract the primary straight lines from images like the one used in the above example that gives a start and and end co-ordinate of the line that can then be used in subsequent calculations please? I'd also prefer to not set the number of lines I'm expecting Mathematica to discover.


title: Overlay graphic3d and ListPlot3d

text: I have this two input for 3d graphics: 1) A `Graphics3D` with a `Polygon` that represents a map 2) A `ListPlot3D` with geo location and prices How can i overlay the `ListPlot3D` over the map?               Mappa = Graphics3D[{Texture[laMappa], Polygon[{{-1, -1, 0}, {1, -1, 0}, {1, 2, 0},               {-1, 2, 0}}, VertexTextureCoordinates -> {{0, 0}, {1, 0}, {1, 1}, {0, 1}}]}]          Grafico = ListPlot3D[{idatinuovi}, Mesh -> None,  ImageSize -> 1000, Filling -> Bottom,                                        PlotRange -> {{41, 43}, {12, 13}, {1000, 15000}},                                        InterpolationOrder -> 0, ColorFunction -> "Rainbow"]      this is a .csv file related to data by:                idatinuovi = idati[[All, 2 ;; 4]]


title: Is it possible to use the global definition of a symbol as part of its new local definition set in a Block?

text: Consider the following expression:               Block[{Print = CopyToClipboard}, doSomething[]]      Now all calls to `Print` from `doSomething[]` (and all other functions it may invoke) will copy values to clipboard rather than printing them. Suppose I want to modify this expression such that all calls to `Print` _both_ copy values to clipboard _and_ print them.               Block[{Print = (CopyToClipboard[#]; Print[#]) &}, doSomething[]]      I need something similar to the above, but which doesn't go into an infinite recursion and where the occurrence of `Print` within the pure function retains its global meaning that it had before the execution entered the `Block`. Is it possible? I understand the difference between `Block` and `Module`, and I need to use `Block` to create a dynamic (rather than lexical) scoping.


title: How do you use ssh-keys instead of a password to run a remote-kernel over ssh?

text: On all of the servers that I regularly interact with, I have ssh-keys setup for passwordless access via ssh. Yet, every time I attempt to start a remote- kernel over ssh, I get asked for a password despite having a ssh-keyagent up and running. What do I need to add to the kernel configuration, such as the launch command, to have it use my ssh-key? I'm using version 8.0.4 on MacOS 10.6.8, if it makes a difference.


title: How to input and output partitioned matrices that show partitions and compute as normal?

text: I want to demonstrate multiplication of partitioned matrices as in the example here. Using the Insert Menu, you can build a matrix and draw lines between rows and columns. However, I want to be able to select a subset of the column and row lines to add. Then I would like to be able to compute the product of two compatible matrices (so the column partition on the left matrix matches the row partition on the right matrix) and have the output show as a partitioned matrix where the row dividers come from the row dividers on the left matrix and the column dividers come from the column dividers on the right matrix. In the example, this matrix multiplication is defined since the column partition on A is, say, {2} and the same for the row partition on B. The resulting matrix C has a row divider after the 2nd row since A had a row divider after the second row (while B had no column dividers).


title: Slideshow export to PDF (landscape)

text: I am trying to export a slideshow presentation that I have created in Mathematica 8 to PDF in case I have trouble connecting my laptop to a projector for a talk I am giving. I have tried Save As, Save Selection As and Print to File after setting the print setting to landscape. The latter produces the enormous file even though I am rasterizing the plots, and the former save my file normally but in portrait and also lose the `Scaled` option of plots. Does anyone have a good workflow for this or know a good reference?


title: Import variables from excel

text: I have an excel file I am trying to import into _Mathematica_. It is in the format of variables names in the first row of each column and then the data underneath. One of the columns specifies hospitals, so each row has observations for each hospital. There are approximately 800 columns and 400 hospitals. Some of the variables are strings, such as hospital names. But I could delete those in Excel if I had to. I imported it and the text filled the screen. I saw there was a program called Link for Excel, but $250 is a lot of money for a student. I then tried ListView and the table was all wrong. The variable names were in the first column. It seems like it was displaying rows as columns. I want to make sure what is imported is a matrix so we can run some linear programming. Thanks so much for your help in advance!


title: Plotting vector field

text: I am trying to plot the gradient vector field of the function: $$sin\sqrt{x^2+y^2}$$ This is my code: `VectorPlot[{(x cos (sqrt (x^2 + y^2)))/ sqrt (x^2 + y^2), (y cos (sqrt (x^2 + y^2)))/ sqrt (x^2 + y^2)}, {x, -3, 3}, {y, -3, 3}]` I am getting a blank plot. Am I entering the code correctly?


title: Why does a variable become real when using the second argument of Dynamic?

text: Version 9, on Windows 7. Please compare these 2 very simple examples               Manipulate[      x,      Row[{Manipulator[Dynamic[x], {0, 10, 1/100}]}],      {{x, 1}, None}     ]      and               Manipulate[      x,      Row[{Manipulator[Dynamic[x, (x = #) &], {0, 10, 1/100}]}],      {{x, 1}, None}     ]      Would you not expect them to work the same way? But in the second case, `x` becomes real, while in the first case it remains rational as expected. ![enter image description here](http://i.stack.imgur.com/90wtw.png) Any one knows the reasoning for this? I read about the second argument of `Dynamic`, but do not see something obvious. I did not read everything about it. But the above behaviour is surprising. The fix is easy:               Manipulate[      x,      Row[{Manipulator[         Dynamic[x, (x = Rationalize[#]) &], {0, 10, 1/100}]}],      {{x, 1}, None}      ]      I am only asking why it happens, so I can learn more.


title: Obtaining the Fourier transform of an operator

text: We know that if we have a function $f(x)$, and we call $g(\omega)$ its Fourier transform, then the Fourier transform of $x f(x)$ is $$\imath \frac{\mathrm{d} g(\omega))}{\mathrm{d}\omega} $$ and viceversa, ${\mathrm{d}f(x)}/{\mathrm{d}x}$ becomes $\imath \omega g(\omega)$ How can I transform an operator, for example $$L(u):= u+x^2 u+\frac{\mathrm{d}u}{\mathrm{d}x}$$ so that it automatically transforms the variable $x$ in a derivative, and $x^n$ in a derivative of order $n$, and viceversa? * * * The specific problem I'm solving is the following:               ds[u_]:= A.{D[u,s],D[u,a]}.{1,0}      da[u_]:= A.{D[u,s],D[u,a]}.{0,1}   (*where A is a constant matrix*)     G1[t_,x_,y_]:= G0[t,x,y] Integrate[K[tt,x,y],{tt,t,T}]      where `G0[t,s,a]` is a function already defined and `K[tt,x,y]` is a polinomial in the variable $x$ and $y$ (so its integral is still a polinomial in $x$ and $y$. I know the inverse Fourier transform of `G0`, that is `F0[t,s,a]`. What I need is to do is to substitute all the the $x^n$ terms in the polinomial with the operator `ds[u]` applied $n$ times and similarly, $y^m$ with the the operator `da[u]` applied $m$ times. thanks


title: Sorting function for non commuting bosons

text: I am trying to write a sorting function which will sort expressions involving products of bosonic objects which do not commute. For example, I can have objects like $a,\ a^\dagger,\ b,\ b^\dagger$ where $[a,a^\dagger] = 1$, $[b,b^\dagger] = 1$ and all other commutators vanish. So the sorting function should for example sort a term like $a a^\dagger b a b^\dagger$ to $a + a^\dagger a^2 + ab^\dagger b+a^\dagger a^2 b^\dagger b$. So I want to sort them first lexicographically and then such that dagger objects will be to the left. I have been using a function from the notebook provided by Simon Tyler at arXiv. So I have created a bosonic object using               Clear[Boson, BosonC, BosonA]     Boson /: MakeBoxes[      Boson[cr : (True | False), sym_], fmt_] :=        With[{sbox = If[StringQ[sym], sym, ToBoxes[sym]]},         With[{abox = If[cr, SuperscriptBox[#, FormBox["\[Dagger]", Bold]], #] &@sbox},         InterpretationBox[abox, Boson[cr, sym]]        ]       ]     BosonA[sym_:"a"] := Boson[False, sym]     BosonC[sym_:"a"] := Boson[True, sym]      Next following the code, I want to turn the atomic expressions into lists and use Mathematica's ordering after that. So I have functions like               NCOrder[Boson[cr : True | False, sym_]] :=      NCOrder[Boson[cr, sym]] = {{sym, Boole[cr]}}          NCSort[expr_] := Module[{temp, NCM},      temp = expr /. NonCommutativeMultiply -> NCM;         temp = temp /. (NCM[a : _Boson ..] :> (NCM[a][[Ordering[#]]] &[NCOrder /@ {a}]));      temp = temp /. NCM[singleArg_] :> singleArg;         temp = temp /. NCM -> NonCommutativeMultiply     ]      This NCSort function right now is just ordering it lexicographically. I want to modify it so that it would sort and order the expressions as I explained above.


title: Export Voronoi Data?

text: I have a Voronoi Diagram (VoronoiDiagram/DiagramPlot) in Mathematica - and I need to export it to something I can import into ArcGIS. I'll settle for any sort of data/vector-based (non-raster) format. Obviously a BMP isn't going to work - I need the actual vector data.