Advanced Search

dharr

Dr. David Harrington

8917 Reputation

22 Badges

21 years, 139 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

Biorthogonality...

dharr 8917
March 13 2024
2 2

@mehdi jafari You asked about the relationship between left and right eigenvectors, and getting a diagonal matrix from them. They are biorthogonal, meaning that left eigenvectors are orthogonal to right eigenvectors corresponding to different eigenvalues. So the product of the matrices of left eigenvectors (as row vectors) and right eigenvectors (as column vectors) is diagonal.

[Edit - 2023 version uploaded] Rerun the worksheet in your version of Maple; the rtable output from Maple 2024 doesn't display correctly in earlier versions (or here on Mapleprimes, SCR submitted).

> 

restart
> 

with(LinearAlgebra)
> 

A := Matrix(3, 3, [1, 0, 1, 1, 1, 0, 2, 0, 1])

Matrix(%id = 36893490995083995364)

Right (column) eigenvectors corresponding to increasing eigenvalues (any order for the eigenvalues will do, so long as it is consistent for left and right eigenvectors)

> 

revs, rvecs1 := Eigenvectors(A); evs, rperm := sort(revs, output = [sorted, permutation], key = evalf); rvecs := rvecs1[() .. (), rperm]

evs, rperm := Vector(3, {(1) = 1-sqrt(2), (2) = 1, (3) = 1+sqrt(2)}), [3, 1, 2]

Matrix(%id = 36893490995083958140)

Left (row) eigenvectors corresponding to increasing eigenvals

> 

levs, lvecs1 := Eigenvectors(LinearAlgebra:-Transpose(A)); evs, lperm := sort(levs, output = [sorted, permutation], key = evalf); lvecs := lvecs1*LinearAlgebra:-Transpose([() .. (), lperm])

evs, lperm := Vector(3, {(1) = 1-sqrt(2), (2) = 1, (3) = 1+sqrt(2)}), [2, 3, 1]

Matrix(%id = 36893490995083988852)

Biorthogonality means that left eigenvectors are orthogonal to right eigenvectors corresponding to different eigenvalues.

e.g, left eigenvector for 1 is orthogonal to right eigenvalues for 1-srqt(2) and 1+sqrt(2)

> 

lvecs[1, () .. ()], rvecs[() .. (), 2]; lvecs[1, () .. ()].rvecs[() .. (), 2]; lvecs[1, () .. ()].rvecs[() .. (), 3]

Vector[row](%id = 36893490995083967172), Vector[column](%id = 36893490995083967292)

0

0

All possibilities

> 

lvecs.rvecs

Matrix(%id = 36893490995083956340)

> 

NULL

Download Biorthogonal.mw

 

normalization...

dharr 8917
March 12 2024
3 0

The eigenvectors from Eigenvectors are arbitrarily scaled. In this case they happen to be normalized but differ in signs. Once you account for this (and a forgotten square root), there is agreement. Notice you are also relying on the orders of the eigenvalues and singular values being consistent.

 

WHY DOES THE TRANSPOSE OF Vt DETERMINED BY THE FUNCTION SingularValues NOT AGREE WITH THE evectors CALCULATED BY THE FUNCTION Eigenvectors?

 

> 

restart:
with(LinearAlgebra):
> 

X := Matrix([[8.79,9.93,9.83,5.45,3.16],
           [6.11,6.91,5.04,-0.27,7.98],
           [-9.15,-7.93,4.86,4.85,3.01],
           [9.57,1.64,8.83,0.74,5.80],
           [-3.49,4.02,9.80,10.00,4.27],
           [9.84,0.15,-8.99,-6.02,-5.31]],
           datatype=float[8],order=Fortran_order);

Matrix(6, 5, {(1, 1) = 8.79, (1, 2) = 9.93, (1, 3) = 9.83, (1, 4) = 5.45, (1, 5) = 3.16, (2, 1) = 6.11, (2, 2) = 6.91, (2, 3) = 5.04, (2, 4) = -.27, (2, 5) = 7.98, (3, 1) = -9.15, (3, 2) = -7.93, (3, 3) = 4.86, (3, 4) = 4.85, (3, 5) = 3.01, (4, 1) = 9.57, (4, 2) = 1.64, (4, 3) = 8.83, (4, 4) = .74, (4, 5) = 5.8, (5, 1) = -3.49, (5, 2) = 4.02, (5, 3) = 9.8, (5, 4) = 10.0, (5, 5) = 4.27, (6, 1) = 9.84, (6, 2) = .15, (6, 3) = -8.99, (6, 4) = -6.02, (6, 5) = -5.31})

 (1)
