Dr. John May

## 2511 Reputation

15 years, 149 days
Maplesoft
Guru

## Social Networks and Content at Maplesoft.com

I am a Senior Developer in the Mathematical Software Group and have been with Maplesoft since 2007. I am also an Adjunct Assistant Professor in the School of Computer Science at the University of Waterloo.

I have a Ph.D in Mathematics from North Carolina State University as well as Masters and Bachelors degrees from the University of Oregon. I have been working on research in computational mathematics since 1997.

My main research interests in are computational linear and polynomial algebra, especially numerical polynomial algebra. I currently work on the exact algebraic solvers as well as other subsystems of Maple.

## MaplePrimes Activity

### These are answers submitted by John May

Lists in Maple are immutable datastructures.  Any changes you make to a list are illusory and involve creating a new list.  Try it!

```A := [1, 2, 3];
B := A;
A[1] := 4;
print(B);
```

vs.

```A := Array([1, 2, 3]);
B := A;
A[1] := 4;
B;
```

The Data Structure sections of the Maple Programming Guide are essential reading on this topic.  Understanding what is going on here is pretty crucial to being able to write high performance in Maple.

## Autoexpand Sections Option...

You can gointo the View menu and turn off the autoexpansion of sections on execution.  By default, it is turned on (shown).

## You might find this handy:   r :=...

You might find this handy:

`r := dismantle['string'](x+sin(x));`

## You can efficiently generate the set of ...

You can efficiently generate the set of equations using the seq command.

```n := 100;
h := 1/(n+1);

# boundary conditions
bdc := {seq([u[0, j] = 0, u[n+1, j] = 0, u[j, 0] = 0, u[j, n+1] = 0][], j = 0 .. n+1)};

# interior points:
eqs := eval(
{seq(seq(
(u[i+1, j]-u[i, j])*(u[i+1, j]-2*u[i, j]+u[i-1, j])+h*(u[i,j+1]-2*u[i, j]+u[i, j-1]) = 0,
j = 1 .. n), i = 1 .. n)},
bdc);
```

However, this is a very large system of quadratic equations and probably isn't condusive to having a general purpose solver applied to it.  fsolve, for example, can solve this for n=6, but gives up after 90 seconds without finding an answer when n=7 (49 equations).  At n=100 you have 10000 equations and probably no hope with fsolve returning at all.

## It looks like you might be able to solve...

It looks like you might be able to solve your equations symbolically in terms of C__A

```a := tau = (C__A0-C__A)/(-r__A);
b := tau = (C__B0-C__B)/(-r__B);
c := tau = (C__C0-C__C)/(-r__C);
d := tau = (C__D0-C__D)/(-r__D);
sol := solve([a, b, c, d], {C__B, C__C, C__D, tau})```

I can't tell for sure since you don't provide the definitions for all the other symbols not solved for.
Even if you can't solve symbolically you might be able to wrap the call to fsolve in a procedure and plot that procedure.

## Refinement using isolve...

@gkokovidis's soltuion does most of the job, but assume, Digits, and assign are not needed here, and the choice of i=0..3 does not need to be made appriori, if you use isolve:

```eqns := {10*cos((6*(1/10))*t)-10*cos((3/10)*t+(1/4)*Pi),
10*sin((6*(1/10))*t)-10*sin((3/10)*t+(1/4)*Pi)};
tsols := solve(eqns, allsolutions);
zsols := isolve({rhs(tsols[]) > 0, rhs(tsols[]) < 70});

sols := {seq(eval(tsols, s), s in zsols)};
```

sols in this case:

## Statistics:-SurfacePlot...

The command you want is probably plots:-surfdata you can then add a rotation annimation with the plot3d/viewpoint option.

If your data is in a Matrix M, something like this should work:

```plots:-surfdata(M, viewpoint="circleleft");
```

## Assignment and evaluation is a very frag...

Assignment and evaluation is a very fragile way to specialize the variables in expressions - and it is especially bad when paired with the wierd evaluation rules for Arrays/Matrix.  I prefer to use eval, which works well here:

```M__11 := Matrix(eval(Mm__11,
{I__0 = 0.102033222998e-2,
I__1 = .80922342232909,
I__2 = 1,
a = .2223212330091,
b = .293939202032292
}), datatype = float[8]);```

## Wrap Solve...

A nice limited way to avoid the message is to wrap solve in a proc.

```wrapsolve := proc () return solve(args) end proc;
solve(sin(x*y) < 0);
#Warning, solutions may have been lost
wrapsolve(sin(x*y) < 0);
#No warning
```

(solve checks if it is being called directly from the top level before issuing the warning - this is an easy way to trick it into not issuing the message)

## There are a lot of problems with the cod...

There are a lot of problems with the code you have posted.  This is my attempt to fix it.

```psi_0 := proc(n,x) # not used
(1/sqrt(sqrt(pi)*2^n*factorial(n))*exp(-x^2/2)*HermiteH(n,x));
end proc;

psi := proc(n, x)
exp(-1/2*x^2)*HermiteH(n, x)/sqrt(sqrt(Pi)*2^n*n!);
end proc;

psi2_0 := proc(a,x) # not used
(1/sqrt(sqrt(Pi)*2^a*factorial(a))*exp(-x^2/2)*HermiteH(a,x));
end proc;

psi2 := proc(a, x)
exp(-1/2*x^2)*HermiteH(a, x)/sqrt(sqrt(Pi)*2^a*a!);
end proc;

result := proc(n,a)
psi(n,x)*psi2(a,x);
end proc;

for n from 0 to 2 do
for a from 0 to 2 do
print(evalf(int(result(n,a),x=0..infinity)));
end do;
end do;```

## The GlobalOptimization toolbox is a sepe...

The GlobalOptimization toolbox is a seperately purchased addon to Maple.  It doesn't look like you have it installed.

## I generally avoid the "File > Ex...

I generally avoid the "File > Export as > LaTeX"  menu since by design it is trying to reproduce the exact formatting of a Maple worksheet in LaTeX. (and, apparently, it appears to have other issues)

If you have mathematical expressions in Maple and you want to export them you can use the latex command and then paste the result of that into your latex document.

latex({solve(a*x^2+b*x+c, x)});

If you are exporting actually code, I would first cut and paste it into a "Code Edit Region" and then paste that code into a preformated (verbatim) region in your latex document.

## Another way, using plot x vs. T...

A key first step is to correct the input: there is a x(3-2*x) that needs to be changed to x*(3-2*x). Then you can solve for T in terms of x and plot:

eq := (1/4)*x*(3-2*x)^2/(1-x)^3 = 18325.73901+exp(36.58420421-10902.09286/T);
expr := rhs(isolate(eq, T));
a := fsolve(expr = 300, x = 0 .. 1);
b := fsolve(expr = 800, x = .9 .. 1.1);
plot(expr, x = a .. b);

## List of Colors...

SpatterPlot requires a list of colors and you've provided just a single color.  The error message does a poor job formatting things, so it's a more cryptic than necessary, unfortunately.

SpatterPlot( [barvy], symbol="box");

will work

John

## SolveTools:-SemiAlgebraic...

The best tool for this sort of polynomial system is probably ?SolveTools:-SemiAlgebraic  it symbolically finds the real solutions to a set of equations and inequalities.

```sols := SolveTools:-SemiAlgebraic({eq1,eq2,r>0,x>0}); # large
evalf[25](sols);```

` `

John

 1 2 3 4 5 6 7 Page 3 of 10
﻿