/*Translated by FOR_C, v3.4.2 (-), on 07/09/115 at 08:33:13 */
/*FOR_C Options SET: ftn=u io=c no=p op=aimnv pf=,p_dtgfi2 s=dbov str=l x=f - prototypes */
#include <math.h>
#include "fcrt.h"
#include <string.h>
#include <stddef.h>
#include <stdio.h>
#include <stdlib.h>
#include "p_dtgfi2.h"
#include "drdtgfi2.h"
/* File: DRDTGFI2.[F|FOR] CONTAINS DRDTGFI2 AND DTGFCN. */
/* PARAMETER translations */
#define MB 16
#define MEVAL 8
#define MP 28
#define MT 336
/* end of PARAMETER translations */
int main( )
{
char title[76];
long int _n, bdry[MB][4], i, ieval, info[3], ip[MP], kf, mode,
ncont, nt, triang[MT];
double dz[MP][2], dzout[2], dztemp[2], dztrue[2], q[2], savwrk[28],
w[MP], z[MP], zout, ztrue;
static double qtab[MEVAL][2]={0.0e0,0.9e0,-2.0e0,-1.1e0,-0.8e0,
-0.7e0,-2.7e0,1.5e0,0.0e0,0.9e0,-2.0e0,-0.5e0,-0.8e0,-0.7e0,-2.0e0,
0.7e0};
static double x[MP]={-0.76059e0,-0.02286e0,-0.44790e0,0.15068e0,
-0.87287e0,-0.23390e0,0.06093e0,-0.84142e0,-0.69173e0,-0.56613e0,
-0.42243e0,-0.17249e0,-0.06484e0,0.48286e0,-0.88784e0,0.10277e0,
0.69087e0,-0.84292e0,0.08784e0,-0.95068e0,0.02496e0,0.94973e0,
0.04588e0,-0.51667e0,-0.77561e0,0.90000e0,-0.70830e0,-0.40500e0};
static double y[MP]={-0.31421e0,0.75657e0,0.14321e0,0.42353e0,
-0.62983e0,0.12326e0,0.34054e0,-0.55144e0,0.34158e0,-0.67143e0,
0.43087e0,0.96081e0,-0.68130e0,0.12095e0,-0.10663e0,0.29219e0,
0.31028e0,0.12934e0,0.10709e0,-0.42307e0,0.49895e0,0.68597e0,
-0.78215e0,-0.12362e0,-0.88827e0,-0.60000e0,-0.88620e0,0.08600e0};
static long np = MP;
static LOGICAL32 wantdz = TRUE;
/* OFFSET Vectors w/subscript range: 1 to dimension */
double *const Dzout = &dzout[0] - 1;
double *const Dztemp = &dztemp[0] - 1;
double *const Dztrue = &dztrue[0] - 1;
long *const Info = &info[0] - 1;
long *const Ip = &ip[0] - 1;
double *const Q = &q[0] - 1;
double *const Savwrk = &savwrk[0] - 1;
long *const Triang = &triang[0] - 1;
double *const W = &w[0] - 1;
double *const X = &x[0] - 1;
double *const Y = &y[0] - 1;
double *const Z = &z[0] - 1;
/* end of OFFSET VECTORS */
/*>> 2001-07-16 DRDTGFI2 Krogh Added exponent 0 to some constants.
*>> 2001-05-22 DRDTGFI2 Krogh Minor change for making .f90 version.
*>> 1997-07-01 DRDTGFI2 Krogh Reversed subscripts in B (CLL suggestion)
*>> 1997-06-20 DRDTGFI2 Krogh Minor changes to get proper C conversion.
*>> 1997-06-18 DRDTGFI2 CLL
*>> 1996-03-04 DRDTGFI2 CLL
*>> 1996-02-02 DRDTGFI2 CLL
*>> 1995-10-31 DRDTGFI2 CLL
* Demo driver for DTGFI, DTGGRD, DTGPD, etc.
* Tests search when extrapolating.
* ------------------------------------------------------------------
*--D replaces "?": DR?TGFI2, ?TGFI, ?TGGRD, ?TGPD, ?TGFCN, ?TGPRG
* ------------------------------------------------------------------ */
/* ------------------------------------------------------------------ */
printf(" Program DRDTGFI2. Demo driver for DTGFI,\n"
" DTGGRD, DTGPD, etc.\n");
dtggrd( x, y, np, ip, w, triang, MT, bdry, MB, &nt, info );
dtgprg( x, y, np, triang, bdry, Info[2], nt );
if (Info[1] != 0)
{
printf(" Error return from DTGGRD. INFO(1) =%5ld\n", Info[1]);
exit(0);
}
/* do 200 ncont = 0,1 */
for (ncont = 1; ncont <= 1; ncont++)
{
/* do 100 KF= 0,4 */
for (kf = 1; kf <= 1; kf++)
{
printf("\n\n New Case. Interpolation with C%1ld continuity.\n", ncont);
for (i = 1; i <= np; i++)
{
dtgfcn( kf, X[i], Y[i], title,76, &Z[i], &Dztemp[1],
&Dztemp[2] );
}
printf("\n%s\n", title);
dtgpd( x, y, z, dz, np, triang, nt, ip );
Savwrk[1] = 0.0e0;
printf("\n Values and partial derivatives Interpolated\n "
" and Extrapolated at selected points.\n\n I X Y"
" Z_INTERP Z_TRUE Z_ERR\n "
" DZ1_INTERP DZ1_TRUE DZ1_ERR\n "
" DZ2_INTERP DZ2_TRUE DZ2_ERR\n\n");
for (ieval = 1; ieval <= MEVAL; ieval++)
{
Q[1] = qtab[ieval - 1][0];
Q[2] = qtab[ieval - 1][1];
dtgfi( x, y, z, dz, triang, nt, bdry, MB, ncont, q,
&zout, wantdz, dzout, &mode, savwrk );
if (mode >= 0)
{
dtgfcn( kf, Q[1], Q[2], title,76, &ztrue, &Dztrue[1],
&Dztrue[2] );
printf(" %4ld", ieval);
for(_n=0L; _n < sizeof(q)/sizeof(double); _n++)
printf("%6.2f", q[_n]);
printf("%6.2f%11.6f%11.6f\n", zout, ztrue, zout - ztrue);
printf(" %11.6f%11.6f%10.2e\n", Dzout[1], Dztrue[1], Dzout[1] - Dztrue[1]);
printf(" %11.6f%11.6f%10.2e\n", Dzout[2], Dztrue[2], Dzout[2] - Dztrue[2]);
}
else
{
printf(" %4ld", ieval);
for(_n=0L; _n < sizeof(q)/sizeof(double); _n++)
printf("%6.2f", q[_n]);
printf(" Error.\n");
}
}
}
}
exit(0);
} /* end of function */
/* ================================================================== */
#include <string.h>
void /*FUNCTION*/ dtgfcn(
long kf,
double x,
double y,
char *title,int title_s,
double *z,
double *zx,
double *zy)
{
static char title1[4-(0)+1][57], title2[4-(0)+1][20];
long int i, _i, _r;
static int _aini = 1;
#define NCHRTMPS 1
CHRTMP _c[NCHRTMPS];
ini_chrtmp(_c,NCHRTMPS);
if( _aini ){ /* Do 1 TIME INITIALIZATIONS! */
{ static char* _itmp0[] = {" CONSTANT FUNCTION Z = 2",
" "," LINEAR FUNCTION Z = ( 1 + 2*X + 3*Y ) / 6"," ",
" QUADRATIC FUNCTION Z = ( -1 + 2*X - 3*Y + 4*X**2 - ",
"X*Y + 9*Y**2 ) / 10"," CUBIC FUNCTION Z = ( 9*X**3 - 2*(X**2)*Y + 3*X*Y**2",
" - 4 * Y**3 ) / 10"," EXPONENTIAL FUNCTION Z = EXP( -2 * (X**2 + Y**2) )",
" "};
for (i = 0, _r = 0; i <= 4; i++)
{
f_strncpy( title1[i], _itmp0[_r++], 56 );
f_strncpy( title2[i], _itmp0[_r++], 19 );
}
}
_aini = 0;
}
/*>> 1995-09-26 CLL Editing for inclusion into MATH77.
* C.L.Lawson, JPL, 1976 Dec 10. Edited comments 1979 Mar 3.
* This subr evaluates a function and its first partial derivs as
* selected by KF. KF can be from 0 to 4.
* Input is KF, X, and Y. Output is TITLE, Z, ZX, and ZY.
* ------------------------------------------------------------------ */
/* ------------------------------------------------------------------ */
f_strncpy( title, f_concat(&_c[0],title1[kf],title2[kf],NULL), title_s -
1 );
if (kf == 0)
{
/* KF=0 CONSTANT FCN. */
*z = 2.0e0;
*zx = 0.0e0;
*zy = 0.0e0;
}
else if (kf == 1)
{
/* KF=1 LINEAR FCN. */
*z = (1.0e0 + 2.0e0*x + 3.0e0*y)/6.0e0;
*zx = 2.0e0/6.0e0;
*zy = 3.0e0/6.0e0;
}
else if (kf == 2)
{
/* KF=2 QUADRATIC FCN. */
*z = (-1.0e0 + 2.0e0*x - 3.0e0*y + 4.0e0*SQ(x) - x*y + 9.0e0*
SQ(y))*0.1e0;
*zx = (2.0e0 + 8.0e0*x - y)*0.1e0;
*zy = (-3.0e0 - x + 18.0e0*y)*0.1e0;
}
else if (kf == 3)
{
/* KF=3 CUBIC FCN. */
*z = (9.0e0*CUBE(x) - 2.0e0*(SQ(x))*y + 3.0e0*x*SQ(y) - 4.0e0*
CUBE(y))*0.1e0;
*zx = (27.0e0*SQ(x) - 4.0e0*x*y + 3.0e0*SQ(y))*0.1e0;
*zy = (-2.0e0*SQ(x) + 6.0e0*x*y - 12.0e0*SQ(y))*0.1e0;
}
else if (kf == 4)
{
/* KF=4 EXPONENTIAL FCN. */
*z = exp( -(SQ(x) + SQ(y))*2.0e0 );
/* NOTE THAT THE INFLECTION POINT OF THIS FCN IN ANY
* RADIAL DIRECTION FROM THE ORIGIN IS AT R = .5
* */
*zx = -4.0e0*x**z;
*zy = -4.0e0*y**z;
}
rel_chrtmp(_c,NCHRTMPS);
return;
#undef NCHRTMPS
} /* end of function */