> 

U,S,Vt:= SingularValues(X, output = ['U', 'S', 'Vt'],thin=true)

U, S, Vt := Matrix(6, 5, {(1, 1) = -.591142376412437, (1, 2) = .263167814714055, (1, 3) = .355430173862827, (1, 4) = .314264362726927, (1, 5) = .229938315364748, (2, 1) = -.397566794202426, (2, 2) = .243799027926330, (2, 3) = -.222390000685446, (2, 4) = -.753466150953458, (2, 5) = -.363589686697497, (3, 1) = -0.334789690624459e-1, (3, 2) = -.600272580693583, (3, 3) = -.450839268922308, (3, 4) = .233449657244714, (3, 5) = -.305475732747932, (4, 1) = -.429706903137018, (4, 2) = .236166806281125, (4, 3) = -.685862863873811, (4, 4) = .331860018200310, (4, 5) = .164927634884511, (5, 1) = -.469747921566658, (5, 2) = -.350891398883703, (5, 3) = .387444603099673, (5, 4) = .158735559582156, (5, 5) = -.518257437353535, (6, 1) = .293358758464403, (6, 2) = .576262119133891, (6, 3) = -0.208529179808710e-1, (6, 4) = .379077667060160, (6, 5) = -.652551600592398}), Vector(6, {(1) = 27.4687324182218, (2) = 22.6431850097747, (3) = 8.55838822848258, (4) = 5.98572320151213, (5) = 2.01489965871576, (6) = 0.}), Matrix(5, 5, {(1, 1) = -.251382792720498, (1, 2) = -.396845551776930, (1, 3) = -.692151007470363, (1, 4) = -.366170444772230, (1, 5) = -.407635238653352, (2, 1) = .814836686086339, (2, 2) = .358661500188002, (2, 3) = -.248888011159285, (2, 4) = -.368593537944618, (2, 5) = -0.979625692668875e-1, (3, 1) = -.260618505584221, (3, 2) = .700768209407253, (3, 3) = -.220811446720437, (3, 4) = .385938483188542, (3, 5) = -.493250142851024, (4, 1) = .396723777130597, (4, 2) = -.450711241216643, (4, 3) = .251321149693754, (4, 4) = .434248601436671, (4, 5) = -.622684072035804, (5, 1) = -.218027763686546, (5, 2) = .140209949871121, (5, 3) = .589119449239943, (5, 4) = -.626528250364817, (5, 5) = -.439551692342332})

 (2)
> 

SDM:= DiagonalMatrix(S[1..5],5,5)

Matrix(5, 5, {(1, 1) = 27.46873241822184, (1, 2) = 0., (1, 3) = 0., (1, 4) = 0., (1, 5) = 0., (2, 1) = 0., (2, 2) = 22.643185009774694, (2, 3) = 0., (2, 4) = 0., (2, 5) = 0., (3, 1) = 0., (3, 2) = 0., (3, 3) = 8.558388228482578, (3, 4) = 0., (3, 5) = 0., (4, 1) = 0., (4, 2) = 0., (4, 3) = 0., (4, 4) = 5.985723201512132, (4, 5) = 0., (5, 1) = 0., (5, 2) = 0., (5, 3) = 0., (5, 4) = 0., (5, 5) = 2.014899658715757})

 (3)

THIS EQUALS TO ORIGINAL X MATRIX

 

> 

U.SDM.Vt

Matrix(6, 5, {(1, 1) = 8.789999999999997, (1, 2) = 9.93, (1, 3) = 9.829999999999995, (1, 4) = 5.449999999999993, (1, 5) = 3.159999999999998, (2, 1) = 6.1099999999999985, (2, 2) = 6.9099999999999975, (2, 3) = 5.0399999999999965, (2, 4) = -.26999999999999996, (2, 5) = 7.980000000000001, (3, 1) = -9.14999999999999, (3, 2) = -7.930000000000001, (3, 3) = 4.859999999999987, (3, 4) = 4.849999999999992, (3, 5) = 3.009999999999995, (4, 1) = 9.569999999999997, (4, 2) = 1.6399999999999977, (4, 3) = 8.82999999999999, (4, 4) = .7399999999999956, (4, 5) = 5.799999999999994, (5, 1) = -3.489999999999992, (5, 2) = 4.019999999999998, (5, 3) = 9.799999999999985, (5, 4) = 9.999999999999988, (5, 5) = 4.269999999999993, (6, 1) = 9.83999999999999, (6, 2) = .15000000000000033, (6, 3) = -8.989999999999982, (6, 4) = -6.0199999999999925, (6, 5) = -5.309999999999993})

 (4)
