Advanced Search

dharr

Dr. David Harrington

7493 Reputation

21 Badges

20 years, 206 days
University of Victoria
Professor or university staff
Victoria, British Columbia, Canada

Social Networks and Content at Maplesoft.com

Maple Application Center
Contact Dr. David Harrington
I am a retired professor of chemistry at the University of Victoria, BC, Canada. 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

Hyperbolas...

dharr 7493
December 19 2024
1 2

No artificial intelligence used.

We consider an ellipse x^2/a^2+y^2/b^2-1=0 and 2 vertices of this ellipse A(a,0) and B(0,b). We imagine a variable equilateral hyperbola passing through the points O, A and B. This curve meets the ellipse at 2 other points A1 and B1. Show that the line A1B1 passes through a fixed point.

https://www.mapleprimes.com/questions/239553-Variable-Equilateral-Hyperbole-Passing

> 

restart
> 

with(Student:-MultivariateCalculus); with(plots); _local(D)

Define ellipse and an arbitrary hyperbola (B^2-4A*C>0) (see Wikipedia)

> 

ellipse := x^2/a^2+y^2/b^2-1; hyperbola := A*x^2+B*x*y+C*y^2+D*x+E*y+F

x^2/a^2+y^2/b^2-1

A*x^2+B*x*y+C*y^2+D*x+E*y+F

The hyperbola is equilateral when A=-C

> 

hyperbola2 := eval(hyperbola, C = -A)

A*x^2-A*y^2+B*x*y+D*x+E*y+F

Intersection points between ellipse and hyperbola are at O (assume this means the origin), A, B. Use these to eliminate D,E,F

> 

eqO := eval(hyperbola2, {x = 0, y = 0}); eqA := eval(hyperbola2, {x = a, y = 0}); eqB := eval(hyperbola2, {x = 0, y = b})

F

A*a^2+D*a+F

-A*b^2+E*b+F

> 

solve({eqA, eqB, eqO}, {D, E, F}); hyperbola3 := eval(hyperbola2, %)

{D = -A*a, E = A*b, F = 0}

-A*a*x+A*b*y+A*x^2-A*y^2+B*x*y

A1 and B1 are two points on both the ellipse and the hyperbola. The 4 solutions are points B, A, A1, B1; A1 and B1 are complicated expressions in A.B,a,b

> 

ans := sort([solve({ellipse, hyperbola3}, {x, y}, explicit)]); nops(%)

4

> 

A1 := eval([x, y], ans[3]); B1 := eval([x, y], ans[4])

If the assertion is true, then the line A1-B1 intersects at the same point for different {A,B} pairs.

> 

L1 := Line(eval([A1, B1][], {A = A__1, B = B__1})); L2 := Line(eval([A1, B1][], {A = A__2, B = B__2}))

And we indeed find a value independent of A__1, A__2, B__1, B__2

> 

X__0 := GetIntersection(L1, L2)

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

Construct a diagram for two different hyperbolas

> 

params0 := {a = 1.5, b = 1.}; params := `union`(params0, {A = 1., B = 1.}); params2 := `union`(remove(has, params, B), {B = 3.})

{A = 1., B = 1., a = 1.5, b = 1.}

{A = 1., B = 3., a = 1.5, b = 1.}

Red points are O,A,B,X; Black and blue points are A1,B1 for the blue and black hyperbolas

> 

ellipseplot := implicitplot(eval(ellipse, params), x = -1.5 .. 1.5, y = -1 .. 1, color = magenta); hyperbola1plot := implicitplot(eval(hyperbola3, params), x = -2 .. 2, y = -1.5 .. 1.5, color = black); hyperbola2plot := implicitplot(eval(hyperbola3, params2), x = -2 .. 2, y = -1.5 .. 1.5, color = blue); pts := pointplot(eval([[0, 0], [a, 0], [0, b], X__0], params), symbol = solidcircle, color = red, symbolsize = 10); pts1 := pointplot(eval([A1, B1], params), symbol = solidcircle, color = black, symbolsize = 10); pts2 := pointplot(eval([A1, B1], params2), symbol = solidcircle, color = blue, symbolsize = 10); pline1 := plottools:-line(eval([A1, B1][], params), color = black, linestyle = dash); pline2 := plottools:-line(eval([A1, B1][], params2), color = blue, linestyle = dash); t1 := textplot(eval([[A1[], A__1], [B1[], B__1]], params), color = black, align = below); t2 := textplot(eval([[A1[], A__1], [B1[], B__1]], params2), color = blue, align = right); t3 := textplot(eval([[0, 0, O], [X__0[], X], [a, 0, A]], params0), color = red, align = right); t4 := textplot(eval([[0, b, B]], params0), color = red, align = above); display(t1, t2, t3, t4, pline1, pline2, ellipseplot, hyperbola1plot, hyperbola2plot, pts, pts1, pts2, axes = none, scaling = constrained)

> 

NULL

Download Ellipse.mw

sort...

dharr 7493
December 17 2024
1 3

I'm not sure exactly what you mean by labelled sets. You can do something like this

ans := sort([solve({x^2 - x + 1}, x)])

sort will sort in a consistent way from session to session, though what exact the order is is nontrivial.

missing multiplication...

dharr 7493
December 17 2024
1 0

You are missing the multiplication between the (x+7) and (x-1), so l45 := (x + 7)*(x - 1) = (1 + x)^2;

You can use a space instead of * if you are using 2D input

 

Octagon...

dharr 7493
December 17 2024
4 0
> 

restart
> 

_local(O); with(Student:-MultivariateCalculus)

Vectors of adjacent sides

> 

A := [a__1, a__2]; B := [b__1, b__2]

Origin and far point C

> 

O := [0, 0]; C := A+B

Midpoints

> 

