My code: Plot3D[{(2*x*y)/(x^2 + y^2)}, {x, -2, 2}, {y, -2, 2}] Output: And now I can export this to STL by using the export STL command, however, when I try to print this on the MakerBot 3D printer there is a problem because the width of the graph is too thin. I need to increase the thickness of the width of the graph, can I do this in Mathematica? Answer Try this: Plot3D[{(2*x*y)/(x^2 + y^2)}, {x, -2, 2}, {y, -2, 2}, PlotStyle -> Thickness[1]] And remember to watch your units - if you print in millimetres, 1 is a bit small... The above no longer works in Mathematica version 10. Instead of Plot3D , use ParametricPlot3D : ParametricPlot3D[{x, y, (2 x y)/(x^2 + y^2)}, {x, -2, 2}, {y, -2, 2}, PlotStyle -> Thickness[1]]

Context In[855]:= D[Abs[x], x] /. x -> 1 Out[855]= Derivative[1][Abs][1] In[856]:= D[x, x] /. x -> 1 Out[856]= 1 Question Why is Derivative[1][Abs][x] not simplifed to 1 when x -> 1 ? Answer $\mathrm{abs}(z)$ defined on the set of complex numbers $\mathbb{C}$ is not a holomorphic function because it violates the Cauchy-Riemann conditions , and the derivative is not well defined. $\mathrm{abs}(x)$ defined on the set of real numbers $\mathbb{R}$ is differentiable everywhere except at $x=0$. Mathematica treats Abs[x] as a function defined on complex numbers and so will not simplify Abs'[x] when it occurs in expressions. However, you can force it to simplify by explicitly stating that the argument is real as follows: FullSimplify[Abs'[x], x ∈ Reals] (* Sign[x] *) Also see this question for a very closely related problem.

Using the code from this answer for a sufficiently large problem setting crashes my kernel. I'll copy/paste the code here: n = 12; m = 4; costs = Range@n/2 // N; vars = Flatten@{ Array[x, {n, m}], Array[y, m], z }; constraints = Flatten@{ Table[Sum[x[i, j], {j, m}] == 1, {i, n}], Table[y[j] == Sum[x[i, j] costs[[i]], {i, n}], {j, m}], Table[z >= y[j], {j, m}] }; bm = CoefficientArrays[Equal @@@ constraints, vars]; solution = LinearProgramming[ Last@CoefficientArrays[z, vars], bm[[2]], Transpose@{-bm[[1]], constraints[[All, 0]] /. {Equal -> 0, GreaterEqual -> 1}}, vars /. {_x -> {0, 1}, (_y | z) -> {0, \[Infinity]}}, vars /. {_x -> Integers, (_y | z) -> Reals} ]; Cases[Pick[vars, solution, 1], x[ij__] :> {ij}]; GroupBy[%, Last -> First] // Values // SortBy[First] Running this for n=40 yields a result almost instantly, whereas n=60 leads to this after a minute or so. I am using Mathematica 10.0.1 on ...

I have a list $\{x,y,z,w,p\}$ and I want to replace $x,y,w$ by $q$. This works {x,y,z,w,p}/.{x->q, y->q, p->q} But is there a way to do this concisely? I tried the following, but doesn't work {x,y,z,w,p}/.{(x || y || p) ->q} Thanks in advance for any help. Answer {x, y, z, w, p} /. Alternatives[x, y, p] -> q (* or *) {x, y, z, w, p} /. Thread[{x, y, p} -> q] {q, q, z, w, q}

One can use $Epilog to do something when the Kernel is quit or put an end.m file next to the init.m . For Wolfram System sessions, $Epilog is conventionally defined to read in a file named end.m. But if $Epilog is set by the user, then end.m is skipped. Question: So what to do if I want something to be done each time Kernel is quit but I also want to be able to play with $Epilog . In that sense that if I set $Epilog I would like my default action to be taken anyway? I need to stress out that I want to establish an action that will not be accidentally overwritten with daily (but advanced) mma usage.

