Dr. Patrick T

## 2083 Reputation

13 years, 259 days

## automatic simplification...

Since this thread is still alive and well, I'm curious about the OP's first question: why doesn't the following simplify automatically?

```
diff(x^(n), x);

n
x  n
----
x

simplify(%);

(n - 1)
x        n
```

## In the above, insert this...

In the above, insert this below the line that starts with A:=, and you will see a plot of the two solutions for A.

```plot( evalf( eval( expr2-T, params ) ), B=-1..1);
```

Perhaps what you want is to assign A1 to the sol in the range -1..0 and A2 to the sol in the range 0..1 and accordingly modify your assignments of u0 and u1 such that u0:=A1 and u1:=A2. Depends on your problem at hand.

## Sometimes it's convenient to...

Sometimes it's convenient to manipulate symbolic expressions only to evaluate them right at the end. This way you can more easily see that your problem comes from the fact that the constant A is missing -- because, as pointed out by Doug the equation does not have a solution, as you wrote it.

I found that A is positive in the range B=0..1 but negative in the range -1..0, but I haven't looked into it at all beyond that.

Some of your code was a little confusing too. After tidying up a little, I got:

```restart;
C1 := u -> u*FresnelC(u)-((1/Pi)*sin((1/2)*Pi*u^2));
S1 := u -> u*FresnelS(u)+((1/Pi)*cos((1/2)*Pi*u^2))-(1/Pi);

params := [R0=0.7, R1=1.5, T=3];

expr1 := (R0+R1)*(((Pi*C1(1))^2+(Pi*S1(1)+1)^2))^(1/2);
eval( expr1, params );
evalf(%);

expr2 := (R0+R1)*(((Pi*C1(B))^2+(Pi*S1(B)+1)^2))^(1/2);
eval( expr2, params );
```
```# Make sure this equation does have a solution, e.g. by changing the value of T.
A := fsolve( evalf( eval( expr2=T, params ) ), B, 0..1);
```
```# In the following I add new constants to the list of parameters
# You don't have to do it this way, but I find it convenient to carry symbolic
# expressions until I need to evaluate them.
a0 = evalf( eval( Pi*R0*A, params ) ): params := [op(params), %]:
a1 = evalf( eval( Pi*R1*A, params ) ): params := [op(params), %]:
D0 = evalf( eval( -Pi*R0*C1(A), params ) ): params := [op(params), %]:
E0 = evalf( eval( -Pi*R0*S1(A)-R0, params ) ): params := [op(params), %]:
D1 = evalf( eval( Pi*R1*C1(A), params ) ): params := [op(params), %]:
E1 = evalf( eval( Pi*R1*S1(A)+R1, params ) ): params := [op(params), %];

with(plots): with(plottools):
c0 := eval( circle([D0, E0], R0), params ):
c1 := eval( circle([D1, E1], R1), params ):
display( [c0,c1], scaling=constrained );

#Clarify what to do with u0 and u1, in your code they are both equal to A
u0 := A:
u1 := -A:
p1 := plot( eval( [-a0*FresnelC(u), -a0*FresnelS(u)], params ), u=u0..u1 ):
p2 := plot( eval( [a1*FresnelC(u), a1*FresnelS(u)], params ), u=u0..u1 ):
display(p1,p2);
```
```#Check Doug's suggestions for more plot options
```

```restart;
ode := diff(y(x), x)-x*y(x)+exp(-y(x));

other ways to gather information about the solution:

```
with(DEtools);
dfieldplot(ode, y(x), x=-10..10, y=-5..5, scene=[x,y(x)]);
DEplot(ode, y(x), x=0..1, [[y(0)=0.1], [y(0)=0.5], [y(0)=0.8]], arrows=medium);

