Dr. Patrick T

## 2113 Reputation

15 years, 223 days

## is it plotsetup?...

Do you mean something along the lines of :

plotsetup(ps,plotoutput=`name.ps`,plotoptions=`nocolor,portrait,noborder,axiswidth=500pt,axisheight=500`):

where the sizes would have to be tweaked to fit your needs.

## option "continuous" in call to int...

The option continuous in the call to int absolutely does the trick.

Thanks a lot Robert for your prompt and, as always, extremely helpful suggestions.

I had tried to make assumptions about the functions, as you mention, but that didn't work, as you note. I had also set _EnvAllSolutions to true, but that didn't seem to help either. The call to continuous in int does the trick! Thanks.

## I get it...

I get it: I wasn't thinking properly about what the surface associated with rhs(qd*xd) actually represents.

I can see things more clearly now (I hope it lasts).

Robert, thank you so much for your help. Very timely too.

Patrick.

## Thanks Robert for your quick...

I think I may be missing something obvious to you here. Please forgive my sloppy understanding of the math involved (I'm no mathematician, just trying to use a tool I do not master).

What I'm hoping to find is the answer to the question:  "How many solutions?" 0 solution / 1 solution / many solutions to the system of 3 equations in 3 unknowns: {xd=0, cd=0, qd=0} and {x,c,q} ?

If the equations defined flat planes I'd expect them, typically, to intersect at a single point. So my first expectation was to find a single point where the three planes intersect.

However, the equation qd=0 contains the term cd/xd which is equal to 0/0. Now you say that qd is undefined when xd=0: is that all there is to it? From that, can we conclude that there is no triplet (x,c,q) that makes {xd=0, cd=0, qd=0} true ?

When I look at the graph, I see things that look like intersections of the three planes: what am I seeing?

1) a continuum of solutions: you can pick all kinds of combinations of (x,c,q) that make {xd=0, cd=0, qd=0} true, and for each there will be different values of cd/xd=0/0. That's how I understand your statement that the limit cd/xd is path-dependent.

2) no solution at all, it's just very close with a hole in the region where it looks like the three planes intersect. But AT xd=0, the system is undefined, and therefore any "limit" (however calculated) will not belong to the domain and is therefore not properly a limit. That's how I understand your statement that qd is undefined at xd=0 together with your statement that the limit exists for some arbitrary path.

3) a combination of 1) and 2): it looks like a continuum because the system gets to within an epsilon of having a continuum of solutions, but none of these can be a solution because the system is undefined at these values.

Sorry for being thick, and thanks a lot for your help, much appreciated,

Patrick.

## Further thoughts on the validity of the ...

As I noted above, there are several candidate solutions for this 3-D system: one solution returned by fsolve and one solution obtained by computing the limit of the 0/0 for in the qd=0 equation, by Robert's procedure. The solutions are pretty close, but when I look at the picture of  { xd=0, cd=0, qd=0 } my mind fills with anxieties and doubts (solving this system is important to me at this time!)

A naive look at the picture would suggest q=b as a candidate solution. That's less than 1 percent away from the other candidate solutions. The logic would be something like: (q-b)*(f'(x)-q) -->0 at the same time as (f(x)-c) --> 0, the ratio converging to a finite limit, i.e. (q-b)*(f'(x)-q)/(f(x)-c) tends to a finite limit at xd=cd=qd=0, rather than (f'(x)-q)/(f(x)-c) tends to a finite limit at xd=cd=qd=0. I don't know, I'm just saying.

Because this is what the picture looks like:

