title: Equation solving with GroebnerBasis

text: I have a system of quadratic equations defined below               nx = 1;     ny = 2;     ne = 1;     nt = 2;     n = nx + ny;     ns = 2;              For[ii = 1, ii <= ns, ii++,           G[ii] = Table[g[ii][i, j], {i, ny}, {j, nx + ne + 1}];           H[ii] = Table[h[ii][i, j], {i, nx}, {j, nx + ne + 1}];           ];              eqns = {-1.1111111111111112` g[1][1,                1] + (-0.9703966865736877` + 1.` g[1][1, 1] + 0.9` g[1][2, 1] +                 0.11111111111111109` g[2][1, 1] +                 0.09999999999999998` g[2][2, 1]) h[1][1, 1],             10.` g[1][1, 1] -              161.` g[1][2,                1] + (143.451` g[1][2, 1] + 15.938999999999995` g[2][2, 1]) h[                1][1, 1],             0.7763173492589502` + 0.6199999999999999` g[1][2, 1] -              0.9703966865736877` h[1][1, 1], -1.1111111111111112` g[2][1,                1] + (-0.9703966865736877` + 0.11111111111111109` g[1][1, 1] +                 0.09999999999999998` g[1][2, 1] + 1.` g[2][1, 1] +                 0.9` g[2][2, 1]) h[2][1, 1],             10.` g[2][1, 1] -              161.` g[2][2,                1] + (15.938999999999995` g[1][2, 1] + 143.451` g[2][2, 1]) h[                2][1, 1],             0.7763173492589502` + 0.13999999999999996` g[2][2, 1] -              0.9703966865736877` h[2][1, 1]};              polys = eqns /. {0. -> 0};      This system has some structure. First, I find the variables               vars = Variables[polys];       Deconstruct the polynomials (I know it is a quadratic system). There are 2 nonzero constant terms, 12 nonzero linear terms, and 12 nonzero 2nd-order terms               {m0, m1, m2} = CoefficientArrays[polys, vars]      Confirm the reconstructed polynomials are the same as the originals (just checking)               polys - (m0 + m1.vars + m2.vars.vars) // Expand // Chop       Of the 12 nonzero 2nd-order terms, of which there are 8 distinct cross- products ...               Flatten[Most[ArrayRules[#]] & /@ m2][[All, 1]] // Union // Length      ... none of which are purely quadratic               Diagonal /@ m2 // Union      So, it is a quadratic system but nothing purely quadratic If I apply `GroebnerBasis`               GBordering =       Flatten[{H[1][[1 ;; nx, 1 ;; nx]], H[2][[1 ;; nx, 1 ;; nx]],         G[1][[1 ;; ny, 1 ;; nx]], G[2][[1 ;; ny, 1 ;; nx]]}]     GB = GroebnerBasis[polys, GBordering];     GB // MatrixForm     soln = Solve[Thread[GB == 0]];     soln // MatrixForm      I get the right solutions. Anyone interested, I will send the code. But you can copy this into your computer and it will give you 9 solutions. Everything OK up to here. But if I make the system a larger (only a little)               nx = 2;     ny = 3;     ne = 1;     nt = 1;     n = nx + ny;     ns = 2;     For[ii = 1, ii <= ns, ii++,       G[ii] = Table[g[ii][i, j], {i, ny}, {j, nx + ne + 1}];       H[ii] = Table[h[ii][i, j], {i, nx}, {j, nx + ne + 1}];       ];     eqns = {9.169301858525786` -          1.` g[1][3,            1] + (-19.453328842995205` + 8.169847955946475` g[1][2, 1] +             0.9077608839940525` g[2][2, 1]) h[1][1,            1] + (0.` + 8.169847955946475` g[1][2, 2] +             0.9077608839940525` g[2][2, 2]) h[1][2, 1],         0.` - 1.` g[1][3,            2] + (-19.453328842995205` + 8.169847955946475` g[1][2, 1] +             0.9077608839940525` g[2][2, 1]) h[1][1,            2] + (0.` + 8.169847955946475` g[1][2, 2] +             0.9077608839940525` g[2][2, 2]) h[1][2, 2], -1.` g[1][3,            1] + (0.` - 1.` g[1][1, 1] + 0.9` g[1][3, 1] -             0.11111111111111109` g[2][1, 1] +             0.09999999999999998` g[2][3, 1]) h[1][1,            1] + (1.0999999999999999` - 1.` g[1][1, 2] + 0.9` g[1][3, 2] -             0.11111111111111109` g[2][1, 2] +             0.09999999999999998` g[2][3, 2]) h[1][2, 1], -1.` g[1][3,            2] + (0.` - 1.` g[1][1, 1] + 0.9` g[1][3, 1] -             0.11111111111111109` g[2][1, 1] +             0.09999999999999998` g[2][3, 1]) h[1][1,            2] + (1.0999999999999999` - 1.` g[1][1, 2] + 0.9` g[1][3, 2] -             0.11111111111111109` g[2][1, 2] +             0.09999999999999998` g[2][3, 2]) h[1][2,            2], -178.88888888888889` g[1][1, 1] -          9.` g[1][3,            1] + (159.39` g[1][1, 1] + 17.709999999999994` g[2][1, 1]) h[1][           1, 1] + (159.39` g[1][1, 2] + 17.709999999999994` g[2][1, 2]) h[            1][2, 1], -178.88888888888889` g[1][1, 2] -          9.` g[1][3,            2] + (159.39` g[1][1, 1] + 17.709999999999994` g[2][1, 1]) h[1][           1, 2] + (159.39` g[1][1, 2] + 17.709999999999994` g[2][1, 2]) h[            1][2, 2],         0.` - 0.21199999999999997` g[1][1, 1] +          0.99` h[1][2, 1], -0.792` - 0.21199999999999997` g[1][1, 2] +          0.99` h[1][2, 2], 0.` + 1.` g[1][2, 1] - 1.` h[1][1, 1],         0.` + 1.` g[1][2, 2] - 1.` h[1][1, 2],         9.169301858525786` -          1.` g[2][3,            1] + (-19.453328842995205` + 0.9077608839940525` g[1][2, 1] +             8.169847955946475` g[2][2, 1]) h[2][1,            1] + (0.` + 0.9077608839940525` g[1][2, 2] +             8.169847955946475` g[2][2, 2]) h[2][2, 1],         0.` - 1.` g[2][3,            2] + (-19.453328842995205` + 0.9077608839940525` g[1][2, 1] +             8.169847955946475` g[2][2, 1]) h[2][1,            2] + (0.` + 0.9077608839940525` g[1][2, 2] +             8.169847955946475` g[2][2, 2]) h[2][2, 2], -1.` g[2][3,            1] + (0.` - 0.11111111111111109` g[1][1, 1] +             0.09999999999999998` g[1][3, 1] - 1.` g[2][1, 1] +             0.9` g[2][3, 1]) h[2][1,            1] + (1.0999999999999999` - 0.11111111111111109` g[1][1, 2] +             0.09999999999999998` g[1][3, 2] - 1.` g[2][1, 2] +             0.9` g[2][3, 2]) h[2][2, 1], -1.` g[2][3,            2] + (0.` - 0.11111111111111109` g[1][1, 1] +             0.09999999999999998` g[1][3, 1] - 1.` g[2][1, 1] +             0.9` g[2][3, 1]) h[2][1,            2] + (1.0999999999999999` - 0.11111111111111109` g[1][1, 2] +             0.09999999999999998` g[1][3, 2] - 1.` g[2][1, 2] +             0.9` g[2][3, 2]) h[2][2, 2], -178.88888888888889` g[2][1, 1] -          9.` g[2][3,            1] + (17.709999999999994` g[1][1, 1] + 159.39` g[2][1, 1]) h[2][           1, 1] + (17.709999999999994` g[1][1, 2] + 159.39` g[2][1, 2]) h[            2][2, 1], -178.88888888888889` g[2][1, 2] -          9.` g[2][3,            2] + (17.709999999999994` g[1][1, 1] + 159.39` g[2][1, 1]) h[2][           1, 2] + (17.709999999999994` g[1][1, 2] + 159.39` g[2][1, 2]) h[            2][2, 2],         0.` - 0.14799999999999996` g[2][1, 1] +          0.99` h[2][2, 1], -0.792` - 0.14799999999999996` g[2][1, 2] +          0.99` h[2][2, 2], 0.` + 1.` g[2][2, 1] - 1.` h[2][1, 1],         0.` + 1.` g[2][2, 2] - 1.` h[2][1, 2]};     polys = eqns /. {0. -> 0};      but with exactly the same structure Find the variables               vars = Variables[polys]      Deconstruct the polynomials (I know it is a quadratic system). There are 10 nonzero constant terms, 228 nonzero linear terms, and 324 nonzero 2nd-order terms               {m0, m1, m2} = CoefficientArrays[polys, vars]      Confirm the reconstructed polynomials are the same as the originals (just checking)               polys - (m0 + m1.vars + m2.vars.vars) // Expand // Chop       Of the 324 nonzero 2nd-order terms, of which there are 180 distinct cross- products ...               Flatten[Most[ArrayRules[#]] & /@ m2][[All, 1]] // Union // Length      ... none of which are purely quadratic               Diagonal /@ m2 // Union      Mathematica cannot handle it. I am not sure whether it is a problem of memory, a problem of speed or something like that. I wounder whether I could take advantage of the quadratic structure without purely quadratic components to simplify the problem and make it possible for Mathematica to handle it.


title: How can I scale the y axis?

text: I have a file I import and read into a table of $2 \times n$. I can plot this data using the following code:               directory = "c:\nameofdirectory";     SetDirectory[directory];     filenames = FileNames[];     readlist = ToString[filenames[[1]]];     table = OpenRead[readlist];     data = ReadList[table, {Number, Number}];     Close[table];          ListPlot[data, Frame -> True, PlotRange -> {All, All}, Joined -> True,        FrameLabel -> {"Voltage (V)", "Current (A)"}];      but the $y$-axis is order of $10^{-6} A$. I'd prefer to have my $y$-data scaled by $10^{-6}$ then label my axis is `microA` instead of `A`. I've tried using `Ticks` but I don't really understand how that works (it never changes my plot at all when I use it). If there's no easy way to do this, can someone tell me how to directly divide my $y$-values from the input? I can separate the $y$-values from "data" and divide them using               xdat = List[];     ydat = List[];     For[j = 1, j < Length[data] + 1, j++,         AppendTo[xdat, data[[j]][[1]]];         AppendTo[ydat, data[[j]][[2]]];        ];          y=ydat/10^-6;      but I don't know how to place them back into the form as before so that it can be easily plotted with `ListPlot`.


title: Using Manipulate to vary the frequency of a sound produced by SoundNote[]

text: I have tried various ways to combine Manipulate and Sound/SoundNote that would allow me to vary the frequency of a sound while it is being produced. The addition of Dynamic didn't help, but did increase the variety of error messages. Is this just not possible?               Manipulate[EmitSound[Sound[SoundNote[n, 3]]], {n, 0, 10}]      My goal to listen to, or observe, the beat frequencies produced as a tone becomes closer to or further away from twice the frequency of a second, steady tone. My underlying issue is to develop insight into how original idea of an octave.


title: Building adjacency matrix from set data

text: I have two sets: games and players. Players pick their games. As such I will have data say g1={p1,p3,p5},g2={p2,p4},g3={p2,p3,p5}. My interest though is to build the connections among players because of the common games they have played. For instance, since p1 played with p3 and p5 in game 1, then they are connected. I can then build an adjacency matrix for the players: {{0,0,1,0,1},{0,0,1,1,1},{1,1,0,0,1},{0,1,0,0,0},{1,1,1,0,0}}.I'll be doing this for hundreds of games and thousands of players. Question is: how do I efficiently get the adjacency matrix given the data? Thanks for the inputs.


title: How can I replace Polygon by Line in the output of ParametricPlot3D while preserving colours?