> 

X -~ U.SDM.Vt

Matrix(6, 5, {(1, 1) = 0.1776356839e-14, (1, 2) = 0., (1, 3) = 0.5329070518e-14, (1, 4) = 0.7105427358e-14, (1, 5) = 0.2220446049e-14, (2, 1) = 0.1776356839e-14, (2, 2) = 0.2664535259e-14, (2, 3) = 0.3552713679e-14, (2, 4) = -0.5551115123e-16, (2, 5) = -0.8881784197e-15, (3, 1) = -0.1065814104e-13, (3, 2) = 0.8881784197e-15, (3, 3) = 0.1332267630e-13, (3, 4) = 0.7993605777e-14, (3, 5) = 0.4884981308e-14, (4, 1) = 0.3552713679e-14, (4, 2) = 0.2220446049e-14, (4, 3) = 0.1065814104e-13, (4, 4) = 0.4440892099e-14, (4, 5) = 0.6217248938e-14, (5, 1) = -0.7993605777e-14, (5, 2) = 0.1776356839e-14, (5, 3) = 0.1598721155e-13, (5, 4) = 0.1243449788e-13, (5, 5) = 0.6217248938e-14, (6, 1) = 0.1065814104e-13, (6, 2) = -0.3330669074e-15, (6, 3) = -0.1776356839e-13, (6, 4) = -0.7105427358e-14, (6, 5) = -0.6217248938e-14})

 (5)

EIGENVALUES

 

> 

S*~S

Vector(6, {(1) = 754.5312606638714, (2) = 512.7138273868854, (3) = 73.24600906942915, (4) = 35.82888224512065, (5) = 4.059820634692874, (6) = 0.})

 (6)

EIGENVECTORS

 

> 

Transpose(Vt)

Matrix(5, 5, {(1, 1) = -.2513827927204978, (1, 2) = .8148366860863387, (1, 3) = -.2606185055842209, (1, 4) = .39672377713059703, (1, 5) = -.21802776368654583, (2, 1) = -.3968455517769299, (2, 2) = .35866150018800186, (2, 3) = .7007682094072526, (2, 4) = -.45071124121664313, (2, 5) = .1402099498711206, (3, 1) = -.6921510074703628, (3, 2) = -.2488880111592855, (3, 3) = -.22081144672043732, (3, 4) = .2513211496937536, (3, 5) = .5891194492399427, (4, 1) = -.3661704447722298, (4, 2) = -.3685935379446182, (4, 3) = .3859384831885419, (4, 4) = .434248601436671, (4, 5) = -.6265282503648171, (5, 1) = -.4076352386533523, (5, 2) = -0.979625692668875e-1, (5, 3) = -.4932501428510237, (5, 4) = -.6226840720358041, (5, 5) = -.4395516923423325})

 (7)

COMPARE LINE (6) AND THE evalues OF  LINE  (8), IN AGREEMENT.

 

COMPARE LINE (7) AND THE evectors OF LINE (8), NOT IN AGREEMENT - differ in sign

 

> 

evalues, evectors:= Eigenvectors(Transpose(X).X)

evalues, evectors := Vector(5, {(1) = 754.531260663872+0.*I, (2) = 512.713827386885+0.*I, (3) = 73.2460090694292+0.*I, (4) = 35.8288822451207+0.*I, (5) = 4.05982063469289+0.*I}), Matrix(5, 5, {(1, 1) = -.251382792720496+0.*I, (1, 2) = -.814836686086340+0.*I, (1, 3) = -.260618505584220+0.*I, (1, 4) = .396723777130597+0.*I, (1, 5) = .218027763686546+0.*I, (2, 1) = -.396845551776929+0.*I, (2, 2) = -.358661500188003+0.*I, (2, 3) = .700768209407252+0.*I, (2, 4) = -.450711241216643+0.*I, (2, 5) = -.140209949871121+0.*I, (3, 1) = -.692151007470364+0.*I, (3, 2) = .248888011159284+0.*I, (3, 3) = -.220811446720437+0.*I, (3, 4) = .251321149693753+0.*I, (3, 5) = -.589119449239943+0.*I, (4, 1) = -.366170444772231+0.*I, (4, 2) = .368593537944617+0.*I, (4, 3) = .385938483188543+0.*I, (4, 4) = .434248601436671+0.*I, (4, 5) = .626528250364817+0.*I, (5, 1) = -.407635238653353+0.*I, (5, 2) = 0.979625692668865e-1+0.*I, (5, 3) = -.493250142851024+0.*I, (5, 4) = -.622684072035804+0.*I, (5, 5) = .439551692342333+0.*I})

 (8)

