{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 Output" 2 20 "" 0 1 0 0 255 1 0 0 0 0 0 0 0 0 0 } {CSTYLE "" -1 256 "" 0 1 0 0 0 0 0 1 0 0 0 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 "Text Output" -1 2 1 {CSTYLE "" -1 -1 "Co urier" 1 10 0 0 255 1 0 0 0 0 0 1 3 0 3 }1 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 8 4 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 8 2 0 0 0 0 0 0 -1 0 }{PSTYLE "" 2 6 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 2 0 0 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "Warning" 2 7 1 {CSTYLE "" -1 -1 "" 0 1 0 0 255 1 0 0 0 0 0 0 1 0 0 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "Maple Output" 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 "" 4 256 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }} {SECT 0 {EXCHG {PARA 256 "" 0 "" {TEXT 256 28 "Exporting Cubicspline t o C++" }}{PARA 0 "" 0 "" {TEXT -1 139 "We want to show one way to expo rt the cubic spline generating algorithm that we have coded in Maple t o C++. Here is the Maple version. " }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 820 "cubicspline := proc(x,y,z1, zn)\n local i,dx,dy,z,d,r,n,a,b,c,t,s,tmp;\n n := nops(x);\n z[1] := z1: z[n]:=zn:\n for i from 1 to n-1 do \n dx[i]:=x[i+1]-x[i]:\n \+ dy[i]:=y[i+1]-y[i] od:\n for i from 1 to n-2 do d[i]:= 2*(dx[i]+dx[ i+1]):\n r[i]:= 6*(dy[i+1]/dx[i+1]-dy[i]/dx[i]) od:\n r[1] := r[ 1]-dx[1]*z[1]: \n for i from 2 to n-2 do \n d[i] := d[i]-dx[i+1]^2/ d[i-1]: \n r[i] := r[i]-dx[i]*r[i-1]/d[i-1] \n od:\n for i from n-1 by -1 to 2 do\n z[i] := (r[i-1]-dx[i]*z[i+1])/d[i-1] od:\n for \+ i from 1 to n-1 do\n a[i] := 1/6*(z[i+1]-z[i])/dx[i]: \n b[i] := z[ i]/2:\n c[i] := dy[i]/dx[i]-a[i]*dx[i]^2-b[i]*dx[i] : \n s[i] := a[ i]*(t-x[i])^3+b[i]*(t-x[i])^2+c[i]*(t-x[i])+y[i]:\n s[i] := unapply(s [i],t) od:\n readlib(piecewise);\n unapply(piecewise(seq(op([t " 0 "" {MPLTEXT 1 0 749 "cubicspline := proc(x, y,z1,zn)\n local i,dx,dy,z,d,r,a,b,c,n,t,cofs;\n n := nops(x);\n \+ z[1]:= z1: z[n]:=zn:\n for i from 1 to n-1 do \n dx[i]:=x[i+1]-x[ i]:\n dy[i]:=y[i+1]-y[i] od:\n for i from 1 to n-2 do d[i]:= 2*(dx[ i]+dx[i+1]):\n r[i]:= 6*(dy[i+1]/dx[i+1]-dy[i]/dx[i]) od:\n r[1] := r[1]-dx[1]*z[1]: \n for i from 2 to n-2 do \n d[i] := d[i]-dx[i +1]^2/d[i-1]: \n r[i] := r[i]-dx[i]*r[i-1]/d[i-1] \n od:\n for \+ i from n-1 by -1 to 2 do\n z[i] := (r[i-1]-dx[i]*z[i+1])/d[i-1] od:\n for i from 1 to n-1 do\n a[i] := 1/6*(z[i+1]-z[i])/dx[i]: \n b[i] := z[i]/2:\n c[i] := dy[i]/dx[i]-a[i]*dx[i]^2-b[i]*dx[i] : \n od:\n for i from 1 to n-1 do\n cofs[i,1]:=a[i]: cofs[i,2]:=b[i]: \n cofs [i,3]:= c[i]: cofs[i,4]:=y[i]: od:\n cofs;\n end:\n" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 14 "with(codegen);" }}{PARA 7 "" 1 "" {TEXT -1 29 "Warning, new definition for C" }}{PARA 12 "" 1 "" {XPPMATH 20 "6#7;%\"CG%%GRADG%)GRADIENTG%(HESSIANG%)JACOBIANG%%costG%( declareG%+dontreturnG%(fortranG%'hornerG%-intrep2mapleG%*joinprocsG%+m akeglobalG%*makeparamG%)makeprocG%)makevoidG%-maple2intrepG%)optimizeG %)packargsG%+packlocalsG%+packparamsG%+prep2transG%*renamevarG%&splitG %)swapargsG" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 598 "Note: There are \+ loads of warning messages about not being able to determine arraybound s for various things. Pay no heed! Just color all of the output an d, using the mouse change the paragraph style (in the context bar abov e) . Then you can modify the exported code and copy it to a separate file for compilation. We will need to do some modification of the \+ exported code. Namely, we need to declare the vectors and matrices us ed in the code using the types VECTOR_double and MATRIX_double defined in Pozo's MV. We also add two include statements for the headers acco mpanying his library." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" } }}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 20 "C(cubicspline,ansi);" }} {PARA 6 "" 1 "" {TEXT -1 69 "C/pretreatment: WARNING: Unable to dete rmine array bounds in z[1]" }}{PARA 6 "" 1 "" {TEXT -1 69 "C/pretrea tment: WARNING: Unable to determine array bounds in z[n]" }}{PARA 6 "" 1 "" {TEXT -1 70 "C/pretreatment: WARNING: Unable to determine \+ array bounds in dx[i]" }}{PARA 6 "" 1 "" {TEXT -1 71 "C/pretreatment : WARNING: Unable to determine array bounds in x[i+1]" }}{PARA 6 " " 1 "" {TEXT -1 69 "C/pretreatment: WARNING: Unable to determine arr ay bounds in x[i]" }}{PARA 6 "" 1 "" {TEXT -1 70 "C/pretreatment: \+ WARNING: Unable to determine array bounds in dy[i]" }}{PARA 6 "" 1 " " {TEXT -1 71 "C/pretreatment: WARNING: Unable to determine array bo unds in y[i+1]" }}{PARA 6 "" 1 "" {TEXT -1 69 "C/pretreatment: WAR NING: Unable to determine array bounds in y[i]" }}{PARA 6 "" 1 "" {TEXT -1 69 "C/pretreatment: WARNING: Unable to determine array boun ds in d[i]" }}{PARA 6 "" 1 "" {TEXT -1 70 "C/pretreatment: WARNING : Unable to determine array bounds in dx[i]" }}{PARA 6 "" 1 "" {TEXT -1 72 "C/pretreatment: WARNING: Unable to determine array boun ds in dx[i+1]" }}{PARA 6 "" 1 "" {TEXT -1 69 "C/pretreatment: WARN ING: Unable to determine array bounds in r[i]" }}{PARA 6 "" 1 "" {TEXT -1 72 "C/pretreatment: WARNING: Unable to determine array boun ds in dy[i+1]" }}{PARA 6 "" 1 "" {TEXT -1 72 "C/pretreatment: WARN ING: Unable to determine array bounds in dx[i+1]" }}{PARA 6 "" 1 "" {TEXT -1 70 "C/pretreatment: WARNING: Unable to determine array boun ds in dy[i]" }}{PARA 6 "" 1 "" {TEXT -1 70 "C/pretreatment: WARNIN G: Unable to determine array bounds in dx[i]" }}{PARA 6 "" 1 "" {TEXT -1 69 "C/pretreatment: WARNING: Unable to determine array boun ds in r[1]" }}{PARA 6 "" 1 "" {TEXT -1 69 "C/pretreatment: WARNING : Unable to determine array bounds in r[1]" }}{PARA 6 "" 1 "" {TEXT -1 70 "C/pretreatment: WARNING: Unable to determine array bounds in \+ dx[1]" }}{PARA 6 "" 1 "" {TEXT -1 69 "C/pretreatment: WARNING: Una ble to determine array bounds in z[1]" }}{PARA 6 "" 1 "" {TEXT -1 69 "C/pretreatment: WARNING: Unable to determine array bounds in d [i]" }}{PARA 6 "" 1 "" {TEXT -1 69 "C/pretreatment: WARNING: Unable \+ to determine array bounds in d[i]" }}{PARA 6 "" 1 "" {TEXT -1 72 "C/ pretreatment: WARNING: Unable to determine array bounds in dx[i+1] " }}{PARA 6 "" 1 "" {TEXT -1 71 "C/pretreatment: WARNING: Unable to \+ determine array bounds in d[i-1]" }}{PARA 6 "" 1 "" {TEXT -1 69 "C/p retreatment: WARNING: Unable to determine array bounds in r[i]" }} {PARA 6 "" 1 "" {TEXT -1 69 "C/pretreatment: WARNING: Unable to dete rmine array bounds in r[i]" }}{PARA 6 "" 1 "" {TEXT -1 70 "C/pretrea tment: WARNING: Unable to determine array bounds in dx[i]" }} {PARA 6 "" 1 "" {TEXT -1 71 "C/pretreatment: WARNING: Unable to dete rmine array bounds in r[i-1]" }}{PARA 6 "" 1 "" {TEXT -1 71 "C/pretr eatment: WARNING: Unable to determine array bounds in d[i-1]" }} {PARA 6 "" 1 "" {TEXT -1 69 "C/pretreatment: WARNING: Unable to dete rmine array bounds in z[i]" }}{PARA 6 "" 1 "" {TEXT -1 71 "C/pretrea tment: WARNING: Unable to determine array bounds in r[i-1]" }} {PARA 6 "" 1 "" {TEXT -1 70 "C/pretreatment: WARNING: Unable to dete rmine array bounds in dx[i]" }}{PARA 6 "" 1 "" {TEXT -1 71 "C/pretre atment: WARNING: Unable to determine array bounds in z[i+1]" }} {PARA 6 "" 1 "" {TEXT -1 71 "C/pretreatment: WARNING: Unable to dete rmine array bounds in d[i-1]" }}{PARA 6 "" 1 "" {TEXT -1 69 "C/pretr eatment: WARNING: Unable to determine array bounds in a[i]" }} {PARA 6 "" 1 "" {TEXT -1 71 "C/pretreatment: WARNING: Unable to dete rmine array bounds in z[i+1]" }}{PARA 6 "" 1 "" {TEXT -1 69 "C/pretr eatment: WARNING: Unable to determine array bounds in z[i]" }} {PARA 6 "" 1 "" {TEXT -1 70 "C/pretreatment: WARNING: Unable to dete rmine array bounds in dx[i]" }}{PARA 6 "" 1 "" {TEXT -1 69 "C/pretre atment: WARNING: Unable to determine array bounds in b[i]" }} {PARA 6 "" 1 "" {TEXT -1 69 "C/pretreatment: WARNING: Unable to dete rmine array bounds in z[i]" }}{PARA 6 "" 1 "" {TEXT -1 69 "C/pretrea tment: WARNING: Unable to determine array bounds in c[i]" }}{PARA 6 "" 1 "" {TEXT -1 70 "C/pretreatment: WARNING: Unable to determine \+ array bounds in dy[i]" }}{PARA 6 "" 1 "" {TEXT -1 70 "C/pretreatment : WARNING: Unable to determine array bounds in dx[i]" }}{PARA 6 " " 1 "" {TEXT -1 69 "C/pretreatment: WARNING: Unable to determine arr ay bounds in a[i]" }}{PARA 6 "" 1 "" {TEXT -1 70 "C/pretreatment: \+ WARNING: Unable to determine array bounds in dx[i]" }}{PARA 6 "" 1 " " {TEXT -1 69 "C/pretreatment: WARNING: Unable to determine array bo unds in b[i]" }}{PARA 6 "" 1 "" {TEXT -1 70 "C/pretreatment: WARNI NG: Unable to determine array bounds in dx[i]" }}{PARA 6 "" 1 "" {TEXT -1 74 "C/pretreatment: WARNING: Unable to determine array boun ds in cofs[i,1]" }}{PARA 6 "" 1 "" {TEXT -1 69 "C/pretreatment: WA RNING: Unable to determine array bounds in a[i]" }}{PARA 6 "" 1 "" {TEXT -1 74 "C/pretreatment: WARNING: Unable to determine array boun ds in cofs[i,2]" }}{PARA 6 "" 1 "" {TEXT -1 69 "C/pretreatment: WA RNING: Unable to determine array bounds in b[i]" }}{PARA 6 "" 1 "" {TEXT -1 74 "C/pretreatment: WARNING: Unable to determine array boun ds in cofs[i,3]" }}{PARA 6 "" 1 "" {TEXT -1 69 "C/pretreatment: WA RNING: Unable to determine array bounds in c[i]" }}{PARA 6 "" 1 "" {TEXT -1 74 "C/pretreatment: WARNING: Unable to determine array boun ds in cofs[i,4]" }}{PARA 6 "" 1 "" {TEXT -1 69 "C/pretreatment: WA RNING: Unable to determine array bounds in y[i]" }}{PARA 6 "" 1 "" {TEXT -1 17 "#include " }}{PARA 6 "" 1 "" {TEXT -1 57 "double \+ cubicspline(double x,double y,double z1,double zn)" }}{PARA 6 "" 1 "" {TEXT -1 1 "\{" }}{PARA 6 "" 1 "" {TEXT -1 8 " int i;" }}{PARA 6 "" 1 "" {TEXT -1 12 " double dx;" }}{PARA 6 "" 1 "" {TEXT -1 12 " doubl e dy;" }}{PARA 6 "" 1 "" {TEXT -1 11 " double z;" }}{PARA 6 "" 1 "" {TEXT -1 11 " double d;" }}{PARA 6 "" 1 "" {TEXT -1 11 " double r;" }}{PARA 6 "" 1 "" {TEXT -1 11 " double a;" }}{PARA 6 "" 1 "" {TEXT -1 11 " double b;" }}{PARA 6 "" 1 "" {TEXT -1 11 " double c;" }} {PARA 6 "" 1 "" {TEXT -1 8 " int n;" }}{PARA 6 "" 1 "" {TEXT -1 11 " \+ double t;" }}{PARA 6 "" 1 "" {TEXT -1 14 " double cofs;" }}{PARA 6 " " 1 "" {TEXT -1 3 " \{" }}{PARA 6 "" 1 "" {TEXT -1 16 " n = nops(x );" }}{PARA 6 "" 1 "" {TEXT -1 14 " z[0] = z1;" }}{PARA 6 "" 1 "" {TEXT -1 16 " z[n-1] = zn;" }}{PARA 6 "" 1 "" {TEXT -1 29 " for( i = 1;i <= n-1.0;i++)" }}{PARA 6 "" 1 "" {TEXT -1 5 " \{" }}{PARA 6 "" 1 "" {TEXT -1 28 " dx[i-1] = x[i]-x[i-1];" }}{PARA 6 "" 1 " " {TEXT -1 28 " dy[i-1] = y[i]-y[i-1];" }}{PARA 6 "" 1 "" {TEXT -1 5 " \}" }}{PARA 6 "" 1 "" {TEXT -1 29 " for(i = 1;i <= n-2.0; i++)" }}{PARA 6 "" 1 "" {TEXT -1 5 " \{" }}{PARA 6 "" 1 "" {TEXT -1 37 " d[i-1] = 2.0*dx[i-1]+2.0*dx[i];" }}{PARA 6 "" 1 "" {TEXT -1 51 " r[i-1] = 6.0*dy[i]/dx[i]-6.0*dy[i-1]/dx[i-1];" }}{PARA 6 "" 1 "" {TEXT -1 5 " \}" }}{PARA 6 "" 1 "" {TEXT -1 24 " r[0] += -dx[0]*z[0];" }}{PARA 6 "" 1 "" {TEXT -1 29 " for(i = 2;i <= n-2.0 ;i++)" }}{PARA 6 "" 1 "" {TEXT -1 5 " \{" }}{PARA 6 "" 1 "" {TEXT -1 36 " d[i-1] += -dx[i]*dx[i]/d[i-2];" }}{PARA 6 "" 1 "" {TEXT -1 39 " r[i-1] += -dx[i-1]*r[i-2]/d[i-2];" }}{PARA 6 "" 1 "" {TEXT -1 5 " \}" }}{PARA 6 "" 1 "" {TEXT -1 29 " for(i = n-1;2.0 <= i;i--)" }}{PARA 6 "" 1 "" {TEXT -1 44 " z[i-1] = (r[i-2]-dx[i -1]*z[i])/d[i-2];" }}{PARA 6 "" 1 "" {TEXT -1 29 " for(i = 1;i <= n -1.0;i++)" }}{PARA 6 "" 1 "" {TEXT -1 5 " \{" }}{PARA 6 "" 1 "" {TEXT -1 41 " a[i-1] = (z[i]-z[i-1])/dx[i-1]/6.0;" }}{PARA 6 "" 1 "" {TEXT -1 26 " b[i-1] = z[i-1]/2.0;" }}{PARA 6 "" 1 "" {TEXT -1 70 " c[i-1] = dy[i-1]/dx[i-1]-a[i-1]*pow(dx[i-1],2.0)-b[i-1]*d x[i-1];" }}{PARA 6 "" 1 "" {TEXT -1 5 " \}" }}{PARA 6 "" 1 "" {TEXT -1 29 " for(i = 1;i <= n-1.0;i++)" }}{PARA 6 "" 1 "" {TEXT -1 5 " \{" }}{PARA 6 "" 1 "" {TEXT -1 28 " cofs[i-1][0] = a[i- 1];" }}{PARA 6 "" 1 "" {TEXT -1 28 " cofs[i-1][1] = b[i-1];" }} {PARA 6 "" 1 "" {TEXT -1 28 " cofs[i-1][2] = c[i-1];" }}{PARA 6 " " 1 "" {TEXT -1 28 " cofs[i-1][3] = y[i-1];" }}{PARA 6 "" 1 "" {TEXT -1 5 " \}" }}{PARA 6 "" 1 "" {TEXT -1 17 " return(cofs);" }}{PARA 6 "" 1 "" {TEXT -1 3 " \}" }}{PARA 6 "" 1 "" {TEXT -1 1 "\}" }}{PARA 6 "" 1 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 116 "Here is the exported code stripped of the error messages, but oth erwise unchanged. Below that is the modified code." }}{PARA 0 "" 0 " " {TEXT -1 0 "" }}{PARA 2 "" 1 "" {TEXT -1 17 "#include " }} {PARA 2 "" 1 "" {TEXT -1 74 "void cubicspline(double x,double y,double z1,double zn,double crea_par[4])" }}{PARA 2 "" 1 "" {TEXT -1 1 "\{" } }{PARA 2 "" 1 "" {TEXT -1 8 " int i;" }}{PARA 2 "" 1 "" {TEXT -1 12 " double dx;" }}{PARA 2 "" 1 "" {TEXT -1 12 " double dy;" }}{PARA 2 " " 1 "" {TEXT -1 11 " double z;" }}{PARA 2 "" 1 "" {TEXT -1 11 " doub le d;" }}{PARA 2 "" 1 "" {TEXT -1 11 " double r;" }}{PARA 2 "" 1 "" {TEXT -1 11 " double a;" }}{PARA 2 "" 1 "" {TEXT -1 11 " double b;" }}{PARA 2 "" 1 "" {TEXT -1 11 " double c;" }}{PARA 2 "" 1 "" {TEXT -1 8 " int n;" }}{PARA 2 "" 1 "" {TEXT -1 11 " double t;" }}{PARA 2 "" 1 "" {TEXT -1 3 " \{" }}{PARA 2 "" 1 "" {TEXT -1 16 " n = nops( x);" }}{PARA 2 "" 1 "" {TEXT -1 14 " z[0] = z1;" }}{PARA 2 "" 1 "" {TEXT -1 16 " z[n-1] = zn;" }}{PARA 2 "" 1 "" {TEXT -1 29 " for( i = 1;i <= n-1.0;i++)" }}{PARA 2 "" 1 "" {TEXT -1 5 " \{" }}{PARA 2 "" 1 "" {TEXT -1 28 " dx[i-1] = x[i]-x[i-1];" }}{PARA 2 "" 1 " " {TEXT -1 28 " dy[i-1] = y[i]-y[i-1];" }}{PARA 2 "" 1 "" {TEXT -1 5 " \}" }}{PARA 2 "" 1 "" {TEXT -1 29 " for(i = 1;i <= n-2.0; i++)" }}{PARA 2 "" 1 "" {TEXT -1 5 " \{" }}{PARA 2 "" 1 "" {TEXT -1 37 " d[i-1] = 2.0*dx[i-1]+2.0*dx[i];" }}{PARA 2 "" 1 "" {TEXT -1 51 " r[i-1] = 6.0*dy[i]/dx[i]-6.0*dy[i-1]/dx[i-1];" }}{PARA 2 "" 1 "" {TEXT -1 5 " \}" }}{PARA 2 "" 1 "" {TEXT -1 24 " r[0] += -dx[0]*z[0];" }}{PARA 2 "" 1 "" {TEXT -1 29 " for(i = 2;i <= n-2.0 ;i++)" }}{PARA 2 "" 1 "" {TEXT -1 5 " \{" }}{PARA 2 "" 1 "" {TEXT -1 36 " d[i-1] += -dx[i]*dx[i]/d[i-2];" }}{PARA 2 "" 1 "" {TEXT -1 39 " r[i-1] += -dx[i-1]*r[i-2]/d[i-2];" }}{PARA 2 "" 1 "" {TEXT -1 5 " \}" }}{PARA 2 "" 1 "" {TEXT -1 29 " for(i = n-1;2.0 <= i;i--)" }}{PARA 2 "" 1 "" {TEXT -1 44 " z[i-1] = (r[i-2]-dx[i -1]*z[i])/d[i-2];" }}{PARA 2 "" 1 "" {TEXT -1 29 " for(i = 1;i <= n -1.0;i++)" }}{PARA 2 "" 1 "" {TEXT -1 5 " \{" }}{PARA 2 "" 1 "" {TEXT -1 41 " a[i-1] = (z[i]-z[i-1])/dx[i-1]/6.0;" }}{PARA 2 "" 1 "" {TEXT -1 26 " b[i-1] = z[i-1]/2.0;" }}{PARA 2 "" 1 "" {TEXT -1 70 " c[i-1] = dy[i-1]/dx[i-1]-a[i-1]*pow(dx[i-1],2.0)-b[i-1]*d x[i-1];" }}{PARA 2 "" 1 "" {TEXT -1 5 " \}" }}{PARA 2 "" 1 "" {TEXT -1 20 " crea_par[0] = a;" }}{PARA 2 "" 1 "" {TEXT -1 20 " \+ crea_par[1] = b;" }}{PARA 2 "" 1 "" {TEXT -1 20 " crea_par[2] = c; " }}{PARA 2 "" 1 "" {TEXT -1 20 " crea_par[3] = y;" }}{PARA 2 "" 1 "" {TEXT -1 11 " return;" }}{PARA 2 "" 1 "" {TEXT -1 3 " \}" }} {PARA 2 "" 1 "" {TEXT -1 1 "\}" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}} {EXCHG {PARA 0 "" 0 "" {TEXT -1 90 "// This C++ word takes data x and y and values z1 and zn for the second derivative at the" }}{PARA 0 " " 0 "" {TEXT -1 89 "// right hand and left hand knots and returns the coefficients cofs of the cubic spline " }}{PARA 0 "" 0 "" {TEXT -1 33 "// which interpolates this data." }}{PARA 0 "" 0 "" {TEXT -1 1134 "MATRIX_double cubicspline(VECTOR_double x,VECTOR_double y,double z1,double zn)\n\{\n int n=x.size(); \n VECTOR_double dx(n,0);\n VE CTOR_double dy(n,0);\n VECTOR_double z(n,0);\n VECTOR_double d(n -1,0);\n VECTOR_double r(n,0);\n VECTOR_double a(n-1,0) ;\n VEC TOR_double b(n-1,0) ;\n VECTOR_double c(n-1,0) ;\n MATRIX_double c ofs(n-1,4);\n double t;\n \n \{\n int i =0;\n z(0) = z1;\n \+ z(n-1) = zn;\n for(i = 1;i <= n-1.0;i++)\n \{\n dx(i-1) = \+ x(i)-x(i-1);\n dy(i-1) = y(i)-y(i-1);\n \}\n for( i = 1;i < = n-2.0;i++)\n \{\n d(i-1) = 2.0*dx(i-1)+2.0*dx(i);\n r(i -1) = 6.0*dy(i)/dx(i)-6.0*dy(i-1)/dx(i-1);\n \}\n r(0) += -dx(0) *z(0);\n for(i = 2;i <= n-2.0;i++)\n \{\n d(i-1) += -dx(i)* dx(i)/d(i-2);\n r(i-1) += -dx(i-1)*r(i-2)/d(i-2);\n \}\n fo r(i = n-1;2.0 <= i;i--)\n z(i-1) = (r(i-2)-dx(i-1)*z(i))/d(i-2); \n for(i = 1;i <= n-1.0;i++)\n \{\n a(i-1) = (z(i)-z(i-1))/ dx(i-1)/6.0;\n b(i-1) = z(i-1)/2.0;\n c(i-1) = dy(i-1)/dx(i- 1)-a(i-1)*pow(dx(i-1),2.0)-b(i-1)*dx(i-1); \n \}\n for(i=0;i " 0 "" {MPLTEXT 1 0 218 "csplfun := proc(x,cofs,t)\n local n,i,j;\n n := nops(x);\n i := n-1;\n for j from 1 to n-2 do\n if t<=x[j+1] then \+ i := j ; break; fi; od; RETURN(((cofs[i,1]*(t-x[i])+cofs[i,2])*(t-x[i] )+cofs[i,3])*(t-x[i])+cofs[i,4]) end: " }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 16 "C(csplfun,a nsi);" }}{PARA 6 "" 1 "" {TEXT -1 71 "C/pretreatment: WARNING: Unabl e to determine array bounds in x[j+1]" }}{PARA 6 "" 1 "" {TEXT -1 74 "C/pretreatment: WARNING: Unable to determine array bounds in c ofs[i,1]" }}{PARA 6 "" 1 "" {TEXT -1 69 "C/pretreatment: WARNING: Un able to determine array bounds in x[i]" }}{PARA 6 "" 1 "" {TEXT -1 74 "C/pretreatment: WARNING: Unable to determine array bounds in c ofs[i,2]" }}{PARA 6 "" 1 "" {TEXT -1 69 "C/pretreatment: WARNING: Un able to determine array bounds in x[i]" }}{PARA 6 "" 1 "" {TEXT -1 74 "C/pretreatment: WARNING: Unable to determine array bounds in c ofs[i,3]" }}{PARA 6 "" 1 "" {TEXT -1 69 "C/pretreatment: WARNING: Un able to determine array bounds in x[i]" }}{PARA 6 "" 1 "" {TEXT -1 74 "C/pretreatment: WARNING: Unable to determine array bounds in c ofs[i,4]" }}{PARA 6 "" 1 "" {TEXT -1 17 "#include " }}{PARA 6 "" 1 "" {TEXT -1 45 "double csplfun(double x,double cofs,double t)" }} {PARA 6 "" 1 "" {TEXT -1 1 "\{" }}{PARA 6 "" 1 "" {TEXT -1 11 " doubl e n;" }}{PARA 6 "" 1 "" {TEXT -1 8 " int i;" }}{PARA 6 "" 1 "" {TEXT -1 8 " int j;" }}{PARA 6 "" 1 "" {TEXT -1 3 " \{" }}{PARA 6 "" 1 "" {TEXT -1 16 " n = nops(x);" }}{PARA 6 "" 1 "" {TEXT -1 12 " i = \+ n-1;" }}{PARA 6 "" 1 "" {TEXT -1 29 " for(j = 1;j <= n-2.0;j++)" }} {PARA 6 "" 1 "" {TEXT -1 21 " if( t <= x[j] )" }}{PARA 6 "" 1 "" {TEXT -1 9 " \{" }}{PARA 6 "" 1 "" {TEXT -1 16 " i = j ;" }}{PARA 6 "" 1 "" {TEXT -1 16 " break;" }}{PARA 6 "" 1 "" {TEXT -1 9 " \}" }}{PARA 6 "" 1 "" {TEXT -1 79 " return(((co fs[i-1][0]*(t-x[i-1])+cofs[i-1][1])*(t-x[i-1])+cofs[i-1][2])*(t-" }} {PARA 6 "" 1 "" {TEXT -1 22 "x[i-1])+cofs[i-1][3]);" }}{PARA 6 "" 1 " " {TEXT -1 3 " \}" }}{PARA 6 "" 1 "" {TEXT -1 1 "\}" }}{PARA 6 "" 1 " " {TEXT -1 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 54 "Now here is a m odification that will compile with gcc." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 410 "#include \"vecdefs.h\"\n#include MA TRIX_H \n#include \n\ndouble csplfun(const VECTOR_double & x,c onst MATRIX_double & cofs,double t)\n\{\n int n;\n int i;\n int j; \n \{\n n = x.size();\n i = n-1;\n for(j = 1;j <= n-2.0;j++) \n if( t <= x(j) )\n \{\n i = j;\n break ;\n \}\n return(((cofs(i-1,0)*(t-x(i-1))+cofs(i-1,1))*(t-x(i -1))+cofs(i-1,2))*(t-\nx(i-1))+cofs(i-1,3));\n \}\n\}\n" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 0 "" }}}}{MARK "17 0 0" 0 }{VIEWOPTS 1 1 0 1 1 1803 }