Memperhalus hasil FindRoot

Question

Saya memiliki satu set kode yang melibatkan pencarian yang sesuai cuntuk masing-masing a(meskipun saya akan memberikan nilai ananti) dan zmenggunakan kendala toroot[a,c,z]dan kemudian menggantikannya ckembali ke ekspresi akhir functionS[a,z]. Saya juga memiliki fungsi lain yang mana ada perubahan variabel di functionSR[l,z]mana a->l-0.01.

d = 3;
zh = 1.5;
toroot[a_, c_?NumericQ, z_] := a - NIntegrate[(c z^(d + 1) x^d)/((1 - ((z x)/zh)^(d + 1)) (1 - c^2 (z x)^(2 d)))^(1/2), {x, 0, 1}, MaxRecursion -> 5, PrecisionGoal -> 4, Method -> "LocalAdaptive"]
cz[a_?NumericQ, z_?NumericQ] := c /. FindRoot[toroot[a, c, z], {c, 0.0009, 0.0000001, 10000}, WorkingPrecision -> 5]
intS[a_?NumericQ, z_?NumericQ] := NIntegrate[With[{b = z/zh}, (((-1)/(d - 1)) cz[a, z]^2 z^(2 d)) x^d ((1 - (b x)^(d + 1))/(1 - cz[a, z]^2 (z x)^(2 d)))^(1/2) - ((b^(d + 1) (d + 1))/(2 (d - 1))) x ((1 - cz[a, z]^2 (z x)^(2 d))/(1 - (b x)^(d + 1)))^(1/2) + (b^(d + 1) x)/((1 - (b x)^(d + 1)) (1 - cz[a, z]^2 (z x)^(2 d)))^(1/2)], {x, 0, 1}, MaxRecursion -> 5, PrecisionGoal -> 4, Method -> "LocalAdaptive"]
functionS[a_, z_] = ((-((1 - cz[a, z]^2 z^(2 d)) (1 - (z/zh)^(d + 1)))^(1/2)/(d - 1)) + intS[a, z] + 1)/(z^(d - 1));
functionSR[l_, z_] = Replace[functionS[a, z], a -> (l - 0.01), Infinity];

Masalah saya adalah ketika saya mencoba menemukan minimum dari functionS[a,z]dan functionSR[l,z]untuk beberapa adan l, katakanlah a=1dan l=1, itu memberi saya kesalahan. Saya pikir itu terkait dengan perilaku ckapan a=1atau l=1.

In[23]:= FindMinimum[functionS[1, z], {z, 1.2, 1.5}] // 
  Quiet // AbsoluteTiming
FindMinimum[functionSR[1, z], {z, 1.2, 1.5}] // Quiet // AbsoluteTiming

During evaluation of In[23]:= NIntegrate::ncvb: NIntegrate failed to converge to prescribed accuracy after 5 recursive bisections in x near {x} = {0.697475}. NIntegrate obtained 0.000944548 -0.00149313 I and 0.0006178735732839699` for the integral and error estimates.

During evaluation of In[23]:= NIntegrate::ncvb: NIntegrate failed to converge to prescribed accuracy after 5 recursive bisections in x near {x} = {0.697475}. NIntegrate obtained 0.000944548 -0.00149313 I and 0.0006178735732839699` for the integral and error estimates.