Find signs of first rows of evectors and Transpose(Vt)

> 

signsev:=convert(map(signum,simplify(fnormal(evectors[1,..]),zero)),list);
signsV :=convert(map(signum,simplify(fnormal(Vt[..,1]),zero)),list);
diffsigns:=zip(`*`,signsev,signsV);

[-1, -1, -1, 1, 1]

  

[-1, 1, -1, 1, -1]

  

[1, -1, 1, 1, -1]

 (9)
> 

Transpose(Vt)-evectors.DiagonalMatrix(diffsigns);

Matrix(5, 5, {(1, 1) = -0.1998401444e-14+0.*I, (1, 2) = -0.8881784197e-15+0.*I, (1, 3) = -0.8326672685e-15+0.*I, (1, 4) = -0.1110223025e-15+0.*I, (1, 5) = -0.8326672685e-16+0.*I, (2, 1) = -0.9436895709e-15+0.*I, (2, 2) = -0.6661338148e-15+0.*I, (2, 3) = 0.5551115123e-15+0.*I, (2, 4) = 0.1665334537e-15+0.*I, (2, 5) = -0.3885780586e-15+0.*I, (3, 1) = 0.1110223025e-14+0.*I, (3, 2) = -0.1221245327e-14+0.*I, (3, 3) = -0.1387778781e-15+0.*I, (3, 4) = 0.5551115123e-15+0.*I, (3, 5) = -0.2220446049e-15+0.*I, (4, 1) = 0.1387778781e-14+0.*I, (4, 2) = -0.1054711873e-14+0.*I, (4, 3) = -0.8881784197e-15+0.*I, (4, 4) = -0.4440892099e-15+0.*I, (4, 5) = -0.2220446049e-15+0.*I, (5, 1) = 0.3330669074e-15+0.*I, (5, 2) = -0.9853229344e-15+0.*I, (5, 3) = 0.5551115123e-15+0.*I, (5, 4) = -0.3330669074e-15+0.*I, (5, 5) = 0.7216449660e-15+0.*I})

 (10)
> 

sdm:= DiagonalMatrix(sqrt~(evalues[1..5]),5,5)

Matrix(5, 5, {(1, 1) = 27.46873241822184+0.*I, (1, 2) = 0.*I, (1, 3) = 0.*I, (1, 4) = 0.*I, (1, 5) = 0.*I, (2, 1) = 0.*I, (2, 2) = 22.64318500977469+0.*I, (2, 3) = 0.*I, (2, 4) = 0.*I, (2, 5) = 0.*I, (3, 1) = 0.*I, (3, 2) = 0.*I, (3, 3) = 8.558388228482581+0.*I, (3, 4) = 0.*I, (3, 5) = 0.*I, (4, 1) = 0.*I, (4, 2) = 0.*I, (4, 3) = 0.*I, (4, 4) = 5.985723201512133+0.*I, (4, 5) = 0.*I, (5, 1) = 0.*I, (5, 2) = 0.*I, (5, 3) = 0.*I, (5, 4) = 0.*I, (5, 5) = 2.0148996587157613+0.*I})

 (11)

THIS SHOULD EQUAL TO THE ORIGINAL X MATRIX

 

> 

U.sdm.Transpose(evectors.DiagonalMatrix(diffsigns));

