# Tag Info

24

Usage Just use this function with any polyhedron in in form: GraphicsComplex[pts_, Polygon[vertices_, ___]]. When I find time and motivation maybe I will add more DownValues so it can be more general. At the moment you can play with solids given by PolyhedronData[... "Faces"]: polyhedronRandomWalk[ PolyhedronData["DuerersSolid", "Faces"] ] ...

23

BernoulliDistribution is a perfect fit for this. RandomVariate[BernoulliDistribution[1 - 0.1], {50}] {1, 1, 1, 0, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1} Also, as kguler states, you can use RandomChoice, but the benefit of BernoulliDistribution is that ...

12

ClearAll[dist] dist [rho_] := CopulaDistribution[{"Binormal", rho}, {UniformDistribution[{0, 1}], UniformDistribution[{0, 1}]}]; data1 = RandomVariate[dist[-.9], 5000]; ListPlot[data1, Frame -> True, AspectRatio -> 1] data2 = RandomVariate[dist[.6], 5000]; ListPlot[data2, Frame -> True, AspectRatio -> 1]

10

In Mathematica it is natural to approach such a task with list operations and pattern matching. dice1 = RandomInteger[{1, 6}, 100]; dice2 = RandomInteger[{1, 6}, 100]; Count[dice1 + dice2, 2 | 3 | 4 | 5 | 6] You seem to be a very new beginner, since you are using x[i] and y[i] as if these are vectors, when they are in fact not, in Mathematica. Mathematica ...

9