During evaluation of In[23]:= NIntegrate::ncvb: NIntegrate failed to converge to prescribed accuracy after 5 recursive bisections in x near {x} = {0.697475}. NIntegrate obtained 0.000949747 -0.00149122 I and 0.000620731102746343` for the integral and error estimates.

During evaluation of In[23]:= General::stop: Further output of NIntegrate::ncvb will be suppressed during this calculation.

During evaluation of In[23]:= FindRoot::reged: The point {1.70561} is at the edge of the search region {1.0000*10^-7,10000.} in coordinate 1 and the computed search direction points outside the region.

During evaluation of In[23]:= FindRoot::reged: The point {1.70561} is at the edge of the search region {1.0000*10^-7,10000.} in coordinate 1 and the computed search direction points outside the region.

During evaluation of In[23]:= FindRoot::reged: The point {1.70561} is at the edge of the search region {1.0000*10^-7,10000.} in coordinate 1 and the computed search direction points outside the region.

During evaluation of In[23]:= General::stop: Further output of FindRoot::reged will be suppressed during this calculation.

During evaluation of In[23]:= FindMinimum::nrnum: The function value 0.436961 -1.38189 I is not a real number at {z} = {1.2}.

During evaluation of In[23]:= FindMinimum::nrnum: The function value 0.436961 -1.38189 I is not a real number at {z} = {1.2}.

Out[23]= {0.760891, FindMinimum[functionS[1, z], {z, 1.2, 1.5}]}

During evaluation of In[23]:= NIntegrate::ncvb: NIntegrate failed to converge to prescribed accuracy after 5 recursive bisections in x near {x} = {0.699811}. NIntegrate obtained 0.00286247 -0.0000971587 I and 0.0005426332486649041` for the integral and error estimates.

During evaluation of In[23]:= NIntegrate::ncvb: NIntegrate failed to converge to prescribed accuracy after 5 recursive bisections in x near {x} = {0.699811}. NIntegrate obtained 0.00286247 -0.0000971587 I and 0.0005426332486649041` for the integral and error estimates.

During evaluation of In[23]:= NIntegrate::ncvb: NIntegrate failed to converge to prescribed accuracy after 5 recursive bisections in x near {x} = {0.699811}. NIntegrate obtained 0.00286812 -0.0000961916 I and 0.0005442259497809905` for the integral and error estimates.

During evaluation of In[23]:= General::stop: Further output of NIntegrate::ncvb will be suppressed during this calculation.

During evaluation of In[23]:= FindRoot::reged: The point {1.68855} is at the edge of the search region {1.0000*10^-7,10000.} in coordinate 1 and the computed search direction points outside the region.

During evaluation of In[23]:= FindRoot::reged: The point {1.68855} is at the edge of the search region {1.0000*10^-7,10000.} in coordinate 1 and the computed search direction points outside the region.

During evaluation of In[23]:= FindRoot::reged: The point {1.68855} is at the edge of the search region {1.0000*10^-7,10000.} in coordinate 1 and the computed search direction points outside the region.

During evaluation of In[23]:= General::stop: Further output of FindRoot::reged will be suppressed during this calculation.

During evaluation of In[23]:= FindMinimum::nrnum: The function value 0.439434 -1.36539 I is not a real number at {z} = {1.2}.

During evaluation of In[23]:= FindMinimum::nrnum: The function value 0.439434 -1.36539 I is not a real number at {z} = {1.2}.

Out[24]= {0.771827, FindMinimum[functionSR[1, z], {z, 1.2, 1.5}]}

Untuk a = 0,1, plotnya jauh lebih mulus

Untuk a = 1, plot berisi lebih banyak tonjolan

Apakah kode saya ditulis dengan buruk untuk diekstrak c? Apakah ada perubahan yang bisa dilakukan? Saya telah membaca di suatu tempat yang Reducejuga dapat digunakan sebagai pengganti FindRoottetapi saya masih memahaminya. Juga, apakah menggunakan LocalAdaptivesebagai metode yang NIntegratecocok di sini?

PEMBARUAN: Harap perhatikan kesalahan ketik, saya telah memperbaikinya. Di plot sebelumnya, saya menulis c=0.1dan c=1tetapi harus a=0.1dan a=1.

Ekspresi masalah saya diberikan oleh,

$$a = c z_s^{d+1}\int_0^1 dx \frac{x^d}{\sqrt{(1-(z_s/z_h)^{d+1} x^{d+1})(1-c^2 z_s^{2d} x^{2d})}} \tag{1}\label{1}$$

