{VERSION 5 0 "IBM INTEL NT" "5.0" } {USTYLETAB {CSTYLE "Maple Input" -1 0 "Courier" 0 1 255 0 0 1 0 1 0 0 1 0 0 0 0 1 }{CSTYLE "2D Math" -1 2 "Times" 0 1 0 0 0 0 0 0 2 0 0 0 0 0 0 1 }{CSTYLE "2D Comment" 2 18 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 } {CSTYLE "" -1 256 "" 1 14 0 0 0 0 0 0 0 0 0 0 0 0 0 1 }{PSTYLE "Normal " -1 0 1 {CSTYLE "" -1 -1 "Times" 1 12 0 0 0 1 2 2 2 2 2 2 1 1 1 1 }1 1 0 0 0 0 1 0 1 0 2 2 0 1 }} {SECT 0 {EXCHG {PARA 0 "" 0 "" {TEXT 256 66 " \+ Lame Polynomials " }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 118 " A Maple worksheet written by Hans Volkmer, F ebruary 10, 2001, updated February 13, 2002, October 15, 2002." }}} {EXCHG {PARA 0 "" 0 "" {TEXT -1 83 " \+ send comments to volkmer@uwm.edu" }}}{EXCHG {PARA 0 " " 0 "" {TEXT -1 63 "In this worksheet we compute the eight types of La me polynomals" }}{PARA 0 "" 0 "" {TEXT -1 17 " " } {XPPEDIT 18 0 "uE[2*n]^m;" "6#)&%#uEG6#*&\"\"#\"\"\"%\"nGF)%\"mG" } {TEXT -1 1 "(" }{XPPEDIT 18 0 "z,k^2;" "6$%\"zG*$%\"kG\"\"#" }{TEXT -1 3 "), " }{XPPEDIT 18 0 "sE[2*n+1]^m;" "6#)&%#sEG6#,&*&\"\"#\"\"\"% \"nGF*F*F*F*%\"mG" }{TEXT -1 1 "(" }{XPPEDIT 18 0 "z,k^2;" "6$%\"zG*$% \"kG\"\"#" }{TEXT -1 2 ")," }{XPPEDIT 18 0 "cE[2*n+1]^m;" "6#)&%#cEG6# ,&*&\"\"#\"\"\"%\"nGF*F*F*F*%\"mG" }{TEXT -1 1 "(" }{XPPEDIT 18 0 "z,k ^2;" "6$%\"zG*$%\"kG\"\"#" }{TEXT -1 2 ")," }{XPPEDIT 18 0 "dE[2*n+1]^ m;" "6#)&%#dEG6#,&*&\"\"#\"\"\"%\"nGF*F*F*F*%\"mG" }{TEXT -1 2 " (" } {XPPEDIT 18 0 "z,k^2;" "6$%\"zG*$%\"kG\"\"#" }{TEXT -1 2 ")," }}{PARA 0 "" 0 "" {XPPEDIT 18 0 "scE[2*n+2]^m;" "6#)&%$scEG6#,&*&\"\"#\"\"\"% \"nGF*F*F)F*%\"mG" }{TEXT -1 1 "(" }{XPPEDIT 18 0 "z,k^2;" "6$%\"zG*$% \"kG\"\"#" }{TEXT -1 3 "), " }{XPPEDIT 18 0 "sdE[2*n+2]^m;" "6#)&%$sdE G6#,&*&\"\"#\"\"\"%\"nGF*F*F)F*%\"mG" }{TEXT -1 1 "(" }{XPPEDIT 18 0 " z,k^2;" "6$%\"zG*$%\"kG\"\"#" }{TEXT -1 3 "), " }{XPPEDIT 18 0 "cdE[2* n+2]^m;" "6#)&%$cdEG6#,&*&\"\"#\"\"\"%\"nGF*F*F)F*%\"mG" }{TEXT -1 1 " (" }{XPPEDIT 18 0 "z,k^2;" "6$%\"zG*$%\"kG\"\"#" }{TEXT -1 3 "), " } {XPPEDIT 18 0 "scdE[2*n+3]^m;" "6#)&%%scdEG6#,&*&\"\"#\"\"\"%\"nGF*F* \"\"$F*%\"mG" }{TEXT -1 2 " (" }{XPPEDIT 18 0 "z,k^2;" "6$%\"zG*$%\"kG \"\"#" }{TEXT -1 1 ")" }}{PARA 0 "" 0 "" {TEXT -1 19 "for given m, n a nd " }{XPPEDIT 18 0 "k^2;" "6#*$%\"kG\"\"#" }{TEXT -1 47 " as polynomi als in Jacobian elliptic functions " }{XPPEDIT 18 0 "sn*z;" "6#*&%#snG \"\"\"%\"zGF%" }{TEXT -1 2 ", " }{XPPEDIT 18 0 "cn*z;" "6#*&%#cnG\"\" \"%\"zGF%" }{TEXT -1 2 ", " }{XPPEDIT 18 0 "dn*z;" "6#*&%#dnG\"\"\"%\" zGF%" }{TEXT -1 35 ". Note that n must be a nonnegative" }}{PARA 0 "" 0 "" {TEXT -1 29 "integer and m = 0,1,2,...,n." }}{PARA 0 "" 0 "" {TEXT -1 54 "The Lame polynomials are normalized as in Chapter LA.." } }}{EXCHG {PARA 0 "" 0 "" {TEXT -1 89 "Warning: Note that Jacobian elli ptic functions are incorrectly implemented in Maple 6.00." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 37 "restart;with(linalg):with(orthopoly ):" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 38 "Digits:=10; # select \+ required precison" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 40 "matrix u:=proc(n,k2) # matrix for type uE" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 13 "local A,r,nu;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 8 "nu:=2*n;" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "A:=matrix(n+1,n+1,0);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 49 "for r from 0 to n do A[r+1,r+1]:=4*r^2*(2-k2) od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 40 "if n>0 then A[1,2]:=(nu-1)* (nu+2)*k2 fi;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 68 "for r from 1 to n- 1 do A[r+1,r+2]:=0.5*(nu-2*r-1)*(nu+2*r+2)*k2 od; " }}{PARA 0 "> " 0 " " {MPLTEXT 1 0 63 "for r from 1 to n do A[r+1,r]:=0.5*(nu-2*r+2)*(nu+2 *r-1)*k2 od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 2 "A;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 4 "end;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 62 "eigenu:=proc(n,k2,q) # q-th normalized eigenvector for matrixu" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 18 "local c,s,i,t,A,K;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 17 "A:=matrixu(n,k2);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "K:=[eigenvectors(A)];" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 45 "K: =sort(K,(a,b)->evalb(Re(a[1]) < Re(b[1])));" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 15 "c:=K[q][3][1]; " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 14 "s:=0.5*c[1]^2;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 38 "for i from 2 to \+ n+1 do s:=s+c[i]^2 od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 14 "s:=1/sqrt (s); " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "t:=0.5*c[1];" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 36 "for i from 2 to n+1 do t:=t+c[i] od;" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 22 "if t<0 then s:=-s fi; " }}{PARA 0 " > " 0 "" {MPLTEXT 1 0 39 "for i from 1 to n+1 do c[i]:=s*c[i] od;" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 2 "c;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 4 "end;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 45 "uE:=proc(m,n,k 2) # gives uE for given m,n,k^2" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 12 " local c,p,E;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 20 "c:=eigenu(n,k2,m+1) ;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "E:=1/2*c[1];" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 65 "for p from 1 to n do E:=E+c[p+1]*T(2*p,JacobiSN( z,sqrt(k2))); od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 2 "E;" }}{PARA 0 " > " 0 "" {MPLTEXT 1 0 4 "end;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 41 "matrixs:=proc(n,k2) # matrix for type sE " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 13 "local A,r,nu;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 10 "n u:=2*n+1;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "A:=matrix(n+1,n+1,0); " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 30 "A[1,1]:=2-k2+0.5*nu*(nu+1)*k2; " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 53 "for r from 1 to n do A[r+1,r+1] :=(2*r+1)^2*(2-k2) od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 68 "for r fro m 0 to n-1 do A[r+1,r+2]:=0.5*(nu-2*r-2)*(nu+2*r+3)*k2 od; " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 61 "for r from 1 to n do A[r+1,r]:=0.5*(nu-2* r+1)*(nu+2*r)*k2 od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 2 "A;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 4 "end;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "eigens:=proc(n,k2,q) " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 18 "local c,s,i,t,A,K;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 19 "A:=ma trixs(n,k2,d);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "K:=[eigenvectors( A)];" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 45 "K:=sort(K,(a,b)->evalb(Re(a [1]) < Re(b[1])));" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 14 "c:=K[q][3][1] ;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 7 "s:=0.0;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 38 "for i from 1 to n+1 do s:=s+c[i]^2 od;" }}{PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 13 "s:=1/sqrt(s);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 7 "t:=0.0;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 36 "for i from 1 to n+1 do t:=t+c[i] od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 22 "if t<0 then s: =-s fi; " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 39 "for i from 1 to n+1 do \+ c[i]:=s*c[i] od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 2 "c;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 4 "end;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 16 "sE:=proc(m,n,k2)" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "local c,p ,E;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 20 "c:=eigens(n,k2,m+1);" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 5 "E:=0;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 66 "for p from 0 to n do E:=E+c[p+1]*T(2*p+1,JacobiSN(z,sqrt(k2))) od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 2 "E;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 4 "end;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 42 "mat rixc:=proc(n,k2) # matrix for type cE " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 13 "local A,r,nu;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 10 "nu:=2*n+1; " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "A:=matrix(n+1,n+1,0);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 30 "A[1,1]:=2-k2-0.5*nu*(nu+1)*k2;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 53 "for r from 1 to n do A[r+1,r+1]:=(2*r+1)^ 2*(2-k2) od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 68 "for r from 0 to n-1 do A[r+1,r+2]:=0.5*(nu-2*r-2)*(nu+2*r+3)*k2 od; " }}{PARA 0 "> " 0 " " {MPLTEXT 1 0 61 "for r from 1 to n do A[r+1,r]:=0.5*(nu-2*r+1)*(nu+2 *r)*k2 od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 2 "A;" }}{PARA 0 "> " 0 " " {MPLTEXT 1 0 4 "end;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "e igenc:=proc(n,k2,q) " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 18 "local c,s,i ,t,A,K;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "A:=matrixc(n,k2,n+1);" } }{PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "K:=[eigenvectors(A)];" }}{PARA 0 " > " 0 "" {MPLTEXT 1 0 45 "K:=sort(K,(a,b)->evalb(Re(a[1]) < Re(b[1]))) ;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 14 "c:=K[q][3][1];" }}{PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 7 "s:=0.0;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 38 "f or i from 1 to n+1 do s:=s+c[i]^2 od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 13 "s:=1/sqrt(s);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 7 "t:=0.0;" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 44 "for i from 1 to n+1 do t:=t+(2*i-1) *c[i] od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 22 "if t<0 then s:=-s fi; \+ " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 39 "for i from 1 to n+1 do c[i]:=s* c[i] od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 2 "c;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 4 "end;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 16 "cE: =proc(m,n,k2)" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "local c,p,E;" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 20 "c:=eigenc(n,k2,m+1);" }}{PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 5 "E:=0;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 64 "for p from 0 to n do E:=E+c[p+1]*U(2*p,JacobiSN(z,sqrt(k2))) od;" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 23 "JacobiCN(z,sqrt(k2))*E;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 4 "end;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 40 "matrixd:=proc(n,k2) # matrix for type dE" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 13 "local A,r,nu;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 10 "n u:=2*n+1;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "A:=matrix(n+1,n+1,0); " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 49 "for r from 0 to n do A[r+1,r+1] :=4*r^2*(2-k2) od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 36 "if n>0 then A [1,2]:=nu*(nu+1)*k2 fi;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 66 "for r fr om 1 to n-1 do A[r+1,r+2]:=0.5*(nu-2*r)*(nu+2*r+1)*k2 od; " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 61 "for r from 1 to n do A[r+1,r]:=0.5*(nu-2*r+ 1)*(nu+2*r)*k2 od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 2 "A;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 4 "end;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "eigend:=proc(n,k2,q) " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 24 "lo cal c,s1,s2,s,i,t,A,K;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "A:=matrix d(n,k2,n+1);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "K:=[eigenvectors(A) ];" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 45 "K:=sort(K,(a,b)->evalb(Re(a[1 ]) < Re(b[1])));" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 14 "c:=K[q][3][1]; " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 15 "s1:=0.5*c[1]^2;" }}{PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 40 "for i from 2 to n+1 do s1:=s1+c[i]^2 od;" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 7 "s2:=0; " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 43 "for i from 1 to n do s2:=s2+c[i]*c[i+1] od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 40 "s:=(1-0.5*k2)*s1-0.5*k2*s2;s:=1/sqrt(s); " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "t:=0.5*c[1];" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 36 "for i from 2 to n+1 do t:=t+c[i] od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 22 "if t<0 then s:=-s fi; " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 39 "for i from 1 to n+1 do c[i]:=s*c[i] od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 2 "c;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 4 "end;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 16 "dE:=proc(m,n,k2)" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "local c,p,E;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 20 "c:=eigend(n,k2,m+1);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "E:=1/2*c[1];" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 64 "for p from 1 \+ to n do E:=E+c[p+1]*T(2*p,JacobiSN(z,sqrt(k2))) od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 26 "E:=JacobiDN(z,sqrt(k2))*E;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 4 "end;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 44 "mat rixsc:=proc(n,k2) # matrix for type scE " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 13 "local A,r,nu;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 10 "n u:=2*n+2;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "A:=matrix(n+1,n+1,0); " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 53 "for r from 0 to n do A[r+1,r+1] :=(2*r+2)^2*(2-k2) od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 68 "for r fro m 0 to n-1 do A[r+1,r+2]:=0.5*(nu-2*r-3)*(nu+2*r+4)*k2 od; " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 61 "for r from 1 to n do A[r+1,r]:=0.5*(nu-2* r)*(nu+2*r+1)*k2 od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 2 "A;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 4 "end;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 22 "eigensc:=proc(n,k2,q) " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 18 "local c,s,i,t,A,K;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 22 "A:=ma trixsc(n,k2,n+1);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "K:=[eigenvecto rs(A)];" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 45 "K:=sort(K,(a,b)->evalb(R e(a[1]) < Re(b[1])));" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 14 "c:=K[q][3] [1];" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 7 "s:=0.0;" }}{PARA 0 "> " 0 " " {MPLTEXT 1 0 38 "for i from 1 to n+1 do s:=s+c[i]^2 od;" }}{PARA 0 " > " 0 "" {MPLTEXT 1 0 13 "s:=1/sqrt(s);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 7 "t:=0.0;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 40 "for i from 1 to n +1 do t:=t+2*i*c[i] od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 22 "if t<0 t hen s:=-s fi; " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 39 "for i from 1 to n +1 do c[i]:=s*c[i] od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 2 "c;" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 4 "end;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 17 "scE:=proc(m,n,k2)" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "local c,p,E;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "c:=eigensc(n,k2 ,m+1);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 5 "E:=0;" }}{PARA 0 "> " 0 " " {MPLTEXT 1 0 66 "for p from 0 to n do E:=E+c[p+1]*U(2*p+1,JacobiSN(z ,sqrt(k2))) od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 23 "JacobiCN(z,sqrt( k2))*E;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 4 "end;" }}}{EXCHG {PARA 0 " > " 0 "" {MPLTEXT 1 0 44 "matrixsd:=proc(n,k2) # matrix for type sdE \+ " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 13 "local A,r,nu;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 10 "nu:=2*n+2;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "A:=matrix(n+1,n+1,0);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 30 "A[1,1]:=2 -k2+0.5*nu*(nu+1)*k2;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 53 "for r from 1 to n do A[r+1,r+1]:=(2*r+1)^2*(2-k2) od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 68 "for r from 0 to n-1 do A[r+1,r+2]:=0.5*(nu-2*r-1)*(nu +2*r+2)*k2 od; " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 61 "for r from 1 to \+ n do A[r+1,r]:=0.5*(nu-2*r)*(nu+2*r+1)*k2 od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 2 "A;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 4 "end;" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 22 "eigensd:=proc(n,k2,q) " }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 24 "local c,s1,s2,s,i,t,A,K;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 18 "A:=matrixsd(n,k2);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "K:=[eigenvectors(A)];" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 45 "K:=sort(K,(a,b)->evalb(Re(a[1]) < Re(b[1])));" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 14 "c:=K[q][3][1];" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 8 "s1:=0.0;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 40 "for i from 1 to n+ 1 do s1:=s1+c[i]^2 od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 16 "s2:=0.5*c [1]^2; " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 43 "for i from 1 to n do s2: =s2+c[i]*c[i+1] od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 40 "s:=(1-0.5*k2 )*s1-0.5*k2*s2;s:=1/sqrt(s);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 7 "t:=0 .0;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 36 "for i from 1 to n+1 do t:=t+ c[i] od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 22 "if t<0 then s:=-s fi; \+ " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 39 "for i from 1 to n+1 do c[i]:=s* c[i] od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 2 "c;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 4 "end;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 17 "sdE :=proc(m,n,k2)" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "local c,p,E;" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "c:=eigensd(n,k2,m+1);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 5 "E:=0;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 66 "fo r p from 0 to n do E:=E+c[p+1]*T(2*p+1,JacobiSN(z,sqrt(k2))) od;" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 26 "E:=JacobiDN(z,sqrt(k2))*E;" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 4 "end;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 44 "matrixcd:=proc(n,k2) # matrix for type cdE " }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 15 "local d,A,r,nu;" }}{PARA 0 "> " 0 " " {MPLTEXT 1 0 10 "nu:=2*n+2;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "A: =matrix(n+1,n+1,0);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 30 "A[1,1]:=2-k2 -0.5*nu*(nu+1)*k2;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 53 "for r from 1 \+ to n do A[r+1,r+1]:=(2*r+1)^2*(2-k2) od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 68 "for r from 0 to n-1 do A[r+1,r+2]:=0.5*(nu-2*r-1)*(nu +2*r+2)*k2 od; " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 61 "for r from 1 to \+ n do A[r+1,r]:=0.5*(nu-2*r)*(nu+2*r+1)*k2 od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 2 "A;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 4 "end;" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 22 "eigencd:=proc(nu,k2,q)" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 24 "local c,s1,s2,s,i,t,A,K;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 18 "A:=matrixcd(n,k2);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "K:=[eigenvectors(A)];" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 45 "K:=sort(K,(a,b)->evalb(Re(a[1]) < Re(b[1])));" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 14 "c:=K[q][3][1];" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 8 "s1:=0.0;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 40 "for i from 1 to n+ 1 do s1:=s1+c[i]^2 od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 16 "s2:=0.5*c [1]^2; " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 43 "for i from 1 to n do s2: =s2-c[i]*c[i+1] od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 40 "s:=(1-0.5*k2 )*s1+0.5*k2*s2;s:=1/sqrt(s);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 7 "t:=0 .0;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 44 "for i from 1 to n+1 do t:=t+ (2*i-1)*c[i] od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 22 "if t<0 then s:= -s fi; " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 39 "for i from 1 to n+1 do c [i]:=s*c[i] od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 2 "c;" }}{PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 4 "end;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 17 "cdE:=proc(m,n,k2)" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "local c,p, E;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "c:=eigencd(n,k2,m+1);" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 5 "E:=0;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 64 "for p from 0 to n do E:=E+c[p+1]*U(2*p,JacobiSN(z,sqrt(k2))) o d;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 47 "E:=JacobiDN(z,sqrt(k2))*Jacob iCN(z,sqrt(k2))*E;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 4 "end;" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 45 "matrixscd:=proc(n,k2) # matr ix for type scd " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 13 "local A,r,nu; " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 10 "nu:=2*n+3;" }}{PARA 0 "> " 0 " " {MPLTEXT 1 0 21 "A:=matrix(n+1,n+1,0);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 53 "for r from 0 to n do A[r+1,r+1]:=(2*r+2)^2*(2-k2) od; " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 68 "for r from 0 to n-1 do A[r+1,r+ 2]:=0.5*(nu-2*r-2)*(nu+2*r+3)*k2 od; " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 63 "for r from 1 to n do A[r+1,r]:=0.5*(nu-2*r-1)*(nu+2*r+2)*k2 od; " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 2 "A;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 4 "end;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 23 "eig enscd:=proc(n,k2,q) " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 24 "local c,s1, s2,s,i,t,A,K;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 19 "A:=matrixscd(n,k2) ;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 21 "K:=[eigenvectors(A)];" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 45 "K:=sort(K,(a,b)->evalb(Re(a[1]) < R e(b[1])));" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 14 "c:=K[q][3][1];" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 8 "s1:=0.0;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 40 "for i from 1 to n+1 do s1:=s1+c[i]^2 od;" }}{PARA 0 " > " 0 "" {MPLTEXT 1 0 9 "s2:=0.0; " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 43 "for i from 1 to n do s2:=s2+c[i]*c[i+1] od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 40 "s:=(1-0.5*k2)*s1-0.5*k2*s2;s:=1/sqrt(s);" }}{PARA 0 " > " 0 "" {MPLTEXT 1 0 7 "t:=0.0;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 40 "for i from 1 to n+1 do t:=t+2*i*c[i] od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 22 "if t<0 then s:=-s fi; " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 39 "for i from 1 to n+1 do c[i]:=s*c[i] od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 2 "c;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 4 "end;" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 18 "scdE:=proc(m,n,k2)" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "local c,p,E;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 22 "c:=eigenscd(n,k2,m+1);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 5 "E:=0;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 66 "for p from 0 to n d o E:=E+c[p+1]*U(2*p+1,JacobiSN(z,sqrt(k2))) od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 47 "E:=JacobiDN(z,sqrt(k2))*JacobiCN(z,sqrt(k2))*E;" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 4 "end;" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 50 "Example 1: Check orthogonality and normalization. " }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 20 "n:=2;k:=0.7;k2:=k*k;" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 16 "K:=EllipticK(k);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 17 "e0:=scdE(0,n,k2);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 17 "e1:=scdE(1,n,k2);" }}}{EXCHG {PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 17 "e2:=scdE(2,n,k2);" }}}{EXCHG {PARA 0 "> " 0 " " {MPLTEXT 1 0 25 "evalf(int(e0*e1,z=0..K));" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 25 "evalf(int(e0*e2,z=0..K));" }}}{EXCHG {PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 25 "evalf(int(e1*e2,z=0..K));" }}}{EXCHG {PARA 0 " > " 0 "" {MPLTEXT 1 0 43 "evalf(4/Pi*int(JacobiDN(z,k)*e0^2,z=0..K)); " }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 43 "evalf(4/Pi*int(JacobiDN (z,k)*e1^2,z=0..K));" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 56 "eva lf(4/Pi*int(JacobiDN(z,k)*e1^2,z=0..K)); # looks good" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 39 "Example 2: Plot these Lame polynomials." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 12 "with(plots):" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 37 "p0:=plot(e0,z=-2*K..2*K,color=black ):" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 35 "p1:=plot(e1,z=-2*K..2 *K,color=red):" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 36 "p2:=plot( e2,z=-2*K..2*K,color=blue):" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 92 "display(p0,p1,p2,title=\"Lame polynomials uE_2n^m for n=2, k=0.7, \+ m=0,1,2 (black,red,blue)\");" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 77 "E xample 3: Plot the magnitude of Lame polynomials as a function of comp lex z." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 15 "k:=0.7;k2:=k^2;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 14 "f:=uE(1,2,k2);" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 16 "K:=EllipticK(k);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 18 "CK:=EllipticCK(k);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 24 "g:=abs(subs(z=x+I*y,f));" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 76 "plot3d(g,x=-3*K...3*K,y=0..2 *CK,axes=boxed,view=0..15,orientation=[-45,45]);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}}{MARK "5 0 0" 2 }{VIEWOPTS 1 1 0 1 1 1803 1 1 1 1 }{PAGENUMBERS 0 1 2 33 1 1 }