dharr

Dr. David Harrington

3020 Reputation

17 Badges

17 years, 303 days
University of Victoria
Professor or university staff
Victoria, British Columbia, Canada

Social Networks and Content at Maplesoft.com

Maple Application Center
I am a professor of chemistry at the University of Victoria, BC, Canada, where my research areas are electrochemistry and surface science. I have been a user of Maple since about 1990.

MaplePrimes Activity


These are answers submitted by dharr

If you use plot with style=point, then you can use the lists directly without pairing.

plot(P,Q,style=point,symbol=solidcircle,symbolsize=15,color=red,view=[-1..1,0..30]);

 

In Maple 2015, the following works. (Not displaying correctly in Mapleprimes - see worksheet).
 

restart;

with(InertForm):

a:=1: b:= 2: c:=1: d:= 6: e:= 2:

P:=`%*`(a,b,c,d,e);

`%*`(1, 2, 1, 6, 2)

Display(P);

`%*`(1, 2, 1, 6, 2)

Display(P,inert=false);

1.2.1.6.2

``

 

Download Inert.mw


 

eq:=2*cosh(xi)^8*alpha*B[2]^2 + 4*A[2]*sinh(xi)^6*cosh(xi)^2 + 4*A[1]*sinh(xi)^5*cosh(xi)^3 + 4*sinh(xi)^4*cosh(xi)^4*A[0];

2*cosh(xi)^8*alpha*B[2]^2+4*A[2]*sinh(xi)^6*cosh(xi)^2+4*A[1]*sinh(xi)^5*cosh(xi)^3+4*sinh(xi)^4*cosh(xi)^4*A[0]

coeff(coeff(eq,cosh(xi)^2),sinh(xi)^6);

4*A[2]

 

Download coeff.mw

This may not be the most efficient way to do it; a lot depends on if you have simple ways to generate the arguments or the functions or function names.
 

restart;

f1:=x->x^2;f2:=x->x^3;f3:=x->sin(x);

proc (x) options operator, arrow; x^2 end proc

proc (x) options operator, arrow; x^3 end proc

proc (x) options operator, arrow; sin(x) end proc

nrows:=4:
A:=Matrix(nrows,3):
for i to nrows do
  A[i,...]:=<f1(i/2.)|f2(i*2.)|f3(i*1.5)>;
end do:

A;

Matrix([[.2500000000, 8., .9974949866], [1.000000000, 64., .1411200081], [2.250000000, 216., -.9775301177], [4.000000000, 512., -.2794154982]])

ExcelTools:-Export(A,cat(currentdir(),"/test.xlsx"));

 

 

Download Excel.mw

For a first order ode, you can specify only one initial/boundary condition. Then your first case returns a solution in terms of Heaviside.

In my version of Maple, I don't get a solution for the second one if I just take the phi(0) condition.

restart;

eq:=diff(theta(t),t)=-(A[1]*Dirac(t-T[1])+A[2]*Dirac(t-T[2]));

diff(theta(t), t) = -A[1]*Dirac(t-T[1])-A[2]*Dirac(t-T[2])

bcon := theta(0)=-v[1];

theta(0) = -v[1]

ans:=dsolve({eq,bcon},theta(t));

theta(t) = -A[1]*Heaviside(t-T[1])-A[2]*Heaviside(t-T[2])-v[1]+A[1]*Heaviside(-T[1])+A[2]*Heaviside(-T[2])

 


 

Download Dirac.mw

Not sure if you already have this in a string or that is the file contents. Either way, this should work.

restart;

str:="[[1. 0.\n]\n[0. 1.]]"

"[[1. 0.
]
[0. 1.]]"

Scan to give a Matrix of size 2,2, ignoring characters "[", "]" and ".".

You could read this in directly from a file using fscanf instead of sscanf.

See ?rtable_scanf for more information

A:=sscanf(str,"%{2,2;d([].)}am")[];

A := Matrix(2, 2, {(1, 1) = 1, (1, 2) = 0, (2, 1) = 0, (2, 2) = 1})

 

Download MatrixInput.mw

I didn't submit this earlier since I had the same problem that @Rouben Rostamian had. But if you assume that it doesn't start at rest, but is given an initial finite y-component of the velocity, then you can derive the semicubical parabola.
 

restart

