dharr

Dr. David Harrington

8270 Reputation

22 Badges

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

You could use string manipulations to automate a file conversion for the simple stuff as in the example file below, and then hand edit the file for the things that don't easily convert. (integrate is a synonym for int, so you could leave that). Of course you could spend a lot of time on a more sophisticated conversion process, but I'm guessing it is easier to do an approximate conversion and then fix things up in Maple. I "read" the modified file in the example, which means that Maple parses the commands, but it might be easier to copy from the converted file to an input region.

MaxToMaple.mw

Spelling error

wiht(LinearAlgebra) should be with(LinearAlgebra)

There is no general formula in terms of radicals for a general 5th degree polynomial. For some cases of coefficients there will be solutions, roots(x^5-1); is a simple example.

As @Carl Love says, walks are enumerated by the adjacency matrix, and can be used to do this.

restart

with(GraphTheory); with(LinearAlgebra)

G := AddVertex(PathGraph(5), [6, 7]); AddEdge(G, {{3, 7}, {4, 6}, {6, 7}}); DrawGraph(G, layout = circle, size = [250, 250])

GRAPHLN(undirected, unweighted, [1, 2, 3, 4, 5, 6, 7], Array(1..7, {(1) = {2}, (2) = {1, 3}, (3) = {2, 4}, (4) = {3, 5}, (5) = {4}, (6) = {}, (7) = {}}), `GRAPHLN/table/2`, 0)

n := NumberOfVertices(G)

7

A := AdjacencyMatrix(G)

B := `~`[`*`](A, Matrix(n, n, symbol = omega))

Matrix(%id = 36893490081594197644)

B^2 enumerates all walks of length 2. Look at the diagonals to find the ones that return. So for example from vertex 3 back to 3 there are three such walks, and the omegas tell where we went.

R := Diagonal(B^2)

Vector[column](%id = 36893490081591008900)

So if we collect the indices in each term into a set, we can see when we have a set of all vertices. Here after 8 steps we aren't done yet for vertex 1.

makesets:=y->evalindets(evalindets[2](expand(y),`^`,1,op),`*`,x->map(u->`if`(type(u,integer),NULL,op(u)),{op(x)})):

a:=(B^8)[1,1];
makesets(a);

(omega[1, 2]^2*omega[2, 1]^2+omega[1, 2]*omega[2, 1]*omega[2, 3]*omega[3, 2])^2+(omega[1, 2]^2*omega[2, 1]*omega[2, 3]+omega[1, 2]*omega[2, 3]*(omega[2, 3]*omega[3, 2]+omega[3, 4]*omega[4, 3]+omega[3, 7]*omega[7, 3]))*(omega[3, 2]*omega[2, 1]^2*omega[1, 2]+(omega[2, 3]*omega[3, 2]+omega[3, 4]*omega[4, 3]+omega[3, 7]*omega[7, 3])*omega[3, 2]*omega[2, 1])+omega[1, 2]*omega[2, 3]*omega[3, 4]*omega[4, 5]*omega[5, 4]*omega[4, 3]*omega[3, 2]*omega[2, 1]+omega[1, 2]*omega[2, 3]*(omega[3, 4]*omega[4, 6]+omega[3, 7]*omega[7, 6])*(omega[4, 3]*omega[6, 4]+omega[6, 7]*omega[7, 3])*omega[3, 2]*omega[2, 1]

{1, 2, 3, 4, 7}+{1, 2, 3, 4, 5}+{1, 2, 3, 4, 6}+2*{1, 2, 3, 4, 6, 7}+{1, 2, 3, 6, 7}+{1, 2, 3, 4}+{1, 2, 3, 7}+{1, 2, 3}+{1, 2}

So here is a procedure for finding the length of the shortest returning walk. For a directed graph it might not terminate.

