Mathematica Stack Exchange is a question and answer site for users of Mathematica. Join them; it only takes a minute:

Sign up
Here's how it works:
  1. Anybody can ask a question
  2. Anybody can answer
  3. The best answers are voted up and rise to the top

I need help with minimizing a functional in Mathematica.

I have a function $V(\xi)=\sum_{i=1}^\infty C_{3_i}J_0(\xi/C_i)+C_{4_i}Y_0(\xi/C_i),~~~\Sigma\leq\xi\leq\xi_1$, and want to find such $\{C_{3_i}\}$ and $\{C_{4_i}\}$ that $V(\xi)\equiv-1$. Here $\{C_i\}$ are known (see below), $J_0(x)$ and $Y_0(x)$ are Bessel functions. To do that, I take $\int\limits_\Sigma^{\xi_1}|V(\xi)+1|d\xi$ and minimize it whit respect to parameters $\{C_{3_i}\}$ and $\{C_{4_i}\}$. I don't know why, but it takes forever for Mathematica to minimize that functional (at least more than 14 hours). Probably that's because integration increases the number of summands adding some with $C_{3_i}C_{4_i}$ coefficients. Am I doing something wrong or is there a method for getting the $\{C_{3_i}\}$ and $\{C_{4_i}\}$ faster?

Important note: I take $V(\xi)\approx\sum_{i=1}^N C_{3_i}J_0(\xi/C_i)+C_{4_i}Y_0(\xi/C_i)$, where $N$ is length of $\{C_i\}$, and $\{C_i\}$ are eigenvalues of Sturm-Liouville problem (see the code for detailed information). $V(\xi)\equiv-1$ is actually obtained from initial condition for the initial function $F(\xi,\tau)=V(\xi)U(\tau)$.

xi1 = 10;
Sigma = 1;
nn = 1;
f3[z_] = BesselJ[1, nn*xi1/z]*BesselY[0, nn*Sigma/z] - 
BesselJ[0, nn*Sigma/z]*BesselY[1, nn*xi1/z];
Plot[f3[z], {z, 0, 10*nn}]
NumberOfRoots = 0;
Timing[roots3 = 
  NDSolve[{x'[z] == f3'[z], x[1] == f3[1]}, {x}, {z, 1000, 0.1}, 
   Method -> {"EventLocator", "Event" -> x[z], 
     "EventAction" :> {Sow[y /. FindRoot[f3[y], {y, z}]], 
       NumberOfRoots++}}, MaxStepSize -> 0.1, 
   MaxSteps -> 1000000]]][[2, 1]]];
Min[Abs[f3 /@ roots3]]
Max[Abs[f3 /@ roots3]]
roots3 = Sort[roots3];

The code above finds $\{C_i\}$ - roots of equation given by boundary conditions $V(\Sigma)=0$ and $\frac{dV}{d\xi}(\xi_1)=0$ $$J_1(\xi_1/C)Y_0(\Sigma/C)-J_0(\Sigma/C)Y_1(\xi_1/C)=0.$$ I use the function $$f_3(C)=J_1(n\xi_1/C)Y_0(n\Sigma/C)-J_0(n\Sigma/C)Y_1(n\xi_1/C),$$ where $n$ is a scale parameter. Bigger $n$ gives more roots. Thus, for $n=1$ we get 29 roots, $n=5$ - 143, $n=10$ - 286, etc. Then I print maximum and minimum values of $f_3(\{C_i\})$ to check if I really get zeros. For $n=1$ $\min f_3(\{C_i\})\sim2*10^{-19}$, $\max f_3(\{C_i\})\sim 2*10^{-15}$, which is really satisfactory.

Array[CCNew, NumberOfRoots]; (*CCNew=C*)
Array[C3, NumberOfRoots];
Array[C4, NumberOfRoots];
For[i = 0, i < NumberOfRoots, CCNew[i] = roots3[[i]]/nn, i++] (*here we take the scale into account*)
InitFunction[xi_] = 
Sum[C3[i]*BesselJ[0, xi/CCNew[i]] + 
C4[i]*BesselY[0, xi/CCNew[i]], {i, 1, NumberOfRoots}]; (*InitFunction=V*)
C5 = Join[Array[C3, NumberOfRoots], Array[C4, NumberOfRoots]];
Timing[Minimize[Integrate[Abs[InitFunction[xi] + 1], {xi, Sigma, xi1}], C5]]

And here it takes forever to find $\{C_5\}=\{C_3\}\cup\{C_4\}$.

Any help with this?

share|improve this question

Your Answer


By posting your answer, you agree to the privacy policy and terms of service.

Browse other questions tagged or ask your own question.