# Tag Info

2

If you only want the pixel coordinates of a line, that can be done much simpler: Clear[pointsOnLine] pointsOnLine[{p1_, p2_}, nPts_] := Array[Round[p1 + # (p2 - p1)] &, nPts, {0., 1.}] pointsOnLine[l_: {p1_, p2_}] := pointsOnLine[l, Round[Max[Abs[p1 - p2]]] + 1] The result is the same as your algorithm: but it's much faster: p1 = {2, 3}; p2 = ...

3

Here's a C++ implementation using LTemplate. I'm using LTemplate because it made it easy enough to write the code that I didn't give up before starting ;-) << LTemplate SetDirectory[$TemporaryDirectory]; (* currently LTemplate writes and reads files to/from the current directory *) code = " #include <cmath> struct Binner { ... 0 I have found one possible way to include all the reflections in one formula: Clear[bresenhamSolve]; bresenhamSolve[p1 : {x1_, y1_}, p2 : {x2_, y2_}] := Block[{dp = Cross[p2 - p1], ddp = Subtract @@ Abs[p2 - p1], x, y}, {x, y} /. Solve[{dp.{x, y} + err Max[Abs[dp]] (2 UnitStep[x y ddp] - 1) == dp.p1, Abs[err] < 1/2 || err == 1/2, {x, y} ... 0 If you want to use a particular maximum five times, there is no need to calculate it five times. c = 1 Manipulate[ max1 = x /. Last[FindMaximum[{(0.8 + x/(0.1 + (a*x))) - (0.8 + (c*x)), 0 <= x <= 1.6}, {x, 0}]]; max2 = x /. Last[FindMaximum[{(0.8 + x/(0.1 + (b*x))) - (0.8 + (c*x)), 0 <= x <= 1.6}, {x, 0}]]; ... 3 Per your example case: (* setup (you'll pay time here, but once done, it's done) *) dist = ParameterMixtureDistribution[PoissonDistribution[m], m \[Distributed] NormalDistribution[1/2, 1/10]]; pdf = PDF[dist, z]; (* pick some end-point with negligible probability *) (* Here, using 30 truncates tail with total p of ... 2 As mentioned in the comments, there is an implementation of Bresenham's line algorithm by halirutan: bresenham[p0_, p1_] := Module[{dx, dy, sx, sy, err, newp}, {dx, dy} = Abs[p1 - p0]; {sx, sy} = Sign[p1 - p0]; err = dx - dy; newp[{x_, y_}] := With[{e2 = 2 err}, {If[e2 > -dy, err -= dy; x + sx, x], If[e2 < dx, err += dx; y + sy, ... 4 Probably the best approach to the problem would be to implement an algorithm which allows to generate pixel positions in the original image along the Line without using the FrontEnd I think ImageTransformation doest just that. Using your definitions: r = Min[id]/2; dir = N@AngleVector[α]; (*modify step distance so we step at least 1 pixel in x or y \ ... 5 As suggested by LLlAMnYP in a comment to this question, this is a humble contribution. The OP has already been answered. This is not answer per se but shows that CompileLength should not always be increased, and should even sometimes be reduced for significant speed gain. Consider the following (stupid) function: x1 = Function[{n, T, t}, (Table[Cos[(Mod[t, ... 6 I just tried to clean up the code a bit. With respect to sortVecs, I used Ordering which is really what you were going for. You were wasting a little bit of time by taking the conjugate of one of the vectors, when the eigenvectors in vecls are all real-valued. Should you move to a different form for H that gives complex eigenvectors, just uncomment the ... 4 Using IndexGraph: g = GridGraph[{250, 250}]; a = DeveloperToPackedArray[ List @@@ EdgeList[IndexGraph[g]]]; // AbsoluteTiming {0.063857, Null} Using AdjacencyMatrix: b = UpperTriangularize[AdjacencyMatrix[g]][ "NonzeroPositions"]; // AbsoluteTiming {0.002584, Null} c = idxEdgeList[g]; // AbsoluteTiming {0.276563, Null} Test ... 2 An order of 10! I can't get you there, but I think I can add some average RT improvement. You have some up-front time costs like fetching values from the research servers- we'll ignore these, as it can't be improved without caching / downloading. Where I see enormous gains are using your simple operations (Fast) before you delve into the expensive ... 0 Pre-process everything into an enormous adjacency matrix data structure, or possibly a graph. It would take a very long time now, but be blazing fast later. 8 I also got interested in this problem and solved it using Quantile Regression. See my blog post Finding local extrema in noisy data using Quantile Regression . The proposed Quantile Regression algorithm is a version of the polynomial fitting solution proposed by Leonid Shifrin above, and has the following advantages: (i) it requires less parameter tweaking, ... 2 fun[lst_, part_] := Flatten[Last@Reap[MapIndexed[Sow[#2[[1]], #1] &, part], lst, #2 &]] I make no comment re: efficiency and have upvoted Leonid Shifrin and am happy to delete if assessed as non-contributory. 12 This seems to be about twice faster on large generic lists, and I made it somewhat faster still on integer lists (about 4-5x faster, per my tests): ClearAll[partitionToVectorLS]; partitionToVectorLS[list : {__Integer}, partitions_, sparsenessThreshold_: 10] := Module[{max = Max[list], min = Min[list], copy = partitions, sparseness, inds, ... 4 one more.. this I think fares well for very large n. result = NestWhile[ Nest[ Complement[#, Rest@Nearest[ # , RandomChoice[#] , { Infinity, .05}]] & , #, Ceiling[(Length@#)/100] ] &, pts, Min[EuclideanDistance @@@ Nearest[#, #, 2]] < .05 & ]; Kind of ugly to double Nest but the convergence test ... 5 The following "solution" has the benefits of: making a very a uniform grid. being fast. It has the (perhaps mortal) drawbacks of: not being automated. being pretty liberal about kicking out points. Nonetheless, I wanted to play a little. Here's my take: generate a square grid of points and use Nearest to pick out the points nearest to the gridpoints: ... 14 pts = Partition[RandomReal[1, 10000], 2]; ListPlot[pts] Use SameTest option with Union pts2 = Union[pts, SameTest -> (Norm[#1 - #2] < 0.05 &)]; Length[pts2] 326 ListPlot[pts2] 15 The following is a much faster, but not optimal, recursive solution: pts = RandomReal[1, {10000, 2}]; f = Nearest[pts]; k[{}, r_] := r k[ptsaux_, r_: {}] := Module[{x = RandomChoice[ptsaux]}, k[Complement[ptsaux, f[x, {Infinity, .05}]], Append[r, x]]] ListPlot@k[pts] Some timings show this is two orders of magnitude faster ... 7 You need to get rid of some extraneous List wrappers A = {{1, 1, 1}, {1, -1, 2}}; B = {4, 0}; X = {x, y, z}; and then thread the two sides of your matrix equation over Equals. ContourPlot3D[Evaluate @ Thread[A.X == B], {x, -10, 10}, {y, -10, 10}, {z, -10, 10}, ContourStyle -> Opacity[.6], MeshStyle -> Gray] The graphics shown above rotate in ... 2 z[x_, y_, t_] = -(1/2)* EllipticThetaPrime[1, (Pi/2)*(x + t*y), Exp[I*Pi*t]]/ EllipticTheta[1, (Pi/2)*(x + t*y), Exp[I*Pi*t]]; t1 = 1/2 + 3 I; In Plot3D use PlotStyle rather than ContourStyle. Instead of using N[_, 50] on argument of the plot, use the option WorkingPrecision. Increase the resolution with the option PlotPoints and set ... 2 How about mv = Compile[{{v, _Real, 1}, {id, _Real, 2}}, id.v, CompilationTarget -> "C" ]; v = RandomReal[{0, 1}, 500]; id = N[IdentityMatrix[500]]; res=mv[v, id] It works without trouble. 0 I think I made it work. Remembering function values is one key, but I needed another one: a Module arc[mn_, mR0_, mlen_, md1_, md2_] := Module[{n = mn, Pa, Pb, Pc, R, M, alpha, R0 = mR0, len = mlen, d1 = md1, d2 = md2}, Pa[n_] := Pa[n] = {d1*(0.5 + n), 0}; Pb[n_] := Pb[n] = Pa[n] + R[n]*{1 - Cos[alpha[n]], Sin[alpha[n]]}; Pc[n_] := Pc[n] = ... 2 The reason your code is slow is that you call PixelValue for each voxel. (You also call ImageDimensions[sh] for each voxel when you only need to compute it once, but that's not a big deal.) Instead, you can use ImageData once on your silhouette, giving you all pixel values at once. Here is my suggestion: {w, h} = ImageDimensions[silb0]; im = ... 5 GaussianRandomField is a special case. More generally, what is required is fast code for the (inverse) Spherical Harmonic Transform (SHT), which will work for any coefficients$a_{l,m}\$. SHTns is a high performance library for Spherical Harmonic Transform written in C and so should be straightforward to link in using MathLink. It would be very nice if this ...

4

A part from efficiency, I noticed that with the definition below, lmax >= 64 works: Clear[field]; field[θ_, ϕ_] := Chop@ Total[ Table[ alms[l, m] SphericalHarmonicY[l, m, θ, ϕ], {l, 0, lmax}, {m, -l, l} ] , 2 ]; nn = 4.; dat = ParallelTable[ field[θ, ϕ], {θ, 0, Pi, Pi/nn}, ...

5

Using the following to sim data of the actual size in your example: length = RandomInteger[{1500, 2000}, 30000]; value = Table[RandomReal[{}, length[[i]]], {i, 30000}]; start = RandomInteger[{1, 118000}, 30000]; This takes a few seconds on a friggin' loungebook so I'd venture under a second on real hardware: totLen = 120000; accum = ConstantArray[0., ...

1

If all you need is the sum at the end, why construct the matrix at all? result = ConstantArray[0, 20]; Do[ result[[start[[i]] - 1 + Range[Length@value[[i]]]]] += value[[i]], {i, 4}]; result For such a sparse operation such a loop might just be the most efficient thing you can do. Edit, a more mathematia-esque approach: Fold[#1 + ...

3

To generate an array of rules without using Table (cf. your comment under yohbs' answer), and to use it to construct a sparse array (see the updates for Method 2): (* Method 1 *) positions = Transpose[{ Join @@ MapThread[ConstantArray, {Range@4, length}], Join @@ MapThread[Range[#1, #1 + #2 - 1] &, {start, length}]}] rules = ...

2

First of all, as a thumb rule - if you're using loops in Mathematica, you're probably not doing it right. See here and here for more details. Second, the whole point of a sparse array is that you don't need to allocate all the memory for all the entries in advance, like you do with the command sparseA = ConstantArray[...]. Here's how it should be done: ...

29

General comments First, if you plan to use multi-dimensional integrals it is better to test with multi-dimensional integrals not with one dimensional ones. One might think that the test in the question is an appropriate one if multi-dimensional integration is done by the integrator in a recursive manner. This seems to be case for scipy.integrate.nquad (see ...

2

For the built-in function ParametricPlot3D, today I found a useful option PerformanceGoal PerformanceGoal is an option for various algorithmic and presentational functions that specifies what aspect of performance to try to optimize with Automatic settings for options. pts = Table[{i, j, (-1)^(i + j)}, {i, 6}, {j, 5}]; CAGDBezierSurface[pts, ...

8

This produces the same result ( about 5x faster ) img1 = Import["ExampleData/lena.tif"]; img2 = ColorNegate[Import["ExampleData/lena.tif"]]; ImageAdd[ ImageMultiply[img1, Image[1 - SparseArray[ ptsToReplace -> 1 , Reverse@ImageDimensions[img2]]]], ImageMultiply[img2, Image[SparseArray[ ptsToReplace -> 1 , ...

5

you can define a function that downsamples in one dimension: downsampleX = Total /@ Partition[#, 10] &; Then call that once on your array and then on each row: downsampleX /@ downsampleX[sbt]; Takes about 0.19 s on my PC

6

You could compile the function. For example, imageBinC = Compile[{{imagedata, _Integer, 2}, {n, _Integer}}, With[{dims = Reverse@Dimensions[imagedata]}, Table[ Total[Flatten[ imagedata[[row ;; row + n - 1, col ;; col + n - 1]]]], {row, 1, dims[[2]] - n + 1, n}, {col, 1, dims[[1]] - n + 1, n} ] ], CompilationTarget -> ...

5

The source of your problem is using Integrate with i[t, freq] in the expression. It doesn't integrate. Using NIntegrate (as suggested by Michael E2) the problem can be solved in a reasonable amount of time. Below are two changes. ReI[freq_] := 2/tmax[freq]*NIntegrate[ReaF[t, freq], {t, 0, tmax[freq]}] ImI[freq_] := 2/tmax[freq]*NIntegrate[ImaF[t, ...

4

You can increase the performance of i considerably simply by making sure that IBV is evaluated only once. i0 = 0.0001; alpha = 0.5; E0 = 0; T = 293; z = 1; Cdl = 0.000001; iLimit = 0.0005; U0 = 0.1; bvVor = z*96485/8.314/T; U[t_, freq_] := U0*Sin[2*π*freq*t] DU[t_, freq_] := 2*freq*π*U0*Cos[2*freq*π*t] iBV[t_, freq_] := i0*(-Exp[(alpha - ...

7

The procedure posted by Xavier can be sped up nicely using some undocumented features, under the same conditions of that code (data comprises sets of positive integers): minSets3[list_] := Module[{u = DeleteDuplicates@list, ju, rl, lsl, x, y, sp, ranger}, ranger[l_, lens_] := With[{al = Accumulate@lens}, Inner[l[[# ;; #2]] &, ...

6

Follow two codes to remove subsets efficiently (courtesy of ciao). Codes minSets minSets[list_List] := Module[{u = Union@list, f, rl}, f[_] = 0; MapThread[f[#] += #2 &, {Join @@ u, Join @@ MapThread[ ConstantArray, {rl = 2^Range[0, Length@u - 1], Length /@ u}]}]; Pick[u, Subtract[BitAnd @@@ Map[f, u, {2}], rl], 0]]; minSets2 ...

2

Why don't you use NFourierSeries? << FourierSeries` f[t_] = Cos[t] + 0.2*Sin[2*t] + Piecewise[{{.1, 1 < Mod[t, 2 Pi] < 2}}, 0]; NFourierSeries[f[t], t, 3]; // AbsoluteTiming {3.39228, Null} This time is measured on my laptop!

Top 50 recent answers are included