Matrix(6, 5, {(1, 1) = 8.789999999999969+0.*I, (1, 2) = 9.929999999999986+0.*I, (1, 3) = 9.830000000000018+0.*I, (1, 4) = 5.450000000000025+0.*I, (1, 5) = 3.1600000000000064+0.*I, (2, 1) = 6.10999999999998+0.*I, (2, 2) = 6.909999999999992+0.*I, (2, 3) = 5.040000000000017+0.*I, (2, 4) = -.2699999999999823+0.*I, (2, 5) = 7.980000000000011+0.*I, (3, 1) = -9.150000000000004+0.*I, (3, 2) = -7.930000000000009+0.*I, (3, 3) = 4.859999999999969+0.*I, (3, 4) = 4.849999999999975+0.*I, (3, 5) = 3.0099999999999865+0.*I, (4, 1) = 9.569999999999972+0.*I, (4, 2) = 1.6399999999999917+0.*I, (4, 3) = 8.830000000000009+0.*I, (4, 4) = .7400000000000133+0.*I, (4, 5) = 5.800000000000008+0.*I, (5, 1) = -3.490000000000021+0.*I, (5, 2) = 4.019999999999979+0.*I, (5, 3) = 9.799999999999986+0.*I, (5, 4) = 10.000000000000002+0.*I, (5, 5) = 4.26999999999999+0.*I, (6, 1) = 9.840000000000016+0.*I, (6, 2) = .1500000000000148+0.*I, (6, 3) = -8.98999999999998+0.*I, (6, 4) = -6.019999999999986+0.*I, (6, 5) = -5.309999999999981+0.*I})

 (12)
> 

X - U.sdm.Transpose(evectors.DiagonalMatrix(diffsigns));;

Matrix(6, 5, {(1, 1) = 0.3019806627e-13+0.*I, (1, 2) = 0.1421085472e-13+0.*I, (1, 3) = -0.1776356839e-13+0.*I, (1, 4) = -0.2486899575e-13+0.*I, (1, 5) = -0.6217248938e-14+0.*I, (2, 1) = 0.2042810365e-13+0.*I, (2, 2) = 0.7993605777e-14+0.*I, (2, 3) = -0.1687538997e-13+0.*I, (2, 4) = -0.1770805724e-13+0.*I, (2, 5) = -0.1065814104e-13+0.*I, (3, 1) = 0.3552713679e-14+0.*I, (3, 2) = 0.8881784197e-14+0.*I, (3, 3) = 0.3108624469e-13+0.*I, (3, 4) = 0.2486899575e-13+0.*I, (3, 5) = 0.1332267630e-13+0.*I, (4, 1) = 0.2842170943e-13+0.*I, (4, 2) = 0.8215650382e-14+0.*I, (4, 3) = -0.8881784197e-14+0.*I, (4, 4) = -0.1332267630e-13+0.*I, (4, 5) = -0.7993605777e-14+0.*I, (5, 1) = 0.2087219286e-13+0.*I, (5, 2) = 0.2042810365e-13+0.*I, (5, 3) = 0.1421085472e-13+0.*I, (5, 4) = -0.1776356839e-14+0.*I, (5, 5) = 0.9769962617e-14+0.*I, (6, 1) = -0.1598721155e-13+0.*I, (6, 2) = -0.1479372180e-13+0.*I, (6, 3) = -0.1953992523e-13+0.*I, (6, 4) = -0.1332267630e-13+0.*I, (6, 5) = -0.1865174681e-13+0.*I})

 (13)
> 

 

Download SVD_test2.mw

 

some examples...

dharr 8917
March 08 2024
1 1

Not sure about documentation or the level of explanation you want, but if examples help, here are some:

> 

restart

infix -> prefix

> 

`and`(a, b, c, d)

a and b and c and d

> 

`+`(a, b, c, d)

a+b+c+d

> 

`=`(a, b)

a = b

> 

`union`({a}, {b, c}, {a, d})

{a, b, c, d}

> 

`$`(2, 3)

2, 2, 2

> 

`$`(2, 3)

2, 2, 2

postfix -> prefix

> 

L := [a, b, c, d]

[a, b, c, d]

> 

L[]

a, b, c, d

> 

`?[]`(L)

a, b, c, d

> 

factorial(3)

6

> 

`!`(3); factorial(3)

`!`(3)

6

more...

> 

`and`(L[])

a and b and c and d

> 

(`@`(`and`, `?[]`))(L)

a and b and c and d

> 

`<,>`(a, b)

Vector[column](%id = 36893490135843830708)

> 

`<,>`(a, b)

Vector[column](%id = 36893490135843825652)

> 

`<|>`(a, b)

Vector[row](%id = 36893490135843820236)

> 

`<|>`(a, b)

Vector[row](%id = 36893490135843822756)

> 

NULL

Download prefixes.mw

spacecurve...

dharr 8917
March 07 2024
1 2

A workaround is to use spacecurve instead of plot3d. Unfortunately, the curves need to be done separately.

Perpendicular_3D_lines.mw

3D PDF from Maple 3D plots...

dharr 8917
March 06 2024
2 2