The brown and blue planes are the xd=cd=0 planes, while the green plane here corresponds to the plane qd=0 (the plane with the 0/0 form). The disturbing sight is the green strip that can be seen between the brown and blue planes (below the little yellow triangle that indicates Maple's fsolve candidate solution), the worrying thing is that it looks like this strip of green indicates a continuum of solution.

Is the visible continuum of solutions an optical illusion or a consequence of the cd/xd term appearing in the qd equation? Indeed, recall that:

qd = (q-b1)*(b2-q) + (q-b)*cd/xd

The picture can be turned around:

Now, let me remind you of what qd=0  (the green plane above) looks like from a different angle:

So, is there a continuum or just the illusion of a continuum? or maybe no solution at all? I'm so confused at this stage, I don't know what to think anymore.

The worksheet for these calculations is here

## The limit and a candidate solution (Part...

> restart: with(plots): plotsetup(default): with(plottools): with(DEtools):

> # Question: Does the 3-D system xd=0,cd=0,qd=0 have a solution in (x,c,q) such that x>0, c>0 and b1<q<b2 ? Is it unique?

> b1:=2/100: b2:=3/100: b:=25/1000: d := 1/10:
f := x-> x^(1/2)-d*x:
xd := f(x)-c;
cd := (diff(f(x),x)-q)*c;
qd := (q-b1)*(b2-q)+(q-b)*(diff(f(x),x)-q)*c/(f(x)-c);

> ss:= fsolve(xd,cd/c,qd, x,c,q, avoid=x=0, c=0, q=0);
xs:=subs(ss,x): cs:=subs(ss,c): qs:=subs(ss,q):

ss := {c = 2.399289639, q = 0.02511089584, x = 15.97164840}

> # Is this solution returned by fsolve valid?

> # Robert Israel's MaplePrimes Post of 8 Dec 2008. Objective: to compute the limit of the cd/xd term that appears in qd, at the point (if any) where xd=cd=qd=0

> solve(f(t^(2))=c,t);

> #Select the first or second solution

> t0 := %[2];

> q0 := D(f)(t0^2) assuming c<5/2;

> simplify ((D(f)(t^2)-q0)/(f(t^2)-c)) assuming t>0 and c<5/2 and 5-(25-10*c)^(1/2)>0;

> series(%,t=t0,2);

> convert(%,polynom);

> # With the first solution

> eq1 := q = 1/(2*((5+(25-10*c)^(1/2))^2)^(1/2))-1/10;

> eq2 := (q-b1)*(b2-q)/(q-b) = c*5/2/(5+(25-10*c)^(1/2))^2/(25-10*c)^(1/2);

> # With the second solution

> eq1 := q = 1/(2*((5-(25-10*c)^(1/2))^2)^(1/2))-1/10;

> p1 := plot(rhs(eq1), c=0..5/2, colour=blue):

> eq2 := (q-b1)*(b2-q)/(q-b) = c*(5/2/(5-(25-10*c)^(1/2))/(-5+(25-10*c)^(1/2))/(25-10*c)^(1/2));

> p2 := implicitplot(eq2, q=b1..b2, c=0..5/2, colour=red):

> display(p2):

> # The command implicitplot failed here, not sure why, but there's a workaround since eq2 is quadratic in q:

> solve(eq2,q):

> eq2b := %:

> p2b := plot(eq2b, c=0..5/2, colour=red):

> display(p1,p2b, view=[2.3..5/2,b1..b2]);

Sorry, couldn't work out how to insert a picture. Please copy-paste the following into your browser:

www.mapleprimes.com/files/9249_2D-System-1.gif

> # The red-line looks very close to horizontal near the value q=b ... but that's not the case for small values of c

> display(p2b, view=[0..0.1,0..250]):

> # I then use the computed limit to solve the 3-D system

> sss := fsolve(eq1,eq2, c=0..5/2, q=b1..b2);
css:=subs(sss,c): qss:=subs(sss,q): xss:=fsolve(f(x)=css,x);

sss := {c = 2.400425833, q = 0.02493342813}

xss := 16.01705604

> # I then compare the solution thus obtained with the one obtained directly by fsolve

> 100*xs/xss; 100*cs/css; 100*qs/qss;

99.71650446

99.95266698

100.7117662

> # The difference is under one percent.

> 100*qss/b;

99.73371252

> # But then, that's also true of the "candidate" solution q=b

> # Why is q=b a "candidate" solution? The graph is very disturbing because it looks like any value of q between b1 and b2 could be a solution.

>

View file details

## The limit and a candidate solution...

Once again Robert, many, many thanks. I don't usually post until I have run out of ideas, and your help is invaluable.

I am still confused. Part of me is inclined to think that there is a solution to the 3-D system; another part of me that there is no solution; yet another part that there is some sort of continuum of solutions. If I had to bet, I would pick the first option, but that's wishful thinking. Let me explain.

First, your calculation of the limit. I think that's what I had in mind. If I understand correctly, you do something like this: (i) you fix c, (ii) you select q=q0 such that both xd=0 and cd=0, so that q0=Q(c), (iii) then you compute the Taylor expansion of the fraction (f'(x)-Q(c))/(f(x)-c) with respect to x=x0, where x0 is the value that puts x on the xd=0 plane, (iv) the first term is the limit of the fraction; obtained by moving along the xd=0 and cd=0. Better said: with respect to the square-root of x, denoted t.