WalkLength:=proc(G::Graph,v)
  uses GraphTheory;
  local x,y,u,vv,A,B,makesets,n,vertset,i,omega;
  if not member(v,Vertices(G),'vv') then error "vertex not in graph" end if;
  if not IsConnected(G) then return infinity end if;
  makesets:=y->evalindets(evalindets[2](expand(y),`^`,1,op),
                      `*`,x->map(u->`if`(type(u,integer),NULL,op(u)),{op(x)})):
  n:=NumberOfVertices(G);
  vertset:={$(1..n)};
  A:=AdjacencyMatrix(G)*~Matrix(n,n,symbol=omega);
  B:=copy(A);
  for i from 2 do
    B:=B.A;
  until has(makesets(B[vv,vv]),[vertset]);
  i;
end proc:

WalkLength(G, 1)

10

NULL

Download Walks3.mw

 

I think this is fairly robust.

get_exponent:=(expr,x::name)->diff(expr,x)/expr*x:
get_exponent(x^2,x);
get_exponent(x^(n-1),x);
get_exponent(3*x^(n+2),x);
get_exponent(3*x[2]^k*x[3]^(n+2),x[3]);

    
                               2

                             n - 1

                             n + 2

                             n + 2

Once you entered floating point values in your Matrix, the implication is that you have given up on symbolic solutions. Since the Matrix wasn't obviously of a structure that has real eigenvalues, the algorithm used was a general one that generated complex values. As @acer pointed out, you can then make them look nicer, but the algorithm used is the same. For a symmetric Matrix, if you use shape=symmetric, then Maple will choose an algorithm that will produce real eigenvalues and eigenvectors (try it on A:=Matrix(2,2,[[0.8,0.3],[0.3,0.7]],shape=symmetric);).

But if you enter exact values in the first place, then Maple will do the whole calculation exactly. In general, this is the power of Maple over a regular programming language. Here I converted the values you entered to rationals to avoid re-entering, but you could have entered these directly.

Download Eigenvectors.mw

Yes, my Maple 2023.0 on Windows 10 shows this behaviour.

Try Insert->Table and choose 1 row, 1 column, Title: "Fig. 1" and tick Show caption. Then Insert->Image->From file and choose the file. After you have it, right-click (cmd-click on Mac) and choose Table -> Properties, Table size 50% and the image will reduce size 50% keeping the aspect ratio. That menu also has a "Caption" tab which gives several ways to specify the caption. As far as I know there isn't automatic renumbering.

