Tag Info

Hot answers tagged

22

First off, your function is very simple without any hard number-crunching, so it will always be hard to get a large speedup for the compiled version. Secondly, your Parallelization option for Compile is useless because it doesn't do any parallelization this way. Let me give slightly changed versions of your examples and explain how you can achieve a large ...


17

Here is a summary of comments (before @ciao's best answer above), with a change in notation. These functions calculate the number of partitions of n into exactly k distinct parts of size at most m. NumberOfWays000[n_, k_, m_] := Count[Map[Length,Map[DeleteDuplicates, IntegerPartitions[n,{k},Range[m]]]], k] NumberOfWays001[n_, k_, m_] := ...


16

This seems pretty quick, particularly on larger cases / larger k, e.g. 451, 29, 101 finishes in a few seconds on the loungebook. N.B. - I have not tested this exhaustively, just thrown together from ideas... If[Min[#3, #1 - Tr@Range@(#2 - 1)] < 0, 0, SeriesCoefficient[QPochhammer[-x y, x, Min[#3, #1 - Tr@Range@(#2 - 1)]], {x , ...


15

Here is a totally different approach based on the fact that successive products forming the generating function are due to multiplication by a binomial $1+t*z^j$. Form a matrix $v$ of zeros with $n+1$ rows and $k+1$ columns. Initialize the top left corner to 1. Iterate $v=v+w$ where $w$ is the matrix $v$ shifted down by $j$ rows and to the right by 1. The ...


9

Chunks of weak compositions Here is slightly modified version of algorithm used in Combinatorica`NextComposition converted to a LibraryFunction. Needs["CCompilerDriver`"] " #include \"WolframLibrary.h\" DLLEXPORT mint WolframLibrary_getVersion() { return WolframLibraryVersion; } DLLEXPORT int WolframLibrary_initialize(WolframLibraryData libData) { ...


9

The problem is that it reevaluates the sum every single time you call it, recomputing every 20^4 term again and again. You just need to compile the function CC2 so that it performs the summation only once. Using the code you have, it takes my machine about 6 seconds to compute a single data point: CC2[0.003] // AbsoluteTiming (* {6.069311, 1.49893} *) ...


7

This is faster for the given array: keepFirstK7[lst_List] := Delete[lst, Flatten[(GatherBy[Position[lst, #], Last] & /@ DeleteDuplicates[Flatten@lst])[[All, All, 2 ;;]], 2]] keepFirstK7@rArray $\ $ {{4, 2, 1, 0}, {2}, {0, 0, 3, 3, 3}, {}, {3, 4}} Another approach keepFirstK[lst_List] := Block[{sowFirst}, sowFirst[int_, count_] := ...


7

Another solution using memoization: keepFirst[lst_List] := Module[{isFirst}, isFirst[elem_, col_] := (isFirst[elem, col] = False; True); Pick[lst, MapIndexed[isFirst[#1, Last[#2]] &, lst, {2}]] ]


7

Well, to be honest, despite I've been using Mathematica for 3 years, I'm getting more and more confused about what's functional programming, but the following solution is at least more elegant and faster than yours: searchSpan2[knots_, u0_] := First@Ordering[UnitStep[u0 - knots], 1] - 1 NonzeroBasis2[p_, u_, u0_] := With[ {i = searchSpan2[u, u0], ...


6

To speed up the CAGDBezierSurface,I have two trials 1, Refactor the pure function AllBasis I don't know why the most efficient method to calculate the Benstein basis is the defintion of Benstein basis. AllBasisNew = Function[{deg, u0}, Bernstein[deg, #, u0] & /@ Range[0, deg]] AllBasisOld = Function[{deg, u0}, ...


5

You're after the Swinnerton-Dyer Polynomials. Take a look and compare with OEIS (which, BTW cites our friend Roman E. Maeder. Programming in Mathematica, Addison-Wesley, 1990, page 105): MinimalPolynomial[Sum[Sqrt[Prime[i]], {i, #}], x] & /@ Range@5 gives the same results shown in the OEIS page, but MinimalPolynomial can't calculate sixth term within ...


5

Compile f and use a memo-ized version of it Since it seems like NIntegrate decides to symbolically evaluate its argument first, I thought I'd force it not to by compiling the function f. This seems to make a significant difference: Clear[f, f1, g] g[x_] = Nest[f[x] + 1./# &, f[x], 500]; f1 = Compile[{x}, Sum[1/100 Erfc[-(x^2/k)], {k, 100}]]; ...


5

Let's assume you have the following points a-f, which form the line segments ab and cd and ef. a = {0, 0}; b = {3, 3}; c = {1, .5}; d = {2, 1.5}; e = {0, 2}; f = {2, 0}; In order to test whether they intersect you can use RegionIntersection. First form lines from your points, then use RegionIntersection on the two lines. Compare the intersection to an ...


5

As this is a special-functions question, I feel justified in using a bit of heavy artillery. Here goes nothing... In effect, what the OP seems to want to do is to evaluate $$\sum_{n=1}^\infty \frac{(q^{n+1};q)_\infty}{(q^n;q)_\infty} q^{n-1}$$ (where $(a;q)_n$ is the $q$-Pochhammer symbol) by approximating it with its partial sums. However, there is a ...


5

For the auxiliary fuction searchSpan[] , which came form the following algorithm of The NURBS Book. where $U=\{\underbrace {a,\cdots ,a}_{p+1},u_{p+1},\cdots u_{m-p-1},\underbrace {b,\cdots,b}_{p+1}\}, \quad n=m-p-1$ Here is a rule-based solution that I implemented according to The Toad's answer NonzeroBasis[{deg_, knots_}, u0_] := Module[{coeff, i, ...


5

I found this question quite interesting, so I thought I would collect the answers contributed in comments for future reference and to have the question appear as answered in search. I generated a slightly bigger matrix to play with, and minimally modified the code to render it independent of the size of the matrix. I also compared timings of each method to ...


4

Evaluate the integral once and for all (cf. cdfc): cdfc[k_] = Integrate[PDF[NormalDistribution[0, 1], y], {y, k, Infinity}]; TCJS[T_, k_] := A/T + c1[T]*d1*T + h1*((d1*T)/2 + k*σ1*Sqrt[T + 1]) + (b1/T)*σ1* Sqrt[T + L1]*(PDF[NormalDistribution[0, 1], k] - k*cdfc[k]); EQ1[T_] := (k*σ1*h1)/(2*Sqrt[T + L1]) - ((b1*σ1)/ ...


4

For a non-symmetric real matrix you can consider using LibraryLink to speed things up. It still won't be as fast as the Total/Tr answer, but it may be useful otherwise (call this C program SumUpperTriangle.c): #include "WolframLibrary.h" DLLEXPORT int SumUpperTriangle(WolframLibraryData libData, mint Argc, MArgument *Args, MArgument Res) { /* Variable ...


4

NSolve with adequate precision works well $Version (* "10.2.0 for Mac OS X x86 (64-bit) (July 7, 2015)" *) eqns = Sin[z + Sin[z + Sin[z]]] == Cos[z + Cos[z + Cos[z]]] && -3 < Re[z] < 3 && -3 < Im[z] < 3; roots = NSolve[eqns, z, WorkingPrecision -> 20]; And @@ (eqns /. roots) (* True *) Note that there are a large ...


4

According to the solution of xzczd that dealing with the calculation of $$N_{i-p,p}(u_0),B_{i-p+1,p}(u_0), \cdots, N_{i,p}(u_0)$$ I mimic this strategy to calculate all the values of Berstern basis function of degree $n$ $$\color{blue}{B_{n,0}(u_0),B_{n,1}(u_0), \cdots, B{n,n}(u_0)}$$ AllBernsteinBasis[n_, u0_] := Nest[MovingAverage[ArrayPad[#, 1], {u0, ...


3

Here's a single-shot adaptation of my earlier (deleted since question changed) ideas. A bit over twice as fast as those for single-shot cases only, e.g. on c=5, 5000 sublists, 30 sublist length, p=0.001 done once. mutateList5[list_, c_, p_] := Module[{fl = Flatten@list, cr = Range@c, cl, lfl = Times @@ Dimensions@list, chg}, cl = Delete[cr, #] & /@ ...


3

This expression, designated exp for convenience, can be simplified substantially as follows. num = Map[FullSimplify[#] &, Numerator[exp]]; Map[FullSimplify[#, θ[_] ∈ Reals && ϕ[_] ∈ Reals] &, Denominator[exp] /. Abs[z_]^2 :> FullSimplify[Abs[z]^2, θ[_] ∈ Reals && ϕ[_] ∈ Reals]; /. Abs[z_]^2 :> z^2] den = ...


2

You can define your tensor contraction routine using the builtins Dot and Transpose. Here is an example: DotAt[T_?TensorQ, U_?TensorQ, m_Integer?Positive, n_Integer?Positive] := With[{dimT = ArrayDepth@T, dimU = ArrayDepth@U}, Dot[Transpose[T, Insert[Range[dimT - 1], dimT, m]], Transpose[U, Insert[Range[2, dimU], 1, n]]]] DotAt[T, ...


2

Uninspiring but thought I'd play: func[lst_] := Module[{res = Thread[{#, Range[Length@#]}] & /@ lst, f = Unique[]}, (f[#] = 1) & /@ Flatten[res, 1]; Cases[#, {x_, 1} :> x] & /@ Map[{#[[1]], f[#]++} &, res, {2}]] so func[rArray] yields: {{4, 2, 1, 0}, {2}, {0, 0, 3, 3, 3}, {}, {3, 4}}


2

This might not be entirely correct or complete but it gives an idea of a recursive descent approach. Clear[coeff]; coeff[dp_Plus, mon_, vars_] := coeff[#, mon, vars] & /@ dp coeff[dp_, Dot[], vars_] /; FreeQ[dp, Dot] && FreeQ[dp, vars] := dp coeff[a_.*dp_, mon_, vars_] /; MatchQ[dp, mon] && FreeQ[a, vars] := a coeff[dp_, _, vars_] /; ...


2

Here you have a way without using neither UnitStep nor PieceWise that improves the performance by 75% wrt your code. It computes 200 functions in a very reasonable time for your toy example. The main trick is to use a numeric (black box) function to be able to take Part[... ] inside it. n = 200; af = Array[f, n]; taf[t_] := Through[af[t]] bb[i_?IntegerQ, ...


2

In the narrowly defined interpretation of the question, one can get an incremental speed gain by flattening the lists, and reshaping the result to the desired dimensions afterwards: ChangeV4[prob_, c_, sizelist_, sizesublist_][list_] := Block[{elemstochange = RandomVariate[BernoulliDistribution[prob], sizelist sizesublist], temp = ...


2

Here's a version allowing MaxFilter to work with windows of even length k. It runs MaxFilter with window radius k/2-1, then corrects the output. MovingMaxEven is slower than Andrew's MovingMax, for example, 1.01 s versus his 0.78 s on 10 million points. MovingMaxEven[s_List, k_?EvenQ] := Block[{r = k/2 - 1, f}, f = MaxFilter[s, r][[r + 1 ;; -r - ...


1

In V10, Solve works, too, and gives 1251 solutions. roots2 = Solve[eqns, z]; // AbsoluteTiming Length@roots2 Solve::incs: Warning: Solve was unable to prove that the solution set found is complete. >> (* {99.1951, Null} 1251 *) Maybe there are more, too. The timing is almost 6 times as long as BobHanlon's NSolve command on my computer. But ...


1

I'm very late here, but why not? First I steal the list of happy numbers below 1000 from wikipedia: happy1000 = {1, 7, 10, 13, 19, 23, 28, 31, 32, 44, 49, 68, 70, 79, 82, 86, 91, 94, 97, 100, 103, 109, 129, 130, 133, 139, 167, 176, 188, 190, 192, 193, 203, 208, 219, 226, 230, 236, 239, 262, 263, 280, 291, 293, 301, 302, 310, 313, 319, 320, 326, ...



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