dharr

Dr. David Harrington

8200 Reputation

22 Badges

20 years, 335 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 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

@C_R has pointed out many of the things I might have mentioned. Here a few more comments along the same lines.

Q1. To me, your formulation uses matrices, but is not very matrixy. That's hard to define and probably largely intuitive on my part. Your matrices are almost all involved in dot products or quadratic forms (a^T.M.x), which produce scalars that you then involve in non-matrix equations. The power of matrices is in setting up matrix equations (I suppose a few vectors are allowed) and solving the matrix equations - often just a linear equation solving of Ax = b, but perhaps AX + XB = C, AXB = C or some nonlinear ones.

Specifically (but vaguely), the entries of r don't look like they have a simple matrix origin - we have square roots of entries, and while the denominators themselves might be matrixy, division is not a matrix thing so the denominators might better be in scalars outside the matrix. For the numerators, there are only diagonal elements of S, so where did the off-diagonals go?

The key point was given by @C_R: "but for that we have to see the structure of the equations in the first place".

Q2. I would say no here. Try to set things up without components, even if later you need to be explicit about them. For example to take a small part of your equations. q[1] = w^T.S.e[1], q[2] = w^T.S.e[2], q[3] = w^T.S.e[3] . So q[1], q[2], q[3] are the first, second and third components of w^T.S, so this is just q^T = w^T.S or q = S^T.w. The w[i]^2 seem to make this type of analysis harder (not matrixy?).