OA := (O+A)*(1/2); OB := (O+B)*(1/2); AC := (A+C)*(1/2); BC := (B+C)*(1/2)

Area of polygon formula. Accepts a sequence of points in cyclic order.

> 

Area := proc () local n, V, i; n := nargs; V := Array(0 .. n-1, [args]); simplify((1/2)*add(-V[`mod`(i+1, n)][1]*V[i][2]+V[`mod`(i+1, n)][2]*V[i][1], i = 0 .. n-1)) end proc

Parallelogram area

> 

p_area := Area(O, B, C, A)

-a__1*b__2+a__2*b__1

Lines forming octagon

> 

lines := [Line(O, BC), Line(OB, C), Line(B, AC), Line(BC, A), Line(C, OA), Line(AC, O), Line(A, OB), Line(OA, B), Line(O, BC)]
> 

intersections := seq(GetIntersection(lines[i], lines[i+1]), i = 1 .. 8)
> 

o_area := Area(intersections)

-(1/6)*b__2*a__1+(1/6)*b__1*a__2

> 

simplify(o_area/p_area)

1/6

> 

NULL

Download Parallelogram.mw

Edit: With plot: Parallelogram2.mw

minor issues...

dharr 7493
December 16 2024
1 0

See changes in red

> 

test := simplify((int((x^2+x)*(x-a+I*b)/abs(x-a-I*b)^2, x = -5 .. 5))/(I*Pi))

(1/2)*((-I*a^2+(-I+2*b)*a+I*b^2+b)*ln(a^2+b^2-10*a+25)+(I*a^2+(I-2*b)*a-I*b^2-b)*ln(a^2+b^2+10*a+25)+(-2*a^2+(-(4*I)*b-2)*a-2*b*(I-b))*arctan((a-5)/b)+(2*a^2+((4*I)*b+2)*a+2*b*(I-b))*arctan((a+5)/b)-(20*I)*a-20*I+20*b)/Pi

Above: the "full" definite integral.  The integral is being taken with respect to x

test1 := simplify(int((x^2+x)*(x-a+I*b)/abs(x-a+I*b)^2, x))NULL

(1/2)*(a^2+(1+(2*I)*b)*a+b*(I-b))*ln(a^2-2*a*x+b^2+x^2)+arctan((-x+a)/b)*(-I*a^2+(-I+2*b)*a+I*b^2+b)+x*(I*b+a+(1/2)*x+1)

The indefinite integralNULL

> 

test2 := simplify(subs(x = 5, test1)-subs(x = -5, test1))

10+(1/2)*(a^2+(1+(2*I)*b)*a+b*(I-b))*ln(a^2+b^2-10*a+25)-(1/2)*(a^2+(1+(2*I)*b)*a+b*(I-b))*ln(a^2+b^2+10*a+25)+arctan((a-5)/b)*(-I*a^2+(-I+2*b)*a+I*b^2+b)+(I*a^2+(I-2*b)*a-I*b^2-b)*arctan((a+5)/b)+(10*I)*b+10*a

"Hard coding" the definite integral

> 

simplify(I*test*Pi-test2)

0

 

The above result is now zero.

Download Complex_valued_integrals2.mw

ListTools:-Classify...

dharr 7493
December 15 2024
1 3

Maybe you only need to classify, and not order the sublists after classification? If so you only need LL and not the sort line in the procedure. I'm not too familiar with these orderings and I'm not sure which one you actually want.

> 

restart;
> 

F := [x^2 - y - 1, y^2 - y, y^3 - y, x^2 - 4, y^3 - 1, z^2 - 4, 2*x^2, x^2*y];

[x^2-y-1, y^2-y, y^3-y, x^2-4, y^3-1, z^2-4, 2*x^2, x^2*y]

> 

monomial_classify := proc(L::list(polynom), ord)
  local LL, Lmon, xx;
  Lmon := xx -> Groebner:-LeadingMonomial(xx, ord);
  LL := map(convert,convert(ListTools:-Classify(Lmon, L),'list'),'list');
  sort(LL, (a, b) -> Groebner:-TestOrder(Lmon(a[1]), Lmon(b[1]), ord));
end proc:
> 

monomial_classify(F, grlex(z,y,x));

[[2*x^2, x^2-4, x^2-y-1], [y^2-y], [z^2-4], [x^2*y], [y^3-1, y^3-y]]

> 

 

Download sort2.mw

Explore...

dharr 7493
December 14 2024
0 0

I'm assuming you want two plots, of Cb(t) and Cm(t) vs time, and not Cb(t) vs Cm(t). You seem to have cut and pasted things, so I set up to solve the problem from live commands. You may have to adjust some things to get exactly what you want. (The plot doesn't show here in Mapleprimes, but does if you download the document.) I had to make a guess about the meaning of c0 and c; they may be the wrong way around.

restart

 

 

 

 

Absorption from stomach

de1 := diff(Cm(t), t) = -Ka*Cm(t)                            

 

Glucose concentration in blood

de2 := diff(Cb(t), t) = Ka*Cm(t)-Ku*Cb(t)

 

Initial conditions

 

ins := Cb(0) = c0, Cm(0) = c 

 

Solve

solution := dsolve({de1, de2, ins})

{Cb(t) = -(-Ka+Ku)*(Ka*c+Ka*c0-Ku*c0)*exp(-Ku*t)/(Ka-Ku)^2-Ka*exp(-Ka*t)*c/(Ka-Ku), Cm(t) = c*exp(-Ka*t)}

Ceqns := eval([Cb(t), Cm(t)], solution)

[-(-Ka+Ku)*(Ka*c+Ka*c0-Ku*c0)*exp(-Ku*t)/(Ka-Ku)^2-Ka*exp(-Ka*t)*c/(Ka-Ku), c*exp(-Ka*t)]

