Hot answers tagged

12

We can treat the variable $y$ as an element $[y_1 \colon y_2]$ of the projective line. In code, this means replacing y[x] by y1[x]/y2[x]. For an IVP $y' = f(x, y), \ y(x_0) = y_0$, we translate the initial condition as $y_1(x_0) = y_0, \ y_2(x_0) = 1$. Since the substitution yields an equation in two variables $y_1$, $y_2$, ...


11

You could try something like the following. bounded[Indeterminate] := $MaxMachineNumber; bounded[x_?NumericQ] := x This gives a very large number instead of Indeterminate so NMinimize keeps searching. NMinimize[bounded[J[{θ0, θ1, θ2}, X1, y1]], {θ0, θ1, θ2}] (* {0.203498, {θ0 -> -25.1613, θ1 -> 0.206231, θ2 -> 0.201471}} *) Incidentally, ...


9

First, in general, I would advise you not to trust numerical algorithms. If there are doubts about the outcomes then solve the same problem with different (numerical or not) methods and see do their results agree. For the integral in the question I assume you can evaluate it with several different invocations of the Monte Carlo method and compare the ...


6

You can use WhenEvent[ ] function, maybe it's not the perfect (not very accurate) solution but it works. sol = With[{eps = 10^5}, First@NDSolve[{y'[x] == y[x]^2 + 1, y[0] == 0, WhenEvent[{y[x] == eps, y[x] == -eps}, y[x] -> -y[x]]}, y, {x, -10, 10}]]; Plot[{y[x] /. sol, Tan[x]}, {x, -6, 6}, PlotLegends -> {"NDSolve", "Tan[x]"}, PlotStyle -> ...


4

Update notice 2: Found starting initial conditions that work in V10. Update notice: bbgodfrey pointed out that the built-in shooting method does not work in V10.0.1 (nor V10.0.2). Set up the shooting method by hand, so you can control the convergence. The problem is that small changes in the initial conditions cause the solutions to blow up before r ...


4

For your first integral, NIntegrate gives warnings or messages about failing to converge when not excluding the singular points. When excluding the singular points I trust the result 13.6216 + 0. I because some alternative methods agree. Monte Carlo sampling: In[76]:= NIntegrate[ 1/Abs[tk], {kx, 0, 2 Sqrt[3] Pi/3}, {ky, 0, 4 Pi/3}, Method -> ...


4

If you mean the vertical line, this is created by only one point. p = Plot[N[Im[isinh]], {e, 0, 2}, PlotRange -> Full, Exclusions -> None, WorkingPrecision -> Infinity]; list = Cases[p, Line[x_] :> x, -1]; ListPlot[list, PlotRange -> All] If you delete this point then: point = Cases[p, {x_ /; Abs[x - 1] < 0.01, y_ /; y < 0.5}, ...


4

There are two problems with your code: 1. Incorrect Usage of WhenEvent As the name indicates WhenEvent is meant to be used when you want to e.g. switch at certain events, what I think you try to do is to set a'[t] to zero for a whole period (0<=t<=0.1), but that's AFAIK not what WhenEvent can be used for directly. Of course you can reformulate your ...


3

I'm not familiar with the physics so I cannot say whether integrating over this set of 2 real dimension, in C^2, is what is wanted. I think the code below will cover the product space of the contours that are requested. ii = z1*z2/(p/z2 + (1 - p) z1 - 1)*Exp[1/z1 + z1 + 1/z2 + z2]; i1 = (ii /. {z1 -> Exp[I*t1], z2 -> Exp[I*t2]})*I*Exp[I*t1]*I* ...


3

You may try something like: raw = FinancialData["GE", All]; fraw = Flatten[raw]; data = Table[fraw[[4*i]], {i, 1, Length[raw]}];(*extracting just the prices*) model = A + B*Abs[c - x]^z; fit = FindFit[data, {model}, {A, B, c, z}, x]; modelf = Function[{t}, Evaluate[model /. fit]] Show[Plot[modelf[x], {x, 0, 12000}], ListPlot@data]


2

Reversing the order of integration produces a solution: ans= Integrate[HeavisideTheta[1 - x - y]/(x 100^2 - y (1 - y) 90^2), {y, 0, 1}, {x, 0, 1}, PrincipalValue -> True] (* Log[(100*10^(38/81))/(81*19^(19/81))]/10000 *) N[ans] (* 0.000060027526501455836 *) Solutions of this sort are what I would expect based on outlining a pencil-and-paper ...


2

I think part of the issue lies with the coordinate singularity at r=0. I made the following transformations, first x1 -> x1/KO2 , x2->x2/KS . Then u1=x1/r and u2=x2/r. I also assumed that you were solving from R0>0 to R because you couldn't solve at R0=0. I set the boundary conditions u1[0]=0, u2[0]=0 . These transformations were used because Laplacian[ ...


2

The use of NumericQ as mentioned by MarcoB and Guess who it is. in the comments seems to be important. Also, estimating the integral using a Total[Table[]] as in the code I posted in the second revision of the question makes this method computational feasible enough to solve my problem.


2

With Limit[] seems to work. Limit[Integrate[1/(t + 1 + \[Epsilon])* DiracDelta[t + 1], {t, -\[Infinity], \[Infinity]}], \[Epsilon] -> 0] $\infty$ Limit[Integrate[1/(t + 1 - \[Epsilon])* DiracDelta[t + 1], {t, -\[Infinity], \[Infinity]}], \[Epsilon] -> 0] $-\infty$ we have at the same point: $(\infty\ \text{and} -\infty) \to ...


2

First, this layer after layer of mandatory manual desktop-publishing followed by manual un-desktop-publishing and screwing around just to get greek characters on the screen and make code acceptably pretty adds layers of uncertainty and potential errors. For example, I suspect that somewhere in the editing process "curly phi" and "phi" got reversed in your ...


2

Here is a way to plot the sampling points of the integration process. (If this is what it is asked.) First using the functions Sow and Reap and the option EvaluationMonitor we gather the sampling points: Block[{x = 0.2, t = 0.2, m = 0.2}, res = Reap[ NIntegrate[ f[y] (2. Sin[(x^2 + y^2)/(8. (t - r)) + \[Pi]/ 4] Cos[(x y)/(4. (t - ...


1

Because boundary conditions are specified at z = 0 and at z = 1, ParametricNDSolve must integrate between those two endpoints. However, this does not in general prevent obtaining a solution. s = ParametricNDSolveValue[{Ω f''[z] - f'[z] - ψ f[z]^3 == 0, f[0] - Ω f'[0] == 1, f'[1] == 0}, f, {z, 0, 1}, {Ω, ψ}]; which then can be evaluated at s[Ω, ...


1

Fast and dirty, try this: Clear[\[CapitalOmega], \[Psi], sol]; n = 3; \[Psi] = 1; lst = Flatten[ Table[sol = NDSolve[{\[CapitalOmega]*f''[z] - f'[z] - \[Psi]*f[z]^n == 0, f[0] - \[CapitalOmega]*f'[0] == 1, f'[1] == 0}, f, {z, 0, 1}]; Table[ {\[CapitalOmega], z, f[z] /. sol}, {z, 0, 1, 0.1}], {\[CapitalOmega], 1, 2, 0.1}] /. ...


1

So as I posted in the Edit, the initial problem was that I was trying to integrate a function of the form 1/x around 0. In this case, as I will describe later, it corresponds to integrating around the singularity given by b(qx^2+qy^2+qz^2) + w ==0. The integral diverges on each side and naively integrating doesn't work. We need some other method, like taking ...


1

Results can be obtained with the following changes. First, delete the two asterisks at the end of the line ik=ik2 (* if evk>phi*)**, which confuse Mathematica. Second replace DSolve by First@NDSolve (and t by { t, 0, T} in NDSolve), because DSolve cannot solve these equations. Third, change T to .5*(1/fre) to avoid stiffness issues. Then, labeling ...


1

One problem is using NumericQ - see What are the most common pitfalls awaiting new users? Another problem is how many times NIntegrate is being called. To do the 2D integral in FF, f4 will be called tens of thousands of times. Computing f4 in turns does three 1D integrations, which will call f1, f2, f3 several hundred times. And there's still one more ...


1

It looks to me like the data extends to the right of the singularity. Is it acceptable to simply remove that part of the data? It doesn't look like it could possibly fit the model. FindFit[data[[;; 9500]], model, {{A, 0}, {B, 10^9}, {c, 10000}, {z, -2}}, x] (* {A -> -0.0521369, B -> 1.54292*10^10, c -> 11142.5, z -> -2.71708} *) ...



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