2

I'm trying to draw a VectorPlot of a function that admittedly probably is not very well suited for such a plot. But at the moment I just get a white sheet and I would expect something more than that.

Some definitions:

G = 6.67 10^-11;
Subscript[m, sun] = 1.988425 10^30; 
Subscript[m, 1] = 0.1 Subscript[m, sun]; Subscript[m, 2] = 0.4 Subscript[m, sun];
Subscript[x, 1] = -7 10^5/(Subscript[m, 1]/Subscript[m, 2] + 1); 
Subscript[x, 2] = -(Subscript[m, 1]/Subscript[m, 2]) Subscript[x, 1];
ω = {0, 0, Sqrt[G (Subscript[m, 1] + Subscript[m, 2])/Abs[Subscript[x, 1] - Subscript[x, 2]]^3]};

I'm sorry about the amount of subscripts (if there is an easy way to convert them I don't know about it, though Hold[expr] will render the subscripts.) Here's the plot:

ϕ[r_] := -(G Subscript[m, 1])/Norm[r - {Subscript[x, 1], 0, 0}] - (G Subscript[m, 2])/Norm[r - {Subscript[x, 2], 0, 0}] - 0.5 (ω\[Cross]r).(ω\[Cross]r)

VectorPlot[Evaluate[-Grad[ϕ[{x, y, 0}], {x, y}]], {x, -10^6, 10^6}, {y, -10^6, 10^6}]

The function has singularities in it so I've considered what's been said in this post and have also been playing with the various options, but haven't been able to get anything other than the blank, white sheet:

blank

In so far that the gradient i orthogonal to level sets, the contour plot can tell us a bit about what kind of arrows we could expect to see. The contour plot:

contourplot

There is something fishy about the Grad which makes me suspect VetorPlot, though it does not complain, does not get any values:

    N[Grad[ϕ[{x, y, 0}], {x, y}] /. {x -> 1, y -> 2}]

{-193.335 + 2.70673*10^9 Derivative[1][Abs][-139999.] + 4.22919*10^7 Derivative[1][Abs][560001.], -386.67 + 38818.9 Derivative[1][Abs][2.]}

It leaves a derivative in there.

Community
  • 1
C. E.
  • 70,533
  • 6
  • 140
  • 264
  • did you check your Grad it is resulting in {0,0}..is it ? – Pankaj Sejwal Sep 05 '13 at 16:09
  • @Blackbird No, I don't think it is but there is something else that I've noticed in my attempt to fix this that might absolutely be important, added it at the end. – C. E. Sep 05 '13 at 16:15
  • I tried to do it in steps and I am ending up with {0,0} every time. So, I thought. – Pankaj Sejwal Sep 05 '13 at 16:16
  • @Blackbird But looking at the contour plot, it seems impossible that the gradient is zero everywhere, right? – C. E. Sep 05 '13 at 16:17
  • To me it is working in only this way, but I don't know if it is relevant or not.Manipulate[ VectorPlot[ Evaluate[FunctionExpand[-Grad[[Phi][{x, y, 0}], {x, y}] /. {x -> n, y -> k}]], {x, -10.^10, 10.^10}, {y, -10.^10, 10.^10}], {n, -10.^10, 10.^10}, {k, -10.^10, 10.^10}]` – Pankaj Sejwal Sep 05 '13 at 17:38
  • @Blackbird You're not plotting the desired function at all if you replace x and y by constants. – Jens Sep 05 '13 at 18:27
  • @Jens: can't disagree, your point is very appropriate. – Pankaj Sejwal Sep 05 '13 at 18:34

1 Answers1

4

The problem is that Norm involves the Abs function whose derivative can't be symbolically simplified:

D[Abs[x], x]

(* ==> Derivative[1][Abs][x] *)

This leads to an indeterminate result that can't be plotted. So you should define the norm more explicitly, for example like this:

norm = Sqrt[#1.#1] &;

Then replace Norm by norm everywhere. This allows the gradient to be evaluated in your plot. You still have to work out the scaling of the vector arrows due to the singularity as in the answer you already linked to.

Jens
  • 97,245
  • 7
  • 213
  • 499
  • Actually on my machine I don't get Derivative[1] in answer and to solve problem of Abs FunctionExpand does the work. But why it is not working, I don't know. – Pankaj Sejwal Sep 05 '13 at 17:24
  • @Blackbird I've tried my solution in version 8 and 9, and it does the trick. In both cases, the Abs'[...] was the culprit. Not sure what exactly you did. – Jens Sep 05 '13 at 17:26
  • what I found is that if replacement is not done, if is still returning expression with unresolved Abs and hence again not valid to be plotted.FunctionExpand will work only after x,y are replaced. I am not sure if my approach is correct,but I thought to contribute. – Pankaj Sejwal Sep 05 '13 at 17:53
  • @Blackbird I appreciate anyone who tries to help me, ty! :) Thank you Jens; also for that other answer. I have a beautiful plot now... – C. E. Sep 05 '13 at 20:10