Skip to main content

Posts

Showing posts from August, 2015

development - How to simplify writing LibraryLink code?

LibraryLink is an API for extending Mathematica through C or C++. It is very fast because it gives direct access to Mathematica's packed array data structure, without even needing to make a copy of the arrays. Unfortunately, working with LibraryLink involves writing a lot of tedious boilerplate code. Mathematica 10 introduced managed library expressions, a way to have C-side data structures automatically destroyed when Mathematica no longer keeps a reference to them. I find this feature useful for almost all non-trivial LibraryLink projects I make, but unfortunately it requires even more boilerplate. For me, writing all this repetitive code (and looking up how to do it in the documentation) has been the single biggest barrier to starting any LibraryLink project. It just takes too much time. Is there an easier way? Answer I wrote a package to automatically generate all the boilerplate needed for LibraryLink and for managed library expressions based on a template that describes a...

geography - Neighboring counties (within the State of Florida) for every county in the State of Florida

When I attempt evaluating: Entity[Alachua County, Florida, United States(administrative division)] [EntityProperty["AdministrativeDivision","BorderingCounties"]] I get: Missing["UnknownEntity", {"AdministrativeDivision", {"AlachuaCounty", "Florida", "UnitedStates"}}] The "BorderingCounties" attribute does NOT work on my computer for some unknown reason. Other attributes (like "BorderingStates") work just fine. I'm attempting to find an alternative way to get the result that I'm seeking. I'm attempting to get a list of neighboring counties (within the State of Florida) for every county in the State of Florida. For example; after evaluating: counties=EntityList[US counties in Florida (administrative divisions)]; I first try to build a list of neighboring counties; surrounding Escambia county; within the State of Florida as follows: counties[[16] gives: "Escambia County, Florida...

plotting - why this call to StreamPlot returns unevaluated? (no error message, no beep)

I am really baffled by this. Call to StreamPlot does not return empty plot, nor an error message, nor a beep. It just return unevaluated. I do not think I've seen something like this before. What would the cause for this? Is this expected behavior? Normally, when a plot can not be generated, either an error is returned or empty plot. ClearAll[x, y, f]; f = (x*y - Sqrt[-1 + x^2 + y^2])/(-1 + x^2); StreamPlot[{1, f}, {x, -2, 2}, {y, -2, 2}] btw, it will plot OK when changing the -1 to 1 in the above, under the square root: So it seems due to value of the function becomes complex over some region. But normally when this happens, an empty plot is returned, right? V 12 on windows 10. Answer I think it should give an error message, perhaps. I agree with the OP's guess that problem stems from the complex values of the vector field. It may be that it's because the field is real around the boundaries and then becomes complex. Here is a workaround: StreamPlot[ Evaluate[Conditio...

Partitioned matrix operations

