diff options
Diffstat (limited to 'src/GeneralizedPolynomial-Uniform/Lagrange/2d.log')
| -rw-r--r-- | src/GeneralizedPolynomial-Uniform/Lagrange/2d.log | 4482 |
1 files changed, 4482 insertions, 0 deletions
diff --git a/src/GeneralizedPolynomial-Uniform/Lagrange/2d.log b/src/GeneralizedPolynomial-Uniform/Lagrange/2d.log new file mode 100644 index 0000000..4fa17e7 --- /dev/null +++ b/src/GeneralizedPolynomial-Uniform/Lagrange/2d.log @@ -0,0 +1,4482 @@ + |\^/| Maple 7 (IBM INTEL LINUX) +._|\| |/|_. Copyright (c) 2001 by Waterloo Maple Inc. + \ MAPLE / All rights reserved. Maple is a registered trademark of + <____ ____> Waterloo Maple Inc. + | Type ? for help. +# util.maple -- misc utility routines +# $Id: util.maple,v 1.3 2002/05/19 13:12:18 jthorn Exp $ +> +# +# fix_rationals - convert numbers to RATIONAL() calls +# nonmatching_names - find names in a list which *don't* have a specified prefix +# sprint_numeric_list - convert a numeric list to a valid C identifier suffix +# print_name_list_dcl - print a C declaration for a list of names +# +# hypercube_points - compute all (integer) points in an N-dimensional hypercube +# +# ftruncate - truncate a file to zero length +# +> +################################################################################ +################################################################################ +################################################################################ +> +# +# This function converts all {integer, rational} subexpressions of its +# input except integer exponents and -1 factors in products, into function +# calls +# RATIONAL(num,den) +# This is useful in conjunction with the C() library function, since +# +# C( (1/3) * foo * bar ) +# t0 = foo*bar/3; +# +# generates a (slow) division (and runs the risk of mixed-mode-arithmetic +# problems), while +# +# C((1.0/3.0) * foo * bar); +# t0 = 0.3333333333*foo*bar; +# +# suffers from roundoff error. With this function, +# +# fix_rationals((1/3) * foo * bar); +# RATIONAL(1,3) foo bar +# C(%); +# t0 = RATIONAL(1.0,3.0)*foo*bar; +# +# which a C preprocessor macro can easily convert to the desired +# +# t0 = (1.0/3.0)*foo*bar; +# +# Additionally, this function can be told to leave certain types of +# subexpressions unconverged. For example, +# fix_rationals(expr, type, specfunc(integer, DATA)); +# will leave all subexpressions of the form DATA(integer arguments) +# unconverted. +# +# Arguments: +# expr = (in) The expression to be converted. +# inert_fn = (optional in) +# If specified, this argument should be a Boolean procedure +# or the name of a Boolean procedure. This procedure should +# take one or more argument, and return true if and only if +# the first argument should *not* be converted, i.e. if we +# should leave this expression unchanged. See the last +# example above. +# ... = (optional in) +# Any further arguments are passed as additional arguments to +# the inert_fn procedure. +# +> fix_rationals := +> proc( +> expr::{ +> algebraic, name = algebraic, +> list({algebraic, name = algebraic}), +> set ({algebraic, name = algebraic}) +> }, +> inert_fn::{name, procedure} +> ) +> local nn, k, +> base, power, fbase, fpower, +> fn, fn_args_list, +> num, den, mult; +> +# do we want to convert this expression? +> if ((nargs >= 2) and inert_fn(expr, args[3..nargs])) +> then return expr; +> end if; +> +# recurse over lists and sets +> if (type(expr, {list,set})) +> then return map(fix_rationals, expr, args[2..nargs]); +> end if; +> +# recurse over equation right hand sides +> if (type(expr, name = algebraic)) +> then return ( lhs(expr) = fix_rationals(rhs(expr), args[2..nargs]) ); +> end if; +> +# recurse over functions other than RATIONAL() +> if (type(expr, function)) +> then +> fn := op(0, expr); +> if (fn <> 'RATIONAL') +> then +> fn_args_list := [op(expr)]; +> fn_args_list := map(fix_rationals, fn_args_list, args[2..nargs]); +> fn; return '%'( op(fn_args_list) ); +> end if; +> end if; +> +> nn := nops(expr); +> +# recurse over sums +> if (type(expr, `+`)) +> then return sum('fix_rationals(op(k,expr), args[2..nargs])', 'k'=1..nn); +> end if; +> +# recurse over products +# ... leaving leading -1 factors intact, i.e. not converted to RATIONAL(-1,1) +> if (type(expr, `*`)) +> then +> if (op(1, expr) = -1) +> then return -1*fix_rationals(remove(type, expr, 'identical(-1)'), +> args[2..nargs]); +> else return product('fix_rationals(op(k,expr), args[2..nargs])', +> 'k'=1..nn); +> end if; +> end if; +> +# recurse over powers +# ... leaving integer exponents intact +> if (type(expr, `^`)) +> then +> base := op(1, expr); +> power := op(2, expr); +> +> fbase := fix_rationals(base, args[2..nargs]); +> if (type(power, integer)) +> then fpower := power; +> else fpower := fix_rationals(power, args[2..nargs]); +> end if; +> return fbase ^ fpower; +> end if; +> +# fix integers and fractions +> if (type(expr, integer)) +> then return 'RATIONAL'(expr, 1); +> end if; +> if (type(expr, fraction)) +> then +> num := op(1, expr); +> den := op(2, expr); +> +> return 'RATIONAL'(num, den); +> end if; +> +# turn Maple floating-point into integer fraction, then recursively fix that +> if (type(expr, float)) +> then +> mult := op(1, expr); +> power := op(2, expr); +> return fix_rationals(mult * 10^power, args[2..nargs]); +> end if; +> +# identity op on names +> if (type(expr, name)) +> then return expr; +> end if; +> +# unknown type +> error "%0", +> "unknown type for expr!", +> " whattype(expr) = ", whattype(expr), +> " expr = ", expr; +> end proc; +fix_rationals := proc(expr::{algebraic, name = algebraic, +list({algebraic, name = algebraic}), set({algebraic, name = algebraic})}, +inert_fn::{procedure, name}) +local nn, k, base, power, fbase, fpower, fn, fn_args_list, num, den, mult; + if 2 <= nargs and inert_fn(expr, args[3 .. nargs]) then return expr + end if; + if type(expr, {set, list}) then + return map(fix_rationals, expr, args[2 .. nargs]) + end if; + if type(expr, name = algebraic) then + return lhs(expr) = fix_rationals(rhs(expr), args[2 .. nargs]) + end if; + if type(expr, function) then + fn := op(0, expr); + if fn <> 'RATIONAL' then + fn_args_list := [op(expr)]; + fn_args_list := + map(fix_rationals, fn_args_list, args[2 .. nargs]); + fn; + return '%'(op(fn_args_list)) + end if + end if; + nn := nops(expr); + if type(expr, `+`) then return + sum('fix_rationals(op(k, expr), args[2 .. nargs])', 'k' = 1 .. nn) + end if; + if type(expr, `*`) then + if op(1, expr) = -1 then return -fix_rationals( + remove(type, expr, 'identical(-1)'), args[2 .. nargs]) + else return product('fix_rationals(op(k, expr), args[2 .. nargs])', + 'k' = 1 .. nn) + end if + end if; + if type(expr, `^`) then + base := op(1, expr); + power := op(2, expr); + fbase := fix_rationals(base, args[2 .. nargs]); + if type(power, integer) then fpower := power + else fpower := fix_rationals(power, args[2 .. nargs]) + end if; + return fbase^fpower + end if; + if type(expr, integer) then return 'RATIONAL'(expr, 1) end if; + if type(expr, fraction) then + num := op(1, expr); den := op(2, expr); return 'RATIONAL'(num, den) + end if; + if type(expr, float) then + mult := op(1, expr); + power := op(2, expr); + return fix_rationals(mult*10^power, args[2 .. nargs]) + end if; + if type(expr, name) then return expr end if; + error "%0", "unknown type for expr!", " whattype(expr) = ", + whattype(expr), " expr = ", expr +end proc + +> +################################################################################ +> +# +# This function finds names in a list which *don't* have a specified prefix. +# +# Arguments: +# name_list = A list of the names. +# prefix = The prefix we want to filter out. +# +# Results: +# This function returns the subset list of names which don't have the +# specified prefix. +# +> nonmatching_names := +> proc( name_list::list({name,string}), prefix::{name,string} ) +> +> select( proc(n) +> evalb(not StringTools[IsPrefix](prefix,n)); +> end proc +> , +> name_list +> ); +> end proc; +nonmatching_names := proc( +name_list::list({name, string}), prefix::{name, string}) + select(proc(n) evalb(not StringTools[IsPrefix](prefix, n)) end proc, + name_list) +end proc + +> +################################################################################ +> +# +# This function converts a numeric list to a string which is a valid +# C identifier suffix: elements are separated by "_", decimal points are +# replaced by "x", and all nonzero values have explicit +/- signs, which +# are replaced by "p"/"m". +# +# For example, [0,-3.5,+4] --> "0_m3x5_p4". +# +> sprint_numeric_list := +> proc(nlist::list(numeric)) +> +# generate preliminary string, eg "+0_-3.5_+4" +> map2(sprintf, "%+a", nlist); +> ListTools[Join](%, "_"); +> cat(op(%)); +> +# fixup bad characters +> StringTools[SubstituteAll](%, "+0", "0"); +> StringTools[CharacterMap](".+-", "xpm", %); +> +> return %; +> end proc; +sprint_numeric_list := proc(nlist::list(numeric)) + map2(sprintf, "%+a", nlist); + ListTools[Join](%, "_"); + cat(op(%)); + StringTools[SubstituteAll](%, "+0", "0"); + StringTools[CharacterMap](".+-", "xpm", %); + return % +end proc + +> +################################################################################ +> +# +# This function prints a C declaration for a list of names. +# +# Argument: +# name_list = A list of the names. +# name_type = The C type of the names, eg. "double". +# file_name = The file name to write the declaration to. This is +# truncated before writing. +# +> print_name_list_dcl := +> proc( name_list::list({name,string}), +> name_type::string, +> file_name::string ) +> local blanks, separator_string; +> +> ftruncate(file_name); +> +# a sequence of blanks with the same length as name_type +> seq(" ", i=1..length(name_type)); +> +# string to separate names +> separator_string := cat(",\n", %, " "); +> +> map(convert, name_list, string); +> ListTools[Join](%, separator_string); +> cat(op(%)); +> +> fprintf(file_name, +> "%s %s;\n", +> name_type, %); +> +> fclose(file_name); +> NULL; +> end proc; +print_name_list_dcl := proc( +name_list::list({name, string}), name_type::string, file_name::string) +local blanks, separator_string; + ftruncate(file_name); + seq(" ", i = 1 .. length(name_type)); + separator_string := cat(",\n", %, " "); + map(convert, name_list, string); + ListTools[Join](%, separator_string); + cat(op(%)); + fprintf(file_name, "%s %s;\n", name_type, %); + fclose(file_name); + NULL +end proc + +> +################################################################################ +################################################################################ +################################################################################ +> +# +# This function computes a list of all the (integer) points in an +# N-dimensional hypercube, in lexicographic order. The present +# implementation requires N <= 4. +# +# Arguments: +# cmin,cmax = N-element lists of cube minimum/maximum coordinates. +# +# Results: +# The function returns a set of d-element lists giving the coordinates. +# For example, +# hypercube([0,0], [2,1] +# returns +# { [0,0], [0,1], [1,0], [1,1], [2,0], [2,1] } +> hypercube_points := +> proc(cmin::list(integer), cmax::list(integer)) +> local N, i,j,k,l; +> +> N := nops(cmin); +> if (nops(cmax) <> N) +> then error +> "must have same number of dimensions for min and max coordinates!"; +> fi; +> +> if (N = 1) +> then return [seq([i], i=cmin[1]..cmax[1])]; +> elif (N = 2) +> then return [ +> seq( +> seq([i,j], j=cmin[2]..cmax[2]), +> i=cmin[1]..cmax[1]) +> ]; +> elif (N = 3) +> then return [ +> seq( +> seq( +> seq([i,j,k], k=cmin[3]..cmax[3]), +> j=cmin[2]..cmax[2] ), +> i=cmin[1]..cmax[1]) +> ]; +> elif (N = 4) +> then return [ +> seq( +> seq( +> seq( +> seq([i,j,k,l], l=cmin[4]..cmax[4]), +> k=cmin[3]..cmax[3] ), +> j=cmin[2]..cmax[2]), +> i=cmin[1]..cmax[1]) +> ]; +> else +> error "implementation restriction: must have N <= 4, got %1!", N; +> fi; +> end proc; +hypercube_points := proc(cmin::list(integer), cmax::list(integer)) +local N, i, j, k, l; + N := nops(cmin); + if nops(cmax) <> N then error + "must have same number of dimensions for min and max coordinates!" + end if; + if N = 1 then return [seq([i], i = cmin[1] .. cmax[1])] + elif N = 2 then return + [seq(seq([i, j], j = cmin[2] .. cmax[2]), i = cmin[1] .. cmax[1])] + elif N = 3 then return [seq( + seq(seq([i, j, k], k = cmin[3] .. cmax[3]), j = cmin[2] .. cmax[2]) + , i = cmin[1] .. cmax[1])] + elif N = 4 then return [seq(seq(seq( + seq([i, j, k, l], l = cmin[4] .. cmax[4]), k = cmin[3] .. cmax[3]), + j = cmin[2] .. cmax[2]), i = cmin[1] .. cmax[1])] + else error "implementation restriction: must have N <= 4, got %1!", N + end if +end proc + +> +################################################################################ +################################################################################ +################################################################################ +> +# +# This function truncates a file to 0 length if it exists, or creates +# it at that length if it doesn't exist. +# +# Arguments: +# file_name = (in) The name of the file. +# +> ftruncate := +> proc(file_name::string) +> fopen(file_name, 'WRITE'); +> fclose(%); +> NULL; +> end proc; +ftruncate := + + proc(file_name::string) fopen(file_name, 'WRITE'); fclose(%); NULL end proc + +# interpolate.maple -- compute generalized interpolation formulas/coefficients +# $Id: interpolate.maple,v 1.4 2002/05/14 15:52:50 jthorn Exp $ +> +# +# <<<representation of numbers, data values, etc>>> +# Lagrange polynomial_interpolant - compute Lagrange polynomial interpolant +# coeff_as_lc_of_data - coefficients of ... (linear combination of data) +# +# print_coeff__lc_of_data - print C code to compute coefficients +# print_data_var_assign - print C code to assign data-value variables +# print_interp_coeff_var_store - print C code to store coeff vars "somewhere" +# print_interp_cmpt__lc_of_data - print C code for computation of interpolant +# +# coeff_name - name of coefficient of data at a given [m] coordinate +# data_var_name - name of variable storing data value at a given [m] coordinate +# +> +################################################################################ +> +# +# ***** representation of numbers, data values, etc ***** +# +# We use RATIONAL(p.0,q.0) to denote the rational number p/q. +# +# We use DATA(...) to represent the data values being interpolated at a +# specified [m] coordinate, where the arguments are the [m] coordinates. +# +# We use COEFF(...) to represent the molecule coefficient at a specified +# [m] coordinate, where the arguments are the [m] coordinates. +# +# For example, the usual 1-D centered 2nd order 1st derivative molecule +# would be written +# RATIONAL(-1.0,2.0)*DATA(-1) + RATIONA(1.0,2.0)*DATA(1) +# and its coefficients as +# COEFF(-1) = RATIONAL(-1.0,2.0) +# COEFF(1) = RATIONAL(1.0,2.0) +# +> +################################################################################ +################################################################################ +################################################################################ +> +# +# This function computes a Lagrange polynomial interpolant in any +# number of dimensions. +# +# Arguments: +# fn = The interpolation function. This should be a procedure in the +# coordinates, having the coefficients as global variables. For +# example, +# proc(x,y) c00 + c10*x + c01*y end proc +# coeff_list = A set of the interpolation coefficients (coefficients in +# the interpolation function), for example [c00, c10, c01]. +# coord_list = A list of the coordinates (independent variables in the +# interpolation function), for example [x,y]. +# posn_list = A list of positions (each a list of numeric values) where the +# interpolant is to use data, for example hypercube([0,0], [1,1]). +# Any positions may be used; if they're redundant (as in the +# example) the least-squares interpolant is computed. +# +# Results: +# This function returns the interpolating polynomial, in the form of +# an algebraic expression in the coordinates and the data values. +# +> Lagrange_polynomial_interpolant := +> proc( +> fn::procedure, coeff_list::list(name), +> coord_list::list(name), posn_list::list(list(numeric)) +> ) +> local posn, data_eqns, coeff_eqns; +> +# coefficients of interpolating polynomial +> data_eqns := { seq( fn(op(posn))='DATA'(op(posn)) , posn=posn_list ) }; +> coeff_eqns := linalg[leastsqrs](data_eqns, {op(coeff_list)}); +> if (has(coeff_eqns, '_t')) +> then error "interpolation coefficients aren't uniquely determined!"; +> end if; +> +# interpolant as a polynomial in the coordinates +> return subs(coeff_eqns, eval(fn))(op(coord_list)); +> end proc; +Lagrange_polynomial_interpolant := proc(fn::procedure, coeff_list::list(name), +coord_list::list(name), posn_list::list(list(numeric))) +local posn, data_eqns, coeff_eqns; + data_eqns := {seq(fn(op(posn)) = 'DATA'(op(posn)), posn = posn_list)}; + coeff_eqns := linalg[leastsqrs](data_eqns, {op(coeff_list)}); + if has(coeff_eqns, '_t') then + error "interpolation coefficients aren't uniquely determined!" + end if; + return subs(coeff_eqns, eval(fn))(op(coord_list)) +end proc + +> +################################################################################ +> +# +# This function takes as input an interpolating polynomial, expresses +# it as a linear combination of the data values, and returns the coefficeints +# of that form. +# +# Arguments: +# interpolant = The interpolating polynomial (an algebraic expression +# in the coordinates and the data values). +# posn_list = The same list of positions as was used to compute the +# interpolating polynomial. +# +# Results: +# This function returns the coefficients, as a list of equations of the +# form COEFF(...) = value , where each value is a polynomial in the +# coordinates. The order of the list matches that of posn_list. +# +> coeff_as_lc_of_data := +> proc( +> interpolant::algebraic, +> posn_list::list(list(numeric)) +> ) +> local data_list, interpolant_as_lc_of_data; +> +# interpolant as a linear combination of the data values +> data_list := [ seq( 'DATA'(op(posn)) , posn=posn_list ) ]; +> interpolant_as_lc_of_data := collect(interpolant, data_list); +> +# coefficients of the data values in the linear combination +> return map( +> proc(posn::list(numeric)) +> coeff(interpolant_as_lc_of_data, DATA(op(posn))); +> 'COEFF'(op(posn)) = %; +> end proc +> , +> posn_list +> ); +> end proc; +coeff_as_lc_of_data := proc( +interpolant::algebraic, posn_list::list(list(numeric))) +local data_list, interpolant_as_lc_of_data; + data_list := [seq('DATA'(op(posn)), posn = posn_list)]; + interpolant_as_lc_of_data := collect(interpolant, data_list); + return map(proc(posn::list(numeric)) + coeff(interpolant_as_lc_of_data, DATA(op(posn))); + 'COEFF'(op(posn)) = % + end proc, posn_list) +end proc + +> +################################################################################ +################################################################################ +################################################################################ +> +# +# This function prints C expressions for the coefficients of an +# interpolating polynomial. (The polynomial is expressed as linear +# combinations of the data values with coefficients which are +# RATIONAL(p,q) calls.) +# +# Arguments: +# coeff_list = A list of the coefficients, as returned from +# coeff_as_lc_of_data() . +# coeff_name_prefix = A prefix string for the coefficient names. +# temp_name_type = The C type to be used for Maple-introduced temporary +# names, eg. "double". +# file_name = The file name to write the coefficients to. This is +# truncated before writing. +# +> print_coeff__lc_of_data := +> proc( coeff_list::list(specfunc(numeric,COEFF) = algebraic), +> coeff_name_prefix::string, +> temp_name_type::string, +> file_name::string ) +> global `codegen/C/function/informed`; +> local coeff_list2, cmpt_list, temp_name_list; +> +# convert LHS of each equation from a COEFF() call (eg COEFF(-1,+1)) +# to a Maple/C variable name (eg coeff_I_m1_p1) +> coeff_list2 := map( +> proc(coeff_eqn::specfunc(numeric,COEFF) = algebraic) +> local posn; +> posn := [op(lhs(coeff_eqn))]; +> coeff_name(posn,coeff_name_prefix); +> convert(%, name); # codegen[C] wants LHS +> # to be an actual Maple *name* +> % = fix_rationals(rhs(coeff_eqn)); +> end proc +> , +> coeff_list +> ); +> +# +# generate the C code +# +> +# tell codegen[C] not to warn about unknown RATIONAL() and DATA() "fn calls" +# via undocumented :( global table +> `codegen/C/function/informed`['RATIONAL'] := true; +> `codegen/C/function/informed`['DATA'] := true; +> +> ftruncate(file_name); +> +# optimized computation sequence for all the coefficients +# (may use local variables t0,t1,t2,...) +> cmpt_list := [codegen[optimize](coeff_list2, tryhard)]; +> +# list of the t0,t1,t2,... local variables +> temp_name_list := nonmatching_names(map(lhs,cmpt_list), coeff_name_prefix); +> +# declare the t0,t1,t2,... local variables (if there are any) +> if (nops(temp_name_list) > 0) +> then print_name_list_dcl(%, temp_name_type, file_name); +> fi; +> +# now print the optimized computation sequence +> codegen[C](cmpt_list, filename=file_name); +> +> fclose(file_name); +> +> NULL; +> end proc; +print_coeff__lc_of_data := proc( +coeff_list::list(specfunc(numeric, COEFF) = algebraic), +coeff_name_prefix::string, temp_name_type::string, file_name::string) +local coeff_list2, cmpt_list, temp_name_list; +global `codegen/C/function/informed`; + coeff_list2 := map(proc( + coeff_eqn::(specfunc(numeric, COEFF) = algebraic)) + local posn; + posn := [op(lhs(coeff_eqn))]; + coeff_name(posn, coeff_name_prefix); + convert(%, name); + % = fix_rationals(rhs(coeff_eqn)) + end proc, coeff_list); + `codegen/C/function/informed`['RATIONAL'] := true; + `codegen/C/function/informed`['DATA'] := true; + ftruncate(file_name); + cmpt_list := [codegen[optimize](coeff_list2, tryhard)]; + temp_name_list := + nonmatching_names(map(lhs, cmpt_list), coeff_name_prefix); + if 0 < nops(temp_name_list) then + print_name_list_dcl(%, temp_name_type, file_name) + end if; + codegen[C](cmpt_list, filename = file_name); + fclose(file_name); + NULL +end proc + +> +################################################################################ +> +# +# This function prints a sequence of C expression to assign the data-value +# variables, eg +# data_m1_p1 = DATA(-1,1); +# +# Arguments: +# posn_list = The same list of positions as was used to compute the +# interpolating polynomial. +# data_var_name_prefix = A prefix string for the data variable names. +# file_name = The file name to write the coefficients to. This is +# truncated before writing. +# +> print_data_var_assign := +> proc( +> posn_list::list(list(numeric)), +> data_var_name_prefix::string, +> file_name::string +> ) +> +> ftruncate(file_name); +> map( +> proc(posn::list(numeric)) +> fprintf(file_name, +> "%s = %a;\n", +> data_var_name(posn,data_var_name_prefix), +> DATA(op(posn))); +> end proc +> , +> posn_list +> ); +> fclose(file_name); +> +> NULL; +> end proc; +print_data_var_assign := proc(posn_list::list(list(numeric)), +data_var_name_prefix::string, file_name::string) + ftruncate(file_name); + map(proc(posn::list(numeric)) + fprintf(file_name, "%s = %a;\n", + data_var_name(posn, data_var_name_prefix), DATA(op(posn))) + end proc, posn_list); + fclose(file_name); + NULL +end proc + +> +################################################################################ +> +# +# This function prints a sequence of C expression to store the interpolation +# coefficients in COEFF(...) expressions, eg +# COEFF(1,-1) = factor * coeff_dx_p1_m1; +# +# Arguments: +# posn_list = The same list of positions as was used to compute the +# interpolating polynomial. +# RHS_factor_name = If this string is non-empty, then the coefficient is +# multiplied by this factor before being stored, eg +# setting this to "factor" would give the example above. +# If this string is empty (""), the multiplication is +# omitted, eg +# COEFF(1,-1) = coeff_dx_p1_m1; +# coeff_name_prefix = A prefix string for the coefficient names. +# file_name = The file name to write the coefficients to. This is +# truncated before writing. +# +> print_interp_coeff_var_store := +> proc( +> posn_list::list(list(numeric)), +> RHS_factor_name::string, +> coeff_name_prefix::string, +> file_name::string +> ) +> +> ftruncate(file_name); +> map( +> proc(posn::list(numeric)) +> if (length(RHS_factor_name) > 0) +> then fprintf(file_name, +> "%a = %s * %s;\n", +> 'COEFF'(op(posn)), +> RHS_factor_name, +> coeff_name(posn,coeff_name_prefix)); +> else fprintf(file_name, +> "%a = %s;\n", +> 'COEFF'(op(posn)), +> coeff_name(posn,coeff_name_prefix)); +> end if; +> end proc +> , +> posn_list +> ); +> fclose(file_name); +> +> NULL; +> end proc; +print_interp_coeff_var_store := proc(posn_list::list(list(numeric)), +RHS_factor_name::string, coeff_name_prefix::string, file_name::string) + ftruncate(file_name); + map(proc(posn::list(numeric)) + if 0 < length(RHS_factor_name) then fprintf(file_name, + "%a = %s * %s;\n", 'COEFF'(op(posn)), RHS_factor_name, + coeff_name(posn, coeff_name_prefix)) + else fprintf(file_name, "%a = %s;\n", 'COEFF'(op(posn)), + coeff_name(posn, coeff_name_prefix)) + end if + end proc, posn_list); + fclose(file_name); + NULL +end proc + +> +################################################################################ +> +# +# This function prints a C expression to compute the interpolant, +# using the coefficients computed by print_coeff__lc_of_data() +# (i.e. expressing the interpolant as a linear combination of the +# data values). +# +# Arguments: +# posn_list = The same list of positions as was used to compute the +# interpolating polynomial. +# result_var_name = The (string) name of the variable to which the +# result is to be assigned. +# coeff_name_prefix = A prefix string for the coefficient names. +# data_var_name_prefix = A prefix string for the data variable names. +# file_name = The file name to write the coefficients to. This is +# truncated before writing. +# +> print_interp_cmpt__lc_of_data := +> proc( +> posn_list::list(list(numeric)), +> result_var_name::string, +> coeff_name_prefix::string, +> data_var_name_prefix::string, +> file_name::string +> ) +> +> ftruncate(file_name); +> +> fprintf(file_name, "%s =\n", result_var_name); +> +# list of "coeff*data_var" terms +> map( +> proc(posn::list(numeric)) +> sprintf("%s*%s", +> coeff_name(posn,coeff_name_prefix), +> data_var_name(posn,data_var_name_prefix)); +> end proc +> , +> posn_list +> ); +> +> ListTools[Join](%, "\n\t+ "); +> cat(op(%)); +> fprintf(file_name, "\t%s;\n", %); +> +> fclose(file_name); +> +> NULL; +> end proc; +print_interp_cmpt__lc_of_data := proc(posn_list::list(list(numeric)), +result_var_name::string, coeff_name_prefix::string, +data_var_name_prefix::string, file_name::string) + ftruncate(file_name); + fprintf(file_name, "%s =\n", result_var_name); + map(proc(posn::list(numeric)) + sprintf("%s*%s", coeff_name(posn, coeff_name_prefix), + data_var_name(posn, data_var_name_prefix)) + end proc, posn_list); + ListTools[Join](%, "\n\t+ "); + cat(op(%)); + fprintf(file_name, "\t%s;\n", %); + fclose(file_name); + NULL +end proc + +> +################################################################################ +################################################################################ +################################################################################ +> +# +# This function computes the name of the coefficient of the data at a +# given [m] position, i.e. it encapsulates our naming convention for this. +# +# Arguments: +# posn = (in) The [m] coordinates. +# name_prefix = A prefix string for the coefficient name. +# +# Results: +# The function returns the coefficient, as a Maple string. +# +> coeff_name := +> proc(posn::list(numeric), name_prefix::string) +> cat(name_prefix, sprint_numeric_list(posn)); +> end proc; +coeff_name := proc(posn::list(numeric), name_prefix::string) + cat(name_prefix, sprint_numeric_list(posn)) +end proc + +> +################################################################################ +> +# +# This function computes the name of the variable in which the C code +# will store the input data at a given [m] position, i.e. it encapsulates +# our naming convention for this. +# +# Arguments: +# posn = (in) The [m] coordinates. +# name_prefix = A prefix string for the variable name. +# +# Results: +# The function returns the variable name, as a Maple string. +# +> data_var_name := +> proc(posn::list(numeric), name_prefix::string) +> cat(name_prefix, sprint_numeric_list(posn)); +> end proc; +data_var_name := proc(posn::list(numeric), name_prefix::string) + cat(name_prefix, sprint_numeric_list(posn)) +end proc + +# Maple code to compute lists of point positions in hypercube-shaped molecules +# $Id: $ +> +################################################################################ +> +# +# 1D interpolation points +# +> posn_list_1d_size2 := hypercube_points([ 0], [+1]); + posn_list_1d_size2 := [[0], [1]] + +> posn_list_1d_size3 := hypercube_points([-1], [+1]); + posn_list_1d_size3 := [[-1], [0], [1]] + +> posn_list_1d_size4 := hypercube_points([-1], [+2]); + posn_list_1d_size4 := [[-1], [0], [1], [2]] + +> posn_list_1d_size5 := hypercube_points([-2], [+2]); + posn_list_1d_size5 := [[-2], [-1], [0], [1], [2]] + +> posn_list_1d_size6 := hypercube_points([-2], [+3]); + posn_list_1d_size6 := [[-2], [-1], [0], [1], [2], [3]] + +> posn_list_1d_size7 := hypercube_points([-3], [+3]); + posn_list_1d_size7 := [[-3], [-2], [-1], [0], [1], [2], [3]] + +> +################################################################################ +> +# +# 2D interpolation points (Fortran ordering) +# +> posn_list_2d_size2 := map(ListTools[Reverse], +> hypercube_points([ 0, 0], [+1,+1])); + posn_list_2d_size2 := [[0, 0], [1, 0], [0, 1], [1, 1]] + +> posn_list_2d_size3 := map(ListTools[Reverse], +> hypercube_points([-1,-1], [+1,+1])); +posn_list_2d_size3 := [[-1, -1], [0, -1], [1, -1], [-1, 0], [0, 0], [1, 0], + + [-1, 1], [0, 1], [1, 1]] + +> posn_list_2d_size4 := map(ListTools[Reverse], +> hypercube_points([-1,-1], [+2,+2])); +posn_list_2d_size4 := [[-1, -1], [0, -1], [1, -1], [2, -1], [-1, 0], [0, 0], + + [1, 0], [2, 0], [-1, 1], [0, 1], [1, 1], [2, 1], [-1, 2], [0, 2], [1, 2], + + [2, 2]] + +> posn_list_2d_size5 := map(ListTools[Reverse], +> hypercube_points([-2,-2], [+2,+2])); +posn_list_2d_size5 := [[-2, -2], [-1, -2], [0, -2], [1, -2], [2, -2], [-2, -1], + + [-1, -1], [0, -1], [1, -1], [2, -1], [-2, 0], [-1, 0], [0, 0], [1, 0], + + [2, 0], [-2, 1], [-1, 1], [0, 1], [1, 1], [2, 1], [-2, 2], [-1, 2], [0, 2], + + [1, 2], [2, 2]] + +> +################################################################################ +> +# +# 3D interpolation points (Fortran ordering) +# +> posn_list_3d_size2 := map(ListTools[Reverse], +> hypercube_points([ 0, 0, 0], [+1,+1,+1])); +posn_list_3d_size2 := [[0, 0, 0], [1, 0, 0], [0, 1, 0], [1, 1, 0], [0, 0, 1], + + [1, 0, 1], [0, 1, 1], [1, 1, 1]] + +> posn_list_3d_size3 := map(ListTools[Reverse], +> hypercube_points([-1,-1,-1], [+1,+1,+1])); +posn_list_3d_size3 := [[-1, -1, -1], [0, -1, -1], [1, -1, -1], [-1, 0, -1], + + [0, 0, -1], [1, 0, -1], [-1, 1, -1], [0, 1, -1], [1, 1, -1], [-1, -1, 0], + + [0, -1, 0], [1, -1, 0], [-1, 0, 0], [0, 0, 0], [1, 0, 0], [-1, 1, 0], + + [0, 1, 0], [1, 1, 0], [-1, -1, 1], [0, -1, 1], [1, -1, 1], [-1, 0, 1], + + [0, 0, 1], [1, 0, 1], [-1, 1, 1], [0, 1, 1], [1, 1, 1]] + +> posn_list_3d_size4 := map(ListTools[Reverse], +> hypercube_points([-1,-1,-1], [+2,+2,+2])); +posn_list_3d_size4 := [[-1, -1, -1], [0, -1, -1], [1, -1, -1], [2, -1, -1], + + [-1, 0, -1], [0, 0, -1], [1, 0, -1], [2, 0, -1], [-1, 1, -1], [0, 1, -1], + + [1, 1, -1], [2, 1, -1], [-1, 2, -1], [0, 2, -1], [1, 2, -1], [2, 2, -1], + + [-1, -1, 0], [0, -1, 0], [1, -1, 0], [2, -1, 0], [-1, 0, 0], [0, 0, 0], + + [1, 0, 0], [2, 0, 0], [-1, 1, 0], [0, 1, 0], [1, 1, 0], [2, 1, 0], + + [-1, 2, 0], [0, 2, 0], [1, 2, 0], [2, 2, 0], [-1, -1, 1], [0, -1, 1], + + [1, -1, 1], [2, -1, 1], [-1, 0, 1], [0, 0, 1], [1, 0, 1], [2, 0, 1], + + [-1, 1, 1], [0, 1, 1], [1, 1, 1], [2, 1, 1], [-1, 2, 1], [0, 2, 1], + + [1, 2, 1], [2, 2, 1], [-1, -1, 2], [0, -1, 2], [1, -1, 2], [2, -1, 2], + + [-1, 0, 2], [0, 0, 2], [1, 0, 2], [2, 0, 2], [-1, 1, 2], [0, 1, 2], + + [1, 1, 2], [2, 1, 2], [-1, 2, 2], [0, 2, 2], [1, 2, 2], [2, 2, 2]] + +> posn_list_3d_size5 := map(ListTools[Reverse], +> hypercube_points([-2,-2,-2], [+2,+2,+2])); +posn_list_3d_size5 := [[-2, -2, -2], [-1, -2, -2], [0, -2, -2], [1, -2, -2], + + [2, -2, -2], [-2, -1, -2], [-1, -1, -2], [0, -1, -2], [1, -1, -2], + + [2, -1, -2], [-2, 0, -2], [-1, 0, -2], [0, 0, -2], [1, 0, -2], [2, 0, -2], + + [-2, 1, -2], [-1, 1, -2], [0, 1, -2], [1, 1, -2], [2, 1, -2], [-2, 2, -2], + + [-1, 2, -2], [0, 2, -2], [1, 2, -2], [2, 2, -2], [-2, -2, -1], [-1, -2, -1], + + [0, -2, -1], [1, -2, -1], [2, -2, -1], [-2, -1, -1], [-1, -1, -1], + + [0, -1, -1], [1, -1, -1], [2, -1, -1], [-2, 0, -1], [-1, 0, -1], [0, 0, -1], + + [1, 0, -1], [2, 0, -1], [-2, 1, -1], [-1, 1, -1], [0, 1, -1], [1, 1, -1], + + [2, 1, -1], [-2, 2, -1], [-1, 2, -1], [0, 2, -1], [1, 2, -1], [2, 2, -1], + + [-2, -2, 0], [-1, -2, 0], [0, -2, 0], [1, -2, 0], [2, -2, 0], [-2, -1, 0], + + [-1, -1, 0], [0, -1, 0], [1, -1, 0], [2, -1, 0], [-2, 0, 0], [-1, 0, 0], + + [0, 0, 0], [1, 0, 0], [2, 0, 0], [-2, 1, 0], [-1, 1, 0], [0, 1, 0], + + [1, 1, 0], [2, 1, 0], [-2, 2, 0], [-1, 2, 0], [0, 2, 0], [1, 2, 0], + + [2, 2, 0], [-2, -2, 1], [-1, -2, 1], [0, -2, 1], [1, -2, 1], [2, -2, 1], + + [-2, -1, 1], [-1, -1, 1], [0, -1, 1], [1, -1, 1], [2, -1, 1], [-2, 0, 1], + + [-1, 0, 1], [0, 0, 1], [1, 0, 1], [2, 0, 1], [-2, 1, 1], [-1, 1, 1], + + [0, 1, 1], [1, 1, 1], [2, 1, 1], [-2, 2, 1], [-1, 2, 1], [0, 2, 1], + + [1, 2, 1], [2, 2, 1], [-2, -2, 2], [-1, -2, 2], [0, -2, 2], [1, -2, 2], + + [2, -2, 2], [-2, -1, 2], [-1, -1, 2], [0, -1, 2], [1, -1, 2], [2, -1, 2], + + [-2, 0, 2], [-1, 0, 2], [0, 0, 2], [1, 0, 2], [2, 0, 2], [-2, 1, 2], + + [-1, 1, 2], [0, 1, 2], [1, 1, 2], [2, 1, 2], [-2, 2, 2], [-1, 2, 2], + + [0, 2, 2], [1, 2, 2], [2, 2, 2]] + +# Maple code to define Lagrange interpolating functions/coords/coeffs +# $Id: 1d.maple,v 1.4 2002/05/14 15:54:01 jthorn Exp $ +> +################################################################################ +> +# +# 1-D interpolating functions +# +> +> fn_1d_order1 := +> proc(x) +> + c0 + c1*x +> end proc; + fn_1d_order1 := proc(x) c0 + c1*x end proc + +> +> fn_1d_order2 := +> proc(x) +> + c0 + c1*x + c2*x^2 +> end proc; + fn_1d_order2 := proc(x) c0 + c1*x + c2*x^2 end proc + +> +> fn_1d_order3 := +> proc(x) +> + c0 + c1*x + c2*x^2 + c3*x^3 +> end proc; + fn_1d_order3 := proc(x) c0 + c1*x + c2*x^2 + c3*x^3 end proc + +> +> fn_1d_order4 := +> proc(x) +> + c0 + c1*x + c2*x^2 + c3*x^3 + c4*x^4 +> end; + fn_1d_order4 := proc(x) c0 + c1*x + c2*x^2 + c3*x^3 + c4*x^4 end proc + +> +> fn_1d_order5 := +> proc(x) +> + c0 + c1*x + c2*x^2 + c3*x^3 + c4*x^4 + c5*x^5 +> end; + fn_1d_order5 := proc(x) c0 + c1*x + c2*x^2 + c3*x^3 + c4*x^4 + c5*x^5 end proc + +> +> fn_1d_order6 := +> proc(x) +> + c0 + c1*x + c2*x^2 + c3*x^3 + c4*x^4 + c5*x^5 + c6*x^6 +> end; +fn_1d_order6 := + + proc(x) c0 + c1*x + c2*x^2 + c3*x^3 + c4*x^4 + c5*x^5 + c6*x^6 end proc + +> +######################################## +> +# coordinates for 1-D interpolating functions +> coord_list_1d := [x]; + coord_list_1d := [x] + +> +######################################## +> +# +# coefficients in 1-D interpolating functions +# +> +> coeff_list_1d_order1 := [c0, c1]; + coeff_list_1d_order1 := [c0, c1] + +> coeff_list_1d_order2 := [c0, c1, c2]; + coeff_list_1d_order2 := [c0, c1, c2] + +> coeff_list_1d_order3 := [c0, c1, c2, c3]; + coeff_list_1d_order3 := [c0, c1, c2, c3] + +> coeff_list_1d_order4 := [c0, c1, c2, c3, c4]; + coeff_list_1d_order4 := [c0, c1, c2, c3, c4] + +> coeff_list_1d_order5 := [c0, c1, c2, c3, c4, c5]; + coeff_list_1d_order5 := [c0, c1, c2, c3, c4, c5] + +> coeff_list_1d_order6 := [c0, c1, c2, c3, c4, c5, c6]; + coeff_list_1d_order6 := [c0, c1, c2, c3, c4, c5, c6] + +> +################################################################################ +> +# +# 2-D interpolating functions +# +> +> fn_2d_order1 := +> proc(x,y) +> + c01*y +> + c00 + c10*x +> end proc; + fn_2d_order1 := proc(x, y) c01*y + c00 + c10*x end proc + +> +> fn_2d_order2 := +> proc(x,y) +> + c02*y^2 +> + c01*y + c11*x*y +> + c00 + c10*x + c20*x^2 +> end proc; +fn_2d_order2 := + + proc(x, y) c02*y^2 + c01*y + c11*x*y + c00 + c10*x + c20*x^2 end proc + +> +> fn_2d_order3 := +> proc(x,y) +> + c03*y^3 +> + c02*y^2 + c12*x*y^2 +> + c01*y + c11*x*y + c21*x^2*y +> + c00 + c10*x + c20*x^2 + c30*x^3 +> end proc; +fn_2d_order3 := proc(x, y) + c03*y^3 + c02*y^2 + c12*x*y^2 + c01*y + c11*x*y + c21*x^2*y + c00 + + c10*x + c20*x^2 + c30*x^3 +end proc + +> +> fn_2d_order4 := +> proc(x,y) +> + c04*y^4 +> + c03*y^3 + c13*x*y^3 +> + c02*y^2 + c12*x*y^2 + c22*x^2*y^2 +> + c01*y + c11*x*y + c21*x^2*y + c31*x^3*y +> + c00 + c10*x + c20*x^2 + c30*x^3 + c40*x^4 +> end; +fn_2d_order4 := proc(x, y) + c04*y^4 + c03*y^3 + c13*x*y^3 + c02*y^2 + c12*x*y^2 + c22*x^2*y^2 + + c01*y + c11*x*y + c21*x^2*y + c31*x^3*y + c00 + c10*x + c20*x^2 + + c30*x^3 + c40*x^4 +end proc + +> +######################################## +> +# coordinates for 2-D interpolating functions +> coord_list_2d := [x,y]; + coord_list_2d := [x, y] + +> +######################################## +> +# +# coefficients in 2-D interpolating functions +# +> +> coeff_list_2d_order1 := [ +> c01, +> c00, c10 +> ]; + coeff_list_2d_order1 := [c01, c00, c10] + +> coeff_list_2d_order2 := [ +> c02, +> c01, c11, +> c00, c10, c20 +> ]; + coeff_list_2d_order2 := [c02, c01, c11, c00, c10, c20] + +> coeff_list_2d_order3 := [ +> c03, +> c02, c12, +> c01, c11, c21, +> c00, c10, c20, c30 +> ]; + coeff_list_2d_order3 := [c03, c02, c12, c01, c11, c21, c00, c10, c20, c30] + +> coeff_list_2d_order4 := [ +> c04, +> c03, c13, +> c02, c12, c22, +> c01, c11, c21, c31, +> c00, c10, c20, c30, c40 +> ]; +coeff_list_2d_order4 := + + [c04, c03, c13, c02, c12, c22, c01, c11, c21, c31, c00, c10, c20, c30, c40] + +> +################################################################################ +> +# +# 3-D interpolating functions +# +> +> fn_3d_order1 := +> proc(x,y,z) +# z^0 ----------- +> + c010*y +> + c000 + c100*x +# z^1 ----------- +> + c001*z +> end proc; + fn_3d_order1 := proc(x, y, z) c010*y + c000 + c100*x + c001*z end proc + +> +> fn_3d_order2 := +> proc(x,y,z) +# z^0 -------------------------- +> + c020*y^2 +> + c010*y + c110*x*y +> + c000 + c100*x + c200*x^2 +# z^1 -------------------------- +> + c011*y*z +> + c001*z + c101*x*z +# z^2 -------------------------- +> + c002*z^2 +> end proc; +fn_3d_order2 := proc(x, y, z) + c020*y^2 + c010*y + c110*x*y + c000 + c100*x + c200*x^2 + c011*y*z + + c001*z + c101*x*z + c002*z^2 +end proc + +> +> fn_3d_order3 := +> proc(x,y,z) +# z^0 ------------------------------------------- +> + c030*y^3 +> + c020*y^2 + c120*x*y^2 +> + c010*y + c110*x*y + c210*x^2*y +> + c000 + c100*x + c200*x^2 + c300*x^3 +# z^1 ------------------------------------------- +> + c021*y^2*z +> + c011*y *z + c111*x*y*z +> + c001 *z + c101*x *z + c201*x^2*z +# z^2 ------------------------------------------- +> + c012*y*z^2 +> + c002 *z^2 + c102*x*z^2 +# z^3 ------------------------------------------- +> + c003 *z^3 +> end proc; +fn_3d_order3 := proc(x, y, z) + c030*y^3 + c020*y^2 + c120*x*y^2 + c010*y + c110*x*y + c210*x^2*y + + c000 + c100*x + c200*x^2 + c300*x^3 + c021*y^2*z + c011*y*z + + c111*x*y*z + c001*z + c101*x*z + c201*x^2*z + c012*y*z^2 + c002*z^2 + + c102*x*z^2 + c003*z^3 +end proc + +> +> fn_3d_order4 := +> proc(x,y,z) +# z^0 -------------------------------------------------------- +> + c040*y^4 +> + c030*y^3 + c130*x*y^3 +> + c020*y^2 + c120*x*y^2 + c220*x^2*y^2 +> + c010*y + c110*x*y + c210*x^2*y + c310*x^3*y +> + c000 + c100*x + c200*x^2 + c300*x^3 + c400*x^4 +# z^1 ------------------------------------------- +> + c031*y^3*z +> + c021*y^2*z + c121*x*y^2*z +> + c011*y *z + c111*x*y *z + c211*x^2*y*z +> + c001 *z + c101*x *z + c201*x^2 *z + c301*x^3*z +# z^2 ------------------------------------------- +> + c022*y^2*z^2 +> + c012*y *z^2 + c112*x*y*z^2 +> + c002 *z^2 + c102*x *z^2 + c202*x^2*z^2 +# z^3 ------------------------------------------- +> + c013*y *z^3 +> + c003 *z^3 + c103*x *z^3 +# z^4 ------------------------------------------- +> + c004 *z^4 +> end; +fn_3d_order4 := proc(x, y, z) + c102*x*z^2 + c012*y*z^2 + c111*x*y*z + c121*x*y^2*z + c211*x^2*y*z + + c112*x*y*z^2 + c010*y + c000 + c100*x + c001*z + c020*y^2 + c110*x*y + + c011*y*z + c101*x*z + c200*x^2 + c002*z^2 + c030*y^3 + c300*x^3 + + c003*z^3 + c040*y^4 + c400*x^4 + c004*z^4 + c120*x*y^2 + c210*x^2*y + + c021*y^2*z + c201*x^2*z + c130*x*y^3 + c220*x^2*y^2 + c310*x^3*y + + c031*y^3*z + c301*x^3*z + c022*y^2*z^2 + c202*x^2*z^2 + c013*y*z^3 + + c103*x*z^3 +end proc + +> +######################################## +> +# coordinates for 3-D interpolating functions +> coord_list_3d := [x,y,z]; + coord_list_3d := [x, y, z] + +> +######################################## +> +# +# coefficients in 3-D interpolating functions +# +> +> coeff_list_3d_order1 := [ +> # z^0 ----- +> c010, +> c000, c100, +> # z^1 ----- +> c001 +> ]; + coeff_list_3d_order1 := [c010, c000, c100, c001] + +> coeff_list_3d_order2 := [ +> # z^0 ----------- +> c020, +> c010, c110, +> c000, c100, c200, +> # z^1 ----------- +> c011, +> c001, c101, +> # z^2 ----------- +> c002 +> ]; +coeff_list_3d_order2 := + + [c020, c010, c110, c000, c100, c200, c011, c001, c101, c002] + +> coeff_list_3d_order3 := [ +> # z^0 ---------------- +> c030, +> c020, c120, +> c010, c110, c210, +> c000, c100, c200, c300, +> # z^1 ---------------- +> c021, +> c011, c111, +> c001, c101, c201, +> # z^2 ---------------- +> c012, +> c002, c102, +> # z^3 ---------------- +> c003 +> ]; +coeff_list_3d_order3 := [c030, c020, c120, c010, c110, c210, c000, c100, c200, + + c300, c021, c011, c111, c001, c101, c201, c012, c002, c102, c003] + +> coeff_list_3d_order4 := [ +> # z^0 ----------------------- +> c040, +> c030, c130, +> c020, c120, c220, +> c010, c110, c210, c310, +> c000, c100, c200, c300, c400, +> # z^1 ----------------------- +> c031, +> c021, c121, +> c011, c111, c211, +> c001, c101, c201, c301, +> # z^2 ----------------------- +> c022, +> c012, c112, +> c002, c102, c202, +> # z^3 ----------------------- +> c013, +> c003, c103, +> # z^4 ----------------------- +> c004 +> ]; +coeff_list_3d_order4 := [c040, c030, c130, c020, c120, c220, c010, c110, c210, + + c310, c000, c100, c200, c300, c400, c031, c021, c121, c011, c111, c211, + + c001, c101, c201, c301, c022, c012, c112, c002, c102, c202, c013, c003, + + c103, c004] + +> +################################################################################ +# 2d.maple -- compute Lagrange interpolation coefficients in 2-D +# $Id: $ +> +################################################################################ +> +# +# 2d, cube, order=1, smoothing=0 (size=2) +# +> +# interpolating polynomial +> interp_2d_cube_order1_smooth0 +> := Lagrange_polynomial_interpolant(fn_2d_order1, coeff_list_2d_order1, +> coord_list_2d, posn_list_2d_size2); +interp_2d_cube_order1_smooth0 := + + (1/2 DATA(1, 1) - 1/2 DATA(0, 0) - 1/2 DATA(1, 0) + 1/2 DATA(0, 1)) y + + - 1/4 DATA(1, 1) + 3/4 DATA(0, 0) + 1/4 DATA(1, 0) + 1/4 DATA(0, 1) + + + (1/2 DATA(1, 1) - 1/2 DATA(0, 0) + 1/2 DATA(1, 0) - 1/2 DATA(0, 1)) x + +> +# I +> coeff_as_lc_of_data(%, posn_list_2d_size2); +[COEFF(0, 0) = 3/4 - 1/2 y - 1/2 x, COEFF(1, 0) = 1/4 - 1/2 y + 1/2 x, + + COEFF(0, 1) = - 1/2 x + 1/4 + 1/2 y, COEFF(1, 1) = 1/2 y - 1/4 + 1/2 x] + +> print_coeff__lc_of_data(%, "coeff_I_", "fp", +> "2d.coeffs/2d.cube.order1.smooth0/coeff-I.compute.c"); +bytes used=1000160, alloc=917336, time=0.07 +> +# d/dx +> simplify( diff(interp_2d_cube_order1_smooth0,x) ); + 1/2 DATA(1, 1) - 1/2 DATA(0, 0) + 1/2 DATA(1, 0) - 1/2 DATA(0, 1) + +> coeff_as_lc_of_data(%, posn_list_2d_size2); + [COEFF(0, 0) = -1/2, COEFF(1, 0) = 1/2, COEFF(0, 1) = -1/2, COEFF(1, 1) = 1/2] + +> print_coeff__lc_of_data(%, "coeff_dx_", "fp", +> "2d.coeffs/2d.cube.order1.smooth0/coeff-dx.compute.c"); +> +# d/dy +> simplify( diff(interp_2d_cube_order1_smooth0,y) ); + 1/2 DATA(1, 1) - 1/2 DATA(0, 0) - 1/2 DATA(1, 0) + 1/2 DATA(0, 1) + +> coeff_as_lc_of_data(%, posn_list_2d_size2); + [COEFF(0, 0) = -1/2, COEFF(1, 0) = -1/2, COEFF(0, 1) = 1/2, COEFF(1, 1) = 1/2] + +> print_coeff__lc_of_data(%, "coeff_dy_", "fp", +> "2d.coeffs/2d.cube.order1.smooth0/coeff-dy.compute.c"); +bytes used=2000328, alloc=1376004, time=0.15 +> +################################################################################ +> +# +# 2d, cube, order=2, smoothing=0 (size=3) +# +> +# interpolating polynomial +> interp_2d_cube_order2_smooth0 +> := Lagrange_polynomial_interpolant(fn_2d_order2, coeff_list_2d_order2, +> coord_list_2d, posn_list_2d_size3); +interp_2d_cube_order2_smooth0 := (1/6 DATA(0, -1) + 1/6 DATA(-1, -1) + + - 1/3 DATA(-1, 0) - 1/3 DATA(1, 0) + 1/6 DATA(1, -1) + 1/6 DATA(-1, 1) + + 2 + + 1/6 DATA(1, 1) - 1/3 DATA(0, 0) + 1/6 DATA(0, 1)) y + (1/6 DATA(-1, 1) + + - 1/6 DATA(1, -1) - 1/6 DATA(0, -1) - 1/6 DATA(-1, -1) + 1/6 DATA(0, 1) + + + 1/6 DATA(1, 1)) y + + + (1/4 DATA(1, 1) - 1/4 DATA(1, -1) + 1/4 DATA(-1, -1) - 1/4 DATA(-1, 1)) x y + + + 2/9 DATA(0, -1) - 1/9 DATA(-1, -1) + 2/9 DATA(-1, 0) + 2/9 DATA(1, 0) + + - 1/9 DATA(1, -1) - 1/9 DATA(-1, 1) - 1/9 DATA(1, 1) + 5/9 DATA(0, 0) + + + 2/9 DATA(0, 1) + (- 1/6 DATA(-1, -1) - 1/6 DATA(-1, 0) + 1/6 DATA(1, 0) + + + 1/6 DATA(1, -1) - 1/6 DATA(-1, 1) + 1/6 DATA(1, 1)) x + ( + + 1/6 DATA(-1, -1) + 1/6 DATA(-1, 0) + 1/6 DATA(1, 0) + 1/6 DATA(1, -1) + + - 1/3 DATA(0, -1) + 1/6 DATA(-1, 1) + 1/6 DATA(1, 1) - 1/3 DATA(0, 0) + + 2 + - 1/3 DATA(0, 1)) x + +> +# I +> coeff_as_lc_of_data(%, posn_list_2d_size3); + 2 2 +[COEFF(-1, -1) = - 1/6 y + 1/6 y - 1/9 - 1/6 x + 1/4 x y + 1/6 x , + + 2 2 + COEFF(0, -1) = - 1/6 y + 2/9 + 1/6 y - 1/3 x , + + 2 2 + COEFF(1, -1) = - 1/6 y - 1/4 x y + 1/6 x + 1/6 y - 1/9 + 1/6 x , + + 2 2 + COEFF(-1, 0) = - 1/3 y + 2/9 + 1/6 x - 1/6 x, + + 2 2 + COEFF(0, 0) = - 1/3 y - 1/3 x + 5/9, + + 2 2 + COEFF(1, 0) = 2/9 - 1/3 y + 1/6 x + 1/6 x, + + 2 2 + COEFF(-1, 1) = 1/6 y + 1/6 y + 1/6 x - 1/4 x y - 1/6 x - 1/9, + + 2 2 + COEFF(0, 1) = 1/6 y - 1/3 x + 2/9 + 1/6 y , + + 2 2 + COEFF(1, 1) = 1/6 y + 1/6 y + 1/6 x - 1/9 + 1/4 x y + 1/6 x] + +> print_coeff__lc_of_data(%, "coeff_I_", "fp", +> "2d.coeffs/2d.cube.order2.smooth0/coeff-I.compute.c"); +bytes used=3000520, alloc=1769148, time=0.22 +bytes used=4002236, alloc=1834672, time=0.30 +> +# d/dx +> simplify( diff(interp_2d_cube_order2_smooth0,x) ); +1/4 y DATA(1, 1) - 1/4 y DATA(1, -1) + 1/4 y DATA(-1, -1) - 1/4 y DATA(-1, 1) + + - 1/6 DATA(-1, -1) - 1/6 DATA(-1, 0) + 1/6 DATA(1, 0) + 1/6 DATA(1, -1) + + - 1/6 DATA(-1, 1) + 1/6 DATA(1, 1) + 1/3 x DATA(-1, -1) + + + 1/3 x DATA(-1, 0) + 1/3 x DATA(1, 0) + 1/3 x DATA(1, -1) + + - 2/3 x DATA(0, -1) + 1/3 x DATA(-1, 1) + 1/3 x DATA(1, 1) + + - 2/3 x DATA(0, 0) - 2/3 x DATA(0, 1) + +> coeff_as_lc_of_data(%, posn_list_2d_size3); +[COEFF(-1, -1) = 1/4 y - 1/6 + 1/3 x, COEFF(0, -1) = - 2/3 x, + + COEFF(1, -1) = 1/6 + 1/3 x - 1/4 y, COEFF(-1, 0) = - 1/6 + 1/3 x, + + COEFF(0, 0) = - 2/3 x, COEFF(1, 0) = 1/6 + 1/3 x, + + COEFF(-1, 1) = - 1/4 y - 1/6 + 1/3 x, COEFF(0, 1) = - 2/3 x, + + COEFF(1, 1) = 1/4 y + 1/3 x + 1/6] + +> print_coeff__lc_of_data(%, "coeff_dx_", "fp", +> "2d.coeffs/2d.cube.order2.smooth0/coeff-dx.compute.c"); +bytes used=5002436, alloc=1834672, time=0.38 +> +# d/dy +> simplify( diff(interp_2d_cube_order2_smooth0,y) ); +1/3 y DATA(0, -1) + 1/3 y DATA(-1, -1) - 2/3 y DATA(-1, 0) - 2/3 y DATA(1, 0) + + + 1/3 y DATA(1, -1) + 1/3 y DATA(-1, 1) + 1/3 y DATA(1, 1) + + - 2/3 y DATA(0, 0) + 1/3 y DATA(0, 1) + 1/6 DATA(-1, 1) - 1/6 DATA(1, -1) + + - 1/6 DATA(0, -1) - 1/6 DATA(-1, -1) + 1/6 DATA(0, 1) + 1/6 DATA(1, 1) + + + 1/4 x DATA(1, 1) - 1/4 x DATA(1, -1) + 1/4 x DATA(-1, -1) + + - 1/4 x DATA(-1, 1) + +> coeff_as_lc_of_data(%, posn_list_2d_size3); +[COEFF(-1, -1) = - 1/6 + 1/4 x + 1/3 y, COEFF(0, -1) = 1/3 y - 1/6, + + COEFF(1, -1) = 1/3 y - 1/6 - 1/4 x, COEFF(-1, 0) = - 2/3 y, + + COEFF(0, 0) = - 2/3 y, COEFF(1, 0) = - 2/3 y, + + COEFF(-1, 1) = 1/3 y - 1/4 x + 1/6, COEFF(0, 1) = 1/6 + 1/3 y, + + COEFF(1, 1) = 1/6 + 1/3 y + 1/4 x] + +> print_coeff__lc_of_data(%, "coeff_dy_", "fp", +> "2d.coeffs/2d.cube.order2.smooth0/coeff-dy.compute.c"); +> +# d^2/dx^2 +> simplify( diff(interp_2d_cube_order2_smooth0,x,x) ); +1/3 DATA(-1, -1) + 1/3 DATA(-1, 0) + 1/3 DATA(1, 0) + 1/3 DATA(1, -1) + + - 2/3 DATA(0, -1) + 1/3 DATA(-1, 1) + 1/3 DATA(1, 1) - 2/3 DATA(0, 0) + + - 2/3 DATA(0, 1) + +> coeff_as_lc_of_data(%, posn_list_2d_size3); +bytes used=6002660, alloc=1834672, time=0.47 +[COEFF(-1, -1) = 1/3, COEFF(0, -1) = -2/3, COEFF(1, -1) = 1/3, + + COEFF(-1, 0) = 1/3, COEFF(0, 0) = -2/3, COEFF(1, 0) = 1/3, + + COEFF(-1, 1) = 1/3, COEFF(0, 1) = -2/3, COEFF(1, 1) = 1/3] + +> print_coeff__lc_of_data(%, "coeff_dxx_", "fp", +> "2d.coeffs/2d.cube.order2.smooth0/coeff-dxx.compute.c"); +> +# d^2/dxdy +> simplify( diff(interp_2d_cube_order2_smooth0,x,y) ); + 1/4 DATA(1, 1) - 1/4 DATA(1, -1) + 1/4 DATA(-1, -1) - 1/4 DATA(-1, 1) + +> coeff_as_lc_of_data(%, posn_list_2d_size3); +[COEFF(-1, -1) = 1/4, COEFF(0, -1) = 0, COEFF(1, -1) = -1/4, COEFF(-1, 0) = 0, + + COEFF(0, 0) = 0, COEFF(1, 0) = 0, COEFF(-1, 1) = -1/4, COEFF(0, 1) = 0, + + COEFF(1, 1) = 1/4] + +> print_coeff__lc_of_data(%, "coeff_dxy_", "fp", +> "2d.coeffs/2d.cube.order2.smooth0/coeff-dxy.compute.c"); +> +# d^2/dy^2 +> simplify( diff(interp_2d_cube_order2_smooth0,y,y) ); +1/3 DATA(0, -1) + 1/3 DATA(-1, -1) - 2/3 DATA(-1, 0) - 2/3 DATA(1, 0) + + + 1/3 DATA(1, -1) + 1/3 DATA(-1, 1) + 1/3 DATA(1, 1) - 2/3 DATA(0, 0) + + + 1/3 DATA(0, 1) + +> coeff_as_lc_of_data(%, posn_list_2d_size3); +[COEFF(-1, -1) = 1/3, COEFF(0, -1) = 1/3, COEFF(1, -1) = 1/3, + + COEFF(-1, 0) = -2/3, COEFF(0, 0) = -2/3, COEFF(1, 0) = -2/3, + + COEFF(-1, 1) = 1/3, COEFF(0, 1) = 1/3, COEFF(1, 1) = 1/3] + +> print_coeff__lc_of_data(%, "coeff_dyy_", "fp", +> "2d.coeffs/2d.cube.order2.smooth0/coeff-dyy.compute.c"); +> +################################################################################ +> +# +# 2d, cube, order=3, smoothing=0 (size=4) +# +> +# interpolating polynomial +> interp_2d_cube_order3_smooth0 +> := Lagrange_polynomial_interpolant(fn_2d_order3, coeff_list_2d_order3, +> coord_list_2d, posn_list_2d_size4); +bytes used=7003016, alloc=1834672, time=0.54 +bytes used=8003280, alloc=1900196, time=0.60 +bytes used=9003556, alloc=1965720, time=0.67 +bytes used=10003800, alloc=1965720, time=0.74 +bytes used=11003984, alloc=1965720, time=0.80 +interp_2d_cube_order3_smooth0 := (1/40 DATA(0, 0) - 1/40 DATA(1, 0) + + + 1/40 DATA(0, 1) - 1/40 DATA(1, 1) + 3/40 DATA(-1, 0) + 1/40 DATA(1, -1) + + - 1/40 DATA(0, -1) - 3/40 DATA(-1, -1) + 3/40 DATA(-1, 1) + + + 1/40 DATA(1, 2) + 3/40 DATA(2, 2) - 3/40 DATA(2, 0) - 3/40 DATA(2, 1) + + 2 / + - 3/40 DATA(-1, 2) - 1/40 DATA(0, 2) + 3/40 DATA(2, -1)) x y + | + \ + + - 1/25 DATA(0, 0) - 1/100 DATA(1, 0) - 1/100 DATA(0, 1) + 3/50 DATA(1, 1) + + 13 + - 1/50 DATA(-1, 0) - --- DATA(1, -1) - 1/50 DATA(0, -1) + 100 + + 13 + + 6/25 DATA(-1, -1) - --- DATA(-1, 1) + 2/25 DATA(1, 2) - 3/50 DATA(2, 2) + 100 + + + 7/100 DATA(2, 0) + 2/25 DATA(2, 1) - 9/100 DATA(-1, 2) + + \ + + 7/100 DATA(0, 2) - 9/100 DATA(2, -1)| x y + (1/40 DATA(0, 0) + / + + + 1/40 DATA(1, 0) - 1/40 DATA(0, 1) - 1/40 DATA(1, 1) - 1/40 DATA(-1, 0) + + + 3/40 DATA(1, -1) + 3/40 DATA(0, -1) - 3/40 DATA(-1, -1) + + + 1/40 DATA(-1, 1) - 3/40 DATA(1, 2) + 3/40 DATA(2, 2) - 1/40 DATA(2, 0) + + + 1/40 DATA(2, 1) + 3/40 DATA(-1, 2) - 3/40 DATA(0, 2) - 3/40 DATA(2, -1)) + + 2 93 37 37 17 + x y + --- DATA(0, 0) + --- DATA(1, 0) + --- DATA(0, 1) - --- DATA(1, 1) + 200 200 200 200 + + 11 11 23 + + -- DATA(-1, 0) - 1/50 DATA(1, -1) + -- DATA(0, -1) - --- DATA(-1, -1) + 50 50 200 + + - 1/50 DATA(-1, 1) + (1/8 DATA(0, 0) + 1/8 DATA(1, 0) - 1/8 DATA(0, 1) + + - 1/8 DATA(1, 1) + 1/8 DATA(-1, 0) - 1/24 DATA(1, -1) - 1/24 DATA(0, -1) + + - 1/24 DATA(-1, -1) - 1/8 DATA(-1, 1) + 1/24 DATA(1, 2) + 1/24 DATA(2, 2) + + + 1/8 DATA(2, 0) - 1/8 DATA(2, 1) + 1/24 DATA(-1, 2) + 1/24 DATA(0, 2) + + 3 / 21 19 + - 1/24 DATA(2, -1)) y + |- -- DATA(0, 0) - -- DATA(1, 0) + \ 80 80 + + 11 23 + + 9/80 DATA(0, 1) + -- DATA(1, 1) - -- DATA(-1, 0) + 9/80 DATA(1, -1) + 80 80 + + 11 13 + + -- DATA(0, -1) + -- DATA(-1, -1) + 7/80 DATA(-1, 1) - 1/80 DATA(1, 2) + 80 80 + + 17 13 + - 3/80 DATA(2, 2) - -- DATA(2, 0) + -- DATA(2, 1) + 3/80 DATA(-1, 2) + 80 80 + + \ 2 / 57 + + 1/80 DATA(0, 2) + 7/80 DATA(2, -1)| y + |- --- DATA(0, 0) + / \ 400 + + 63 117 103 31 + - --- DATA(1, 0) + --- DATA(0, 1) + --- DATA(1, 1) - --- DATA(-1, 0) + 400 400 400 400 + + 157 223 109 + - ---- DATA(1, -1) - ---- DATA(0, -1) - ---- DATA(-1, -1) + 1200 1200 1200 + + 111 37 149 49 + + --- DATA(-1, 1) + ---- DATA(1, 2) - ---- DATA(2, 2) - --- DATA(2, 0) + 400 1200 1200 400 + + 69 131 43 89 \ + + --- DATA(2, 1) - ---- DATA(-1, 2) + ---- DATA(0, 2) + ---- DATA(2, -1)| + 400 1200 1200 1200 / + + / 57 117 63 103 + y + |- --- DATA(0, 0) + --- DATA(1, 0) - --- DATA(0, 1) + --- DATA(1, 1) + \ 400 400 400 400 + + 223 111 31 109 + - ---- DATA(-1, 0) + --- DATA(1, -1) - --- DATA(0, -1) - ---- DATA(-1, -1) + 1200 400 400 1200 + + 157 69 149 43 + - ---- DATA(-1, 1) + --- DATA(1, 2) - ---- DATA(2, 2) + ---- DATA(2, 0) + 1200 400 1200 1200 + + 37 89 49 131 \ + + ---- DATA(2, 1) + ---- DATA(-1, 2) - --- DATA(0, 2) - ---- DATA(2, -1)| + 1200 1200 400 1200 / + + / 21 19 11 + x + |- -- DATA(0, 0) + 9/80 DATA(1, 0) - -- DATA(0, 1) + -- DATA(1, 1) + \ 80 80 80 + + 11 23 13 + + -- DATA(-1, 0) + 7/80 DATA(1, -1) - -- DATA(0, -1) + -- DATA(-1, -1) + 80 80 80 + + 13 + + 9/80 DATA(-1, 1) + -- DATA(1, 2) - 3/80 DATA(2, 2) + 1/80 DATA(2, 0) + 80 + + 17 \ + - 1/80 DATA(2, 1) + 7/80 DATA(-1, 2) - -- DATA(0, 2) + 3/80 DATA(2, -1)| + 80 / + + 2 + x + (1/8 DATA(0, 0) - 1/8 DATA(1, 0) + 1/8 DATA(0, 1) - 1/8 DATA(1, 1) + + - 1/24 DATA(-1, 0) - 1/8 DATA(1, -1) + 1/8 DATA(0, -1) - 1/24 DATA(-1, -1) + + - 1/24 DATA(-1, 1) + 1/24 DATA(2, 1) - 1/8 DATA(1, 2) + 1/24 DATA(2, 2) + + + 1/24 DATA(2, 0) - 1/24 DATA(-1, 2) + 1/8 DATA(0, 2) + 1/24 DATA(2, -1)) + + 3 13 + x - 2/25 DATA(1, 2) + 7/200 DATA(2, 2) + --- DATA(2, 0) - 2/25 DATA(2, 1) + 100 + + 17 13 17 + - --- DATA(-1, 2) + --- DATA(0, 2) - --- DATA(2, -1) + 200 100 200 + +> +# I +> coeff_as_lc_of_data(%, posn_list_2d_size4); +bytes used=12006624, alloc=1965720, time=0.86 + 13 2 3 2 3 13 2 +[COEFF(-1, -1) = -- y - 1/24 y - 3/40 x y - 1/24 x + 6/25 x y + -- x + 80 80 + + 2 109 23 109 2 3 + - 3/40 x y - ---- x - --- - ---- y, COEFF(0, -1) = - 1/40 x y + 1/8 x + 1200 200 1200 + + 2 23 2 11 3 31 11 2 223 + - 1/50 x y + 3/40 x y - -- x + -- - 1/24 y - --- x + -- y - ---- y, + 80 50 400 80 1200 + + 2 2 2 3 111 + COEFF(1, -1) = 9/80 y + 3/40 x y + 1/40 x y - 1/24 y - 1/50 + --- x + 400 + + 157 2 13 3 + - ---- y + 7/80 x - --- x y - 1/8 x , COEFF(2, -1) = - 9/100 x y + 1200 100 + + 2 2 3 89 17 2 3 + + 3/40 x y + 7/80 y - 1/24 y + ---- y - --- - 3/40 x y + 1/24 x + 1200 200 + + 131 2 11 31 2 + - ---- x + 3/80 x , COEFF(-1, 0) = -- - 1/50 x y - --- y + 3/40 x y + 1200 50 400 + + 223 3 23 2 2 11 2 3 + - ---- x - 1/24 x - -- y - 1/40 x y + -- x + 1/8 y , COEFF(0, 0) = + 1200 80 80 + + 21 2 21 2 3 3 57 57 2 + - 1/25 x y - -- y - -- x + 1/8 x + 1/8 y - --- x - --- y + 1/40 x y + 80 80 400 400 + + 2 93 19 2 63 2 + + 1/40 x y + ---, COEFF(1, 0) = - -- y - --- y - 1/40 x y - 1/100 x y + 200 80 400 + + 3 2 117 37 3 2 + + 1/8 y + 9/80 x + --- x + --- - 1/8 x + 1/40 x y, COEFF(2, 0) = + 400 200 + + 17 2 2 49 3 2 3 2 + - -- y + 1/80 x - --- y + 1/8 y - 1/40 x y + 1/24 x - 3/40 x y + 80 400 + + 13 43 2 157 + + 7/100 x y + --- + ---- x, COEFF(-1, 1) = 1/40 x y - 1/50 - ---- x + 100 1200 1200 + + 111 3 3 2 13 2 2 + + --- y - 1/24 x - 1/8 y + 3/40 x y - --- x y + 7/80 y + 9/80 x , + 400 100 + + 2 37 3 63 19 2 3 + COEFF(0, 1) = - 1/40 x y + --- + 1/8 x - --- x - -- x - 1/8 y + 200 400 80 + + 2 117 2 103 2 + + 1/40 x y + --- y + 9/80 y - 1/100 x y, COEFF(1, 1) = --- y - 1/40 x y + 400 400 + + 11 2 11 2 103 3 2 17 3 + + -- x + -- y + 3/50 x y + --- x - 1/8 x - 1/40 x y - --- - 1/8 y , + 80 80 400 200 + + 3 2 2 37 3 2 + COEFF(2, 1) = - 1/8 y + 1/40 x y - 3/40 x y + ---- x + 1/24 x - 1/80 x + 1200 + + 13 2 69 2 89 + - 2/25 + -- y + --- y + 2/25 x y, COEFF(-1, 2) = 3/80 y + ---- x + 80 400 1200 + + 2 131 2 2 3 3 + + 3/40 x y - ---- y - 3/40 x y + 7/80 x + 1/24 y - 9/100 x y - 1/24 x + 1200 + + 17 49 3 43 2 2 + - ---, COEFF(0, 2) = - --- x + 1/24 y + ---- y - 3/40 x y + 1/80 y + 200 400 1200 + + 17 2 2 13 3 69 + - -- x - 1/40 x y + --- + 7/100 x y + 1/8 x , COEFF(1, 2) = --- x + 80 100 400 + + 2 2 2 37 3 13 2 + - 3/40 x y - 1/80 y + 1/40 x y + ---- y - 2/25 + 1/24 y + -- x + 1200 80 + + 3 2 3 3 2 + + 2/25 x y - 1/8 x , COEFF(2, 2) = 3/40 x y + 1/24 x + 1/24 y - 3/80 x + + 149 149 2 2 + - 3/50 x y - ---- y - ---- x + 3/40 x y - 3/80 y + 7/200] + 1200 1200 + +> print_coeff__lc_of_data(%, "coeff_I_", "fp", +> "2d.coeffs/2d.cube.order3.smooth0/coeff-I.compute.c"); +bytes used=13006948, alloc=2031244, time=0.93 +bytes used=14011252, alloc=2031244, time=1.01 +bytes used=15011472, alloc=2031244, time=1.11 +bytes used=16011652, alloc=2031244, time=1.19 +bytes used=17011812, alloc=2031244, time=1.29 +bytes used=18012000, alloc=2031244, time=1.39 +bytes used=19012216, alloc=2031244, time=1.50 +bytes used=20012436, alloc=2031244, time=1.62 +> +# d/dx +> simplify( diff(interp_2d_cube_order3_smooth0,x) ); +bytes used=21013024, alloc=2031244, time=1.73 +3/20 x y DATA(-1, 2) + 1/20 x y DATA(0, 0) + 1/20 x y DATA(1, 0) + + - 1/20 x y DATA(0, 1) - 1/20 x y DATA(1, 1) - 1/20 x y DATA(-1, 0) + + + 3/20 x y DATA(1, -1) + 3/20 x y DATA(0, -1) - 3/20 x y DATA(-1, -1) + + 57 117 63 103 + - --- DATA(0, 0) + --- DATA(1, 0) - --- DATA(0, 1) + --- DATA(1, 1) + 400 400 400 400 + + 2 + + 1/20 x y DATA(-1, 1) - 3/20 x y DATA(1, 2) - 1/8 x DATA(-1, 2) + + 2 2 2 + + 3/8 x DATA(0, 2) + 1/8 x DATA(2, -1) + 1/8 x DATA(2, 1) + + 2 + + 1/8 x DATA(2, 0) - 3/20 x y DATA(2, -1) - 3/20 x y DATA(0, 2) + + + 1/20 x y DATA(2, 1) - 1/20 x y DATA(2, 0) + 3/20 x y DATA(2, 2) + + + 2/25 y DATA(2, 1) + 7/100 y DATA(2, 0) + 2/25 y DATA(1, 2) + + - 9/100 y DATA(2, -1) + 7/100 y DATA(0, 2) - 9/100 y DATA(-1, 2) + + 13 + - 1/50 y DATA(0, -1) - --- y DATA(1, -1) - 1/50 y DATA(-1, 0) + 100 + + 2 13 + - 3/50 y DATA(2, 2) + 3/40 y DATA(2, -1) - --- y DATA(-1, 1) + 100 + + 2 + + 6/25 y DATA(-1, -1) + 1/40 y DATA(1, 2) - 1/100 y DATA(0, 1) + + - 1/100 y DATA(1, 0) - 1/25 y DATA(0, 0) + 3/50 y DATA(1, 1) + + 2 2 2 + - 1/40 y DATA(0, 2) - 3/40 y DATA(2, 1) - 3/40 y DATA(2, 0) + + 2 2 2 + + 3/40 y DATA(2, 2) - 3/40 y DATA(-1, 2) + 3/40 y DATA(-1, 1) + + 2 2 2 + - 1/40 y DATA(0, -1) + 1/40 y DATA(1, -1) + 3/40 y DATA(-1, 0) + + 2 2 2 + - 3/40 y DATA(-1, -1) + 1/40 y DATA(0, 1) - 1/40 y DATA(1, 0) + + 2 11 + + 1/40 y DATA(0, 0) + -- x DATA(-1, 0) + 9/40 x DATA(1, 0) + 40 + + 19 2 21 13 + - -- x DATA(0, 1) - 1/40 y DATA(1, 1) - -- x DATA(0, 0) + -- x DATA(1, 2) + 40 40 40 + + 23 13 11 + - -- x DATA(0, -1) + -- x DATA(-1, -1) + -- x DATA(1, 1) + 40 40 40 + + + 7/40 x DATA(1, -1) + 9/40 x DATA(-1, 1) - 3/40 x DATA(2, 2) + + + 1/40 x DATA(2, 0) - 1/40 x DATA(2, 1) + 7/40 x DATA(-1, 2) + + 2 2 17 + - 3/8 x DATA(1, 0) + 3/8 x DATA(0, 1) - -- x DATA(0, 2) + 40 + + 2 2 + + 3/40 x DATA(2, -1) - 3/8 x DATA(1, -1) + 3/8 x DATA(0, -1) + + 2 2 2 + + 3/8 x DATA(0, 0) - 3/8 x DATA(1, 1) - 1/8 x DATA(-1, 0) + + 2 2 2 + - 3/8 x DATA(1, 2) - 1/8 x DATA(-1, -1) - 1/8 x DATA(-1, 1) + + 2 223 111 31 + + 1/8 x DATA(2, 2) - ---- DATA(-1, 0) + --- DATA(1, -1) - --- DATA(0, -1) + 1200 400 400 + + 109 157 69 149 + - ---- DATA(-1, -1) - ---- DATA(-1, 1) + --- DATA(1, 2) - ---- DATA(2, 2) + 1200 1200 400 1200 + + 43 37 89 49 + + ---- DATA(2, 0) + ---- DATA(2, 1) + ---- DATA(-1, 2) - --- DATA(0, 2) + 1200 1200 1200 400 + + 131 + - ---- DATA(2, -1) + 1200 + +> coeff_as_lc_of_data(%, posn_list_2d_size4); +bytes used=22015264, alloc=2031244, time=1.80 + 13 2 109 2 +[COEFF(-1, -1) = -- x - 1/8 x - ---- + 6/25 y - 3/20 x y - 3/40 y , + 40 1200 + + 31 2 2 23 + COEFF(0, -1) = - --- + 3/8 x - 1/40 y - -- x - 1/50 y + 3/20 x y, + 400 40 + + 2 111 2 13 + COEFF(1, -1) = 3/20 x y + 1/40 y + 7/40 x + --- - 3/8 x - --- y, + 400 100 + + 2 131 2 + COEFF(2, -1) = 1/8 x - 9/100 y - 3/20 x y - ---- + 3/40 y + 3/40 x, + 1200 + + 223 2 11 2 + COEFF(-1, 0) = - ---- + 3/40 y - 1/20 x y - 1/50 y + -- x - 1/8 x , + 1200 40 + + 57 21 2 2 + COEFF(0, 0) = - --- - -- x + 3/8 x - 1/25 y + 1/40 y + 1/20 x y, + 400 40 + + 117 2 2 + COEFF(1, 0) = 9/40 x + --- - 1/40 y - 1/100 y - 3/8 x + 1/20 x y, + 400 + + 2 43 2 + COEFF(2, 0) = 1/8 x + ---- - 3/40 y + 1/40 x - 1/20 x y + 7/100 y, + 1200 + + 2 157 13 2 + COEFF(-1, 1) = 1/20 x y - 1/8 x + 9/40 x - ---- - --- y + 3/40 y , + 1200 100 + + 63 2 19 2 + COEFF(0, 1) = - --- - 1/100 y + 3/8 x - 1/20 x y - -- x + 1/40 y , + 400 40 + + 11 2 2 103 + COEFF(1, 1) = - 1/20 x y + -- x - 1/40 y - 3/8 x + 3/50 y + ---, + 40 400 + + 2 2 37 + COEFF(2, 1) = - 1/40 x + 1/20 x y + 1/8 x - 3/40 y + ---- + 2/25 y, + 1200 + + 89 2 2 + COEFF(-1, 2) = - 9/100 y + 3/20 x y + ---- - 3/40 y + 7/40 x - 1/8 x , + 1200 + + 2 2 17 49 + COEFF(0, 2) = - 1/40 y - 3/20 x y + 7/100 y + 3/8 x - -- x - ---, + 40 400 + + 2 13 69 2 + COEFF(1, 2) = - 3/8 x + -- x + --- - 3/20 x y + 2/25 y + 1/40 y , + 40 400 + + 149 2 2 + COEFF(2, 2) = 3/20 x y - ---- + 1/8 x - 3/40 x + 3/40 y - 3/50 y] + 1200 + +> print_coeff__lc_of_data(%, "coeff_dx_", "fp", +> "2d.coeffs/2d.cube.order3.smooth0/coeff-dx.compute.c"); +bytes used=23015460, alloc=2031244, time=1.87 +bytes used=24015620, alloc=2031244, time=1.96 +bytes used=25016144, alloc=2096768, time=2.07 +bytes used=26016296, alloc=2096768, time=2.18 +> +# d/dy +> simplify( diff(interp_2d_cube_order3_smooth0,y) ); +bytes used=27016596, alloc=2096768, time=2.28 +- 3/20 x y DATA(-1, 2) + 1/20 x y DATA(0, 0) - 1/20 x y DATA(1, 0) + + + 1/20 x y DATA(0, 1) - 1/20 x y DATA(1, 1) + 3/20 x y DATA(-1, 0) + + + 1/20 x y DATA(1, -1) - 1/20 x y DATA(0, -1) - 3/20 x y DATA(-1, -1) + + 57 63 117 103 + - --- DATA(0, 0) - --- DATA(1, 0) + --- DATA(0, 1) + --- DATA(1, 1) + 400 400 400 400 + + 2 + + 3/20 x y DATA(-1, 1) + 1/20 x y DATA(1, 2) + 3/40 x DATA(-1, 2) + + 2 2 2 + - 3/40 x DATA(0, 2) - 3/40 x DATA(2, -1) + 1/40 x DATA(2, 1) + + 2 + - 1/40 x DATA(2, 0) + 3/20 x y DATA(2, -1) - 1/20 x y DATA(0, 2) + + - 3/20 x y DATA(2, 1) - 3/20 x y DATA(2, 0) + 3/20 x y DATA(2, 2) + + 13 17 + + -- y DATA(2, 1) - -- y DATA(2, 0) - 1/40 y DATA(1, 2) + 40 40 + + + 7/40 y DATA(2, -1) + 1/40 y DATA(0, 2) + 3/40 y DATA(-1, 2) + + 11 23 + + -- y DATA(0, -1) + 9/40 y DATA(1, -1) - -- y DATA(-1, 0) + 40 40 + + 2 + - 3/40 y DATA(2, 2) - 1/8 y DATA(2, -1) + 7/40 y DATA(-1, 1) + + 13 2 + + -- y DATA(-1, -1) + 1/8 y DATA(1, 2) + 9/40 y DATA(0, 1) + 40 + + 19 21 11 2 + - -- y DATA(1, 0) - -- y DATA(0, 0) + -- y DATA(1, 1) + 1/8 y DATA(0, 2) + 40 40 40 + + 2 2 2 + - 3/8 y DATA(2, 1) + 3/8 y DATA(2, 0) + 1/8 y DATA(2, 2) + + 2 2 2 + + 1/8 y DATA(-1, 2) - 3/8 y DATA(-1, 1) - 1/8 y DATA(0, -1) + + 2 2 2 + - 1/8 y DATA(1, -1) + 3/8 y DATA(-1, 0) - 1/8 y DATA(-1, -1) + + 2 2 2 + - 3/8 y DATA(0, 1) + 3/8 y DATA(1, 0) + 3/8 y DATA(0, 0) + + - 1/50 x DATA(-1, 0) - 1/100 x DATA(1, 0) - 1/100 x DATA(0, 1) + + 2 + - 3/8 y DATA(1, 1) - 1/25 x DATA(0, 0) + 2/25 x DATA(1, 2) + + - 1/50 x DATA(0, -1) + 6/25 x DATA(-1, -1) + 3/50 x DATA(1, 1) + + 13 13 + - --- x DATA(1, -1) - --- x DATA(-1, 1) - 3/50 x DATA(2, 2) + 100 100 + + + 7/100 x DATA(2, 0) + 2/25 x DATA(2, 1) - 9/100 x DATA(-1, 2) + + 2 2 + + 1/40 x DATA(1, 0) - 1/40 x DATA(0, 1) + 7/100 x DATA(0, 2) + + 2 2 + - 9/100 x DATA(2, -1) + 3/40 x DATA(1, -1) + 3/40 x DATA(0, -1) + + 2 2 2 + + 1/40 x DATA(0, 0) - 1/40 x DATA(1, 1) - 1/40 x DATA(-1, 0) + + 2 2 2 + - 3/40 x DATA(1, 2) - 3/40 x DATA(-1, -1) + 1/40 x DATA(-1, 1) + + 2 31 157 + + 3/40 x DATA(2, 2) - --- DATA(-1, 0) - ---- DATA(1, -1) + 400 1200 + + 223 109 111 37 + - ---- DATA(0, -1) - ---- DATA(-1, -1) + --- DATA(-1, 1) + ---- DATA(1, 2) + 1200 1200 400 1200 + + 149 49 69 131 + - ---- DATA(2, 2) - --- DATA(2, 0) + --- DATA(2, 1) - ---- DATA(-1, 2) + 1200 400 400 1200 + + 43 89 + + ---- DATA(0, 2) + ---- DATA(2, -1) + 1200 1200 + +> coeff_as_lc_of_data(%, posn_list_2d_size4); + 2 109 2 13 +[COEFF(-1, -1) = - 1/8 y - ---- - 3/40 x + -- y + 6/25 x - 3/20 x y, + 1200 40 + + 2 2 223 11 + COEFF(0, -1) = - 1/8 y + 3/40 x - 1/50 x - ---- - 1/20 x y + -- y, + 1200 40 + + 13 2 2 157 + COEFF(1, -1) = - --- x - 1/8 y + 9/40 y + 3/40 x - ---- + 1/20 x y, + 100 1200 + + 2 2 89 + COEFF(2, -1) = - 1/8 y + 3/20 x y + 7/40 y - 3/40 x - 9/100 x + ----, + 1200 + + 23 2 31 2 + COEFF(-1, 0) = - -- y - 1/50 x + 3/8 y - --- - 1/40 x + 3/20 x y, + 40 400 + + 21 57 2 2 + COEFF(0, 0) = - -- y - --- + 1/20 x y - 1/25 x + 3/8 y + 1/40 x , + 40 400 + + 2 63 2 19 + COEFF(1, 0) = - 1/100 x + 3/8 y - --- - 1/20 x y + 1/40 x - -- y, + 400 40 + + 2 49 17 2 + COEFF(2, 0) = 3/8 y + 7/100 x - --- - -- y - 1/40 x - 3/20 x y, + 400 40 + + 2 2 111 13 + COEFF(-1, 1) = - 3/8 y + 1/40 x + 3/20 x y + --- + 7/40 y - --- x, + 400 100 + + 117 2 2 + COEFF(0, 1) = --- - 1/100 x - 1/40 x + 9/40 y + 1/20 x y - 3/8 y , + 400 + + 11 103 2 2 + COEFF(1, 1) = -- y + --- - 3/8 y + 3/50 x - 1/40 x - 1/20 x y, + 40 400 + + 2 2 69 13 + COEFF(2, 1) = - 3/8 y + 1/40 x + 2/25 x - 3/20 x y + --- + -- y, + 400 40 + + 2 2 131 + COEFF(-1, 2) = 1/8 y + 3/40 x + 3/40 y - 9/100 x - ---- - 3/20 x y, + 1200 + + 2 43 2 + COEFF(0, 2) = 7/100 x - 1/20 x y - 3/40 x + 1/40 y + ---- + 1/8 y , + 1200 + + 2 2 37 + COEFF(1, 2) = 1/8 y + 2/25 x - 3/40 x - 1/40 y + ---- + 1/20 x y, + 1200 + + 2 2 149 + COEFF(2, 2) = 1/8 y + 3/20 x y - 3/50 x + 3/40 x - ---- - 3/40 y] + 1200 + +> print_coeff__lc_of_data(%, "coeff_dy_", "fp", +> "2d.coeffs/2d.cube.order3.smooth0/coeff-dy.compute.c"); +bytes used=28017044, alloc=2096768, time=2.36 +bytes used=29017220, alloc=2096768, time=2.45 +bytes used=30017984, alloc=2096768, time=2.55 +bytes used=31018144, alloc=2096768, time=2.67 +> +# d^2/dx^2 +> simplify( diff(interp_2d_cube_order3_smooth0,x,x) ); + 21 19 11 +- -- DATA(0, 0) + 9/40 DATA(1, 0) - -- DATA(0, 1) + -- DATA(1, 1) + 40 40 40 + + + 1/20 y DATA(2, 1) - 1/20 y DATA(2, 0) - 3/20 y DATA(1, 2) + + - 3/20 y DATA(2, -1) - 3/20 y DATA(0, 2) + 3/20 y DATA(-1, 2) + + + 3/20 y DATA(0, -1) + 3/20 y DATA(1, -1) - 1/20 y DATA(-1, 0) + + + 3/20 y DATA(2, 2) + 1/20 y DATA(-1, 1) - 3/20 y DATA(-1, -1) + + - 1/20 y DATA(0, 1) + 1/20 y DATA(1, 0) + 1/20 y DATA(0, 0) + + - 1/20 y DATA(1, 1) - 1/4 x DATA(-1, 0) - 3/4 x DATA(1, 0) + + + 3/4 x DATA(0, 1) + 3/4 x DATA(0, 0) - 3/4 x DATA(1, 2) + + + 3/4 x DATA(0, -1) - 1/4 x DATA(-1, -1) - 3/4 x DATA(1, 1) + + - 3/4 x DATA(1, -1) - 1/4 x DATA(-1, 1) + 1/4 x DATA(2, 2) + + + 1/4 x DATA(2, 0) + 1/4 x DATA(2, 1) - 1/4 x DATA(-1, 2) + + 11 + + 3/4 x DATA(0, 2) + 1/4 x DATA(2, -1) + -- DATA(-1, 0) + 7/40 DATA(1, -1) + 40 + + 23 13 13 + - -- DATA(0, -1) + -- DATA(-1, -1) + 9/40 DATA(-1, 1) + -- DATA(1, 2) + 40 40 40 + + - 3/40 DATA(2, 2) + 1/40 DATA(2, 0) - 1/40 DATA(2, 1) + 7/40 DATA(-1, 2) + + 17 + - -- DATA(0, 2) + 3/40 DATA(2, -1) + 40 + +> coeff_as_lc_of_data(%, posn_list_2d_size4); +bytes used=32018340, alloc=2096768, time=2.81 + 13 23 +[COEFF(-1, -1) = - 3/20 y - 1/4 x + --, COEFF(0, -1) = - -- + 3/20 y + 3/4 x, + 40 40 + + COEFF(1, -1) = 7/40 + 3/20 y - 3/4 x, + + COEFF(2, -1) = - 3/20 y + 1/4 x + 3/40, + + 11 21 + COEFF(-1, 0) = - 1/4 x + -- - 1/20 y, COEFF(0, 0) = 1/20 y - -- + 3/4 x, + 40 40 + + COEFF(1, 0) = 9/40 + 1/20 y - 3/4 x, COEFF(2, 0) = 1/40 - 1/20 y + 1/4 x, + + 19 + COEFF(-1, 1) = - 1/4 x + 1/20 y + 9/40, COEFF(0, 1) = - 1/20 y - -- + 3/4 x, + 40 + + 11 + COEFF(1, 1) = - 1/20 y + -- - 3/4 x, COEFF(2, 1) = 1/20 y + 1/4 x - 1/40, + 40 + + 17 + COEFF(-1, 2) = 3/20 y - 1/4 x + 7/40, COEFF(0, 2) = - -- - 3/20 y + 3/4 x, + 40 + + 13 + COEFF(1, 2) = -- - 3/20 y - 3/4 x, COEFF(2, 2) = - 3/40 + 1/4 x + 3/20 y] + 40 + +> print_coeff__lc_of_data(%, "coeff_dxx_", "fp", +> "2d.coeffs/2d.cube.order3.smooth0/coeff-dxx.compute.c"); +bytes used=33018592, alloc=2096768, time=2.91 +> +# d^2/dxdy +> simplify( diff(interp_2d_cube_order3_smooth0,x,y) ); +bytes used=34018772, alloc=2096768, time=3.02 +- 1/25 DATA(0, 0) - 1/100 DATA(1, 0) - 1/100 DATA(0, 1) + 3/50 DATA(1, 1) + + - 3/20 y DATA(2, 1) - 3/20 y DATA(2, 0) + 1/20 y DATA(1, 2) + + + 3/20 y DATA(2, -1) - 1/20 y DATA(0, 2) - 3/20 y DATA(-1, 2) + + - 1/20 y DATA(0, -1) + 1/20 y DATA(1, -1) + 3/20 y DATA(-1, 0) + + + 3/20 y DATA(2, 2) + 3/20 y DATA(-1, 1) - 3/20 y DATA(-1, -1) + + + 1/20 y DATA(0, 1) - 1/20 y DATA(1, 0) + 1/20 y DATA(0, 0) + + - 1/20 y DATA(1, 1) - 1/20 x DATA(-1, 0) + 1/20 x DATA(1, 0) + + - 1/20 x DATA(0, 1) + 1/20 x DATA(0, 0) - 3/20 x DATA(1, 2) + + + 3/20 x DATA(0, -1) - 3/20 x DATA(-1, -1) - 1/20 x DATA(1, 1) + + + 3/20 x DATA(1, -1) + 1/20 x DATA(-1, 1) + 3/20 x DATA(2, 2) + + - 1/20 x DATA(2, 0) + 1/20 x DATA(2, 1) + 3/20 x DATA(-1, 2) + + - 3/20 x DATA(0, 2) - 3/20 x DATA(2, -1) - 1/50 DATA(-1, 0) + + 13 13 + - --- DATA(1, -1) - 1/50 DATA(0, -1) + 6/25 DATA(-1, -1) - --- DATA(-1, 1) + 100 100 + + + 2/25 DATA(1, 2) - 3/50 DATA(2, 2) + 7/100 DATA(2, 0) + 2/25 DATA(2, 1) + + - 9/100 DATA(-1, 2) + 7/100 DATA(0, 2) - 9/100 DATA(2, -1) + +> coeff_as_lc_of_data(%, posn_list_2d_size4); +[COEFF(-1, -1) = - 3/20 y - 3/20 x + 6/25, + + 13 + COEFF(0, -1) = 3/20 x - 1/20 y - 1/50, COEFF(1, -1) = 3/20 x + 1/20 y - ---, + 100 + + COEFF(2, -1) = - 9/100 - 3/20 x + 3/20 y, + + COEFF(-1, 0) = - 1/50 - 1/20 x + 3/20 y, + + COEFF(0, 0) = - 1/25 + 1/20 y + 1/20 x, + + COEFF(1, 0) = - 1/100 + 1/20 x - 1/20 y, + + COEFF(2, 0) = - 1/20 x + 7/100 - 3/20 y, + + 13 + COEFF(-1, 1) = - --- + 1/20 x + 3/20 y, + 100 + + COEFF(0, 1) = - 1/20 x - 1/100 + 1/20 y, + + COEFF(1, 1) = 3/50 - 1/20 x - 1/20 y, + + COEFF(2, 1) = - 3/20 y + 1/20 x + 2/25, + + COEFF(-1, 2) = - 3/20 y + 3/20 x - 9/100, + + COEFF(0, 2) = 7/100 - 3/20 x - 1/20 y, + + COEFF(1, 2) = - 3/20 x + 1/20 y + 2/25, + + COEFF(2, 2) = 3/20 x + 3/20 y - 3/50] + +> print_coeff__lc_of_data(%, "coeff_dxy_", "fp", +> "2d.coeffs/2d.cube.order3.smooth0/coeff-dxy.compute.c"); +bytes used=35018968, alloc=2096768, time=3.11 +> +# d^2/dy^2 +> simplify( diff(interp_2d_cube_order3_smooth0,y,y) ); + 21 19 11 +- -- DATA(0, 0) - -- DATA(1, 0) + 9/40 DATA(0, 1) + -- DATA(1, 1) + 40 40 40 + + - 3/4 y DATA(2, 1) + 3/4 y DATA(2, 0) + 1/4 y DATA(1, 2) + + - 1/4 y DATA(2, -1) + 1/4 y DATA(0, 2) + 1/4 y DATA(-1, 2) + + - 1/4 y DATA(0, -1) - 1/4 y DATA(1, -1) + 3/4 y DATA(-1, 0) + + + 1/4 y DATA(2, 2) - 3/4 y DATA(-1, 1) - 1/4 y DATA(-1, -1) + + - 3/4 y DATA(0, 1) + 3/4 y DATA(1, 0) + 3/4 y DATA(0, 0) + + - 3/4 y DATA(1, 1) + 3/20 x DATA(-1, 0) - 1/20 x DATA(1, 0) + + + 1/20 x DATA(0, 1) + 1/20 x DATA(0, 0) + 1/20 x DATA(1, 2) + + - 1/20 x DATA(0, -1) - 3/20 x DATA(-1, -1) - 1/20 x DATA(1, 1) + + + 1/20 x DATA(1, -1) + 3/20 x DATA(-1, 1) + 3/20 x DATA(2, 2) + + - 3/20 x DATA(2, 0) - 3/20 x DATA(2, 1) - 3/20 x DATA(-1, 2) + + 23 + - 1/20 x DATA(0, 2) + 3/20 x DATA(2, -1) - -- DATA(-1, 0) + 40 + + 11 13 + + 9/40 DATA(1, -1) + -- DATA(0, -1) + -- DATA(-1, -1) + 7/40 DATA(-1, 1) + 40 40 + + 17 13 + - 1/40 DATA(1, 2) - 3/40 DATA(2, 2) - -- DATA(2, 0) + -- DATA(2, 1) + 40 40 + + + 3/40 DATA(-1, 2) + 1/40 DATA(0, 2) + 7/40 DATA(2, -1) + +> coeff_as_lc_of_data(%, posn_list_2d_size4); +bytes used=36019260, alloc=2096768, time=3.25 + 13 11 +[COEFF(-1, -1) = - 1/4 y + -- - 3/20 x, COEFF(0, -1) = - 1/20 x - 1/4 y + --, + 40 40 + + COEFF(1, -1) = - 1/4 y + 1/20 x + 9/40, + + 23 + COEFF(2, -1) = 3/20 x - 1/4 y + 7/40, COEFF(-1, 0) = 3/4 y - -- + 3/20 x, + 40 + + 21 19 + COEFF(0, 0) = - -- + 1/20 x + 3/4 y, COEFF(1, 0) = - -- + 3/4 y - 1/20 x, + 40 40 + + 17 + COEFF(2, 0) = - 3/20 x + 3/4 y - --, COEFF(-1, 1) = 3/20 x - 3/4 y + 7/40, + 40 + + 11 + COEFF(0, 1) = 1/20 x + 9/40 - 3/4 y, COEFF(1, 1) = - 1/20 x + -- - 3/4 y, + 40 + + 13 + COEFF(2, 1) = - 3/20 x + -- - 3/4 y, COEFF(-1, 2) = - 3/20 x + 1/4 y + 3/40, + 40 + + COEFF(0, 2) = 1/4 y + 1/40 - 1/20 x, COEFF(1, 2) = 1/20 x - 1/40 + 1/4 y, + + COEFF(2, 2) = 3/20 x + 1/4 y - 3/40] + +> print_coeff__lc_of_data(%, "coeff_dyy_", "fp", +> "2d.coeffs/2d.cube.order3.smooth0/coeff-dyy.compute.c"); +bytes used=37019440, alloc=2096768, time=3.33 +> +################################################################################ +> +# +# 2d, cube, order=4, smoothing=0 (size=5) +# +> +# interpolating polynomial +> interp_2d_cube_order4_smooth0 +> := Lagrange_polynomial_interpolant(fn_2d_order4, coeff_list_2d_order4, +> coord_list_2d, posn_list_2d_size5); +bytes used=38019764, alloc=2096768, time=3.45 +bytes used=39020548, alloc=2096768, time=3.52 +bytes used=40020708, alloc=2162292, time=3.58 +bytes used=41022132, alloc=2162292, time=3.65 +bytes used=42022568, alloc=2162292, time=3.72 +bytes used=43024680, alloc=2162292, time=3.78 + 541 246 +interp_2d_cube_order4_smooth0 := ---- DATA(0, 0) + ---- DATA(1, 0) + 1225 1225 + + 246 24 246 24 + + ---- DATA(0, 1) - ---- DATA(1, 1) + ---- DATA(-1, 0) - ---- DATA(1, -1) + 1225 1225 1225 1225 + + 246 24 24 + + ---- DATA(0, -1) - ---- DATA(-1, -1) - ---- DATA(-1, 1) + ( + 1225 1225 1225 + + - 1/30 DATA(1, 0) - 1/30 DATA(1, 1) + 1/30 DATA(-1, 0) - 1/30 DATA(1, -1) + + + 1/30 DATA(-1, -1) + 1/30 DATA(-1, 1) - 1/30 DATA(1, 2) + 1/60 DATA(2, 2) + + + 1/60 DATA(2, 0) + 1/60 DATA(2, 1) + 1/30 DATA(-1, 2) + 1/60 DATA(2, -1) + + - 1/60 DATA(-2, 2) - 1/60 DATA(-2, 1) - 1/60 DATA(-2, -1) + + - 1/60 DATA(-2, 0) + 1/60 DATA(2, -2) + 1/30 DATA(-1, -2) + + 3 + - 1/30 DATA(1, -2) - 1/60 DATA(-2, -2)) x + (1/20 DATA(0, 0) + + - 1/30 DATA(1, 0) + 1/20 DATA(0, 1) - 1/30 DATA(1, 1) - 1/30 DATA(-1, 0) + + - 1/30 DATA(1, -1) + 1/20 DATA(0, -1) - 1/30 DATA(-1, -1) + + - 1/30 DATA(-1, 1) - 1/30 DATA(1, 2) + 1/120 DATA(2, 2) + 1/120 DATA(2, 0) + + + 1/120 DATA(2, 1) - 1/30 DATA(-1, 2) + 1/20 DATA(0, 2) + + + 1/120 DATA(2, -1) + 1/120 DATA(-2, 2) + 1/120 DATA(-2, 1) + + + 1/120 DATA(-2, -1) + 1/120 DATA(-2, 0) + 1/120 DATA(2, -2) + + - 1/30 DATA(-1, -2) + 1/20 DATA(0, -2) - 1/30 DATA(1, -2) + + 4 + + 1/120 DATA(-2, -2)) x + (1/20 DATA(0, 0) + 1/20 DATA(1, 0) + + - 1/30 DATA(0, 1) - 1/30 DATA(1, 1) + 1/20 DATA(-1, 0) - 1/30 DATA(1, -1) + + - 1/30 DATA(0, -1) - 1/30 DATA(-1, -1) - 1/30 DATA(-1, 1) + + + 1/120 DATA(1, 2) + 1/120 DATA(2, 2) + 1/20 DATA(2, 0) - 1/30 DATA(2, 1) + + + 1/120 DATA(-1, 2) + 1/120 DATA(0, 2) - 1/30 DATA(2, -1) + + + 1/120 DATA(-2, 2) - 1/30 DATA(-2, 1) - 1/30 DATA(-2, -1) + + + 1/20 DATA(-2, 0) + 1/120 DATA(2, -2) + 1/120 DATA(-1, -2) + + 4 / + + 1/120 DATA(0, -2) + 1/120 DATA(1, -2) + 1/120 DATA(-2, -2)) y + | + \ + + 17 31 31 17 + --- DATA(0, 1) + --- DATA(1, 1) - --- DATA(1, -1) - --- DATA(0, -1) + 105 210 210 105 + + 31 31 31 + - --- DATA(-1, -1) + --- DATA(-1, 1) + 1/84 DATA(1, 2) - --- DATA(2, 2) + 210 210 420 + + 11 17 11 + + --- DATA(2, 1) + 1/84 DATA(-1, 2) + --- DATA(0, 2) - --- DATA(2, -1) + 105 420 105 + + 31 11 11 31 + - --- DATA(-2, 2) + --- DATA(-2, 1) - --- DATA(-2, -1) + --- DATA(2, -2) + 420 105 105 420 + + 17 + - 1/84 DATA(-1, -2) - --- DATA(0, -2) - 1/84 DATA(1, -2) + 420 + + 31 \ /17 31 17 + + --- DATA(-2, -2)| y + |--- DATA(1, 0) + --- DATA(1, 1) - --- DATA(-1, 0) + 420 / \105 210 105 + + 31 31 31 11 + + --- DATA(1, -1) - --- DATA(-1, -1) - --- DATA(-1, 1) + --- DATA(1, 2) + 210 210 210 105 + + 31 17 11 + - --- DATA(2, 2) + --- DATA(2, 0) + 1/84 DATA(2, 1) - --- DATA(-1, 2) + 420 420 105 + + 31 + + 1/84 DATA(2, -1) + --- DATA(-2, 2) - 1/84 DATA(-2, 1) + 420 + + 17 31 11 + - 1/84 DATA(-2, -1) - --- DATA(-2, 0) - --- DATA(2, -2) - --- DATA(-1, -2) + 420 420 105 + + 11 31 \ / 57 + + --- DATA(1, -2) + --- DATA(-2, -2)| x + |- --- DATA(0, 0) + 105 420 / \ 196 + + 83 53 181 83 + + --- DATA(1, 0) - --- DATA(0, 1) + ---- DATA(1, 1) + --- DATA(-1, 0) + 735 196 1470 735 + + 181 53 181 + + ---- DATA(1, -1) - --- DATA(0, -1) + ---- DATA(-1, -1) + 1470 196 1470 + + 181 113 289 191 + + ---- DATA(-1, 1) + --- DATA(1, 2) - ---- DATA(2, 2) + ---- DATA(2, 0) + 1470 735 5880 5880 + + 71 113 41 71 + + ---- DATA(2, 1) + --- DATA(-1, 2) - --- DATA(0, 2) + ---- DATA(2, -1) + 5880 735 196 5880 + + 289 71 71 + - ---- DATA(-2, 2) + ---- DATA(-2, 1) + ---- DATA(-2, -1) + 5880 5880 5880 + + 191 289 113 41 + + ---- DATA(-2, 0) - ---- DATA(2, -2) + --- DATA(-1, -2) - --- DATA(0, -2) + 5880 5880 735 196 + + 113 289 \ 2 + + --- DATA(1, -2) - ---- DATA(-2, -2)| x + (- 1/70 DATA(1, 0) + 735 5880 / + + - 1/140 DATA(1, 1) + 1/70 DATA(-1, 0) - 1/140 DATA(1, -1) + + + 1/140 DATA(-1, -1) + 1/140 DATA(-1, 1) + 1/70 DATA(1, 2) + + + 1/35 DATA(2, 2) - 1/35 DATA(2, 0) - 1/70 DATA(2, 1) - 1/70 DATA(-1, 2) + + - 1/70 DATA(2, -1) - 1/35 DATA(-2, 2) + 1/70 DATA(-2, 1) + + + 1/70 DATA(-2, -1) + 1/35 DATA(-2, 0) + 1/35 DATA(2, -2) + + 2 + - 1/70 DATA(-1, -2) + 1/70 DATA(1, -2) - 1/35 DATA(-2, -2)) x y + ( + + - 1/60 DATA(1, 1) + 1/60 DATA(1, -1) - 1/60 DATA(-1, -1) + 1/60 DATA(-1, 1) + + + 1/120 DATA(1, 2) + 1/60 DATA(2, 2) - 1/30 DATA(2, 1) - 1/120 DATA(-1, 2) + + + 1/30 DATA(2, -1) - 1/60 DATA(-2, 2) + 1/30 DATA(-2, 1) + + - 1/30 DATA(-2, -1) - 1/60 DATA(2, -2) + 1/120 DATA(-1, -2) + + 3 / 57 + - 1/120 DATA(1, -2) + 1/60 DATA(-2, -2)) x y + |- --- DATA(0, 0) + \ 196 + + 53 83 181 53 + - --- DATA(1, 0) + --- DATA(0, 1) + ---- DATA(1, 1) - --- DATA(-1, 0) + 196 735 1470 196 + + 181 83 181 + + ---- DATA(1, -1) + --- DATA(0, -1) + ---- DATA(-1, -1) + 1470 735 1470 + + 181 71 289 41 + + ---- DATA(-1, 1) + ---- DATA(1, 2) - ---- DATA(2, 2) - --- DATA(2, 0) + 1470 5880 5880 196 + + 113 71 191 113 + + --- DATA(2, 1) + ---- DATA(-1, 2) + ---- DATA(0, 2) + --- DATA(2, -1) + 735 5880 5880 735 + + 289 113 113 41 + - ---- DATA(-2, 2) + --- DATA(-2, 1) + --- DATA(-2, -1) - --- DATA(-2, 0) + 5880 735 735 196 + + 289 71 191 + - ---- DATA(2, -2) + ---- DATA(-1, -2) + ---- DATA(0, -2) + 5880 5880 5880 + + 71 289 \ 2 + + ---- DATA(1, -2) - ---- DATA(-2, -2)| y + (- 1/30 DATA(0, 1) + 5880 5880 / + + - 1/30 DATA(1, 1) + 1/30 DATA(1, -1) + 1/30 DATA(0, -1) + + + 1/30 DATA(-1, -1) - 1/30 DATA(-1, 1) + 1/60 DATA(1, 2) + 1/60 DATA(2, 2) + + - 1/30 DATA(2, 1) + 1/60 DATA(-1, 2) + 1/60 DATA(0, 2) + 1/30 DATA(2, -1) + + + 1/60 DATA(-2, 2) - 1/30 DATA(-2, 1) + 1/30 DATA(-2, -1) + + - 1/60 DATA(2, -2) - 1/60 DATA(-1, -2) - 1/60 DATA(0, -2) + + 3 + - 1/60 DATA(1, -2) - 1/60 DATA(-2, -2)) y + (1/49 DATA(0, 0) + + + 1/98 DATA(1, 0) + 1/98 DATA(0, 1) + 1/196 DATA(1, 1) + 1/98 DATA(-1, 0) + + + 1/196 DATA(1, -1) + 1/98 DATA(0, -1) + 1/196 DATA(-1, -1) + + + 1/196 DATA(-1, 1) - 1/98 DATA(1, 2) + 1/49 DATA(2, 2) - 1/49 DATA(2, 0) + + - 1/98 DATA(2, 1) - 1/98 DATA(-1, 2) - 1/49 DATA(0, 2) - 1/98 DATA(2, -1) + + + 1/49 DATA(-2, 2) - 1/98 DATA(-2, 1) - 1/98 DATA(-2, -1) + + - 1/49 DATA(-2, 0) + 1/49 DATA(2, -2) - 1/98 DATA(-1, -2) + + 2 2 / + - 1/49 DATA(0, -2) - 1/98 DATA(1, -2) + 1/49 DATA(-2, -2)) x y + | + \ + + 37 37 37 37 + --- DATA(1, 1) - --- DATA(1, -1) + --- DATA(-1, -1) - --- DATA(-1, 1) + 300 300 300 300 + + 21 11 21 21 + + --- DATA(1, 2) - --- DATA(2, 2) + --- DATA(2, 1) - --- DATA(-1, 2) + 200 150 200 200 + + 21 11 21 21 + - --- DATA(2, -1) + --- DATA(-2, 2) - --- DATA(-2, 1) + --- DATA(-2, -1) + 200 150 200 200 + + 11 21 21 11 \ + + --- DATA(2, -2) + --- DATA(-1, -2) - --- DATA(1, -2) - --- DATA(-2, -2)| + 150 200 200 150 / + + x y + (- 1/70 DATA(0, 1) - 1/140 DATA(1, 1) + 1/140 DATA(1, -1) + + + 1/70 DATA(0, -1) + 1/140 DATA(-1, -1) - 1/140 DATA(-1, 1) + + - 1/70 DATA(1, 2) + 1/35 DATA(2, 2) + 1/70 DATA(2, 1) - 1/70 DATA(-1, 2) + + - 1/35 DATA(0, 2) - 1/70 DATA(2, -1) + 1/35 DATA(-2, 2) + 1/70 DATA(-2, 1) + + - 1/70 DATA(-2, -1) - 1/35 DATA(2, -2) + 1/70 DATA(-1, -2) + + 2 + + 1/35 DATA(0, -2) + 1/70 DATA(1, -2) - 1/35 DATA(-2, -2)) x y + ( + + - 1/60 DATA(1, 1) + 1/60 DATA(1, -1) - 1/60 DATA(-1, -1) + 1/60 DATA(-1, 1) + + - 1/30 DATA(1, 2) + 1/60 DATA(2, 2) + 1/120 DATA(2, 1) + 1/30 DATA(-1, 2) + + - 1/120 DATA(2, -1) - 1/60 DATA(-2, 2) - 1/120 DATA(-2, 1) + + + 1/120 DATA(-2, -1) - 1/60 DATA(2, -2) - 1/30 DATA(-1, -2) + + 3 99 + + 1/30 DATA(1, -2) + 1/60 DATA(-2, -2)) x y - ---- DATA(1, 2) + 1225 + + 51 96 99 99 + + ---- DATA(2, 2) + ---- DATA(2, 0) - ---- DATA(2, 1) - ---- DATA(-1, 2) + 1225 1225 1225 1225 + + 96 99 51 99 + + ---- DATA(0, 2) - ---- DATA(2, -1) + ---- DATA(-2, 2) - ---- DATA(-2, 1) + 1225 1225 1225 1225 + + 99 96 51 + - ---- DATA(-2, -1) + ---- DATA(-2, 0) + ---- DATA(2, -2) + 1225 1225 1225 + + 99 96 99 + - ---- DATA(-1, -2) + ---- DATA(0, -2) - ---- DATA(1, -2) + 1225 1225 1225 + + 51 + + ---- DATA(-2, -2) + 1225 + +> +# I +> coeff_as_lc_of_data(%, posn_list_2d_size5); +bytes used=44033940, alloc=2162292, time=3.85 +bytes used=45039820, alloc=2162292, time=3.91 + 51 3 4 3 2 2 +[COEFF(-2, -2) = ---- + 1/60 x y + 1/120 x - 1/60 y - 1/35 x y - 1/35 x y + 1225 + + 4 289 2 2 2 31 3 31 11 + + 1/120 y - ---- x + 1/49 x y + --- y - 1/60 x + --- x - --- x y + 5880 420 420 150 + + 3 289 2 21 2 4 + + 1/60 x y - ---- y , COEFF(-1, -2) = --- x y - 1/70 x y + 1/120 y + 5880 200 + + 3 3 99 113 2 71 2 11 3 + + 1/120 x y + 1/30 x - ---- + --- x + ---- y - --- x - 1/60 y + 1225 735 5880 105 + + 4 3 2 2 2 + - 1/30 x - 1/30 x y - 1/84 y - 1/98 x y + 1/70 x y, COEFF(0, -2) = + + 2 2 2 4 41 2 17 96 4 + 1/35 x y - 1/49 x y + 1/20 x - --- x - --- y + ---- + 1/120 y + 196 420 1225 + + 3 191 2 2 2 99 + - 1/60 y + ---- y , COEFF(1, -2) = - 1/98 x y - ---- - 1/84 y + 5880 1225 + + 2 113 2 3 4 3 2 + + 1/70 x y + --- x - 1/120 x y + 1/120 y + 1/30 x y + 1/70 x y + 735 + + 21 4 3 71 2 11 3 + - --- x y - 1/30 x - 1/30 x + ---- y + --- x - 1/60 y , COEFF(2, -2) = + 200 5880 105 + + 2 2 4 3 4 31 31 3 + 1/49 x y + 1/120 y - 1/60 x y + 1/120 x + --- y - --- x + 1/60 x + 420 420 + + 2 51 3 3 11 289 2 2 + - 1/35 x y + ---- - 1/60 x y - 1/60 y + --- x y - ---- y + 1/35 x y + 1225 150 5880 + + 289 2 2 2 2 4 + - ---- x , COEFF(-2, -1) = - 1/98 x y + 1/70 x y + 1/120 x + 5880 + + 3 11 21 71 2 2 3 + + 1/120 x y - --- y + --- x y + ---- x - 1/70 x y + 1/30 y - 1/84 x + 105 200 5880 + + 3 3 113 2 4 99 24 + - 1/60 x - 1/30 x y + --- y - 1/30 y - ----, COEFF(-1, -1) = - ---- + 735 1225 1225 + + 31 4 181 2 2 31 3 3 + - --- x - 1/30 y + ---- x + 1/140 x y - --- y - 1/60 x y + 1/30 y + 210 1470 210 + + 4 3 3 181 2 37 2 2 + - 1/30 x + 1/30 x - 1/60 x y + ---- y + --- x y + 1/196 x y + 1470 300 + + 2 246 83 2 3 53 2 17 + + 1/140 x y , COEFF(0, -1) = ---- + --- y + 1/30 y - --- x - --- y + 1225 735 196 105 + + 2 4 2 2 4 181 2 + + 1/70 x y + 1/20 x + 1/98 x y - 1/30 y , COEFF(1, -1) = ---- y + 1470 + + 4 4 3 3 24 181 2 2 + - 1/30 x - 1/30 y - 1/30 x + 1/60 x y - ---- + ---- x + 1/140 x y + 1225 1470 + + 37 31 2 3 3 2 2 31 + - --- x y - --- y - 1/140 x y + 1/30 y + 1/60 x y + 1/196 x y + --- x + 300 210 210 + + 3 21 4 2 2 2 + , COEFF(2, -1) = 1/30 x y - --- x y - 1/30 y - 1/70 x y - 1/98 x y + 200 + + 99 11 71 2 113 2 3 2 + - ---- + 1/84 x - --- y + ---- x + --- y + 1/30 y - 1/70 x y + 1225 105 5880 735 + + 4 3 3 4 3 + + 1/120 x - 1/120 x y + 1/60 x , COEFF(-2, 0) = 1/120 x - 1/60 x + + 191 2 17 2 2 4 2 41 2 96 + + ---- x - --- x - 1/49 x y + 1/20 y + 1/35 x y - --- y + ----, + 5880 420 196 1225 + + 4 17 83 2 246 53 2 4 + COEFF(-1, 0) = 1/20 y - --- x + --- x + ---- - --- y - 1/30 x + 105 735 1225 196 + + 2 3 2 2 + + 1/70 x y + 1/30 x + 1/98 x y , + + 57 2 4 4 57 2 2 2 541 + COEFF(0, 0) = - --- x + 1/20 y + 1/20 x - --- y + 1/49 x y + ----, + 196 196 1225 + + 246 17 2 53 2 4 4 + COEFF(1, 0) = ---- + --- x - 1/70 x y - --- y - 1/30 x + 1/20 y + 1225 105 196 + + 3 83 2 2 2 2 4 + - 1/30 x + --- x + 1/98 x y , COEFF(2, 0) = - 1/35 x y + 1/20 y + 735 + + 4 96 3 17 2 2 41 2 191 2 + + 1/120 x + ---- + 1/60 x + --- x - 1/49 x y - --- y + ---- x , + 1225 420 196 5880 + + 2 2 3 99 3 2 11 + COEFF(-2, 1) = - 1/98 x y - 1/60 x - ---- - 1/30 y + 1/70 x y + --- y + 1225 105 + + 2 4 113 2 3 4 3 + + 1/70 x y - 1/30 y + --- y + 1/30 x y + 1/120 x - 1/120 x y + 735 + + 21 71 2 2 2 3 + - 1/84 x - --- x y + ---- x , COEFF(-1, 1) = 1/196 x y + 1/30 x + 200 5880 + + 3 2 4 3 37 31 31 + + 1/60 x y - 1/140 x y - 1/30 x + 1/60 x y - --- x y - --- x + --- y + 300 210 210 + + 181 2 2 3 24 181 2 4 + + ---- y + 1/140 x y - 1/30 y - ---- + ---- x - 1/30 y , COEFF(0, 1) + 1470 1225 1470 + + 3 2 2 53 2 246 4 17 2 + = - 1/30 y + 1/98 x y - --- x + ---- - 1/30 y + --- y - 1/70 x y + 196 1225 105 + + 4 83 2 3 24 3 + + 1/20 x + --- y , COEFF(1, 1) = - 1/30 y - ---- - 1/60 x y + 735 1225 + + 2 2 37 3 181 2 31 4 31 + + 1/196 x y + --- x y - 1/60 x y + ---- x + --- y - 1/30 y + --- x + 300 1470 210 210 + + 4 2 3 181 2 2 + - 1/30 x - 1/140 x y - 1/30 x + ---- y - 1/140 x y, COEFF(2, 1) = + 1470 + + 3 113 2 3 2 2 21 + - 1/30 x y + 1/84 x + --- y + 1/120 x y - 1/98 x y + --- x y + 735 200 + + 2 99 4 71 2 3 11 4 + - 1/70 x y - ---- - 1/30 y + ---- x + 1/60 x + --- y + 1/120 x + 1225 5880 105 + + 3 2 2 3 11 + - 1/30 y + 1/70 x y, COEFF(-2, 2) = 1/35 x y - 1/60 x y + --- x y + 150 + + 2 31 289 2 4 2 2 289 2 31 + - 1/35 x y - --- y - ---- y + 1/120 y + 1/49 x y - ---- x + --- x + 420 5880 5880 420 + + 3 4 3 51 3 99 + - 1/60 x y + 1/120 x + 1/60 y + ---- - 1/60 x , COEFF(-1, 2) = - ---- + 1225 1225 + + 3 21 71 2 2 2 2 + - 1/120 x y - --- x y + 1/84 y + ---- y - 1/98 x y - 1/70 x y + 200 5880 + + 113 2 4 3 4 3 11 3 + + --- x - 1/30 x + 1/30 x y + 1/120 y + 1/30 x - --- x + 1/60 y + 735 105 + + 2 4 2 17 41 2 191 2 + - 1/70 x y, COEFF(0, 2) = 1/20 x - 1/35 x y + --- y - --- x + ---- y + 420 196 5880 + + 96 4 2 2 3 99 + + ---- + 1/120 y - 1/49 x y + 1/60 y , COEFF(1, 2) = - ---- + 1225 1225 + + 2 2 2 3 71 2 3 + - 1/98 x y - 1/70 x y - 1/30 x y + 1/84 y + ---- y + 1/120 x y + 5880 + + 3 4 4 11 113 2 2 21 + - 1/30 x - 1/30 x + 1/120 y + --- x + --- x + 1/70 x y + --- x y + 105 735 200 + + 3 2 4 2 2 289 2 + + 1/60 y , COEFF(2, 2) = 1/35 x y + 1/120 x + 1/49 x y - ---- y + 5880 + + 3 3 31 4 289 2 3 2 + + 1/60 x y + 1/60 x y - --- y + 1/120 y - ---- x + 1/60 x + 1/35 x y + 420 5880 + + 11 31 51 3 + - --- x y - --- x + ---- + 1/60 y ] + 150 420 1225 + +> print_coeff__lc_of_data(%, "coeff_I_", "fp", +> "2d.coeffs/2d.cube.order4.smooth0/coeff-I.compute.c"); +bytes used=46039980, alloc=2162292, time=4.00 +bytes used=47040144, alloc=2162292, time=4.08 +bytes used=48048756, alloc=2162292, time=4.16 +bytes used=49050036, alloc=2227816, time=4.23 +bytes used=50050424, alloc=2227816, time=4.33 +bytes used=51050604, alloc=2227816, time=4.48 +bytes used=52050764, alloc=2227816, time=4.55 +bytes used=53052936, alloc=2227816, time=4.64 +bytes used=54053412, alloc=2227816, time=4.74 +bytes used=55053604, alloc=2227816, time=4.89 +bytes used=56053828, alloc=2227816, time=4.95 +bytes used=57054040, alloc=2293340, time=5.07 +bytes used=58054312, alloc=2293340, time=5.26 +bytes used=59054740, alloc=2293340, time=5.33 +bytes used=60054904, alloc=2293340, time=5.43 +bytes used=61055148, alloc=2293340, time=5.58 +bytes used=62055404, alloc=2293340, time=5.65 +bytes used=63058408, alloc=2293340, time=5.83 +bytes used=64058672, alloc=2293340, time=5.89 +bytes used=65058956, alloc=2293340, time=6.08 +bytes used=66059172, alloc=2293340, time=6.16 +> +# d/dx +> simplify( diff(interp_2d_cube_order4_smooth0,x) ); +bytes used=67059496, alloc=2293340, time=6.26 +bytes used=68059676, alloc=2293340, time=6.33 +bytes used=69060256, alloc=2358864, time=6.40 +17 31 17 31 +--- DATA(1, 0) + --- DATA(1, 1) - --- DATA(-1, 0) + --- DATA(1, -1) +105 210 105 210 + + 31 31 11 31 + - --- DATA(-1, -1) - --- DATA(-1, 1) + --- DATA(1, 2) - --- DATA(2, 2) + 210 210 105 420 + + 17 11 + + --- DATA(2, 0) + 1/84 DATA(2, 1) - --- DATA(-1, 2) + 1/84 DATA(2, -1) + 420 105 + + 31 2 + + --- DATA(-2, 2) - 1/84 DATA(-2, 1) - 1/70 y DATA(2, 1) + 420 + + 2 2 2 + - 1/70 y DATA(-1, 2) - 1/70 y DATA(2, -1) - 1/35 y DATA(-2, 2) + + 2 2 2 + + 1/70 y DATA(-2, 1) + 1/70 y DATA(-2, -1) + 1/35 y DATA(-2, 0) + + 2 2 2 + + 1/35 y DATA(2, -2) - 1/70 y DATA(-1, -2) + 1/35 y DATA(2, 2) + + 2 2 2 + - 1/70 y DATA(1, 0) - 1/140 y DATA(1, 1) + 1/70 y DATA(-1, 0) + + 2 2 2 + - 1/140 y DATA(1, -1) + 1/140 y DATA(-1, -1) + 1/140 y DATA(-1, 1) + + - 1/35 x y DATA(-1, 2) - 2/35 x y DATA(0, 2) - 1/35 x y DATA(0, 1) + + 2 2 + - 1/20 x y DATA(-1, -1) - 1/20 x y DATA(1, 1) - 1/35 x y DATA(-2, -1) + + - 2/35 x y DATA(2, -2) + 2/35 x y DATA(2, 2) + 1/35 x y DATA(2, 1) + + 37 2 2 + + --- y DATA(1, 1) - 1/49 x y DATA(-2, -1) + 2/49 x y DATA(2, -2) + 300 + + 2 2 + - 1/49 x y DATA(-1, -2) - 1/49 x y DATA(1, -2) - 1/35 x y DATA(2, -1) + + 21 + + 2/35 x y DATA(-2, 2) + 1/35 x y DATA(-2, 1) - --- y DATA(-1, 2) + 200 + + 21 11 21 + - --- y DATA(2, -1) + --- y DATA(-2, 2) - --- y DATA(-2, 1) + 200 150 200 + + 21 11 21 + + --- y DATA(-2, -1) + --- y DATA(2, -2) + --- y DATA(-1, -2) + 200 150 200 + + 21 11 37 + - --- y DATA(1, -2) - --- y DATA(-2, -2) - --- y DATA(1, -1) + 200 150 300 + + 2 2 2 + + 1/70 y DATA(1, -2) - 1/35 y DATA(-2, -2) + 1/70 y DATA(1, 2) + + 2 2 + - 1/35 y DATA(2, 0) + 1/35 x y DATA(-1, -2) - 2/49 x y DATA(0, 2) + + 2 2 2 + - 2/49 x y DATA(0, -2) + 1/49 x y DATA(-1, 0) - 2/49 x y DATA(2, 0) + + 2 + - 2/49 x y DATA(-2, 0) + 2/35 x y DATA(0, -2) + 1/35 x y DATA(1, -2) + + 2 2 + - 2/35 x y DATA(-2, -2) - 1/49 x y DATA(2, 1) - 1/49 x y DATA(2, -1) + + 2 2 37 + + 2/49 x y DATA(-2, 2) - 1/49 x y DATA(-2, 1) + --- y DATA(-1, -1) + 300 + + 3 3 3 + + 1/30 y DATA(-2, 1) + 1/60 y DATA(-1, 1) - 1/60 y DATA(1, 1) + + 3 37 21 + + 1/60 y DATA(1, -1) - --- y DATA(-1, 1) + --- y DATA(1, 2) + 300 200 + + 11 21 3 + - --- y DATA(2, 2) + --- y DATA(2, 1) + 1/60 y DATA(-2, -2) + 150 200 + + 71 289 71 + + ---- x DATA(2, -1) - ---- x DATA(-2, 2) + ---- x DATA(-2, 1) + 2940 2940 2940 + + 71 289 226 + + ---- x DATA(-2, -1) - ---- x DATA(2, -2) + --- x DATA(-1, -2) + 2940 2940 735 + + 226 289 53 + + --- x DATA(1, -2) - ---- x DATA(-2, -2) - -- x DATA(0, 1) + 735 2940 98 + + 53 41 41 + - -- x DATA(0, -1) - -- x DATA(0, 2) - -- x DATA(0, -2) + 98 98 98 + + 166 3 3 + + --- x DATA(-1, 0) - 1/60 y DATA(2, -2) - 1/60 y DATA(-1, -1) + 735 + + 3 3 3 + + 1/120 y DATA(1, 2) + 1/60 y DATA(2, 2) - 1/30 y DATA(2, 1) + + 3 3 3 + - 1/120 y DATA(-1, 2) + 1/30 y DATA(2, -1) - 1/60 y DATA(-2, 2) + + 71 191 3 + + ---- x DATA(2, 1) + ---- x DATA(-2, 0) - 1/30 y DATA(-2, -1) + 2940 2940 + + 3 3 191 + + 1/120 y DATA(-1, -2) - 1/120 y DATA(1, -2) + ---- x DATA(2, 0) + 2940 + + 17 3 + - 1/84 DATA(-2, -1) - --- DATA(-2, 0) - 2/15 x DATA(1, -2) + 420 + + 3 3 3 + + 1/30 x DATA(-2, -2) + 1/5 x DATA(0, 1) + 1/5 x DATA(0, -1) + + 3 3 3 + + 1/5 x DATA(0, 2) + 1/5 x DATA(0, -2) - 2/15 x DATA(-1, 0) + + 3 3 57 + + 1/30 x DATA(2, 0) + 1/30 x DATA(-2, 0) - -- x DATA(0, 0) + 98 + + 226 181 181 + + --- x DATA(-1, 2) + --- x DATA(1, 1) + --- x DATA(1, -1) + 735 735 735 + + 181 181 226 + + --- x DATA(-1, -1) + --- x DATA(-1, 1) + --- x DATA(1, 2) + 735 735 735 + + 289 166 3 + - ---- x DATA(2, 2) + --- x DATA(1, 0) - 2/15 x DATA(-1, 2) + 2940 735 + + 3 3 3 + - 2/15 x DATA(1, 1) - 2/15 x DATA(1, -1) - 2/15 x DATA(-1, -1) + + 3 3 3 + - 2/15 x DATA(-1, 1) + 1/30 x DATA(2, 1) + 1/30 x DATA(2, -1) + + 3 3 3 + + 1/30 x DATA(-2, 2) + 1/30 x DATA(-2, 1) + 1/30 x DATA(-2, -1) + + 3 3 2 + + 1/30 x DATA(2, -2) - 2/15 x DATA(-1, -2) - 1/20 x DATA(-2, -2) + + 3 3 3 + - 2/15 x DATA(1, 2) + 1/30 x DATA(2, 2) - 2/15 x DATA(1, 0) + + 2 2 2 + - 1/10 x DATA(1, 0) - 1/10 x DATA(1, 1) + 1/10 x DATA(-1, 0) + + 2 2 2 + - 1/10 x DATA(1, -1) + 1/10 x DATA(-1, -1) + 1/10 x DATA(-1, 1) + + 3 2 2 + + 1/5 x DATA(0, 0) - 1/10 x y DATA(-1, -2) + 1/20 x y DATA(1, -1) + + 2 2 2 + + 1/40 x y DATA(-2, -1) - 1/49 x y DATA(1, 2) + 2/49 x y DATA(2, 2) + + 2 2 2 + + 1/49 x y DATA(1, 0) + 2/49 x y DATA(0, 0) - 1/49 x y DATA(-1, 2) + + 2 2 2 + + 1/98 x y DATA(1, 1) + 1/98 x y DATA(1, -1) + 1/98 x y DATA(-1, -1) + + 2 2 2 + + 1/98 x y DATA(-1, 1) + 2/49 x y DATA(-2, -2) + 1/49 x y DATA(0, 1) + + 2 + + 1/49 x y DATA(0, -1) - 1/35 x y DATA(1, 2) + 1/70 x y DATA(1, -1) + + 2 2 + + 1/35 x y DATA(0, -1) + 1/10 x y DATA(1, -2) + 1/10 x y DATA(-1, 2) + + 2 31 + - 1/40 x y DATA(-2, 1) - 1/70 x y DATA(1, 1) - --- DATA(2, -2) + 420 + + 11 11 2 + - --- DATA(-1, -2) + --- DATA(1, -2) - 1/40 x y DATA(2, -1) + 105 105 + + 2 2 2 + + 1/20 x y DATA(-2, -2) - 1/10 x y DATA(1, 2) + 1/40 x y DATA(2, 1) + + 2 2 2 + + 1/20 x y DATA(2, 2) - 1/20 x y DATA(2, -2) - 1/20 x y DATA(-2, 2) + + 2 + + 1/20 x y DATA(-1, 1) + 1/70 x y DATA(-1, -1) - 1/70 x y DATA(-1, 1) + + 2 2 2 + - 1/10 x DATA(1, 2) + 1/20 x DATA(2, 0) + 1/20 x DATA(2, 1) + + 2 2 2 + + 1/10 x DATA(-1, 2) + 1/20 x DATA(2, -1) + 1/20 x DATA(2, 2) + + 2 2 2 + - 1/20 x DATA(-2, 2) - 1/20 x DATA(-2, 1) - 1/20 x DATA(-2, -1) + + 2 2 2 + - 1/20 x DATA(-2, 0) + 1/20 x DATA(2, -2) + 1/10 x DATA(-1, -2) + + 2 31 + - 1/10 x DATA(1, -2) + --- DATA(-2, -2) + 420 + +> coeff_as_lc_of_data(%, posn_list_2d_size5); +bytes used=70063472, alloc=2424388, time=6.47 + 31 3 2 11 3 +[COEFF(-2, -2) = --- + 1/60 y - 2/35 x y + 2/49 x y - --- y + 1/30 x + 420 150 + + 2 2 2 289 2 + - 1/20 x + 1/20 x y - 1/35 y - ---- x, COEFF(-1, -2) = - 1/10 x y + 2940 + + 2 2 3 2 226 11 + - 1/70 y - 1/49 x y + 1/120 y + 1/10 x + --- x - --- + 1/35 x y + 735 105 + + 21 3 2 41 3 + + --- y - 2/15 x , COEFF(0, -2) = - 2/49 x y - -- x + 1/5 x + 2/35 x y, + 200 98 + + 3 2 3 2 + COEFF(1, -2) = - 2/15 x - 1/49 x y - 1/120 y - 1/10 x + 1/35 x y + + 21 226 11 2 2 11 + - --- y + --- x + --- + 1/10 x y + 1/70 y , COEFF(2, -2) = --- y + 200 735 105 150 + + 2 2 2 31 289 3 3 + + 2/49 x y + 1/20 x + 1/35 y - --- - ---- x - 1/60 y + 1/30 x + 420 2940 + + 2 3 21 3 + - 2/35 x y - 1/20 x y, COEFF(-2, -1) = 1/30 x + --- y - 1/30 y + 200 + + 2 2 2 71 2 + - 1/49 x y - 1/35 x y + 1/40 x y + 1/70 y + ---- x - 1/20 x - 1/84, + 2940 + + 3 31 2 3 2 + COEFF(-1, -1) = - 2/15 x - --- + 1/140 y + 1/70 x y - 1/60 y - 1/20 x y + 210 + + 37 2 2 181 + + --- y + 1/10 x + 1/98 x y + --- x, + 300 735 + + 3 2 53 3 + COEFF(0, -1) = 1/35 x y + 1/5 x + 1/49 x y - -- x, COEFF(1, -1) = 1/60 y + 98 + + 3 2 2 31 2 181 + - 2/15 x + 1/20 x y - 1/10 x + --- + 1/98 x y + --- x + 1/70 x y + 210 735 + + 2 37 2 2 + - 1/140 y - --- y, COEFF(2, -1) = - 1/70 y - 1/49 x y - 1/35 x y + 300 + + 21 2 3 3 2 71 + - --- y - 1/40 x y + 1/30 x + 1/30 y + 1/84 + 1/20 x + ---- x, + 200 2940 + + 2 2 191 3 2 17 + COEFF(-2, 0) = - 2/49 x y + 1/35 y + ---- x + 1/30 x - 1/20 x - ---, + 2940 420 + + 2 17 166 2 2 3 + COEFF(-1, 0) = 1/49 x y - --- + --- x + 1/10 x + 1/70 y - 2/15 x , + 105 735 + + 57 3 2 + COEFF(0, 0) = - -- x + 1/5 x + 2/49 x y , + 98 + + 2 2 17 2 3 166 + COEFF(1, 0) = 1/49 x y - 1/70 y + --- - 1/10 x - 2/15 x + --- x, + 105 735 + + 2 2 17 2 3 191 + COEFF(2, 0) = - 1/35 y - 2/49 x y + --- + 1/20 x + 1/30 x + ---- x, + 420 2940 + + 3 3 2 2 71 + COEFF(-2, 1) = 1/30 y + 1/30 x - 1/49 x y + 1/35 x y + 1/70 y + ---- x + 2940 + + 2 2 21 2 + - 1/20 x - 1/40 x y - --- y - 1/84, COEFF(-1, 1) = 1/20 x y - 1/70 x y + 200 + + 2 2 2 181 3 3 31 37 + + 1/10 x + 1/98 x y + 1/140 y + --- x - 2/15 x + 1/60 y - --- - --- y + 735 210 300 + + 53 2 3 181 + , COEFF(0, 1) = - -- x + 1/49 x y + 1/5 x - 1/35 x y, COEFF(1, 1) = --- x + 98 735 + + 2 3 3 31 37 2 + - 1/10 x - 1/60 y - 2/15 x + --- + --- y - 1/140 y - 1/70 x y + 210 300 + + 2 2 21 3 + + 1/98 x y - 1/20 x y, COEFF(2, 1) = --- y - 1/30 y + 1/35 x y + 200 + + 2 71 2 3 2 2 + + 1/20 x + ---- x - 1/70 y + 1/30 x + 1/84 + 1/40 x y - 1/49 x y , + 2940 + + 2 2 2 11 2 + COEFF(-2, 2) = - 1/20 x - 1/20 x y - 1/35 y + --- y + 2/49 x y + 150 + + 3 3 289 31 3 + - 1/60 y + 1/30 x + 2/35 x y - ---- x + ---, COEFF(-1, 2) = - 2/15 x + 2940 420 + + 2 3 2 21 2 226 + - 1/35 x y + 1/10 x - 1/120 y + 1/10 x y - --- y - 1/49 x y + --- x + 200 735 + + 2 11 2 41 3 + - 1/70 y - ---, COEFF(0, 2) = - 2/35 x y - 2/49 x y - -- x + 1/5 x , + 105 98 + + 11 2 2 3 2 + COEFF(1, 2) = --- - 1/10 x - 1/49 x y - 2/15 x - 1/35 x y + 1/70 y + 105 + + 21 2 3 226 2 2 + + --- y - 1/10 x y + 1/120 y + --- x, COEFF(2, 2) = 1/20 x y + 1/35 y + 200 735 + + 31 2 11 3 289 2 + + 2/35 x y - --- + 1/20 x - --- y + 1/30 x - ---- x + 2/49 x y + 420 150 2940 + + 3 + + 1/60 y ] + +> print_coeff__lc_of_data(%, "coeff_dx_", "fp", +> "2d.coeffs/2d.cube.order4.smooth0/coeff-dx.compute.c"); +bytes used=71064040, alloc=2424388, time=6.54 +bytes used=72064272, alloc=2424388, time=6.62 +bytes used=73075772, alloc=2424388, time=6.70 +bytes used=74075968, alloc=2424388, time=6.88 +bytes used=75076136, alloc=2424388, time=6.95 +bytes used=76076304, alloc=2424388, time=7.04 +bytes used=77076556, alloc=2424388, time=7.20 +bytes used=78076920, alloc=2424388, time=7.29 +bytes used=79077292, alloc=2424388, time=7.46 +bytes used=80077492, alloc=2424388, time=7.59 +bytes used=81077684, alloc=2424388, time=7.72 +> +# d/dy +> simplify( diff(interp_2d_cube_order4_smooth0,y) ); +bytes used=82108412, alloc=2424388, time=7.80 +bytes used=83108600, alloc=2424388, time=7.88 + 17 31 2 +- 2/35 x y DATA(2, 0) + --- DATA(0, 1) + --- DATA(1, 1) + 1/35 x DATA(0, -2) + 105 210 + + 53 166 166 57 + - -- y DATA(1, 0) + --- y DATA(0, 1) + --- y DATA(0, -1) - -- y DATA(0, 0) + 98 735 735 98 + + 191 2 2 + + ---- y DATA(0, -2) + 1/20 y DATA(0, 2) - 1/20 y DATA(0, -2) + 2940 + + 3 3 3 + + 1/5 y DATA(-1, 0) + 1/5 y DATA(2, 0) + 1/5 y DATA(-2, 0) + + 53 41 41 + - -- y DATA(-1, 0) - -- y DATA(2, 0) - -- y DATA(-2, 0) + 98 98 98 + + 191 2 2 + + ---- y DATA(0, 2) - 1/10 y DATA(0, 1) + 1/10 y DATA(0, -1) + 2940 + + 2 2 2 + - 2/49 x y DATA(0, -2) + 1/70 x DATA(0, -1) - 1/35 x DATA(0, 2) + + 3 3 3 + + 1/30 y DATA(0, 2) + 1/30 y DATA(0, -2) - 2/15 y DATA(0, -1) + + 3 3 3 + + 1/5 y DATA(0, 0) + 1/5 y DATA(1, 0) - 2/15 y DATA(0, 1) + + 31 17 31 31 + - --- DATA(1, -1) - --- DATA(0, -1) - --- DATA(-1, -1) + --- DATA(-1, 1) + 210 105 210 210 + + 31 11 + + 1/84 DATA(1, 2) - --- DATA(2, 2) + --- DATA(2, 1) + 1/84 DATA(-1, 2) + 420 105 + + 17 11 31 11 + + --- DATA(0, 2) - --- DATA(2, -1) - --- DATA(-2, 2) + --- DATA(-2, 1) + 420 105 420 105 + + 2 2 2 + - 1/10 y DATA(2, 1) + 1/20 y DATA(-1, 2) + 1/10 y DATA(2, -1) + + 2 2 2 + + 1/20 y DATA(-2, 2) - 1/10 y DATA(-2, 1) + 1/10 y DATA(-2, -1) + + 2 2 2 + - 1/20 y DATA(2, -2) - 1/20 y DATA(-1, -2) + 1/20 y DATA(2, 2) + + 2 2 2 + - 1/10 y DATA(1, 1) + 1/10 y DATA(1, -1) + 1/10 y DATA(-1, -1) + + 2 2 + - 1/10 y DATA(-1, 1) - 1/35 x y DATA(-1, 2) + 1/98 x y DATA(-1, -1) + + 2 + + 1/98 x y DATA(1, 1) + 1/35 x y DATA(-2, -1) + 2/35 x y DATA(2, -2) + + 181 + + 2/35 x y DATA(2, 2) - 1/35 x y DATA(2, 1) + --- y DATA(1, 1) + 735 + + 2 2 2 + - 1/10 x y DATA(-2, -1) - 1/20 x y DATA(2, -2) + 1/40 x y DATA(-1, -2) + + 2 + - 1/40 x y DATA(1, -2) - 1/35 x y DATA(2, -1) - 2/35 x y DATA(-2, 2) + + 71 226 + + 1/35 x y DATA(-2, 1) + ---- y DATA(-1, 2) + --- y DATA(2, -1) + 2940 735 + + 289 226 226 + - ---- y DATA(-2, 2) + --- y DATA(-2, 1) + --- y DATA(-2, -1) + 2940 735 735 + + 289 71 71 + - ---- y DATA(2, -2) + ---- y DATA(-1, -2) + ---- y DATA(1, -2) + 2940 2940 2940 + + 289 181 2 + - ---- y DATA(-2, -2) + --- y DATA(1, -1) - 1/20 y DATA(1, -2) + 2940 735 + + 2 2 + - 1/20 y DATA(-2, -2) + 1/20 y DATA(1, 2) - 1/35 x y DATA(-1, -2) + + 2 + + 1/35 x y DATA(1, -2) - 2/35 x y DATA(-2, -2) - 1/10 x y DATA(2, 1) + + 2 2 2 + + 1/10 x y DATA(2, -1) - 1/20 x y DATA(-2, 2) + 1/10 x y DATA(-2, 1) + + 181 3 3 + + --- y DATA(-1, -1) - 2/15 y DATA(-2, 1) - 2/15 y DATA(-1, 1) + 735 + + 3 3 181 + - 2/15 y DATA(1, 1) - 2/15 y DATA(1, -1) + --- y DATA(-1, 1) + 735 + + 71 289 226 + + ---- y DATA(1, 2) - ---- y DATA(2, 2) + --- y DATA(2, 1) + 2940 2940 735 + + 3 21 11 + + 1/30 y DATA(-2, -2) - --- x DATA(2, -1) + --- x DATA(-2, 2) + 200 150 + + 21 21 11 + - --- x DATA(-2, 1) + --- x DATA(-2, -1) + --- x DATA(2, -2) + 200 200 150 + + 21 21 11 + + --- x DATA(-1, -2) - --- x DATA(1, -2) - --- x DATA(-2, -2) + 200 200 150 + + 3 3 3 + + 1/30 y DATA(2, -2) - 2/15 y DATA(-1, -1) + 1/30 y DATA(1, 2) + + 3 3 3 + + 1/30 y DATA(2, 2) - 2/15 y DATA(2, 1) + 1/30 y DATA(-1, 2) + + 3 3 21 + - 2/15 y DATA(2, -1) + 1/30 y DATA(-2, 2) + --- x DATA(2, 1) + 200 + + 3 3 3 + - 2/15 y DATA(-2, -1) + 1/30 y DATA(-1, -2) + 1/30 y DATA(1, -2) + + 11 3 3 + - --- DATA(-2, -1) + 1/30 x DATA(1, -2) + 1/60 x DATA(-2, -2) + 105 + + 21 37 37 + - --- x DATA(-1, 2) + --- x DATA(1, 1) - --- x DATA(1, -1) + 200 300 300 + + 37 37 21 + + --- x DATA(-1, -1) - --- x DATA(-1, 1) + --- x DATA(1, 2) + 300 300 200 + + 11 3 3 + - --- x DATA(2, 2) + 1/30 x DATA(-1, 2) - 1/60 x DATA(1, 1) + 150 + + 3 3 3 + + 1/60 x DATA(1, -1) - 1/60 x DATA(-1, -1) + 1/60 x DATA(-1, 1) + + 3 3 3 + + 1/120 x DATA(2, 1) - 1/120 x DATA(2, -1) - 1/60 x DATA(-2, 2) + + 3 3 3 + - 1/120 x DATA(-2, 1) + 1/120 x DATA(-2, -1) - 1/60 x DATA(2, -2) + + 3 2 3 + - 1/30 x DATA(-1, -2) - 1/35 x DATA(-2, -2) - 1/30 x DATA(1, 2) + + 3 2 2 + + 1/60 x DATA(2, 2) - 1/140 x DATA(1, 1) + 1/140 x DATA(1, -1) + + 2 2 2 + + 1/140 x DATA(-1, -1) - 1/140 x DATA(-1, 1) - 1/49 x y DATA(-1, -2) + + 2 2 2 + + 1/98 x y DATA(1, -1) - 1/49 x y DATA(-2, -1) + 1/40 x y DATA(1, 2) + + 2 2 2 + + 1/20 x y DATA(2, 2) - 1/40 x y DATA(-1, 2) - 1/20 x y DATA(1, 1) + + 2 2 2 + + 1/20 x y DATA(1, -1) - 1/20 x y DATA(-1, -1) + 1/20 x y DATA(-1, 1) + + 2 + + 1/20 x y DATA(-2, -2) + 1/35 x y DATA(1, 2) - 1/70 x y DATA(1, -1) + + 2 2 2 + - 1/49 x y DATA(1, -2) - 1/49 x y DATA(-1, 2) - 1/49 x y DATA(-2, 1) + + 31 + - 1/70 x y DATA(1, 1) + --- DATA(2, -2) - 1/84 DATA(-1, -2) + 420 + + 17 2 + - --- DATA(0, -2) - 1/84 DATA(1, -2) - 1/49 x y DATA(2, -1) + 420 + + 2 2 2 + + 2/49 x y DATA(-2, -2) - 1/49 x y DATA(1, 2) - 1/49 x y DATA(2, 1) + + 2 2 2 + + 2/49 x y DATA(2, 2) + 2/49 x y DATA(2, -2) + 2/49 x y DATA(-2, 2) + + 2 + + 1/98 x y DATA(-1, 1) + 1/70 x y DATA(-1, -1) + 1/70 x y DATA(-1, 1) + + 2 2 2 + - 1/70 x DATA(1, 2) + 1/70 x DATA(2, 1) - 1/70 x DATA(-1, 2) + + 2 2 2 + - 1/70 x DATA(2, -1) + 1/35 x DATA(2, 2) + 1/35 x DATA(-2, 2) + + 2 2 2 + + 1/70 x DATA(-2, 1) - 1/70 x DATA(-2, -1) - 1/35 x DATA(2, -2) + + 2 2 31 + + 1/70 x DATA(-1, -2) + 1/70 x DATA(1, -2) + --- DATA(-2, -2) + 420 + + 2 2 + + 2/35 x y DATA(-2, 0) + 1/49 x y DATA(-1, 0) - 2/49 x y DATA(2, 0) + + 2 2 + - 2/49 x y DATA(-2, 0) - 2/49 x y DATA(0, 2) - 1/35 x y DATA(1, 0) + + 2 2 + + 1/35 x y DATA(-1, 0) + 1/49 x y DATA(0, -1) + 2/49 x y DATA(0, 0) + + 2 2 2 + + 1/49 x y DATA(1, 0) + 1/49 x y DATA(0, 1) - 1/70 x DATA(0, 1) + +> coeff_as_lc_of_data(%, posn_list_2d_size5); +bytes used=84109804, alloc=2424388, time=7.95 + 31 3 2 2 11 3 +[COEFF(-2, -2) = --- + 1/30 y - 1/35 x - 1/20 y - --- x + 1/60 x + 420 150 + + 2 2 289 2 + + 2/49 x y + 1/20 x y - ---- y - 2/35 x y, COEFF(-1, -2) = - 1/20 y + 2940 + + 21 2 2 2 71 3 + + --- x - 1/49 x y + 1/40 x y + 1/70 x - 1/84 + ---- y - 1/30 x + 200 2940 + + 3 + - 1/35 x y + 1/30 y , + + 191 17 2 2 3 2 + COEFF(0, -2) = ---- y - --- - 2/49 x y + 1/35 x + 1/30 y - 1/20 y , + 2940 420 + + 3 2 2 21 + COEFF(1, -2) = 1/35 x y + 1/30 x - 1/49 x y - 1/84 - 1/20 y - --- x + 200 + + 2 2 71 3 3 + + 1/70 x - 1/40 x y + ---- y + 1/30 y , COEFF(2, -2) = - 1/60 x + 2940 + + 2 31 11 289 2 2 + + 2/35 x y - 1/35 x + --- + --- x - ---- y - 1/20 y + 2/49 x y + 420 150 2940 + + 2 3 21 2 3 226 + - 1/20 x y + 1/30 y , COEFF(-2, -1) = --- x - 1/10 x y - 2/15 y + --- y + 200 735 + + 3 2 2 11 2 + + 1/120 x + 1/35 x y + 1/10 y - 1/70 x - --- - 1/49 x y, COEFF(-1, -1) + 105 + + 2 37 2 181 2 2 + = - 1/20 x y + 1/70 x y + --- x + 1/10 y + --- y + 1/98 x y + 1/140 x + 300 735 + + 3 3 31 + - 2/15 y - 1/60 x - ---, + 210 + + 17 2 166 2 2 3 + COEFF(0, -1) = - --- + 1/49 x y + --- y + 1/10 y + 1/70 x - 2/15 y , + 105 735 + + 181 2 2 2 31 3 + COEFF(1, -1) = --- y + 1/140 x + 1/98 x y + 1/10 y - --- - 2/15 y + 735 210 + + 37 3 2 11 + - --- x + 1/60 x - 1/70 x y + 1/20 x y , COEFF(2, -1) = - 1/35 x y - --- + 300 105 + + 2 21 226 3 3 2 2 + + 1/10 x y - --- x + --- y - 1/120 x - 2/15 y - 1/49 x y + 1/10 y + 200 735 + + 2 3 41 2 + - 1/70 x , COEFF(-2, 0) = 1/5 y - -- y + 2/35 x y - 2/49 x y, + 98 + + 3 53 2 + COEFF(-1, 0) = 1/5 y - -- y + 1/49 x y + 1/35 x y, + 98 + + 2 57 3 + COEFF(0, 0) = 2/49 x y - -- y + 1/5 y , + 98 + + 3 53 2 + COEFF(1, 0) = 1/5 y - -- y + 1/49 x y - 1/35 x y, + 98 + + 3 2 41 + COEFF(2, 0) = 1/5 y - 2/35 x y - 2/49 x y - -- y, COEFF(-2, 1) = + 98 + + 3 3 11 2 2 226 + - 1/120 x - 2/15 y + --- - 1/49 x y + 1/70 x + 1/35 x y + --- y + 105 735 + + 2 2 21 37 2 + + 1/10 x y - 1/10 y - --- x, COEFF(-1, 1) = - --- x + 1/20 x y + 200 300 + + 3 31 2 2 3 2 + + 1/70 x y + 1/60 x + --- - 1/140 x - 1/10 y - 2/15 y + 1/98 x y + 210 + + 181 + + --- y, + 735 + + 17 3 2 166 2 2 + COEFF(0, 1) = --- - 2/15 y - 1/70 x + --- y - 1/10 y + 1/49 x y, + 105 735 + + 31 3 3 2 2 + COEFF(1, 1) = --- - 1/60 x - 2/15 y - 1/70 x y - 1/20 x y - 1/10 y + 210 + + 181 2 37 2 2 21 + + --- y + 1/98 x y + --- x - 1/140 x , COEFF(2, 1) = - 1/49 x y + --- x + 735 300 200 + + 226 2 2 11 2 3 + + --- y - 1/10 y + 1/70 x - 1/35 x y + --- - 1/10 x y + 1/120 x + 735 105 + + 3 3 2 289 31 3 + - 2/15 y , COEFF(-2, 2) = - 1/60 x + 1/35 x - ---- y - --- + 1/30 y + 2940 420 + + 11 2 2 2 + + --- x + 2/49 x y - 2/35 x y - 1/20 x y + 1/20 y , COEFF(-1, 2) = + 150 + + 21 2 2 71 2 2 + - --- x - 1/49 x y + 1/20 y + ---- y - 1/40 x y - 1/70 x + 1/84 + 200 2940 + + 3 3 + - 1/35 x y + 1/30 x + 1/30 y , + + 2 17 191 2 3 2 + COEFF(0, 2) = 1/20 y + --- + ---- y - 1/35 x + 1/30 y - 2/49 x y, + 420 2940 + + 2 2 3 2 21 + COEFF(1, 2) = - 1/70 x + 1/35 x y - 1/49 x y - 1/30 x + 1/20 y + --- x + 200 + + 3 71 2 3 289 + + 1/30 y + ---- y + 1/84 + 1/40 x y , COEFF(2, 2) = 1/60 x - ---- y + 2940 2940 + + 2 2 3 31 2 2 11 + + 1/35 x + 1/20 y + 1/30 y - --- + 2/49 x y + 1/20 x y - --- x + 420 150 + + + 2/35 x y] + +> print_coeff__lc_of_data(%, "coeff_dy_", "fp", +> "2d.coeffs/2d.cube.order4.smooth0/coeff-dy.compute.c"); +bytes used=85110104, alloc=2424388, time=8.01 +bytes used=86110296, alloc=2424388, time=8.09 +bytes used=87116120, alloc=2424388, time=8.16 +bytes used=88116280, alloc=2424388, time=8.30 +bytes used=89116760, alloc=2424388, time=8.38 +bytes used=90116912, alloc=2424388, time=8.48 +bytes used=91117288, alloc=2424388, time=8.64 +bytes used=92117484, alloc=2424388, time=8.73 +bytes used=93117892, alloc=2424388, time=8.88 +bytes used=94118116, alloc=2424388, time=8.97 +bytes used=95118268, alloc=2424388, time=9.12 +> +# d^2/dx^2 +> simplify( diff(interp_2d_cube_order4_smooth0,x,x) ); +bytes used=96118472, alloc=2424388, time=9.19 + 57 166 53 181 +- -- DATA(0, 0) + --- DATA(1, 0) - -- DATA(0, 1) + --- DATA(1, 1) + 98 735 98 735 + + 2 + + 3/5 x DATA(0, -2) - 1/35 y DATA(0, 1) + 1/35 y DATA(0, -1) + + 2 2 + + 2/35 y DATA(0, -2) - 2/49 y DATA(0, 2) - 2/49 y DATA(0, -2) + + 2 2 + - 2/35 y DATA(0, 2) + 1/49 y DATA(0, 1) + 1/49 y DATA(0, -1) + + 2 2 166 + + 3/5 x DATA(0, -1) + 3/5 x DATA(0, 2) + --- DATA(-1, 0) + 735 + + 181 53 181 181 + + --- DATA(1, -1) - -- DATA(0, -1) + --- DATA(-1, -1) + --- DATA(-1, 1) + 735 98 735 735 + + 2 2 + - 2/49 y DATA(-2, 0) - 1/10 x DATA(-2, 0) - 2/49 y DATA(2, 0) + + 2 2 2 + + 2/49 y DATA(0, 0) + 1/49 y DATA(1, 0) + 1/49 y DATA(-1, 0) + + 2 2 + + 1/10 x DATA(2, 0) + 3/5 x DATA(0, 0) - 2/5 x DATA(1, 0) + + 2 2 2 + + 1/10 x DATA(-2, 0) - 2/5 x DATA(-1, 0) + 1/10 x DATA(2, 0) + + 226 289 + - 1/5 x DATA(1, 0) + 1/5 x DATA(-1, 0) + --- DATA(1, 2) - ---- DATA(2, 2) + 735 2940 + + 191 71 226 41 + + ---- DATA(2, 0) + ---- DATA(2, 1) + --- DATA(-1, 2) - -- DATA(0, 2) + 2940 2940 735 98 + + 71 289 71 + + ---- DATA(2, -1) - ---- DATA(-2, 2) + ---- DATA(-2, 1) + 2940 2940 2940 + + 2 2 2 + - 1/49 y DATA(2, 1) - 1/49 y DATA(-1, 2) - 1/49 y DATA(2, -1) + + 2 2 2 + + 2/49 y DATA(-2, 2) - 1/49 y DATA(-2, 1) - 1/49 y DATA(-2, -1) + + 2 2 2 + + 2/49 y DATA(2, -2) - 1/49 y DATA(-1, -2) + 2/49 y DATA(2, 2) + + 2 2 2 + + 1/98 y DATA(1, 1) + 1/98 y DATA(1, -1) + 1/98 y DATA(-1, -1) + + 2 + + 1/98 y DATA(-1, 1) + 1/5 x y DATA(-1, 2) + 1/20 x y DATA(-2, -1) + + - 1/10 x y DATA(2, -2) + 1/10 x y DATA(2, 2) + 1/20 x y DATA(2, 1) + + - 1/70 y DATA(1, 1) - 1/20 x y DATA(2, -1) - 1/10 x y DATA(-2, 2) + + - 1/20 x y DATA(-2, 1) - 1/35 y DATA(-1, 2) - 1/35 y DATA(2, -1) + + + 2/35 y DATA(-2, 2) + 1/35 y DATA(-2, 1) - 1/35 y DATA(-2, -1) + + - 2/35 y DATA(2, -2) + 1/35 y DATA(-1, -2) + 1/35 y DATA(1, -2) + + 2 + - 2/35 y DATA(-2, -2) + 1/70 y DATA(1, -1) - 1/49 y DATA(1, -2) + + 2 2 + + 2/49 y DATA(-2, -2) - 1/49 y DATA(1, 2) - 1/5 x y DATA(-1, -2) + + + 1/5 x y DATA(1, -2) + 1/10 x y DATA(-2, -2) + 1/70 y DATA(-1, -1) + + - 1/70 y DATA(-1, 1) - 1/35 y DATA(1, 2) + 2/35 y DATA(2, 2) + + + 1/35 y DATA(2, 1) + 1/10 x DATA(2, -1) - 1/10 x DATA(-2, 2) + + - 1/10 x DATA(-2, 1) - 1/10 x DATA(-2, -1) + 1/10 x DATA(2, -2) + + + 1/5 x DATA(-1, -2) - 1/5 x DATA(1, -2) - 1/10 x DATA(-2, -2) + + 71 191 + + 1/10 x DATA(2, 1) + ---- DATA(-2, -1) + ---- DATA(-2, 0) + 2940 2940 + + + 1/5 x DATA(-1, 2) - 1/5 x DATA(1, 1) - 1/5 x DATA(1, -1) + + + 1/5 x DATA(-1, -1) + 1/5 x DATA(-1, 1) - 1/5 x DATA(1, 2) + + 2 2 + + 1/10 x DATA(2, 2) + 1/10 x DATA(-2, -2) - 2/5 x DATA(1, 1) + + 2 2 2 + - 2/5 x DATA(1, -1) - 2/5 x DATA(-1, -1) - 2/5 x DATA(-1, 1) + + - 1/5 x y DATA(1, 2) + 1/10 x y DATA(1, -1) - 1/10 x y DATA(1, 1) + + 289 226 41 226 + - ---- DATA(2, -2) + --- DATA(-1, -2) - -- DATA(0, -2) + --- DATA(1, -2) + 2940 735 98 735 + + 2 + - 1/10 x y DATA(-1, -1) + 1/10 x y DATA(-1, 1) - 2/5 x DATA(1, 2) + + 2 2 2 + + 1/10 x DATA(2, 1) - 2/5 x DATA(-1, 2) + 1/10 x DATA(2, -1) + + 2 2 2 + + 1/10 x DATA(2, 2) + 1/10 x DATA(-2, 2) + 1/10 x DATA(-2, 1) + + 2 2 2 + + 1/10 x DATA(-2, -1) + 1/10 x DATA(2, -2) - 2/5 x DATA(-1, -2) + + 2 289 2 + - 2/5 x DATA(1, -2) - ---- DATA(-2, -2) + 3/5 x DATA(0, 1) + 2940 + +> coeff_as_lc_of_data(%, posn_list_2d_size5); +bytes used=97118648, alloc=2424388, time=9.26 + 2 289 2 +[COEFF(-2, -2) = 1/10 x y - 2/35 y + 2/49 y - ---- + 1/10 x - 1/10 x, + 2940 + + 2 226 2 + COEFF(-1, -2) = - 1/5 x y - 2/5 x + 1/35 y + --- + 1/5 x - 1/49 y , + 735 + + 2 2 41 + COEFF(0, -2) = 2/35 y - 2/49 y + 3/5 x - --, + 98 + + 226 2 2 + COEFF(1, -2) = --- - 2/5 x + 1/35 y + 1/5 x y - 1/5 x - 1/49 y , + 735 + + 2 289 2 + COEFF(2, -2) = 2/49 y - ---- - 2/35 y + 1/10 x + 1/10 x - 1/10 x y, + 2940 + + 71 2 2 + COEFF(-2, -1) = - 1/10 x + ---- - 1/49 y + 1/10 x - 1/35 y + 1/20 x y, + 2940 + + 2 2 181 + COEFF(-1, -1) = - 2/5 x + 1/98 y + 1/5 x + --- - 1/10 x y + 1/70 y, + 735 + + 53 2 2 + COEFF(0, -1) = - -- + 1/35 y + 1/49 y + 3/5 x , + 98 + + 2 2 181 + COEFF(1, -1) = 1/70 y + 1/98 y - 2/5 x - 1/5 x + 1/10 x y + ---, + 735 + + 2 2 71 + COEFF(2, -1) = 1/10 x - 1/35 y - 1/49 y + 1/10 x + ---- - 1/20 x y, + 2940 + + 2 2 191 + COEFF(-2, 0) = - 1/10 x + 1/10 x - 2/49 y + ----, + 2940 + + 2 2 166 + COEFF(-1, 0) = - 2/5 x + 1/49 y + --- + 1/5 x, + 735 + + 2 57 2 + COEFF(0, 0) = 2/49 y - -- + 3/5 x , + 98 + + 2 2 166 + COEFF(1, 0) = - 1/5 x + 1/49 y - 2/5 x + ---, + 735 + + 2 191 2 + COEFF(2, 0) = 1/10 x - 2/49 y + ---- + 1/10 x , + 2940 + + 71 2 2 + COEFF(-2, 1) = ---- - 1/49 y + 1/10 x + 1/35 y - 1/20 x y - 1/10 x, + 2940 + + 181 2 2 + COEFF(-1, 1) = 1/5 x + --- + 1/10 x y - 2/5 x + 1/98 y - 1/70 y, + 735 + + 53 2 2 + COEFF(0, 1) = - -- + 3/5 x + 1/49 y - 1/35 y, + 98 + + 181 2 2 + COEFF(1, 1) = --- + 1/98 y - 2/5 x - 1/70 y - 1/5 x - 1/10 x y, + 735 + + 2 71 2 + COEFF(2, 1) = 1/10 x + 1/35 y + 1/10 x + ---- + 1/20 x y - 1/49 y , + 2940 + + 2 289 2 + COEFF(-2, 2) = 2/49 y + 2/35 y - ---- - 1/10 x y - 1/10 x + 1/10 x , + 2940 + + 226 2 2 + COEFF(-1, 2) = --- - 2/5 x + 1/5 x + 1/5 x y - 1/49 y - 1/35 y, + 735 + + 41 2 2 + COEFF(0, 2) = - -- - 2/35 y - 2/49 y + 3/5 x , + 98 + + 2 226 2 + COEFF(1, 2) = - 1/5 x - 2/5 x - 1/5 x y + --- - 1/49 y - 1/35 y, + 735 + + 2 289 2 + COEFF(2, 2) = 1/10 x + 2/49 y - ---- + 1/10 x + 2/35 y + 1/10 x y] + 2940 + +> print_coeff__lc_of_data(%, "coeff_dxx_", "fp", +> "2d.coeffs/2d.cube.order4.smooth0/coeff-dxx.compute.c"); +bytes used=98118800, alloc=2424388, time=9.34 +bytes used=99118996, alloc=2424388, time=9.43 +bytes used=100119256, alloc=2424388, time=9.58 +bytes used=101122084, alloc=2424388, time=9.74 +bytes used=102122280, alloc=2424388, time=9.85 +bytes used=103122400, alloc=2424388, time=9.97 +> +# d^2/dxdy +> simplify( diff(interp_2d_cube_order4_smooth0,x,y) ); +bytes used=104122552, alloc=2424388, time=10.13 +bytes used=105123136, alloc=2424388, time=10.20 +37 +--- DATA(1, 1) - 4/49 x y DATA(-2, 0) - 4/49 x y DATA(0, -2) +300 + + + 4/49 x y DATA(0, 0) + 2/49 x y DATA(0, -1) - 4/49 x y DATA(0, 2) + + + 2/49 x y DATA(0, 1) + 2/49 x y DATA(1, 0) - 4/49 x y DATA(2, 0) + + + 2/49 x y DATA(-1, 0) + 1/35 y DATA(-1, 0) + 2/35 y DATA(-2, 0) + + 37 37 37 + - --- DATA(1, -1) + --- DATA(-1, -1) - --- DATA(-1, 1) - 2/35 y DATA(2, 0) + 300 300 300 + + 21 11 + - 1/35 y DATA(1, 0) - 1/35 x DATA(0, 1) + --- DATA(1, 2) - --- DATA(2, 2) + 200 150 + + 21 21 21 11 + + --- DATA(2, 1) - --- DATA(-1, 2) - --- DATA(2, -1) + --- DATA(-2, 2) + 200 200 200 150 + + 21 2 2 + - --- DATA(-2, 1) - 1/10 y DATA(2, 1) - 1/40 y DATA(-1, 2) + 200 + + 2 2 2 + + 1/10 y DATA(2, -1) - 1/20 y DATA(-2, 2) + 1/10 y DATA(-2, 1) + + 2 2 2 + - 1/10 y DATA(-2, -1) - 1/20 y DATA(2, -2) + 1/40 y DATA(-1, -2) + + 2 2 2 + + 1/20 y DATA(2, 2) - 1/20 y DATA(1, 1) + 1/20 y DATA(1, -1) + + 2 2 + - 1/20 y DATA(-1, -1) + 1/20 y DATA(-1, 1) - 2/49 x y DATA(-1, 2) + + - 2/49 x y DATA(-2, -1) + 4/49 x y DATA(2, -2) + 4/49 x y DATA(2, 2) + + - 2/49 x y DATA(2, 1) - 1/70 y DATA(1, 1) - 2/49 x y DATA(2, -1) + + + 4/49 x y DATA(-2, 2) - 2/49 x y DATA(-2, 1) - 1/35 y DATA(-1, 2) + + - 1/35 y DATA(2, -1) - 2/35 y DATA(-2, 2) + 1/35 y DATA(-2, 1) + + + 1/35 y DATA(-2, -1) + 2/35 y DATA(2, -2) - 1/35 y DATA(-1, -2) + + + 1/35 y DATA(1, -2) - 2/35 y DATA(-2, -2) - 1/70 y DATA(1, -1) + + 2 2 2 + - 1/40 y DATA(1, -2) + 1/20 y DATA(-2, -2) + 1/40 y DATA(1, 2) + + - 2/49 x y DATA(-1, -2) - 2/49 x y DATA(1, -2) + 4/49 x y DATA(-2, -2) + + + 1/70 y DATA(-1, -1) + 1/70 y DATA(-1, 1) + 1/35 y DATA(1, 2) + + + 2/35 y DATA(2, 2) - 1/35 y DATA(2, 1) - 1/35 x DATA(2, -1) + + + 2/35 x DATA(-2, 2) + 1/35 x DATA(-2, 1) - 1/35 x DATA(-2, -1) + + - 2/35 x DATA(2, -2) + 1/35 x DATA(-1, -2) + 1/35 x DATA(1, -2) + + 21 + - 2/35 x DATA(-2, -2) + 1/35 x DATA(2, 1) + --- DATA(-2, -1) + 200 + + - 1/35 x DATA(-1, 2) - 1/70 x DATA(1, 1) + 1/70 x DATA(1, -1) + + + 1/70 x DATA(-1, -1) - 1/70 x DATA(-1, 1) - 1/35 x DATA(1, 2) + + 2 2 + + 2/35 x DATA(2, 2) + 1/20 x DATA(-2, -2) - 1/20 x DATA(1, 1) + + 2 2 2 + + 1/20 x DATA(1, -1) - 1/20 x DATA(-1, -1) + 1/20 x DATA(-1, 1) + + - 2/49 x y DATA(1, 2) + 1/49 x y DATA(1, -1) + 1/49 x y DATA(1, 1) + + 11 21 21 + + --- DATA(2, -2) + --- DATA(-1, -2) - --- DATA(1, -2) + 150 200 200 + + 2 + + 1/49 x y DATA(-1, -1) + 1/49 x y DATA(-1, 1) - 1/10 x DATA(1, 2) + + 2 2 2 + + 1/40 x DATA(2, 1) + 1/10 x DATA(-1, 2) - 1/40 x DATA(2, -1) + + 2 2 2 + + 1/20 x DATA(2, 2) - 1/20 x DATA(-2, 2) - 1/40 x DATA(-2, 1) + + 2 2 2 + + 1/40 x DATA(-2, -1) - 1/20 x DATA(2, -2) - 1/10 x DATA(-1, -2) + + 2 11 + + 1/10 x DATA(1, -2) - --- DATA(-2, -2) - 2/35 x DATA(0, 2) + 150 + + + 1/35 x DATA(0, -1) + 2/35 x DATA(0, -2) + +> coeff_as_lc_of_data(%, posn_list_2d_size5); + 2 2 11 +[COEFF(-2, -2) = 4/49 x y + 1/20 x + 1/20 y - 2/35 y - 2/35 x - ---, + 150 + + 2 2 21 + COEFF(-1, -2) = - 2/49 x y - 1/35 y + 1/35 x + 1/40 y - 1/10 x + ---, + 200 + + COEFF(0, -2) = 2/35 x - 4/49 x y, + + 2 21 2 + COEFF(1, -2) = - 2/49 x y + 1/35 x + 1/10 x + 1/35 y - --- - 1/40 y , + 200 + + 11 2 2 + COEFF(2, -2) = --- - 1/20 y - 1/20 x - 2/35 x + 4/49 x y + 2/35 y, + 150 + + 21 2 2 + COEFF(-2, -1) = --- + 1/35 y - 1/10 y - 2/49 x y - 1/35 x + 1/40 x , + 200 + + 2 2 37 + COEFF(-1, -1) = 1/49 x y - 1/20 y - 1/20 x + 1/70 x + 1/70 y + ---, + 300 + + COEFF(0, -1) = 1/35 x + 2/49 x y, + + 37 2 2 + COEFF(1, -1) = - --- + 1/49 x y + 1/20 x - 1/70 y + 1/70 x + 1/20 y , + 300 + + 2 2 21 + COEFF(2, -1) = - 1/35 x + 1/10 y - 2/49 x y - 1/35 y - 1/40 x - ---, + 200 + + COEFF(-2, 0) = 2/35 y - 4/49 x y, COEFF(-1, 0) = 1/35 y + 2/49 x y, + + COEFF(0, 0) = 4/49 x y, COEFF(1, 0) = - 1/35 y + 2/49 x y, + + COEFF(2, 0) = - 4/49 x y - 2/35 y, + + 2 2 21 + COEFF(-2, 1) = 1/35 y + 1/35 x - 2/49 x y + 1/10 y - 1/40 x - ---, + 200 + + 2 37 2 + COEFF(-1, 1) = 1/20 y - --- - 1/70 x + 1/20 x + 1/70 y + 1/49 x y, + 300 + + COEFF(0, 1) = - 1/35 x + 2/49 x y, + + 37 2 2 + COEFF(1, 1) = 1/49 x y + --- - 1/20 y - 1/70 x - 1/20 x - 1/70 y, + 300 + + 21 2 2 + COEFF(2, 1) = --- + 1/40 x - 1/10 y + 1/35 x - 1/35 y - 2/49 x y, + 200 + + 2 2 11 + COEFF(-2, 2) = - 1/20 x + 4/49 x y + 2/35 x - 1/20 y - 2/35 y + ---, + 150 + + 2 2 21 + COEFF(-1, 2) = 1/10 x - 1/35 y - 1/35 x - 1/40 y - --- - 2/49 x y, + 200 + + COEFF(0, 2) = - 2/35 x - 4/49 x y, + + 2 2 21 + COEFF(1, 2) = - 1/35 x + 1/40 y - 2/49 x y + 1/35 y - 1/10 x + ---, + 200 + + 2 11 2 + COEFF(2, 2) = 2/35 y + 4/49 x y + 1/20 x - --- + 2/35 x + 1/20 y ] + 150 + +> print_coeff__lc_of_data(%, "coeff_dxy_", "fp", +> "2d.coeffs/2d.cube.order4.smooth0/coeff-dxy.compute.c"); +bytes used=106123844, alloc=2424388, time=10.27 +bytes used=107128472, alloc=2424388, time=10.35 +bytes used=108128708, alloc=2424388, time=10.53 +bytes used=109128920, alloc=2424388, time=10.69 +bytes used=110129208, alloc=2424388, time=10.83 +bytes used=111129472, alloc=2424388, time=11.03 +bytes used=112129776, alloc=2424388, time=11.21 +> +# d^2/dy^2 +> simplify( diff(interp_2d_cube_order4_smooth0,y,y) ); +bytes used=113130192, alloc=2424388, time=11.31 +bytes used=114130376, alloc=2424388, time=11.38 + 57 53 166 181 +- -- DATA(0, 0) - -- DATA(1, 0) + --- DATA(0, 1) + --- DATA(1, 1) + 98 98 735 735 + + 53 181 166 181 + - -- DATA(-1, 0) + --- DATA(1, -1) + --- DATA(0, -1) + --- DATA(-1, -1) + 98 735 735 735 + + 181 71 289 41 + + --- DATA(-1, 1) + ---- DATA(1, 2) - ---- DATA(2, 2) - -- DATA(2, 0) + 735 2940 2940 98 + + 226 71 191 226 + + --- DATA(2, 1) + ---- DATA(-1, 2) + ---- DATA(0, 2) + --- DATA(2, -1) + 735 2940 2940 735 + + 289 226 2 + - ---- DATA(-2, 2) + --- DATA(-2, 1) - 2/5 y DATA(2, 1) + 2940 735 + + 2 2 2 + + 1/10 y DATA(-1, 2) - 2/5 y DATA(2, -1) + 1/10 y DATA(-2, 2) + + 2 2 2 + - 2/5 y DATA(-2, 1) - 2/5 y DATA(-2, -1) + 1/10 y DATA(2, -2) + + 2 2 2 + + 1/10 y DATA(-1, -2) + 1/10 y DATA(2, 2) - 2/5 y DATA(1, 1) + + 2 2 2 + - 2/5 y DATA(1, -1) - 2/5 y DATA(-1, -1) - 2/5 y DATA(-1, 1) + + - 1/20 x y DATA(-1, 2) - 1/5 x y DATA(-2, -1) - 1/10 x y DATA(2, -2) + + + 1/10 x y DATA(2, 2) - 1/5 x y DATA(2, 1) - 1/5 y DATA(1, 1) + + + 1/5 x y DATA(2, -1) - 1/10 x y DATA(-2, 2) + 1/5 x y DATA(-2, 1) + + + 1/10 y DATA(-1, 2) + 1/5 y DATA(2, -1) + 1/10 y DATA(-2, 2) + + - 1/5 y DATA(-2, 1) + 1/5 y DATA(-2, -1) - 1/10 y DATA(2, -2) + + - 1/10 y DATA(-1, -2) - 1/10 y DATA(1, -2) - 1/10 y DATA(-2, -2) + + 2 2 + + 1/5 y DATA(1, -1) + 1/10 y DATA(1, -2) + 1/10 y DATA(-2, -2) + + 2 + + 1/10 y DATA(1, 2) + 1/20 x y DATA(-1, -2) - 1/20 x y DATA(1, -2) + + + 1/10 x y DATA(-2, -2) + 1/5 y DATA(-1, -1) - 1/5 y DATA(-1, 1) + + + 1/10 y DATA(1, 2) + 1/10 y DATA(2, 2) - 1/5 y DATA(2, 1) + + - 1/35 x DATA(2, -1) - 2/35 x DATA(-2, 2) + 1/35 x DATA(-2, 1) + + + 1/35 x DATA(-2, -1) + 2/35 x DATA(2, -2) - 1/35 x DATA(-1, -2) + + + 1/35 x DATA(1, -2) - 2/35 x DATA(-2, -2) - 1/35 x DATA(2, 1) + + 226 41 + + --- DATA(-2, -1) - -- DATA(-2, 0) - 1/35 x DATA(-1, 2) + 735 98 + + - 1/70 x DATA(1, 1) - 1/70 x DATA(1, -1) + 1/70 x DATA(-1, -1) + + + 1/70 x DATA(-1, 1) + 1/35 x DATA(1, 2) + 2/35 x DATA(2, 2) + + 2 2 2 + + 2/49 x DATA(-2, -2) + 1/98 x DATA(1, 1) + 1/98 x DATA(1, -1) + + 2 2 + + 1/98 x DATA(-1, -1) + 1/98 x DATA(-1, 1) + 1/20 x y DATA(1, 2) + + 289 + + 1/10 x y DATA(1, -1) - 1/10 x y DATA(1, 1) - ---- DATA(2, -2) + 2940 + + 71 191 71 + + ---- DATA(-1, -2) + ---- DATA(0, -2) + ---- DATA(1, -2) + 2940 2940 2940 + + 2 + - 1/10 x y DATA(-1, -1) + 1/10 x y DATA(-1, 1) - 1/49 x DATA(1, 2) + + 2 2 2 + - 1/49 x DATA(2, 1) - 1/49 x DATA(-1, 2) - 1/49 x DATA(2, -1) + + 2 2 2 + + 2/49 x DATA(2, 2) + 2/49 x DATA(-2, 2) - 1/49 x DATA(-2, 1) + + 2 2 2 + - 1/49 x DATA(-2, -1) + 2/49 x DATA(2, -2) - 1/49 x DATA(-1, -2) + + 2 289 + - 1/49 x DATA(1, -2) - ---- DATA(-2, -2) - 1/35 x DATA(1, 0) + 2940 + + + 1/35 x DATA(-1, 0) + 1/5 y DATA(0, -1) - 2/35 x DATA(2, 0) + + 2 + + 1/10 y DATA(0, 2) - 1/10 y DATA(0, -2) + 3/5 y DATA(-2, 0) + + 2 2 + + 3/5 y DATA(-1, 0) + 3/5 y DATA(1, 0) + 2/35 x DATA(-2, 0) + + 2 2 2 + + 3/5 y DATA(0, 0) - 2/5 y DATA(0, 1) - 2/5 y DATA(0, -1) + + 2 2 2 + + 1/10 y DATA(0, 2) + 1/10 y DATA(0, -2) + 1/49 x DATA(1, 0) + + 2 2 2 + + 3/5 y DATA(2, 0) - 2/49 x DATA(-2, 0) + 1/49 x DATA(-1, 0) + + 2 2 + - 1/5 y DATA(0, 1) + 2/49 x DATA(0, 0) + 1/49 x DATA(0, 1) + + 2 2 2 + + 1/49 x DATA(0, -1) - 2/49 x DATA(0, 2) - 2/49 x DATA(0, -2) + + 2 + - 2/49 x DATA(2, 0) + +> coeff_as_lc_of_data(%, posn_list_2d_size5); + 2 2 289 +[COEFF(-2, -2) = 1/10 y + 1/10 x y - 2/35 x + 2/49 x - 1/10 y - ----, + 2940 + + 71 2 2 + COEFF(-1, -2) = ---- + 1/10 y + 1/20 x y - 1/10 y - 1/35 x - 1/49 x , + 2940 + + 2 191 2 + COEFF(0, -2) = 1/10 y + ---- - 2/49 x - 1/10 y, + 2940 + + 2 2 71 + COEFF(1, -2) = - 1/10 y + 1/35 x - 1/20 x y - 1/49 x + 1/10 y + ----, + 2940 + + 2 289 2 + COEFF(2, -2) = 2/49 x - ---- + 1/10 y - 1/10 x y + 2/35 x - 1/10 y, + 2940 + + 2 2 226 + COEFF(-2, -1) = 1/5 y - 2/5 y - 1/49 x + 1/35 x + --- - 1/5 x y, + 735 + + 2 181 2 + COEFF(-1, -1) = - 2/5 y - 1/10 x y + 1/5 y + 1/70 x + --- + 1/98 x , + 735 + + 2 2 166 + COEFF(0, -1) = 1/49 x - 2/5 y + --- + 1/5 y, + 735 + + 2 2 181 + COEFF(1, -1) = - 1/70 x + 1/5 y - 2/5 y + 1/10 x y + 1/98 x + ---, + 735 + + 2 226 2 + COEFF(2, -1) = - 1/35 x - 2/5 y + --- - 1/49 x + 1/5 x y + 1/5 y, + 735 + + 41 2 2 + COEFF(-2, 0) = - -- + 2/35 x - 2/49 x + 3/5 y , + 98 + + 2 2 53 + COEFF(-1, 0) = 3/5 y + 1/35 x + 1/49 x - --, + 98 + + 2 57 2 + COEFF(0, 0) = 3/5 y - -- + 2/49 x , + 98 + + 2 2 53 + COEFF(1, 0) = 3/5 y + 1/49 x - 1/35 x - --, + 98 + + 41 2 2 + COEFF(2, 0) = - -- - 2/49 x - 2/35 x + 3/5 y , + 98 + + 226 2 2 + COEFF(-2, 1) = 1/35 x - 1/5 y + --- + 1/5 x y - 2/5 y - 1/49 x , + 735 + + 2 181 2 + COEFF(-1, 1) = - 1/5 y + 1/70 x + 1/10 x y - 2/5 y + --- + 1/98 x , + 735 + + 2 2 166 + COEFF(0, 1) = - 2/5 y + 1/49 x - 1/5 y + ---, + 735 + + 181 2 2 + COEFF(1, 1) = - 1/70 x + --- - 1/10 x y + 1/98 x - 2/5 y - 1/5 y, + 735 + + 226 2 2 + COEFF(2, 1) = --- - 1/5 y - 1/49 x - 1/35 x - 2/5 y - 1/5 x y, + 735 + + 2 289 2 + COEFF(-2, 2) = - 2/35 x + 2/49 x - 1/10 x y - ---- + 1/10 y + 1/10 y , + 2940 + + 71 2 2 + COEFF(-1, 2) = ---- - 1/49 x + 1/10 y - 1/35 x + 1/10 y - 1/20 x y, + 2940 + + 2 2 191 + COEFF(0, 2) = - 2/49 x + 1/10 y + 1/10 y + ----, + 2940 + + 71 2 2 + COEFF(1, 2) = ---- + 1/10 y + 1/20 x y + 1/10 y + 1/35 x - 1/49 x , + 2940 + + 289 2 2 + COEFF(2, 2) = - ---- + 1/10 y + 2/49 x + 2/35 x + 1/10 x y + 1/10 y ] + 2940 + +> print_coeff__lc_of_data(%, "coeff_dyy_", "fp", +> "2d.coeffs/2d.cube.order4.smooth0/coeff-dyy.compute.c"); +bytes used=115130604, alloc=2424388, time=11.46 +bytes used=116138968, alloc=2424388, time=11.55 +bytes used=117139208, alloc=2424388, time=11.72 +bytes used=118140008, alloc=2424388, time=11.80 +bytes used=119140212, alloc=2424388, time=12.00 +bytes used=120140404, alloc=2424388, time=12.21 +bytes used=121141040, alloc=2424388, time=12.38 +> +################################################################################ +> quit +bytes used=121567580, alloc=2424388, time=12.41 |