Skip to main content

Posts

Showing posts from June, 2016

programming - Underlying Algorithms for List Manipulation Functions

Does anyone know, or know of a link for the underlying algorithms used for list manipulation functions? Specifically: This question came about as there has been numerous questions that have asked about removing duplicates or some other list manipulation function. Some of the best answers often included usage of GatherBy due to the performance it offers (a few hundred times over DeleteDuplicates in some cases). Related: usages of GatherBy and Union : Checking for duplicates in sublists GatherBy Ordering function with recognition of duplicates GatherBy Gather list elements by labels GatherBy Retrieving duplicates from nested list GatherBy Delete duplicates in a list, depending on the sequence of numbers Union How to delete duplicate solutions of this system of equations? Union Answer This is an incomplete answer; I will continue it tomorrow. Work In Progress: errors may abound. Preamble hat-tip to Leonid For the variations with custom test or ordering functions we can snoop on...

sorting - Sort lists according to the order of another

I have three parallel lists (i.e. the elements in position i of each list are related). I want to sort the first list using the function Sort and make the same changes to the other lists so that I still have parallel lists when finished. How can I do this? As an example: Given the lists {2, 3, 1} , {a, b, c} , and {alpha, beta, gamma} , sorting all lists according the first one gives {1, 2, 3} , {c, a, b} , and {gamma, alpha, beta} . Answer lists = {list1, list2, list3} = {{1, 3, 2}, {a, b, c}, {x, y, z}}; Another option SortBy[lists\[Transpose], First]\[Transpose] {{1, 2, 3}, {a, c, b}, {x, z, y}}

filtering - How to define patterns to filter lists of replacement rules?

I have a list of replacement rules like list˘of˘replacement˘rules = {a -> 1, b -> 2, c -> 3} and I want to filter all the elements where ‘b’ is replaced by using Cases[] and a search pattern. The following code does not work: Cases[list˘of˘replacement˘rules, b -> _] The result is {}. How should the pattern be defined in order to get the result {b->2}? Answer rules = {a -> 1, b -> 2, c -> 3}; FilterRules[rules, b] (* {b -> 2} )* or Cases[{a -> 1, b -> 2, c -> 3}, PatternSequence[b -> _]] (* {b -> 2} )*

data - Calculate mean of values in bins

Problem description : I have a list of $\{x,y\}$ pairs. I'd like to divide $x$ into equal[*] bins, say $bx_1, bx_2, \ldots$, calculate $\left $ for every bin and then plot the bin values versus the means. I.e., plot $bx_i$ versus $\left _i$ with ListPlot[] . Question : I've already done it manually, but I was wondering whether: exists a builtin function in Mathematica that does what I want. exists a builtin function that I could use in my own implementation (e.g., HistogramList[] ?). [*] Bins with the same interval. Equally spaced intervals. EDIT : Largely off-topic, but in R it turns out to be very easy with the fields package: library('fields') df off-axis distance, V2 -> energy st df2 # Plot mean energy for every distance bin # EDIT: Actually I should plot against `centers` of `st`, but anyway. names(df2) plot(df2$mean.energy, type="s", xlab="Off-axis distance (mm)", ylab="Mean Energy (MeV)")

undocumented - What's purpose of function System`Private`LookupCodeByName and System`Private`LookupNameByCode

These are useful function from their name.And they are kernel function since: System`Private`HasAnyCodesQ /@ {System`Private`LookupCodeByName, System`Private`LookupNameByCode} {True, True} But how to use them? Answer I tried random inputs until I got something to work. It looks like it has to do with character codes: System`Private`LookupNameByCode[200] "CapitalEGrave" And then we can plug in the reverse: System`Private`LookupCodeByName["CapitalEGrave"] 200 Edit by yode as J.M.'s comment Grid[DeleteCases[{#, FromCharacterCode[#], System`Private`LookupNameByCode[#]} & /@ Range[100, 174], {_, _, $Failed}], Frame -> All]

plotting - Exporting the solutions of two ODEs sampled over same set of domain values