Interesting; I didn't know about such a thing. I looked into it and played around a bit, and various of the exported plot formats can be used to make a u3d file, which can then be imported into Acrobat to make a 3D PDF or used in pdflatex to make a 3D PDF.

Here are the steps I tried

1. Make a 3d plot of a red cube, and export it to a file

p:=plots:-display(plottools:-cuboid([0,0,0],[1,1,1],color=red));
Export("redcube.off", p, base=currentdir);

I tried .off, .ply, .stl via the Export command and .x3d geometry and COLLADA (.dae) via the plot context menu. (These formats are all known to MeshLab.)

2. Convert to .u3d

I used MeshLab from www.meshlab.net (Windows 64 bit but other OSs are available). File -> Import mesh -> choose file. .ply gave an error, .stl asked about removing duplicate vertices, which I accepted. Only .off showed as a red cube; the others had no color.

3a. Open Acrobat Pro and create a blank page .pdf (using file->create). Find the rich media tool and open it. Select add 3d and drag a rectangle where you want the object to be. Choose the .u3d file. You will have to overide a media safety warning, then you can click on the image to activate it. Then click and you can move around. There was no color, and some faces of the cube were transparent. (A sphere from complexplot3d also lost color but was faithfully reproduced.)

3b. The .u3d file can instead be incorporated into .pdf via the LaTeX package media9 and pdfLaTeX. I used overleaf (www.overleaf.com), which uses TeX Live. The .zip file here contains the .tex file, the .u3d file and the generated .pdf file.

3D_PDF.zip

Comments - only the 3D objects are exported, e.g. from POLYGONS PLOT3D structures. These can be colored face by face, but by default the cuboid command sets them all at once to red. I tried coloring different faces differently, but Maple was not exporting the colors to the .obj file (and maybe others).

no value...

dharr 8917
March 06 2024
1 0

The default value is shown when there is no argument passed corresponding to that parameter. (Using NULL achieves the same thing because the NULL is removed from the sequence.)

> 

restart
> 

test := proc (a := 1, b := 2) a+b end proc

proc (a := 1, b := 2) a+b end proc

> 

test(3, 1)

4

> 

test(3)

5

> 

test(3, NULL)

5

> 

3, NULL

3

> 

NULL

Download test.mw

arrows palette...

dharr 8917
March 05 2024
0 0

The arrows palette has diagonal arrows, which can be inserted in 2D, e.g. in a Table:

x

-5

 

-2

 

0

 

2

 

3

 

4

diff(f(x), x)

|

-

0

+

0

-

0

-

0

+

|

f(x)

|

a

lok. min

b

lok. max

a

van v.tg

a

lok. min

b

|
 

Download arrows.mw

 

One way...

dharr 8917
March 05 2024
0 0

Here's one way.

> 

restart
> 

_local(D, gamma)
> 

eqP1 := P1 = (D-C)/(1+gamma); eqP2 := P2 = (D*gamma+C)/(1+gamma)
> 

Theta_exact := -Y*(-D+C)/(1+gamma)+(D*gamma+C)/(1+gamma)

-Y*(-D+C)/(1+gamma)+(D*gamma+C)/(1+gamma)

> 

eliminate({eqP1, eqP2, Theta_exact}, {C, D}); Theta_exact2 := %[2][]

[{C = -P1*gamma+P2, D = P1+P2}, {P1*Y+P2}]

P1*Y+P2

NULL

Download substitution_problem.mw

Mathcontainer...

dharr 8917
March 03 2024
2 1

The font size in a TextArea component cannot be changed as far as I can see. But if you change to a MathContainer component, you can use all the typesetting capabilities. Here is one way to do it:

S4ThalèsAnimation.mw

Zeckendorf's theorem...

dharr 8917
March 03 2024
2 2

Zeckendorf's theorem in the title of your post is about a unique sum of Fibonacci's, not all sums. @mmcdara has just given a version of all sums; here is an efficient implementation of the unique version.

Zeckendorf's theorem. There is a unique sequence of distinct fibonacci numbers adding to F if we exclude consecutive fibonacci numbers.

> 

restart

Golden ratio

> 

phi := (1+sqrt(5))*(1/2)

1/2+(1/2)*5^(1/2)

n for largest fibonacci not greater than F - see wikipedia

> 

nlarge := unapply(floor(log[phi](sqrt(5)*(F+1/2))), F)