Explore(plot(Ceqns, t = 0 .. 5, 0 .. 30), Ka = 1 .. 10, Ku = 1 .. 10, c0 = 0 .. 20, c = 0 .. 20, initialvalues = [Ka = 2, Ku = 1, c0 = 10, c = 10])

 

NULL

NULL

``

NULL

NULL

Download Explore_plot.mw

sfloat-hfloat conversion...

dharr 7493
December 11 2024
3 0

See the ?evalhf help page, There is overhead converting from software floats to hardware floats and back again. The conversion to hfloats is presumably still present below, but it doesn't have to convert back.

> 

restart;
> 

UseHardwareFloats:=false;

false

> 

CodeTools:-Usage( for i from 1 to 100 by 0.001 do exp(i) end do):
CodeTools:-Usage( for i from 1 to 100 by 0.001 do evalhf[hfloat](exp(i)) end do):

memory used=0.87GiB, alloc change=109.00MiB, cpu time=2.31s, real time=2.42s, gc time=93.75ms
memory used=12.84MiB, alloc change=0 bytes, cpu time=78.00ms, real time=68.00ms, gc time=0ns

 

 

Download evalhf.mw

General method...

dharr 7493
December 11 2024
2 1

Since all trig and exp functions of interest are expressible in terms of exponentials, this method may be more general; at least it doesn't require the generation of different side relations for each case.

The solved parameters are sometimes complex, which has led to psi having cosh rather than cos. This may be OK or maybe solving for an appropriate subset of the variables solves this problem; which subset to use wasn't clear to me.

> 

restart
> 

with(PDEtools); _local(gamma)
> 

undeclare(prime)

`There is no more prime differentiation variable; all derivatives will be displayed as indexed functions`

> 

declare(phi(x, t)); declare(psi(x, t))

phi(x, t)*`will now be displayed as`*phi

psi(x, t)*`will now be displayed as`*psi

> 

pde := diff(phi(x, t), `$`(t, 2))+diff(phi(x, t), `$`(x, 3), t)-3*(diff((diff(phi(x, t), x))*(diff(phi(x, t), t)), x))

diff(diff(phi(x, t), t), t)+diff(diff(diff(diff(phi(x, t), t), x), x), x)-3*(diff(diff(phi(x, t), x), x))*(diff(phi(x, t), t))-3*(diff(phi(x, t), x))*(diff(diff(phi(x, t), t), x))

> 

eq2 := phi(x, t) = 3*(diff(ln(psi(x, t)), x))

phi(x, t) = 3*(diff(psi(x, t), x))/psi(x, t)

> 

eq3 := (1/3)*numer(normal(eval(pde, eq2)))

(diff(diff(diff(psi(x, t), t), t), x))*psi(x, t)^4+(diff(diff(diff(diff(diff(psi(x, t), t), x), x), x), x))*psi(x, t)^4-2*(diff(diff(psi(x, t), t), x))*(diff(psi(x, t), t))*psi(x, t)^3-13*(diff(diff(psi(x, t), t), x))*(diff(diff(diff(psi(x, t), x), x), x))*psi(x, t)^3-(diff(diff(diff(diff(psi(x, t), x), x), x), x))*(diff(psi(x, t), t))*psi(x, t)^3-(diff(psi(x, t), x))*(diff(diff(psi(x, t), t), t))*psi(x, t)^3-4*(diff(diff(diff(diff(psi(x, t), t), x), x), x))*(diff(psi(x, t), x))*psi(x, t)^3-15*(diff(diff(diff(psi(x, t), t), x), x))*(diff(diff(psi(x, t), x), x))*psi(x, t)^3+69*(diff(diff(psi(x, t), t), x))*(diff(psi(x, t), x))*(diff(diff(psi(x, t), x), x))*psi(x, t)^2+2*(diff(psi(x, t), x))*(diff(psi(x, t), t))^2*psi(x, t)^2+17*(diff(diff(diff(psi(x, t), x), x), x))*(diff(psi(x, t), t))*(diff(psi(x, t), x))*psi(x, t)^2+15*(diff(diff(psi(x, t), x), x))^2*(diff(psi(x, t), t))*psi(x, t)^2+21*(diff(diff(diff(psi(x, t), t), x), x))*(diff(psi(x, t), x))^2*psi(x, t)^2-60*(diff(diff(psi(x, t), t), x))*(diff(psi(x, t), x))^3*psi(x, t)-90*(diff(diff(psi(x, t), x), x))*(diff(psi(x, t), t))*(diff(psi(x, t), x))^2*psi(x, t)+60*(diff(psi(x, t), x))^4*(diff(psi(x, t), t))

We seek a solution of this form

> 

eq4 := psi(x, t) = gamma[1]*cos(theta[0]*(t*omega[0]+x))+gamma[2]*exp(theta[1]*(t*epsilon[0]+x))+exp(-theta[1]*(t*epsilon[0]+x))

psi(x, t) = gamma[1]*cos(theta[0]*(t*omega[0]+x))+gamma[2]*exp(theta[1]*(t*epsilon[0]+x))+exp(-theta[1]*(t*epsilon[0]+x))

After substitution into the ode, we want to find the coefficients of the trig functions, which we can set to zero since we want it to be true for all x and t. In principle solve/identity with two identity variables would work, but Maple allows only one.
Notice that substituting in leads us not only to have cos, but also sin, and if we have sin^2 and cos^2 the we can simplify using sin^2 + cos^2 =1 as a side relation. But this means tailoring the answer to the particular ode, since another may have other trig relationships.

> 

eq5 := normal(eval(eq3, eq4)); indets(eq5)

{t, x, epsilon[0], gamma[1], gamma[2], omega[0], theta[0], theta[1], cos(theta[0]*(t*omega[0]+x)), exp(theta[1]*(t*epsilon[0]+x)), exp(-theta[1]*(t*epsilon[0]+x)), sin(theta[0]*(t*omega[0]+x))}