I have two functions h[r] and c[r] defined by two differentials equation. operh[h_, c_] := h D[D[h, r, r] + 1/r D[h, r] - 1/lc^2 h, r] + m D[c, r] operc[h_, c_] := cmt h + ξ ((ϕ - H) (1 - ϕ) )/Sqrt[R^2 - r^2] - 1/r D[r h (1 + e^2/1680 h^4 (D[c, r])^2 ) D[c, r], r] I numerically solve them with NDSolve (with solution called "s" here under). Then I plot them with no problem: Plot[h[r] /. s, {r, rleft, rright}, PlotRange -> All, AxesOrigin -> {0, 0}] Plot[c[r] /. s, {r, rleft, rright}, PlotRange -> All, AxesOrigin -> {0, 0}] With same rright and rleft. What I want is a .txt or .csv file with 3 columns, containing r, h[r] and c[r] . From what I have read on this site I tried: ploth = Plot[h[r] /. s, {r, rleft, rright}, PlotRange -> All, AxesOrigin -> {0, 0}] pts = Cases[ploth, Line[pts_] :> pts, Infinity][[1]] Export[SystemDialogInput["FileSave", "filename.csv"], pts] It gave me a list of 157 pairs {r, h[r]} . Then I did the ...

input forms - Explicit digit-count (precision) of real number

I'm trying to do something which seems like it should be simpler than it is (at least in my attempts at it). I have some code where I read in a list of numbers generated for me by a coworker. These numbers have varying degrees of precision, so the list might be something like hisData = {0.05467`, 12.34230`, 4.69`, 9.3452`, ...} I want to compare this to a list that I generated myself, where each number may have different degrees of precision. So my list might read myData = {0.0547`, 12.34231`, 4.6877`, 9.345`, ...} The comparison must be done in a particular way: all decimal digits that appear explicitly in any pair of numbers must agree exactly, except possibly when rounding is needed to make the number of digits match, but the number of digits given need not be the same. In the example lists above, all but the second pair of elements "agree". This would be easy enough to do if I had a function called, say, explicitPrecision, that counted the number of digits before and ...

security - List of dangerous functions

I no longer believe it is a good idea to work with a list of dangerous functions. I have tried to edit the question is such a way that it respects my earlier perspective, but also so that it does not give users the impression that certain things are safe whereas they are not. Introduction Running Mathematica code from an untrusted source is very dangerous. Luckily, Mathematica warns you about dynamic content, so that it is probably safe to open any notebook, as long as you don't run it. I feel it would be very nice to have a list of functions that can harm a system. Using this, we may be able to make a function that checks if a notebook is safe. A function to find dangerous functions A naive function to see if dangerous expressions are present in a NotebookObject , could be the following. Through[{Unprotect, ClearAll}[dangerousFunctionsQ]] dangerousFunctionsQ[nb_] := ! FreeQ[ToExpression[Unevaluated[#], StandardForm, HoldComplete] & /@ (NotebookRead@Cells[nb])[[All,...

Find all roots of an interpolating function (solution to a differential equation)

