{VERSION 3 0 "IBM INTEL NT" "3.0" } {USTYLETAB {CSTYLE "Maple Input" -1 0 "Courier" 0 1 255 0 0 1 0 1 0 0 1 0 0 0 0 }{CSTYLE "2D Math" -1 2 "Times" 0 1 0 0 0 0 0 0 2 0 0 0 0 0 0 }{CSTYLE "2D Comment" 2 18 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 } {CSTYLE "2D Input" 2 19 "" 0 1 255 0 0 1 0 0 0 0 0 0 0 0 0 }{CSTYLE "2 D Output" 2 20 "" 0 1 0 0 255 1 0 0 0 0 0 0 0 0 0 }{CSTYLE "" 0 21 "" 0 1 0 0 0 1 0 0 0 0 2 0 0 0 0 }{PSTYLE "Normal" -1 0 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 } {PSTYLE "Heading 1" 0 3 1 {CSTYLE "" -1 -1 "" 1 18 0 0 0 0 0 1 0 0 0 0 0 0 0 }1 0 0 0 6 6 0 0 0 0 0 0 -1 0 }{PSTYLE "Heading 2" 3 4 1 {CSTYLE "" -1 -1 "" 1 14 0 0 0 0 0 0 0 0 0 0 0 0 0 }0 0 0 -1 4 4 0 0 0 0 0 0 -1 0 }{PSTYLE "Heading 3" 4 5 1 {CSTYLE "" -1 -1 "" 1 12 0 0 0 0 1 0 0 0 0 0 0 0 0 }0 0 0 -1 0 0 0 0 0 0 0 0 -1 0 }{PSTYLE "Maple O utput" 0 11 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }3 3 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 11 12 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }1 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "Title" 0 18 1 {CSTYLE "" -1 -1 "" 1 18 0 0 0 0 0 1 1 0 0 0 0 0 0 }3 0 0 -1 12 12 0 0 0 0 0 0 19 0 }{PSTYLE "Author" 0 19 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }3 0 0 -1 8 8 0 0 0 0 0 0 -1 0 }} {SECT 0 {PARA 18 "" 0 "" {TEXT -1 47 "Convergence Rates of Matrix Subd ivision Schemes" }}{PARA 19 "" 0 "" {TEXT -1 22 "Denis Zorin, 1997-200 0" }}{SECT 0 {PARA 3 "" 0 "" {TEXT -1 12 "Introduction" }}{PARA 0 "" 0 "" {TEXT -1 424 "The procedures defined in this worksheet can be us ed to estimate the convergence rates of scalar and matrix subdivision \+ schemes. It is based on the generalization to the matrix case of sev eral theorems in \"Stationary Subdivision\" by Cavaretta, Dahmen and M icchelli [CDM], and to the multivariate case of several theorems fro m \"Matrix subdivision\" by Cohen, Dyn and Levin [CDL]. The goal is \+ to compute four constants " }{XPPEDIT 18 0 "gamma, c, gamma[D], c[D] " "6&%&gammaG%\"cG&F#6#%\"DG&F$6#F'" }{TEXT -1 58 " , such that a give n (possibly matrix) subdivision scheme " }{XPPEDIT 18 0 "S" "6#%\"SG" }{TEXT -1 28 " and its difference scheme " }{XPPEDIT 18 0 "S[D]" "6#& %\"SG6#%\"DG" }{TEXT -1 12 " satisfy \n " }{XPPEDIT 18 0 "D[](S^n*p) \+ " 0 "" {XPPEDIT 19 1 "with(linalg): kernelopts(ASSERT = tru e):" "6#C$-%%withG6#%'linalgG-%+kerneloptsG6#/%'ASSERTG%%trueG" }}} {PARA 0 "" 0 "" {TEXT -1 0 "" }}{SECT 0 {PARA 0 "" 0 "" {MPLTEXT 0 21 20 "IsLaurentPolynom(A) " }{TEXT -1 92 "Check if the argument is a s calar Laurent polynomial in z1,z2 with arbitrary coefficients." }} {EXCHG {PARA 0 "> " 0 "" {XPPEDIT 19 1 "IsLaurentPolynom := proc(A::ra tpoly(anything,[z1,z2]))\nglobal z1,z2; type(denom(factor(A)), monomi al(anything,[z1,z2]));\nend:" "6#>%1IsLaurentPolynomGR6#'%\"AG-%(ratpo lyG6$%)anythingG7$%#z1G%#z2G7\"6\"F1-%%typeG6$-%&denomG6#-%'factorG6#F (-%)monomialG6$F,7$F.F/F16$F.F/F1" }}}}{SECT 0 {PARA 0 "" 0 "" {MPLTEXT 0 21 19 "IsLaurentMatrix(A) " }{TEXT -1 92 "Check if the arg ument is a matrix Laurent polynomial in z1,z2 with arbitrary coeffici ents." }}{EXCHG {PARA 0 "> " 0 "" {XPPEDIT 19 1 "IsLaurentMatrix := pr oc(A::matrix(ratpoly(anything,[z1,z2]))) local res, ind; \n# has to do this manually, 'matrix'(LaurentPolynom) does not work\n res := true ; \n for ind in indices(A) do \n res := res and IsLaurentPolynom (A[op(ind)]); \n od;\n res;\nend:" "6#>%0IsLaurentMatrixGR6#'%\"AG-% 'matrixG6#-%(ratpolyG6$%)anythingG7$%#z1G%#z2G7$%$resG%$indG6\"F6C%>F4 %%trueG?&F5-%(indicesG6#F(F9>F43F4-%1IsLaurentPolynomG6#&F(6#-%#opG6#F 5F4F6F6F6" }}}}{SECT 0 {PARA 0 "" 0 "" {MPLTEXT 0 21 24 "IsIntegerVect or2List(W) " }{TEXT -1 253 "Check if the argument is a list of vectors with integer components, each vector of length 2 each component in th e range 0..3, and if one of the components is 0, the other is also 0. \+ Used to check validity of arguments representing contraction functions ." }}{EXCHG {PARA 0 "> " 0 "" {XPPEDIT 19 1 "IsIntegerVector2List := p roc(W::list(list(integer))) local i,res;\n res := true;\n for i from 1 to vectdim(W) do \n res := res and evalb(vectdim(op(i,W)) = 2); \n res := res and evalb(W[i][1] >= 0 and W[i][1] <= 3);\n res := res and evalb(W[i][2] >= 0 and W[i][2] <= 3);\n res := res and eva lb( W[i][1] <> W[i][2]) or ( ( W[i][1] = 0) and ( W[i][2] = 0) ); \n od;\n res;\nend:" "6#>%5IsIntegerVector2ListGR6#'%\"WG-%%listG6#-F *6#%(integerG7$%\"iG%$resG6\"F2C%>F1%%trueG?(F0\"\"\"\"\"\"-%(vectdimG 6#F(F5C&>F13F1-%&evalbG6#/-F:6#-%#opG6$F0F(\"\"#>F13F1-F@6#31\"\"!&&F( 6#F06#\"\"\"1&&F(6#F06#\"\"\"\"\"$>F13F1-F@6#31FO&&F(6#F06#\"\"#1&&F(6 #F06#\"\"#\"\"$>F153F1-F@6#0&&F(6#F06#\"\"\"&&F(6#F06#\"\"#3/&&F(6#F06 #\"\"\"FO/&&F(6#F06#\"\"#FOF1F2F2F2" }}}}{SECT 0 {PARA 4 "" 0 "" {TEXT -1 17 " Type definitions" }}{EXCHG {PARA 0 "> " 0 "" {XPPEDIT 19 1 "`type/LaurentPolynom` := eval(IsLaurentPolynom):" "6#>%4type/Lau rentPolynomG-%%evalG6#%1IsLaurentPolynomG" }}}{EXCHG {PARA 0 "> " 0 " " {XPPEDIT 19 1 "`type/LaurentMatrix` := eval(IsLaurentMatrix):" "6#> %3type/LaurentMatrixG-%%evalG6#%0IsLaurentMatrixG" }}}{EXCHG {PARA 0 " > " 0 "" {XPPEDIT 19 1 "`type/IntegerVector2List` := eval(IsIntegerVec tor2List):" "6#>%8type/IntegerVector2ListG-%%evalG6#%5IsIntegerVector2 ListG" }}}}{SECT 0 {PARA 0 "" 0 "" {MPLTEXT 0 21 25 "MakeContractionMa trix(W) " }{TEXT -1 70 "Make a matrix Laurent polynomial out of a list of vectors W such that " }{MPLTEXT 0 21 28 "IsIntegerVector2List(W,1, 3) " }{TEXT -1 34 "is true. The matrix is defined as " }{XPPEDIT 18 0 "diag(p[1]..p[n])" "6#-%%diagG6#;&%\"pG6#\"\"\"&F(6#%\"nG" }{TEXT -1 9 " , where " }{XPPEDIT 18 0 "p[i]" "6#&%\"pG6#%\"iG" }{TEXT -1 76 " \+ are vectors corresponding to vectors of W, obtained by replacing 1 w ith " }{XPPEDIT 18 0 "z[1]^2-1" "6#,&*$&%\"zG6#\"\"\"\"\"#\"\"\"\"\"\" !\"\"" }{TEXT -1 10 ", 2 with " }{XPPEDIT 18 0 " z[2]^2-1" "6#,&*$&% \"zG6#\"\"#\"\"#\"\"\"\"\"\"!\"\"" }{TEXT -1 13 " and 3 with " } {XPPEDIT 18 0 " z[1]^2*z[2]^2 -1" "6#,&*&&%\"zG6#\"\"\"\"\"#&F&6#\"\"# \"\"#\"\"\"\"\"\"!\"\"" }{TEXT -1 2 ". " }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{EXCHG {PARA 0 "> " 0 "" {XPPEDIT 19 1 "MakeContractionMatrix := \+ proc(W::IntegerVector2List) local p, i,Wmatrix; global z1,z2; \n p := [z1-1,z2-1,z1*z2-1];\n Wmatrix := matrix(2*vectdim(W),vectdim(W),0); \n for i from 1 to vectdim(W) do \n Wmatrix[2*i-1,i] := p[W[i][1] ];\n Wmatrix[2*i,i] := p[W[i][2]];\n od;\n evalm(Wmatrix);\n en d:" "6#>%6MakeContractionMatrixGR6#'%\"WG%3IntegerVector2ListG7%%\"pG% \"iG%(WmatrixG6\"F.C&>F+7%,&%#z1G\"\"\"\"\"\"!\"\",&%#z2GF4\"\"\"F6,&* &F3F4F8F4F4\"\"\"F6>F--%'matrixG6%*&\"\"#F4-%(vectdimG6#F(F4-FD6#F(\" \"!?(F,\"\"\"F4-FD6#F(%%trueGC$>&F-6$,&*&\"\"#F4F,F4F4\"\"\"F6F,&F+6#& &F(6#F,6#\"\"\">&F-6$*&\"\"#F4F,F4F,&F+6#&&F(6#F,6#\"\"#-%&evalmG6#F-F .6$F3F8F." }}}}{SECT 0 {PARA 0 "" 0 "" {TEXT -1 1 " " }{MPLTEXT 0 21 17 "MakePolynom(coef)" }{TEXT -1 77 " make a polynomial in z1,z2 out \+ of a two-dimesional array of coefficients. " }}{EXCHG {PARA 0 "> " 0 "" {XPPEDIT 19 1 "MakePolynom := proc(polycoef::array(2)) local ind, p ;\nglobal z1,z2; p := 0;\n for ind in indices(polycoef) do\n p := p + z1^op(1,ind)*z2^op(2,ind)*polycoef[op(ind)];\n od;\n eval(p);\n end:" "6#>%,MakePolynomGR6#'%)polycoefG-%&arrayG6#\"\"#7$%$indG%\"pG6 \"F0C%>F/\"\"!?&F.-%(indicesG6#F(%%trueG>F/,&F/\"\"\"*()%#z1G-%#opG6$ \"\"\"F.F;)%#z2G-F@6$\"\"#F.F;&F(6#-F@6#F.F;F;-%%evalG6#F/F06$F>FDF0" }}}}{SECT 0 {PARA 0 "" 0 "" {MPLTEXT 0 21 11 "CoefList(A)" }{TEXT -1 88 " Make a 2d table of matrix coefficients of a matrix Laurent poly nomial in 2 variables." }}{EXCHG {PARA 0 "" 0 "" {TEXT -1 68 " con vert polynomials to tables of coefficient indexed by powers " }}{PARA 0 "> " 0 "" {XPPEDIT 19 1 "extractCoeffs := proc(x::LaurentPolynom) lo cal i, dummy, polycoeffs,terms, powers, extrPowers; global z1,z2; \n \+ polycoeffs := [coeffs(expand(factor(x)), [z1, z2], 'terms')]; \n \+ extrPowers := unapply( [ 'degree'(dummy, z1), 'degree'(dummy, z2)], \+ dummy); \n powers := map(extrPowers, [terms]); \n table( [ seq( ( op( op(i, powers) ) = op(i, polycoeffs ), i=1..vectdim(polycoef fs) ))] ); \n end: \n \n " "6#>%.extractCoeffsGR6#'%\"xG%/LaurentPolyn omG7(%\"iG%&dummyG%+polycoeffsG%&termsG%'powersG%+extrPowersG6\"F1C&>F -7#-%'coeffsG6%-%'expandG6#-%'factorG6#F(7$%#z1G%#z2G.F.>F0-%(unapplyG 6$7$-.%'degreeG6$F,F?-.FI6$F,F@F,>F/-%$mapG6$F07#F.-%&tableG6#7#-%$seq G6$/-%#opG6#-Ffn6$F+F/-Ffn6$F+F-/F+;\"\"\"-%(vectdimG6#F-F16$F?F@F1" } }{PARA 0 "" 0 "" {TEXT -1 34 "given a coefficient, make a matrix" }} {PARA 0 "> " 0 "" {XPPEDIT 19 1 "makeCoeffMatrix := proc(A::LaurentMat rix, i::integer,j::integer) \n local l,m, matrcoeff;\n matrcoeff : = matrix( rowdim(A), coldim(A),0 );\n for l from 1 to rowdim(A) do\n for m from 1 to coldim(A) do\n if member( [i,j], \{indices( A[l,m])\}) then\n matrcoeff[l,m] := A[l,m][i,j];\n fi;\n od;\n od;\n eval(matrcoeff);\n end:\n" "6#>%0makeCoeffMatrixG R6%'%\"AG%.LaurentMatrixG'%\"iG%(integerG'%\"jGF,7%%\"lG%\"mG%*matrcoe ffG6\"F3C%>F2-%'matrixG6%-%'rowdimG6#F(-%'coldimG6#F(\"\"!?(F0\"\"\"\" \"\"-F:6#F(%%trueG?(F1\"\"\"FB-F=6#F(FE@$-%'memberG6$7$F+F.<#-%(indice sG6#&F(6$F0F1>&F26$F0F1&&F(6$F0F16$F+F.-%%evalG6#F2F3F3F3" }}{PARA 0 " > " 0 "" {XPPEDIT 19 1 "CoeffList := proc(A::LaurentMatrix) local x, m atrcoeff,n, zeroMatrix,\n ld1, ld2, ud1, ud2,\n l,m, pcoeff, allp owers; \n global extractCoeffs, makeCoeffMatrix;\n pcoeff := map( \+ extractCoeffs, A );\n allpowers := map( unapply( \{'indices'(x)\},x) , pcoeff);\n allpowers := [op( `union`(op( map( op, [entries( allpow ers ) ] ))))];\n ld1 := min( op( map( unapply( 'op'(1,x), x ), allpo wers ) ));\n ld2 := min( op( map( unapply( 'op'(2,x), x ), allpowers ) ));\n ud1 := max( op( map( unapply( 'op'(1,x), x ), allpowers ) ) );\n ud2 := max( op( map( unapply( 'op'(2,x), x ), allpowers ) ));\n matrcoeff := array(ld1..ud1, ld2..ud2); \n for n in allpowers do \n matrcoeff[op(n)] := makeCoeffMatrix( pcoeff, op(n) );\n od; \n zeroMatrix := matrix(rowdim(A), coldim(A),0);\n \n for l \+ from ld1 to ud1 do\n for m from ld2 to ud2 do\n if not membe r( [l,m], allpowers) then\n matrcoeff[l,m] := evalm(zeroMatrix ); \n fi;\n od;\n od;\n eval(matrcoeff);\nend:" "6#>%*Co effListGR6#'%\"AG%.LaurentMatrixG7.%\"xG%*matrcoeffG%\"nG%+zeroMatrixG %$ld1G%$ld2G%$ud1G%$ud2G%\"lG%\"mG%'pcoeffG%*allpowersG6\"F7C.>F5-%$ma pG6$%.extractCoeffsGF(>F6-F;6$-%(unapplyG6$<#-.%(indicesG6#F+F+F5>F67# -%#opG6#-%&unionG6#-FL6#-F;6$FL7#-%(entriesG6#F6>F/-%$minG6#-FL6#-F;6$ -FB6$-.FL6$\"\"\"F+F+F6>F0-Fen6#-FL6#-F;6$-FB6$-.FL6$\"\"#F+F+F6>F1-%$ maxG6#-FL6#-F;6$-FB6$-.FL6$\"\"\"F+F+F6>F2-F`p6#-FL6#-F;6$-FB6$-.FL6$ \"\"#F+F+F6>F,-%&arrayG6$;F/F1;F0F2?&F-F6%%trueG>&F,6#-FL6#F--%0makeCo effMatrixG6$F5-FL6#F->F.-%'matrixG6%-%'rowdimG6#F(-%'coldimG6#F(\"\"!? (F3F/\"\"\"F1F`r?(F4F0FgsF2F`r@$4-%'memberG6$7$F3F4F6>&F,6$F3F4-%&eval mG6#F.-%%evalG6#F,F76$F=FgrF7" }}{PARA 0 "" 0 "" {TEXT -1 1 "\010" } {MPLTEXT 1 0 0 "" }}}}{SECT 0 {PARA 0 "" 0 "" {MPLTEXT 0 21 71 "Genera teCode(OutputFile,SchemeName,Consts,Mask,ControlSize,triangular) " } {TEXT -1 138 "Generate code for a function that creates an instance o f class RegScheme, defining a scheme. See file regscheme.h for class \+ definition. " }{MPLTEXT 0 21 6 "Consts" }{TEXT -1 82 " is a list of t ables of convergence constants, each table generted by a call to " } {MPLTEXT 0 21 21 "ConvergenceEstimates " }{TEXT -1 23 "or a similar fu nction. " }{MPLTEXT 0 21 5 "Mask " }{TEXT -1 48 "is the standard coeff icient mask of the scheme,\n" }{MPLTEXT 0 21 12 "ControlSize " }{TEXT -1 330 "is the number of layers of control vertices outside a patch re quired to define the surface on the patch (clearly ths number can be c omputed from the mask; however, computing the minimal number in the g eneral case requires some care; typically, in each specific case this number is obvious, and we simply specify it explicitly); " }{MPLTEXT 0 21 11 "triangular " }{TEXT -1 129 " indicates the symmetry type: if \+ it is true, than 3-directional symmetry is assumed; otherwise 2-direct ional\nsymmetry is assumed." }}{EXCHG {PARA 0 "> " 0 "" {XPPEDIT 19 1 "GenerateCode := proc( OutputFile::string, SchemeName::string, Consts: :list(table), Mask::array(2), ControlSize::integer, triangular::boolea n) \n local i, g, dg, mask, mmin1, mmin2, mmax1, mmax2;\n global S; \n assert( ControlSize >= 0, `ControlSize must be nonnegative`);\n f printf(OutputFile, `RegScheme* %s() \{\\n`, SchemeName);\n mask := ev al(Mask);\n for i in indices( Mask ) do\n if Mask[op(i)] = 0 then \+ mask[op(i)] := BV; else mask[op(i)] := Mask[op(i)]; fi;\n od;\n mask := map( ConvertToFloat, eval(mask) ); \n mmin1 := min( op( map( unap ply( 'op'(1, x), x), [indices(Mask)]))); \n mmin2 := min( op( map( un apply( 'op'(2, x), x), [indices(Mask)]))); \n mmax1 := max( op( map( \+ unapply( 'op'(1, x), x), [indices(Mask)]))); \n mmax2 := max( op( map ( unapply( 'op'(2, x), x), [indices(Mask)]))); \n fprintf(OutputFile, ` Float* g = new Float[%d];\\n`, vectdim(Consts) ); \n fprintf(Output File,` Float* dg = new Float[%d];\\n`, vectdim(Consts) ); \n fprintf( OutputFile,` Float** mask = new pFloat[%d];\\n`, mmax1 - mmin1 + 1 ); \+ \n fprintf(OutputFile,` for( int i = 0; i < %d; i++)\\n mask[i] = n ew Float[%d];\\n`, \n mmax1 - mmin1 + 1,mmax2 - mmin2 + 1 ); \n \+ g := array(1..vectdim(Consts));\n dg := array(1..vectdim(Consts));\n \+ for i from 1 to vectdim(Consts) do \n g[i] := Consts[i][`gamma`]; \n dg[i] := Consts[i][`gammaD`]; \n od;\n g := map(Convert ToFloat, eval(g)); \n dg := map(ConvertToFloat,eval(dg)); \n S := [ ]; \n C(mask); C(g); C(dg);\n for i from 1 to vectdim(S) do \n fp rintf(OutputFile, `%s`, S[i]);\n od;\n \n fprintf( OutputFile, `ret urn new RegScheme( %A, %d, g, %A, dg, new Grid( %d, %d, \+ %d, %d, mask), %d, %s); \\n\}\\n\\n`, \n ConvertToFloat(Consts[1][ `c`]), vectdim(Consts), ConvertToFloat(Consts[1][`cD`]), mmin1, mmax1, mmin2, mmax2, ControlSize, triangular ); \n NULL;\n\nend: " "6#>%-G enerateCodeGR6('%+OutputFileG%'stringG'%+SchemeNameGF)'%'ConstsG-%%lis tG6#%&tableG'%%MaskG-%&arrayG6#\"\"#'%,ControlSizeG%(integerG'%+triang ularG%(booleanG7*%\"iG%\"gG%#dgG%%maskG%&mmin1G%&mmin2G%&mmax1G%&mmax2 G6\"FGC;-%'assertG6$1\"\"!F9%@ControlSize~must~be~nonnegativeG-%(fprin tfG6%F(%3RegScheme*~%s()~|fr|+GF+>FB-%%evalG6#F3?&F?-%(indicesG6#F3%%t rueG@%/&F36#-%#opG6#F?FM>&FB6#-F[o6#F?%#BVG>&FB6#-F[o6#F?&F36#-F[o6#F? >FB-%$mapG6$%/ConvertToFloatG-FU6#FB>FC-%$minG6#-F[o6#-F^p6$-%(unapply G6$-.F[o6$\"\"\"%\"xGFbq7#-FY6#F3>FD-Fep6#-F[o6#-F^p6$-F\\q6$-.F[o6$\" \"#FbqFbq7#-FY6#F3>FE-%$maxG6#-F[o6#-F^p6$-F\\q6$-.F[o6$\"\"\"FbqFbq7# -FY6#F3>FF-Fhr6#-F[o6#-F^p6$-F\\q6$-.F[o6$\"\"#FbqFbq7#-FY6#F3-FP6%F(% <~Float*~g~=~new~Float[%d];|+G-%(vectdimG6#F--FP6%F(%=~Float*~dg~=~new ~Float[%d];|+G-F[u6#F--FP6%F(%A~Float**~mask~=~new~pFloat[%d];|+G,(FE \"\"\"FC!\"\"\"\"\"Ffu-FP6&F(%en~for(~int~i~=~0;~i~<~%d;~i++)|+~~~mask [i]~=~new~Float[%d];|+G,(FEFfuFCFgu\"\"\"Ffu,(FFFfuFDFgu\"\"\"Ffu>F@-F 56#;\"\"\"-F[u6#F->FA-F56#;\"\"\"-F[u6#F-?(F?\"\"\"Ffu-F[u6#F-FenC$>&F @6#F?&&F-6#F?6#%&gammaG>&FA6#F?&&F-6#F?6#%'gammaDG>F@-F^p6$F`p-FU6#F@> FA-F^p6$F`p-FU6#FA>%\"SG7\"-%\"CG6#FB-Fay6#F@-Fay6#FA?(F?\"\"\"Ffu-F[u 6#F^yFen-FP6%F(%#%sG&F^y6#F?-FP6-F(%bqreturn~new~RegScheme(~%A,~~%d,~g ,~%A,~dg,~new~Grid(~%d,~%d,~%d,~%d,~mask),~%d,~~~%s);~|+| hr|+|+G-F`p6#&&F-6#\"\"\"6#%\"cG-F[u6#F--F`p6#&&F-6#\"\"\"6#%#cDGFCFEF DFFF9F<%%NULLGFG6#F^yFG" }}}}}{SECT 0 {PARA 3 "" 0 "" {TEXT -1 50 "Con vergence estimates for bivariate matrix schemes" }}{PARA 0 "" 0 "" {TEXT -1 143 "Note: all polynomials are required to be matrix polynom ials; to use the functions for scalar polynomials, convert it to a 1 b y 1 matrix first." }}{SECT 0 {PARA 4 "" 0 "" {TEXT -1 48 "Decompose a \+ bivariate matrix Laurent polynomial " }}{SECT 0 {PARA 0 "" 0 "" {MPLTEXT 0 21 32 "IdealDecompose(A,i,j,symmetrize)" }{TEXT -1 32 " Giv en a Laurent polynomial " }{XPPEDIT 18 0 "A(z[1],z[2])" "6#-%\"AG6 $&%\"zG6#\"\"\"&F'6#\"\"#" }{TEXT -1 40 " known to be in the ideal ge nerated by " }{XPPEDIT 18 0 "z[1]^2 - 1" "6#,&*$&%\"zG6#\"\"\"\"\"#\" \"\"\"\"\"!\"\"" }{TEXT -1 6 " and " }{XPPEDIT 18 0 " z[2]^2 - 1" "6# ,&*$&%\"zG6#\"\"#\"\"#\"\"\"\"\"\"!\"\"" }{TEXT -1 19 " find polynomi als " }{XPPEDIT 18 0 "TA[1]" "6#&%#TAG6#\"\"\"" }{TEXT -1 4 "and " } {XPPEDIT 18 0 "TA[2]" "6#&%#TAG6#\"\"#" }{TEXT -1 14 " such that " }{XPPEDIT 18 0 "A = (z[1]^2 -1)*TA[1] + (z[2]^2-1)*TA[2]" "6#/%\"AG,&* &,&*$&%\"zG6#\"\"\"\"\"#\"\"\"\"\"\"!\"\"F.&%#TAG6#\"\"\"F.F.*&,&*$&F* 6#\"\"#\"\"#F.\"\"\"F0F.&F26#\"\"#F.F." }{TEXT -1 7 " , or " } {XPPEDIT 18 0 "A = (z[1]^2*z[2]^2-1)*TA[1]+(z[2]^2-1)*TA[2]" "6#/%\"AG ,&*&,&*&&%\"zG6#\"\"\"\"\"#&F*6#\"\"#\"\"#\"\"\"\"\"\"!\"\"F2&%#TAG6# \"\"\"F2F2*&,&*$&F*6#\"\"#\"\"#F2\"\"\"F4F2&F66#\"\"#F2F2" }{TEXT -1 6 " or " }{XPPEDIT 18 0 "A = (z[1]^2-1)*TA[1]+(z[1]^2*z[2]^2-1)*TA[2 ]" "6#/%\"AG,&*&,&*$&%\"zG6#\"\"\"\"\"#\"\"\"\"\"\"!\"\"F.&%#TAG6#\"\" \"F.F.*&,&*&&F*6#\"\"\"\"\"#&F*6#\"\"#\"\"#F.\"\"\"F0F.&F26#\"\"#F.F. " }{TEXT -1 254 " ; the expressions are derived from the expressions \+ of Lemma 2.3 [CDM] . The decomposition is not unique; the expressions are asymmetric even in the first case; the choice of the pair of pol ynomials is given by the arguments i and j; i corresponds to " } {XPPEDIT 18 0 "z[1]^2-1" "6#,&*$&%\"zG6#\"\"\"\"\"#\"\"\"\"\"\"!\"\"" }{TEXT -1 8 ", 2 to " }{XPPEDIT 18 0 " z[2]^2-1" "6#,&*$&%\"zG6#\"\"# \"\"#\"\"\"\"\"\"!\"\"" }{TEXT -1 11 " and 3 to " }{XPPEDIT 18 0 " z[ 1]^2*z[2]^2 -1" "6#,&*&&%\"zG6#\"\"\"\"\"#&F&6#\"\"#\"\"#\"\"\"\"\"\"! \"\"" }{TEXT -1 111 ". The order matters; the last argument is used \+ is used to request symmetrized decomposition with respect to " } {XPPEDIT 18 0 "z[1]^2-1" "6#,&*$&%\"zG6#\"\"\"\"\"#\"\"\"\"\"\"!\"\"" }{TEXT -1 8 ", and " }{XPPEDIT 18 0 " z[2]^2-1" "6#,&*$&%\"zG6#\"\"# \"\"#\"\"\"\"\"\"!\"\"" }{TEXT -1 99 " ; if one of i,j is 3, the last argument is ignored. The function returns a list of 2 polynomial" } {MPLTEXT 0 21 2 "s." }}{EXCHG {PARA 0 "> " 0 "" {XPPEDIT 19 1 "IdealDe compose := proc(A::LaurentPolynom, i::integer,j::integer,symmetrize::b oolean) \n local TA1, TA2, Q,p; global z1,z2; \n ASSERT( ( i >= 1 a nd i <= 3 and j >= 1 and j <= 3 and i <> j), \n `second and third a rgument cannot coincide and should be in the range 1..3`);\n ASSERT( \+ (simplify(subs( \{ z1 = -1, z2 = -1 \} , A )) = 0 ) and \n (s implify(subs( \{ z1 = 1, z2 = -1 \} , A )) = 0 ) and \n (sim plify(subs( \{ z1 = -1, z2 = 1 \} , A )) = 0 ) and \n (simpl ify(subs( \{ z1 = 1, z2 = 1 \} , A )) = 0 ) ,\n `Laurent polynom ial does not have all roots [(-1)^e, (-1)^g], e,g = 0,1` );\n\n p[1] \+ := z1^2 - 1; p[2] := z2^2-1; p[3] := z1^2*z2^2-1;\n\n Q[1,2] := ( (1 -z1)*subs( z1 = -1, A) + (1+z1)*subs( z1 = 1, A));\n Q[1,3] := ( (1 -z1)*subs( \{z2 = -z1*z2, z1 = -1\}, A) + \n (1+z1)*sub s( \{z2 = z1*z2, z1 = 1\}, A) );\n Q[2,1] := ( (1-z2)*subs( z2 = - 1, A) + (1+z2)*subs( z2 = 1, A));\n Q[2,3] := ( (1-z2)*subs( \{z1 = -z1*z2, z2 = -1\}, A) + \n (1+z2)*subs( \{z1 = z1*z2, z2 = 1\}, A));\n Q[3,1] := ( (1-z1*z2)*subs( z2 = -1/z1, A) + (1+ z1*z2)*subs( z2 = 1/z1, A));\n Q[3,2] := ( (1-z1*z2)*subs( z1 = -1/z 2, A) + (1+z1*z2)*subs( z1 = 1/z2, A)); \n if symmetrize and i = 1 and j = 2 then\n TA1 := factor( (1/(4*(z1^2-1)) )*(2*A + Q[2,1] - \+ Q[1,2]));\n TA2 := factor( (1/(4*(z2^2-1)) )*(2*A + Q[1,2] - Q[2,1] ));\n elif symmetrize and i = 2 and j = 1 then\n TA2 := factor( (1 /(4*(z1^2-1)) )*(2*A + Q[2,1] - Q[1,2]));\n TA1 := factor( (1/(4*(z 2^2-1)) )*(2*A + Q[1,2] - Q[2,1]));\n else \n TA1 := factor( Q[j,i ]/(2*p[i]) ) ;\n TA2 := factor( (2*A - Q[j,i])/(2*p[j]) ) ;\n fi; \n \n if simplify( p[i]*TA1 + p[j]*TA2 - A) <> 0 then \n \+ print( `decomposition error`);\n fi;\n [TA1, TA2];\nend: " "6#>%/I dealDecomposeGR6&'%\"AG%/LaurentPolynomG'%\"iG%(integerG'%\"jGF,'%+sym metrizeG%(booleanG7&%$TA1G%$TA2G%\"QG%\"pG6\"F7C0-%'ASSERTG6$33331\"\" \"F+1F+\"\"$1\"\"\"F.1F.\"\"$0F+F.%dosecond~and~third~argument~cannot~ coincide~and~should~be~in~the~range~1..3G-F:6$333/-%)simplifyG6#-%%sub sG6$<$/%#z1G,$\"\"\"!\"\"/%#z2G,$\"\"\"FenF(\"\"!/-FQ6#-FT6$<$/FX\"\" \"/Fgn,$\"\"\"FenF(Fjn/-FQ6#-FT6$<$/FX,$\"\"\"Fen/Fgn\"\"\"F(Fjn/-FQ6# -FT6$<$/FX\"\"\"/Fgn\"\"\"F(Fjn%aoLaurent~polynomial~does~not~have~all ~roots~[(-1)^e,~(-1)^g],~e,g~=~0,1G>&F66#\"\"\",&*$FX\"\"#\"\"\"\"\"\" Fen>&F66#\"\"#,&*$Fgn\"\"#Fcq\"\"\"Fen>&F66#\"\"$,&*&FX\"\"#Fgn\"\"#Fc q\"\"\"Fen>&F56$\"\"\"\"\"#,&*&,&\"\"\"FcqFXFenFcq-FT6$/FX,$\"\"\"FenF (FcqFcq*&,&\"\"\"FcqFXFcqFcq-FT6$/FX\"\"\"F(FcqFcq>&F56$\"\"\"\"\"$,&* &,&\"\"\"FcqFXFenFcq-FT6$<$/Fgn,$*&FXFcqFgnFcqFen/FX,$\"\"\"FenF(FcqFc q*&,&\"\"\"FcqFXFcqFcq-FT6$<$/Fgn*&FXFcqFgnFcq/FX\"\"\"F(FcqFcq>&F56$ \"\"#\"\"\",&*&,&\"\"\"FcqFgnFenFcq-FT6$/Fgn,$\"\"\"FenF(FcqFcq*&,&\" \"\"FcqFgnFcqFcq-FT6$/Fgn\"\"\"F(FcqFcq>&F56$\"\"#\"\"$,&*&,&\"\"\"Fcq FgnFenFcq-FT6$<$/FX,$*&FXFcqFgnFcqFen/Fgn,$\"\"\"FenF(FcqFcq*&,&\"\"\" FcqFgnFcqFcq-FT6$<$/FX*&FXFcqFgnFcq/Fgn\"\"\"F(FcqFcq>&F56$\"\"$\"\"\" ,&*&,&\"\"\"Fcq*&FXFcqFgnFcqFenFcq-FT6$/Fgn,$*&\"\"\"FcqFXFenFenF(FcqF cq*&,&\"\"\"Fcq*&FXFcqFgnFcqFcqFcq-FT6$/Fgn*&\"\"\"FcqFXFenF(FcqFcq>&F 56$\"\"$\"\"#,&*&,&\"\"\"Fcq*&FXFcqFgnFcqFenFcq-FT6$/FX,$*&\"\"\"FcqFg nFenFenF(FcqFcq*&,&\"\"\"Fcq*&FXFcqFgnFcqFcqFcq-FT6$/FX*&\"\"\"FcqFgnF enF(FcqFcq@'33F0/F+\"\"\"/F.\"\"#C$>F3-%'factorG6#*(\"\"\"Fcq*&\"\"%Fc q,&*$FX\"\"#Fcq\"\"\"FenFcqFen,(*&\"\"#FcqF(FcqFcq&F56$\"\"#\"\"\"Fcq& F56$\"\"\"\"\"#FenFcq>F4-Fd\\l6#*(\"\"\"Fcq*&\"\"%Fcq,&*$Fgn\"\"#Fcq\" \"\"FenFcqFen,(*&\"\"#FcqF(FcqFcq&F56$\"\"\"\"\"#Fcq&F56$\"\"#\"\"\"Fe nFcq33F0/F+\"\"#/F.\"\"\"C$>F4-Fd\\l6#*(\"\"\"Fcq*&\"\"%Fcq,&*$FX\"\"# Fcq\"\"\"FenFcqFen,(*&\"\"#FcqF(FcqFcq&F56$\"\"#\"\"\"Fcq&F56$\"\"\"\" \"#FenFcq>F3-Fd\\l6#*(\"\"\"Fcq*&\"\"%Fcq,&*$Fgn\"\"#Fcq\"\"\"FenFcqFe n,(*&\"\"#FcqF(FcqFcq&F56$\"\"\"\"\"#Fcq&F56$\"\"#\"\"\"FenFcqC$>F3-Fd \\l6#*&&F56$F.F+Fcq*&\"\"#Fcq&F66#F+FcqFen>F4-Fd\\l6#*&,&*&\"\"#FcqF(F cqFcq&F56$F.F+FenFcq*&\"\"#Fcq&F66#F.FcqFen@$0-FQ6#,(*&&F66#F+FcqF3Fcq Fcq*&&F66#F.FcqF4FcqFcqF(FenFjn-%&printG6#%4decomposition~errorG7$F3F4 F76$FXFgnF7" }}{PARA 12 "" 0 "" {TEXT -1 0 "" }}}}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{SECT 0 {PARA 0 "" 0 "" {MPLTEXT 0 21 36 "IdealDecompo seMatrix(A,W,symmetrize)" }{TEXT -1 38 " Given a matrix Laurent polyno mial " }{XPPEDIT 18 0 "A(z[1],z[2])" "6#-%\"AG6$&%\"zG6#\"\"\"&F'6# \"\"#" }{TEXT -1 40 " known to be in the ideal generated by " } {XPPEDIT 18 0 "z[1]^2 - 1" "6#,&*$&%\"zG6#\"\"\"\"\"#\"\"\"\"\"\"!\"\" " }{TEXT -1 6 " and " }{XPPEDIT 18 0 " z[2]^2 - 1" "6#,&*$&%\"zG6#\" \"#\"\"#\"\"\"\"\"\"!\"\"" }{TEXT -1 49 " find a matrix Laurent polyn omial Q such that " }{XPPEDIT 18 0 "A = QW(z[1]^2,z[2]^2)" "6#/%\"AG -%#QWG6$*$&%\"zG6#\"\"\"\"\"#*$&F*6#\"\"#\"\"#" }{TEXT -1 20 " ; comp onentwise, " }{XPPEDIT 18 0 " A[i,j] (z[1],z[2])= Q[i,2*j-1](z [1],z[2])*W[j,1](z[1]^2,z[2]^2) + Q[i,2*j](z[1],z[2])*W[i,2](z[1]^2,z[ 2]^2)\n" "6#/-&%\"AG6$%\"iG%\"jG6$&%\"zG6#\"\"\"&F,6#\"\"#,&*&-&%\"QG6 $F(,&*&\"\"#\"\"\"F)F;F;\"\"\"!\"\"6$&F,6#\"\"\"&F,6#\"\"#F;-&%\"WG6$F )\"\"\"6$*$&F,6#\"\"\"\"\"#*$&F,6#\"\"#\"\"#F;F;*&-&F66$F(*&\"\"#F;F)F ;6$&F,6#\"\"\"&F,6#\"\"#F;-&FG6$F(\"\"#6$*$&F,6#\"\"\"\"\"#*$&F,6#\"\" #\"\"#F;F;" }{TEXT -1 97 "\n Q is twice as wide as A and has the same number of rows. W is a list of pairs of integers " }{XPPEDIT 18 0 "p[i] " "6#&%\"pG6#%\"iG" }{TEXT -1 85 " of length n equal to the wi dth of A. Each integer is one of 1,2,3 corresponding to " }{XPPEDIT 18 0 "z[1]-1" "6#,&&%\"zG6#\"\"\"\"\"\"\"\"\"!\"\"" }{TEXT -1 3 ", " }{XPPEDIT 18 0 "z[2]-1" "6#,&&%\"zG6#\"\"#\"\"\"\"\"\"!\"\"" }{TEXT -1 6 " and " }{XPPEDIT 18 0 "z[1]*z[2]-1" "6#,&*&&%\"zG6#\"\"\"\"\"\" &F&6#\"\"#F)F)\"\"\"!\"\"" }{TEXT -1 107 " respectively. In each pair , the integers cannot coincide. The matrix Laurent polynomial W is de fined as " }{XPPEDIT 18 0 "diag(p[1]..p[n])" "6#-%%diagG6#;&%\"pG6#\" \"\"&F(6#%\"nG" }{TEXT -1 56 " . The last argument is used for the sam e purpose as in " }{MPLTEXT 0 21 14 "IdealDecompose" }{TEXT -1 1 "." } }{EXCHG {PARA 0 "> " 0 "" {XPPEDIT 19 1 "IdealDecomposeMatrix := proc( A::LaurentMatrix,W::IntegerVector2List,symmetrize::boolean) \n local \+ Q, i,j, Wmatrix,testA,Qdecomp; global z1,z2;\n ASSERT( vectdim(W) = c oldim(A), `the length of the 2nd argument (list) should be the same as width of the first argument (matrix)`);\n Q := matrix(rowdim(A),2*co ldim(A));\n for i from 1 to rowdim(A) do \n for j from 1 to coldim (A) do \n # this conditional is a tweak to make difference m atrices diagonal for \n # factorizable schemes; \n if( i m od 2 = 1 ) then\n Qdecomp := IdealDecompose( A[i,j],W[j][2] , W[j][1], symmetrize);\n Q[i,2*j-1] := op(2,Qdecomp); Q[i,2 *j] := op(1,Qdecomp);\n else \n Qdecomp := IdealDecomp ose( A[i,j],W[j][1], W[j][2], symmetrize);\n ASSERT( IsLaurent Polynom( op(1, Qdecomp)) and IsLaurentPolynom( op(2, Qdecomp)),`decomp osition failed` );\n Q[i,2*j-1] := op(1,Qdecomp); Q[i,2*j] := \+ op(2,Qdecomp);\n fi;\n od;\n od;\n Wmatrix := MakeContra ctionMatrix(W);\n testA := evalm( A - Q &* subs( \{z1 = z1^2, z2 = z 2^2\}, evalm(Wmatrix) ) ); \n ASSERT(norm(testA) = 0, `scheme wa s not decomposed correctly`);\n evalm(Q);\n end:" "6#>%5IdealDecompo seMatrixGR6%'%\"AG%.LaurentMatrixG'%\"WG%3IntegerVector2ListG'%+symmet rizeG%(booleanG7(%\"QG%\"iG%\"jG%(WmatrixG%&testAG%(QdecompG6\"F7C)-%' ASSERTG6$/-%(vectdimG6#F+-%'coldimG6#F(%[qthe~length~of~the~2nd~argume nt~(list)~should~be~the~same~as~width~of~the~first~argument~(matrix)G> F1-%'matrixG6$-%'rowdimG6#F(*&\"\"#\"\"\"-FA6#F(FM?(F2\"\"\"FM-FI6#F(% %trueG?(F3\"\"\"FM-FA6#F(FT@%/-%$modG6$F2\"\"#\"\"\"C%>F6-%/IdealDecom poseG6&&F(6$F2F3&&F+6#F36#\"\"#&&F+6#F36#\"\"\"F.>&F16$F2,&*&\"\"#FMF3 FMFM\"\"\"!\"\"-%#opG6$\"\"#F6>&F16$F2*&\"\"#FMF3FM-Fdp6$\"\"\"F6C&>F6 -F]o6&&F(6$F2F3&&F+6#F36#\"\"\"&&F+6#F36#\"\"#F.-F:6$3-%1IsLaurentPoly nomG6#-Fdp6$\"\"\"F6-Fcr6#-Fdp6$\"\"#F6%5decomposition~failedG>&F16$F2 ,&*&\"\"#FMF3FMFM\"\"\"Fbp-Fdp6$\"\"\"F6>&F16$F2*&\"\"#FMF3FM-Fdp6$\" \"#F6>F4-%6MakeContractionMatrixG6#F+>F5-%&evalmG6#,&F(FM-%#&*G6$F1-%% subsG6$<$/%#z1G*$Fau\"\"#/%#z2G*$Feu\"\"#-Fft6#F4Fbp-F:6$/-%%normG6#F5 \"\"!%Dscheme~was~not~decomposed~correctlyG-Fft6#F1F76$FauFeuF7" }}} {PARA 5 "" 0 "" {TEXT -1 0 "" }}}}{SECT 0 {PARA 4 "" 0 "" {TEXT -1 47 "Compute a difference scheme for a matrix scheme" }}{SECT 0 {PARA 0 " " 0 "" {MPLTEXT 0 21 15 "DiffScheme(A,W)" }{TEXT -1 141 " Compute the \+ difference matrix scheme Q for a bivariate matrix scheme; the choice o f directions for differences is given by W, defined as in " }{MPLTEXT 0 21 20 "IdealDecomposeMatrix" }{TEXT -1 73 ". We assume that scheme is specified by its matrix Laurent polynomial " }{XPPEDIT 18 0 "A(z[ 1],z[2])" "6#-%\"AG6$&%\"zG6#\"\"\"&F'6#\"\"#" }{TEXT -1 71 ", and rep roduces constant vectors; the calculation solves the equation" } {XPPEDIT 18 0 "W(z[1],z[2])*A (z[1],z[2]) = Q(z[1],z[2]) * W( z[1]^ 2 , z[2]^2 );" "6#/*&-%\"WG6$&%\"zG6#\"\"\"&F)6#\"\"#\"\"\"-%\"AG6$&F) 6#\"\"\"&F)6#\"\"#F/*&-%\"QG6$&F)6#\"\"\"&F)6#\"\"#F/-F&6$*$&F)6#\"\" \"\"\"#*$&F)6#\"\"#\"\"#F/" }{TEXT -1 30 " componentwise, \+ \n" }{XPPEDIT 18 0 " W[i,1](z[1],z[2])*A[i,j] (z[1],z[2])= Q[2 *i-1,2*j-1](z[1],z[2])*W[j,1](z[1]^2,z[2]^2) + Q[2*i-1,2*j](z[1],z[2]) *W[i,2](z[1]^2,z[2]^2)\n" "6#/*&-&%\"WG6$%\"iG\"\"\"6$&%\"zG6#\"\"\"&F -6#\"\"#\"\"\"-&%\"AG6$F)%\"jG6$&F-6#\"\"\"&F-6#\"\"#F3,&*&-&%\"QG6$,& *&\"\"#F3F)F3F3\"\"\"!\"\",&*&\"\"#F3F8F3F3\"\"\"FJ6$&F-6#\"\"\"&F-6# \"\"#F3-&F'6$F8\"\"\"6$*$&F-6#\"\"\"\"\"#*$&F-6#\"\"#\"\"#F3F3*&-&FD6$ ,&*&\"\"#F3F)F3F3\"\"\"FJ*&\"\"#F3F8F36$&F-6#\"\"\"&F-6#\"\"#F3-&F'6$F )\"\"#6$*$&F-6#\"\"\"\"\"#*$&F-6#\"\"#\"\"#F3F3" }{XPPEDIT 18 0 "W[i,2 ](z[1],z[2])*A[i,j] (z[1],z[2])= Q[2*i ,2*j-1](z[1],z[2])*W[i,1](z[1 ]^2,z[2]^2) + Q[2*i ,2*j](z[1],z[2])*W[i,2](z[1]^2,z[2]^2)\n" "6#/* &-&%\"WG6$%\"iG\"\"#6$&%\"zG6#\"\"\"&F-6#\"\"#\"\"\"-&%\"AG6$F)%\"jG6$ &F-6#\"\"\"&F-6#\"\"#F3,&*&-&%\"QG6$*&\"\"#F3F)F3,&*&\"\"#F3F8F3F3\"\" \"!\"\"6$&F-6#\"\"\"&F-6#\"\"#F3-&F'6$F)\"\"\"6$*$&F-6#\"\"\"\"\"#*$&F -6#\"\"#\"\"#F3F3*&-&FD6$*&\"\"#F3F)F3*&\"\"#F3F8F36$&F-6#\"\"\"&F-6# \"\"#F3-&F'6$F)\"\"#6$*$&F-6#\"\"\"\"\"#*$&F-6#\"\"#\"\"#F3F3" }} {EXCHG {PARA 0 "> " 0 "" {XPPEDIT 19 1 "DiffScheme := proc (A::Laurent Matrix,W::IntegerVector2List) \n local WA, Wmatrix, Q, testA; global z1,z2;\n ASSERT(coldim(A) = rowdim(A),`1st argument should be a squa re matrix`);\n ASSERT( vectdim(W) = coldim(A), `the length of the 2nd argument (list) should be the same as either dimension of the 1st arg ument (matrix)`);\n Wmatrix := MakeContractionMatrix(W); \n WA := ev alm( Wmatrix &* A ); \n Q := IdealDecomposeMatrix(WA,W,false);\n tes tA := evalm( WA - Q &* subs( \{z1 = z1^2, z2 = z2^2\}, evalm(Wmatrix) \+ ) ); \n ASSERT(norm(testA) = 0, `difference scheme was not found corr ectly`);\n evalm(Q);\nend:" "6#>%+DiffSchemeGR6$'%\"AG%.LaurentMatrix G'%\"WG%3IntegerVector2ListG7&%#WAG%(WmatrixG%\"QG%&testAG6\"F2C*-%'AS SERTG6$/-%'coldimG6#F(-%'rowdimG6#F(%G1st~argument~should~be~a~square~ matrixG-F56$/-%(vectdimG6#F+-F96#F(%dqthe~length~of~the~2nd~argument~( list)~should~be~the~same~as~either~dimension~of~the~1st~argument~(matr ix)G>F/-%6MakeContractionMatrixG6#F+>F.-%&evalmG6#-%#&*G6$F/F(>F0-%5Id ealDecomposeMatrixG6%F.F+%&falseG>F1-FN6#,&F.\"\"\"-FQ6$F0-%%subsG6$<$ /%#z1G*$F^o\"\"#/%#z2G*$Fbo\"\"#-FN6#F/!\"\"-F56$/-%%normG6#F1\"\"!%Jd ifference~scheme~was~not~found~correctlyG-FN6#F0F26$F^oFboF2" }}}} {PARA 0 "" 0 "" {TEXT -1 0 "" }}}{SECT 0 {PARA 4 "" 0 "" {TEXT -1 37 " Infinity norm for matrix polynomials " }}{SECT 0 {PARA 0 "" 0 "" {MPLTEXT 0 21 18 "InfNnormSums(A,n) " }{TEXT -1 21 "compute the list o f " }{XPPEDIT 18 0 " 4^n" "6#)\"\"%%\"nG" }{TEXT -1 73 " vectors of \+ sums of magnitudes of entries of matrix coefficients of " } {XPPEDIT 18 0 "A(z)*A(z^2)*A(z^4)..A(z^(2^n))" "6#;*(-%\"AG6#%\"zG\"\" \"-F&6#*$F(\"\"#F)-F&6#*$F(\"\"%F)-F&6#)F()\"\"#%\"nG" }{TEXT -1 8 " w here " }{XPPEDIT 18 0 "z = (z[1],z[2])" "6#/%\"zG6$&F$6#\"\"\"&F$6#\" \"#" }{TEXT -1 141 "; each component in each vector corresponds to \+ a row in the matrix coefficients sums[i,j] includes all matrix coeff icients of monomials " }{XPPEDIT 18 0 "z1^l * z2^k " "6#*&)%#z1G%\"lG \"\"\")%#z2G%\"kGF'" }{TEXT -1 13 " such that " }{XPPEDIT 18 0 " l m od 2^n = i-1 " "6#/-%$modG6$%\"lG)\"\"#%\"nG,&%\"iG\"\"\"\"\"\"!\"\"" }{TEXT -1 6 " and " }{XPPEDIT 18 0 " k mod 2^n = j-1" "6#/-%$modG6$% \"kG)\"\"#%\"nG,&%\"jG\"\"\"\"\"\"!\"\"" }}{PARA 0 "" 0 "" {TEXT -1 1 " " }}{EXCHG {PARA 0 "> " 0 "" {XPPEDIT 19 1 "InfNnormSums := proc(A:: LaurentMatrix, n::integer) \n local x, pn, cfs, i,j,k,l,m, coeffsums, \n ld1, ud1, ld2, ud2, ind1, ind2; \n global z1,z2; \n ASSERT( n \+ = 1 or rowdim(A) = coldim(A),`the matrix should be square if n > 1`); \+ \n pn := A;\n # if n = 1 this is never executed\n for i from 1 to n -1 do \n pn := evalm( pn &* subs( \{ z1 = z1^(2^i), z2 = z2^(2^i)\} , evalm(A)));\n od;\n pn := map( expand, evalm(pn)); \n cfs := Coef fList(pn);\n coeffsums := array(0..(2^n-1),0..(2^n-1));\n for i from 0 to 2^n-1 do\n for j from 0 to 2^n-1 do \n coeffsums[i,j] \+ := vector(rowdim(A),0);\n od;\n od;\n ind1 := map( unapply( 'op '(1, x ), x ), [indices( cfs )]);\n ind2 := map( unapply( 'op'(2, x ) , x ), [indices( cfs )]);\n ld1 := min(op(ind1)); ud1 := max(op(ind1) ); \n ld2 := min(op(ind2)); ud2 := max(op(ind2));\n for l from ld1 t o ud1 do\n for k from ld2 to ud2 do \n for m from 1 to coldim( A) do \n coeffsums[ l mod 2^n, k mod 2^n ] := evalm( coeffsu ms[ l mod 2^n, k mod 2^n ] + map( abs, col(cfs[l,k],m) ));\n od; \n od;\n od;\n eval(coeffsums);\nend:" "6#>%-InfNnormSumsGR6$'%\" AG%.LaurentMatrixG'%\"nG%(integerG71%\"xG%#pnG%$cfsG%\"iG%\"jG%\"kG%\" lG%\"mG%*coeffsumsG%$ld1G%$ud1G%$ld2G%$ud2G%%ind1G%%ind2G6\"F=C1-%'ASS ERTG6$5/F+\"\"\"/-%'rowdimG6#F(-%'coldimG6#F(%Ethe~matrix~should~be~sq uare~if~n~>~1G>F/F(?(F1\"\"\"\"\"\",&F+FP\"\"\"!\"\"%%trueG>F/-%&evalm G6#-%#&*G6$F/-%%subsG6$<$/%#z1G)F[o)\"\"#F1/%#z2G)F`o)\"\"#F1-FW6#F(>F /-%$mapG6$%'expandG-FW6#F/>F0-%*CoeffListG6#F/>F6-%&arrayG6$;\"\"!,&) \"\"#F+FP\"\"\"FS;Ffp,&)\"\"#F+FP\"\"\"FS?(F1FfpFP,&)\"\"#F+FP\"\"\"FS FT?(F2FfpFP,&)\"\"#F+FP\"\"\"FSFT>&F66$F1F2-%'vectorG6$-FG6#F(Ffp>F;-F ho6$-%(unapplyG6$-.%#opG6$\"\"\"F.F.7#-%(indicesG6#F0>F<-Fho6$-Ffr6$-. Fjr6$\"\"#F.F.7#-F_s6#F0>F7-%$minG6#-Fjr6#F;>F8-%$maxG6#-Fjr6#F;>F9-F_ t6#-Fjr6#F<>F:-Fet6#-Fjr6#F&F66$-%$modG6$F4)\"\"#F+-F]v6$F3)\"\"#F+-FW6#,&&F66$-F]v6$F4)\"\"# F+-F]v6$F3)\"\"#F+FP-Fho6$%$absG-%$colG6$&F06$F4F3F5FP-%%evalG6#F6F=6$ F[oF`oF=" }}}{PARA 5 "" 0 "" {TEXT -1 0 "" }}}{SECT 0 {PARA 0 "" 0 "" {MPLTEXT 0 21 14 "InfNnormV(A,n)" }{TEXT -1 47 " get the vector of co mponentwise inf norms of " }{XPPEDIT 18 0 "A(z)*A(z^2)*A(z^4)..A(z^(2^ n))" "6#;*(-%\"AG6#%\"zG\"\"\"-F&6#*$F(\"\"#F)-F&6#*$F(\"\"%F)-F&6#)F( )\"\"#%\"nG" }{TEXT -1 8 " where " }{XPPEDIT 18 0 "z = (z[1],z[2])" " 6#/%\"zG6$&F$6#\"\"\"&F$6#\"\"#" }{TEXT -1 77 ", as the componentwise maximum of the vector sums computed by InfNnormSums" }}{EXCHG {PARA 0 "> " 0 "" {XPPEDIT 19 1 "InfNnormV := proc(A::LaurentMatrix, n ::integer) local x,i,j,m, sums, vsumlist;\n for i from 1 to rowdim(A) do \n vsumlist[i] := [];\n od;\n sums := InfNnormSums(A,n);\n f or i from 0 to 2^n-1 do\n for j from 0 to 2^n-1 do\n for m fro m 1 to rowdim(A) do \n vsumlist[m] := [op(vsumlist[m]),sums[i,j ][m]];\n od;\n od; \n od;\n map( unapply( 'max'('op'(x)),x) , vsumlist );\n end:" "6#>%*InfNnormVGR6$'%\"AG%.LaurentMatrixG'%\"nG% (integerG7(%\"xG%\"iG%\"jG%\"mG%%sumsG%)vsumlistG6\"F4C&?(F/\"\"\"\"\" \"-%'rowdimG6#F(%%trueG>&F36#F/7\">F2-%-InfNnormSumsG6$F(F+?(F/\"\"!F8 ,&)\"\"#F+F8\"\"\"!\"\"F&F36#F17$-%#opG6#&F36#F1&&F26$F/F06#F1-%$mapG6$-%(unapplyG6$-.% $maxG6#-.FZ6#F.F.F3F4F4F4" }}{PARA 5 "> " 0 "" {MPLTEXT 1 0 0 "" }}}} {SECT 0 {PARA 0 "" 0 "" {MPLTEXT 0 21 13 "InfNnorm(A,n)" }{TEXT -1 22 " get the inf norm of " }{XPPEDIT 18 0 "A(z)*A(z^2)*A(z^4)..A(z^(2^n) )" "6#;*(-%\"AG6#%\"zG\"\"\"-F&6#*$F(\"\"#F)-F&6#*$F(\"\"%F)-F&6#)F() \"\"#%\"nG" }{TEXT -1 8 " where " }{XPPEDIT 18 0 "z = (z[1],z[2])" "6 #/%\"zG6$&F$6#\"\"\"&F$6#\"\"#" }{TEXT -1 75 ", as the maximum of all \+ components of vector sums computed by InfNnormSums." }}{EXCHG {PARA 0 "> " 0 "" {XPPEDIT 19 1 "InfNnorm := proc(A::LaurentMatrix, n::integer ) \n max(op(map(op,[entries(InfNnormV(A,n))])));\nend:" "6#>%)InfNnorm GR6$'%\"AG%.LaurentMatrixG'%\"nG%(integerG7\"6\"F.-%$maxG6#-%#opG6#-%$ mapG6$F37#-%(entriesG6#-%*InfNnormVG6$F(F+F.F.F." }}}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{PARA 5 "" 0 "" {TEXT -1 0 "" }}}{SECT 0 {PARA 4 "" 0 "" {TEXT -1 72 "Convergence rate estimates for a matrix scheme and i ts difference scheme" }}{SECT 0 {PARA 0 "" 0 "" {TEXT -1 1 " " } {MPLTEXT 0 21 38 "ConvergenceEstimates(A,n,nD, W, WD, B)" }{TEXT -1 31 " Compute convergence constants " }{XPPEDIT 18 0 "gamma, c, gamma[D ], c[D]" "6&%&gammaG%\"cG&F#6#%\"DG&F$6#F'" }{TEXT -1 143 " , describ ed in the introduction. The function returns a table indexed by names. The argument W allows one to define a contraction function " } {XPPEDIT 18 0 "abs(W)[infinity] " "6#&-%$absG6#%\"WG6#%)infinityG" } {TEXT -1 147 " , similarly, WD defines the contraction function for th e difference scheme. both W and WD are represented by lists of pairs of integers, as for " }{MPLTEXT 0 21 10 "DiffScheme" }{TEXT -1 194 ". The comparison scheme for the difference scheme is taken to be the \+ difference scheme of B . This is not necessary in general, but conveni ent for our purposes. We assume that B has factors " }{XPPEDIT 18 0 "w(z[1]^2,z[2]^2)/w(z[1],z[2])" "6#*&-%\"wG6$*$&%\"zG6#\"\"\"\"\"#*$&F )6#\"\"#\"\"#\"\"\"-F%6$&F)6#\"\"\"&F)6#\"\"#!\"\"" }{TEXT -1 39 ", w here w is any of the polynomials " }{XPPEDIT 18 0 "z[1]-1, z[2]-1, \+ z[1]*z[2] -1" "6%,&&%\"zG6#\"\"\"\"\"\"\"\"\"!\"\",&&F%6#\"\"#F(\"\"\" F*,&*&&F%6#\"\"\"F(&F%6#\"\"#F(F(\"\"\"F*" }{TEXT -1 14 " used in W. \+ " }}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 1330 "ConvergenceEstimates \+ := proc(A::LaurentMatrix, n::posint, nD::posint, \n W::IntegerVector 2List, WD::IntegerVector2List, B::LaurentPolynom) \n local BD, Q, Q2, DA, DQ, BM, i, p, res; global z1,z2, c, cD, gamma, gammaD; \n ASSER T(rowdim(A) = coldim(A),`1st argument should be a square matrix`);\n \+ ASSERT( vectdim(W) = coldim(A), `the length of the 4nd argument (list) should be the same as either dimension of the first argument (matrix) `);\n ASSERT( vectdim(WD) = 2*coldim(A), `the length of the 5nd argum ent (list) should be twice either dimension of the first argument (m atrix)`);\n p := [z1-1,z2-1,z1*z2-1]; \n BM:= evalm(B*diag( seq(1, i = 1..rowdim(A))));\n BD := map( factor, diag( seq( op( [\n 2* BM[i,i]*p[W[i][1]]/subs( \{z1 = z1^2, z2 = z2^2\}, p[W[i][1]]), \n \+ 2*BM[i,i]*p[W[i][2]]/subs( \{z1 = z1^2, z2 = z2^2\}, p[W[i][2]]) ]) ,\n i = 1..rowdim(A)))); \n ASSERT(IsLaurentMatrix(BD),`the co mparison scheme does not have necessary factors`); \n Q := DiffSchem e(A,W); \n Q2 := DiffScheme(evalm( 2*Q ),WD);\n DA := map(expand,map (factor, evalm(A - BM )));\n DQ := map(expand,map(factor, evalm(2*Q - BD)));\n res[gamma] := InfNnorm(Q,n);\n res[c] := InfNnorm(IdealDe composeMatrix(DA,W,true),1);\n res[gammaD] := InfNnorm(Q2,nD);\n res [cD] := InfNnorm(IdealDecomposeMatrix(DQ,WD,true),1); \n eval(res); \nend:" }}}}}}{SECT 0 {PARA 3 "" 0 "" {TEXT -1 46 "Estimates for schem es with a comparison factor" }}{PARA 0 "" 0 "" {TEXT -1 93 "Another ty pe of simplification can be made for schemes with Laurent polynomials of the form " }{XPPEDIT 18 0 "B(z[1],z[2])*Q(z[1],z[2])" "6#*&-%\"BG6 $&%\"zG6#\"\"\"&F(6#\"\"#\"\"\"-%\"QG6$&F(6#\"\"\"&F(6#\"\"#F." } {TEXT -1 34 " (see complete conditions below)." }}{SECT 0 {PARA 0 "" 0 "" {MPLTEXT 0 21 22 "LinearDecompose(A,i,j)" }{TEXT -1 32 " Given a \+ Laurent polynomial " }{XPPEDIT 18 0 "A(z[1],z[2])" "6#-%\"AG6$&%\" zG6#\"\"\"&F'6#\"\"#" }{TEXT -1 40 " known to be in the ideal generat ed by " }{XPPEDIT 18 0 "z[1]- 1" "6#,&&%\"zG6#\"\"\"\"\"\"\"\"\"!\"\" " }{TEXT -1 6 " and " }{XPPEDIT 18 0 " z[2] - 1" "6#,&&%\"zG6#\"\"#\" \"\"\"\"\"!\"\"" }{TEXT -1 19 " find polynomials " }{XPPEDIT 18 0 "TA [1]" "6#&%#TAG6#\"\"\"" }{TEXT -1 4 "and " }{XPPEDIT 18 0 "TA[2]" "6#& %#TAG6#\"\"#" }{TEXT -1 14 " such that " }{XPPEDIT 18 0 "A = (z[1] \+ -1)*TA[1] + (z[2]-1)*TA[2]" "6#/%\"AG,&*&,&&%\"zG6#\"\"\"\"\"\"\"\"\"! \"\"F,&%#TAG6#\"\"\"F,F,*&,&&F)6#\"\"#F,\"\"\"F.F,&F06#\"\"#F,F," } {TEXT -1 7 " , or " }{XPPEDIT 18 0 "A = (z[1]*z[2]-1)*TA[1]+(z[2]-1)* TA[2]" "6#/%\"AG,&*&,&*&&%\"zG6#\"\"\"\"\"\"&F*6#\"\"#F-F-\"\"\"!\"\"F -&%#TAG6#\"\"\"F-F-*&,&&F*6#\"\"#F-\"\"\"F2F-&F46#\"\"#F-F-" }{TEXT -1 6 " or " }{XPPEDIT 18 0 "A = (z[1]-1)*TA[1]+(z[1]*z[2]-1)*TA[2]" "6#/%\"AG,&*&,&&%\"zG6#\"\"\"\"\"\"\"\"\"!\"\"F,&%#TAG6#\"\"\"F,F,*&,& *&&F)6#\"\"\"F,&F)6#\"\"#F,F,\"\"\"F.F,&F06#\"\"#F,F," }{TEXT -1 254 " ; the expressions are derived from the expressions of Lemma 2.3 [CD M] . The decomposition is not unique; the expressions are asymmetric e ven in the first case; the choice of the pair of polynomials is given by the arguments i and j; i corresponds to " }{XPPEDIT 18 0 "z[1]-1" "6#,&&%\"zG6#\"\"\"\"\"\"\"\"\"!\"\"" }{TEXT -1 8 ", 2 to " } {XPPEDIT 18 0 " z[2]-1" "6#,&&%\"zG6#\"\"#\"\"\"\"\"\"!\"\"" }{TEXT -1 11 " and 3 to " }{XPPEDIT 18 0 " z[1]*z[2] -1" "6#,&*&&%\"zG6#\"\" \"\"\"\"&F&6#\"\"#F)F)\"\"\"!\"\"" }{TEXT -1 103 ". The order matters ; the last argument is used to request symmetrized decomposition with respect to " }{XPPEDIT 18 0 "z[1]-1" "6#,&&%\"zG6#\"\"\"\"\"\"\"\"\" !\"\"" }{TEXT -1 9 ", and " }{XPPEDIT 18 0 " z[2]-1" "6#,&&%\"zG6# \"\"#\"\"\"\"\"\"!\"\"" }{TEXT -1 98 " ; if one of i,j is 3, the last \+ argument is ignored . The function returns a list of 2 polynomial" } {MPLTEXT 0 21 2 "s." }}{EXCHG {PARA 0 "> " 0 "" {XPPEDIT 19 1 "LinearD ecompose := proc(A::LaurentPolynom, i::integer,j::integer,symmetrize:: boolean) \n local TA1, TA2, Q,p; \n global z1,z2; \n ASSERT( ( i >= 1 and i <= 3 and j >= 1 and j <= 3 and i <> j), \n `second and thi rd argument cannot coincide and should be in the range 1..3`);\n ASSE RT( (simplify(subs( \{ z1 = 1, z2 = 1 \} , A )) = 0 ),`Laurent polynom ial does not have root [1,1]` );\n\n p[1] := z1 - 1; p[2] := z2-1; p[ 3] := z1*z2-1;\n Q[1,2] := subs( z1 = 1, A);\n Q[1,3] := subs( \{z 2 = z1*z2, z1 = 1\}, A); \n Q[2,1] := subs( z2 = 1, A);\n Q[2,3] : = subs( \{z1 = z1*z2, z2 = 1\}, A); \n Q[3,1] := subs( z2 = 1/z1, \+ A);\n Q[3,2] := subs( z1 = 1/z2, A); \n if symmetrize and i = 1 a nd j = 2 then\n TA1 := factor( (1/(2*(z1-1)) )*(A + Q[2,1] - Q[1,2] ));\n TA2 := factor( (1/(2*(z2-1)) )*(A + Q[1,2] - Q[2,1]));\n eli f symmetrize and i = 2 and j = 1 then\n TA2 := factor( (1/(2*(z1-1) ) )*(A + Q[2,1] - Q[1,2]));\n TA1 := factor( (1/(2*(z2-1)) )*(A + Q [1,2] - Q[2,1]));\n else \n TA1 := factor( Q[j,i]/p[i] ) ;\n TA 2 := factor( (A - Q[j,i])/(p[j]) ) ;\n fi;\n ASSERT(IsLaurentPolynom (TA1) and IsLaurentPolynom(TA2), [TA1,TA2]); \n if simplify( p[ i]*TA1 + p[j]*TA2 - A) <> 0 then \n print( `decomposition error`); \n fi;\n [TA1, TA2];\nend: " "6#>%0LinearDecomposeGR6&'%\"AG%/Laur entPolynomG'%\"iG%(integerG'%\"jGF,'%+symmetrizeG%(booleanG7&%$TA1G%$T A2G%\"QG%\"pG6\"F7C1-%'ASSERTG6$33331\"\"\"F+1F+\"\"$1\"\"\"F.1F.\"\"$ 0F+F.%dosecond~and~third~argument~cannot~coincide~and~should~be~in~the ~range~1..3G-F:6$/-%)simplifyG6#-%%subsG6$<$/%#z1G\"\"\"/%#z2G\"\"\"F( \"\"!%LLaurent~polynomial~does~not~have~root~[1,1]G>&F66#\"\"\",&FU\" \"\"\"\"\"!\"\">&F66#\"\"#,&FXF[o\"\"\"F]o>&F66#\"\"$,&*&FUF[oFXF[oF[o \"\"\"F]o>&F56$\"\"\"\"\"#-FQ6$/FU\"\"\"F(>&F56$\"\"\"\"\"$-FQ6$<$/FX* &FUF[oFXF[o/FU\"\"\"F(>&F56$\"\"#\"\"\"-FQ6$/FX\"\"\"F(>&F56$\"\"#\"\" $-FQ6$<$/FU*&FUF[oFXF[o/FX\"\"\"F(>&F56$\"\"$\"\"\"-FQ6$/FX*&\"\"\"F[o FUF]oF(>&F56$\"\"$\"\"#-FQ6$/FU*&\"\"\"F[oFXF]oF(@'33F0/F+\"\"\"/F.\" \"#C$>F3-%'factorG6#*(\"\"\"F[o*&\"\"#F[o,&FUF[o\"\"\"F]oF[oF]o,(F(F[o &F56$\"\"#\"\"\"F[o&F56$\"\"\"\"\"#F]oF[o>F4-Fct6#*(\"\"\"F[o*&\"\"#F[ o,&FXF[o\"\"\"F]oF[oF]o,(F(F[o&F56$\"\"\"\"\"#F[o&F56$\"\"#\"\"\"F]oF[ o33F0/F+\"\"#/F.\"\"\"C$>F4-Fct6#*(\"\"\"F[o*&\"\"#F[o,&FUF[o\"\"\"F]o F[oF]o,(F(F[o&F56$\"\"#\"\"\"F[o&F56$\"\"\"\"\"#F]oF[o>F3-Fct6#*(\"\" \"F[o*&\"\"#F[o,&FXF[o\"\"\"F]oF[oF]o,(F(F[o&F56$\"\"\"\"\"#F[o&F56$\" \"#\"\"\"F]oF[oC$>F3-Fct6#*&&F56$F.F+F[o&F66#F+F]o>F4-Fct6#*&,&F(F[o&F 56$F.F+F]oF[o&F66#F.F]o-F:6$3-%1IsLaurentPolynomG6#F3-Fgz6#F47$F3F4@$0 -FN6#,(*&&F66#F+F[oF3F[oF[o*&&F66#F.F[oF4F[oF[oF(F]oFZ-%&printG6#%4dec omposition~errorG7$F3F4F76$FUFXF7" }}{PARA 12 "" 0 "" {TEXT -1 0 "" }} }}{PARA 0 "" 0 "" {TEXT -1 1 " " }}{SECT 0 {PARA 0 "" 0 "" {MPLTEXT 0 21 48 " FactorizableConvergenceEstimates(A,n,nD,W,WD,B)" }{TEXT -1 13 ". Similar to " }{MPLTEXT 0 21 17 "ComputeEstimates " }{TEXT -1 166 " \+ but with additional assumptions on A and B. A should be a Laurent poly nomial (not matrix) divisible by B; W and WD should have lengths 1 and 2 respectively, and if " }{XPPEDIT 18 0 "W[i]" "6#&%\"WG6#%\"iG" } {TEXT -1 72 " are the Laurent polynomials corresponding to the compon ents of W, and " }{XPPEDIT 18 0 "WD[i]" "6#&%#WDG6#%\"iG" }{TEXT -1 119 " are the Laurent polynomials of the differences corresponding to the components of WD, then B should be divisible by " }{XPPEDIT 18 0 "W[i](z[1]^2,z[2]^2)*WD[i,j](z[1]^2,z[2]^2)/(W[i](z[1],z[2])*WD[i,j] (z[1],z[2]))" "6#*(-&%\"WG6#%\"iG6$*$&%\"zG6#\"\"\"\"\"#*$&F,6#\"\"#\" \"#\"\"\"-&%#WDG6$F(%\"jG6$*$&F,6#\"\"\"\"\"#*$&F,6#\"\"#\"\"#F5*&-&F& 6#F(6$&F,6#\"\"\"&F,6#\"\"#F5-&F86$F(F:6$&F,6#\"\"\"&F,6#\"\"#F5!\"\" " }{TEXT -1 21 " ; for example, if " }{XPPEDIT 18 0 "W[i] = z[1] -1 " "6#/&%\"WG6#%\"iG,&&%\"zG6#\"\"\"\"\"\"\"\"\"!\"\"" }{TEXT -1 7 " a nd " }{XPPEDIT 18 0 "W[i,j] = z[1]*z[2]-1" "6#/&%\"WG6$%\"iG%\"jG,&*& &%\"zG6#\"\"\"\"\"\"&F,6#\"\"#F/F/\"\"\"!\"\"" }{TEXT -1 33 ", then B \+ should be divisible by " }{XPPEDIT 18 0 "(z[1]+1)*(z[1]*z[2]+1)" "6#* &,&&%\"zG6#\"\"\"\"\"\"\"\"\"F)F),&*&&F&6#\"\"\"F)&F&6#\"\"#F)F)\"\"\" F)F)" }{TEXT -1 22 ". The estimates for " }{XPPEDIT 18 0 " gamma" "6 #%&gammaG" }{TEXT -1 6 " and " }{XPPEDIT 18 0 "gamma[D]" "6#&%&gammaG 6#%\"DG" }{TEXT -1 75 " calculated by this function are the same as \+ the estimates calculated by " }{MPLTEXT 0 21 16 "ComputeEstimates" } {TEXT -1 56 "; however, the estimates for c and cD tend to be better. " }}{EXCHG {PARA 0 "> " 0 "" {XPPEDIT 19 1 "FactorizableConvergenceEst imates := proc(A::LaurentPolynom, n::posint, nD::posint, \n W::Intege rVector2List, WD::IntegerVector2List, B::LaurentPolynom) \n local q, pplus, psq, \n a1,a2,a11,a12,a21,a22,b1,b2,b11,b12,b21,b22,\n tq,t q1,tq2, res;\n global z1, z2, gamma, gammaD, c, cD;\n ASSERT( vectdi m(W) = 1, `the length of the 4nd argument (list) should be 1`);\n ASS ERT( vectdim(WD) = 2, `the length of the 5nd argument (list) should be 2`);\n\n q := factor(A/B);\n ASSERT(IsLaurentPolynom(q),`A should \+ be divisible by B`);\n pplus := [z1+1,z2+1,z1*z2+1]; \n psq := [z1^ 2-1,z2^2-1,z1^2*z2^2-1]; \n\n b1 := factor(B/pplus[W[1][1]]);\n b2 \+ := factor(B/pplus[W[1][2]]);\n b11 := factor(b1/pplus[WD[1][1]]); b12 := factor(b1/pplus[WD[1][2]]);\n b21 := factor(b2/pplus[WD[2][1]]); \+ b22 := factor(b2/pplus[WD[2][2]]);\n \n ASSERT( IsLaurentPolynom(b11 ) and IsLaurentPolynom(b12) and\n IsLaurentPolynom(b21) and IsL aurentPolynom(b22),\n `the comparson scheme Laurent polynomial is not divisible by some of the difference polynomials`);\n tq := LinearDe compose(q-1, W[1][1], W[1][2],true); \n tq1 := LinearDecompose(q-1, W D[1][1], WD[1][2],true); \n tq2 := LinearDecompose(q-1, WD[2][1], WD[ 2][2],true); \n a1 := factor(A/pplus[ W[1][1] ]); a2 := factor(A/pplu s[W[1][2]]); \n res[gamma] := max( InfNnorm(matrix([[a1]]),n),InfNno rm(matrix([[a2]]),n));\n\n ASSERT( simplify( A - B - b1*tq[1]*psq[W[1 ][1]] - b2*tq[2]*psq[W[1][2]]) = 0, `bad decomp`); \n res[c] :=InfNno rm( matrix([[b1*tq[1]]]),1) + InfNnorm( matrix([[b2*tq[2]]]),1);\n \n \+ a11 := factor(2*a1/pplus[WD[1][1]]); a12 := factor(2*a1/pplus[WD[1][2 ]]); \n a21 := factor(2*a2/pplus[WD[2][1]]); a22 := factor(2*a2/pplu s[WD[2][2]]); \n res[gammaD] := max( InfNnorm(matrix([[a11]]),nD),In fNnorm(matrix([[a12]]),nD),\n InfNnorm(matrix([[a21]]) ,nD),InfNnorm(matrix([[a22]]),nD));\n \n res[cD] := max( InfNnorm(ma trix([[2*b11*tq1[1]]]),1) + InfNnorm( matrix([[2*b12*tq1[2]]]),1),\n \+ InfNnorm(matrix([[2*b21*tq2[1]]]),1) + InfNnorm( matrix([[2 *b22*tq2[2]]]),1));\n eval(res);\nend:" "6#>%AFactorizableConvergence EstimatesGR6('%\"AG%/LaurentPolynomG'%\"nG%'posintG'%#nDGF,'%\"WG%3Int egerVector2ListG'%#WDGF1'%\"BGF)75%\"qG%&pplusG%$psqG%#a1G%#a2G%$a11G% $a12G%$a21G%$a22G%#b1G%#b2G%$b11G%$b12G%$b21G%$b22G%#tqG%$tq1G%$tq2G%$ resG6\"FJC>-%'ASSERTG6$/-%(vectdimG6#F0\"\"\"%Rthe~length~of~the~4nd~a rgument~(list)~should~be~1G-FM6$/-FQ6#F3\"\"#%Rthe~length~of~the~5nd~a rgument~(list)~should~be~2G>F7-%'factorG6#*&F(\"\"\"F5!\"\"-FM6$-%1IsL aurentPolynomG6#F7%;A~should~be~divisible~by~BG>F87%,&%#z1GF[o\"\"\"F[ o,&%#z2GF[o\"\"\"F[o,&*&FfoF[oFioF[oF[o\"\"\"F[o>F97%,&*$Ffo\"\"#F[o\" \"\"F\\o,&*$Fio\"\"#F[o\"\"\"F\\o,&*&Ffo\"\"#Fio\"\"#F[o\"\"\"F\\o>F@- Fhn6#*&F5F[o&F86#&&F06#\"\"\"6#\"\"\"F\\o>FA-Fhn6#*&F5F[o&F86#&&F06#\" \"\"6#\"\"#F\\o>FB-Fhn6#*&F@F[o&F86#&&F36#\"\"\"6#\"\"\"F\\o>FC-Fhn6#* &F@F[o&F86#&&F36#\"\"\"6#\"\"#F\\o>FD-Fhn6#*&FAF[o&F86#&&F36#\"\"#6#\" \"\"F\\o>FE-Fhn6#*&FAF[o&F86#&&F36#\"\"#6#\"\"#F\\o-FM6$333-F`o6#FB-F` o6#FC-F`o6#FD-F`o6#FE%ipthe~comparson~scheme~Laurent~polynomial~is~not ~divisible~by~some~of~the~difference~polynomialsG>FF-%0LinearDecompose G6&,&F7F[o\"\"\"F\\o&&F06#\"\"\"6#\"\"\"&&F06#\"\"\"6#\"\"#%%trueG>FG- Fev6&,&F7F[o\"\"\"F\\o&&F36#\"\"\"6#\"\"\"&&F36#\"\"\"6#\"\"#Few>FH-Fe v6&,&F7F[o\"\"\"F\\o&&F36#\"\"#6#\"\"\"&&F36#\"\"#6#\"\"#Few>F:-Fhn6#* &F(F[o&F86#&&F06#\"\"\"6#\"\"\"F\\o>F;-Fhn6#*&F(F[o&F86#&&F06#\"\"\"6# \"\"#F\\o>&FI6#%&gammaG-%$maxG6$-%)InfNnormG6$-%'matrixG6#7#7#F:F+-Fh[ l6$-F[\\l6#7#7#F;F+-FM6$/-%)simplifyG6#,*F(F[oF5F\\o*(F@F[o&FF6#\"\"\" F[o&F96#&&F06#\"\"\"6#\"\"\"F[oF\\o*(FAF[o&FF6#\"\"#F[o&F96#&&F06#\"\" \"6#\"\"#F[oF\\o\"\"!%+bad~decompG>&FI6#%\"cG,&-Fh[l6$-F[\\l6#7#7#*&F@ F[o&FF6#\"\"\"F[o\"\"\"F[o-Fh[l6$-F[\\l6#7#7#*&FAF[o&FF6#\"\"#F[o\"\" \"F[o>F<-Fhn6#*(\"\"#F[oF:F[o&F86#&&F36#\"\"\"6#\"\"\"F\\o>F=-Fhn6#*( \"\"#F[oF:F[o&F86#&&F36#\"\"\"6#\"\"#F\\o>F>-Fhn6#*(\"\"#F[oF;F[o&F86# &&F36#\"\"#6#\"\"\"F\\o>F?-Fhn6#*(\"\"#F[oF;F[o&F86#&&F36#\"\"#6#\"\"# F\\o>&FI6#%'gammaDG-Fe[l6&-Fh[l6$-F[\\l6#7#7#FF.-Fh[l6$-F[\\l6#7#7#F?F.>&FI6#%#cDG-Fe[l6$,&-F h[l6$-F[\\l6#7#7#*(\"\"#F[oFBF[o&FG6#\"\"\"F[o\"\"\"F[o-Fh[l6$-F[\\l6# 7#7#*(\"\"#F[oFCF[o&FG6#\"\"#F[o\"\"\"F[o,&-Fh[l6$-F[\\l6#7#7#*(\"\"#F [oFDF[o&FH6#\"\"\"F[o\"\"\"F[o-Fh[l6$-F[\\l6#7#7#*(\"\"#F[oFEF[o&FH6# \"\"#F[o\"\"\"F[o-%%evalG6#FIFJ6(FfoFioFc[lFhclFi^lFfelFJ" }}}}}{SECT 0 {PARA 3 "" 0 "" {TEXT -1 28 "Results for specific schemes" }}{SECT 0 {PARA 4 "" 0 "" {TEXT -1 22 "Bilinear and Trilinear" }}{PARA 5 "" 0 "" {TEXT -1 0 "" }}{EXCHG {PARA 0 "" 0 "" {TEXT -1 371 "This is just a sanity check; note that because we use difference schemes of a linear scheme as comparison scheme for derivatives,\nwhich is not continuous , we cannot infer C1-continuity from the fact that cD = 0. In general , we assume that C1-continuity is known, and we \nare interested in th e rate of approximation. To establish C1-continuity, one has to look \+ only at " }{XPPEDIT 18 0 "gamma[D]" "6#&%&gammaG6#%\"DG" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 154 "Linear := (1+z)^2/2: Bilinear := evalm(diag( 1)*(1/4)*(1 + z1^(-1))^2*(1 + z2^(-1))^2*z1*z2):\nTrilinear := (1/2)*( 1+z1)*(1+z2)*(1+z1*z2)*z1^(-1)*z2^(-1): \n" }}}{EXCHG {PARA 0 "> " 0 " " {MPLTEXT 1 0 80 "op(ConvergenceEstimates(Bilinear, 1,1,[[1,2]], [[1, 2],[2,1]] ,Bilinear[1,1])); " }}{PARA 11 "" 1 "" {XPPMATH 20 "6#7&/%& gammaG#\"\"\"\"\"#/%'gammaDGF'/%\"cG\"\"!/%#cDGF-" }}}{EXCHG {PARA 0 " > " 0 "" {MPLTEXT 1 0 96 "op(FactorizableConvergenceEstimates(Bilinear [1,1],1,1,[[1,2]], [[1,2],[2,1]] ,Bilinear[1,1])); " }}{PARA 11 "" 1 "" {XPPMATH 20 "6#7&/%&gammaG#\"\"\"\"\"#/%#cDG\"\"!/%'gammaDGF'/%\"cG F+" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 139 "TrilinearEstimate : = ConvergenceEstimates(evalm( Trilinear * diag(1)), 1,1,[[1,2]], [[2,3 ],[3,1]] ,Trilinear) : op(op(TrilinearEstimate));" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#7&/%&gammaG#\"\"\"\"\"#/%#cDG\"\"!/%'gammaDGF'/%\"cGF+ " }}}}{SECT 0 {PARA 4 "" 0 "" {TEXT -1 31 "3-directional box spline (L oop)" }}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 278 "Trilinear := (1/2)* (1+z1)*(1+z2)*(1+z1*z2)*z1^(-1)*z2^(-1): \nThreeDirCoeffs := array( -2 ..2,-2..2, [\n [ 1/16, 1/8, 1/16, 0, 0 ],\n [ 1/8, 3/8, 3/8, 1/8, 0 ],\n [ 1/16, 3/8, 5/8, 3/8, 1/16 ],\n [ 0, 1/8, 3/ 8, 3/8, 1/8 ],\n [ 0, 0, 1/16, 1/8, 1/16 ]\n]):" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 65 "ThreeDir := evalm( diag(1)* factor( MakePolynom(ThreeDirCoeffs)));" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#>%)T hreeDirG-%'matrixG6#7#7#,$*&*(),&\"\"\"F/%#z2GF/\"\"#\"\"\"),&F/F/%#z1 GF/F1F2),&F/F/*&F5F/F0F/F/F1F2F2*&)F5\"\"#F2)F0\"\"#F2!\"\"#F/\"#;" }} }{EXCHG {PARA 0 "" 0 "" {TEXT -1 244 "Because the Laurent polynomial \+ has all three difference factors, and each is squared, we can use a va riety of contraction functions.\nAs it is typically the case, the spec ial function for factorizable polynomials gives better estimates for \+ cD." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 74 "op(ConvergenceEstimates(ThreeDir, 1,1,[[1,2]], [[3,2],[1,3]] ,Tril inear));" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#7&/%&gammaG#\"\"\"\"\"#/%# cDGF'/%'gammaDGF&/%\"cGF&" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 78 "op(ConvergenceEstimates(ThreeDir, 1,1,[[1,2]], [[1,2],[2,1]] ,Bili near[1,1]));" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#7&/%&gammaG#\"\"\"\"\" #/%#cDG#\"\"&\"\")/%'gammaDGF&/%\"cGF&" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 78 "op(ConvergenceEstimates(ThreeDir, 1,1,[[1,2]], [[1,3] ,[1,3]] ,Bilinear[1,1]));" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#7&/%&gamm aG#\"\"\"\"\"#/%#cDG#\"\"$\"\"%/%'gammaDGF&/%\"cGF&" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 91 "op(FactorizableConvergenceEstimates(Three Dir[1,1], 1,1,[[1,2]], [[3,2],[1,3]] ,Trilinear));" }}{PARA 11 "" 1 " " {XPPMATH 20 "6#7&/%&gammaG#\"\"\"\"\"#/%#cDGF&/%'gammaDGF&/%\"cGF&" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 137 "ThreeDirEstimate := Fact orizableConvergenceEstimates(ThreeDir[1,1], 1,1,[[1,2]], [[1,2],[1,2]] ,Bilinear[1,1]): op(op(ThreeDirEstimate));" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#7&/%&gammaG#\"\"\"\"\"#/%#cDGF&/%'gammaDGF&/%\"cGF&" }} }}{SECT 0 {PARA 4 "" 0 "" {TEXT -1 9 "Butterfly" }}{EXCHG {PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 282 "ButterflyCoeffs := proc(w) array( -3..3,-3..3 , [\n [ 0, -w, -w, 0, 0, 0, 0 ],\n [ -w, 0,2*w, 0, -w, 0, 0 ] ,\n [ -w,2*w,1/2,1/2,2*w, -w, 0 ],\n [ 0, 0,1/2, 1,1/2, 0, 0 ], \n [ 0, -w,2*w,1/2,1/2,2*w, -w ],\n [ 0, 0, -w, 0,2*w, 0, -w ],\n [ 0, 0, 0, 0, -w, -w, 0 ]\n]);" }}{PARA 0 "" 0 "" {TEXT -1 1 " \010" }{MPLTEXT 1 0 4 "end:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 62 "Butterfly := evalm( diag(1)* MakePolynom(ButterflyCoeffs(w))):" }} }{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 47 "Butterfly := subs( w = 1/16 , evalm(Butterfly));" }}{PARA 12 "" 1 "" {XPPMATH 20 "6#>%*ButterflyG- %'matrixG6#7#7#,T\"\"\"F+*&)%#z1G\"\"$\"\"\"%#z2GF+#!\"\"\"#;F.#F+\"\" #F1F5*&F.F0*$)F1\"\"#F0!\"\"F2*&*$)F.F6F0F0F1F;F2*&F.F+F1F0F5*&F0F0F.F ;F5*&F0F0F1F;F5*&F.F0)F1F6F0#F+\"\")*&F>F0F1F0FD*&F1F0F.F;FD*&F0F0*&F. \"\"\")F1\"\"#F0F;FD*&F0F0*&)F.\"\"#F0F1\"\"\"F;FD*&F0F0*&F.\"\"\"F1\" \"\"F;F5*&F.F0F1F;FD*&F>F0)F1F/F0F2*&F-F0FCF0F2*&F0F0*&)F.\"\"$F0F1\" \"\"F;F2*&F0F0*&)F.\"\"$F0)F1\"\"#F0F;F2*&*$FCF0F0F.F;F2*&F0F0*&)F.\" \"#F0)F1\"\"$F0F;F2*&F1F0*$)F.\"\"#F0F;F2*&F.F0FXF0F2*&F0F0*&F.\"\"\") F1\"\"$F0F;F2" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 123 "ButterflyEstimate1 := Conve rgenceEstimates(Butterfly, 1,1,[[1,2]], [[2,3],[3,1]] ,Trilinear): op( eval(ButterflyEstimate1));" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#7&/%&gam maG#\"\"(\"\")/%#cDG\"\"\"/%'gammaDGF+/%\"cG#F+\"\"#" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 60 "map(simplify, op(subs(w = 1/16, eval(Butt erflyEstimate1))));" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#7&/%&gammaG#\" \"(\"\")/%#cDG\"\"\"/%'gammaDGF+/%\"cG#F+\"\"#" }}}{EXCHG {PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 93 "ButterflyEstimate2 := ConvergenceEstimates(But terfly, 2,2,[[1,2]], [[2,3],[3,1]] ,Trilinear):" }}}{EXCHG {PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 109 "assume(W > 0); additionally(W < 1/8); \nsolve ( simplify( subs( w = W, eval(ButterflyEstimate2[gammaD]))) < 1);" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 57 "op(map( eval, subs(w = 1/16, eval(ButterflyEstimate2))));" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#7&/%& gammaG#\"#J\"#k/%#cDG\"\"\"/%'gammaDG#\"\"(\"\")/%\"cG#F+\"\"#" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 93 "ButterflyEstimate3 := Conver genceEstimates(Butterfly, 3,3,[[1,2]], [[2,3],[3,1]] ,Trilinear):" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 100 "assume(W > 1/12); additiona lly(W < 1/8); \nsimplify( subs( w = W, eval(ButterflyEstimate3[gammaD] )));" }}{PARA 11 "" 1 "" {XPPMATH 20 "6##\"#6\"#;" }}}{EXCHG {PARA 0 " > " 0 "" {MPLTEXT 1 0 57 "op(map( eval, subs(w = 1/16, eval(ButterflyE stimate3))));" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#7&/%&gammaG#\"$h#\"%C 5/%#cDG\"\"\"/%'gammaDG#\"#6\"#;/%\"cG#F+\"\"#" }}}{EXCHG {PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 321 "LevinCoeffs := proc(w) array( -3..3,-3..3, [ \n [ 0, 0, -w/2, -w, -w/2, 0, 0 ],\n [ 0, 0, 0, 0, 0, 0, 0 ],\n [ -w/2, 0, 5*w, 9*w, 5*w, 0, -w/2],\n [ -w, 0, 9*w, 1, 9*w, 0, -w],\n [ -w/2, 0, 5*w, 9*w, 5*w, 0, -w/2],\n [ 0, \+ 0, 0, 0, 0, 0, 0 ],\n [ 0, 0, -w/2, -w, -w/2, 0, 0 ] \n]);" }}{PARA 0 "" 0 "" {TEXT -1 1 "\010" }{MPLTEXT 1 0 4 "end:" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 54 "Levin := evalm( diag(1)* Mak ePolynom(LevinCoeffs(w))):" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 39 "Levin := subs( w = 1/16, evalm(Levin));" }}{PARA 12 "" 1 "" {XPPMATH 20 "6#>%&LevinG-%'matrixG6#7#7#,L*&\"\"\"F,*$)%#z1G\"\"$F,!\" \"#!\"\"\"#;*&F,F,*&F/\"\"\"%#z2G\"\"\"F1#\"\"&F4F8#\"\"*F4*&F,F,F/F1F <*&F/\"\"\"F8F@F:*&F8F,*$)F/\"\"$F,F1#F3\"#K*&F,F,F8F1F<*&*$)F/\"\"$F, F,F8F1FE*&F,F,*$)F8\"\"$F,F1F2*$FJF,F2*&F/F,*$)F8\"\"$F,F1FEF@F@*&FJF, F8F,FE*&F/F,)F8FKF,FE*&F/F,F8F1F:*&F,F,*&F/\"\"\")F8\"\"$F,F1FEF/F<*&* $FWF,F,F/F1FE*&F,F,*&)F/\"\"$F,F8\"\"\"F1FE*$FWF,F2*&F8F,F/F1F:" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 19 "factor(Levin[1,1]);" }} {PARA 11 "" 1 "" {XPPMATH 20 "6#,$*&*(),&\"\"\"F(%#z2GF(\"\"#\"\"\"),& F(F(%#z1GF(F*F+,4*&)F.\"\"%F+)F)F*F+F(*&)F.\"\"$F+F3F+!\"#*&)F)F2F+)F. F*F+F(*$F:F+F(*&F)F(F:F+F7*&F:F+)F)F6F+F7*&F:F+F3F+!\"%*&F3F+F.F(F7*$F 3F+F(F(F+*&)F.\"\"$F+)F)\"\"$F+!\"\"#!\"\"\"#K" }}}{EXCHG {PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 121 "cfs := [coeffs( -(z1^4*z2^2-2*z1^3*z2^2+z2^4* z1^2+z1^2-2*z2*z1^2-2*z1^2*z2^3-4*z1^2*z2^2-2*z2^2*z1+z2^2), [z1,z2], \+ 't')];" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#>%$cfsG7+\"\"%!\"\"F'\"\"#F( F(F'F'F(" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 9 "t := [t];" }} {PARA 11 "" 1 "" {XPPMATH 20 "6#>%\"tG7+*&)%#z1G\"\"#\"\"\")%#z2GF)F** $F'F**$F+F**&F,\"\"\"F'F**&F+F*F(F0*&)F(\"\"$F*F+F**&)F(\"\"%F*F+F**&) F,F7F*F'F**&F'F*)F,F4F*" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 " " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 19 "M := matrix(5,5,0):" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 85 "for j from 1 to vectdim(cfs) do M[d egree(t[j],z1)+1, degree(t[j],z2)+1] := cfs[j];od:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 9 "evalm(M);" }}{PARA 11 "" 1 "" {XPPMATH 20 "6 #-%'matrixG6#7'7'\"\"!F(!\"\"F(F(7'F(F(\"\"#F(F(7'F)F+\"\"%F+F)F*F'" } }}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 116 "LevinEstimate1 := Converg enceEstimates(Levin, 1,1,[[1,2]], [[1,2],[2,1]] ,Bilinear[1,1]): op(e val(LevinEstimate1));" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#7&/%&gammaG# \"\"$\"\"%/%'gammaDG#\"\"&F(/%\"cG#F'\"\")/%#cDG#\"\"(\"#;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 91 "LevinEstimate2 := ConvergenceEstimates(Levin, 2,2,[[ 1,2]], [[1,2],[2,1]] ,Bilinear[1,1]): " }}{PARA 12 "" 1 "" {XPPMATH 20 "6#>%/LevinEstimate2G-%&TABLEG6#7&/%&gammaG#\"#D\"#k/%'gammaDG#\"#: \"#;/%\"cG#\"\"$\"\")/%#cDG#\"\"(F2" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 25 "op(eval(LevinEstimate2));" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#7&/%&gammaG#\"#D\"#k/%'gammaDG#\"#:\"#;/%\"cG#\"\"$\"\" )/%#cDG#\"\"(F-" }}}{EXCHG {PARA 11 "" 1 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 91 "LevinEstimate3 := ConvergenceEstima tes(Levin, 3,3,[[1,2]], [[1,2],[2,1]] ,Bilinear[1,1]): " }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 23 "op(op(LevinEstimate3));" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#7&/%&gammaG#\"$0\"\"$7&/%'gammaDG#\"%48\"%[? /%\"cG#\"\"$\"\")/%#cDG#\"\"(\"#;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 56 "[gamma = 425/2048, cD = 31/64, gammaD = 5/8, c = 13/3 2];" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 17 "425/2048-105/512;" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#7&/%&gammaG#\"$D%\"%[?/%#cDG#\"#J\"#k/%'gamm aDG#\"\"&\"\")/%\"cG#\"#8\"#K" }}{PARA 11 "" 1 "" {XPPMATH 20 "6##\"\" &\"%[?" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}{PARA 11 "" 1 "" {XPPMATH 20 "6##\"\"$ \"#k" }}}{PARA 0 "" 0 "" {TEXT -1 0 "" }}}{SECT 0 {PARA 4 "" 0 "" {TEXT -1 45 "Tensor product of quadric splines (Doo-Sabin)" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 127 "Qua dric := expand((1/4)*(z^(-1) + 1)^3*z); QuadricTensorSpline := evalm(d iag(1)*(1/16)*(z1^(-1) + 1)^3*(z2^(-1) + 1)^3*z1*z2 );" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#>%(QuadricG,**&\"\"\"F'*$)%\"zG\"\"#F'!\"\"#\"\"\" \"\"%*&F'F'F*F,#\"\"$F/F1F.F*F-" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#>%4 QuadricTensorSplineG-%'matrixG6#7#7#,$**),&\"\"\"F.*&\"\"\"F0%#z1G!\" \"F.\"\"$F0),&F.F.*&F0F0%#z2GF2F.F3F0F1F.F7F.#F.\"#;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 89 "op(ConvergenceEstimates(QuadricTensorSpli ne, 1,1,[[1,2]], [[1,2],[2,1]] ,Bilinear[1,1]));" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#7&/%&gammaG#\"\"\"\"\"#/%#cDG#\"\"$\"\"%/%'gammaDGF&/% \"cGF&" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 158 "QuadricTensorEst imate := FactorizableConvergenceEstimates(QuadricTensorSpline[1,1], 1, 1,[[1,2]], [[1,2],[2,1]] ,Bilinear[1,1]): op(op(QuadricTensorEstimate) );" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#7&/%&gammaG#\"\"\"\"\"#/%#cDG#\" \"$\"\"%/%'gammaDGF&/%\"cGF&" }}}}{SECT 0 {PARA 4 "" 0 "" {TEXT -1 46 "Tensor product of cubic splines (Camull-Clark)" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 119 "Cubic := (1 /8)*(1+z^(-1))^4*z^2; TensorCubicSpline := evalm(diag(1)*(1/64)*(z1^(- 1) + 1)^4*(z2^(-1) + 1)^4*z1^2*z2^2 );" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#>%&CubicG,$*&),&\"\"\"F)*&\"\"\"F+%\"zG!\"\"F)\"\"%F+)F,\"\"#F+# F)\"\")" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#>%2TensorCubicSplineG-%'mat rixG6#7#7#,$**),&\"\"\"F.*&\"\"\"F0%#z1G!\"\"F.\"\"%F0),&F.F.*&F0F0%#z 2GF2F.F3F0)F1\"\"#F0)F7F9F0#F.\"#k" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 87 "op(ConvergenceEstimates(TensorCubicSpline, 1,1,[[1,2] ], [[2,1],[1,2]] ,Bilinear[1,1]));" }}{PARA 11 "" 1 "" {XPPMATH 20 "6# 7&/%&gammaG#\"\"\"\"\"#/%#cDG#\"#6\"#;/%'gammaDGF&/%\"cGF&" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 152 "TensorCubicEstimate:= Factorizable ConvergenceEstimates(TensorCubicSpline[1,1], 1,1,[[1,2]], [[2,1],[1,2] ] ,Bilinear[1,1]): op(op(TensorCubicEstimate)); " }}{PARA 11 "" 1 "" {XPPMATH 20 "6#7&/%&gammaG#\"\"\"\"\"#/%#cDGF&/%'gammaDGF&/%\"cGF&" }} }{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}}{SECT 0 {PARA 4 "" 0 " " {TEXT -1 34 "Tensor product of quartic splines " }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 134 "Quartic := \+ expand((1/4)*(z^(-1) + 1)^5*z^2); QuarticTensorSpline := evalm(diag(1) *(1/256)*(z1^(-1) + 1)^5*(z2^(-1) + 1)^5*z1^2*z2^2 );" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#>%(QuarticG,.*$)%\"zG\"\"#\"\"\"#\"\"\"\"\"%F(#\" \"&F-#F/F)F,*&F*F*F(!\"\"F0*&F*F**$)F(\"\"#F*F2F.*&F*F**$)F(\"\"$F*F2F +" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#>%4QuarticTensorSplineG-%'matrixG 6#7#7#,$**),&\"\"\"F.*&\"\"\"F0%#z1G!\"\"F.\"\"&F0),&F.F.*&F0F0%#z2GF2 F.F3F0)F1\"\"#F0)F7F9F0#F.\"$c#" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 89 "op(ConvergenceEstimates(QuarticTensorSpline, 1,1,[[1,2]], [[1, 2],[2,1]] ,Bilinear[1,1]));" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#7&/%&ga mmaG#\"\"\"\"\"#/%#cDG#\"#h\"#k/%'gammaDGF&/%\"cG#\"\"$\"\"%" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 158 "QuarticTensorEstimate := Fa ctorizableConvergenceEstimates(QuarticTensorSpline[1,1], 1,1,[[1,2]], \+ [[1,2],[2,1]] ,Bilinear[1,1]): op(op(QuarticTensorEstimate));" }} {PARA 11 "" 1 "" {XPPMATH 20 "6#7&/%&gammaG#\"\"\"\"\"#/%#cDG#\"\"(\" \")/%'gammaDGF&/%\"cG#\"\"$\"\"%" }}}}{SECT 0 {PARA 4 "" 0 "" {TEXT -1 43 "Tensor product of 4-point schemes (Kobbelt)" }}{EXCHG {PARA 0 " > " 0 "" {MPLTEXT 1 0 223 "FourPt := factor((-w*z^(-3) + (1/2+w)*z^(-1 ) + 1 + (1/2+w)*z - w*z^3)); TensorFourPt := matrix([[factor((-w*z1^(- 3) + (1/2+w)*z1^(-1) + 1 + (1/2+w)*z1 - w*z1^3)*( -w*z2^(-3) + (1/2+w) *z2^(-1) + 1 + (1/2+w)*z2 -w*z2^3 ))]]);" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#>%'FourPtG,$*&*&),&\"\"\"F*%\"zGF*\"\"#\"\"\",.*&)F+\"\"%F-%\"wG F*F,*&F2F-)F+\"\"$F-!\"%*$)F+F,F-!\"\"*&F8F-F2F-F1*&F2F-F+F*F6F2F,F*F- *$)F+\"\"$F-!\"\"#F9F," }}{PARA 12 "" 1 "" {XPPMATH 20 "6#>%-TensorFou rPtG-%'matrixG6#7#7#,$*&**),&\"\"\"F/%#z1GF/\"\"#\"\"\",.*&)F0\"\"%F2% \"wGF/F1*&F7F2)F0\"\"$F2!\"%*$)F0F1F2!\"\"*&F=F2F7F2F6*&F7F2F0F/F;F7F1 F/),&F/F/%#z2GF/F1F2,.*&)FCF6F2F7F2F1*&F7F2)FCF:F2F;*$)FCF1F2F>*&FJF2F 7F2F6*&F7F2FCF/F;F7F1F/F2*&)FC\"\"$F2)F0\"\"$F2!\"\"#F/F6" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 54 "TensorFourPt := subs( w = 1/16, eva lm( TensorFourPt));" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#>%-TensorFourPt G-%'matrixG6#7#7#,$*&**),&\"\"\"F/%#z1GF/\"\"#\"\"\",,*$)F0\"\"%F2#F/ \"\")*$)F0\"\"$F2#!\"\"F6*$)F0F1F2#!\"$F6F0F " 0 "" {MPLTEXT 1 0 103 "TensorFourPtE stimate1 := ConvergenceEstimates(TensorFourPt, 1,1,[[1,2]], [[1,2],[2, 1]] ,Bilinear[1,1]);" }}{PARA 12 "" 1 "" {XPPMATH 20 "6#>%6TensorFourP tEstimate1G-%&TABLEG6#7&/%&gammaG#\"#D\"#K/%#cDG#\"#J\"#k/%'gammaDG#\" \"&\"\"%/%\"cG#\"#8F-" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 62 "op ( map( eval, subs( w = 1/16, eval(TensorFourPtEstimate1))));" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#7&/%&gammaG#\"#D\"#K/%#cDG#\"#J\"#k/%'gammaD G#\"\"&\"\"%/%\"cG#\"#8F(" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 103 "TensorFourPtEstimate2 := ConvergenceEstimates(TensorFourPt, 2,2,[ [1,2]], [[1,2],[2,1]] ,Bilinear[1,1]):" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 61 "op(map( eval, subs( w = 1/16, eval(TensorFourPtEstima te2))));" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#7&/%&gammaG#\"$0\"\"$c#/%# cDG#\"#J\"#k/%'gammaDG#\"#:\"#;/%\"cG#\"#8\"#K" }}}{EXCHG {PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 165 "TensorFourPtEstimate3 := ConvergenceEstimates (TensorFourPt, 3,3,[[1,2]], [[1,2],[2,1]],Bilinear[1,1]): \nop(map( ev al, subs( w = 1/16, eval(TensorFourPtEstimate3))));" }}{PARA 11 "" 1 " " {XPPMATH 20 "6#7&/%&gammaG#\"$D%\"%[?/%#cDG#\"#J\"#k/%'gammaDG#\"\"& \"\")/%\"cG#\"#8\"#K" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 125 "Te nsorFourPtEstimate1f := op(FactorizableConvergenceEstimates(TensorFour Pt[1,1], 1,1,[[1,2]], [[2,1],[1,2]] ,Bilinear[1,1]));" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#>%7TensorFourPtEstimate1fG7&/%&gammaG#\"#D\"#K/%#c DG#\"#F\"#k/%'gammaDG#\"\"&\"\"%/%\"cG#\"\"*F*" }}}{EXCHG {PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 62 "op(map( eval, subs( w = 1/16, eval(TensorFourP tEstimate1f))));" }}{PARA 11 "" 1 "" {XPPMATH 20 "6&/%&gammaG#\"#D\"#K /%#cDG#\"#F\"#k/%'gammaDG#\"\"&\"\"%/%\"cG#\"\"*F'" }}}}}{SECT 1 {PARA 3 "" 0 "" {TEXT -1 15 "Code generation" }}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 81 "currentdir(\"D:\\\\users\\\\dzorin\\\\maple\"); re adlib(C); read `subdivmatrix-util.mpl`;" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#Q6D:\\users\\dzorin\\maple6\"" }}{PARA 11 "" 1 "" {XPPMATH 20 "6 #R6\"60%%locsG%%parsG%%glosG%%declG%\"kG%\"xG%&errorG%+old_DigitsG%*ol d_QuietG%,output_fileG%(statseqG%$optG%\"AG%(unknownG6#%aoCopyright~(c )~1992~by~the~University~of~Waterloo.~All~rights~reserved.GF$C9@$3-%%t ypeG6$&9\"6#\"\"\"%&arrayG4-F:6$F<%%nameGC$>81F<-%'RETURNG6#-%\"CG6$FG &F=6#;\"\"#9#>8+%'DigitsG>FU\"#;@%4-%)assignedG6#%*precisionG>%,C/prec isionG.%'doubleG>FinFgn@$4-Fen6#&%*infolevelG6#.FL>&Fbo6#FL\"\"!>%+C/f ilenameG%)terminalG>8$%%NULLG>8&F^p>8%F^p>8'F^p>%'C/ansiG%&falseG>8/Fg p>8-.F[p?&8(7#FN%%trueG@7/F^q.%*optimizedG>FipF`q/F^q.%%ansiG>FfpF`q/F ^q.%&xansiG>FfpFgp-F:6$F^q./-%*identicalG6#%'digitsG%'posintG>FU-%$rhs G6#F^q-F:6$F^q./-Fcr6#%)filenameG<$%'symbolG%'stringGC$>F[qFhr>FjoFhr- F:6$F^q./-FcrFfnFcs>Fin-%#opG6$FQF^q-F:6$F^q./-Fcr6#%(globalsG-%%listG 6#Fcs>F`pF^q-F:6$F^q./-Fcr6#%'localsGFht>F]pF^q-F:6$F^q./-Fcr6#%+param etersGFht>FbpF^q-F:6$F^q./-Fcr6#%-declarationsG%)anythingGC%>FdpFhr@%- F:6$Fdp<%%#::G-.Fit6#Fjv-.%$setGF]wF$-%&ERRORG6$%1bad~declarationsGFdp >FdpF^q-Fbw6$%+bad~optionGF^q@'/Fin.%'singleG>FU\"\"(/FinFjn>FUFW-Fbw6 $%:bad~setting~for~precisionGFgn@%0F[qF\\q-%*traperrorG6#-%&fopenG6$F[ q.%'APPENDG>F[q.F[q@'3FB-F:6$F<.F@C%>8)--%(readlibG6#%-tools/gensymG6# F<-%'assignG6$Fgy-%/C/arrayreindexGF]z@%Fip>8.7#--Fjy6#.%)optimizeG6#F gy>Fez7#-%.C/arrayexpandGF\\[l3-F:6$F<.Fcs-F:6$F<.%*procedureGC(-%)use rinfoG6%FQFL%Bgenerating~output~for~a~procedureG@%FipC$>Fip-FhzF]z>8*- Fgx6#--Fjy6#%,C/procedureG6%FipFFc\\l-Fgx6#-Fg\\l6$FFUFT@%/Fc\\l%*lasterrorG-Fbw6#Fc\\l-FI 6#F^pC$>FgyF<@)Fip>Fez7#-Fhz6&FgyF`pF]pFbp-F:6$Fgy-Fit6#/FD%*algebraic G>FezFgy-F:6$FgyF[_l>Fez7#/%#t0GFgy-Fbw6$%1invalid~argumentGFgy>Fez-%2 C/ConvertAndCheckG6#Fez@$Fc]l>FjoF[q>Fc\\l-Fgx6#-%6C/ComputationSequen ceGFi_l@$Fc]lFe]l>FUFT@%Fi]lF[^lF^pF$6(%)C/bufferGFfpFgnFinFjoFboF$" } }}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 340 "OutputFile := \"rschemes. cpp\":\nfprintf(OutputFile , `// This file is automatically generated; do not change it!\\n`):\nfprintf(OutputFile , `// Change the Maple wo rksheet polynomial.mws\\n`):\nfprintf(OutputFile , `#include \"regsche me.h\"\\n`):\nfprintf(OutputFile , `typedef Float* pFloat;\\n\\n`):\nf printf(OutputFile , `#define BV BadValue;\\n\\n`):\n" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 143 "GenerateCode( OutputFile, \"Trilinear\", [ eval(TrilinearEstimate)], map( unapply( x[1,1], x), CoeffList(evalm (Trilinear * diag(1)))), 0, true); " }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 127 "GenerateCode( OutputFile, \"TriBox\",[ eval(ThreeDir Estimate)], map( unapply( x[1,1], x), CoeffList(eval(ThreeDir))), 1, \+ true); " }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 141 "GenerateCode( O utputFile, \"Bicubic\",[ eval(TensorCubicEstimate)], map( unapply( x[ 1,1], x), CoeffList(eval(TensorCubicSpline))), 1, false); " }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 224 "GenerateCode( OutputFile, \"RegBut terfly\", subs( w = 1/16, [ eval(ButterflyEstimate1), eval(ButterflyEs timate2), eval(ButterflyEstimate3)]), map( unapply( x[1,1], x), Coeff List(subs( w = 1/16, eval(Butterfly) ))), 2, true); " }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 237 "GenerateCode( OutputFile, \"TensorFourPt \", subs( w = 1/16, [ eval(TensorFourPtEstimate1), eval(TensorFourPtEs timate2), eval(TensorFourPtEstimate3)]), map( unapply( x[1,1], x), Co effList(subs( w = 1/16, eval(TensorFourPt) ))), 2, false); " }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 147 "GenerateCode( OutputFile, \+ \"Biquadric\",[ eval(QuadricTensorEstimate)], map( unapply( x[1,1], x ), CoeffList(eval(QuadricTensorSpline))), 1, false); " }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 19 "fclose(OutputFile);" }}}}{PARA 3 "" 0 "" {TEXT -1 0 " " }}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}}{MARK "6 7 1 2 0" 0 }{VIEWOPTS 0 0 0 1 1 1803 }