Solve isochrone problem (solution is Neile's semicubical parabola).

We start a mass m at positive coordinates (x__0, y__0) with x component of velocity zero and y component of velocity -v__y and let it move toward the origin under the action of gravity along some curve for which the y (height) coordinate changes by equal amounts in equal times. What is the curve?

 

Solution. Run the problem backwards. At time zero, let the mass leave the origin with velocity components v__x(0) = v__0 and v__y. The problem statement means that v__y is constant.

Let t be the time the particle arrives at coordinates (x,y), so we know y = v__y*t

yt := v__y*t;

v__y*t

At time zero, there is only kinetic energy

E__0 := (1/2)*m*(v__0^2+v__y^2);

(1/2)*m*(v__0^2+v__y^2)

At time t, the energy must be the same

eq := E__0 = (1/2)*m*(v__x^2+v__y^2)+m*g*v__y*t

(1/2)*m*(v__0^2+v__y^2) = (1/2)*m*(v__x^2+v__y^2)+m*g*v__y*t

So the x component of the velocity as a function of time is

vt__x := solve(eq, v__x)[1];

(-2*g*t*v__y+v__0^2)^(1/2)

And this will be zero at time t__0, after which the problem doesn't make sense

t__0 := solve(vt__x, t);

(1/2)*v__0^2/(g*v__y)

So the x coordinate is

xt := `assuming`([int(vt__x, t = 0 .. t)], [v__0::positive]);

-(1/3)*((-2*g*t*v__y+v__0^2)^(3/2)-v__0^3)/(g*v__y)

x__0 := eval(xt, t = t__0);

(1/3)*v__0^3/(g*v__y)

And the y coordinate at the end is

y__0 := v__y*t__0;

(1/2)*v__0^2/g

Eliminate t to get y(x)

ans := solve({x = xt, y = yt}, {t, y}, explicit)[1];

{t = -(1/2)*((-3*g*v__y*x+v__0^3)^(2/3)-v__0^2)/(g*v__y), y = -(1/2)*((-3*g*v__y*x+v__0^3)^(2/3)-v__0^2)/g}

This is a shifted form of the semicubical parabola.

simplify(eval((y-y__0)^3/(x-x__0)^2, ans));

-(9/8)*v__y^2/g

NULL

 

Download Niele.mw

Edit: The conclusion is that to get the conventional form, take (x0,y0) at the origin, which means to take v0 =0, as deduced by Rouben. A more straightforward way rather than going backwards is here:

Niele2.mw

cos(t) goes to zero, meaning 2/cos(t) goes to infinity. So restrict the range with coordinateview. The answer is a vertical line, as you said.
 

plots:-polarplot(2/cos(t), t = 0 .. 2*Pi,coordinateview=[0..6,0..2*Pi],color=red)

``

 

Download polarplot.mw

Not clear exactly how efficient it needs to be. The straightforward way is just to tell sort how to order them.

Download norm.mw

Edit: A one pass using selectremove is more efficient since you are not sorting them all

Download selectremove.mw

PDE:=diff(u(x,t),t)+diff(conjugate(u(x,t)),x)

is correct Maple syntax.

(ans on Mapleprimes doesn't show all the content in the worksheet.)

restart;

f := 2*y^3*z - y^2 -2*m;

2*y^3*z-y^2-2*m

ans:=solve(f,y,parametric,real);

ans := piecewise(And(z = 0, m = 0), [[y = 0]], And(z = 0, m < 0), [[y = sqrt(-2*m)], [y = -sqrt(-2*m)]], And(0 < z, m = 0), [[y = 0], [y = 1/(2*z)]], And(0 < z, m = -1/(54*z^2)), [[y = -1/(6*z)], [y = 1/(3*z)]], And(0 < z, 0 < m), [[y = RootOf(2*_Z^3*z-_Z^2-2*m, index = real[1])]], And(0 < z, m < -1/(54*z^2)), [[y = RootOf(2*_Z^3*z-_Z^2-2*m, index = real[1])]], And(z < 0, m = 0), [[y = 0], [y = 1/(2*z)]], And(z < 0, m = -1/(54*z^2)), [[y = -1/(6*z)], [y = 1/(3*z)]], And(z < 0, 0 < m), [[y = RootOf(2*_Z^3*z-_Z^2-2*m, index = real[1])]], And(z < 0, m < -1/(54*z^2)), [[y = RootOf(2*_Z^3*z-_Z^2-2*m, index = real[1])]], And(0 < z, -1/(54*z^2) < m, m < 0), [[y = RootOf(2*_Z^3*z-_Z^2-2*m, index = real[1])], [y = RootOf(2*_Z^3*z-_Z^2-2*m, index = real[2])], [y = RootOf(2*_Z^3*z-_Z^2-2*m, index = real[3])]], And(z < 0, -1/(54*z^2) < m, m < 0), [[y = RootOf(2*_Z^3*z-_Z^2-2*m, index = real[1])], [y = RootOf(2*_Z^3*z-_Z^2-2*m, index = real[2])], [y = RootOf(2*_Z^3*z-_Z^2-2*m, index = real[3])]], [])

So there are two cases that give three real roots. Let's try the case with And(z < 0, -1/(54*z^2) < m, m < 0)

And(z < 0, -(1/54)/z^2 < m, m < 0)

mz:={m=-1/2,z=-1/20};

{m = -1/2, z = -1/20}

is(eval((And(z<0, -1/(54*z^2) < m, m < 0)),mz));

true

ans1:=eval(ans,mz);
evalf(ans1);

[[y = RootOf(_Z^3+10*_Z^2-10, index = real[1])], [y = RootOf(_Z^3+10*_Z^2-10, index = real[2])], [y = RootOf(_Z^3+10*_Z^2-10, index = real[3])]]

[[y = -9.897926849], [y = -1.057474507], [y = .9554013566]]

Symbolic

map(x->convert(value(rhs(x[])),radical),ans1);

[-(1/6)*(-865+(15*I)*1119^(1/2))^(1/3)-(50/3)/(-865+(15*I)*1119^(1/2))^(1/3)-10/3+((1/2)*I)*3^(1/2)*((1/3)*(-865+(15*I)*1119^(1/2))^(1/3)-(100/3)/(-865+(15*I)*1119^(1/2))^(1/3)), -(1/6)*(-865+(15*I)*1119^(1/2))^(1/3)-(50/3)/(-865+(15*I)*1119^(1/2))^(1/3)-10/3-((1/2)*I)*3^(1/2)*((1/3)*(-865+(15*I)*1119^(1/2))^(1/3)-(100/3)/(-865+(15*I)*1119^(1/2))^(1/3)), (1/3)*(-865+(15*I)*1119^(1/2))^(1/3)+(100/3)/(-865+(15*I)*1119^(1/2))^(1/3)-10/3]

NULL

 

Download CubicSolve.mw

The plot structure is hard to modify, but if you get the graph from SubgroupLattice, then you can use RelabelVertices from the GraphTheory package. The problem is that the vertex labels need to be unique. One way is to add the GroupOrder in parentheses after the original label.

 graphs.mw

Here is another way. [Edit: LinearAlgebra is for RandomMatrix and Copy only, but copy (not in LinearAlgebra) will also work.]

restart

with(LinearAlgebra):

A := RandomMatrix(3, 3); B := RandomMatrix(3, 3)

A := Matrix(3, 3, {(1, 1) = 27, (1, 2) = 99, (1, 3) = 92, (2, 1) = 8, (2, 2) = 29, (2, 3) = -31, (3, 1) = 57, (3, 2) = -76, (3, 3) = -32})

B := Matrix(3, 3, {(1, 1) = 57, (1, 2) = -76, (1, 3) = -32, (2, 1) = 27, (2, 2) = -72, (2, 3) = -74, (3, 1) = -93, (3, 2) = -2, (3, 3) = -4})

Replace last row of A with first row of B (A is changed)

A[3, () .. ()] := B[1, () .. ()]; A

Matrix([[27, 99, 92], [8, 29, -31], [57, -76, -32]])

Replace the last row of A with second row of B with result in C

C := Copy(A); C[3, () .. ()] := B[2, () .. ()]; C

Matrix([[27, 99, 92], [8, 29, -31], [27, -72, -74]])

``

 

Download Rowsubs.mw

plot(F(R), R = 0 .. 100) is correct and works for me. You seemed to have some extra character in that line.

How_to_plot_IR.mw

The discriminant is negative for arbitrary values of the a[i,j], so the roots are usually complex (edit: for real a[i,j]). But for values that make the disciminant zero, you will get two equal real roots. (@tomleslie left off a term of your expression).
 

z := k^2*a[1, 1]^2+k^2*b[1, 1]^2+4*k*a[0, 2]*a[1, 1]+4*k*b[0, 2]*b[1, 1]+4*a[0, 2]^2+4*b[0, 2]^2

k^2*a[1, 1]^2+k^2*b[1, 1]^2+4*k*a[0, 2]*a[1, 1]+4*k*b[0, 2]*b[1, 1]+4*a[0, 2]^2+4*b[0, 2]^2

ans := [solve(z, k)];

[2*(I*a[0, 2]*b[1, 1]-I*a[1, 1]*b[0, 2]-a[0, 2]*a[1, 1]-b[0, 2]*b[1, 1])/(a[1, 1]^2+b[1, 1]^2), -2*(I*a[0, 2]*b[1, 1]-I*a[1, 1]*b[0, 2]+a[0, 2]*a[1, 1]+b[0, 2]*b[1, 1])/(a[1, 1]^2+b[1, 1]^2)]

The discriminant "b^2-4*a*c" can be factored, so the square root disappears. But it is negative, so the roots are complex, unless the discriminant becomes zero

discrim(z, k);

-16*a[0, 2]^2*b[1, 1]^2+32*a[0, 2]*a[1, 1]*b[0, 2]*b[1, 1]-16*a[1, 1]^2*b[0, 2]^2

-16*(a[0, 2]*b[1, 1]-a[1, 1]*b[0, 2])^2

realvals := {a[0, 2] = 1, a[1, 1] = 1, b[0, 2] = 2, b[1, 1] = 2}

{a[0, 2] = 1, a[1, 1] = 1, b[0, 2] = 2, b[1, 1] = 2}

eval(discrim(z, k), realvals);

0

eval(ans, realvals);

[-2, -2]

``


 

Download discrim.mw

1 2 3 4 5 6 7 Last Page 1 of 39