proc (F) options operator, arrow; floor(ln(5^(1/2)*(F+1/2))/ln(1/2+(1/2)*5^(1/2))) end proc

Largest fibonacci not greater than 100.

> 

combinat:-fibonacci(nlarge(100))

89

Greedy algorithm - successively find largest fibonacci's f not greater than the remainder r.
Example: Adding to 100

> 

"r:=100: do f:=combinat:-fibonacci(nlarge(r)); print(f); r:=r-f; until r=0:"

89

8

3

> 

zeck := proc(n::posint)
  local nlarge, r, f, i;
  nlarge := F -> floor(log[(1 + sqrt(5))/2](sqrt(5)*(F + 1/2)));
  r := n;
  f := table();
  for i do
    f[i] := combinat:-fibonacci(nlarge(r));
    r := r - f[i];
  until r = 0;
  convert(f, list);
end proc:
> 

zeck(100)

[89, 8, 3]

> 

zeck(10^12); add(%)

[956722026041, 32951280099, 7778742049, 1836311903, 701408733, 9227465, 832040, 121393, 46368, 2584, 987, 233, 89, 13, 3]

1000000000000

> 

NULL

Download Zeckendorf.mw

zero outside range...

dharr 8917
February 29 2024
1 6

Using trace and showstat shows that this doesn't go through GAMMA, despite the help page. Hope I've traced this correctly; I haven't given all the details. (For positive integers, binomial(a,b) = 0 for b>a is consistent with no ways to choose b out of a.)

> 

restart;
> 

trace(binomial);

[binomial:-ModuleApply]

Consider first integers outside the range 0<=k<=n; the result is just zero

> 

binomial(4, 6);

execute binomial:-ModuleApply, args = 4, 6

{--> enter binomial:-ModuleApply, args = 4, 6

6

-2

4

0

<-- exit binomial:-ModuleApply (now at top level) = 0}

0

The code is just:
elif type(A,'integer') and type(B,'integer') then
            if 0 < A then
                if B < 0 or A < B then
                    res := 0 (result)

Now consider

> 

binomial(n, n + 2) assuming n::posint;

