/**************************************************************
!*      Purpose: This program computes modified Struve        *
!*               function L0(x) using subroutine STVL0        *
!*      Input :  x   --- Argument of L0(x) ( x ò 0 )          *
!*      Output:  SL0 --- L0(x)                                *
!*      Example:                                              *
!*                  x        L0(x)                            *
!*              ------------------------                      *
!*                 0.0   .00000000D+00                        *
!*                 5.0   .27105917D+02                        *
!*                10.0   .28156522D+04                        *
!*                15.0   .33964933D+06                        *
!*                20.0   .43558283D+08                        *
!*                30.0   .78167230D+12                        *
!*                40.0   .14894775D+17                        *
!*                50.0   .29325538D+21                        *
!* ---------------------------------------------------------- *
!* REFERENCE: "Fortran Routines for Computation of Special    *
!*             Functions, jin.ece.uiuc.edu/routines/routines  *
!*             .html".                                        *
!*                                                            *
!*                          C++ Release By J-P Moreau, Paris. *
!*                                  (www.jpmoreau.fr)         *
!*************************************************************/  
#include <stdio.h>
#include <math.h>

double X, SL0;

double Sqr(double X) {
  return X*X;
}

void STVL0(double X, double *SL0) {
/*      ================================================
!       Purpose: Compute modified Struve function L0(x)
!       Input :  x   --- Argument of L0(x) ( x ò 0 )
!       Output:  SL0 --- L0(x)
!       ================================================ */
        double A0,A1,BI0,PI,R,S;
        int K, KM;
	    
	PI=3.141592653589793;
        S=1.0;
        R=1.0;
        if (X <= 20.0) {
           A0=2.0*X/PI;
           for (K=1; K<61; K++) {
              R=R*Sqr(X/(2.0*K+1.0));
              S=S+R;
              if (fabs(R/S) < 1.0e-12) goto e15;
           }
e15:       *SL0=A0*S;
        }
        else {
           KM=int(0.5*(X+1.0));
           if (X >= 50.0) KM=25;
           for (K=1; K<=KM; K++) {
              R=R*Sqr((2.0*K-1.0)/X);
              S=S+R;
              if (fabs(R/S) < 1.0e-12) goto e25;
           }
e25:       A1=exp(X)/sqrt(2.0*PI*X);
           R=1.0;
           BI0=1.0;
           for (K=1; K<17; K++) {
              R=0.125*R*Sqr(2.0*K-1.0)/(K*X);
              BI0=BI0+R;
              if (fabs(R/BI0) < 1.0e-12) goto e35;
           }
e35:       BI0=A1*BI0;
           *SL0=-2.0/(PI*X)*S+BI0;
        }
}


void main() {

        printf("\n Please enter x: ");
        scanf("%lf", &X);
        printf("\n   x        L0(x)    \n");
        printf("-----------------------\n");
        
	STVL0(X, &SL0);
       
	printf(" %5.1f  %e\n\n", X, SL0);
}

// end of file mstvl0.cpp