0

I have a (somewhat long) script which runs cleanly until the very last stage. Here, I'm using FindRoot to attempt to parameterize implicit curves in the zero set of a two variable function. The FindRoot is giving me the familiar error

"is neither a list of replacement rules nor a valid dispatch table..."

and

"is not a list of numbers with dimensions..."

Now, I know using ?NumericQ is a common solution here, and I think I've done that correctly below, but I'm still getting errors. I'm hoping someone can perhaps parse through what I'm overlooking here:

M = 1; tau = (0.5) + (0.4)*I; w1 = Pi/2; w2 = Pi*(tau)/2; inv = WeierstrassInvariants[{w1, w2}]; E2[t_] = 1 - 24*Sum[(n*Exp[2*Pi*I*(t)*n])/(1 - Exp[2*Pi*I*(t)*n]), {n, 1, 300}]; z[u_?NumericQ] = (I* M/2)*(WeierstrassZeta[u, inv] - ((1/3)*N[E2[tau], 50]*(u))); WP[x_, y_] = WeierstrassP[w1*x + w2*y, inv]; L = -(1/3)*N[E2[tau], 50]; f[x_, y_] = Re[WP[x, y] - L]; g[x_, y_] = Im[WP[x, y] - L]; ContourPlot[{f[x, y] == 0, g[x, y] == 0}, {x, 0, 2}, {y, -2, 2}] V1 = Quiet[ FindRoot[{f[x, y] == 0, g[x, y] == 0}, {x, 0.8}, {y, 0.5}, WorkingPrecision -> 50]]; V2 = Quiet[ FindRoot[{f[x, y] == 0, g[x, y] == 0}, {x, 1.2}, {y, 1.5}, WorkingPrecision -> 50]]; V3 = Quiet[ FindRoot[{f[x, y] == 0, g[x, y] == 0}, {x, 0.8}, {y, -1.5}, WorkingPrecision -> 50]]; V4 = Quiet[ FindRoot[{f[x, y] == 0, g[x, y] == 0}, {x, 1.2}, {y, -0.5}, WorkingPrecision -> 50]]; A1 = x /. V1; B1 = y /. V1; A2 = x /. V2; B2 = y /. V2; A3 = x /. V3; B3 = y /. V3; A4 = x /. V4; B4 = y /. V4; Z1 = Quiet[N[z[w1*A1 + w2*B1], 50]] Z2 = Quiet[N[z[w1*A2 + w2*B2], 50]] Z3 = Quiet[N[z[w1*A3 + w2*B3], 50]] Z4 = Quiet[N[z[w1*A4 + w2*B4], 50]] m = (Im[Z1] - Im[Z2])/((Re[Z1] - Re[Z2])); ListPlot[{{Re[Z1], Im[Z1]}, {Re[Z2], Im[Z2]}, {Re[Z3], Im[Z3]}, {Re[Z4], Im[Z4]}}, PlotRange -> {{-5, 5}, {-5, 5}}] Zed[x_?NumericQ, y_?NumericQ] = z[w1*x + w2*y]; Quiet[ContourPlot[{Im[ N[Zed[x, y]]] - (m*(Re[N[Zed[x, y]]] - Re[(M/2)]) - Im[M/2]), 0}, {x, 0, 2}, {y, -4, 4}, MeshFunctions -> {#3 &}, MeshStyle -> {{Thick, Blue}}, Mesh -> {{0}}, PlotStyle -> Directive[Orange, Opacity[0.1], Specularity[White, 30]], PlotPoints -> 75, WorkingPrecision -> 50, ClippingStyle -> None]] deriv[expr_, var_] := D[expr, var] //. {Re'[e_] :> Re[D[e, var]]/D[e, var], Im'[e_] :> Im[D[e, var]]/D[e, var], Arg'[e_] :> Im[D[e, var]/e]/D[e, var], Abs'[e_] :> (Re[e] Re[D[e, var]] + Im[e] Im[D[e, var]])/(Abs[e] D[e, var])}; Fun[x_?NumericQ, y_?NumericQ] := Im[N[Zed[x, y]]] - (m*(Re[N[Zed[x, y]]] - Re[(M/2)]) - Im[M/2]); TGrad[x_?NumericQ, y_?NumericQ] := {deriv[Fun[x, y], x], deriv[Fun[x, y], y]}; rot90[{x_, y_}] := {y, -x} step[f_, pt : {x_, y_}, pt0 : {x0_, y0_}, resolution_] := Module[{grad, grad0, t, contourPoint}, grad = TGrad[x, y] /. {x -> pt[1], y -> pt[2]}; grad0 = grad /. Thread[pt -> pt0]; contourPoint = grad0 t + pt0 /. First@FindRoot[f /. Thread[pt -> grad0 t + pt0], {t, 0}]; Sow[contourPoint]; grad = Normalize[grad /. Thread[pt -> contourPoint]]; contourPoint + rot90[grad] resolution]; result = Reap@NestList[step[Fun[x, y], {x, y}, #, .08] &, {0, 1}, 100];

Believe it or not, all that really does need to be there in some capacity. I should point out that the portion of the above code with deriv[expr_, var_] came from this very helpful post Derivative of real functions including Re and Im and the portion of the code at the end, came from this post How to find an approximate solution in a region when Solve or NSolve does not work?

Can anyone see where my syntax is going wrong leading to the errors I'm getting above?

Community
  • 1
Benighted
  • 1,327
  • 1
  • 10
  • 14
  • 2
    Try this first: 1) use inv = N[Weierstrass...];, 2) use := instead of = everywhere you have patterns on the left hand side. Now you will get some error messages that clearly indicate that Fun and TGrad are not being called in a way that works. – Marius Ladegård Meyer Jun 17 '16 at 08:10
  • What do you mean by patterns? I thought I used := in all functions. I'm not totally sure I see how this will give more transparent errors – Benighted Jun 17 '16 at 20:24
  • 1
    Patterns are the things with underscores, like x_. When your function depends on an input, you usually want to evaluate the expression for the function for each given input. Hence one uses :=. There are lots of such functions above where you use =. Just try it on a fresh kernel and see! – Marius Ladegård Meyer Jun 17 '16 at 22:03
  • @MariusLadegårdMeyer So making all those changes, I got precisely the same errors, except the code immediately stopped, as opposed to running indefinitely. Not sure if that signifies anything helpful. – Benighted Jun 18 '16 at 00:51

1 Answers1

2

There seem to be other oddities in this code. TGrad[x, y] /. {x -> pt[1], y -> pt[2]} doesn't make sense to me. You probably meant to use pt[[1]] which is shorthand for Part[pt, 1]. But even this shouldn't be needed as values of x and y should be populated from the destructuring pattern pt : {x_, y_}. More importantly deriv would appear to need a Symbol for its second argument, as it is passed to D, but TGrad is expressly defined for a numeric parameter x_?NumericQ and this is passed to deriv[Fun[x, y], x].

Mr.Wizard
  • 271,378
  • 34
  • 587
  • 1,371
  • 1
    @spietro, what Mr.Wizard observes is what I would deduce from the error messages that appear with the changes I gave you. Often, when you define functions such as TGrad to do something, but you still see TGrad explicitly in the error message, it means that it haven't evaluated the way you would expect. – Marius Ladegård Meyer Jun 18 '16 at 10:51
  • @MariusLadegårdMeyer I appreciate both comments above, thanks. So I got rid of the junk next to TGrad (I didn't know how the destructuring pattern worked) and I also added ?NumericQ to the second argument of deriv. Running the code, there are similar, but different errors saying "0 (and 1) are not valid variables." The functions which show up in the errors are Fun, TGrad, and D. So apparently I'm not calling these properly. Specifically, the error claim that my output of Fun is not a list of numbers with dimension {1}. Can you see where I'm calling this improperly? – Benighted Jun 18 '16 at 18:46