Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

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 = 
Reap[First[
  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]]];
NumberOfRoots
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*)
Clear[xi];
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
add comment

Your Answer

 
discard

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

Browse other questions tagged or ask your own question.