According to the help page, rsolve requires A(n), not A[n]. (But see @acer's answer.)

restart

sol := rsolve({A(n+1)*(5*A(n)+2) = A(n), A(1) = 1}, A)

1/(3*2^n-5)

eval(sol, n = 3)

1/19

seq(sol, n = 1 .. 10)

1, 1/7, 1/19, 1/43, 1/91, 1/187, 1/379, 1/763, 1/1531, 1/3067

``

NULL

Download ask_recrusive_sequence_of_A[n].mw

It seems to me that case (I) covers all sign combinations of lambda and mu. Case (II) looks incorrect to me - just putting the negative in doesn't work unless lambda and mu are redefined.

restart

de1 := diff(F(eta), eta) = lambda*G(eta); de2 := diff(G(eta), eta) = mu*F(eta)

diff(F(eta), eta) = lambda*G(eta)

diff(G(eta), eta) = mu*F(eta)

ans := dsolve({de1, de2}, {F(eta), G(eta)})

{F(eta) = c__1*exp(lambda^(1/2)*mu^(1/2)*eta)+c__2*exp(-lambda^(1/2)*mu^(1/2)*eta), G(eta) = mu^(1/2)*(c__1*exp(lambda^(1/2)*mu^(1/2)*eta)-c__2*exp(-lambda^(1/2)*mu^(1/2)*eta))/lambda^(1/2)}

ans2 := convert(ans, trigh)

{F(eta) = (c__1+c__2)*cosh(lambda^(1/2)*mu^(1/2)*eta)+(c__1-c__2)*sinh(lambda^(1/2)*mu^(1/2)*eta), G(eta) = mu^(1/2)*(c__1-c__2)*cosh(lambda^(1/2)*mu^(1/2)*eta)/lambda^(1/2)+mu^(1/2)*(c__1+c__2)*sinh(lambda^(1/2)*mu^(1/2)*eta)/lambda^(1/2)}

odetest(ans2, [de1, de2])

[0, 0]

ans3 := subs(lambda = -lambda, mu = -mu, ans2)

{F(eta) = (c__1+c__2)*cosh((-lambda)^(1/2)*(-mu)^(1/2)*eta)+(c__1-c__2)*sinh((-lambda)^(1/2)*(-mu)^(1/2)*eta), G(eta) = (-mu)^(1/2)*(c__1-c__2)*cosh((-lambda)^(1/2)*(-mu)^(1/2)*eta)/(-lambda)^(1/2)+(-mu)^(1/2)*(c__1+c__2)*sinh((-lambda)^(1/2)*(-mu)^(1/2)*eta)/(-lambda)^(1/2)}

simplify(odetest(ans3, [de1, de2]))

[2*(-lambda)^(1/2)*((c__1-c__2)*cosh((-lambda)^(1/2)*(-mu)^(1/2)*eta)+sinh((-lambda)^(1/2)*(-mu)^(1/2)*eta)*(c__1+c__2))*(-mu)^(1/2), -2*mu*((c__1+c__2)*cosh((-lambda)^(1/2)*(-mu)^(1/2)*eta)+(c__1-c__2)*sinh((-lambda)^(1/2)*(-mu)^(1/2)*eta))]

Take a specific numerical example

params := {c__1 = 1, c__2 = 1, lambda = -1, mu = -1}

{c__1 = 1, c__2 = 1, lambda = -1, mu = -1}

ans3num := eval(ans3, params); desnum := eval({de1, de2}, params)

{F(eta) = 2*cosh(eta), G(eta) = 2*sinh(eta)}

{diff(F(eta), eta) = -G(eta), diff(G(eta), eta) = -F(eta)}

simplify(eval(desnum, ans3num))

{2*cosh(eta) = -2*cosh(eta), 2*sinh(eta) = -2*sinh(eta)}

NULL

Download dsolve.mw

I don't understand the seed issue, but have some other answers.

My understanding is that Maple does some preprocessing and then calls the actual NAG compiled code.

Numerical_integration_Using_MC_methods.mw

 

I didn't have the patience to wait for your code to finish to diagnose the results. (I'd guess Equal doesn't try too hard for very complicated entries.) The main issue here is that the inverse of a matrix with symbolic entries is too complicated to be useful, except for very small matrices. But if you have numeric entries, you can easily verify the relationship works. Here's how I would do that. I used integer entries to avoid numerical roundoff issues, but for floating point entries it can be faster.

restart

with(LinearAlgebra)

`⊗` := KroneckerProduct

A := RandomMatrix(10, 10); B := RandomMatrix(7, 7)

Matrix(10, 10, {(1, 1) = -38, (1, 2) = 12, (1, 3) = -82, (1, 4) = 82, (1, 5) = 22, (1, 6) = 76, (1, 7) = 31, (1, 8) = -16, (1, 9) = -98, (1, 10) = -4, (2, 1) = 91, (2, 2) = 45, (2, 3) = -70, (2, 4) = 72, (2, 5) = 14, (2, 6) = -44, (2, 7) = -50, (2, 8) = -9, (2, 9) = -77, (2, 10) = 27, (3, 1) = -1, (3, 2) = -14, (3, 3) = 41, (3, 4) = 42, (3, 5) = 16, (3, 6) = 24, (3, 7) = -80, (3, 8) = -50, (3, 9) = 57, (3, 10) = 8, (4, 1) = 63, (4, 2) = 60, (4, 3) = 91, (4, 4) = 18, (4, 5) = 9, (4, 6) = 65, (4, 7) = 43, (4, 8) = -22, (4, 9) = 27, (4, 10) = 69, (5, 1) = -23, (5, 2) = -35, (5, 3) = 29, (5, 4) = -59, (5, 5) = 99, (5, 6) = 86, (5, 7) = 25, (5, 8) = 45, (5, 9) = -93, (5, 10) = 99, (6, 1) = -63, (6, 2) = 21, (6, 3) = 70, (6, 4) = 12, (6, 5) = 60, (6, 6) = 20, (6, 7) = 94, (6, 8) = -81, (6, 9) = -76, (6, 10) = 29, (7, 1) = -26, (7, 2) = 90, (7, 3) = -32, (7, 4) = -62, (7, 5) = -95, (7, 6) = -61, (7, 7) = 12, (7, 8) = -38, (7, 9) = -72, (7, 10) = 44, (8, 1) = 30, (8, 2) = 80, (8, 3) = -1, (8, 4) = -33, (8, 5) = -20, (8, 6) = -48, (8, 7) = -2, (8, 8) = -18, (8, 9) = -2, (8, 10) = 92, (9, 1) = 10, (9, 2) = 19, (9, 3) = 52, (9, 4) = -68, (9, 5) = -25, (9, 6) = 77, (9, 7) = 50, (9, 8) = 87, (9, 9) = -32, (9, 10) = -31, (10, 1) = 22, (10, 2) = 88, (10, 3) = -13, (10, 4) = -67, (10, 5) = 51, (10, 6) = 9, (10, 7) = 10, (10, 8) = 33, (10, 9) = -74, (10, 10) = 67})

Matrix(%id = 36893491254531674588)

Equal(`⊗`(1/A, 1/A), 1/`⊗`(A, A)); Equal(`⊗`(1/A, 1/B), 1/`⊗`(A, B))

true

true

NULL

Download Kronecker.mw

Maple uses indexed RootOfs to distinguish the roots - these are in a specific order - see the ?RootOf/indexed help page. Solve gives a generic RootOf, but you can use allvalues to get all the values, and it tracks the different possibilities by default (see the dependent option). An example for a polynomial system is given below.

Edit: For commands like irreduc, the default field for the coefficients is the rationals, but Maple can handle extension fields through evala and its subcommands. As @Carl Love notes, the indexed RootOfs are for all the complex roots.

I'm not sure I have fully understand your question. If you give a more specific example of what you want to achieve or what didn't work as you expected, I'll try to clarify.

Download Rootof.mw

Maple has some algorithm to find the best set of points to use for plotting a function (not just at equal intervals across the x-axis), and here it seems to be getting confused. You can fix this by using adaptive=false or numpoints=2000.

I'm not sure how integration by parts would apply here. You can just exchange the sum and integration:

restart

inte_eq := int(-(sum(B[i]*beta[i]* exp(-beta[i] * (t-tau)),i=1..n)),tau=0..t);

int(-(sum(B[i]*beta[i]*exp(-beta[i]*(t-tau)), i = 1 .. n)), tau = 0 .. t)

summand := Student:-Calculus1:-Summand(inte_eq)[];

B[i]*beta[i]*exp(-beta[i]*(t-tau))

integ := int(summand, tau=0..t);

-B[i]*exp(-beta[i]*t)+B[i]

ans := simplify(-sum(integ, i=1..n));

sum(B[i]*(-1+exp(-beta[i]*t)), i = 1 .. n)

Download integral.mw

@one man @vv I found plex wasn't much slower than tdeg and put all the pieces together. RationalUnivariateRepresentation should allow fsolve to work on a single polynomial but there are numerical and speed issues.

(I originally tried these polynomial systems with SolveTools:-PolynomialSystem, but it was much slower and I gave up.)

nonlinear2.mw

Edit: I wasn't substituting back into the correct original equations; now fixed.

Well, this still uses the tickmarks option, but I think it achieves what you want. 

restart; with(plots); with(DEtools)

Ode1 := diff(y(t), t, t) = cos(t)

diff(diff(y(t), t), t) = cos(t)

IC := (D(y))(0) = 0, y(0) = 1

(D(y))(0) = 0, y(0) = 1

gsol := dsolve(Ode1, y(t))

y(t) = -cos(t)+c__1*t+c__2

exactsol := dsolve({IC, Ode1}, y(t), numeric, output = listprocedure)

odeplot(exactsol, [t, y(t)], 0 .. 4*Pi, color = blue, title = "Solution to the Differential Equation", labels = ["t", "y(t)"], tickmarks = [piticks, default])

 

NULL

Download directly_integrable.mw

Depends on the exact type of question. But Maple's geometry package might be able to help. Or you can set up relevant equations in Maple and then solve them.

You can't use B[1] or B[2] as parameter names - use B1, B2 instead. Note: if you were using 1-D input, you would receive a more intelligible error message about unexpected [ and the cursor would point you to the offending [

In your first approach notice that there is an infinity at the end, so Maple is telling you the result is +/- infinity. Within the signum function, there two terms that might have different signs. Assuming n>1 would not be enough to resolve this. To check the infinity part, you could put in some numerical values for all but gamma, and plot vs gamma.

So for the second part, you will get the same result, but now the X and lambda are more complicated, and deciding the signs of things is also more complicated, so it will take much longer to decide if it can.

Your question is not very specific about what you want to do, so maybe this doesn't answer your question. The LinearAlgebra:-Generic package allows you to define generic field operations, but for extensions of Q you can do many basic things just with the LinearAlgebra package.

restart

alias(alpha = RootOf(x^2-2))

alpha

M := Matrix(2, 2, [1+alpha, alpha, 2, 2])

Matrix(%id = 36893490439769270860)

evala(1/M)

Matrix(%id = 36893490439769228340)

M.M; evala(%)

Matrix(2, 2, {(1, 1) = (1+RootOf(_Z^2-2))^2+2*RootOf(_Z^2-2), (1, 2) = (1+RootOf(_Z^2-2))*RootOf(_Z^2-2)+2*RootOf(_Z^2-2), (2, 1) = 6+2*RootOf(_Z^2-2), (2, 2) = 2*RootOf(_Z^2-2)+4})

Matrix(%id = 36893490439840085884)

NULL

Download Field.mw

"I tried converting an outupt to a string which removed all the backslashes in the file path. Maybe there is a dedicated command that I overlooked."

SplitPath does this, in a way that doesn't care whether the separator is \\ or /

Other related manipulations are in the FileTools package

f:=currentdir(); #some path

"C:\Users\dharr\OneDrive - University of Victoria\Desktop"

FileTools:-SplitPath(f);

["C:", "Users", "dharr", "OneDrive - University of Victoria", "Desktop"]

NULL

Download SpliPath.mw

First 10 11 12 13 14 15 16 Last Page 12 of 81