(As you probably already found, you can insert images in a text area but resizing with the handles doesn't maintain the aspect ratio.)

(Notice that _B2 stands for a binary variable with values 0 or 1, but _Z1 and _Z2 are (different) integers.)

restart

sol := solve({-sin(alpha)-5*cos(beta)+5*tan(lambda) = 0, -sin(alpha)+5*cos(beta)+3*tan(lambda) = 0, 2*sin(alpha)+2*cos(beta)+3*tan(lambda) = 0, 0 <= alpha and alpha <= 2*Pi, 0 <= beta and beta <= 2*Pi, 0 <= lambda and lambda <= 2*Pi}, {alpha, beta, lambda}, allsolutions, explicit)

{alpha = Pi*_Z1, beta = (1/2)*Pi-Pi*_Z1-_B2*Pi+2*_B2*Pi*_Z1+2*Pi*_Z2, lambda = -arctan((2/3)*sin(Pi*_Z1)+(2/3)*cos((1/2)*Pi-Pi*_Z1-_B2*Pi+2*_B2*Pi*_Z1+2*Pi*_Z2))+Pi*_Z4}

undervars := indets(sol, suffixed({_B, _Z}))

{_B2, _Z1, _Z2, _Z4}

subs(undervars[1] = n, undervars[2] = n, undervars[3] = n, sol)

{alpha = n*Pi, beta = (1/2)*Pi+2*n^2*Pi, lambda = arctan((2/3)*sin(2*n^2*Pi))+Pi*_Z4}

NULL

Download underline.mw

Your code has a list M := [0, 0.5, 1.0, 1.5, 5] which you insert into other expressions as though it were a single number. So even for F[3] that isn't sensible. This ultimately leads to the error message.

mul(map(mul,ListTools:-Collect(A)))

CN and N1 are 

CN := (G, u) -> add(map2(GraphTheory:-Degree, G, GraphTheory:-Neighborhood(G, u, 'closed')));
N1 := G -> add(map(e->(CN(G,e[1])+CN(G,e[2]))/2, convert(GraphTheory:-Edges(G),list)));

The others are straightforward modifications.

This gives the answer in term of four algebraic numbers, each of which can be written in terms of the first root, which I think is more in the spirit of Galois theory than radicals. I also show the minimal polynomials that show the Galois group is order 72, But I'm not using the Galois group for much of this, as the OP asked for. Sorry if this is already obvious to the OP; perhaps we can have some clarification if the answers here are off track.

restart

_local(gamma); poly := x^6-3*x+3

gamma

x^6-3*x+3

Doesn't factor in Q

irreduc(poly)

true

Let alpha be the first root (index=1 by default)

alias(alpha = RootOf(poly, x))

alpha

alpha is the only root that can be written in terms of alpha. Notice that the rest factors, i.e. the quintic is reducible in  Q(alpha)

p6alpha := evala(Factor(poly, {alpha}))

(x^3+((6/17)*alpha-(6/17)*alpha^5-(14/17)*alpha^4-(10/17)*alpha^3-(12/17)*alpha^2+15/17)*x^2+((6/17)*alpha^2-33/17+(3/17)*alpha^5+(7/17)*alpha^4+(5/17)*alpha^3-(3/17)*alpha)*x+(5/17)*alpha^3+18/17+(7/17)*alpha^4+(6/17)*alpha^2+(3/17)*alpha^5-(3/17)*alpha)*((11/17)*alpha^5+(3/17)*alpha^4+(7/17)*alpha^3+(5/17)*alpha^2+(6/17)*alpha-36/17+((6/17)*alpha^5+(14/17)*alpha^4+(10/17)*alpha^3+(12/17)*alpha^2+(11/17)*alpha-15/17)*x+x^2)*(x-alpha)

Collect the factors

p3alpha, p2alpha, p1alpha := op(p6alpha)

x^3+((6/17)*alpha-(6/17)*alpha^5-(14/17)*alpha^4-(10/17)*alpha^3-(12/17)*alpha^2+15/17)*x^2+((6/17)*alpha^2-33/17+(3/17)*alpha^5+(7/17)*alpha^4+(5/17)*alpha^3-(3/17)*alpha)*x+(5/17)*alpha^3+18/17+(7/17)*alpha^4+(6/17)*alpha^2+(3/17)*alpha^5-(3/17)*alpha, (11/17)*alpha^5+(3/17)*alpha^4+(7/17)*alpha^3+(5/17)*alpha^2+(6/17)*alpha-36/17+((6/17)*alpha^5+(14/17)*alpha^4+(10/17)*alpha^3+(12/17)*alpha^2+(11/17)*alpha-15/17)*x+x^2, x-alpha

Let beta be the first root of the quadratic; the other can be expressed in terms of alpha and beta

p2alpha; alias(beta = RootOf(p2alpha, x))

(11/17)*alpha^5+(3/17)*alpha^4+(7/17)*alpha^3+(5/17)*alpha^2+(6/17)*alpha-36/17+((6/17)*alpha^5+(14/17)*alpha^4+(10/17)*alpha^3+(12/17)*alpha^2+(11/17)*alpha-15/17)*x+x^2

alpha, beta

p2alphabeta := evala(Factor(p2alpha, {alpha, beta}))

(x+(6/17)*alpha^5+(14/17)*alpha^4+(10/17)*alpha^3+(12/17)*alpha^2+(11/17)*alpha-15/17+beta)*(x-beta)

The cubic can't be done in terms of alpha and beta, so let gamma be its first root.

p3alpha; irreduc(p3alpha, {alpha, beta}); alias(gamma = RootOf(p3alpha, x))

x^3+((6/17)*alpha-(6/17)*alpha^5-(14/17)*alpha^4-(10/17)*alpha^3-(12/17)*alpha^2+15/17)*x^2+((6/17)*alpha^2-33/17+(3/17)*alpha^5+(7/17)*alpha^4+(5/17)*alpha^3-(3/17)*alpha)*x+(5/17)*alpha^3+18/17+(7/17)*alpha^4+(6/17)*alpha^2+(3/17)*alpha^5-(3/17)*alpha

true

alpha, beta, gamma

It factors into a quadratic and the linear term, so we need another algebraic number as the first root of the quadratic

p3alphabetagamma := evala(Factor(p3alpha, {alpha, beta, gamma}))

(x^2+(-(6/17)*alpha^5-(14/17)*alpha^4-(10/17)*alpha^3-(12/17)*alpha^2+(6/17)*alpha+15/17+gamma)*x+(3/17)*alpha^5+(7/17)*alpha^4+(5/17)*alpha^3+(6/17)*alpha^2-(3/17)*alpha-33/17-(6/17)*gamma*alpha^5-(14/17)*gamma*alpha^4-(10/17)*gamma*alpha^3-(12/17)*gamma*alpha^2+(6/17)*gamma*alpha+(15/17)*gamma+gamma^2)*(x-gamma)

p2alphabetagamma, p1alphabetagamma := op(p3alphabetagamma)

x^2+(-(6/17)*alpha^5-(14/17)*alpha^4-(10/17)*alpha^3-(12/17)*alpha^2+(6/17)*alpha+15/17+gamma)*x+(3/17)*alpha^5+(7/17)*alpha^4+(5/17)*alpha^3+(6/17)*alpha^2-(3/17)*alpha-33/17-(6/17)*gamma*alpha^5-(14/17)*gamma*alpha^4-(10/17)*gamma*alpha^3-(12/17)*gamma*alpha^2+(6/17)*gamma*alpha+(15/17)*gamma+gamma^2, x-gamma

alias(delta = RootOf(p2alphabetagamma, x))

alpha, beta, gamma, delta

evala(Factor(p2alphabetagamma, {alpha, beta, delta, gamma}))

(x-(6/17)*alpha^5-(14/17)*alpha^4-(10/17)*alpha^3-(12/17)*alpha^2+(6/17)*alpha+15/17+gamma+delta)*(x-delta)

Can put it all together

evala(AFactor(poly))

(x-delta)*(x+(6/17)*alpha^5+(14/17)*alpha^4+(10/17)*alpha^3+(12/17)*alpha^2+(11/17)*alpha-15/17+beta)*(x-beta)*(x-alpha)*(x-gamma)*(x-(6/17)*alpha^5-(14/17)*alpha^4-(10/17)*alpha^3-(12/17)*alpha^2+(6/17)*alpha+15/17+gamma+delta)

This can be done all in one step with PolynomialTools:-Split - the splitting field K with four algebraic numbers is returned (but figuring out the nested RootOf's is painful).

PolynomialTools:-Split(poly, x, 'K'); K

(x-delta)*(x+(6/17)*alpha^5+(14/17)*alpha^4+(10/17)*alpha^3+(12/17)*alpha^2+(11/17)*alpha-15/17+beta)*(x-beta)*(x-alpha)*(x-gamma)*(x-(6/17)*alpha^5-(14/17)*alpha^4-(10/17)*alpha^3-(12/17)*alpha^2+(6/17)*alpha+15/17+gamma+delta)

{alpha, beta, delta, gamma}

evala(Minpoly(alpha, x)); dalpha := degree(%, x)

x^6-3*x+3

6

evala(Minpoly(beta, x, {alpha})); dbeta := degree(%, x)

x^2+(6/17)*x*alpha^5+(14/17)*alpha^4*x+(10/17)*alpha^3*x+(12/17)*alpha^2*x+(11/17)*x*alpha-(15/17)*x+(11/17)*alpha^5+(3/17)*alpha^4+(7/17)*alpha^3+(5/17)*alpha^2+(6/17)*alpha-36/17

2

evala(Minpoly(gamma, x, {alpha, beta})); dgamma := degree(%, x)

x^3-(6/17)*alpha^5*x^2-(14/17)*alpha^4*x^2-(10/17)*alpha^3*x^2-(12/17)*alpha^2*x^2+(6/17)*alpha*x^2+(15/17)*x^2+(3/17)*x*alpha^5+(7/17)*alpha^4*x+(5/17)*alpha^3*x+(6/17)*alpha^2*x-(3/17)*x*alpha-(33/17)*x+(3/17)*alpha^5+(7/17)*alpha^4+(5/17)*alpha^3+(6/17)*alpha^2-(3/17)*alpha+18/17

3

evala(Minpoly(delta, x, {alpha, beta, gamma})); ddelta := degree(%, x)

x^2-(6/17)*x*alpha^5-(14/17)*alpha^4*x-(10/17)*alpha^3*x-(12/17)*alpha^2*x+(6/17)*x*alpha+(15/17)*x+x*gamma+(3/17)*alpha^5+(7/17)*alpha^4+(5/17)*alpha^3+(6/17)*alpha^2-(3/17)*alpha-33/17-(6/17)*gamma*alpha^5-(14/17)*gamma*alpha^4-(10/17)*gamma*alpha^3-(12/17)*gamma*alpha^2+(6/17)*gamma*alpha+(15/17)*gamma+gamma^2

2

Order of Galois group.

dalpha*dbeta*dgamma*ddelta

72

NULL

Download Splits.mw

That is a surprising result, and sufficiently far enough from the correct answer that I would call it a bug, even though it is just a numerical issue. Certainly since the matrix is "nice" (symmetric, diagonally dominant). As @vv says it depends on the precision, but certainly you don't need to go to 40 to get OK results. At Digits=14, M_MP.M_MP - C is recognizably close to zero and at Digits=16, the entries are all about 1e-6 or 1e-7.

In the general case, the algorithm for matrix functions is nontrivial, but Maple often has special cases for shape=symmetric or shape = hermitian, so then I think just applying the function to the eigenvalues as you did is straightforward.

@vv's comment about the determinant might be true (I agree), but it seems to imply that the reason that <0,1;0,0> doesn't have a square root is because its determinant is zero. That has more to do with the fact that it is non-diagonalizable; for example <1,0;0,0> (positive semidefinite) and <1,1;0,0> both have square roots and have determinant zero. 

For the edge separator, Maple's IsCutSet can be used - it just checks for an increase in the number of components. Checking for a Cut is more complicated since just checking for 2 components can fail to detect extra removed edges that are in one of the components. Here is one way to do it.

[Edit - this actually tests for a cutset, which is a minimal Cut, so I've renamed it IsMinCut; see below for the test for a Cut]

restart;

IsMinCut returns true if the edges are a minimal cut (cutset), i.e., removal of the edges cuts the graph into two components, and no edge can be added back without reconnecting the graph.

This version doesn't handle loops for efficiency, but could be modified to use the underlying graph and ignore loops (provided edges doesn't contain a loop).

IsMinCut:=proc(G::GRAPHLN,edges::set)
       local partition,G1,G2;
       uses GraphTheory;
       if not IsConnected(G) or IsDirected(G) then error "Graph must be connected and undirected" end if;
       partition:=ConnectedComponents(DeleteEdge(G,edges,'inplace'=false));
       if nops(partition) <> 2 then return false end if;
       G1:=InducedSubgraph(G,partition[1]);
       G2:=InducedSubgraph(G,partition[2]);
       if (Edges(G) minus Edges(G1) minus Edges(G2)) = edges then
          true
       else
          false
       end if;
end proc:

with(GraphTheory):
with(combinat):
G:=CompleteGraph(3,3):
edge:=choose(Edges(G), 7):
`and`(seq(IsCutSet(G,s),s in edge)); # any 7 edges cuts the graph into 2 or more components

true

`or`(seq(IsMinCut(G,s),s in edge)); # no set of 7 edges is a cut

false

DrawGraph(G,size=[250,250]);

edges:={{1,4},{1,5},{1,6},{3,6}}; # {3,6} edge is extra
IsMinCut(G,edges);
 

{{1, 4}, {1, 5}, {1, 6}, {3, 6}}

false

edges:={{1,4},{1,5},{1,6}};
IsMinCut(G,edges);

{{1, 4}, {1, 5}, {1, 6}}

true

edges:={{1,4},{1,5},{1,6},{3,6},{3,5},{3,4}}; #3 components after cut
IsMinCut(G,edges);

{{1, 4}, {1, 5}, {1, 6}, {3, 4}, {3, 5}, {3, 6}}

false

NULL

Download IsMinCut.mw

First 32 33 34 35 36 37 38 Last Page 34 of 82