\begin{align} S &= \frac{1}{4 z_s^{d-1}}\Bigg(1 -\frac{\sqrt{(1-c^2 z_s^{2d})(1-b^{d+1})}}{d-1} - \frac{1}{d-1} c^2 z_s^{2d} \int^1_0 dx x^d \sqrt{\frac{(1-(b x)^{d+1})}{(1-c^2(z_s x)^{2d})}}\\ & -\frac{b^{d+1}(d+1)}{2(d-1)} \int^1_0 dx x \sqrt{\frac{(1-c^2(z_s x)^{2d})}{(1-(b x)^{d+1})}}\\ & + b^{d+1}\int^1_0 dx \frac{x}{\sqrt{(1-(b x)^{d+1})(1-c^2(z_s x)^{2d})}}\Bigg) \tag{2}\label{2} \end{align}

dimana $b=\frac{z_s}{z_h}$ dan catat itu $c=c(z_s)$( c=c[z]) meskipun dalam kode c=c[a,z],$c$ seharusnya hanya bergantung pada $z_s$( z) sejak$a$ akan ditentukan pada akhirnya.

Juga, mungkin ada cara yang lebih baik untuk menemukan desain $c$. Sebenarnya, saya bisa mendapatkan kendala lain di mana$\frac{dS}{dz_s} = 0$ (Itu karena pada akhirnya saya harus meminimalkan $S$ dengan hormat $z_s$) dan mungkin turunan dari $\eqref{1}$ dengan hormat $z_s$, agar ini bisa digunakan untuk mencari $c$?

Answer 1

Sumber NIntegratepesan kesalahan dapat dilihat dari faktor integrand x^d/Sqrt[1-c x^d z^d],, dari toroot. Karena c > z^-3, integrand adalah singular untuk beberapa titik dalam domain {x, 0, 1},. Selain itu, jika NIntegratedapat diintegrasikan melalui singularitas (dan, dengan bantuan, dapat), hasilnya adalah bilangan kompleks, yang (mungkin) tidak diinginkan. Untuk melanjutkan, ubah variabel integrasi ke xd = x^(d+1)dan terapkan yang sesuai Methoddari sini .

toroot[a_, c_?NumericQ, z_] := a - NIntegrate[((1 - xd (z /zh)^(d + 1)) 
   (1 - c^2 xd^(2 d/(d + 1))  z^(2 d)))^(-1/2), {xd, 0, 1}, Method -> {"GlobalAdaptive",
    "SingularityHandler" -> "DoubleExponential"}] (c z^(d + 1))/4

Selain itu, tentukan ulang czuntuk menggunakan garis potong Methoddan ikat pencarian cantara 0dan z^-3.

cz[a_?NumericQ, z_?NumericQ] := c /. 
    FindRoot[toroot[a, c, z], {c, .5 z^-3, .6 z^-3/2, 0, z^-3}]

(Tebakan awal, .5 z^-3dan .6 z^-3, dipilih secara sewenang-wenang.) Dengan definisi ini, czmengembalikan nilai yang benar dari c, jika ada, dan z^-3bersama dengan FindRoot::regedpesan kesalahan sebaliknya. Dengan definisi tersebut maka kedua plot pada soal dapat diperoleh dengan benar sebagai berikut. Untuk a = 1,

Plot[Check[cz[1, z], Null], {z, 1.42, zh}, AxesLabel -> {z, c}, 
    ImageSize -> Large, LabelStyle -> {15, Bold, Black}]

Checkmencegah plotting dari jarak dekat c = 1.42, di mana tidak ada solusi, meskipun tidak menghilangkan pesan kesalahan yang sesuai. Plot kedua, untuk a = .1, adalah

LogPlot[Check[cz[.1, z], Null], {z, .2, zh}, AxesLabel -> {z, c}, 
    ImageSize -> Large, LabelStyle -> {15, Bold, Black}]