/*Translated by FOR_C, v3.4.2 (-), on 07/09/115 at 08:33:08 */
/*FOR_C Options SET: ftn=u io=c no=p op=aimnv pf=,p_dbesjn s=dbov str=l x=f - prototypes */
#include <math.h>
#include "fcrt.h"
#include <stdio.h>
#include <stdlib.h>
#include "p_dbesjn.h"
/*     DRDBESJN
 *>> 1999-01-07 DRDBESJN Krogh Added external statement.
 *>> 1996-05-31 DRDBESJN Krogh  Changes to use M77CON
 *>> 1994-09-01 DRDBESJN WVS Moved formats to top for C conversion
 *>> 1992-04-29 DRDBESJN CAO Replaced '1' in format.
 *>> 1987-12-09 DRDBESJN Lawson  Initial Code.
 *--D replaces "?": DR?BESJN, ?BESJN, ?BESYN
 *     DEMONSTRATION PROGRAM FOR BESSEL function.
 * */
 
int main( )
{
	long int i, n;
	double bj[2], by[2], z;
	static double x[3]={0.5e0,1.5e0,3.2e0};
	static double alpha[3]={1.5e0,3.0e0,7.8e0};
	static double pi2 = 1.5707963267948966192313216e0;
		/* OFFSET Vectors w/subscript range: 1 to dimension */
	double *const Alpha = &alpha[0] - 1;
	double *const Bj = &bj[0] - 1;
	double *const By = &by[0] - 1;
	double *const X = &x[0] - 1;
		/* end of OFFSET VECTORS */
 
 
 
 
	printf("     X         NU           J(NU,X)           Y(NU,X)            Z\n");
	printf("                           J(NU+1,X)         Y(NU+1,X)\n \n");
 
	for (i = 1; i <= 3; i++)
	{
		n = 2;
		dbesjn( X[i], Alpha[i], n, bj );
		dbesyn( X[i], Alpha[i], n, by );
		z = pi2*X[i]*(Bj[2]*By[1] - Bj[1]*By[2]) - 1.e0;
		printf(" %6.2f     %6.2f    %15.8g     %15.8g%13.2g\n", X[i], Alpha[i], Bj[1], By[1],
		   z);
		printf("                      %15.8g     %15.8g     \n \n", Bj[2], By[2]);
	}
 
	exit(0);
} /* end of function */