We have two distinct arguments, call them X and Y

> 

eq6 := eval(eq5, {theta[0]*(t*omega[0]+x) = X, theta[1]*(t*`ε`[0]+x) = Y}); indets(eq6)

{X, Y, epsilon[0], gamma[1], gamma[2], omega[0], theta[0], theta[1], cos(X), exp(Y), exp(-Y), sin(X)}

Try converting all the trig functions to exp.

> 

eq7 := convert(eq6, exp); indets(eq7)

{X, Y, epsilon[0], gamma[1], gamma[2], omega[0], theta[0], theta[1], exp(Y), exp(-(2*I)*X), exp(-I*X), exp(I*X), exp((2*I)*X), exp(-3*Y), exp(-2*Y), exp(-Y), exp(3*Y), exp(Y-(2*I)*X), exp(Y+(2*I)*X), exp(2*Y-I*X), exp(2*Y+I*X), exp(3*Y-(2*I)*X), exp(3*Y+(2*I)*X), exp(4*Y-I*X), exp(4*Y+I*X)}

All the exps are products of exp(Y) and exp(IX) - convert these to eY and eX

> 

tr := solve({eX = exp(I*X), eY = exp(Y)}, {X, Y}); eq8 := simplify(eval(eq7, tr)); indets(eq8)

{X = -I*ln(eX), Y = ln(eY)}

{eX, eY, epsilon[0], gamma[1], gamma[2], omega[0], theta[0], theta[1]}

So now collecting wrt eX and eY should give only independent relationships? From which the required equations are the coefficients.

> 

eqs := {coeffs(collect(numer(normal(eq8)), {eX, eY}, 'distributed'), {eX, eY})}; nops(%); indets(eqs)

32

{epsilon[0], gamma[1], gamma[2], omega[0], theta[0], theta[1]}

Next we solve for some subset of the variables in terms of the others. Here I just solve for all of them

> 

vars := indets(eqs); ans := solve(eqs, vars)

{epsilon[0], gamma[1], gamma[2], omega[0], theta[0], theta[1]}

{epsilon[0] = 0, gamma[1] = 0, gamma[2] = gamma[2], omega[0] = omega[0], theta[0] = theta[0], theta[1] = theta[1]}, {epsilon[0] = epsilon[0], gamma[1] = 0, gamma[2] = 0, omega[0] = omega[0], theta[0] = theta[0], theta[1] = theta[1]}, {epsilon[0] = epsilon[0], gamma[1] = 0, gamma[2] = gamma[2], omega[0] = omega[0], theta[0] = theta[0], theta[1] = 0}, {epsilon[0] = 0, gamma[1] = gamma[1], gamma[2] = gamma[2], omega[0] = omega[0], theta[0] = 0, theta[1] = theta[1]}, {epsilon[0] = epsilon[0], gamma[1] = gamma[1], gamma[2] = gamma[2], omega[0] = omega[0], theta[0] = 0, theta[1] = 0}, {epsilon[0] = epsilon[0], gamma[1] = gamma[1], gamma[2] = gamma[2], omega[0] = 0, theta[0] = theta[0], theta[1] = 0}, {epsilon[0] = 0, gamma[1] = gamma[1], gamma[2] = gamma[2], omega[0] = 0, theta[0] = theta[0], theta[1] = theta[1]}, {epsilon[0] = epsilon[0], gamma[1] = gamma[1], gamma[2] = (1/4)*gamma[1]^2, omega[0] = -epsilon[0], theta[0] = -I*theta[1], theta[1] = theta[1]}, {epsilon[0] = epsilon[0], gamma[1] = gamma[1], gamma[2] = (1/4)*gamma[1]^2, omega[0] = -epsilon[0], theta[0] = I*theta[1], theta[1] = theta[1]}

> 

eqI := ans[9]

{epsilon[0] = epsilon[0], gamma[1] = gamma[1], gamma[2] = (1/4)*gamma[1]^2, omega[0] = -epsilon[0], theta[0] = I*theta[1], theta[1] = theta[1]}

> 

eqpsi := eval(eq4, eqI)

psi(x, t) = gamma[1]*cosh(theta[1]*(-t*epsilon[0]+x))+(1/4)*gamma[1]^2*exp(theta[1]*(t*epsilon[0]+x))+exp(-theta[1]*(t*epsilon[0]+x))

> 

eqphi := eval(eq2, eqpsi)

phi(x, t) = 3*(gamma[1]*theta[1]*sinh(theta[1]*(-t*epsilon[0]+x))+(1/4)*gamma[1]^2*theta[1]*exp(theta[1]*(t*epsilon[0]+x))-theta[1]*exp(-theta[1]*(t*epsilon[0]+x)))/(gamma[1]*cosh(theta[1]*(-t*epsilon[0]+x))+(1/4)*gamma[1]^2*exp(theta[1]*(t*epsilon[0]+x))+exp(-theta[1]*(t*epsilon[0]+x)))

> 

simplify(eval(pde, eqphi))

0

> 

NULL

Download PDE-Hard.mw

hand calculation...

dharr 7493
December 11 2024
1 4

Your hand calculation of ode1 does not appear to be correct. See attached.

problem.mw

CP function...

dharr 7493
December 10 2024
2 1

The problem is that your function CP can't take T1(t) as an argument since that is not what CoolProp wants. But dsolve will complain if you don't have it. So you can make the procedure just return unevaluated if given a symbolic input to keep dsolve happy, but pass a numeric value to CoolProp. See attached.

Download Coolprop.mw

Solution...

dharr 7493
December 09 2024
2 2
> 

restart
> 

with(PDEtools); _local(gamma)
> 