```

## criteria...

Let me offer a rather uninformed opinion here: what makes the specific parameter combination special in the picture you have linked to? Is it aesthetic beauty or a specific mathematical property? If the former, I doubt there's any rigorous method to get them, but if the latter then there may be hope.

## do I understand?...

Am I summarizing your problem correctly below? It looks like the maximum is at a corner of your range, both g and h at their smallest values. Correct my input if I misunderstood:

restart;
with(plots);
Br:=1.43;
L:=5E-3;
First:=Int(Int(((2*L)/(x^2+(L/2)^2+(h-z)^2)^(3/2)),z=(-L/2)..(L/2)),x=(g/2)..((g/2)+L));
Second:=Int(Int(((2*h+L)/(x^2+y^2+(h+(L/2))^2)^(3/2)),y=(-L/2)..(L/2)),x=(g/2)..((g/2)+L));
Third:=Int(Int(((-2*h+L)/(x^2+y^2+(h-(L/2))^2)^(3/2)),y=(-L/2)..(L/2)),x=(g/2)..((g/2)+L));
HorBx:=(Br/(4*Pi))*(First+Second+Third);
p := plot3d(HorBx, g=0.1e-3..5e-3, h=0..10e-3):
display(p, axes=box, style=patchcontour);

## explicit plot...

The reference to x,y,z, as opposed to x,y,f(x,y) made me think the OP wanted implicitplot, but my assumption seems to have been wrong.

To change the range, color and styles, you can combine Christopher's and Doug's suggestions and do something like:

p1 := plot3d( 0.5, x=-2*Pi..2*Pi, y=-2*Pi..2*Pi, color=black):
p2 := plot3d( sin(x*y), x=-Pi..Pi, y=-Pi..Pi, style=patchcontour):
display([p1,p2],axes=box, orientation=[50,50]);

## e.g. use implicitplot3d...

```restart;
with(plots);
p1 := z-0.5;
implicitplot3d({p1},x=-1..1,y=-1..1,z=-1..1, axes=box);

```

## same using phaseportrait...

with the phaseportrait command you can add particular trajectories to the phase plane:

```
phaseportrait([diff(x(t),t) = -2*x(t)-y(t),
diff(y(t),t) = y(t)-x(t)^2+x(t)],
[x(t),y(t)],
t=0..5,
[[x(0)=8,y(0)=8],[x(0)=-5,y(0)=10]],
x=-10..10, y=-10..10,
scene=[x(t),y(t)],
color=black);```
```
```

## you can also add labels afterwards...

if your purpose is to produce plots for publication or presentations, one approach is to add labels and titles afterwards. I use LaTeX and pstricks to do that, once you've got a routine going it's very fast.

## dfieldplot is one way, phaseportrait ano...

```with(DEtools);
```
```dfieldplot([diff(x(t),t) = -2*x(t)-y(t),diff(y(t),t) = y(t)-x(t)^2+x(t)], [x(t),y(t)], t=-5..5, x=-10..10, y=-10..10, scene=[x(t),y(t)]);
```

## student assignment...

This should probably have been posted in the Student forum (scroll down on the main page).

What have you done so far? Post your code so far, it'll help you get feedback. It is usually thought not very ethical to do students' homework from scratch.

## singularity at x=-2...

I don't understand your problem, but just copying your equations and running it with Maple 13/Classic, I get a warning about a singularity, see below.

Note that the exponential function is exp() and the irrational number pi is Pi with a capital letter.

I think you'd benefit from explaining your task in more detail, assuming we have no idea what you're after.

```
```
```ode := diff(t(x), x, x)+2*x*(diff(t(x), x))
+140*(1/10)+((1-t(x))*(1/14)+.5-.5*erf(x))*((1-t(x))*(1/14)
+.5+.5*erf(x))*exp(-1/t(x)) = 0;
bc := t(-2) = 1, (D(t))(-2) = 2*(7.5-1)*exp(-(-2)^2)/sqrt(Pi);
sol := dsolve({bc, ode}, numeric, t(x));
```
```plots[odeplot](sol, [[x, t(x)]], -2 .. 1);
Warning, cannot evaluate the solution further right of -1.6827785, probably a singularity

```

## parfrac, x...

```> convert(1/(z^2*(z^2+9)), parfrac);

1         1
- ---------- + ----
2           2
9 (z  + 9)   9 z

```

Now beware that if your expression contains several symbols you need to tell Maple over which one to perform the partial fraction representation:

```> convert(1/(z^2*(z^2+a)), parfrac,z);

1         1
- ---------- + ----
2           2
a (z  + a)   a z

> convert(1/(z^2*(z^2+a)), parfrac);
Error, (in convert/parfrac) the variable name (for conversion to partial fractions) must be provided

```

Maybe not so useful, but you should also know about fullparfrac:

```> convert(1/(z^2*(z^2+9)), fullparfrac);

/   -----                      \
|    \       /     _alpha     \|    1
|     )      |----------------|| + ----
|    /       \162 (z - _alpha)/|      2
|   -----                      |   9 z
\_alpha = %1                   /

2
%1 := RootOf(_Z  + 9)

```

## 1D input...

and you want to enter the equations as follows:

eq1 := 3*x+5*y+7*x=25;

﻿