Finding roots of nontrivial equation

Question

I am trying to solve the following nontrivial equation:

    fun[z_] := (Log[Gamma[n]] + n Log[Gamma[z + 1]] + z n Log[n] - Log[Gamma[n (z + 1)]])/n
    fun[z] /. n -> 8
    m = Table[i, {i, -1, 1, 0.1}] 
    sols = Table[FindRoot[fun'[z] + 0.5 == m[[i]], {z, 1}], {i, Length[m]}]
    sols2 = Table[FindRoot[fun'[z] + 0.5 == m[[i]], {z, sols[[j]]}], {i, Length[m]}, 
        {j,   Length[sols]}]

Basically, I get roots from the equation that are stored in sols. However, I also need all values from sols as starting values for the next iteration that will be executed in sols2. However, this doesn't work and produce the following errors:

    {{FindRoot[\!\(\*SuperscriptBox["fun", "\[Prime]",MultilineFunction->None]\)[z] +0.5 == m[[i]], {z, sols[[j]]}], 
    FindRoot[\!\(\*SuperscriptBox["fun", "\[Prime]",MultilineFunction->None]\)[z] + 0.5 == m[[i]], {z, sols[[j]]}], 
    FindRoot[\!\(\*SuperscriptBox["fun", "\[Prime]", MultilineFunction->None]\)[z] + 0.5 == m[[i]], {z, sols[[j]]}], 
    ....

Any suggestion on how to fix this?

Answer 1

Let's polish up your code a bit...

With[{n = 8},
 fun[z_] := (LogGamma[n] + n LogGamma[z + 1] + z n Log[n] - LogGamma[n (z + 1)])/n;
 m = Range[-1, 1, 1/10];
 Map[(z /. First[FindRoot[fun'[z] + 1/2 == #, {z, 1}]]) &, m]
 ]

Some notes:

• You don't really need or want the last step in your question here; a far better route would have been to tweak the options of FindRoot[], as well as finding better initial estimates for your roots.

• I have elected to use the built-in LogGamma[] function here, which avoids all the problems of a construction like Log[Gamma[z]]. (There's something strange about having to deal with large numbers just to get a reasonably-sized result...)

• Range[] is as useful as Table[] for generating a list of numbers in arithmetic progression, and should be always kept in mind.

• Map[] is useful if you want to run some function through the elements of a list.

Answer 2

Your fun depends upon n and, though you plugged in n->8 at line 2, your subsequent code doesn't refer to that line. Thus, FindRoot receives input that doesn't return numerical output. That's the source of the error.

I don't see the purpose behind the two steps, as in both steps you're solving equations of the form fun'[z]+0.5==c, where c is constant. These are easily solved.

fun[n_, z_] := (Log[Gamma[n]] + n Log[Gamma[z + 1]] + 
  z n Log[n] - Log[Gamma[n (z + 1)]])/n;
f[z_] = fun[8, z];
FindRoot[f'[z] + 1/2 == -3/10, {z, 0}]

{z -> -0.321733}

Answer 3

Syntactically correct version of sols2 is

 sols2 = Table[FindRoot[fun'[z] + 0.5 == m[[i]], {z, sols[[j, 1, 2]]}], 
 {i, Length[m]}, {j, Length@sols}]

using fun[z] as in Mark's answer.

However, this does not serve your purpose for two reasons:

First, the solution of equation i from sols is used as starting value for the search for the root of equation j (for different i and j) in the second round which is not likely to improve the results. Perhaps, you meant to use a single index as in:

 sols2 = Table[FindRoot[fun'[z] + 0.5 == m[[i]], {z, sols[[i, 1, 2]]}], 
 {i, Length[m]}]

Second, considering the error message you get for both sols and sols2

FindRoot::lstol: The line search decreased the step size to within tolerance specified by AccuracyGoal and PrecisionGoal but was unable to find a sufficient decrease in the merit function. You may need more than MachinePrecision digits of working precision to meet these tolerances.

you need to consider playing with the options for FindRoot as J.M. suggested in his answer and comments.

EDIT: Inpecting the plots could suggest an explanation of the numerical results.

Plotting three examples from the 21 equations:

 Row[Plot[fun'[z] + 0.5 - m[[#]], {z, -3, 3}] & /@ {1, 15, 21}]

So, as the region z<=-1 seems hopeless, one needs to restrict root search to the appropriate region (z>-1).

Next, plotting all 21 functions

 Plot[Evaluate[(fun'[z] + 0.5 - m[[#]]) & /@ Range@Length@m], {z, -1, 5}]

we get a hint as to why the first 14 equations are solved without any issue while the last 7 are not. In particular, equation 15 simplifies to

  FindRoot[fun'[z] == 0, {z, 1}] 

Note that Limit[fun'[z] ,z->Infinity]=0 and fun'[z] is strictly increasing on (0,Infinity) (This reference and citations therein on complete monotonicity properties of Gamma and related functions). These observations together with the fact that

 fun'[0] = - 0.513416 < 0

explain why the equation fun'[z]==0 does not have a solution in [0,Infinity]. By the same reasoning, equations 15-21 (which differ from equation 15 by a constant) do not have a solution in [0,Infinity].