undeclare(prime)

`There is no more prime differentiation variable; all derivatives will be displayed as indexed functions`

> 

declare(phi(x, t)); declare(psi(x, t))

phi(x, t)*`will now be displayed as`*phi

psi(x, t)*`will now be displayed as`*psi

> 

pde := diff(phi(x, t), t)+diff(phi(x, t), x)+(3/2)*phi(x, t)*(diff(phi(x, t), x, t))-(3/2)*(diff(phi(x, t), x))*(diff(phi(x, t), t))+(1/2)*phi(x, t)^2*(diff(phi(x, t), t))

diff(phi(x, t), t)+diff(phi(x, t), x)+(3/2)*phi(x, t)*(diff(diff(phi(x, t), t), x))-(3/2)*(diff(phi(x, t), x))*(diff(phi(x, t), t))+(1/2)*phi(x, t)^2*(diff(phi(x, t), t))

> 

eq2 := phi(x, t) = 3*(diff(ln(psi(x, t)), x))

phi(x, t) = 3*(diff(psi(x, t), x))/psi(x, t)

> 

eq3 := (1/3)*numer(normal(eval(pde, eq2)))

2*(diff(diff(psi(x, t), t), x))*psi(x, t)-9*(diff(diff(psi(x, t), t), x))*(diff(diff(psi(x, t), x), x))+2*(diff(diff(psi(x, t), x), x))*psi(x, t)-2*(diff(psi(x, t), x))^2-2*(diff(psi(x, t), x))*(diff(psi(x, t), t))+9*(diff(psi(x, t), x))*(diff(diff(diff(psi(x, t), t), x), x))

I'm following the notation in the paper - there were missing multiplications in this line

> 

eq4 := psi(x, t) = gamma[1]*cos(theta[0]*(t*omega[0]+x))+gamma[2]*exp(theta[1]*(t*epsilon[0]+x))+exp(-theta[1]*(t*epsilon[0]+x))

psi(x, t) = gamma[1]*cos(theta[0]*(t*omega[0]+x))+gamma[2]*exp(theta[1]*(t*epsilon[0]+x))+exp(-theta[1]*(t*epsilon[0]+x))

After the following substitution, we are supposed to find 5 equations from the coefficients, i.e., independent of x and t.

> 

eq5a := eval(eq3, eq4)

2*(-gamma[1]*theta[0]^2*omega[0]*cos(theta[0]*(t*omega[0]+x))+gamma[2]*theta[1]^2*epsilon[0]*exp(theta[1]*(t*epsilon[0]+x))+theta[1]^2*epsilon[0]*exp(-theta[1]*(t*epsilon[0]+x)))*(gamma[1]*cos(theta[0]*(t*omega[0]+x))+gamma[2]*exp(theta[1]*(t*epsilon[0]+x))+exp(-theta[1]*(t*epsilon[0]+x)))-9*(-gamma[1]*theta[0]^2*omega[0]*cos(theta[0]*(t*omega[0]+x))+gamma[2]*theta[1]^2*epsilon[0]*exp(theta[1]*(t*epsilon[0]+x))+theta[1]^2*epsilon[0]*exp(-theta[1]*(t*epsilon[0]+x)))*(-gamma[1]*theta[0]^2*cos(theta[0]*(t*omega[0]+x))+gamma[2]*theta[1]^2*exp(theta[1]*(t*epsilon[0]+x))+theta[1]^2*exp(-theta[1]*(t*epsilon[0]+x)))+2*(-gamma[1]*theta[0]^2*cos(theta[0]*(t*omega[0]+x))+gamma[2]*theta[1]^2*exp(theta[1]*(t*epsilon[0]+x))+theta[1]^2*exp(-theta[1]*(t*epsilon[0]+x)))*(gamma[1]*cos(theta[0]*(t*omega[0]+x))+gamma[2]*exp(theta[1]*(t*epsilon[0]+x))+exp(-theta[1]*(t*epsilon[0]+x)))-2*(-gamma[1]*theta[0]*sin(theta[0]*(t*omega[0]+x))+gamma[2]*theta[1]*exp(theta[1]*(t*epsilon[0]+x))-theta[1]*exp(-theta[1]*(t*epsilon[0]+x)))^2-2*(-gamma[1]*theta[0]*sin(theta[0]*(t*omega[0]+x))+gamma[2]*theta[1]*exp(theta[1]*(t*epsilon[0]+x))-theta[1]*exp(-theta[1]*(t*epsilon[0]+x)))*(-gamma[1]*theta[0]*omega[0]*sin(theta[0]*(t*omega[0]+x))+gamma[2]*theta[1]*epsilon[0]*exp(theta[1]*(t*epsilon[0]+x))-theta[1]*epsilon[0]*exp(-theta[1]*(t*epsilon[0]+x)))+9*(-gamma[1]*theta[0]*sin(theta[0]*(t*omega[0]+x))+gamma[2]*theta[1]*exp(theta[1]*(t*epsilon[0]+x))-theta[1]*exp(-theta[1]*(t*epsilon[0]+x)))*(gamma[1]*theta[0]^3*omega[0]*sin(theta[0]*(t*omega[0]+x))+gamma[2]*theta[1]^3*epsilon[0]*exp(theta[1]*(t*epsilon[0]+x))-theta[1]^3*epsilon[0]*exp(-theta[1]*(t*epsilon[0]+x)))

We have sin^2 and cos^1 with the same coefficients

> 

eq5b := collect(normal(eq5a), {cos, sin})