I need to solve a discontinuous equation which is typical in theory of plasticity. For a simple case I get the following equation system (reformulated for numerical implementation): $$\begin{align*} s(t) &= \frac{\sigma(t)}{C_1} + s_{ep}\\ s_{ep}'(t) &= \begin{cases}\frac{C_1}{C_1+C_2}s'(t) & \text{for } |\sigma(t)-C_2 s_{ep}(t)| \ge \sigma_{gr} \land \sigma(t)s'(t)>0\\0 & \text{otherwise}\end{cases} \end{align*}$$ with "zero" initial conditions. I'd like to get the solution for $s(t)$ for given parameters $C_1$ , $C_2$ , $\sigma_{gr}$ and a known function $\sigma(t)$ . I assumed: $\sigma(t) = 40000\sin(0.02t)$ , $C_1=80000$ , $C_2 = 20000$ , $\sigma_{gr} = 15000$ . This should give a hysteresis loop on a plane $\sigma(t)-s(t)$ . So in Mathematica I tried to use automatic discontinuity handling by defining the second equation using a Piecewise function: σ[t_] := 40000*Sin[0.02*t]; eq1 = s[t] == σ[t]/C1 + sep[t]; eq2 = sep'[t] == Piece...

I have strings of the following form: string = "ABC123DEFG456HI89UZXX1"; Letter keys of variable lengths are followed by (positive) integers. I want to get this transformation: {{"ABC", 123}, {"DEFG", 456}, {"HI", 89}, {"UZXX", 1}} I have written: chars = Characters @ string; runs = Length /@ Split[IntegerQ @ ToExpression @ # & /@ chars] {3, 3, 4, 3, 2, 2, 4, 1} takes = Transpose[{# - runs + 1, #}] & [Accumulate @ runs] {{1, 3}, {4, 6}, {7, 10}, {11, 13}, {14, 15}, {16, 17}, {18, 21}, {22, 22}} result = Partition[StringJoin /@ Map[Take[chars, #] &, takes], 2] /. {a_String, b_String} :> {a, ToExpression @ b} {{"ABC", 123}, {"DEFG", 456}, {"HI", 89}, {"UZXX", 1}} I have two questions: (1) How could the above coding be shortened / improved? (it seems to be too long for such a trivial problem) (2) How would a "direct" method ( StringCases , StringSplit ...) loo...

I really miss having something like a struct in Mathematica. I know of (and regularly use) a couple of programming techniques which feel like a struct (e.g., using downvalues ), but are ultimately unsatisfactory (perhaps I'm using downvalues incorrectly). What programming approaches are available which provide similar functionality to a struct ? Here's an abbreviated (and hopefully not too obtuse) example of how I use downvalues to emulate a struct. In this case, I'm distinguishing between TLC and TEC (these are sets of parameters for two different phases of a Moon mission, trans-lunar cruise and trans-earth cruise): deadBandWidth[X][TLC] ^= 10. \[Degree]; deadBandWidth[Y][TLC] ^= 10. \[Degree]; deadBandWidth[Z][TLC] ^= 20. \[Degree]; sunSearchAngle[Z][TLC] ^= 230. \[Degree]; sunSearchRate[Z][TLC] ^= 1. \[Degree]/Second; sunSearchAngle[X][TLC] ^= 75. \[Degree]; sunSearchRate[X][TLC] ^= 1. \[Degree]/Second; safingSpinRate[TLC] ^= (360. \[Degree])/Day; sunVector[TLC] ^= {-C...

I'm using the words "Feynman diagrams" for indexing purposes, but the question is in fact purely graph-theoretic. Given a list n={n1,n2,...} of non-negative integers, I want a function that generates all graphs with n1 1-valent vertices, n2 2-valent vertices, etc., with their corresponding symmetry factors. (These are, by definition, the number of automorphisms of the corresponding graph). I'm aware of the existence of some packages that generate such graphs, but I was looking for something that doesn't require any installation. Just to copy+paste code. This is a self-answered question, but any alternative answer is obviously welcome and appreciated. Note: by graph I mean pseudograph , i.e., multiple edges and self-loops are allowed. Note also: the $j$-th element in the list n is the cardinality of $j$ in the degree sequence of the graphs. I'm not sure if there is an accepted name for the list n , although it would be nice to know. This list is important...

I want to be able to define n shapes and find the intersection between them and a line. I'm treating the shapes as "one" shape, i.e. I do not want any internal intersections (hence the volume equations). I can do this by typing out each shape and equation as follows: line = InfiniteLine[{{0.5, 0.5, 0.5}, {0.5, 0.5, 1}}]; (*shapes - I can define these*) shape1 = Ball[]; shape2 = Cone[{{0, 0, 0}, {1, 1, 1}}, 1/2]; (*------------ I want from here automated ---------*) (*surface equations *) surfaceequations1 = RegionMember[RegionBoundary[shape1], {x, y, z}]; surfaceequations2 = RegionMember[RegionBoundary[shape2], {x, y, z}]; volumeequation1 = RegionMember[shape1, {x, y, z}]; volumeequation2 = RegionMember[shape2, {x, y, z}]; intersection = NSolve[{x, y, z} \[Element] line && (surfaceequations1 || surfaceequations2) && ! (volumeequation1 && volumeequation2), {x, y, z}]; (* ---------- automation can stop here ---------- *) points = Poin...

I am in a resizing mood these days. This code works quite well, up to a point. Grid[Partition[ Table[Pane[Style[FromCharacterCode[96 + i], FontSize -> 150], ImageSizeAction -> "ResizeToFit", ImageSize -> {Scaled[1], Scaled[1]}, Alignment -> {Center, Center}], {i, 1, 25}], 5], Dividers -> All, ItemSize -> {Scaled[.2], Scaled[.2]}] It makes a 5x5 grid of 25 letters with boxes extending and shrinking according to the Notebook's width, and the strings centered in each box. The letters also shrink and expand accordingly. That's the minimal example of what I wanted. But when I shrink the window a little too much, the result is clipped instead of shrinked: And when I enlarge the window a little too much, depending on the fontsize used, the letters stop to enlarge as well, as you could probably guess, this is not what I want. Moreover, I would like the code not to be dependent on an heuristic font size. I have tried various options but without...

I'm trying to diagonalize a certain 2 x 2 matrix, but Mathematica refuses to find the eigenvectors. In particular, when I input Eigenvectors[{{1.8741*10^7 + 1.40161*10^6 B, 2.79374*10^7}, {2.79374*10^7, -3.1235*10^7 - 1.40161*10^6 B}}] (B is a parameter). I receive the error Eigenvectors::eivec0: Unable to find all eigenvectors and the bogus output {{0, 0}, {0, 0}} On the other hand, Mathematica has no problem with the eigenvalues, for some reason. Can someone help me figure out what's going on here? Answer Although the documentation could be clearer on this point, Eigenvectors doesn't like to work symbolically on a matrix containing elements with approximate numbers. The solution is to Rationalize the matrix. data = {{1.8741*10^7 + 1.40161*10^6 B, 2.79374*10^7}, {2.79374*10^7, -3.1235*10^7 - 1.40161*10^6 B}}; Eigenvectors[Rationalize[data]] // Column

EDIT: The answer posted by @Jens works for built-in functions only. However I'm starting to wonder if this is related to my Mathematica install somehow. I have version 7. Here is an example in the Mathematica documentation center that implies that if I ask InverseFunction[] to find the inverse of a non-built-in pure function, it will be able to. However, as a test I tried to duplicate the result: The documentation center gives In:= InverseFunction[(a # + b)/(c # + d) &] Out:= (-b+d#1)/(a-c#1) & . Using the exact same input, I get Out:= InverseFunction[(a # + b)/(c # + d) &] , which is just the inputted command. :\ The documentation center is talking about version 9. Is there some fundamental difference in how versions 7 and 9 use InverseFunction ? Is there any workaround? I've also tried to use other definitions other than InverseFunction[] , such as Solve[] ...but in that case Solve[] uses inverse functions to find the solution (as evinced by the message I get: ...

I am a new user of Mathematica and have some questions about the simplifications of calculated expressions. I am unable to attach an image of the session, but my Mathematica commands are: Element[{x,y,z},Reals] Element[{x0,y0,z0},Reals] rhatV={x-x0,y-y0,z-z0} rhat=Norm[rhatV] In the expression for rhat, I am unable to get rid of the Abs functions, despite the Reals declarations. phi=1/rhat D[phi,x] In the evaluated derivative is there a way to have x-x0 in the numerator recognized as rhatV[[1]] and the denominator as rhat^3, such that it can be used in additional operations? Answer The formulation of the assumptions are one problem, and the fact that you're not using them is another: $Assumptions = {Element[{x, y, z}, Reals], Element[{x0, y0, z0}, Reals] }; rhatV = {x - x0, y - y0, z - z0}; rhat = Simplify[Norm[rhatV]]; phi = 1/rhat; D[phi, x] $-\frac{x-\text{x0}}{\left((x-\text{x0})^2+(y-\text{y0 })^2+(z-\text{z0})^2\right)^{3/2}}$ Here I put the assumptions in a special...

Take a nonlinear equation such as exp = (x + 3)/4 + Exp[x] + 1 + (c + x) Note that this is not a polynomial. Now, I want to extract the coefficients of this. A few are easy: Coefficient[exp, x ] (* Correctly gives 5/4*) Coefficient[exp, Exp[x] ] (* Correctly gives 1*) But how can I extract the coefficient on the constant term? I can't figure out how to write the "form" for the coefficient function to extract it. Note that treating it as a polyomial and asking for the 0th order will not work (e.g. Coefficient[exp, x,0] is not correct) Answer You can use CoefficientList : exp = (x + 3)/4 + Exp[x] + 1 + (c + x); CoefficientList[exp, {x, Exp[x]}] (* {{7/4 + c, 1}, {5/4, 0}} *) Whether that is a convenient way depends on what you want to do with it. Following extracts the parts: {{const, ExpC}, {xC, xExpC}} = %

(This is a follow-up question to the answer found here .) I would like to define a clip plane dynamically from the current ViewPoint of my Graphics3D . Based on the solution found in the link above, I would like to do something like this (where the ClipPlanes option depends on the ViewPoint variable vp ): DynamicModule[ {vp}, {vp} = Options[Graphics3D, ViewPoint][[All, 2]]; Graphics3D[{FaceForm[Red, Blue], Sphere[]}, Axes -> True, ViewPoint -> Dynamic[vp], ClipPlanes -> {{Sequence@@vp,0}}, ClipPlanesStyle -> Directive[Opacity[.2], Blue] ] ] Unfortunately this doesn't work, and I couldn't figure out how to do it. Answer I think this may be the behavior you desire. I am borrowing Kuba's modified ClipPanes specification. One should not need RawBoxes here I think but it does serve the purpose to get our Dynamic expression into the Front End box form. DynamicModule[{vp}, {vp} = Options[Graphics3D, ViewPoint][[All, 2]]; Graphics3D[{FaceFor...

In fact My problem is this $$\frac{\partial u}{\partial t}+\ sin(y)\frac{\partial u}{\partial x}=\nu(\frac{\partial^2 u}{\partial x^2} +\frac{\partial^2 u}{\partial y^2})$$ But I wanted to test the method first to the heat equation and check if the L^2 norm of the solution behaves like this $$|u|_{L^2} =(\int_{-\pi}^{\pi} \int_{-\pi}^{\pi} u^2 dx dy)^{1/2} \leq e^{-\nu t}$$ Given that $$\frac{\partial u}{\partial t}=\nu\Bigl(\frac{\partial^2u}{\partial x^2}+\frac{\partial^2u}{\partial y^2}\Bigr)$$ With the following periodic boundary conditions: $$u(-\pi,y,t)=u(\pi,y,t) \\ u(x,-\pi,t)=u(x,\pi,t) \\u_x(-\pi,y,t)=u_x(\pi,y,t)\\ u_y(x,-\pi,t)=u_y(x,\pi,t)\\ u(x,y,0)=\sin(x)$$ I have tried to solve this using fourier collocation method in mathematica And then using NDSolve to solve the system of ODe. n = 11; ν = 1; T = 100; u[x_, y_, t_] := \!\( \*UnderoverscriptBox[\(∑\), \(k = 0\), \(n - 1\)]\( \*UnderoverscriptBox[\(∑\), \(l = 0\), \(n - 1\)]\(a[k, l]\)[t]* EXP[I*k*x]*EXP[I*l*y]...

I am working with noncommuting objects and basically I am using Mathematica to sort large expressions by normal ordering and simplifying them. To this end, I first create an object called Boson (These satisfy $[a,a^\dagger]=1$): 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_: String "a"] := Boson[False, sym] BosonC[sym_: String "a"] := Boson[True, sym] Next I alias the noncommutative product with CenterDot as follows: Unprotect[NonCommutativeMultiply]; Clear[NonCommutativeMultiply, CenterDot]; CenterDot[a__] := NonCommutativeMultiply[a]; NonCommutativeMultiply /: MakeBoxes[NonCommutativeMultiply[a__], fmt_] := With[{cbox = ToBoxes[HoldForm[CenterDot[a]]]}, InterpretationBox[cbox, NonCommutativeMultiply[a]]] Protect[NonCommuta...

I wish to color a 2-dimensional lattice grid according to the value of a function at each lattice-node. More specifically, if I have 9 angles in a 3x3 array, angles={{0, π, π}, {0, 0, π/2}, {π/2, 0, 3 π/3}} then one can plot these angles on a lattice-grid by the following code in Mathematica: angles = {{0, π, π}, {0, 0, π/2}, {π/2, 0, 3 π/3}}; GraphicsGrid[ Map[Graphics[{ LightGray, Circle[{0, 0}, 1], Hue[#/(2 π), .6, .8], Thick, Arrowheads[Medium], Arrow[{{0, 0}, {Cos[#], Sin[#]}}]}] &, angles, {2}]] Here, the colouring is done using the Hue command according to the value of each angle. However, now I want to compute a function f[i,j] ; to be more specific, f[i_,j_]:=Cos[angles[[i + 1, j]] - angles[[i, j]]] + Cos[angles[[i, j+1]] - angles[[i, j]]]; with angles[[n+1,i_]]:=angles[[1,i]]; angles[[i_,n+1]]:=angles[[i,1]]; i.e., the boundary conditions. In the first code, Hue is used with the angles[[i,j]] (through # ), which is probably straightforward. But is it possible ...

As we know, for most Gaussian integrals, we can get the analytical result. Now I have many Gaussian integrals to treat, which have the following general form, Integrate[ x1^n1 x2^n2 x3^n3 ... xd^ nd Exp[-{x1, x2, x3 ... xd}.matrix.{x1, x2, x3 ... xd} + c1 x1 + c2 x2 ... cd xd], {x1, -Infinity, Infinity}, {x2, -Infinity, Infinity} ... {xd, -Infinity, Infinity}]. The matrix is a symmetrical positive definite matrix. We can get the analytical integral result for this general form. However, it is a difficult task for Mathematica to deal with the similiar expressions, which takes so much time. I have written a package for treating this problem. But, it is still not sufficient because I need to deal with lots of these kinds of integral expressions. So the computation time is important. So have you ever seen the package for doing this kind of Gaussian integral? Or can you give me some suggestions about this problem? Thank you very much! Thanks all of you. For the concept ...

The following problem is an exploration of what causes the error "Input is not an ordinary differential equation" in Mathematica as it seems to have changed from version 8 to version 9. Specifically I have found that in Mathematica 9 I have come across the issue that NDSolve[{x'[t] == 5, C'[t] == 5, x[0] == 1, C[0] == 1}, {x, C}, {t, 0, 100}] gives the error, whereas in version 8 it does not. So something is up, but my question has more to do with if this is a bug in v8 or a bug in v9, that is, I get that C is a built in symbol: so is this error a correct catching of an input mistake? If so then NDSolve[{C'[t] == 5, C[0] == 1}, C, {t, 0, 100}] Should also give an error and it does not. I really want to understand what is possibly going on, since as an ecologist I have to deal with equations that have C in them for "consumer" and I am forced to troubleshoot any issues that might come up with people not understanding the shadowing issue. Clarification It...

What is the quickest and best way to export a plot image from Mathematica to a $\LaTeX$ Editor? Is there something quicker than exporting it to an image and then including it in the $\LaTeX$ document? Answer I usually Export to hi-res PNG bitmaps for ease of use (there are a number of discussions on how best to export high-quality images on this forum. Take a peek at the right column of this page under "Related"). Personally I like notebooks that do not need any mouse-clicking or any other user interaction to produce output which makes reruns that much more comfortable and and most of all reproducible. For ease of use when dealing with many different graphics I prepare a notebook for each of these and also generate the LaTeX code for insertion in my document. The graphics get their name from the notebook so you can easily spawn different versions by renaming the notebook. plot = Plot[Sin[x], {x, 0, 10}] (*gfxname = StringTake[FileNameSplit[NotebookFileName[]][[-1]], {1, -4...

I want to flatten a series of fisheye images by remapping them to a rectinlinear projection. To achieve this, I need to be able to remap the pixels of the image using fisheye correction formulas for the x-, and y-coordinates. How can I achieve this? I have found this question and this question but I wonder how to use this for fisheye correction. So far I've tried to use ImageTransformation for this, but I can't get the function to work properly. f[pt_] := With[{s = {.5, .5}}, Module[{rd, polarcoor, ru, newcoor}, rd = Norm[pt - s]^2/Norm[s]; polarcoor = CoordinateTransform[{"Cartesian" -> "Polar", 2}, (pt - s)]; ru = 1*Tan[2*ArcSin[((polarcoor[[1]])/(2*1))]]; CoordinateTransform["Polar" -> "Cartesian", {ru, polarcoor[[2]]}] ] ] ImageTransformation[image,f] This should first translate the image coordinates to polar coordinates, then calculate the new r (ru=r undistorted), and than translate these back to car...

UPDATE Today I was having these WaveletScalogram problems again. I think this is ultimately some kind of interaction issue with Dynamic . I have narrowed it down to this: Graphics[Raster[{{0}}, ColorFunction -> "AvocadoColors"]] Run this line, and then disable Dynamic Updating ("Evaluation" menu) and run it again. On my computer it comes out blank when Dynamic Updating is disabled. If you remove the ColorFunction option, it will work fine with or without Dynamic Updating. So new question: Why does this happen? Is it a bug? I wonder if others can even reproduce this. Previous version of this question: I was trying to export WaveletScalogram plots (to image formats) in parallel, and some of the images were coming out blank. I've managed to reduce the issue to what appears to be quirky Rasterize behavior on WaveletScalogram plots specifically, in parallel computations. I'm on version 8.0.4.0 on Win7 64. With the following: ParallelTable[ Rasterize@ Wa...

Related question: How can I convert a complex number into an exponent form Mathematica insists on displaying complex number in form a+I b when a or b are not exact: Clear["Global`*"] z = 3 + 4 I; Abs[z] Exp[I Arg[z]] z = 3.0 + 4 I; Abs[z] Exp[I Arg[z]] What I'd like is a polarForm wrapper that keeps the polar form even when a or b are not exact, like this: Clear["Global`*"] (z = 3.0 + 4 I) // polarForm (z = 3 + 4 I) // polarForm In the above, polarForm is the wrapper needed. Answer Specify the display format of something using MakeBoxes , like so: MakeBoxes[polarForm[z_Complex], form_] := With[{r = Abs[z], ϕ = Arg[z]}, RowBox[{If[r == 1, Sequence @@ {}, MakeBoxes[r, form]], If[ϕ == 0, Sequence @@ {}, SuperscriptBox[MakeBoxes[E, form], RowBox[{MakeBoxes["\[ImaginaryI]", form], If[ϕ == 1, Sequence @@ {}, MakeBoxes[ϕ, form]]}]]]}]] Sqrt[5] E^(I ArcTan[2]) // N // polarForm (* 2.23607E^(I1.10715) *)

Trying to unflatten an array that was part of ragged array of sub matrices. As input we have flattened array: myflatarray={7.7056, 4.20225, 3.02775, 7.60807, 9.77169, 6.18476, 4.80993, 6.32965, 2.69882, 0.268505, 1.78048, 9.67702, 0.875699, 7.18504, 1.62811, 0.0313174, 5.66048, 1.61613, 5.71987, 0.971484, 3.87503, 2.8537, 2.68776, 3.47433, 2.81483, 4.0387, 1.74939, 3.12037, 4.18016, 2.35264, 0.853583, 1.22412, 0.382633, 3.67874, 8.68059, 8.02419, 6.9547, 7.11112}; and the original dimensions of the submatrices in the ragged array. mydimensions={{4, 5}, {3, 4}, {2, 3}}; In other words the original object looked like original={ {{7.7056, 4.20225, 3.02775, 7.60807, 9.77169}, {6.18476, 4.80993, 6.32965, 2.69882, 0.268505}, {1.78048, 9.67702, 0.875699, 7.18504, 1.62811}, {0.0313174, 5.66048, 1.61613, 5.71987, 0.971484}}, {{3.87503, 2.8537, 2.68776, 3.47433}, {2.81483, 4.0387,1.74939, 3.12037}, {4.18016, 2.35264, 0.853583, 1.22412}}, {{0.382633, 3.67...

The documentation on RuleDelayed states that lhs :> rhs represents a rule that transforms lhs to rhs , evaluating rhs only after the rule is used. Consider the following rule: a/;condition :> (condition=False; RandomReal[]) When this rule is used, it sets the variable condition to False , and therefore the lefthand side of the rule cannot be matched any more. So I expected that in the following command only the first symbol a would be replaced by a real number. condition=True; {a,a,a} /. a/;condition :> (condition=False; RandomReal[]) (* {0.420057,0.183578,0.751358} *) But all three are replaced, so the rule is used three times. What is going on here? Answer There seems to be a subtlety in the way delayed rules are used. Have a look at the following: {a,a,a} /. a/;(Print["lhs evaluated"];True) :>(Print["rhs valuated"]; RandomReal[]) (*lhs evaluated lhs evaluated lhs evaluated rhs evaluated rhs evaluated rhs evaluated *) (* {0.797753,0....