text: {m0, m1, m2, m3, m4, m5, m6, m7} = {4, 3, 2, 3, 6, 2, 6, 4};     r = Sin[m0*phi]^m1 + Cos[m2*phi]^m3 + Sin[m4*theta]^m5 + Cos[m6*theta]^m7;     g = ParametricPlot3D[{r*Sin[phi]*Cos[theta], r*Cos[phi],          r*Sin[phi]*Sin[theta]}, {phi, 0, Pi}, {theta, 0, 2 Pi},         ColorFunction -> (Blend[{Red, Yellow, Green, Cyan, Blue}, #] &),         Mesh -> None, Boxed -> False, Axes -> 0];     g /. Polygon -> Line      gives this: ![enter image description here](http://i.stack.imgur.com/BAssf.jpg) The colour information is gone. I want to get the output below: ![enter image description here](http://i.stack.imgur.com/ukhBI.jpg) How can I get it?


title: FindRoot: Minimum Execution Time

text: Is there a minimum execution time for the command FindRoot? If so, can I get rid of it? For example, if I compute the time to solve the fundtion x=0, it takes around 0.27ms per evaluation:               In[581]:= AbsoluteTiming[Do[FindRoot[x, {x, 300}], {try, 1, 100}]]          Out[581]= {0.0270027, Null}      However, this is virtually the same time as it takes to solve a much more complex equation system {x=0,y=0,z=0}               In[580]:= AbsoluteTiming[Do[FindRoot[{x, y, z}, {{x, 300}, {y, 300}, {z, 300}}], {try, 1, 100}]]          Out[580]= {0.0280028, Null}      Since I am running findroot many times (the actual system is more complicated), I care about this "fixed" running cost. Is there a way to get rid of it?


title: Mapping multiple functions over an array

text: I am looking for an elegant way to Map multiple functions over an array. My problem looks something like this: Consider a list of functions:               f = {f1,f2,f3,f4} (* example, can be larger *)      And an array:               A = {{a1, b1, c1, d1, e1}, {a2, b2, c2, d2, e2}, {a3, b3, c3, d3, e3}}; (*example, can be larger*)      **Q1:** I’d like to map the functions in f over the array A so that the final result looks like:               B= {{1 - a1, f1[a1]*b1, f2[a1]*c1, f3[a1]*d1, f4[a1]*e1}, {1 - a2,      f1[a2] b2, f2[a2] c2, f3[a2] d2, f4[a2] e2}, {1 - a3, f1[a3] b3,      f2[a3] c3, f3[a3] d3, f4[a3] e3}}      ![Sample input/output](http://i.stack.imgur.com/WdrAF.png) **Q2:** How do I map the functions a little differently so that the output looks like:               c= {{a1, f1[a1] b1, f2[b1] c1, f3[c1] d1, f4[d1] e1}, {a2, f1[a2] b2,      f2[b2] c2, f3[c2] d2, f4[d2] e2}, {a3, f1[a3] b3, f2[b3] c3,      f3[c3] d3, f4[d3] e3}} // tf      ![Alternate mapping](http://i.stack.imgur.com/McXL7.png)


title: Merging two lists by prepending elements

text: Being a newcomer to _Mathematica_ , I have the following question. Given               x = Range@9~Partition~3;      and               y = {a, b, c};      I want to combine `x` and `y` so that               result = {{a, 1, 2, 3}, {b, 4, 5, 6}, {c, 7, 8, 9}};      I have written               MapPrepend[a_, b_] := Prepend[a[[#]], b[[#]]] & /@ Range@Length@b      which does the job, but I wonder if there are other, more functional solutions? Probably this question has already been answered. In this case it would be very nice if you just would direct me to the links.


title: how to create Dataset after importing txt file

text: I would like to translate a imported txt file to Dataset.               test = Import["E:\\mathematica\\example.txt", "table",             {"FieldSeparators" -> ",", "CharacterEncoding" -> "ASCII"}]      The result is: > {{"naam A", "naam B", "naam C", "naam D"}, {1, 2, 1.2, "a"}, {2, 3, 2.3, > "b"}, {3, 4, 2.3, "c"}, {4, 5, 4.5, "d"}, {5, 6, 5.6, "e"}, {6, 7, 6.7, > "f"}, {7, 8, 7.8, "g"}, {8, 9, 8.9, "h"}, {9, 10, 9.1, "i"}, {10, 11, 10.1, > "j"}} Now I would like to translate it to `Dataset` which is new in _Methematica 10.0_.               dataset = Dataset[{        <|"naam A" -> 1, "naam B" -> "2", "naam C" -> {1.1}, "naam D" -> {a}|>,        <|"naam A" -> 2, "naam B" -> "3", "naam C" -> {1.2}, "naam D" -> {b}|>,        <|"naam A" -> 3, "naam B" -> "4", "naam C" -> {1.3}, "naam D" -> {c}|>,        ....        }]      I tried it with the next code:               data = Dataset[<|"key1" -> test[[All, 1]], "key2" -> test[[All, 2]],         "key3" -> test[[All, 3]], "key4" -> test[[All, 4]]|>]      This code gives the wrong output. Anyone a suggestion how to do it?


title: Dealing with numbers too large for machine precision in Graphics

text: `Graphics` only supports machine precision numbers (i.e. number that can be converted to machine precision). Take for example               Graphics[Point[{0, Exp[1000.]}]]      which produces an error because `Exp[1000.]` is too large to be representable as a `double`. What is the _simplest_ way to deal with situations like this? Replacing all too large numbers with just some arbitrary very large machine precision number is fine, as they'll be outside the plot range anyway. I'm sure there must be some built-in function that can fix this easily.


title: Why is my plot cut off?

text: I am trying to use `Show` to combine two plots so I do the following:               g1 = Plot[x, {0, 100}, AspectRatio -> 1]     g2 = ListContourPlot[data, FrameLabel -> {"label1", "label2"},        PlotLabel -> "title", AspectRatio -> 1,        DataRange -> {{0, 100}, {0, 100}}]          g3 = Show[{g1, g2}]      Everything looks great on the computer. Now, I want to export.               Export[g3, "test.pdf", ImageResolution -> 600]      This creates a problem. The y-axis label is completely missing. The image is cut off on the left-hand side. It's not white space or anything like that. It's just not there. The whole left side of the image is missing. The title of the plot is partially cut off too. How do I fix this?


title: Sort per column

text: I have some data:               data={{"a", 1, 1}, {"a", 1, 2}, {"a", 1, 3},            {"c", 2, 1}, {"b", 2, 2}, {"b", 2, 3},            {"c", 3, 1}, {"a", 3, 2}, {"a", 3, 3}}      When I use:               Sort[data]      I get the list ordered in incremental order first with column 1 then by 2 and last by the third. Suppose I want the first column in increasing order but the second in decreasing order and the third again in increasing order. How do I do that. I tried with Sort and SortBy but that doesn't work. Output should be: > a 3 2 > > a 3 3 > > a 1 1 > > a 1 2 > > a 1 3 > > b 2 2 > > b 2 3 > > c 3 1 > > c 2 1 It's probably a simple question but apparently not obvious to me.


title: Simplify function with normalized vector

text: I made a very simple case:               testF[x_, y_, z_] = x^2 + y^2 + z^2     Simplify[testF[x, y, z]]      If {x, y, z} is normalized, the simplified result should just be 1. How could I specify this condition/assumption?


title: "Off-label use" of Mathematica functions for performance

text: Assisting an acquaintance with a simulation where a very large number of binomial probability calculations are required, I pointed them to a slightly obscure function to dramatically speed up their algorithm. Specifically, for               Probability[x == b, x \[Distributed] BinomialDistribution[n, p]]      the same result can be had with               BernsteinBasis[n, b, p]      but the latter is 15-20X faster (and surprisingly generally more accurate with inexact _p_ ) than the former, and 1.5-2X faster than _PDF_ calculation. Only explicit calculation using the formula for binomial probability is faster. I'm curious what other "off-label" uses of built-in functions have been noted to have performance advantages.


title: Numerically solving a 2D oscillating integral

text: I'm having trouble solving this integral numerically:               Integrand2 :=  I/λ Holo[ξ, η] E^(I (2 π)/λ r)/r     ResultsTable2 = Table[NIntegrate[Integrand2, {ξ, -20, 20}, {η, -20, 20},           Method -> "MultidimensionalRule"], {x, -15, 15, 5}, {y, -15, 15, 5}]      The function called `Holo` is modelling an aperture that looks like this: ![enter image description here](http://i.stack.imgur.com/gE7Zo.jpg) To put it simply: I'm trying to do _Fresnel_ diffraction with this as the aperture. My goal is to get a 2D inverse _Fourier transform_ by summing over the entire `Holo` term (a function of `ξ` and `η`). As you can see from the image, it oscillates rapidly, so unless I integrate over a very small area of `ξ` and `η`, I get the following error: > NIntegrate::slwcon: Numerical integration converging too slowly; suspect one > of the following: singularity, value of the integration is 0, highly > oscillatory integrand, or WorkingPrecision too small. To get a good 2D inverse _Fourier transform_ , I need to have a larger area. Ideally, integrating over 100 by 100 would be good. Is there some kind of method or work around that can make this happen? I can provide more details if needed. I'm new to _Mathematica_ and to this site, so any answer is appreciated. I can learn from any comments or suggestions you have :) EDIT: I guess to make my question a bit clearer: _how can I numerically integrate a 2D oscillating function?_


title: Point belonging to a Disk in a Cartesian Coordinate System

text: The image below represents a human subject fixations while observing this abstract pattern for 3 seconds. I would like to know how much time they spent looking at the actual disk. That is a red point within a Black Disk. ![enter image description here](http://i.stack.imgur.com/GcxFx.png) Each fixations (each sublist of "fixations" below is coded as such : {{X,Y},Duration in ms}                fixations=      {{{20.3899, 14.8931}, 238}, {{27.0063, 18.8899},428},        {{25.8113, 24.8679}, 377}, {{24.2579, 22.022},106},        {{25.3208, 24.022}, 130}, {{21.739, 12.1792},175},        {{29.2673, 8.88994}, 295}, {{30.3868, 17.6572},160},        {{31.217, 22.6761}, 145}, {{22.9686, 20.6918},155},        {{19.6321, 20.2704}, 145}}      Each Disk is coded as {X,Y,R}                Disks=      {{22.8176, 19.9696, 0.974938}, {29.5314, 10.7197, 0.974938},        {17.5112,19.7207, 0.974938}, {30.8997, 23.2454, 0.974938},       {28.0588,6.09759, 0.974938}, {30.8524, 17.0661, 1.53205},        {21.0393, 10.7137,1.53205}, {25.451, 25.1336, 1.53205}}      How could I : -Count the number of fixations that are within a disk. -Then how much actual time they spent there.  The only solution I thought about was to compute the EuclideanDistance of each fixation to each disk but I feel there is a smarter faster way to do this. THank You for your attention.


title: Throwing away all expression on a loopback link

text: What is the simplest way to throw away all expressions waiting to be read on a loopback link, effectively resetting the link? I am building an expression on a loopback link using the C MathLink interface. If the operation succeeds, I want to send the expression to `stdlink`. If it fails, I need to _clear_ the loopback link. What is the recommended way to do this?


title: Problem with general assumptions

text: I hope I am not to redundant with my questions. As part of solving a physical problem, I want to assume once and for all that all my variables and parameters are real, positive and finite. Is it possible make such assumptions for all _Mathematica_ functions? I tried this piece of code I found on the forum:               $Assumptions = Element[_, Reals] && {_ ∈ Reals} && {_} > 0;      But still, `Sign` does not aways return 1, `Minimize` does not treat all parameters as negative, etc.


title: Define a mathematical set

text: I have a terribly short question: Is there a way to define a pure mathematical set? I tried a bit of googling, but the word "set" has two meanings and the other one (as in setting) as by far more dominant and I couldn't find any reference to mathematical sets in _Mathematica_. As for the application, I would like to be able to define two sets `A`, `B` and test whether the two sets are the same or not, namely have the same elements.


title: Using ParametricNDSolveValue with initial conditions of another ParametricNDSolveValue

text: I have a system of parameter-dependent ODEs, which I can solve using `ParametricNDSolveValue`. Now, I have another ODE, which has as initial condition the solution of one component of the first system at a specific time point. Consider the example               pfun = ParametricNDSolveValue[{x'[t] == a*x[t], x[0] == 1}, {x}, {t, 0, 10}, {a}];     qfun = ParametricNDSolveValue[{y'[s] == a^2*y[s], y[0] == First[pfun[a]][2]}, {y}, {s, 0, 10}, {a}];      Then,               qfun[4]      gives the output: _ParametricNDSolveValue::ndinnt: "Initial condition 4.[2.] is not a number or a rectangular array of numbers"_. This is odd, since `First[pfun[.5]][2]` produces $2.71828$ as output, as one expects. However, omitting the curly brackets of `x` in the first line gives the correct result:               pfun = ParametricNDSolveValue[{x'[t] == a*x[t], x[0] == 1}, x, {t, 0, 10}, {a}];     qfun = ParametricNDSolveValue[{y'[s] == a^2*y[s], y[0] == pfun[a][2]}, {y}, {s, 0, 10}, {a}];     First[qfun[4]][2]      produces the output $2.35385*10^{17}$, which is $\exp(2a)\exp(a^2s)$ for $a=4$ and $s=2$. As mentioned, `pfun` is a _system_ of equations in the real problem, so omitting curly brackets is no option for me. Any ideas? **UPDATE** I found out that the problem with the first two lines of code is that the fragment `First[pfun[a]]` is first evaluated with `a` having no value, thus `First[pfun[a]]` reduces to `a`. How can I specify that `First` is evaluated _after_ `a` is given a specific value?


title: How to use NDSolve to track a deviation from equilibrium?

text: i have a task to find a event, to achieve this i think need to integrate until a constant function suddenly breaks the equilibrium and give a sharp peak. I think that EventLocator should work, thus tried the proposals from How to use NDSolve to track equilibrium?, with a threshold, but the numerical calculator is still complaining that the value before and after the step are not real numbers, but this would not be the case in theory, the code:               psi=Pi;     For[i=1,i<141,i++,w[i]=1/2];     fext=1;     k = 35.529;     ci = Table[y[i][0] == RandomReal[{-Pi, Pi}], {i, 1, 141}];     AppendTo[ci,y[0][0] == y[141][0] - psi];          lckm = {y[0]'[t] == fext, y[141]'[t] == fext,Table[y[i]'[t] ==               w[i] +                k (Sin[y[i + 1][t] - y[i][t]] +                   Sin[y[i - 1][t] - y[i][t]]), {i, 1,               140}]};          sol = NDSolve[{lckm, ci},        Table[y[i], {i, 0, 141}], {t, 0.0, 800.0},       WorkingPrecision -> 30, AccuracyGoal -> 15, PrecisionGoal -> 15,        MaxSteps -> Infinity];          v1 = y[140][t] /. sol[[1]];     v2 = y[141][t] /. sol[[1]];     norm = Sin[v1 - v2] /. t -> 800;     k = 35.526;         ci = Table[            y[i][0] ==              Part[{y[i][t] /. sol} /. t -> 800.0, 1], {i, 1, 141}];         AppendTo[ci, y[0][0] == y[141][0] - psi];               lckm = {y[0]'[t] == fext, y[141]'[t] == fext,Table[y[i]'[t] ==           w[i] +            k (Sin[y[i + 1][t] - y[i][t]] +               Sin[y[i - 1][t] - y[i][t]]),{i, 1,           140}]};     sol = NDSolve[{lckm, ci},        Table[y[i], {i, 0, 141}], {t, 0.0, 800.0},        Method -> {"EventLocator",          "Event" ->           Boole[Abs[             Sin[y[140][t] - y[141][t]]/norm - 1] >             1/10], "EventAction" -> Print["t:", t]},        WorkingPrecision -> 30, AccuracyGoal -> 15, PrecisionGoal -> 15,        MaxSteps -> Infinity]      All the constant have a value well ascertained value, including norm that is a normalization responsable for the deviation detect. The answer to this, in mathematica 7, is:               NDSolve::precw: The precision of the differential equation ({((y[0])^\[Prime])[t]==1,((y[141])^\[Prime])[t]==1,((y[1])^\[Prime])[t]==1/2+66.439 (Sin[Subscript[\[Null], <<2>>][<<1>>]+Times[<<2>>]]-Sin[Plus[<<2>>]]),<<6>>,((y[8])^\[Prime])[t]==1/2+66.439 (Sin[\[Null], <<2>>][<<1>>]+Times[<<2>>]]-Sin[Plus[<<2>>]]),<<274>>}) is less than WorkingPrecision (30.`). >>          NDSolve::evre: The value of the event function at t = 2.79284074034703538067228005241`30.*^-10 was not a real number.  The event will be ignored in steps where it does not evaluate to real numbers at both ends. >>      Please note that this logic is set to work, i.e. this error is not dependent in the value of norm in the sense that it needs it to detect the deviation. The point is: why Mathematica is outputing an error saying that something is not a real number ? where is the thing that is not a real number ? Iam sorry if the example is too extense, i think that the problem poses it in this form or iam too newbie in mathematica and dont know rewrite. Does anyone has any suggestion in order to circumvent it ? Thank you !!!


title: Colours in BarChart[]

text: This question is a follow up from this one. I want to be able to specify a list of colours for my bars in `BarChart`. In responses to the question I linked to, @kguler produced the following, which sets the bars to the left of x=8 as a different colour:               cedF[{{xmin_, xmax_}, {ymin_, ymax_}}, ___] :={If[xmax <= 8,         RGBColor[1, 0, 0], Sequence[]],       Dynamic@EdgeForm[Directive[Thickness[.015], Lighter@CurrentValue["Color"]]],      Rectangle[{xmin, ymin}, {xmax, ymax}, RoundingRadius -> 5]};          Histogram[RandomVariate[NormalDistribution[10, 2], 500],ChartStyle -> "Pastel",        ChartElementFunction -> cedF]      ![kguler example](http://i.stack.imgur.com/lHIGV.png) My question is how to generalise this: I would like to be able to have a list that specifies the colours of every bin uniquely. How would I go about this? If[] statements seem a bit clunky for this application.


title: How do I create a histogram with unequal bin widths?

text: I would like to create a histogram for human population by age, but the $2010$ census data I have available is for five uneven cohorts: ages $0$-$17$, $18$-$24$, $25$-$44$, $45$-$64$, and $65+$. My county has the following data: `{86031, 26671, 91927, 93983, 32232}`, for a total population of $330844$. The corresponding histogram bins should have widths `{17, 7, 20, 20, 20}` and heights `{86031/17, 26671/7, 91927/20, 93983/20, 32232/20}`. I will eventually want to make this pretty by adding labels, etc., but I don't see how to even get started.


title: Is the System Information that Mathematica displays not completely correct?

text: Or maybe I did not understand the meaning of Resolution. This is the code I used -               SystemInformation["Devices", "ScreenInformation"]      I get 120 as my system resolution, and Full Screen Area as 1366x768 (~14" diagonal screen). I checked with two of my colleagues who have mathematica installed. The one with a ~15.4" diagonal screen got Resolution as 96 & Full Screen Area as 1366x768. This made sense, sort of, a bigger screen is lesser dense pixel-wise for the same pixel screen area. But when I checked with another ~14" diagonal screen laptop, mathematica returned Resolution as 96 & Full Screen Area as 1366x768. This doesn't make sense, two physically equivalent screens with same pixel area but different resolutions. So I cross-checked in these two websites - http://pediddle.net/dpi-ruler.html and http://concentriclivers.com/screen_density.html Both my laptop and the other 14" laptop got similar dpi's of 112, and the 15" laptop returned a dpi value of about 100. I did the math and these numbers fit it just well. So is Mathematica wrong here? Or am I barking up the wrong tree for the dpi value? I am new here, so thanks for any help extended, in advance!


title: Solving differential equations with sums (power series)

text: I have sets of 10 differential equations, but for this purpose I'll demonstrate what I need on one example that can be solved by hand. My equation is this: $$4GJ\Omega(\theta)\Omega'(\theta)\xi^\theta(\tau,r,\theta,\varphi)+4GJ\Omega^2(\theta)\frac{\partial \xi^\theta}{\partial\theta}=\mathcal{O}(r^{-1})$$ I need to assume that the solution of $\xi^\theta$ is in the form $$\xi^\theta=\sum\limits_n\xi_n^\theta(\tau,\theta,\varphi)r^n$$ and it truncates for large n, where n is negative. If I put that in my equation, and say that all powers greater or equal of 0 will be zero, I get $$\Omega'\xi^\theta_n+\Omega \xi_{n,\theta}^\theta=0,\quad for\ n\geq0$$ where $\xi_{n,\theta}^\theta=\frac{\partial \xi_n^\theta}{\partial\theta}$ for shorthand notation simplification. Now, this basically means that my general solution for $\xi^\theta$ is $$\xi^\theta=\xi^\theta_{-1}r^{-1}+\mathcal{O}(r^{-2})$$ which is what I need to get. Now, my question is: is it possible to do this for all other equations in Mathematica. I got my equations, and I made substitution:               ξθ=Sum[ξ[τ, θ, φ][m]*r^m, m]      and got the result               4 G J Ω[θ]((Sum[D[ξ[τ, θ, φ][m],θ]*r^m, m])Ω[θ]+(Sum[ξ[τ, θ, φ][m]*r^m, m])Ω'[θ])      As I should. Is it possible now, to somehow get the most general result, like I got by solving by hand? I'm asking because I have 10 equations like this, and some have other terms ($\xi^r,\xi^\tau,\xi^\varphi$) that I need to find the general solution to. Do I actually expand it to some order and somehow truncate higher terms? EDIT: The power series solution is in the negative expansion.


title: Attaching a file to a notebook

text: Is there a way to attach a file into a notebook and open it later with, for instance, a double click (or another action button)? In Mathematica windows version, the insert menu has the "object..." entry and I can indeed embed an object into my notebook. But then, I don't know how to exploit it in an attachment kind of way. For instance, in MS Word, if we "insert an object", we can then open it by double clicking its icon: windows opens the file with its predefined application. For the purpose of the discussion, let’s suppose I want to attach a pdf file. If this functionality isn't native on Mathematica, probably a technique based on this post can be used, with an intermediate step of export to a temporary folder, and then a `Run` kind of command... I can also imagine a `Dynamic` interface, with a list of attached files (whose internal data/content is “kept” by a `DynamicModule` internal variable), and four buttons: add, delete, export and open. Since a `Dynamic` cell can be easily copied from a notebook to another, I could easily use this small embedding app on different notebooks. There’s probably another button that would be useful: import. This would make the file content available, as a string, on a global context variable, or at least export the file to a temporary folder, and make its path available on a global context variable. (if you go the "dynamic app" way, please consider, from the beginning, app conflicts if two apps are added to the same notebook; and a more tricky, or probably impossible task, to make some of its functionality work on the Player) **EDIT - 2014-01-12** **(there was an error with the OPEN, that has now been corrected)** Current version: ![enter image description here](http://i.stack.imgur.com/PP3Rv.png) Code and example file: download here Any help on making it better is appreciated, specially in the domain of the compression ratio and safety (currently using GZIP). Future developments (if anyone can help):   * attach the notebook itself without passing through files (poor man self-included versioning)   * attach a Save of a specific variable/definition or a list of them, without passing through files.   * attach a Save of the entire session, to record current state of sessions (kind of persistent memory)   * load stored Saves (with the option of overlaying the existing memory, substitute completely the existing memory (clearing other definitions), adding definitions with an indexed definition name or or an indexed context, so to allow comparisons...)   * export recorded files to variables, so to work on them (ImportString, etc) without passing through files   * and any other crazy infeasible idea... **EDIT - 2014-02-09** It is now possible to save the definitions of the current context. It is very rough, and I would greatly appreciate some help to improve it. Things that are not great:   * it stores its own definition, which could probably be avoidable   * something better should be thought for the problem of creating a new context registry when an old one is opened   * why only the current context, and not a checkbox tree list where one checks contexts and/or definitions that would like to be stored   * no warning for very large content being saved The file can be downloaded here ![enter image description here](http://i.stack.imgur.com/5yLw0.png)


title: “automated” image embedding app

text: Based on the following image attach procedure, can a `Dynamic` app be developed that:   1. at start, only shows a small import button   2. once clicked, opens a choose image file dialog   3. after the image file chosen by the user, asks for representation dimension (pixels or percentage of the screen, etc), some other stuff possible: type and degree of compression (native, JPG, etc), frame, tear down effect etc   4. after options chosen, removes the button and replaces it with the low footprint image following the chosen options   5. right click on the image, and we can choose two "empty" the space (to restart the process) Actually, there should be two import buttons: one for a file dialog, and another for importing from the clipboard (very useful when used with a screen snag tool). An output cell with this dynamic app could then be easily duplicated and past to different notebooks, in order to facilitate the procedure of the above referenced post.


title: Provided an array of 2D point clouds, can we quickly replace distance- and index-based cliques of points with their median or mean values?

text: Imagine I have $N \approx 10^5$ clouds of points in an array of the form: UPDATED               numClouds = 10^3;     cloudList = Table[Table[{RandomInteger[{0, 2}] + RandomReal[{-0.1, 0.1}], RandomInteger[{0, 2}] + RandomReal[{-0.1, 0.1}]}, {j, 1, RandomInteger[{0, 50}]}], {i, 1, numClouds}];      Notice that we generate points here by first randomly snapping them to a coordinate in an integer lattice, and then adding some real number valued noise:               ListPlot[Flatten[cloudList,1]]      I'd like to move a sliding window over `cloudList` of size $k$ (e.g. $k = 10$) where we take all of the points in `{cloudList[[i]], cloudList[[i+1]], ..., cloudList[[i+k]]}`, establish cliques for points within a Euclidean distance $d$ of one-another (i.e. if we can place a circle of diameter $d$ over a set of points, it's a clique), and then replace all of the points within each clique with the median or mean point value of the given clique. I would like to do this without disturbing the array structure of `cloudList`. In terms of finding a fast way to grab all the points within a critical radius of $\frac{d}{2}$, we can use a poorly documented version of `Nearest`: > `Nearest[data, x, {n, r}]`   >  give up to the `n` nearest elements to `x` within a radius `r` The above was from an answer by the user "ssch" (http://mathematica.stackexchange.com/users/1517/ssch) to my question: Selecting for 2D points that are within a threshold distance of an upper- and lower-bound number of points * * * Update - I neglected to specify what should occur when points can belong to more than one clique. The user belisarius (http://mathematica.stackexchange.com/users/193/belisarius), for example, asked about the cliques for a point set `Tuples[{1, 0}, 2]` where `r = 3/2`. Here, all four points belong to the same clique:               list = Tuples[{1, 0}, 2];     nf = Nearest[list];     nf[list[[1]], {Infinity, N[3/2]}]          out :: {{1, 1}, {1, 0}, {0, 1}, {0, 0}}      So we'd set each point in this clique to `{1/2,1/2}`. However, if we set `r = 1`, we have `Binomial[4, 3] = 4` seperate cliques with overlapping points:               {{1, 1}, {1, 0}, {0, 1}}     {{1, 0}, {1, 1}, {0, 0}}     {{0, 1}, {1, 1}, {0, 0}}     {{0, 0}, {1, 0}, {0, 1}}      To resolve this, I'd like to generate each clique as above, compute a mean for each clique, and where points fall in multiple cliques, to randomly assign these points the mean or median value selected with uniform probability from each clique. However, I'd be open to simpler solutions that are easier to implement. * * * (2nd) Update - I have changed the method of generating the random test set s.t. we now first snap a randomly generated test point to a bounded integer lattice coordinate, then add noise. * * * (3rd) Update - This is more of a note, but I would be fine employing a simple greedy algorithm to partition points into clusters, or using `FindClusters` with a `EuclideanDistance` distance function (this is slightly less comfortable since I'm not exactly sure what the underlying algorithm is).


title: How do I remove the subscript after calling BaseForm[]?

text: Say I want to convert some number $N$ to base $k$. Here, I might use the command `BaseForm[N,`k]`, e.g.,`BaseForm[10,`2]` which gives the output 1010[2] (where [2] is subscript). How can I drop the subscript and manipulate the number as if it is a normal integer in decimal form? I'd ideally like to do this without any string processing?


title: Extracting a single layer from a SHP file

text: I have an ESRI ShapeFile (SHP) that I've loaded. It has a number of layers which I can see the names of by using:               Import["file.shp", "LayerNames"]      I get something like > {"Africover_Kenya_Spatial_Aggregation", "ke-landcover", "ke-forests", \ "ke- > rangeland", "ke-agriculture", "ke-waterbodies", "ke-bareareas", \ "ke- > urban", "ke-politicalboundary", "theme1"} I would like to import just one layer ("ke-landcover" for example) from the data but I can't figure this out. If I ask for the "Graphics" element, I think I get the layers individually, but I don't know how to verify whether this is true or not. I'm new to mathematica, and documentation for GIS applications is limited at best. Can anyone help me? Thanks!


title: Maximizing Large Multivariable Function over Region with Parameter Constraints

text: Brand new Mathematica user here. I have a complicated function to maximize, and thought I would try software like Mathematica instead of taking derivatives by hand. I found some Mathematica documentation and believe I've entered the code correctly. Here is my entire code:               Maximize[{(1/2)*(1/(1+2*y-y^2-x^2-2*x^2*y+2*x*y^2))*(a-c*(y^2/(-2*x^2*y^2+2*x*y^2-y^2+2*y+1)+((1+2*y-2*y^2-x^2-2*x^2*y+2*x*y^2)/y^2)*(1-1/(-2*x^2*y^2+2*x*y^2-y^2+2*y+1)))), 0<=x<=1, 0<=y<=1, 0<c<a}, {x,y}]      I believe my syntax is correct, however I'm in doubt over whether I've put my constraints in the right place. The variables $x$ and $y$ are the probabilities with which I want to maximize my function, so I "constrained" them between $0$ and $1$. And the parameters $a$ and $c$ are such that $c>0$, and $a>c$ (these are not variables -- I only want to maximize over the probabilities). Running the script causes Mathematica to say "Running..." and never give a result. Does Mathematica have a limit on how complicated a function can be? Thanks!


title: Mapping points on the unit circle

text: I have a set of points P that are evenly spaced on a square grid of length 1... I would like to map these points to the circle of radius 1 as nicely as possible, while preserving relative positions as much as possible. That is, the points in the square grid should map to the points in the unit circle. How can this be done in Mathematica?


title: Find names similar to a given name in a large data set

text: I'm trying to find a 'similar' name to every item in a list                nearName[name_, namesList_] :=         Module[{lista, pos},          lista = EditDistance[ToLowerCase[name], ToLowerCase[ToString[#]]] & /@ (namesList);          pos = Flatten[Position[lista, Min[lista]]];          If[Length[pos] > 1 || Min[lista] > 3 , Missing["NotAvailable"], namesList[[pos]]]]      The problem is that with a very large input (`nameList` contains about 400k items and I make about 1k calls to `nearName`), `nearName` is using a lot of CPU. Is there any way I can make this faster?


title: How to create a spliced survival statistical distribution?

text: I am trying to use Mathematica for more of my reliability engineering work and am having trouble using SplicedDistribution. I have a Weibull distribution (dist1) associated with a low stress condition, β1 = 2.59 and η1 = 25889. I have a different Weibull distribution (dist2) associated with a high stress condition with the characteristic life being 1/3 that for dist1, β2 = 2.59 and η2 = 8630. The need is to create a spliced distribution such that for the first 5000 hours only dist1 applies and thereafter only dist2 applies. Think of a device that operates in a relatively clean, stress-free environment for the first 5000 hours of its life, then it is abruptly placed in a stressful, corrosive atmosphere thereafter. My attempt in using SplicedDistribution does not seem to work as shown in the comparison plot below. SplicedDistribution appears to be “mixing” the two distributions for all time, t as opposed to letting only dist1 govern for the first 5000 hours then allowing only dist2 to govern thereafter. I would expect that dist3 would be equivalent to dist1 for the first 5000 hours then its survival (reliability) curve would depart from dist1. Can someone show me how to accomplish what I need using Mathematica ? ![enter image description here](http://i.stack.imgur.com/BvJRy.png) Here is my version 9 code:               {β1 = 2.59, η1 = 25889.};     dist1 = WeibullDistribution[β1, η1]     plot1 = Plot[SurvivalFunction[dist1][t], {t, 0, 18000},        PlotStyle -> Black, Frame -> True,        FrameLabel -> {"Time (hours)", "Reliabilty,  R[t]"},        PlotRange -> {0, 1.05},        PlotLabel -> "Distribution 1 Reliability Plot"]          {β2 = 2.59, η2 = 8630.};     dist2 = WeibullDistribution[β2, η2]     plot2 = Plot[SurvivalFunction[dist2][t], {t, 0, 18000},        PlotStyle -> Blue, Frame -> True,        FrameLabel -> {"Time (hours)", "Reliabilty,  R[t]"},        PlotRange -> {0, 1.05},        PlotLabel -> "Distribution 2 Reliability Plot"]          plot3 = Show[plot1, plot2, PlotRange -> All,        PlotLabel -> "Distribution 1 & 2 Comparison Plot"]          dist3 = SplicedDistribution[{1, 1}, {0, 5000, ∞}, {dist1, dist2}]     plot4 = Plot[SurvivalFunction[dist3][t], {t, 0, 18000},        PlotStyle -> Red, Frame -> True,        FrameLabel -> {"Time (hours)", "Reliabilty,  R[t]"},        PlotLabel -> "Distribution 3 (Spliced) Reliability Plot"]          plot5 = Show[plot3, plot4, PlotRange -> All, AxesOrigin -> {0, 0},        PlotLabel -> "Comparison Plot",        Epilog -> {Text["Distribution 1 (low stress)", {15000, .9}],         Text["Distribution 2 (high stress)", {10500, .6}],         Text["Distribution 3 (spliced)", {6000, .2}]}]


title: Downvalues vs. Scoping for Functions

text: Regarding my recent question on using a default value for a function argument when a pattern was not met yielded some interesting answers, but the general consensus was "Yes this can be done, but there are better options." Assuming the same conditions that certain function arguments (`rt` and `pt` in this case) should default to zero when less than zero, the three main options are as follows: **Downvalues:**               func[qi_, dei_, b_, dmin_, rt_, pt_, t_] :=           func[{qi, dei, b, dmin, Max[0, rt], Max[0, pt], t}]          func[{qi_, dei_, b_, dmin_, rt_, pt_, t_}] := (*do stuff*)      _Or..._               func[qi_, dei_, b_, dmin_, rt_?Negative, pt_, t_] := func[qi, dei, b, dmin, 0, pt, t]     func[qi_, dei_, b_, dmin_, rt_, pt_?Negative, t_] := func[qi, dei, b, dmin, rt, 0, t]          func[qi_, dei_, b_, dmin_, rt_, pt_, t_] :=(*do stuff*)      **Scoping:**               func[qi_, dei_, b_, dmin_, rt_, pt_, t_] := Module[{rtt = Max[0,rt], ptt = Max[0,pt]},        (*do stuff with rtt and ptt instead*) ]      What are the pros and cons of each function definition? Which is the most efficient?


title: How can this type of optical illusion be created in Mathematica?

text: I see these around the web and would like to make them in Mathematica.   Combining them in an array is actually quite mesmerizing! ![Moving Illusion](http://i.stack.imgur.com/zHC6P.jpg)


title: large matrix eigenvalue problem

text: I need solve a very large complex matrix (not sparse and not symmetry) eigenvalue problem, e.g., 1e4*1e4 or even 1e6*1e6. How large dimensions of the matrix can Mathematica support? And, how about the run time? Or, any one have good suggestions for this?


title: Vectors with a certain magnitude in Mathematica

text: For a user specified magnitude, is there a way to have Mathematica produce any 3D vector that fits that magnitude?


title: Solving a system of equations with variable coefficients

text: I am trying to solve a system of equation which seems to be simple but much as I tried I couldn't solve it in _Mathematica_ and it says the coefficients need to be exact. Now, I was wondering if someone walk me through the step by step procedure of solving such equations. I have a variable, `phi`, which is a function of `P`. `P` is a real number and `0 < P < 1`.               phi = 45 beta^P;       where `beta` is a constant. Now the system of equations:               x^(-1/2) y^(-P) - x^(-1/2) y^P phi^2 = 3 phi     2 x y^(-P) - x y^(1 - P) = 4/phi


title: How do I use a while-loop to find the gain of a transfer function's limit?

text: At the moment I can find the limit of my transfer function as follows:               (5 4001.6`)/(5 4001.6` + (4 + s) (200 + s)) /. s -> 0      This works as I would expect and provides a final value of 0.961553 In this example my value of gain, k, is 4001.6` so re-writing the expression for arbitrary k would give:                  (5 k)/(5k + (4 + s) (200 + s)) /. s -> 0      I want to find a value of k that makes the result 0.95. To do this I was thinking about using a while loop. I thought the while loop could start with k=1, compute the limit, test to see if less than 0.95, and if not increment k by one, and perform the test again, until the correct value of k is found: This is what I have tried:                   k = 1; While[(5 k)/(5 k + (4 + s) (200 + s)) /. s < 0.95, k++]      However I get the error: "ReplaceAll::reps: {s<=0.95} is neither a list of replacement rules nor a valid dispatch table, and so cannot be used for replacing. >>" Any help would be most appreciated.


title: compile function and floating points

text: The compile function I have written computes a table of relative frequencies of some data. Thus, the output is a list of increase elements, p, in [0,1]. Since I have say 1000 data points the elements p are equal to x/1000 and x in [0,1000]. Hence, each p should have a maximum of three positions after the decimal point. However, if I run the compile function then values having many more positions after the decimal point appear. However, this appears to be true not for all elements of the resulting list. For example, let the true p equal 0.894, then the compile output reads p = 0.8949999999999. How can I force the result to be 0.894? Unfortunately, `Round[p,1./1000.]`does not work within the compile function. It does not work either if I apply `Round`to the resulting `list`. What works, but seem inefficient, is `Round[list,1/1000]//N`. Thanks for help.


title: Understanding "IntensityData" coordinates

text: Update - Playing around a bit, is it true that we need to process the list of "IntensityData" values $r$ as follows:               r = Partition[r,k];     r = Transpose@r;      And then apply:               vmapped= c + {v[[1]] - k/2 - 0.5, k/2 - v[[2]] + 0.5}      It looks like this properly maps the point $v$ to the larger image containing the morphological component with "IntensityData" $r$, however, it's difficult to tell if there's a small error here? ** I would like to ask a follow-up to a previous question - Understanding the computation of IntensityCentroid? \- since this question was answered well, but I am still a bit confused about the coordinates for "IntensityData". Let $m$ be a square morphological component with $k \times k$ pixels. Let $c$ be the output of `ComponentMeasurements[{m,img},"Centroid"]`, let $r$ be the output of `ComponentMeasurements[{m,img},"IntensityData"]`. For an example, we can set:               k = 7;     c = {157.5, 92.5};     r = {0.509752, 0.958211, 0.177492, 0.679586, 0.651634, 0.159438, 0.249294,                0.840233, 0.403427, 0.706836, 0.0924188, 0.523653, 0.907517, 0.219175,           0.134596, 0.599349, 0.320155, 0.417359, 0.767355, 0.280942, 0.960531,           0.332421, 0.756754, 0.186532, 0.188323, 0.457144, 0.80254, 0.252036, 0.66885,           0.0437862, 0.674661, 0.125856, 0.134035, 0.62525, 0.328863, 0.84693,           0.644811, 0.877901, 0.0685183, 0.307649, 0.481993, 0.695101};      We then partition $r$ in the following way to structure the data in the form of ImageData:               r = Partition[r,k];      I then do "something" to calculate a coordinate with respect to the coordinate frame in `Image[r]`. Call this coordinate $v$, where we might have, for example, `v = {4.232,7.827}` (perhaps this is an intensity centroid in the the coordinate system for $r$). Using image coordinates rather than indices, i.e. the coordinate system used when outputting a value for `ComponentMeasurements[{m,img},"IntensityCentroid"]`, how can I use only the component centroid $c$ to map $v$ to my original image, `img`? Note that I only provide you with $k$, $c$, and $v$ for this mapping, not the original image. Ok, let's provide a small example using ten morphological components. Here, the first set of coordinates in each three-element array represents $v$, which is something similar to an intensity centroid computed with the output of `r = Partition[r,7];` where $r$ is the `"IntensityData"` output for some box-like morphological component with a 7-pixel edge-length. The second set of coordinates in each array ({72.5,250.5}, {208.5,246.5}, etc.) represents the output of `ComponentMeasurements[{m,img},"Centroid"]` for the morphological component used to grab $r$. The final set of coordinates is the output of `ComponentMeasurements[{m,img},"IntensityCentroid"]` for the particular morphological component.               componentOneArray={{4.1209,2.39237},{72.5,250.5},{73.1307,250.376}}     componentTwoArray={{4.10345,2.51433},{208.5,246.5},{208.993,246.415}}     componentThreeArray={{2.74777,11.3037},{18.5,229.5},{16.6019,229.926}}     componentFourArray={{3.75292,2.4059},{92.5,228.5},{93.1028,228.757}}     componentFiveArray={{3.74452,3.1833},{109.5,224.5},{109.302,224.697}}     componentSixArray={{3.76436,2.40024},{137.5,223.5},{138.098,223.743}}     componentSevenArray={{3.52801,2.6132},{221.5,221.5},{221.89,221.992}}     componentEightArray={{4.10626,2.43366},{64.5,209.5},{65.0603,209.391}}     componentNineArray={{4.08676,3.06379},{98.5,208.5},{98.4326,208.419}}     componentTenArray={{3.90224,3.10008},{30.5,205.5},{30.3848,205.606}}      How can we use use the second element of each list, $c$, to map the first element of each list, $v$, to the coordinate system used for the computed intensity centroid of the morphological component?


title: Hide both input cell and output cell

text: You can double-click the right gutter of a input cell to hide a corresponding output cell or double-click the gutter of the output cell to hide the input cell. I would like to hide both the input and the output cell. Is this possible to do within the GUI?


title: Calculate proportion of values for dates

text: At work, we were discussing when is it the best time to change to winter tires for bikes and/or cars. Using `WeatherData[]` and `DateListPlot[]`, it was fairly straightforward for me to create the diagram below: ![Mean temperature per day in Stockholm](http://i.stack.imgur.com/Y5wNT.png) _Fig 1 Mean temperature per day in Stockholm, red is negative and a risk without winter tires._ The code for this is               cityTemp = WeatherData["Stockholm", "MeanTemperature", {{1977, 1, 1}, {2011, 12, 31}, "Day"}];     iceRiskDays = Select[cityTemp, Last[#] < 0 &];          yearStrip[ dataItem_] := {{0, Part[First[dataItem], 2], Last[First[dataItem]]}, Last[dataItem]}          DateListPlot[{yearStrip[#] & /@ cityTemp, yearStrip[#] & /@ iceRiskDays} ]      My question is: **How do I calculate for each day the proportion of values for that day that are below 0 Celsius?** (e.g. for a date at the end of November, the proportion is likely to be bigger than 0.5)? My attempts to do this ended with trying to create separate lists for each date, but I felt that this was a less elegant way and also creates less "fit" with DateListPlot.


title: How to normalize parts of matrix elements?

text: To normalize the last 2 or more columns of a matrix,               m = {{0, 0, 135, 276}, {0, 1, 762, 852}, {0, 2, 762, 852}, {1, 0, 162,     852}, {1, 1, 362, 852}, {1, 2, 262, 852}}      I can do that by               mm = Transpose@      m /. {m[[All, 3]] -> m[[All, 3]]/(m[[All, 3]] + m[[All, 4]]),       m[[All, 4]] -> m[[All, 4]]/(m[[All, 3]] + m[[All, 4]])} // Transpose // MatrixForm      Or step by step               mt = Transpose[ m]          msum=mt[[3]]+mt[[4]]          m3=mt[[3]]/msum          m4=mt[[4]]/msum          mlast = Transpose[mt /. {mt[[3]] -> m3, mt[[4]] -> m4}]      are there any more elegant method in Mathematica?


title: Apply rule only if applicable? Check rule applicability? Custom unapplicable rule?

text: How to avoid recursion in the following example?               Clear[MyRule];     Clear[MyHead1];          MyRule1 = MyHead1[x_MyHead2] | MyHead1[x_MyHead3] -> x;          MyHead1[x_] := Replace[MyHead1[x], MyRule1];          MyHead1[a]          MyHead1[MyHead2[b]]      How to make first rule in example below treated by `Replace` as not applicable?               In[17]:= ComplexRule1 = a -> Block[{},         12         ];          In[18]:= ComplexRule2 = a -> Block[{},         13         ];          In[19]:= Replace[a, {ComplexRule1, ComplexRule2}]          Out[19]= 12


title: How to write a MathLink function in C that returns nothing

text: What is the proper way to write a C MathLink function that returns nothing at all? The following works, but I'm not sure if all this work is actually necessary: The template file:               :Begin:     :Function:      fun     :Pattern:       fun[]     :Arguments:     {}     :ArgumentTypes: {}     :ReturnType:    Manual     :End:      The C code:               void fun() {         MLPutSymbol(stdlink, "Null");     }      In short, here I set the return type to `Manual` and explicitly returned `Null`. Is this all necessary or is there a shorter way?


title: Why does the mean in my coin tossing simulation not approach 1/2?

text: I simulate a simple experiment, a coin flip. What I do is accumulate the mean of the results up to the i-th experiment. What I can't figure out is why the computed means do not asymptotically approach 1/2. What did I do wrong? The code:               n = 10^4;     c = Table[{i, RandomInteger[]}, {i, 1, n}];     s = 0;     m = {};     For[i = 1, i <= n, i++,       s += c[[i, 2]];       AppendTo[m, {i, N[s/i]}];     ]


title: Convert boolean test function to pattern?

text: What is a "reverse" to `MatchQ`? `MatchQ` converts pattern into boolean function. So how to convert back? Suppose I want to use boolean function in `Switch`.


title: Vlookup function as Excel in Mathematica

text: How can we create a Vlookup function, similar to what we have in Excel? So if we have a table `a` and another `b`, we could add columns to table `a`, based on a common key to both, in order to create a table `c`. ![Tables](http://i.stack.imgur.com/89Bad.png)


title: Custom attributes?

text: Is it possible to assign custom attributes to symbols and check them later?               SetAttributes[a, b]      says               Attributes::attnf: b is not a known attribute.


title: How do I overlay a histogram with a plot of CDF?

text: data = Flatten[Import["speedtest1B.csv"]];     distrib = HistogramDistribution[data];     histogram = Histogram[data, ChartStyle -> {Blue}];     Export["speedtest1B.pdf",          Overlay[{          Histogram[             data,              ChartStyle -> {Blue},             ImageSize -> {500,500},             AxesLabel -> {"start-up time [ms]", "number of samples"}             ],         DiscretePlot[             CDF[distrib, x],              {x, 0, Max[data]},             ImageSize -> {500,500},             Axes -> False             ]             }         ]     ];      From the above I get: ![enter image description here](http://i.stack.imgur.com/G4mqh.png) I would like the plots to be properly aligned. Displaying the CDF x-axis on the right side of the plot would be a big plus. Also I don't realy get color- theaming so could someone please explain how to easily change the color od CDF to red.


title: Upsampling matrix by inserting zeros inbetween elements

text: I want to essentially take an image with $n \times n$ pixels and expand it into a sparse array, spacing out the pixels by a factor $2$ or $3$. Each pixel in the input image occupies one of the corners of a cell of $2 \times 2$ or $3 \times 3$ new pixels in the expanded image, with the other $n \times (n - 1)$ pixels (in each cell) set to $0$. So, I want to create a padded array where: If one line of data is               WXYZ      then the output should be               W0 X0 Y0 Z0       00 00 00 00      I think _Mathematica_ 's `SparseArray` might be the way to go, but I'm not sure how to go about implementing this.


title: Export list elements each to it's own page in a multi-paged PDF document

text: I have a list, `myList`, of let's say 20 (different) objects               myList = Table["object", {20}]      I know would like to export every single element of `myList` to its own PDF page in a 20-paged PDF document; i.e., `myList[[1]]` comes on page 1, ..., myList[[20]] comes on page 20. Is this possible in _Mathematica_? I have heard of the `Export` function, of course, but have not seen any way of combining 20 exports in one single document.


title: How to constrain the generation of all possible orderings?

text: Here is code from Simon Woods' answer for getting all possible weak (equal ranks allowed) orderings for $N=3$ objects:                ClearAll[f]; SetAttributes[f, Orderless];      ReplaceList[f[a, b, c], f[a___, b___, c___] :> {{a}, {b}, {c}}] //      DeleteCases[#, {}, -1] & // Union // Column      It gives $13$ such orderings:               {{a, b, c}}     {{a}, {b, c}}     {{b}, {a, c}}     {{c}, {a, b}}     {{a, b}, {c}}     {{a, c}, {b}}     {{b, c}, {a}}     {{a}, {b}, {c}}     {{a}, {c}, {b}}     {{b}, {a}, {c}}     {{b}, {c}, {a}}     {{c}, {a}, {b}}     {{c}, {b}, {a}}      How can I modify this code for the case when not more than $2$ subsets are allowed? The desired output is:               {{a, b, c}}     {{a}, {b, c}}     {{b}, {a, c}}     {{c}, {a, b}}     {{a, b}, {c}}     {{a, c}, {b}}     {{b, c}, {a}}      I am trying to find way for doing such reductions in general $N$ and for any number of subsets-restriction.


title: Setting OperatorSubstitution to false causes Plot's labels of sufficiently small size not to appear in the specified typeface. Is this a bug?

text: If I set `OperatorSubstitution` to `False`, the labels of `Plot` will not use the typeface specified by `FontFamily` unless `FontSize` is more than certain number. In this test scenario, that certain number is `14` and the typeface is `Ubuntu` (Regular) (TTF fonts available from font.ubuntu.com).               SetOptions[Plot, AspectRatio -> 1, ImageSize -> {200, 200}];      Row[{       Plot[x, {x, 0, 1},          LabelStyle -> {FontFamily -> "Ubuntu", FontSize -> 14},          BaseStyle -> {PrivateFontOptions -> {"OperatorSubstitution" -> True}}],        Plot[x, {x, 0, 1},          LabelStyle -> {FontFamily -> "Ubuntu", FontSize -> 14},          BaseStyle -> {PrivateFontOptions -> {"OperatorSubstitution" -> False}}],        Plot[x, {x, 0, 1},          LabelStyle -> {FontFamily -> "Ubuntu", FontSize -> 15},          BaseStyle -> {PrivateFontOptions -> {"OperatorSubstitution" -> False}}]     }]      ![15 is the threshold.](http://i.stack.imgur.com/w4rMJ.png) Only some typefaces are affected (eg Inconsolata is fine). Is this a bug? I'm using Mathematica 9 for Mac.


title: Block Matrix Algebra with Mathematica

text: I have come up with some BlockMatrix Algebra for Mathematica to make notations easier. I have the following:               Needs["Notation`"];     Notation[ParsedBoxWrapper[     RowBox[{"A_", "\[CenterDot]", "B_"}]] \[DoubleLongLeftRightArrow]         ParsedBoxWrapper[     RowBox[{"BlockMultiply", "[",      RowBox[{"A_", ",", "B_"}], "]"}]]];     Notation[ParsedBoxWrapper[     SuperscriptBox["A_", "Inv"]] \[DoubleLongLeftRightArrow]         ParsedBoxWrapper[     RowBox[{"BlockInv", "[", "A_", "]"}]]];     Notation[ParsedBoxWrapper[     SuperscriptBox[     RowBox[{"(", "A_", ")"}], "T"]] \[DoubleLongLeftRightArrow]         ParsedBoxWrapper[     RowBox[{"BlockTransp", "[", "A_", "]"}]]];      And the actual matrix algebra functions:               ClearAll[BlockTransp]     (*SetAttributes[BlockTransp,OneIdentity]*)          BlockTransp[M_List] := Map[BlockTransp[#] &, Transpose[M], {2}]      BlockTransp[M__ + B__ ] := BlockTransp[M] + BlockTransp[B]     BlockTransp[a_ /; NumberQ[a]] := a;     BlockTransp[a_Scalar ] := a;     BlockTransp[a_   M__ /; NumberQ[a]] := a BlockTransp[M];     BlockTransp[a_Scalar   M__] := a BlockTransp[M];          BlockTransp[BlockMultiply[x_, y_] ] :=       BlockMultiply[BlockTransp[y], BlockTransp[x]]     BlockTransp[BlockTransp[x_]] := x          ClearAll[BlockMultiply]     BlockMultiply[mats1 : {{_ ..} ..}, mats2 : {{_ ..} ..}] :=       Inner[BlockMultiply, mats1, mats2]     (*For numbers*)          BlockMultiply[a_ A_ /; NumberQ[a], B_] := a BlockMultiply[A, B]     BlockMultiply[A_ , a_ B_ /; NumberQ[a]] := a BlockMultiply[A, B]     BlockMultiply[A__, a_ /; NumberQ[a]] := a A     BlockMultiply[a_ /; NumberQ[a], A__ ] := a A     (*For scalar symbols*)          BlockMultiply[a_Scalar A_, B_] := a BlockMultiply[A, B]     BlockMultiply[A_ , a_Scalar B_ ] := a BlockMultiply[A, B]     BlockMultiply[A__, a_Scalar ] := a A     BlockMultiply[a_Scalar, A__ ] := a A                    ClearAll[BlockInv]     BlockInv[mats1 : {{_ ..} ..}] := Map[BlockInv[#] &, mats1, {2}];     BlockInv[ a_ /; NumberQ[a]] := a          (*Make some formatting expression*)     Scalar /: MakeBoxes[Scalar[a_], StandardForm] :=        MakeBoxes[a, StandardForm];      I can do now some funny cool thing like:   **1\. Example:**   ![First Example](http://i.stack.imgur.com/LEmOc.png)   **2\. Example:**   ![Second Example](http://i.stack.imgur.com/s2Btb.png)   **3\. Example:**   ![Third Example](http://i.stack.imgur.com/o5P5z.png)   **4\. Example:**   ![Fourth Example](http://i.stack.imgur.com/42Yfj.png) **Remark:** Of course the `BlockInverse` Function is sensless in the sense that it does not really invert correctly. The input to `BlockInverse` should only be a **diagonal block matrix**! There is a general formula for block matrices Block Inversion, but it is not implemented here for the sake of simplicity. I have two questions which I could not come up so far:   1. How can I make the 3. example such that the `scalar` variable is treated like a scalar? Can I overload `NumberQ[scalar]:= ...` ? Or is there a better way? **[SOLVED]** See the above definitions with added `_Scalar` that matches any pattern with head `Scalar`. See example 4!   2. How can I pretty print these things such that superfluous brackets are neglected? Is there some elegant way for the function definitions such that, associativity and non-commutativity is handled well?


title: A good choice to find solutions of ODE

text: In order to rephrase the question. I would ask you to take a look at the following problem:   1. I have this system of ODE: $$ \begin{cases} \dot{x}(t)= x(t)[r_1(t) -b_1(t)x(t- \tau(t)) - c_1(t)y^m(t)] & t\neq t_k\\\ \dot{y}(t)= y(t)[- r_2(t)- b_2(t) y(t- \sigma(t)) +c_2(t) x(t) y^{m-1}(t)] & \\\ \Delta x(t_k)= x(t_k^+)- x(t_k)=d_k x(t_k) & t=t_k\\\ \Delta y(t_k)= y(t_k^+)- y(t_k) = f_k y(t_k) & \end{cases}. $$   2. $x(t)$ and $y(t)$ are discontinuous function at points $t_k$.   3. All these functions: $r_1(t), \, r_2(t), \, b_1(t), \, b_2(t), \, c_1(t),\mbox{ and } c_2(t)$ are positive and periodic with the same period $\omega$.   4. $\tau(t)$ and $\sigma(t)$ are continuously differentiable positive periodic functions with the same period $\omega$ such that $\tau'(t) <1 \mbox{ and } \sigma'(t)<1$.   5. $m\in (0,1)$   6. $-1 < d_k \leq 0, \, -1 < f_k \leq 0$ are constants for $k=1,2, \cdots,$ such that there exists a positive integer $q$ $t_{k+q}=t_k + \omega$. I am trying to make the right choice for these parameters to solve this with _Mathematica_.


title: Problems with manipulating plot in a dialog

text: G'day, guys. Here is a silly problem I am struggling with.               (*file name: 1.nb*)Button["Load the program", NotebookEvaluate[SystemDialogInput["FileOpen"]], Method->"Queued"]          (*file name: 2.nb*)eff = 1 - y Abs[1 + x]/(2 + Abs[1 + x]);CreateWindow[DialogNotebook[{TabView[{"eff"->Panel[With[{expr = eff}, Manipulate[Plot[expr, {x, -5, 5}], {y, 0, 0.1, 0.001}]]]}]       The 2 lines of codes are saved in separate file. The idea is that I want to use the first file to open the 2nd file and show the manipulate output in a dialog. The 2nd file runs fine with itself from notebook (i.e. output shows fine in notbook and a dialog window). However, if I choose the 2nd file by running the 1st file, the graph doesn't update dynamically. Anyone got an idea for that? Thanks in advance!!


title: Percentage axis ticks

text: To create a percentage axe in some plot I use this `Ticks` configuration:               ListPlot[          Range[0,1,0.1]         ,Ticks->{Automatic, {#,If[IntegerPart@Mod[100#,10]==0,Row[{100#,"%"}],""]}&/@Range[0,1,0.05]}     ]      ![enter image description here](http://i.stack.imgur.com/HieCK.png) That has 3 disadvantage:   1. Non automatic range scale   2. No subticks with different sizes   3. It's nor elegant nor simple I miss a Mathematica `TicksFunction` or equivalent, so I could have better control over it. I also tried something like: `ListPlot[Quantity[Range[0, 1, 0.1], "Percent"]]` without succes. Some clue into a better approach?


title: What is the real default CellEvaluationFunction?

text: The default `CellEvaluationFunction` is `Identity` but one can quickly see that this is handled specially as an equivalent function, `(# &)`, does not work at all the same. The sole example in the documentation uses `ToExpression` but this also does not work the same, returning only the output of the last line in the Input Cell:               SetOptions[EvaluationNotebook[], CellEvaluationFunction -> ToExpression]      Now in one cell:               2 + 2     5 + 7      > >     12 >   My question is then: What is the _real_ default function that correctly handles multi-line input/output including the `Out[xx]=` cell labels?


title: Any code to add elements of a list to make a new list?

text: Is there any code built in to add the numbers in a list sequentially? For example, lets say I have a list B. B has 6 numbers like below.               B={0,6,5,4,3,2}      I want to produce a new list C. C will add the first two numbers in the list, so 0+6 = 6. That will be put into the new list C. Then the next term in the list C will be 0+6+5 = 11. And then 0+6+5+4 = 15, and so on. The result of the new list C will be as below.               C={6,11,15,18,20}      Any help is appreciated. Thanks.


I'm getting this error in Wolfram Workbench, and I don't know where it comes from: "Syntax error: newline character interpreted as multiplication" This error happens around an If[] function, in a custom function, defined as               myFunction[x_,y_] :=          Module[{variable},          (* some code *)     variable = 1     If[x==2,x=3,x=4]          (* some code *)          ]      No matter where I try to put my If[] function inside myFunction[], I get this error, but I don't get it if I put it outside myFunction[]. Maybe the problem comes from the Module[] function, as when I take it out the error disappears. Any idea what's going on?


title: Constructing a list that includes a leading zero (01,02,03 ... 55, 56, etc.)

text: I am constructing a list for use in file operations. Say I want to create `list` which contains values from 01 to 87. The kicker is that the first nine integers need to keep a value of `01, 02, 03, 04, 05, 06, 07, 08, 09, 10, 11.. etc.` So if I create `list=Range[87]` I generate a list of integers. Whenever I'm trying to then replace values with '01', etc. they become strings. If they are integers, they revert back to the single value. Ideally, I need a list that I can feed through `Do[]`. Is there a quick way to create a list that uses the two digit form the whole time? FWIW, the specific application of this that I have a list of files i.e. `{1997-01-page100.jpeg, 1997-01-page101,jpeg, <<33 379>>, 1997-87-page9.jpeg, 1997-87-page10,jpeg}`. I'm using the following command to generate specific lists of the sub-file using `importlist=FileNames["1997-"<>ToString[x]<>"-*.jpeg"]` but am being thrown off by the zero in front of the single-digit file number.


title: What changes in Dynamic have occurred from Mathematica 7 to Mathematica 8 that might induce incompatibilities?

text: I have a function that outputs a Manipulate (say                FAnimated[params_]:= Manipulate[F[params,t], {Time, Animator[t]}]      where `F` is some function that takes `params` and `t` and makes a plot. Suppose I make some sliders that control various of the params:               Slider[Dynamic[p1], {3, 7}]      and call               Dynamic[FAnimated[{p1,p2,etc}]].      The hope is that the sliders should let me vary parameters in my plots while they run their animation. In Mathematica 7, this works great, and runs very smoothly. In Mathematica 8, the exact same code no longer lets me drag the sliders while the animation is running, and also doesn't let me drag the sliders smoothly even when the animation is paused! So the question is, "What changes (either in `Dynamic` functionality or in some options that might possibly be set different by default) might cause this code to break?"


title: Efficient conditional division

text: I was curious about the following little problem. Suppose $n \in \mathbb{N}$ is some large number and $p$ is a (possibly large) prime, and I want to update $n$ according to $$n := \begin{cases} n/p, & (n,p) = p; \\\ n, & (n,p) = 1.\end{cases}$$ In other words, I want to divide by $p$ if $p$ divides $n$, and not divide by $p$ if $p$ does not divide $n$. Also, I would like to know if a division took place. Of course I could implement this with something like               If[Divisible[n,p],        n = Quotient[n,p];        Print["A division took place!"];     ];      However, this feels inefficient, since to check divisibility the machine may have already computed $n/p$. Would there be a more efficient way to do this in Mathematica? i.e. A way that uses only $1$ division if $p$ divides $n$ and at most $1$ division if $p$ does not divide $n$? Also, I'm curious which Mathematica functions would be most suited for these purposes in general. Is there any difference in performance between e.g. `Divisible[n,p]`, `Mod[n,p] == 0`, `IntegerQ[n/p]`, `FractionalPart[n/p] == 0`, `CoprimeQ[n,p]`, `GCD[n,p] == p`? Thanks in advance.


title: Can graph interpolating function, but not use it. Any thoughts?

text: When solving a system of differential equations numerically, I get two interpolating functions. I can graph both of these, but when I try to evaluate either of the functions at a point, it just spits out:               {x1[5]}      `x1[t]` is the interpolating function, and 5 is the number I chose to try to evaluate it at. I set `prsol=NDSolve[blah blah]` To evaluate it, I've tried:               Evaluate[x1[5]/.prsol]      and               one[t_]:=x1[t]/.prsol     one[5]      Both of these spit out `{x1[5]}` which is contrary to what the F1 Help says it should do.               P0h = 1.856*10^7;     V0h = .826063;     T0 = 299.817     g = 9.81;     Aac=.16;     Abe=.025;     Vr=1.067;     Vst=.006306;     Vh[x_] := Vr - x*Aac;     Vl[x_] := Vst + x*Abe;     m1 = 290000;     m4 = 9021;     k3 = 8*10^6;     R = 8.3144621;     Tc = 126.1;     Pc = 3.394*10^6;     ω = .040;     b = .07779607*(R*Tc)/Pc;     Trr[T_] := T/Tc;     κ = .37464 + 1.54226 ω - .26993 ω^2     α[T_] := (1 + κ (1 - Sqrt[Trr[T]]))^2     ac = .45723553*(R^2*Tc^2)/Pc;     a[T_] := ac*α[T]     nh = (P0h*V0h)/(R*T0);     nl = (P0l*V0l)/(R*T0);     Vmolh[x_] := Vh[x]/nh     Vmoll[x_] := Vl[x]/nl          Ppr[x_] := (R*T0)/(Vmolh[x] - b) -        a[T0]/(Vmolh[x]^2 + 2*b*Vmolh[x] - b^2)          Fpr[x_] :=       Aac ((R*T0)/(Vmolh[x] - b) -           a[T0]/(Vmolh[x]^2 + 2*b*Vmolh[x] - b^2)) -        Abe ((R*T0)/(Vmoll[x] - b) -           a[T0]/(Vmoll[x]^2 + 2*b*Vmoll[x] - b^2))          A1 = .1;     ω2 = π/4;     c = 10000;     d = 300000;          ctipmot[t_] := A1 Cos[ω2 t]          prsol = NDSolve[{-c x4'[t] -           k3 ((m4 g + Fpr[54/39.3701])/k3 + ctipmot[t] + x4[t]) +           Fpr[x1[t] - x4[t] + 54/39.3701] + m4 g ==          m4 x4''[t], -c x1'[t] - d x1'[t]*Abs[x1'[t]] -           Fpr[x1[t] - x4[t] + 54/39.3701] + m1 g == m1 x1''[t], x1[0] == 1,         x4[0] == 0, x1'[0] == 0, x4'[0] == 0}, {x1[t], x4[t]}, {t, 0,         100}]      > {{x1[t] -> InterpolatingFunction[{{0.,100.}},<>][t],x4[t] -> > InterpolatingFunction[{{0.,100.}},<>][t]}}               Plot[Evaluate[x1[t] /. prsol], {t, 0, 32}]      ^That gives me a graph               Evaluate[x1[5] /. prsol]      ^That gives me:               {x1[5]}      instead of a value.


title: How to batch run Mathematica code from Python?

text: I have this notebook:               arg1 = "input.txt";     arg2 = "output.txt";     u1 = Import[arg1, "Table"];     u2 = u1(*extremely complicated and long chain of operations*);     Export[arg2, u2]      and I would like to run it in a Python program like this:               def run(arg1,arg2):         my_mathematica_code(arg1,arg2)          for i in range(100):         run("input{}".format(i), "output{}".format(i))      How can I do this? Preferably without GUI and several threads in parallel (I can use 4 cores in Mathematica).


title: How can I draw an asymptotic bode plot?

text: When I use _Mathematica_ to draw a bode plot, it always gives the precise version, like this:               BodePlot[(10 (s + 3))/(s (s + 2) (s^2 + s + 2))]      ![enter image description here](http://i.stack.imgur.com/ZjVkf.png) How can I draw the asymptotic version, the kind where `1/(1 + jw)` is represented as a -20dB straight line but not a curve?


title: Matrix and string variable

text: There are $100$ `.dat` data files which contain numerical parameters.               m={1,2,.....,34};     n={8,10,12,.....19};     l={67,68,69,...96}       and so on. I constructed a $100*100$ matrix with the code below:               xmatrixx = Table[0, {100}, {100}];     mflatten=Flatten[m];          Do[xmatrixx[[i, j]] = mflatten[[j]], {i, 1, 1}, {j, 1, 34}]];      But I want to put `n` and `l` and the other files in the other rows of the same matrix `xmatrixx`. How can I change the `mflatten` for obtaining `n` and `l`? All the Best


title: How can I get the solution of complicated implicit function?

text: The question is what is the method to solve the implicit function has real and imaginary number. For example, The function is $$F(x,y)=(-I*x + 2*y^2)^2 + x^2 - 4*y^4*Sqrt[1 - I*x/(y^2)]$$ Although I used a `NSolve` and `Solve` I can't get the solution. (It is too complicated to solve by those command i think.) What is the value y? (y may be a complex number when x is given) And is it possible to draw a graph of `Real[y] - x` and `Imaginary[y] - x`?


title: ParametricNDSolve and Constraint satisfaction problem

text: I am trying to find the values of two parameters that allow for a specific result of my differential equations, under different initial conditions. Given a system of equations:               TestSystem = {H'[t] == p - H[t] (k1 sh[t]), H[0] == H0,        sh'[t] == -k1 H[t] sh[t] + k2 ss[t], sh[0] == sh0,             st == sh[t] + ss[t]};      I use `ParametricNDSolve` to compute the values of the dependent variables, while leaving `k1` and `k2` unspecified:               TestSolver[H0e_, st0_] :=      ParametricNDSolve[       TestSystem /. p -> 0 /. {H0 -> H0e}           /. sh0 -> st0 /. st -> st0         /. k2 -> k2estim /. k1 -> k1estim, {h, sh, st}, {t, 0,         300}, {k1estim, k2estim}];      The conditions that I want to evaluate are then specified by `m,d,hh`, at `t=300`:               m = 100 sh[k1estim, k2estim][300]/st;     d = 100 ss[k1estim, k2estim][300]/st;     hh = h[k1estim, k2estim][300];          cond1a[k1estim_, k2estim_] :=       Evaluate[m /. st -> 5.*^-6 /. TestSolver[1.*^-5, 5.*^-6]]     cond1b[k1estim_, k2estim_] :=       Evaluate[d /. st -> 5.*^-6 /. TestSolver[1.*^-5, 5.*^-6]]     cond1c[k1estim_, k2estim_] :=       Evaluate[hh /. TestSolver[1.*^-5, 5.*^-6]]          cond2a[k1estim_, k2estim_] :=       Evaluate[m /. st -> 5.*^-6 /. TestSolver[5.*^-5, 5.*^-6]]     cond2b[k1estim_, k2estim_] :=       Evaluate[d /. st -> 5.*^-6 /. TestSolver[5.*^-5, 5.*^-6]]     cond2c[k1estim_, k2estim_] :=       Evaluate[hh /. TestSolver[5.*^-5, 5.*^-6]]      Now, say that I wanted to find the following               1 < cond1a < 100,     0 < cond1b < 1,     1.*^-8 < cond1c < 1*^-5,     1 < cond2a < 50,     50 < cond2b < 100,     1.*^-8 < cond2c < 1*^-5,      How could I find the values of `k1` and `k2` that would satisfy such a thing? I was considering using `NMinimize` like this:               NMinimize[{1/cond1a+cond1b+cond1c+cond2a+1/cond2b+cond2c,cons},{k1,k2}]      (where `cons` would be the required constrains). Yet, this does not seem like an efficient and plausible way to solve this problem. Does anyone have any ideas on how I could approach this problem? Given that this is a multi-objective optimization problem, how would you suggest that I could specify the objective function? Thank you so much!


title: How to substitute $x$ in a expression $f(x)$ but not calculate the value of the $f(x)$ at the point $x$?

text: I want to subtitute $x = 1$ into the expression $$f(x)=x^2 + 3 (x^2 - m^2) + \sqrt{x^2 + 1}$$ but I don't calculate the value $f(1)$, that is mean, I want $$1^2 + 3(1^2-m^2)+\sqrt{1^2 + 1}.$$ If I tried               f[x_] := x^2 + 3 (x^2 - m^2) + Sqrt[x^2 + 1];     f[1]      I receive `1 + Sqrt[2] + 3 (1 - m^2)` How do I tell _Mathematica_ to do that?


title: Comparing Positions of Elements in a List

text: How can I check if `Position[i] > Position[j]` for two arbitrary elements in a given list? The problem is that `Position[i]` returns `{{k}}` instead of just `k`, and I don't know how to compare `{{k}}` and `{{l}}`.


title: Finding the maximum order derivative of an expression

text: I have an expression that is the the sum and product of undefined derivatives and their derivatives. e.g. `fg+ Dt[f,t]g +Dt[g,{t,2}]f+` What I would like to find is the maximum order of the derivatives, and so 2 would be returned for the above example. My first attempts use `Cases` to make lists of the order of derivatives. This would send `Dt[x_,t]->1` and `Dt[x_,{t,n_}]->n`. However for some reason `Dt[f,t]` and `Dt[f,{t,2}]` fail to match with their respective patterns and don't trigger the rule. If, in the rules, the `x_` is replaced with `f` the rules will trigger but only for `f` and no other function. Is there a way to do this search without a strict list of the available functions.


title: Programming a numerical method in the functional style

text: I am new to _Mathematica_ and I would like to learn a bit more about functional programming. At the moment I have assignments like programming different numerical methods (for integration: trapezoidal rule, simpson's rule ,.. , for differential equations: euler's method, midpoint Runge-Kutta, ..). I've been implementing all of the methods with procedural programing, with `While` loops and `For` loops, but now with Eluer's method I've started using functional programming and have used `Table` and `FoldList` instead of a `For` loop. It comes together quite nicely. Here is an example:               euler[f_Symbol, y0_, {tz_, tk_}, h0_: 0.1] :=       Module[{t0 = tz, h = h0, n = Abs[tz - tk]/h0, Y, T, y},         y = FoldList[#1 + h*f[#1, #2] &, y0, Table[t0 + i*h, {i, 0, n - 1}]];         T = Table[t0 + i*h, {i, 0, n + 1}];         Y = Table[N[{T[[j]], y[[j]]}], {j, 1, n + 1}];         ListPlot[Y, Joined -> True, PlotRange -> All]       ]      I would like to do something similar with the midpoint Runge-Kutta method, but I don't know how. Here is the procedural code:               midpointRK[f_, y0_, ta_, tb_, h0_: 0.1] :=       Module[{tA = ta, tB = tb, h = h0, n, Y, g, y},         n = Abs[tB - tA]/h;         y[0] = y0;         For[i = 0, i <= n, i++,           y[i + 0.5] = y[i] + (h/2)*f[y[i], tA + i*h];           g = f[y[i + 0.5], tA + i*h*0.5];           y[i + 1] = y[i] + h*g;           Y = Table[{N[tA + j*h], y[j]}, {j, 0, n}]         ];         ListPlot[Y, Joined -> True, PlotRange -> All]          ]      Can someone please help me in this regard and explain to me a little, what each command does?


title: How can I detect and modify one area of an image while preserving the rest

text: I have an image with a big black disk bounded by a white boundary in its center. The image contains some other stuff outside of the spot and its boundary. This stuff has colors which include both black and white. I wonder whether _Mathematica_ can do either of the following:   1. Find the position of all the points in the central disk   2. Change the color of the central disk to red while preserving all the other parts of the image.  ![enter image description here](http://i.stack.imgur.com/kJVmP.png)


title: How to convert rgb to YCbCr?

text: I looked at `ColorConvert`, but haven't found a relevant option. Is there a built-in function do this?


title: Improving a function that filters each table from a list of tables

text: I have written the following function which can take in a list of tables as an arguement, and then filters each one. The function I wrote is below, and the commands I combined are listed below the function in case it is not clear what the steps are (sometimes I get confused if things are all strung together). Can this be cleaned up to be more concise, and readable?               Function[myTable, {tunableORRData =      Replace[myTable , a_ /; Abs[a] >= 10 -> 0, {2}] - ORR /.         targets; d = Sign@Diagonal@tunableORRData;         returnValue[row_, col_, myTable] := myTable[[row, col]] ;         tunablePairs =          Position[(tunableORRData\[Transpose] d)\[Transpose], _?          Negative]  /. Reverse /@ elementList;         pairEnergy =          returnValue[#1, #2] & @@ {#1, #2} & @@@           Position[(tunableORRData\[Transpose] d)\[Transpose], _?            Negative] ;         tempVar = Partition[Riffle[tunablePairs, pairEnergy], 2];          filterImposs = Replace[tempVar, a_ /; Abs[a] >= 5 -> 0 , {2}] ;          filterImposs = Replace[filterImposs, a_ /; a > 1 -> 0 , {2}] ;          Cases[           filterImposs, {x_, y_} /; y != 0]}][#] & /@ {smallBindingTable,            bigBindingTable}      Targets is used to set any value to the table equal to zero that is exactly the target. element list is used to map elements onto rows and columns of the tables.               targets = {ORR -> -1.51}     elementList = {"Sc" -> 1, "Ti" -> 2, "V" -> 3, "Cr" -> 4, "Mn" -> 5,       "Fe" -> 6, "Co" -> 7, "Ni" -> 8, "Cu" -> 9, "Zn" -> 10, "Y" -> 11,       "Zr" -> 12, "Nb" -> 13, "Mo" -> 14, "Tc" -> 15, "Ru" -> 16,       "Rh" -> 17, "Pd" -> 18, "Ag" -> 19, "Hf" -> 20, "Ta" -> 21,       "W" -> 22, "Re" -> 23, "Os" -> 24, "Ir" -> 25, "Pt" -> 26,       "Au" -> 27};      Leave only those elements that are within 0.05 of the target non zero               idealORRSmall =      Replace[smallBindingTable - (ORR) /. targets,       a_ /; Abs[a] >= 0.05 -> 0, {2}] ;      Get those elements and their positions. Then map a rule replacement to put the element names in, and the table name (to make the output human readable) Then filter the same table broadly for a different search.               tunableORRData =        Replace[smallBindingTable , a_ /; Abs[a] >= 10 -> 0, {2}] - ORR /.         targets; (*set any spurious values to zero *)      Here we filter to see if the diagonal element is on opposite sides of the target. IE (opposite in sign from any member of the column)               d = Sign@Diagonal@tunableORRData;     returnValue[row_, col_] :=       smallBindingTable[[row,         col]] ; (* lets us take row, col data from theidealORR data and \     get the energy from the small binding table*)          tunablePairs =       Position[(tunableORRData\[Transpose] d)\[Transpose], _?Negative]  /.        Reverse /@         elementList; (* get tunable pairs, map on element names *)          pairEnergy =       returnValue[#1, #2] & @@ {#1, #2} & @@@        Position[(tunableORRData\[Transpose] d)\[Transpose], _?Negative] ;     tempVar = Partition[Riffle[tunablePairs, pairEnergy], 2];          filterImposs = Replace[tempVar, a_ /; Abs[a] >= 5 -> 0 , {2}] ;     filterImposs = Replace[filterImposs, a_ /; a > 0 -> 0 , {2}] ;          Cases[filterImposs, {x_, y_} /; y != 0]      I have included a sample data table needed to make this work in a txt file. Here is the link:smallBindingTable Below is the output from filtering the smallBindingTable. I would like to do the same thing with a list of such tables.               {{{"Sc", "Pd"}, -1.47187}, {{"Sc", "Ag"}, -1.69071}, {{"Ti", "Mn"},       0.860106}, {{"Ti", "Au"}, -1.19295}, {{"V",       "Au"}, -1.10513}, {{"Cr", "Au"}, -0.131714}, {{"Mn",       "Ag"}, -2.51796}, {{"Mn", "Pt"}, -1.29537}, {{"Mn",       "Au"}, -1.21843}, {{"Fe", "Ag"}, -1.81939}, {{"Fe",       "Au"}, -1.44088}, {{"Co", "Au"}, -1.30395}, {{"Ni",       "Au"}, -0.671956}, {{"Zn", "Ag"}, -2.13397}, {{"Y",       "W"}, -0.189288}, {{"Nb", "Au"}, -1.087}, {{"Mo",       "Au"}, -1.1764}, {{"Tc", "Au"}, -0.922638}, {{"Ru",       "Au"}, -0.787955}, {{"Rh", "Au"}, -0.482705}, {{"Pd",       "Au"}, -0.379127}, {{"Ta", "Ag"}, -2.3501}, {{"Ta",       "Au"}, -1.12029}, {{"W", "Au"}, -1.16135}, {{"Re",       "Pt"}, -1.22974}, {{"Re", "Au"}, -1.04512}, {{"Os",       "Au"}, -0.76104}, {{"Ir", "Au"}, -0.578992}, {{"Pt",       "Au"}, -0.593863}}


title: Multi-table SQLSelects

text: I'm trying to figure out how these statements really work, and am having some trouble with the documentation (and the software!). Using the example publishers database (Mathematica 8), the following works fine (as documented):               SQLSelect[conn,{"TITLES","ROYSCHED"},{{"TITLES","TITLE"},         {"TITLES","TITLE_ID"},{"ROYSCHED","TITLE_ID"},{"ROYSCHED","ROYALTY"}},       "MaxRows"->5]      Now, there is no column `ROYALTY` in the `TITLES` table, so perhaps I should be able to write               SQLSelect[conn,{"TITLES","ROYSCHED"},{{"TITLES","TITLE"},         {"TITLES","TITLE_ID"},{"ROYSCHED","TITLE_ID"},"ROYALTY"},       "MaxRows"->5]      But this results in the error               DatabaseLink`SQL`Private`SQLValue::illegalvalue:          The value {{TITLES,TITLE},{TITLES,TITLE_ID},{ROYSCHED,TITLE_ID},ROYALTY}          cannot be converted to a value in an SQL statement.      So apparently either none or all of the entries in the columns list must themselves be lists. Is this right, or a bug, or what? So assuming that this is correct, and all the columns must be lists in any multi-table inner join, how can I then apply a function such as `COALESCE` to the results to cast NULLs to zero? I tried `{"TITLES","COALESCE(TITLE_ID,0) AS TITLE_ID"}`, which failed because the string passed to the JDBC driver was `TITLES.COALESCE(TITLE_ID,0) AS TITLE_ID`, which is of course incorrect. At the risk of asking too many questions in one question, let me add: I'm perfectly OK with simply returning the rows I want, with whatever NULLs may be in them. But can someone tell me the right way to turn all those NULLs to zeros once the result set is a Mathematica list?


title: PlotLegends from a List

text: I want to generate plots showing dependence of a function on another parameter, and automatically add correct PlotLegend that says which graph is which parameter.               f[x_, a_] := x + a ; Params = {1, 2};     Plot[Evaluate[f[x,#] & /@ Params], {x, 0, 1},           PlotLegends -> {StringForm["a = ``", #] & /@ Params }]      But this doesnt work, Mathematica complains that functions and legend are not the same shape.


title: Web deployment of document containing URLSave function

text: Is it possible to deploy a _Mathematica_ notebook which contains a dynamic `URLSave` function? My code has the following structure:               web = URLSave["http://www.wolfram.com/mathematica/", Close[OpenTemporary[]]]     web = Import[web, {"GZIP", "Text"}]     Dynamic[web, UpdateInterval -> 20, TrackedSymbols -> {}, SynchronousUpdating -> False]      The problem is when I start to deploy it, in the last step after I indicate in which folder to save the CDF file, then this error appears: > Select::normal: Nonatomic expression expected at position 1 in > Select[C:\Users\Morry\Dropbox\MEPS\ALISalehi\HTML\146\Untitled.cdf,StringMatchQ[#1,*.png]&].


title: $\LaTeX$ and Mathematica

text: I quite often would like to draw graphics in my $\LaTeX$ documents using Mathematica. I have encountered three problems. I would like to know if there are any workarounds to these problems I would like to make my graphics homogeneous with my document. That means that I would like to use the same font in the graphics (labels for axis etc) as the main text. Mathematica does not support Computer Modern. I found a workaround using PSFrag, saving graphics as EPS. It is possible using PSfrag to rename the text in the graphic into $\LaTeX$ code. A big downside is that this method does not allow me to use pdflatex. Many other packages (hyperlink) therefore do not work. `Graphics3D` objects are extremely big. If I save it using a bitmap, the picture usually becomes horrible. I often would like to use transparency. If I use `Opacity` to make some part of the graphic transparent, the exported file in Mathematica is horrible.


After executing MessageDialog[] how can I wait for it to close before continuing? I want to popup up a message and wait for the user to hit OK before continuing. Is there a mechanism for this?


title: Choosing a position in an array with uniform random probability conditioned on it having a specific element

text: I have an array `exArray` with two types of integer elements: `1` and `2`. I'd like to choose a position $k$ in the array with uniform random probability conditioned on it having an element with the value `1`. Is there a one-liner to do this? Right now I'm using `RandomInteger[{1,arrayLength}]`, which is a ridiculous way to proceed.


title: Solve function doesn't work

text: I have the following system of equations that I need to solve:               ArrayOfEquations={-(u/10) + 0.4 v[1] - 0.16 v[2] + 0.064 v[3] - 0.0256 v[4] +      0.01024 v[5] - 0.004096 v[6] + 0.0016384 v[7] - 0.00065536 v[8] +      0.000262144 v[9] == g[1],      -0.4 v[1] + 0.32 v[2] - 0.192 v[3] + 0.1024 v[4] -      0.0512 v[5] + 0.024576 v[6] - 0.0114688 v[7] + 0.00524288 v[8] -      0.0023593 v[9] == g[2],      -0.16 v[2] + 0.192 v[3] - 0.1536 v[4] + 0.1024 v[5] -      0.06144 v[6] + 0.0344064 v[7] - 0.0183501 v[8] + 0.00943718 v[9] ==g[3],      -0.064 v[3] + 0.1024 v[4] - 0.1024 v[5] + 0.08192 v[6] -      0.057344 v[7] + 0.0367002 v[8] - 0.0220201 v[9] == g[4],      -0.0256 v[4] + 0.0512 v[5] - 0.06144 v[6] + 0.057344 v[7] -      0.0458752 v[8] + 0.0330301 v[9] == g[5],      -0.01024 v[5] + 0.024576 v[6] - 0.0344064 v[7] +      0.0367002 v[8] - 0.0330301 v[9] == g[6],      -0.004096 v[6] + 0.0114688 v[7] - 0.0183501 v[8] +      0.0220201 v[9] == g[7],      -0.0016384 v[7] + 0.00524288 v[8] - 0.00943718 v[9] == g[8],      -0.00065536 v[8] + 0.0023593 v[9] == g[9],      -0.000262144 v[9] == g[10]}      I want to solve it for variables $u$ and $v[i]$, $i=1,..9$. When I call               Solve[ArrayOfEquations, Join[{u}, Array[v, 9]]]      I get {} as an output. Why doesn't it work? I understand that {} in general means no solution. However, I tried it with some numbers substituted for $g[i]$ and got a numerical output. Can't understand what goes wrong if I keep the $g[i]$ symbolically Thanks for any help!


title: Band Intensity Quantification

text: I would like to measure the density of my western blot band by measuring the area under the peak from this plot. More specifically, my image will have several lanes with a pattern of bands in each lane (see link below). I would like to select an area and then measure the intensity of the bands as scanned from top to bottom. I would do this for each lane by selecting same area (essentially move my rectangle over to the new lane) and then plot the density vs length of gel (down direction). I would then like to integrate the area under each curve to get a value which can then be compared to the pattern and density (increased or decreased) of the neighboring lane. Is this possible to do in _Mathematica_ and, if so, how can I do it. I'm very new to _Mathematica_ and any help would greatly be appreciated. ![http://i2.wp.com/homemadeprojects.org/wp- content/uploads/2014/02/gel.png](http://i.stack.imgur.com/HnRBX.png) * * * _Keywords_ : line, lines, spectrum, spectral


title: How can I format a table to look nice?

text: I have this code, but i would like to see the results well formated, How can I obtain a well formatted data with aligment?, similar to this one ![http://img854.imageshack.us/img854/9938/snap2002.png](http://i.stack.imgur.com/W62Pm.png) insted what i obatin, you can see here ![http://img856.imageshack.us/img856/2972/codeandtalbamalformatte.png](http://i.stack.imgur.com/GBQ7W.png) * * * Here is the code, where i want to obtain a well formated results:               itmax = 10;     tolerancia = 0.2;     alfa = (-1 + Sqrt[5])/2;     ttheta[Lambda_] = Lambda^2 + 2*Lambda;     a = -3;     b = 5;     La = a + (1 - alfa)*(b - a) // N;     Mu = a + alfa*(b - a) // N;     tthetaLa = ttheta[La];     tLacalc = "*";     tthetaMu = ttheta[Mu];     tMucalc = "*";     Print["iteración_k   |   a_k    b_k   |   Lambda_k     Mu_k   |   \     theta(Lambda_k)  theta(Mu_k)"];     Do[       L = b - a;       Print[i, "   |   ", a, "     ", b, "   |  ", La, "     ", Mu,         "   |  ", tthetaLa, tLacalc "      ", tthetaMu, tMucalc ];       If[tthetaLa > tthetaMu,        If[b - a < tolerancia, Break[]];             a = La;        La = Mu;        Mu = a + alfa (b - a);        tthetaLa = tthetaMu;        tLacalc = " ";        tthetaMu = ttheta[Mu];        tMucalc = " ";,             b = Mu;        Mu = La;        La = a + (1 - alfa) (b - a);        tthetaMu = tthetaLa;        tMucalc = " ";        tthetaLa = ttheta[La];        tLacalc = "*";        ]       , {i, 1, itmax}       ];     solucionestimada = (Mu + La)/2;     valorestimado = ttheta[solucionestimada];     Print["The stimated solution is ", solucionestimada,        " where functions value is ", valorestimado];