(-9*gamma[1]^2*omega[0]*theta[0]^4-2*gamma[1]^2*omega[0]*theta[0]^2-2*gamma[1]^2*theta[0]^2)*cos(theta[0]*(t*omega[0]+x))^2+(9*exp(theta[1]*(t*epsilon[0]+x))*epsilon[0]*gamma[1]*gamma[2]*theta[0]^2*theta[1]^2+9*exp(theta[1]*(t*epsilon[0]+x))*gamma[1]*gamma[2]*omega[0]*theta[0]^2*theta[1]^2+9*exp(-theta[1]*(t*epsilon[0]+x))*epsilon[0]*gamma[1]*theta[0]^2*theta[1]^2+9*exp(-theta[1]*(t*epsilon[0]+x))*gamma[1]*omega[0]*theta[0]^2*theta[1]^2+2*exp(theta[1]*(t*epsilon[0]+x))*epsilon[0]*gamma[1]*gamma[2]*theta[1]^2-2*exp(theta[1]*(t*epsilon[0]+x))*gamma[1]*gamma[2]*omega[0]*theta[0]^2-2*exp(theta[1]*(t*epsilon[0]+x))*gamma[1]*gamma[2]*theta[0]^2+2*exp(theta[1]*(t*epsilon[0]+x))*gamma[1]*gamma[2]*theta[1]^2+2*exp(-theta[1]*(t*epsilon[0]+x))*epsilon[0]*gamma[1]*theta[1]^2-2*exp(-theta[1]*(t*epsilon[0]+x))*gamma[1]*omega[0]*theta[0]^2-2*exp(-theta[1]*(t*epsilon[0]+x))*gamma[1]*theta[0]^2+2*exp(-theta[1]*(t*epsilon[0]+x))*gamma[1]*theta[1]^2)*cos(theta[0]*(t*omega[0]+x))+(-9*gamma[1]^2*omega[0]*theta[0]^4-2*gamma[1]^2*omega[0]*theta[0]^2-2*gamma[1]^2*theta[0]^2)*sin(theta[0]*(t*omega[0]+x))^2+(-9*exp(theta[1]*(t*epsilon[0]+x))*epsilon[0]*gamma[1]*gamma[2]*theta[0]*theta[1]^3+9*exp(theta[1]*(t*epsilon[0]+x))*gamma[1]*gamma[2]*omega[0]*theta[0]^3*theta[1]+9*exp(-theta[1]*(t*epsilon[0]+x))*epsilon[0]*gamma[1]*theta[0]*theta[1]^3-9*exp(-theta[1]*(t*epsilon[0]+x))*gamma[1]*omega[0]*theta[0]^3*theta[1]+2*exp(theta[1]*(t*epsilon[0]+x))*epsilon[0]*gamma[1]*gamma[2]*theta[0]*theta[1]+2*exp(theta[1]*(t*epsilon[0]+x))*gamma[1]*gamma[2]*omega[0]*theta[0]*theta[1]+4*exp(theta[1]*(t*epsilon[0]+x))*gamma[1]*gamma[2]*theta[0]*theta[1]-2*exp(-theta[1]*(t*epsilon[0]+x))*epsilon[0]*gamma[1]*theta[0]*theta[1]-2*exp(-theta[1]*(t*epsilon[0]+x))*gamma[1]*omega[0]*theta[0]*theta[1]-4*exp(-theta[1]*(t*epsilon[0]+x))*gamma[1]*theta[0]*theta[1])*sin(theta[0]*(t*omega[0]+x))-36*exp(theta[1]*(t*epsilon[0]+x))*exp(-theta[1]*(t*epsilon[0]+x))*epsilon[0]*gamma[2]*theta[1]^4+8*exp(theta[1]*(t*epsilon[0]+x))*exp(-theta[1]*(t*epsilon[0]+x))*epsilon[0]*gamma[2]*theta[1]^2+8*exp(theta[1]*(t*epsilon[0]+x))*exp(-theta[1]*(t*epsilon[0]+x))*gamma[2]*theta[1]^2

> 

eq5c := simplify(eq5b, {cos(theta[0]*(t*omega[0]+x))^2+sin(theta[0]*(t*omega[0]+x))^2 = 1})

(-36*(epsilon[0]*theta[1]^2-(2/9)*epsilon[0]-2/9)*theta[1]^2*gamma[2]*exp(theta[1]*(t*epsilon[0]+x))-9*gamma[1]*((((-omega[0]-epsilon[0])*theta[0]^2-(2/9)*epsilon[0]-2/9)*theta[1]^2+(2/9)*theta[0]^2*(omega[0]+1))*cos(theta[0]*(t*omega[0]+x))+sin(theta[0]*(t*omega[0]+x))*theta[0]*theta[1]*(omega[0]*theta[0]^2-epsilon[0]*theta[1]^2+(2/9)*omega[0]+(2/9)*epsilon[0]+4/9)))*exp(-theta[1]*(t*epsilon[0]+x))+9*gamma[1]*(((((omega[0]+epsilon[0])*theta[0]^2+(2/9)*epsilon[0]+2/9)*theta[1]^2-(2/9)*theta[0]^2*(omega[0]+1))*cos(theta[0]*(t*omega[0]+x))+sin(theta[0]*(t*omega[0]+x))*theta[0]*theta[1]*(omega[0]*theta[0]^2-epsilon[0]*theta[1]^2+(2/9)*omega[0]+(2/9)*epsilon[0]+4/9))*gamma[2]*exp(theta[1]*(t*epsilon[0]+x))-gamma[1]*theta[0]^2*(omega[0]*theta[0]^2+(2/9)*omega[0]+2/9))

> 

eq5d := collect(eq5c, {cos, sin})

