126 Chapter 3. Interpolation and Extrapolation y(x1,x2) Citi-lw-1 y.1(1,2) 入i-1e-2/ (3.6.6) y.2(x1,T2）= ∑∑0-1e,ti-1-2 du/drz i=1j=1 ,12(x1,t2） ∑i-1)6-1)c,t-2-2(dt/da1)(du/d2) 83g i=1j=】 granted for (including this one) 19881992 where t and u are again given by equation(3.6.4). void bcucof(float y[],float y1[],float y2[],float y12[],float d1,float d2, 111-800872 f1oat米*c) Given arrays y[1..4].y1[1..4],y2[1..4],and y12[1..4].containing the function.gra- dients,and cross derivative at the four grid points of a rectangular grid cell (numbered coun- from NUMERICAL RECIPES IN C: terclockwise from the lower left),and given d1 and d2,the length of the grid cell in the 1-and 2-directions,this routine returns the table c[1..4][1..4]that is used by routine bcuint for bicubic interpolation. (North America to any server computer,is strictly prohibited. tusers to make one paper by Cambridge University Press.Programs Copyright(C) THE static int wt[16][16]= 1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0, -3,0,0,3,0,0,0,0,-2,0,0,-1,0,0,0,0 onty),or 2,0,0,-2,0,0,0,0,1,0,0,1,0,0,0,0 0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0, copy for their 0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0, 0,0,0,0,-3,0,0,3,0,0,0,0,-2,0,0,-1, 0,0,0,0,2,0,0,-2,0,0,0,0,1,0,0,1 -3,3,0,0,-2,-1,0,0,0,0,0,0,0,0,0,0 0,0,0,0,0,0,0,0,-3,3,0,0,-2,-1,0,0 9,-9,9,-9,6,3,-3,-6,6,-6,-3,3,4,2,1,2 -6,6,-6,6,-4,-2,2,4,-3,3,3,-3,-2,-1,-1,-2 2,-2,0,0,1,1,0,0,0,0,0,0,0,0,0,0 0,0,0,0,0,0,0,0,2,-2,0,0,1,1,0,0 -6,6,-6,6,-3,-3,3,3,-4,4,2,-2,-2,-2,-1,-1 4,-4,4,-4,2,2,-2,-2,2,-2,-2,2,1,1,1,1]； ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5) int 1,k,j,i; f1 oat xx,d1d2,c1[16],x[16] personal use.Further reproduction,or 1988-1992 by Numerical Recipes d1d2=d1*d2: for(i=1:1<=4;i++) Pack a temporary vector x. x[i-1]y[1]; x[i+3]=y1[i]*d1; rsend email to directcustserv@cambridge.org(outside North America) Software. x[i+7]=y2[i]*d2; x[i+11]=y12[i]*d1d2; for(1=0;1<=15;1++)[ Matrix multiply by the stored table. Xx=0.0; for(k=0;k<=15;k++)xx+=wt[1][k]*x[k]; c1[i]=x; 1=0; for(1=1;1<=4;1++) Unpack the result into the output table. for(j=1;j<=4;j+)c[1][j]=c1[1++];126 Chapter 3. Interpolation and Extrapolation Permission is granted for internet users to make one paper copy for their own personal use. Further reproduction, or any copyin Copyright (C) 1988-1992 by Cambridge University Press. Programs Copyright (C) 1988-1992 by Numerical Recipes Software. Sample page from NUMERICAL RECIPES IN C: THE ART OF SCIENTIFIC COMPUTING (ISBN 0-521-43108-5) g of machine￾readable files (including this one) to any server computer, is strictly prohibited. To order Numerical Recipes books or CDROMs, visit website http://www.nr.com or call 1-800-872-7423 (North America only), or send email to directcustserv@cambridge.org (outside North America). y(x1, x2) =
4 i=1
4 j=1 cij t i−1uj−1 y,1(x1, x2) =
4 i=1
4 j=1 (i − 1)cij t i−2uj−1(dt/dx1) y,2(x1, x2) =
4 i=1
4 j=1 (j − 1)cij t i−1uj−2(du/dx2) y,12(x1, x2) =
4 i=1
4 j=1 (i − 1)(j − 1)cij t i−2uj−2(dt/dx1)(du/dx2) (3.6.6) where t and u are again given by equation (3.6.4). void bcucof(float y[], float y1[], float y2[], float y12[], float d1, float d2, float **c) Given arrays y[1..4], y1[1..4], y2[1..4], and y12[1..4], containing the function, gra￾dients, and cross derivative at the four grid points of a rectangular grid cell (numbered coun￾terclockwise from the lower left), and given d1 and d2, the length of the grid cell in the 1- and 2-directions, this routine returns the table c[1..4][1..4] that is used by routine bcuint for bicubic interpolation. { static int wt[16][16]= { 1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0, 0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0, -3,0,0,3,0,0,0,0,-2,0,0,-1,0,0,0,0, 2,0,0,-2,0,0,0,0,1,0,0,1,0,0,0,0, 0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0, 0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0, 0,0,0,0,-3,0,0,3,0,0,0,0,-2,0,0,-1, 0,0,0,0,2,0,0,-2,0,0,0,0,1,0,0,1, -3,3,0,0,-2,-1,0,0,0,0,0,0,0,0,0,0, 0,0,0,0,0,0,0,0,-3,3,0,0,-2,-1,0,0, 9,-9,9,-9,6,3,-3,-6,6,-6,-3,3,4,2,1,2, -6,6,-6,6,-4,-2,2,4,-3,3,3,-3,-2,-1,-1,-2, 2,-2,0,0,1,1,0,0,0,0,0,0,0,0,0,0, 0,0,0,0,0,0,0,0,2,-2,0,0,1,1,0,0, -6,6,-6,6,-3,-3,3,3,-4,4,2,-2,-2,-2,-1,-1, 4,-4,4,-4,2,2,-2,-2,2,-2,-2,2,1,1,1,1}; int l,k,j,i; float xx,d1d2,cl[16],x[16]; d1d2=d1*d2; for (i=1;i<=4;i++) { Pack a temporary vector x. x[i-1]=y[i]; x[i+3]=y1[i]*d1; x[i+7]=y2[i]*d2; x[i+11]=y12[i]*d1d2; } for (i=0;i<=15;i++) { Matrix multiply by the stored table. xx=0.0; for (k=0;k<=15;k++) xx += wt[i][k]*x[k]; cl[i]=xx; } l=0; for (i=1;i<=4;i++) Unpack the result into the output table. for (j=1;j<=4;j++) c[i][j]=cl[l++]; }