I'm trying to find all the roots of the solution to a differential equation. Using NSolve or Reduce I don't get the roots, so I'm using an iterative method which I found in physicsforums.com. This method solves my problem, but you have to choose the increment and thus in some cases it might give headaches. I wonder if there is any more general approach. Following is a sample code: data = NDSolve[{1.09 x''[t] - 0.05 x'[t] + 1.1759 Sin[x[t]] == 0, x[0] == Pi/3, x'[0] == 0}, x, {t, 0, 50}] {{x->InterpolatingFunction[{{0.,50.}}, ]}} First attempt with no success: sol = NSolve[x'[t] == 0 /. data , t] NSolve::ifun: Inverse functions are being used by NSolve, so some solutions may not be found; use Reduce for complete solution information. >> {{t->InverseFunction[InterpolatingFunction[{{0.,50.}}, ],1,1][0.]}} Second attempt with no success: sol = Reduce[x'[t] == 0 /. data , t] Reduce::inex: Reduce was unable to solve the system with inex...

timing - Throwaway memoization makes built-ins faster?

Yesterday I got into an argument with @UnchartedWorks over in the comment thread here . At first glance, he posted a duplicate of Marius' answer , but with some unnecessary memoization: unitize[x_] := unitize[x] = Unitize[x] pick[xs_, sel_, patt_] := pick[xs] = Pick[xs, sel, patt] and proposed the following test to justify his claim that his approach is faster: RandomSeed[1]; n = -1; data = RandomChoice[Range[0, 10], {10^8, 3}]; AbsoluteTiming[Pick[data, Unitize@data[[All, n]], 1] // Length] AbsoluteTiming[pick[data, unitize@data[[All, n]], 1] // Length] (* {7.3081, 90913401} {5.87919, 90913401} *) A significant difference. Naturally, I was skeptical. The evaluation queue for his pick is (I believe) as follows: pick is inert, so evaluate the arguments. data is just a list, 1 is inert, data[[All, n]] quickly evaluates to a list unitize@data[[All, n]] writes a large DownValue ... ...calling Unitize@data[[All, n]] in the process, returning the unitized list. Another large DownV...

differential equations - NDSolve's output ignores multiple valid solutions

I'm looking for solutions to a boundary problem involving a non-linear Hamiltonian $$ H(q,p) = \frac{1}{4}\left(q^{2}+p^{2}\right)^{2}, $$ whose solutions are oscillatory but have a complex time dependence. I'm interested in all possible solutions $\left(q(t),p(t)\right)$ that satisfy the following boundary conditions: $$\begin{cases} q(0)&=-1 \\ q(\pi)&=1 \end{cases}$$ and I am absolutely sure there are a lot (maybe an infinity) of them. When I ask Mathematica to solve the boundary problem NDSolve[{q'[t] == p[t] (p[t]^2 + q[t]^2), p'[t] == -q[t] (p[t]^2 + q[t]^2), q[0] == -1, q[Pi] == 1}, q[t], {t, 0, Pi}] I get only one solution, which looks like and satisfies the boundary problem. What I can't figure out is that, manually, I found another solution: NDSolve[{q'[t] == p[t] (p[t]^2 + q[t]^2), p'[t] == -q[t] (p[t]^2 + q[t]^2), q[0] == -1, p[0] == 1.200859}, q[t], {t, 0, Pi}] whose graph is . How can I manipulate NDSolve such that it display...

plotting - How can I get the interpolation line from Listplot?

My discrete data is: data={{1., 5827.}, {3., 6242.}, {15., 6066.}, {60., 5972.}} I want to get the interpolation line from ListPlot: pl=ListPlot[data, Joined -> True,InterpolationOrder -> 2, PlotStyle -> Dashed, PlotRange -> {5900, 6300}] But I can't get the same dash line from InterPolation: Plot[Interpolation[data, InterpolationOrder -> 2][x], {x, 1, 60}, PlotStyle -> Dashed] How can I get the curve like 1st one? Or can I get the dash line points automatically? Answer More of an extended comment here, but it is clear that ListPlot must use a different interpolation function than is available via Interpolation . We can get a very close approximation by Janus's answer here and amending the Method option. The idea is that ListPlot performs an parametric interpolation on the x and y axes separately. data = {{1., 5827.}, {3., 6242.}, {15., 6066.}, {60., 5972.}}; pl = ListPlot[data, Joined -> True, InterpolationOrder -> 2, PlotStyle -> Dashe...

finite element method - Interface points of NDEigensystem

When solving an eigenvalue problem with "NDEigensystem", e.g. a 1D Eigenvalue problem with the interval composed of different materials, which should be solved by the pure numerical method such as FEM, is it necessary to assign the locations of the interface points? In fact, for FEM we know that the interface points should be the element boundary nodes. For example, the code assigning the interface positions is (*This seems to be incorrect.*) Needs["NDSolve`FEM`"] a1 = 0.1; a2 = 0.2; a3 = 0.3; a4 = 0.4; b = 0.7; sigma0 = 37*10^6; omega = 2*10^3*Pi; mu0 = 4*10^-7*Pi; sigma = If[a1 bm = ToBoundaryMesh[ "Coordinates" -> {{0}, {a1}, {a2}, {a3}, {a4}, {b}}, "BoundaryElements" -> {PointElement[{{1}, {2}, {3}, {4}, {5}, {6}}]}]; bm["Wireframe"] rm = ToElementMesh[bm, "MeshOrder" -> 2, MaxCellMeasure -> 0.0001]; L = psi''[r] + psi'[r]/r - (1/r^2 + I*omega*sigma*mu0)*psi[r]; B = DirichletCondition[psi[...

graphics - ListPlot: Plotting large data fast

Mathematica produces fantastic-looking graphics, but it can be slow on large data sets. Here is an example for a (random) time series: rv = RandomVariate[ExponentialDistribution[2], 10^5]; Plotting this takes quite some time: t = AbsoluteTime[]; ListLinePlot[rv, PlotRange -> All] AbsoluteTime[] - t (* Put this line into the NEXT cell, and evaluate both cells together*) The option PerformanceGoal->"Speed" has no effect. Turning off antialiasing makes it much faster (but if you increase the data size to 10^6 instead of 10^5, it is still VERY slow. A time series of 10^6 points is quite reasonable in my applications): t = AbsoluteTime[]; Style[ListLinePlot[rv, PlotRange -> All], Antialiasing -> False] AbsoluteTime[] - t (* This line in separate cell! *) Reducing the number of MaxPlotPoints makes it much faster, but completely distorts the shape of the data: t = AbsoluteTime[]; ListLinePlot[rv, PlotRange -> All, MaxPlotPoints -> 1000] AbsoluteTime[] - t (* This...

calculus and analysis - Maximize violating constraints

I have Maximize[{(h*10)/(300*(100 - (l^.5 + d^.4 + H^.6))), (l + d + H + h) == 669, l > 0, d > 0, H > 0, h > 0}, {h, l, d, H}] I believe this should maximize my formula, but when I run it, I get The function value 12073.4 -0.284149 I is not a real number at {d,h,H,l} = {-28.781, 682.337, 2.43779, 13.006} Can anyone help me figure out how to enter this problem so I get the correct answer? I thought I'd explicitly made it not evaluate the formula for any variable less than 0, but it seems to be trying to do that. Answer I suspect you are correct in your assessment. Since there are approximate numbers in the input, Maximize punts to NMaximize , which uses penalty methods to enforce some constraints (not sure why it needs them here for linear constraints; I need to check into that). You can get better behavior by forcing real values. NMaximize[{Re[(h*10)/(300*(100 - (l^.5 + d^.4 + H^.6)))], (l + d + H + h) == 669, l > 0, d > 0, H > 0, h > 0}, {h,...

programming - when is f@g not the same as f[g]?

I have always thought that f@g will give the same result as f[g] in all cases, and it is just a matter of style which one to use and that g will always evaluates first, and then f will evaluate using the result of g evaluation. I never thought that there can be any precedence issue here, since no one ever mentioned it in all the times I have been using Mathematica. So I was really surprised when I found one case where this was not so. So my question is: How does one know when f@g is not the same as f[g] ? The help says nothing about this (thanks to chat room for giving me the link to this, I searched and could not find it) http://reference.wolfram.com/mathematica/ref/Prefix.html Even though one can see the word precedence and grouping but no explanation of where these are talked about and no more links to follow Prefix[expr, h, precedence, grouping] can be used to specify how the output form should be parenthesized. clearly this is a precedence issue. But I have never seen thi...

core language - Custom operator form for system functions?

There are some functions for which an operator form would make sense but isn't implemented. For example, I'd like to make this work: RandomSample[3] @ Range[10] You can do something like this: Unprotect[RandomSample]; Off[RandomSample::lrwl] RandomSample[n_][x_] := RandomSample[x, n] But I wonder if there is a standard technique to define operator forms for system functions?

output formatting - Format something Inactive

ops = Map[Inactive, Tuples[{Times, Plus, Subtract, Divide}, {4 - 1}], {2}]; rules = Thread[{a, b, c} -> #] & /@ ops; matchQ[list_] := Module[{}, res = Fold[Replace[#1, List -> #2[[1]], {#2[[2]]}, Heads -> True] &, list, Transpose@{{a, b, c}, {1, 2, 3}}] /. rules //. Inactive[Subtract][x_, y_] :> Inactive[Plus][x, -y]; Select[res, Activate[#] == 24 &]] list = {6, 4, 2, 3}; pattern = {{#[[1 ;; 2]], #[[3 ;; 4]]}, {#[[ 1]], {#[[2]], #[[3 ;; 4]]}}, Reverse@{#[[1]], {#[[2]], #[[3 ;; 4]]}}, {{{#[[1]], #[[2]]}, #[[ 3]]}, #[[4]]}, Reverse@{{{#[[1]], #[[2]]}, #[[3]]}, #[[4]]}} &; parts = Flatten[(pattern) /@ Permutations[list, {4}], 1]; res1 = matchQ /@ parts // Quiet // Flatten; When I'm playing with a calculating 24, I used some function like Inactive . My question is how can I make the output more reasonble in human writing habbit? For example, here I replace Subtract with Plus . I want to change ((6+4)+-2)*3 to ((6+4)-2)*3 StringForm is...

web access - More complete "MutipartData" POSTs using URLFetch

I'm trying to upload a file to a wiki using Mathematica 's URLFetch . Some preliminary work is required - one needs to log in to the wiki, obtain an edit token, and exchange some cookies with the wiki server. URLFetch can handle all of that. But then for the actual file upload, the MediaWiki API (e.g. http://www.mediawiki.org/wiki/API:Upload ) specifies that one has to use the POST method with Content-Type=multipart/form-data , but the Mathematica implementation of multipart within URLFetch seems to be incomplete. The Mathematica documentation of the option "MultipartData" of URLFetch says only: to upload multipart data, each part must be of the form {name, mimeType, {bytes}} , where {bytes} is a list of bytes This suggests that in the "Content-Disposition" header of each part of a multipart request, a name can be specified as in: Content-Disposition: form-data; name="file" but there seems to be no way to also specify a filename, as in Conte...

plotting - How do I create a triangulated surface from points?

I have a set of points in a nx3 matrix and I would like to convert them into a surface, so that I may calculate its surface area. The function ListSurfacePlot3D creates the surface how I want it. How do I compute the surface area from there? I have tried other functions, such as DelaunayMesh. However, DelaunayMesh creates a closed region defining a volume, causing the surface area to be greater than expected. Which function do I use? If this surface is created correctly, I would use RegionMeasure to get the area. Edit: Here is my set of data points https://www.dropbox.com/s/db79ewahx24it3r/Data.txt?dl=0 Using Geomagic Control allows me to get the area of the wrapped points, at about 11522 mm2. I don't need the visual, I only want the area from these points. I'll be running the area calculation through a routine of about 6000 samples. Answer It turns out ListSurfacePlot3D does a terribly poor job of approximating the surface in the OP, otherwise one will just apply Discretize...

plotting - Animation: presenting events that persist in time

I have a list of events that occurred on specific moments in time. In the example below, the three events occurred at 23:09, 23:13 and 23:17. I was able to build an animation that flows from time 23:00 to 0:00 and shows at the right time the event by plotting a disc of a specified magnitude at the specified coordinates. data = {{{2012, 3, 19, 23, 9}, 44.891, 12.202, 2.6, 2.5, 3541187340}, {{2012, 3, 19, 23, 13}, 44.898, 11.258, 6.2, 4.1, 3541187580}, {{2012, 3, 19, 23, 17}, 44.186, 11.496, 4.2, 2.2, 3541187820}}; Each event is specified by its time of occurrence, x and y coordinates, magnitude 1, magnitude 2, absolute time. First, I specified the start and end time of the animation in absolute time and build a list of times from start to end in steps of one minute: start = AbsoluteTime[{2012, 3, 19, 23, 0, 0}]; end = AbsoluteTime[{2012, 3, 20, 0, 0, 0}]; time = Table[start + 60*i, {i, 0, (end - start)/60}]; For my convenience I reshuffled my data: d3 = Table[{data[[i, 6]], dat...

plotting - adjust the alignment of different legends

I used Show to combine a ListPlot and Plot , and each plot has its own legend. The code is as follows: list1 = {{0.1, 0.3}, {0.3, 0.42}, {0.5, 0.48}, {0.7, 0.57}, {0.9, 0.72}}; list2 = {{0, 0.2}, {0.2, 0.25}, {0.4, 0.31}, {0.6, 0.34}, {0.8, 0.39}}; curve1[x_] := 0.25 + 0.5 x; curve2[x_] := 0.2 + 0.25 x; p1 = ListPlot[{list1, list2}, PlotRange -> {{0, 1}, {0, 0.8}}, PlotStyle -> {{Red}, {Orange}}, Frame -> True, PlotLegends -> Placed[PointLegend[{"Experimental data 1", "Experimental data 2"}], {0.8, 0.3}], ImageSize -> Large]; p2 = Plot[{curve1[x], curve2[x]}, {x, 0, 1}, PlotRange -> {{0, 1}, {0, 0.8}}, PlotStyle -> {Blue, Black}, Frame -> True, PlotLegends -> Placed[LineLegend[{"Theoretical line for data 1", "Theoretical line for data 2"}], {0.8, 0.2}], ImageSize -> Large]; Show[p1, p2] Now I want the four legends to have the same vertical spacing between ...

counting - Count Column Elements in a Matrix

Let a square matrix be: mat={{1,2,3},{4,5,6},{0,0,0}}; Get a total for each column: In: Total[mat] Out: {5,7,9} I would like to obtain the number of non-zero elements of each column (and the result to be in the form {a,b,c}). Thanks Answer I think the simplest version is: Total @ Unitize @ mat {2, 2, 2}

files and directories - DumpSave'ing while lengthy program runs

tldr: is there a way to load all files following a structure/append to DumpSave file i would like to run a script following the structure: Get['lib.wl'] data = initializeDataSet[] DumpSave['save0', data] result1 = compute1[] DumpSave['save1', result1] Clear[result1] . . . resultn = computen[] DumpSave['saven', resultn] Clear[resultn] is there a good way to load all saves without going >>save0 to >>saven. Or is this manual work needed? A Get[name] where name allows for joker would suffice, but I don't think it exists. Putting them into a separate folder seems sensible. Just saving at the end doesn't seem viable/safe to me in case of an unexpected early exit. Edit: just discovered FileNames[] is a thing. the solution would be like: names = FileNames[pattern,{directory}] Do[Get[fileName], {fileName, names}] Answer You can use Get/@ FileNames["save*"] .

performance tuning - Happy 2K prime question

This being the Q number 2K in the site, and this being the day we got the confirmation of mathematica.se graduating soon, I think a celebration question is in order. So... What is the fastest way to compute the happy prime number 2000 in Mathematica? Edit Here are the Timing[ ] results so far: {"JM", {5.610, 137653}} {"Leonid", {5.109, {12814, 137653}}} {"wxffles", {4.11, {12814, 137653}}} {"Rojo", {0.765, {12814, 137653}}} {"Rojo1", {0.547, {12814, 137653}}} Answer This answer should be read upside down, since the last edit has the fastest, neatest and shortest answer Module[{$guard = True}, happyQ[i_] /; $guard := Block[{$guard = False, appeared}, appeared[_] = False; happyQ[i] ] ] e : happyQ[_?appeared] := e = False; happyQ[1] = True; e : happyQ[i_] := e = (appeared[i] = True; happyQ[#.#&@IntegerDigits[i]]) Now, taking this from @LeonidShiffrin happyPrimeN[n_] := Module[{m = 0, pctr = 0}, While[m ...

differential equations - (NDSolve) Non-linear 2nd order ODE, regular singular point (looking for good methods for this problem)

I am solving this set of non-linear 2nd order ODE by NDSolve, $$r^2\frac{d^2f}{dr^2} = 2f(1-f)(1-2f)+\frac{r^2}{4}h^2(f-1)$$ $$\frac{d}{dr}\left[r^2\frac{dh}{dr} \right]=2h(1-f)^2+\lambda r^2(h^2-1)h$$ $$ f(0)=h(0)=0, \quad f(\infty)=h(\infty)=1$$ From the series solution, I also know the behavior at $r\approx0$, $$ f\approx \alpha r^2, \quad h\approx \beta r$$ Therefore the singularity at $r=0$ is a regular singular point. I tried using NDSolve, \[Lambda] = 625/2048; eqn = {r^2*f''[r] == 2 f[r] (1 - f[r]) (1 - 2 f[r]) - r^2/4 (h[r])^2 (1 - f[r]), D[r^2*h'[r], r] == 2 h[r] (1 - f[r])^2 + \[Lambda]*r^2 (h[r]^2 - 1) h[r]}; bc = {f[0.000001] == 0, f[20] == 1, h[0.000001] == 0, h[20] == 1}; sol = NDSolve[{eqn, bc}, {h, f}, {r, 0.000001, 20}] However, it doesn't solve it says At r =...., step size is effectively zero; singularity or stiff system suspected. I also tried, NDSolve[{eqn, bc} /. {f -> (#^2 g[#] &), h -> (# j[#] &)}, {g, j}, {r, r1, r2}] with boundar...

graphics - How to draw rectangle in mathematica with arrow

I have a basic question but i cant seem to find an answer for it. How do i draw a rectangle in mathematica with either clockwise or counterclockwise arrow heads. In this picture i know two points i.e. (xmin,ymin) and (xmax,ymax). I want to stepwise quench my system and show it in this picture and hence many rectangles between min and max. Rectangle command doesnt seem to provide an answer. And i want to do it with very many points, hence couple of answers that are available in forum are not so useful. (-1+X+1/(1/3+(-3.+3 X)/(3 (1.*10^-18+4.5*10^-12 (2+X)+2.59808*10^-9 Sqrt[4.*10^-12-4 (-1+X)^3+1.*10^-6 (8-(-20+X) X)])^(1/3))+333333. (1.*10^-18+4.5*10^-12 (2+X)+2.59808*10^-9 Sqrt[4.*10^-12-4 (-1+X)^3+1.*10^-6 (8-(-20+X) X)])^(1/3)))/X This is the plotted function f. Now i want to select two X values say X=1.1 and X=3 and divide interval into N parts X = {X1, X2, X3, X4, X5} and directed rectangle between each of the two consecutive values. Because i cant do it in mathematica i drew ...

plotting - PlotLegends "Expressions" default for ListPlot?

gwr pointed out an apparent problem with setting a default Option value of PlotLegends for DateListPlot . I traced the problem to the existence of PlotLegends -> "Expressions" in the default PlotTheme . This is true for ListPlot as well. Since the explicit options within the Plot Theme override options set with SetOptions it appears that our option is ignored. Actually "Expressions" is used across many Themes for ListPlot : PlotLegends /. Charting`ResolvePlotTheme["Default", ListPlot] & /@ {"Default", "Web", "Classic", "Detailed", "Scientific"} {"Expressions", "Expressions", "Expressions", "Expressions", "Expressions"} Why would this be? This option value doesn't seem to even work with these plot types, i.e. dat = Table[{k, PDF[BinomialDistribution[50, p], k]}, {p, {0.3, 0.5, 0.8}}, {k, 0, 50}]; ListPlot[dat, Filling -> Axis, PlotLege...

equation solving - Long Execution Time of Reduce

I have a set of linear homogeneous equations with 8580 variables. I want to estimate the time Mathematica takes to solve the system, by solving smaller sets of equations. For example, I can take three equations, apply Reduce and see the AbsoluteTiming . Then, I repeat the same thing with four, five, ten etc.. equations and I can fit the execution time of Reduce with a polynomial or exponential function. However, the system is very big, only three equations store about 2MB in a text-file and indeed Reduce takes too much time when applied to only one equation (which instead should be very easy, since it should just resolve on one unknown ). Instead, Solve is pretty fast but I have some experiences that Solve does not find all the solutions, so I would like to rely on the algorithms used by Reduce . For example, you can find one equation here (I am forced to refer to external links, as the output of the system is a messy and wouldn't fit here) that you can save in a file and im...