(-9*gamma[1]*(((-omega[0]-epsilon[0])*theta[0]^2-(2/9)*epsilon[0]-2/9)*theta[1]^2+(2/9)*theta[0]^2*(omega[0]+1))*exp(-theta[1]*(t*epsilon[0]+x))+9*gamma[1]*(((omega[0]+epsilon[0])*theta[0]^2+(2/9)*epsilon[0]+2/9)*theta[1]^2-(2/9)*theta[0]^2*(omega[0]+1))*gamma[2]*exp(theta[1]*(t*epsilon[0]+x)))*cos(theta[0]*(t*omega[0]+x))+(-9*gamma[1]*theta[0]*theta[1]*(omega[0]*theta[0]^2-epsilon[0]*theta[1]^2+(2/9)*omega[0]+(2/9)*epsilon[0]+4/9)*exp(-theta[1]*(t*epsilon[0]+x))+9*gamma[1]*theta[0]*theta[1]*(omega[0]*theta[0]^2-epsilon[0]*theta[1]^2+(2/9)*omega[0]+(2/9)*epsilon[0]+4/9)*gamma[2]*exp(theta[1]*(t*epsilon[0]+x)))*sin(theta[0]*(t*omega[0]+x))-36*(epsilon[0]*theta[1]^2-(2/9)*epsilon[0]-2/9)*theta[1]^2*gamma[2]*exp(theta[1]*(t*epsilon[0]+x))*exp(-theta[1]*(t*epsilon[0]+x))-9*gamma[1]^2*theta[0]^2*(omega[0]*theta[0]^2+(2/9)*omega[0]+2/9)

coeffs of cos, sin and constant parts

> 

cospart, rest := selectremove(has, eq5d, cos); cospart := cospart/cos(theta[0]*(t*omega[0]+x)); sinpart, rest := selectremove(has, rest, sin); sinpart := sinpart/sin(theta[0]*(t*omega[0]+x)); rest

-9*gamma[1]*(((-omega[0]-epsilon[0])*theta[0]^2-(2/9)*epsilon[0]-2/9)*theta[1]^2+(2/9)*theta[0]^2*(omega[0]+1))*exp(-theta[1]*(t*epsilon[0]+x))+9*gamma[1]*(((omega[0]+epsilon[0])*theta[0]^2+(2/9)*epsilon[0]+2/9)*theta[1]^2-(2/9)*theta[0]^2*(omega[0]+1))*gamma[2]*exp(theta[1]*(t*epsilon[0]+x))

-9*gamma[1]*theta[0]*theta[1]*(omega[0]*theta[0]^2-epsilon[0]*theta[1]^2+(2/9)*omega[0]+(2/9)*epsilon[0]+4/9)*exp(-theta[1]*(t*epsilon[0]+x))+9*gamma[1]*theta[0]*theta[1]*(omega[0]*theta[0]^2-epsilon[0]*theta[1]^2+(2/9)*omega[0]+(2/9)*epsilon[0]+4/9)*gamma[2]*exp(theta[1]*(t*epsilon[0]+x))

-36*(epsilon[0]*theta[1]^2-(2/9)*epsilon[0]-2/9)*theta[1]^2*gamma[2]*exp(theta[1]*(t*epsilon[0]+x))*exp(-theta[1]*(t*epsilon[0]+x))-9*gamma[1]^2*theta[0]^2*(omega[0]*theta[0]^2+(2/9)*omega[0]+2/9)

The exponentials in "rest" cancel leaving one of the equations

> 

eq51 := expand(simplify(rest))

-9*gamma[1]^2*omega[0]*theta[0]^4-36*epsilon[0]*gamma[2]*theta[1]^4-2*gamma[1]^2*omega[0]*theta[0]^2+8*epsilon[0]*gamma[2]*theta[1]^2-2*gamma[1]^2*theta[0]^2+8*gamma[2]*theta[1]^2

We get the others from the coefficients of the exps in cospart and sinpart

> 

eq53 := expand(coeff(cospart, exp(-theta[1]*(t*`ε`[0]+x)))); eq55 := expand(coeff(cospart, exp(theta[1]*(t*`ε`[0]+x)))); eq54 := expand(coeff(sinpart, exp(theta[1]*(t*`ε`[0]+x)))); eq52 := expand(coeff(sinpart, exp(-theta[1]*(t*`ε`[0]+x))))

9*epsilon[0]*gamma[1]*theta[0]^2*theta[1]^2+9*gamma[1]*omega[0]*theta[0]^2*theta[1]^2+2*epsilon[0]*gamma[1]*theta[1]^2-2*gamma[1]*omega[0]*theta[0]^2-2*gamma[1]*theta[0]^2+2*gamma[1]*theta[1]^2

9*epsilon[0]*gamma[1]*gamma[2]*theta[0]^2*theta[1]^2+9*gamma[1]*gamma[2]*omega[0]*theta[0]^2*theta[1]^2+2*epsilon[0]*gamma[1]*gamma[2]*theta[1]^2-2*gamma[1]*gamma[2]*omega[0]*theta[0]^2-2*gamma[1]*gamma[2]*theta[0]^2+2*gamma[1]*gamma[2]*theta[1]^2

-9*epsilon[0]*gamma[1]*gamma[2]*theta[0]*theta[1]^3+9*gamma[1]*gamma[2]*omega[0]*theta[0]^3*theta[1]+2*epsilon[0]*gamma[1]*gamma[2]*theta[0]*theta[1]+2*gamma[1]*gamma[2]*omega[0]*theta[0]*theta[1]+4*gamma[1]*gamma[2]*theta[0]*theta[1]

9*epsilon[0]*gamma[1]*theta[0]*theta[1]^3-9*gamma[1]*omega[0]*theta[0]^3*theta[1]-2*epsilon[0]*gamma[1]*theta[0]*theta[1]-2*gamma[1]*omega[0]*theta[0]*theta[1]-4*gamma[1]*theta[0]*theta[1]

> 

ans := solve({eq51, eq52, eq53, eq54, eq55})