Now, if I compute the indeterminate limit in the way you suggest, I end up with a 3-D system with 3 unknowns (x,c,q). I can then ask Maple to solve it. If a solution is found, a solution exists, if not it doesn't, right?

So I went ahead with this logic in mind. At one point in your calculations, multiple solutions arise; you selected the first one; this first solution does not yield a solution to the 3-D system; however, the second solution that presents itself does. It's not too far from the solution of the original 3-D system found by fsolve.

In my next post, I will show the steps (I have so far failed to mix normal text, html, and maple code, so I won't try again)

## How to delete a post on mapleprimes...

Thanks a lot Doug, I'll flag the unwanted postings. I'll keep yours as it may serve as a useful tip in the future. Cheers, Patrick.

## my maple worksheet with comments...

something went wrong with my attempt to post. I can edit this post but do not know how to delete it, sorry.

## Does Maple's fsolve return a valid solut...

The little yellow triangle (that's hard to see on the pictures I uploaded above) indicates the position of the solution that fsolve returns. It looks possible to me... is it a serious candidate solution or some fluke?

## And now the html version of the same wor...

> restart: with(plots): plotsetup(default): with(plottools): with(DEtools):

> Digits:=100: interface(displayprecision=3):

> #Question: Does the 3-D system {xd=0,cd=0,qd=0} have a solution in (x,c,q) such that x>0, c>0 and 0.02<0.03 ?

> f := x-> x^.5-.1*x:
xd := f(x)-c;
cd := (diff(f(x),x)-q)*c;
qd := -(q-.2e-1)*(q-.3e-1)+(q-.25e-1)*(diff(f(x),x)-q)*c/(f(x)-c);

> Digits:=10: interface(displayprecision=10):
ss:= fsolve({xd,cd/c,qd}, {x,c,q}, avoid={x=0, c=0, q=0});
xs:=subs(ss,x): cs:=subs(ss,c): qs:=subs(ss,q):

> #Question: Is the solution Maple's fsolve return valid, precise, unique ? The system has an indeterminate form 0/0, which makes it difficult to check the alleged solution directly:

> eval(xd,{x=xs,c=cs,q=qs}); eval(cd,{x=xs,c=cs,q=qs}); eval(qd,{x=xs,c=cs,q=qs});

> #Let's visualize the 3-D system:

px:=implicitplot3d({xd}, x=0..25,c=0..5,q=0.02..0.05, numpoints=10000, plotopts, colour=brown): pc:=implicitplot3d({cd/c}, x=0..25,c=0..5,q=0.02..0.05, numpoints=10000, plotopts, colour=green):
pq:=implicitplot3d({qd}, x=0..25,c=0..5,q=0.02..0.05, numpoints=10000, plotopts, colour=blue):
display3d({pq}, orientation=[-25,70]);

> #The plot is ragged, probably because of the 0/0 form. If I multiply by xd by (f(x)-c) and plot, I eliminate some of the raggedness and what looks like a vertical plane caused by the 0/0 form.

> pq2:=implicitplot3d({qd*((x^.5-.1*x-c))}, x=0..25,c=0..5,q=0.02..0.05, numpoints=10000, plotopts, colour=black):
display3d({pq,pq2}, orientation=[155,35]);

> #It looks like ploting xd*(f(x)-c) goes some way towards making the plot cleaner without eliminating (too much or any) information.

> #I plot the 3-D system together with the alleged solution.

> points:=[[xs,cs,qs]]:

> ps:=pointplot3d(points,colour=yellow,symbol=circle,symbolsize=50):

> display3d({px,pc,pq,ps}, orientation=[-60,-100]);
display3d({px,pc,pq2,ps}, orientation=[-60,-100]);

> #Based on the above 3-D plots, Maple's fsolve solution looks like what I'm looking for. But I can't be sure, because, frankly, the plot is not all that clean. Choosing a different orientation introduces some doubt into my mind.

> display3d({px,pc,pq2,ps}, orientation=[170,95]);

>

This post was generated using the MaplePrimes File Manager