# Tag Info

14

The mesh seems to be fine and you can see that it is by doing: region = ImplicitRegion[! (Norm[{x, y, z}] < 1), {{x, -5, 5}, {y, -5, 5}, {z, 0, 5}}]; m = DiscretizeRegion[region, {{-2, 2}, {-2, 2}, {0, 1}}] To view as wireframe you can do: Needs["NDSolveFEM"] mesh = ToElementMesh[m] // Quiet; Then: Show[mesh["Wireframe"]] If you want to ...

10

Note that sin is not a function in Mathematica, instead use Sin with a capital S and functions use square brackets to hold their arguments. So Sin[x] instead of sin(x) or even sinx as you have it. You can use NestWhile or NestWhileList (for the list of values) also for this problem: f[x_] := x - 0.8 - 0.2 Sin[x] newtonsMethod[foo_, k_, s_: 0.0001] := ...

10

Here are a few additions to @RunnyKine suggestions. If you are ever in doubt about the quality of a mesh (an ElementMesh to be exact) you can query the mesh. Needs["NDSolveFEM"] region = ImplicitRegion[! (Norm[{x, y, z}] < 1), {{x, -5, 5}, {y, -5, 5}, {z, 0, 5}}]; mesh = ToElementMesh[region]; Min[mesh["Quality"]] 0.004439742441262357 So the ...

7