{epsilon[0] = epsilon[0], gamma[1] = 0, gamma[2] = 0, omega[0] = omega[0], theta[0] = theta[0], theta[1] = theta[1]}, {epsilon[0] = 2/(9*theta[1]^2-2), gamma[1] = 0, gamma[2] = gamma[2], omega[0] = omega[0], theta[0] = theta[0], theta[1] = theta[1]}, {epsilon[0] = epsilon[0], gamma[1] = 0, gamma[2] = gamma[2], omega[0] = omega[0], theta[0] = theta[0], theta[1] = 0}, {epsilon[0] = epsilon[0], gamma[1] = gamma[1], gamma[2] = gamma[2], omega[0] = omega[0], theta[0] = 0, theta[1] = 0}, {epsilon[0] = 2*(9*theta[0]^2-2)/(81*theta[0]^2*theta[1]^2+4), gamma[1] = gamma[1], gamma[2] = -(1/4)*gamma[1]^2*theta[0]^2/theta[1]^2, omega[0] = -2*(9*theta[1]^2+2)/(81*theta[0]^2*theta[1]^2+4), theta[0] = theta[0], theta[1] = theta[1]}, {epsilon[0] = 2/(9*theta[1]^2-2), gamma[1] = gamma[1], gamma[2] = gamma[2], omega[0] = 2/(9*theta[1]^2-2), theta[0] = RootOf(_Z^2+1)*theta[1], theta[1] = theta[1]}, {epsilon[0] = -1, gamma[1] = gamma[1], gamma[2] = 0, omega[0] = omega[0], theta[0] = 0, theta[1] = theta[1]}, {epsilon[0] = epsilon[0], gamma[1] = 2*RootOf((epsilon[0]+1)*_Z^2-gamma[2]*epsilon[0]+gamma[2]), gamma[2] = gamma[2], omega[0] = epsilon[0]-1, theta[0] = RootOf((9*epsilon[0]-9)*_Z^2+2*epsilon[0]+2), theta[1] = (1/3)*RootOf(_Z^2-2)}

> 

eqI := ans[5]

{epsilon[0] = 2*(9*theta[0]^2-2)/(81*theta[0]^2*theta[1]^2+4), gamma[1] = gamma[1], gamma[2] = -(1/4)*gamma[1]^2*theta[0]^2/theta[1]^2, omega[0] = -2*(9*theta[1]^2+2)/(81*theta[0]^2*theta[1]^2+4), theta[0] = theta[0], theta[1] = theta[1]}

> 

eqpsi := eval(eq4, eqI)

psi(x, t) = gamma[1]*cos(theta[0]*(-2*t*(9*theta[1]^2+2)/(81*theta[0]^2*theta[1]^2+4)+x))-(1/4)*gamma[1]^2*theta[0]^2*exp(theta[1]*(2*t*(9*theta[0]^2-2)/(81*theta[0]^2*theta[1]^2+4)+x))/theta[1]^2+exp(-theta[1]*(2*t*(9*theta[0]^2-2)/(81*theta[0]^2*theta[1]^2+4)+x))

> 

eqphi := eval(eq2, eqpsi)

phi(x, t) = 3*(-gamma[1]*theta[0]*sin(theta[0]*(-2*t*(9*theta[1]^2+2)/(81*theta[0]^2*theta[1]^2+4)+x))-(1/4)*gamma[1]^2*theta[0]^2*exp(theta[1]*(2*t*(9*theta[0]^2-2)/(81*theta[0]^2*theta[1]^2+4)+x))/theta[1]-theta[1]*exp(-theta[1]*(2*t*(9*theta[0]^2-2)/(81*theta[0]^2*theta[1]^2+4)+x)))/(gamma[1]*cos(theta[0]*(-2*t*(9*theta[1]^2+2)/(81*theta[0]^2*theta[1]^2+4)+x))-(1/4)*gamma[1]^2*theta[0]^2*exp(theta[1]*(2*t*(9*theta[0]^2-2)/(81*theta[0]^2*theta[1]^2+4)+x))/theta[1]^2+exp(-theta[1]*(2*t*(9*theta[0]^2-2)/(81*theta[0]^2*theta[1]^2+4)+x)))

> 

simplify(eval(pde, eqphi))

0

> 

NULL

Download Sol1.mw

subsop...

dharr 7493
December 09 2024
1 6

I don't know what happened - it disappeared just after I prepared my answer. 

> 

restart;
> 

myproc := proc(x) sin(x) end proc;

proc (x) sin(x) end proc

> 

myproc2 := subsop(3=hfloat,eval(myproc))

proc (x) option hfloat; sin(x) end proc

 

 

Download proc.mw

MathematicalFunctions:-Get(plot, Zeta)...

dharr 7493
December 09 2024
2 1

The easiest way to extract the plot is

p := MathematicalFunctions:-Get(plot, Zeta)

Then plottools:-getdata(p[2][1]) extracts data from the plot, or use lprint(p[2][1]) to see the plot structure. As usual right-clicking on this plot can let you see various options, but I'm not aware of anything that can reconstruct the plot command. (Perhaps you were thinking about a plot component for that.)

In this case the plot can be generated by 

plots:-complexplot3d(Zeta(z), z = -4 - 10*I .. 4 + 40*I,
      view = [default, default, 0 .. 3], orientation = [-32, 68]);

though I didn't try to get all options exactly right.

(When I right-clicked on the plot from FunctionAdvisor, Maple disappeared completely, though not reproducibly!)

confirm on other worksheet...

dharr 7493
December 08 2024
1 1

I found that with one of that OPs worksheets (one I loaded in that thread) with a few second delays each keystroke (not 20 s, and not always), but it did have one large plot and multiple large expressions. In the past I have attributed those delays to working from a remote disk, but that wasn't the case today.

3 4 5 6 7 8 9 Last Page 5 of 76