There are many ways you can do this, e.g. ri = RandomInteger[{1, 6}, {100, 2}]; SortBy[Normal@GroupBy[ri, Total, Length@#/100. &], First] yielding: {2 -> 0.01, 3 -> 0.05, 4 -> 0.11, 5 -> 0.08, 6 -> 0.13, 7 -> 0.12, 8 -> 0.17, 9 -> 0.14, 10 -> 0.11, 11 -> 0.07, 12 -> 0.01} rules linking sum to frequency. You can also exploit ...

7

While the answers so far have covered a lot of ground already I have not seen EmpiricalDistribution. I would like to build upon this observation by providing a couple of general considerations that I have found to be useful when doing statistical experiments using Mathematica. What users of Mathematica may take for granted may surprise newcomers: You can ...

5

First, you're not using a fixed step method. (An Euler scheme may be applied to any step size and to one that varies.) To get a true fixed step method you have to turn off "DiscontinuityProcessing" when you have a discontinuous ODE; otherwise, NDSolve will try to adapt the steps to account for the discontinuity. The "DiscontinuityProcessing" stage resets ...

5

For loop is always slow. You may try this: f = Compile[{{x, _Real}, {y, _Real}}, If[y >= 12. Cos[x] && y >= 10 + x^3, 1, 0]]; vol[n_Integer /; n <= 10^6] := 3.* Total[f@@@Transpose@{RandomReal[{0, 1}, n], RandomReal[{10, 13}, n]}]/n; The calculation of 1000000 samples takes 1.1 s on my i5-3210M.

5

Let me give a model solution which can easily be adapted. 2 dimensions Consider a photon which moves in the x-y-plane, starting at time t = 0 in the origin {0,0} and moving towards the positive x-axis. At each tick of the clock, corresponding to a constant distance 1 travelled by the photon, the photon will experience a scattering event which leads to a ...

5

My interpretation is: the photon will be scattered by an angle $\alpha$ (given by getScatterAngle), and the deviation will occur with equal probability in every direction. (For example, a photon initially going along the $z$ axis will be rotated by an angle $\alpha$ about a randomly chosen axis that lies in the $x,y$ plane.) When I've written Monte Carlo ...

5

Take a uniform random distribution and check if it is above some threshold (0.9 in your case). For example: dist[] := If[RandomReal[] > 0.9, 0, 1]; Table[dist[], {i, 100}]

4

RandomChoice[{.9, .1} -> {1, 0}, 10] (* {0, 1, 1, 1, 1, 1, 1, 0, 1, 1} *) Timing results Timing[RandomVariate[BernoulliDistribution[.9], {10^8}];] (* {3.38014, Null} *) Timing[RandomChoice[{.9, .1} -> {1, 0}, 10^8];](* {5.64937, Null} *) dist[] := If[RandomReal[] > 0.9, 0, 1]; Timing[Table[dist[], {i, 10^8}];] (* {13.9842, Null} *) One nice ...

4

The reason it's not finishing is that you set MaxSteps -> Infinity with a high AccuracyGoal. We can actually solve this system analytically, by replacing NDSolve with DSolve: qq = DSolve[{q'[t] == v[t]/r - q[t]/r*(1/c), q1'[t] == (((q[t] + q1[t])/(2*a*ϵ0*k2))*δ), q[0] == 10*10^-9, q1[0] == 10*10^-9}, {q, q1}, t]; We can inspect the solution: ...

4

I misunderstood and commented only about computing the Cesaro means as per the question In s I have all partial sums, but I do not know how to divide them by corresponding n. The desired scatter plot of 500 Monte Carlo attempts with samples of increasing length could be obtained with something like ticks = Range[-0.06, 0.06, 0.02]; s = Table[u = ...

4

In your example plot the Monte-Carlo integral is computed afresh for each new amount of sampling points: s[n_] := Total[Sqrt[1 - #^2] & /@ RandomReal[1, n]]/n Then take one random point, find the mean, take two new random points, find their mean, and so on until 500. The result is: ListPlot[Table[s[n], {n, 500}], PlotRange -> {0.7, 0.9}] The ...

4

You can use RandomVariate to sample from a DiscreteUniformDistribution and then add up the pairs, calculate the probability of the sums observed, and then extract the probabilities of interest. (#/100. & /@ Counts[Plus @@@ RandomVariate[ DiscreteUniformDistribution[{1, 6}], {100, 2}]] )[#] & /@ Range[2, 6] Hope this helps.

4

dist = TransformedDistribution[x + y, { Distributed[x, DiscreteUniformDistribution[{1, 6}]], Distributed[y, DiscreteUniformDistribution[{1, 6}]]}]; SeedRandom[1] For small sample sizes, the match to the theoretical values is poor. data = Total /@ RandomInteger[{1, 6}, {100, 2}]; Show[ Histogram[data, {1.5, 12.5, 1}, "PDF", ...

3

nextGen[n_] := Total@RandomChoice[{11/32, 3/8, 3/16, 3/32} -> {0, 1, 2, 3}, n] simulate[n0_, nrOfGenerations_] := Total@NestList[nextGen, n0, nrOfGenerations] Now we can simulate six generations a hundred times and compute the mean value. The initial number of organisms is 10 in this example. Table[simulate[10, 6], {100}] // Mean // N (* Out: 75.42 *) ...

3

data10 = RandomVariate[f[13, 0.5], {10, 25}]; (* 10 data sets from distribution f *) lls= LogLikelihood[EstimatedDistribution[#, f[a, b], {{a, 1}, {b, 1}}], #] & /@data10 (*{-32.4994, -25.2268, -21.9671, -26.8963, -25.9164, -22.8958, -26.5247, -24.9622, -33.9319, -28.6512}*) maxll=Block[{k=1}, MaximalBy[Last][ With[{dist = ...

3

I'll try it as an answer though it's short... getTrajectory[startV_, steps_] := Accumulate@NestWhileList[stepVector, startV, RandomReal[] >= 1/steps &];

3

You can also use RandomVariate with DiscreteUniformDistribution: rW[a_, b_, n_] := Accumulate[Prepend[ RandomVariate[DiscreteUniformDistribution[{{-a, a}, {-b, b}}], n]], {0,0}] dt = rW[10, 20, 100]; Graphics[{PointSize[Large], Red, Point@#, Thick, Blue, Line@#} &@dt, Frame -> True, Axes->True, AspectRatio -> 1/GoldenRatio] We get the ...

3

f[r_?NumberQ,n_Integer]:={First[#],#.{r,Sqrt[1-r^2]}}&/@RandomReal[NormalDistribution[0,1],{n,2}]; Produces n pairs of numbers with the correlation r.

2

Maybe you can try a matrix approach. 1/ The idea is to generate a matrix like this one : mat = {{a, 1, 0}, {b, 0, 1}, {c, 0, 0}}; then you can see that: {x1, x2, x3}.{{a, 1, 0}, {b, 0, 1}, {c, 0, 0}} {a x1 + b x2 + c x3, x1, x2} gives you the format you want. 2/ Then you can use that matrix directly in NestList (without the need to define a ...

2

This is not a direct answer to your problem, but rather a generalization of @Jens' code from double to n-tuple pendulums. Meaning you can also use it for double pendulums if you like. I'm providing it due to popular demand. Needs["VariationalMethods`"] n = 10; (* number of pendulum segments *) rate = 10; (* animation frame rate *) Clear[s, ϕ, t, g, m]; ...

2

Using Event-Handling and Euler method Expanding the solution given by Michael E2 one might figure that using EventSeries is coming closest to the "real thing" which after all is a series of discrete events: eventTimes = Range[ 0, 10, 1 ]; (* example *) isEventTimeQ[ t_?NumericQ ] := Piecewise[ Table[ { 1 , t == eventTime }, {eventTime, eventTimes}], ...

2

What you are missing, I believe, is sufficient experience of Mathematica's core language at the functional level that you experimenting with. I give you credit for making a good try at formulating your code in a functional way, but I'm afraid you gone somewhat wide of the mark. I have put your code into a form that works and which I think preserves your ...

2

To help get you started I attempted to translate equation (5) from your snapshot into Mathematica. Step 1 I(H) from the snapshot becomes: i[h_] := Piecewise[{ {Gamma[1 - 2 h]/h Sin[π/2 (1 - 2 h)], 0 < h < 1/2}, {Gamma[2 (1 - h)]/(h (2 h - 1)) Sin[π/2 (2 h - 1)], 1/2 < h < 1}, {π, h == 1/2} }] This is a function that can be ...

2

The Finite Element solver in Mathematica does run in parallel, both element computation and the linear solve process are spread over the CPU cores available. Additionally, the option "MeshElementBlocks" for ToElementMesh splits the mesh elements in blocks which could be used for a domain decomposition. To get a more detailed answer you'd need to clarify ...

2

This is not an answer, but an extended comment on @beliarius's answer which I much admire. Notwithstanding my admiration, I have some nits to pick. The simulation should be parameterized by the number of steps to run, the number of dice, and the delay between histograms. oneCycle should not flatten news twice. For a million steps, as proposed by the OP, ...

2

(* First we calc the ways to add up n with two dice *) alts[x_] := Union[Join @@ Map[{#, Reverse@#} &, IntegerPartitions[x, {2}, Range@6]]] ways = alts /@ Range[1, 12]; (* an initial dice config *) m1 = RandomVariate[DiscreteUniformDistribution[{1, 6}], 300]; (* the "collision result" function*) oneCycle[m1_] := Module[{p1, s1, news}, (* form the ...

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