execute binomial:-ModuleApply, args = \`n~\`, n+2

{--> enter binomial:-ModuleApply, args = \`n~\`, n+2

n+2

-2

n

1

1

0

<-- exit binomial:-ModuleApply (now in \`assuming\`) = 0}

0

So what happens here:

B:=normal(b) -> n+2;
Work out AmB: = a-B -> -2
A:=normal(a) -> n;
s: = signum(A)  -> 1

If s = 1 or -1 and AmB::negint (yes)

then

i := signum(0,B,0) -> 1 (sign of B treating 0 as 0)
if i is 1 or -1 then
               if s = 1 or s = i then
                      res := 0

So result is zero.

 

In other words, since A and B are positive and A-B is a negative integer the result is zero

 

binomial(n, n+1) uses  `tool/type` which makes it a bit more obscure but seems to be a similar argument.

> 

 

Download binomial.mw

ad hoc method for your second case...

dharr 8917
February 27 2024
1 13

This is not as systematic as @mmcdara's method, but answers the question for your second example containing sigma__d3. You could use the approximate form for Lambda instead.

> 

restart
> 

local gamma:with(plots):

 # Define the assumptions - or Maple wouldn't know

> 

assume(0 < gamma, 0 < sigma__d, 0 < sigma__d3, 0 < sigma__v);
interface(showassumed=0);

1

# Input lambda parameters

> 

Lambda := RootOf(8*_Z^4 + 12*Gamma*_Z^3 + (5*Gamma^2 - 4)*_Z^2 - 4*Gamma*_Z - Gamma^2);
f__1 := Lambda*(sigma__v/sigma__d):
f__2 := f__1:
f__3 := sqrt(2)*sigma__v/sigma__d3:

RootOf(8*_Z^4+12*Gamma*_Z^3+(5*Gamma^2-4)*_Z^2-4*Gamma*_Z-Gamma^2)

 Let's scale these. If we divide f__1 by sqrt(2)*(sigma__v/sigma__d) we get Lambda/sqrt(2); we divide f__2 and f__3 by the same.

> 

f__1sc := f__1/(sqrt(2)*(sigma__v/sigma__d));
f__2sc := f__2/(sqrt(2)*(sigma__v/sigma__d));
f__3sc := f__3/(sqrt(2)*(sigma__v/sigma__d));

 

(1/2)*RootOf(8*_Z^4+12*Gamma*_Z^3+(5*Gamma^2-4)*_Z^2-4*Gamma*_Z-Gamma^2)*2^(1/2)

(1/2)*RootOf(8*_Z^4+12*Gamma*_Z^3+(5*Gamma^2-4)*_Z^2-4*Gamma*_Z-Gamma^2)*2^(1/2)

sigma__d/sigma__d3

Now Gamma = gamma*sigma__v*sigma__d = gamma*sigma__v*sigma__d3*f3sc = a*f3sc, where a = gamma*sigma__v*sigma__d3.

We want to know when f3sc = f1sc, or Gamma/a = Lambda/sqrt(2), which will only depend on the combined parameter a = gamma* sigma__v*sigma__d3

> 

Gamma/a=f__1sc;
aval:=solve(%,a);

Gamma/a = (1/2)*RootOf(8*_Z^4+12*Gamma*_Z^3+(5*Gamma^2-4)*_Z^2-4*Gamma*_Z-Gamma^2)*2^(1/2)

Gamma*2^(1/2)/RootOf(8*_Z^4+12*Gamma*_Z^3+(5*Gamma^2-4)*_Z^2-4*Gamma*_Z-Gamma^2)

Plot f__1sc = f__3sc = sigma__d/sigma__d3 vs a = gamma*sigma__v*sigma__d3

> 

pp:=plot([aval,f__1sc,Gamma=0..2],labels=[gamma*sigma__v*sigma__d3,sigma__d/sigma__d3]):
> 

display(pp,textplot([3,0.3,typeset('f__3' < 'f__1')]),textplot([3,0.45,typeset('f__3' > 'f__1')]));

Test for some parameters

> 

params:={gamma = 1, sigma__v = 1, sigma__d3 =3, sigma__d = 1};
evalf(eval([gamma*sigma__v*sigma__d3, sigma__d/sigma__d3], params));

{gamma = 1, sigma__d3 = 3, sigma__d = 1, sigma__v = 1}

[3., .3333333333]

This point is below the line, so we expect f3 < f1

> 

evalf(eval(f__1sc, Gamma=eval(gamma*sigma__v*sigma__d,params)));
evalf(eval(f__3sc, params));

.3975956314

.3333333333

> 

NULL

Download complicated_comparison2.mw

results can be different...

dharr 8917
February 26 2024
0 0

The documentation for IsFrobeniusGroup gives an example where IsFrobeniusGroup is true but IsFrobeniusPermGroup is false, so I think this result is OK. The documentation you referred to (for FrobeniusGroup) refers to the situation where  IsFrobeniusPermGroup is true, and so IsFrobeniusGroup will be true.

different perm sets...

dharr 8917
February 25 2024
1 2

The small group (8,3) is expressed as a permutation group on 1,2,..8, but DihedralGroup(4) is expressed as a permutation group on 1,2,3,4. So the domain [1,2,3,4] is only half the points for (8,3) but it is all the points for D4. ReducedDegreePermGroup will convert (8,3) to the form you need.

> 

restart
> 

with(GroupTheory)
> 

G1 := DihedralGroup(4)

_m2189854788256

> 

G2 := SmallGroup(IdentifySmallGroup(G1))

_m2189934002208

> 

AreIsomorphic(G1, G2)

true

> 

IsTransitive(G1); IsTransitive(G2)

true

true

> 

Elements(G1)

{_m2189853538880, _m2189854757216, _m2189854758688, _m2189854759104, _m2189854759616, _m2189854760000, _m2189854787040, _m2189854788096}

> 

Elements(G2)

{_m2189840564448, _m2189840564832, _m2189840565216, _m2189840565600, _m2189854757216, _m2189934009088, _m2189934009312, _m2189934009504}

> 

MinimumPermutationRepresentationDegree(G2)

4

> 

G3 := ReducedDegreePermGroup(G2)

_m2189866943712

> 

AreIsomorphic(G1, G3)

true

> 

IsTransitive(G3, [1, 2, 3, 4])

true

> 

NULL

Download Groups.mw

congruent numbers...

dharr 8917
February 24 2024
1 1

You have missed the fact that the algorithm is only for squarefree integers. For example 20 should give true, but to decide this you need to divide by the largest square factor 4 and apply the algorithm to the quotient 5. You also missed 6, which is squarefree, for some reason that wasn't clear to me. 

(I used @Ronan's version replacing print("True") with return true etc, but have my own working version if you are interested.)

First 24 25 26 27 28 29 30 Last Page 26 of 86