With eqn[{k_, r_, H0_, P0_}] := {H'[t] == r (1 - H[t]/k) - d H[t] P[t], P'[t] == -s P[t] + e H[t] P[t], H[0] == H0, P[0] == P0} d = 0.01; s = 0.3; e = 0.02; I would define one simulation as sim := Module[ {k = RandomVariate[NormalDistribution[150, 20]], r = RandomVariate[NormalDistribution[0.4, 0.003]], H0 = RandomVariate[UniformDistribution[{50, ...

6

There is a difficulty with the statement of the problem. Generally the problem can be solved as shown below. In this case there is a stipulation that $1 < x < 5$ and $1 < y < 5$. Unfortunately the solution to the system does not satisfy these constraints (also shown below). If we agree to use only numerical techniques and pretend that Solve ...

6

Here is another way to look at the OP's plot, which came to me after reading Daniel Lichtblau's comment under the question. His comment is worth emphasizing, which I think the following will show. (Mr.Wizard alludes to the issue in his remark on validity near the end of his answer.) Let's say that machine precision is binary64, which has 53 bits of ...

6

The reason is that you're using machine-precision input, so the result you get will always be machine precision. N[LegendreP[5, 0.1], 20] 0.17882875 Precision /@ {0.1, 0.17882875} {MachinePrecision, MachinePrecision} To get result in arbitrary precision you can use exact input or non-machine precision input (which is what happened when you ...

6

I am not getting the same results in version 10. Also I see no problem with the results already indicated, and I'll say a bit about that as we proceed. Here is what I obtain in version 10. smatrix = {{1 - 2.96392/u2, 0.0000196744/u2}, {1. - 2.96392/u2, -2.11737*10^-10 + 0.0000196746/u2}}; det = Det[smatrix]; sols = Solve[det == 0, u2] During ...

5

T[t_] := {(-Sin[t])/(Sqrt[1 + Cos[t/2]^2]), (Cos[t])/(Sqrt[1 + Cos[t/2]^2]), (Cos[t/2])/(Sqrt[1 + Cos[t/2]^2])} B[t_] := Rationalize@{(-0.5 Cos[t] Sin[t/2] + Cos[t/2] Sin[t]) Sqrt[1.625 + 0.375 Cos[t]], (-Cos[t/2] Cos[t] - 0.5 Sin[t/2] Sin[t])/Sqrt[1.625 + 0.375 Cos[t]], 1/Sqrt[1.625 + ...

5

The sorting (ordering) done by Union is different for different forms of expressions, e.g., analytic versus numeric expressions for a number. Union[{2., (Sqrt[5] + 1)/2}] {2., 1/2 (1 + Sqrt[5])} % // N {2., 1.61803} Union[{2., (Sqrt[5] + 1.0)/2}] {1.61803, 2.} SortBy[{2., (Sqrt[5] + 1)/2}, N] {1/2 (1 + Sqrt[5]), 2.}

4

eldo already posed and deleted something like this, but I think it works well in many cases: poly = Expand @ FromDigits[RandomComplex[2 + 2 I, 6], x]; poly /. n_?NumberQ :> Round[n] (1 + I) + (1 + I) x + I x^2 + (1 + I) x^3 + x^4 + (1 + I) x^5 Note that this will round exponents as well. More robust is the method of rhermans, which might also be ...

4

You don't give an example so I create my own: poly = Expand@FromDigits[RandomComplex[2 + 2 I, 6], x] (0.117797 + 0.674094 I) + (0.980296 + 1.90575 I) x + (0.190167 + 1.68039 I) x^2 + (1.65725 + 1.83193 I) x^3 + (1.07084 + 1.19757 I) x^4 + (0.473445 + 1.37764 I) x^5 Probably you can extract your coefficients first, and then ...

3

Amplifying on Chenminqi's answer g = 1/10; func = (p^2*Sqrt[1 - (2*p)/(-g + p*(1 + p))])/(-1 + (p/g)*(p - 1)) - p^2/((1 + p/g)*(1 + p))*Sqrt[1 - (2*p)/(g + p*(1 + p))] // Simplify; fd = FunctionDomain[func, p] This is equivalent to requiring that the arguments of Sqrt be positive fd == Reduce[ Thread[ Cases[func, Sqrt[x_] -> x, Infinity] ...

3

Guard against premature evaluation of the inside minimization by putting it inside a function which won't evaluate until a numeric argument is supplied: f[x_, y_] = 1 + x^2 - y^2; fminx[y_?NumericQ] := NMinValue[f[x, y], {x}]; FindMaximum[fminx[y], {y, 1.}] // AbsoluteTiming (* {0.515133, {1., {y -> -7.45058*10^-9}}} *) NMaximize[fminx[y], {y}] // ...

2

To get the list that you want f[x_] = x - .8 - .2* Sin[x] // Simplify; NewtonsMethodList[f_, x0_, n : _Integer : 20] := NestList[# - f[#]/f'[#] &, x0, n] (nml1 = NewtonsMethodList[f, Pi // N]) // InputForm {3.141592653589793, 1.1902654422649652, 0.9692779750858744, 0.9643361576782574, 0.9643338876957006, 0.9643338876952227, ...

2

If you want to see the convergents, I recommendFixedPointList. f[x_] := x - 0.8 - 0.2 Sin[x] newtonsMethodList[f_, x0_, n_] := With[{iter = # - f[#]/D[f[#], #]}, FixedPointList[iter &, x0, n]] newtonsMethodList[f, N @ Pi, 10] {3.14159, 1.19027, 0.969278, 0.964336, 0.964334, 0.964334, 0.964334} Notice that, although I set a limit of ten ...

2

fun = x - 0.8 - 0.2 Sin[x] newton1[fun_, n_] := With[{f = fun/D[fun, x]}, Nest[# - f /. x -> # &, 2., n]] newton1[fun, 10] 0.964334 newton2[fun_, n_] := With[{f = fun/D[fun, x]}, NestList[# - f /. x -> # &, 2., n]] ListLinePlot[newton2[fun, 10], AxesOrigin -> {0, 0}, Mesh -> All, MeshStyle -> ...

2

The finite element method can be used on this problem if we make a change of variables to convert the domain $[0, \infty)$ to a finite interval. I believe only MachinePrecision is available in FEM. Since AiryAi vanishes so rapidly, it will make a precise result for a large argument difficult to obtain. Another difficulty in obtaining a precise solution is ...

2

Second answer -- OK, the first answer was hogwash (the curious can inspect the edit history). It pays sometimes to write out the equation and think about it first. The code for this one looks so complicated, but the equations basically have the form (here d = 80000) rc'[t] == 2.05594*10^-10/rc[t] - 8189.14 rc[t] + 80000. rm[t], rm'[t] == -80000. rc[t] + ...

2

@Kellen Myers comment is useful. Since your coefficients are reals you have an approximate solution Whenever a number carries a decimal point as your solution output, it is an approximate real. 0.//Head (* out *) Real As mentioned in documentation center for Real you may change an approximate real number in an exact rational number by ...

2

My first observation is that {{x -> 2.36147 - 1.11022*10^-16 I}, {x -> -2.52892 + 0. I}, {x -> 0.167449 + 0. I}} is not a set of solutions for x^3 - 5 x^2 - x + 1 == 0 This can be seen by plotting the polynomial Plot[x^3 - 5 x^2 - x + 1, {x, -1., 6.}] However, the problem of imaginary fuzz in the roots remains. Solve[x^3 - 5 x^2 - x ...

1

Thanks for updating your Question. With the new, clearer example I believe I can see the issue. Analysis The first method uses evenper on Symbolic values that are in canonical order: r1 = evenper[{a, b, c, d}] {{a, b, c, d}, {a, c, d, b}, {a, d, b, c}, {b, a, d, c}, {b, c, a, d}, {b, d, c, a}, {c, a, b, d}, {c, b, d, a}, {c, d, a, b}, {d, a, c, b}, ...

1

This is really totally normal. Your characteristic polynomial has coefficients ~10^(-10) and your solutions sometimes seem to evaluate to f(sol)= 10^(-20) instead of 0. The precision with which Mathematica will NSolve an equation may be different than the care it takes when you plug in various floating-point reals back into the function they might satisfy. ...

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