I have a $(2 \times 2)$ matrix whose elements are themselves $(2 \times 2)$ matrices (i.e., a partitioned matrix ), e.g: $$ M = \begin{pmatrix} A & B \\ C & D \end{pmatrix}, \quad \quad \text{with } A,B,C,D \in \text{Mat}_{2 \times 2}. $$ Problem: I want to compute some things with this matrix, in particular an expression for $M^{n}$, I have tried using MatrixPower[M, n] in Mathematica but this does not work, (the error message says that the argument at position 1 is not a non-empty square matrix). If I pretend that $A,B,C,D$ are scalars then Mathematica will automatically assume multiplicative commutativity when computing powers of $M$. Additional info : Furthermore, matrix multiplication does not work, if I try $M.M$ then I expect to obtain a result which looks like: $$ \begin{pmatrix} A.A + B.C & A.B+B.D \\ C.A+D.C &C.B+D.D \end{pmatrix}, $$ however, instead, I obtain a result in which the entries of $M.M$ are $(2 \times 2)$ matrices of $(2 \times 2)$ matrices (a...

bugs - Dynamic panel crashes Mathematica 10.0.0.0

This answer by @paw have revealed a bug that causes Mathematica v 10.0.0.0 to crash the GUI and kernel. So far the offending code is this: Needs["QuantityUnits`"] table = Keys[QuantityUnits`Private`$UnitReplacementRules]; Panel[DynamicModule[{f = ""}, Column[{Text[Style["Mathematica Unit Search:", Bold]], InputField[Dynamic[f], String, ContinuousAction -> True], Dynamic[Union@Flatten[StringCases[#, f ~~ ___] & /@ table] // TableForm]}]]] This creates a dynamic panel and the crash happens when pressing "panel format" in the suggestion bar. The bug has been confirmed so far only on M10 under Windows 7 64 bits. My skills are not enough to reduce that code to a minimum working example that reproduces the crash, nor I have available Mathematica under other OS. So I'm both reporting the bug and asking: What would be a minimum working example of the crash? Does this crash also in other platforms? The Wolfram Technical Support...

Print out all functions in Mathematica

Sometimes I only remember parts of a function name and I want find that function name quickly, and if I search in the document center, it will give too much informations related to that and difficult to find what I'm looking for. So is there a way to print out all the MMA functions in alphabetical order in a notebook, so that I can search only in the function names? Thanks. Answer There is Documentation`HelpLookup["guide/AlphabeticalListing"] If you want to tweak it yourself, you can start with SetAttributes[makeSearchable, Listable]; (* Thanks @rm-rf *) makeSearchable[s_String] := Hyperlink[s, "paclet:ref/" s] firstLetter = StringTake[#, 1] &; CreatePalette[ OpenerView@{firstLetter@First@#, Column@makeSearchable@#} & /@ SplitBy[ Pick[#, StringMatchQ[#, WordCharacter ..]] &@Names["System`*"], firstLetter] // Column, WindowTitle -> "Function list", WindowElements -> {"VerticalScrollBar"}] However...

graphics - "The function DiscretizeGraphics is not implemented for GraphicsComplex" confusion

Why does DiscretizeGraphics seems to work on one GraphicsComplex and not the other? Here is an example that works: v = {{1, 0}, {0, 1}, {-1, 0}, {0, -1}}; p1 = Graphics[GraphicsComplex[v, Polygon[{1, 2, 3, 4}]]]; DiscretizeGraphics[p1] But this does not p2 = Graphics3D[First@ParametricPlot3D[{Cos[t], Sin[u], c Sin[t]}, {u, 0, 2 Pi}, {t, 0, 2 Pi}]]; DiscretizeGraphics[p2] (*The function DiscretizeGraphics is not implemented for \ GraphicsComplex[{{0.9999999999998993`,4.487989505128125`*^-7,*) But p2 is a GraphicsComplex? Looking at FullForm[p2] Here is the FullForm for p1 Are not p1 and p2 both GraphicsComplex ? p1 is 2D and p2 is 3D , but are they not both considered GraphicsComplex ? It will good to know exactly what can and what can not be discretized. I tried to find this, but could not. All what I see are examples of usages so far. reference: http://www.wolfram.com/mathematica/new-in-10/data-and-mesh-regions/discretizing-graphics.html http://reference.wolfram.com/lan...

graphics - Generate randomly sized non-overlapping disks?

I need to generate an image of $n$ randomly sized but non-overlapping blobs in a fixed rectangular region ; for example, a square of 300 pixels. The blobs could be disks to keep things simple. The non-overlapping part is tricky; this is what I have so far: Clear @ pair; pair[n_] := Module[{pts=RandomReal[1,{n,2}]}, Image @ Rasterize[Graphics[{{PointSize@RandomReal[{0,.5}],Point[#]}&/@pts}, PlotRange->{{0,1},{0,1}},PlotRangePadding->Scaled[.1]], ImageSize->300]->n ] As you can see the ten disks are overlapping: Answer Just a quick modification of the code here , distinctDisks[n_, range_:{0, 1}, radiusRange_:{0.03, 0.15}] := Module[ {d, f, p, r}, d = {Disk[RandomReal[range, 2], RandomReal[radiusRange]]}; Do[f = RegionDistance[RegionUnion @@ d]; While[ r = RandomReal[radiusRange]; p = RandomReal[range, 2]; f[p] d = Append[d, Disk[p, r]], {n - 1}]; d] distinctDisks[25, {0, 5}, {0, 2}] // Gra...

Why does N not upgrade precision?

Precision[N[1.0, 20]] Precision[N[1, 20]] MachinePrecision 20. It would be so much more intuitive and less error prone, if Precision[N[1.0, 20]] would be 20 and not MachinePrecision . Why do I have to explicitly use N[Rationalized[1.0], 20] to upgrade the precision? If it is about the warning messages during the calculations, Mathematica could return a warning already at the Precision[N[1.0, 20]] step. Edit I am not attempting to use N for rounding. I also understand the difference between MachinePrecision and arbitrary precision, which is very well described in this answer . I want to know, why does Mathematica not consider the description of the number existing in the underling binary representation as exact if I apply N , and why I need to wrap it with Rationalize . Is it performance? Is it some deep semantic meaning of N vs SetPrecision ? Is set SetPrecision any different from N@Rationalize@ ? Answer N can only lower precision, it cannot raise precision. N[1.3`4, 10] ...

programming - Fastest way to measure Hamming distance of integers

I am looking for a fast and robust way to calculate the Hamming distance of integers. The Hamming distance of two integers is the number of matching bits in their binary representations. I expect that clever methods can easily outpace HammingDistance as it works on vectors instead of integers and on any vector not just binary. My naive bitwise method is faster than HammingDistance but I'm pretty sure that it can be further optimized. While compilation would help, it won't work on big integers ($\ge 10^{19}$). Nevertheless, I am interested in compiled solutions! max = 10^10; n = Length@IntegerDigits[max, 2]; data = RandomInteger[{0, max}, {100000, 2}]; m1 = Map[HammingDistance[IntegerDigits[First@#, 2, n], IntegerDigits[Last@#, 2, n]] &, data]; // AbsoluteTiming m2 = Map[Total@IntegerDigits[BitXor @@ #, 2] &, data]; // AbsoluteTiming m1 === m2 {0.967202, Null} {0.624001, Null} True It would be nice to work entirely on the binary representations, and I thought t...

plotting - Probability density histogram with unequal bin widths

I am confused by output of Histogram and HistogramList for probability density ( "PDF" ) when bin widths are not equal. According to this and this pages and other sources, density histograms are computed by dividing counts by bar widths and total number of observations. But Histogram obviously uses another algorithm as one can see from the following example: SeedRandom[1]; data = RandomVariate[NormalDistribution[0, 1], 15]; HistogramList[data, {{-2, 0, 1}}, "PDF"][[2]] BinCounts[data, {{-2, 0, 1}}]/((Length[data] - 3) Differences[{-2, 0, 1}]) {1/6, 2/3} {1/6, 2/3} The outputs are identical when bin counts are divided by total number of observations minus 3 in this case. Why is it? What algorithm Histogram uses for determining this difference (in other cases I got other numbers)? An addition Rod has answered the original question but there is another issue: if one gives upper bin boundary equal to the upper datapoint value then this value will be excluded fro...

options - How can I work out which functions work with SetOptions?

Not all functions seem to work with SetOptions . e.g. SetOptions[Grid, BaseStyle -> Directive[Red]]; Grid[{{"hello", "world"}}] hello world the font is not red. SetOptions[Row, BaseStyle -> Directive[Red]]; Row[{"hello", "world"}] hello world ...and the font is red. SetOptions[InputField, FieldSize -> 5]; InputField[Dynamic[x]] the input field size is much larger than 5. But on the other hand InputField[Dynamic[x], Sequence @@ Options[InputField]] yields an input field with field size 5. ...and so on. What is the easiest way to work out (i.e. make a list of ...) which functions can't be used with SetOptions ? Answer As noted in the question, when you set an option to a function which appears not to work with SetOptions the options do get set, e.g. from the question: InputField[Dynamic[x], Sequence @@ Options[InputField]] but for whatever reason the global setting does not get used locally by default. Another interesting case i...

graphics3d - How to make this 3D animation

How to make an animation of following gif by using Mathematica ? Answer You can use the new RandomPoint - function to sample points from a sphere and use the second argument of Text to position some text in 3D. Text will automatically make the string face towards the viewer/camera for you. pts = RandomPoint[Sphere[{0, 0, 0}], 100]; Text["test", #] & /@ pts // Graphics3D[#, SphericalRegion -> True, Boxed -> False] & Another approach is instead of rotating the entire graphic (as done above) to rotate the points themselves via RotationTransform . This enables the use of coordinate dependent color and size. An appropriate transformation could be r[angle_?NumericQ, pivot : {_, _}] := RotationTransform[angle Degree, pivot]; and can be used as follows to archive something closer to what you are looking for with color and size scaling done in respect to the y -coordinate Graphics3D[(Style[Text["test", #], FontSize -> 14 - 5*(#[[2]] + 1), FontColor ->...

manipulate - How to express Log[b,x] in traditional format in Epilog without changing to base e

I'm having difficulty understanding how to express text in Epilog that is dynamically updated using Log[b, x] . Mathematica changes this to base $e$, but I would like it to be Log[b, x] in traditional format with base $b$, and I can't seem to make it work. I'm guessing I need to break the $\log$ apart using boxes or something, but don't know how to make a subscript that is a dynamically updated b value. Any ideas? Manipulate[ Plot[{b^x, x, Log[b, x]}, {x, -10, 10}, PlotRange -> {{-5, 5}, {-10, 10}}, PerformanceGoal -> "Quality", ImageSize -> All, Epilog -> {Text[b^x, {-3, 4}], Text[Log[b, x], {3, -5}]}, GridLines -> {Range[-10, 10, 1], Range[-10, 10, 1]}, GridLinesStyle -> Opacity[.04]], {{b, 2, "Choose a base"}, 0.01, 4}]

Palette with editable InputField

I need a Palette with editable InputField . When I try the code CreatePalette[Pane[InputField["Enter a string"]]] it makes the pallete but all input goes into an open notebook, not into the InputField . I just need a Palette with editable InputField . How it can be done?

plotting - ListPlot coloring by density of points

I have a set of points on a circumference: points = Map[Normalize, RandomReal[{-1, 1}, {3000, 2}]] And I want to plot them and color the plot according to the density of points. You see, ListPlot[points] will just give me that, a simple plot, on which I can't see the distribution of the points on the circle. Maybe there is an option? I hope you guys can understand what I'm asking. Edit/Follow-up I used SmoothDensityHistogram[points, ColorFunction -> "Rainbow"] Now, it's not really what I'm looking for. I kind of wanted it to look like the original plot, just the circumference. Answer Here is what I think you wanted: points = Map[Normalize, RandomReal[{-1, 1}, {3000, 2}]]; d = SmoothKernelDistribution[points]; colors = Hue /@ Rescale[PDF[d, #] & /@ points]; Graphics[Transpose[{colors, Point /@ points}]] Here the SmoothKernelDistribution is evaluated in the plane, giving you a two-dimensional interpretation of density. One could also understand your q...

interoperability - Executing Maple Code Inside Mathematica

Is it possible to create a C++ program using MathLink that will invoke the Maple kernel and execute a very basic Maple command using OpenMaple ? I'm envisioning the following: Install["MapleLink"] starts the external program MapleLink and installs the Mathematica definition RunMaple. ?RunMaple RunMaple[expr] returns a string which is the output when Maple executes expression expr A sample output: RunMaple["seq(i, i=0..3);"] "0, 1, 2, 3" I'm also open to using J/Link .

Trouble with differential equation

I tried to solve this differential equation: $$\epsilon y''(x)+xy'(x)=-\epsilon \pi^2 \cos(\pi x)-\pi x\sin(\pi x)$$ with boundary conditions: $y(-1)=-2, \space y(1)=0$. If we take $\epsilon=0.1$, Mathematica can solve it without any trouble Block[{e = 0.1, min = -1, max = 1}, Plot[Evaluate[ y[x] /. NDSolve[{e y''[x] + y'[x] x == -e Pi^2 Cos[Pi x] - Pi x Sin[Pi x], y[min] == -2, y[max] == 0}, y, {x, min, max}]], {x, min, max}] ] But if we want a smaller $\epsilon$, let say 0.01, Mathematica seems unable to handle it. Is there any options to invoke or methods to employ to get the desired result? Anyway, this is the solution for $\epsilon=0.0001$. Thank you.

plotting - ErrorListPlot ignores some exact data points

Bug introduced in 10.0 and persisting through 11.0.1 or later Possible bug? Or am I missing something / expecting too much? Forgive the seemingly obscure data points but they're where I've noticed a problem. Needs["ErrorBarPlots`"] data = { {{(10 π)/8, 1}, ErrorBar[0.2]}, {{(11 π)/8, 1}, ErrorBar[0.2]}, {{(12 π)/8, 1}, ErrorBar[0.2]} }; ErrorListPlot won't plot the error bar at x = 11 Pi/8 unless the datapoint is first converted to a numerical value: ErrorListPlot[data, PlotRange -> {{3, 5}, All}] ErrorListPlot[N@data, PlotRange -> {{3, 5}, All}] What gives?

calculus and analysis - How to solve system of integral equations

I have two integral equations -Sin[v]*Integrate[ (x[v, z]*b*Sin[v])/Sqrt[(a^2*(Cos[v])^2 + b^2*(Sin[v])^2 + z^2)^3], {v, 0, 2*Pi}, {z, -l, l}] == 0 -Cos[v]*Integrate[ (x[v, z]*a*Cos[v])/Sqrt[(a^2*(Cos[v])^2 + b^2*(Sin[v])^2 + z^2)^3], {v, 0, 2*Pi}, {z, -l, l}] == 0 I would like to find a solution for x[v,z] , that is going to solve bough equations. Than I would like to plot 3D graph with axes v,z,x[v,z] . Can anyone help me in finding a solution?

plotting - Vertical alignment of FrameTicks Text

The text box of a FrameTick/Tick is, by default, vertically centered on the tick mark. Is there a way to adjust the vertical alignment of FrameTick text? Specifically, along the lines of Baseline, Top, etc. FrameTicksStyle in combination with TextAlignment/Alignment does nothing. Answer Using Pane as in Szabolcs's answer or Framed (with FrameStyle->None ) with a combination of settings for BaselinePosition , ImageMargins and FrameMargins : Graphics[Circle[], Frame -> True, FrameTicks -> {{ {{0, Style[Pane[Style[0.4, FontSize -> 48], BaselinePosition -> (Top -> Bottom), FrameMargins -> {{0, 0}, {0, 40}}, ImageMargins -> {{0, 0}, {0, -40}}], Background -> Green]}}, {{0, Style[Pane[Style[0.4, FontSize -> 48], BaselinePosition -> (Bottom -> Top), FrameMargins -> {{0, 0}, {40, 0}}, ImageMargins -> {{0, 0}, {-40, 0}}], Background -> Green]}}}, {None, None}}] you get

List of parameters to manipulate from an external list

I am trying to achieve the functionality of the following command: Manipulate[ z, {z, 0, 1}] but with passing the list of parameters into Manipulate from and external list, e.g. list = {z, 0, 1}; Manipulate[ z, list ] but with this syntax I recieve the following error: Manipulate::vsform: "Manipulate argument \!\(ToExpression[{\"z\", 0, 1}]\) does not have the correct form for a variable specification" How can I overcome this?

functions - What Method options does FindInstance accept?

While answering Trouble getting FindInstance to return multiple results for certain constraints it became clear that for this problem FindInstance would be better starting with small numbers, but it does not. The function has the Method option: Options[FindInstance, Method] (* Out: {Method -> Automatic} *) Yet, I cannot find documentation for the methods (and methods sub-options) that it accepts. What Method values does FindInstance accept? Will any of them scan an Integer search space small-values first? Answer Here is a partial answer. I believe for Method -> {opts} , opts may be any of the following, and they will have whatever effect they have: Internal`InequalitySolverOptions[] Internal`ReduceOptions[] Internal`NSolveOptions[] (* {"ARSDecision" -> False, "BrownProjection" -> True, "CAD" -> True, "CADAlgebraicCoefficients" -> True, "CADBacksubstitution" -> Automatic, "CADCombineCells" -...

export - HTMLEncode a string for webMathematica

Is there a native Mathematica way to convert a string like " " into "<select></select>" and similarly for other reserved characters in html (or another encoding that would work in a webpage)?. I'm working in webMathematica , but I can't find a good way to print the literal string " " into a page - it either gets converted into an node, or it gets wrapped in a tag. Answer The different HTML entities are stored in System`Convert`MLStringDataDump`$HTMLEntities on version 9 and from here, it's a simple StringReplace : StringReplace[" ", System`Convert`MLStringDataDump`$HTMLEntities] (* "<select></select>" *)

computational geometry - How to generate higher dimensional convex hull?

I like the function ConvexHullMesh very much, since I often need to take a look at the convex hull of points in 2D and 3D, which as an example can be generated and tested with the following code (*Dimension and number of points*) d = 2; np = 10^1; (*Generate data*) data = RandomReal[{-1, 1}, {np, d}]; (*Construct convex hull*) ch = ConvexHullMesh@data; (*Test new point*) new = RandomReal[{-1, 1}, d]; Element[new, ch] True But just out of fun, I wanted to do this with higher dimensional points, e.g., with 5D points. But how do you construct in Mathematica the convex hull then out of the data and create a set for testing new points? Answer As others mentioned, Qhull can do this. There are multiple ways to access Qhull from Mathematica. One way is through the mPower package, now part of xCellerator . An answer discussing how to install the package is here: https://mathematica.stackexchange.com/a/18909/12 Another way is the qh-math package, described in this answer: https://mathematic...

plotting - Contour plot with points and line together

I have drawn a contour plot like attached in the figure. Now, what I want is points along with the contour lines. For example: I have a data set like: 3.34 4.74 0.0024 3.34 4.76 0.0023 ... ... .... I need the points (3.34, 4.74), (3.34, 4.76) etc.. on the curve. My code for this curve is: abc = Import[ "/home/users/chakrtdm/Desktop/test_25_alpha/ratio-a-b.dat", "Table"] ListContourPlot[abc, Contours -> 20, AspectRatio -> Automatic, PlotRange -> {{3.24, 3.38}, {4.70, 4.84}, {0.01, 0.001}}, ContourLabels -> All, TextStyle -> {FontFamily -> "Helvetica", FontSize -> 11, FontWeight -> Bold}, FrameLabel -> {"a", "b"}, Background -> White, LabelStyle -> Directive[Black, Bold]] Please suggest me how can I improve my code to get the points in the contour lines. Thanks! Answer data = Table[ Sin[i + j^2], {i, 0, 3, 0.1}, {j, 0, 3, 0.1}]; ListContourPlot[data, Epilog -> {PointSize[Large], Red,...

numerics - How to make the computer consider two numbers equal up to a certain precision

My problem is that I have a matrix A and the computer says is not Hermitian (self-adjoint). Then I check which elements make A be not Hermitian by calculating: B=A-ConjugateTranspose[A] In case A is Hermitian, this difference should be zero. What I get is that A is not Hermitian, but the terms I get in matrix B all are smaller than 10^(-15) . So how can I fix it at the beggining of the code, Mathematica to consider two numbers equal, up to 10^(-10) precision for example, that is it to consider: 0.00000000010 = 0.00000000012 or saying it in an alternative way, I want Mathematica to consider the following matrix Hermitian: C={{0,I*1.1^(-10)},{-I*1.2^(-10),0}} throughout the whole code. Answer You can check to see if the difference is less than the desired amount, and if so, to force it to be symmetric: c = {{0, I 1.1*10^(-10)}, {-I 1.2*10^(-10), 0}}; cOut = If[Abs[Total[c - ConjugateTranspose[c], 2]] 1/2 (c + ConjugateTranspose[c]),c] This returns cOut which is a sy...

performance tuning - Transferring a large amount of data in parallel calculations

Bug fixed in version 11.1 Functions like MemberQ , FreeQ , etc. no longer unpack. Yay! This question is inspired by one of @whuber's answers Consider the following code: μ = RandomReal[{0, 1}, 100]; Σ = DiagonalMatrix[Exp[RandomReal[{0, 1}, 100]]]; AbsoluteTiming[ RandomVariate[MultinormalDistribution[μ, Σ], 400000];] It runs in 3.5 seconds here. Now let's parallelize it: LaunchKernels[] AbsoluteTiming[ Join @@ ParallelTable[ RandomVariate[MultinormalDistribution[μ, Σ], 200000], {2}];] This runs in 6.3 seconds on a 2-core machine---much slower . It also uses a lot of memory (which I checked using a process monitor). Now let's suppress returning the results from the subkernels by including a semicolon: AbsoluteTiming[ Join @@ ParallelTable[ RandomVariate[MultinormalDistribution[μ, Σ], 200000];, {2}];] This one runs in 2.6 seconds---a speedup. What is happening here? Why is the calculation that returns the result so much slower? Is it a general rule with parallel calcu...

programming - Morphing a "sheet of paper" into a torus

How can I visualize the standard topological "rubber-sheet" construction of a torus, that is, morphing a square into a torus? How can I start or are there any examples in the Mathematica documentation respectively? Answer Edit I had some time so I've added full surface torus. Old code in edit history. DynamicModule[{x = 2., l = 100., x2 = 2., l2 = 100., grid, fast, slow}, Grid[{{ Graphics3D[{ Dynamic[Map[{Blue, Polygon[#[[{1, 2, 4, 3}]]]} &, Join @@@ (Join @@ Partition[#, {2, 2}, 1]) ]&[ ControlActive[fast[l, l2], slow[l, l2]]] ] }, PlotRange -> {{-7, 7}, {-7, 7}, {-1, 2}}, ImageSize -> 600, Axes -> True, BaseStyle -> 18] , Column[{ Slider[Dynamic[x, (l = 10.^#; x = #) &], {.0001, 2.}], Slider[Dynamic[x2, (l2 = 10.^#; x2 = #) &], {.0001, 2.}] }] }}] , Initialization :> ( grid[l_, l2_, n_, ...

Computing kernels of a matrix in distinct cases

I have the following matrix: {{c1 d1 - e1 f1, c1 d2 - e1 f2, c1 d3 - e1 f3, c1 d4 - e1 f4}, {c2 d1 - e2 f1, c2 d2 - e2 f2, c2 d3 - e2 f3, c2 d4 - e2 f4}, {c3 d1 - e3 f1, c3 d2 - e3 f2, c3 d3 - e3 f3, c3 d4 - e3 f4}, {c4 d1 - e4 f1, c4 d2 - e4 f2, c4 d3 - e4 f3, c4 d4 - e4 f4}, {c5 d1 - e5 f1, c5 d2 - e5 f2, c5 d3 - e5 f3, c5 d4 - e5 f4}, {c6 d1 - e6 f1, c6 d2 - e6 f2, c6 d3 - e6 f3, c6 d4 - e6 f4}} From the way it's constructed, I know that the rank can be at most two. If I ask Mathematica to row reduce the matrix, I end up with the following: {{1, 0, (d3 f2 - d2 f3)/(-d2 f1 + d1 f2), (d4 f2 - d2 f4)/(-d2 f1 + d1 f2)}, {0, 1, (d3 f1 - d1 f3)/(d2 f1 - d1 f2), (d4 f1 - d1 f4)/(d2 f1 - d1 f2)}, {0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0}} As can be seen, this is only valid when $d_1f_2\neq d_2f_1$. Furthermore, if I look at the left kernel, then Mathematica gives me the following basis: {{-((-c6 e2 + c2 e6)/(c2 e1 - c1 e2)), -((c6 e1 - c1 e6)/(c2 e1 - c1 e2)), 0...

calculus and analysis - Generating a polynomial that's accurate to within an error of no more than 1/10^5

I'm currently stuck on a question for class that asks... "Find a polynomial p[x] that you can use to calculate 6 ArcTan[x] to within an error of no more than 10^(-5) for all the x's with -(1/Sqrt[3]) ." I used a series expansion below. Clear[x]; approx6arctan[x_] = Normal[Series[6Tan[x], {x, 0, 200}]] However, this can only generate a function that's accurate only to the fourth decimal, no matter how much I expand the series (200 is already huge). Any hints on how to generate a polynomial that's accurate to the fifth decimal? Thanks in advance.

numerical integration - Gillespie Stochastic Simulation Algorithm

The Gillespie SSA is a Monte Carlo stochastic simulation algorithm to find the trajectory of a dynamic system described by a reaction (or interaction) network, e.g. chemical reactions or ecological problems. It was introduced by Dan Gillespie in 1977 ( see paper here ). It is used in case of small molecular numbers (or species abundance) where numerical integration of the related differential equation system is not appropriate due to hard stochastic effects (i.e. the death of a single individual might make a large impact on the population). Can you do it with Mathematica ? Is it general enough? Is it fast? Answer Yes you can. Below is a fairly general, Mathematica -compiled, fast and robust version. Examples 1. Michaelis-Menten kinetics Michaelis-Menten kinetics for enzyme-directed substrate conversion. The enzyme ( e ) converts the susbtrate ( s ) through an enzyme-substrate complex ( c ) to the product ( p ). For comparison, I've included the deterministic ODE system solved b...