************80
! CV0 computes the initial characteristic value of Mathieu functions.
Discussion:
This procedure computes the initial characteristic value of Mathieu
functions for m <= 12 or q <= 300 or q <= m*m.
Licensing:
This routine is copyrighted by Shanjie Zhang and Jianming Jin. However,
they give permission to incorporate this routine into a user program
provided that the copyright is acknowledged.
Modified:
03 August 2012
Author:
Shanjie Zhang, Jianming Jin
Reference:
Shanjie Zhang, Jianming Jin,
Computation of Special Functions,
Wiley, 1996,
ISBN: 0-471-11963-6,
LC: QA351.C45.
Parameters:
Input, integer ( kind = 4 ) KD, the case code:
1, for cem(x,q) ( m = 0,2,4,...)
2, for cem(x,q) ( m = 1,3,5,...)
3, for sem(x,q) ( m = 1,3,5,...)
4, for sem(x,q) ( m = 2,4,6,...)
Input, integer ( kind = 4 ) M, the order of the functions.
Input, real ( kind = 8 ) Q, the parameter of the functions.
Output, real ( kind = 8 ) A0, the characteristic value.
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| integer(kind=4) | :: | kd | ||||
| integer(kind=4) | :: | m | ||||
| real(kind=8) | :: | q | ||||
| real(kind=8) | :: | a0 |
subroutine cv0 ( kd, m, q, a0 ) !*****************************************************************************80 ! !! CV0 computes the initial characteristic value of Mathieu functions. ! ! Discussion: ! ! This procedure computes the initial characteristic value of Mathieu ! functions for m <= 12 or q <= 300 or q <= m*m. ! ! Licensing: ! ! This routine is copyrighted by Shanjie Zhang and Jianming Jin. However, ! they give permission to incorporate this routine into a user program ! provided that the copyright is acknowledged. ! ! Modified: ! ! 03 August 2012 ! ! Author: ! ! Shanjie Zhang, Jianming Jin ! ! Reference: ! ! Shanjie Zhang, Jianming Jin, ! Computation of Special Functions, ! Wiley, 1996, ! ISBN: 0-471-11963-6, ! LC: QA351.C45. ! ! Parameters: ! ! Input, integer ( kind = 4 ) KD, the case code: ! 1, for cem(x,q) ( m = 0,2,4,...) ! 2, for cem(x,q) ( m = 1,3,5,...) ! 3, for sem(x,q) ( m = 1,3,5,...) ! 4, for sem(x,q) ( m = 2,4,6,...) ! ! Input, integer ( kind = 4 ) M, the order of the functions. ! ! Input, real ( kind = 8 ) Q, the parameter of the functions. ! ! Output, real ( kind = 8 ) A0, the characteristic value. ! implicit none real ( kind = 8 ) a0 integer ( kind = 4 ) kd integer ( kind = 4 ) m real ( kind = 8 ) q real ( kind = 8 ) q2 q2 = q * q if ( m == 0 ) then if ( q <= 1.0D+00 ) then a0 = ((( & 0.0036392D+00 * q2 & - 0.0125868D+00 ) * q2 & + 0.0546875D+00 ) * q2 & - 0.5D+00 ) * q2 else if ( q <= 10.0D+00 ) then a0 = (( & 3.999267D-03 * q & - 9.638957D-02 ) * q & - 0.88297D+00 ) * q & + 0.5542818D+00 else call cvql ( kd, m, q, a0 ) end if else if ( m == 1 ) then if ( q <= 1.0D+00 .and. kd == 2 ) then a0 = ((( & - 6.51D-04 * q & - 0.015625D+00 ) * q & - 0.125D+00 ) * q & + 1.0D+00 ) * q & + 1.0D+00 else if ( q <= 1.0D+00 .and. kd == 3 ) then a0 = ((( & - 6.51D-04 * q & + 0.015625D+00 ) * q & - 0.125D+00 ) * q & - 1.0D+00 ) * q & + 1.0D+00 else if ( q <= 10.0D+00 .and. kd == 2 ) then a0 = ((( & - 4.94603D-04 * q & + 1.92917D-02 ) * q & - 0.3089229D+00 ) * q & + 1.33372D+00 ) * q & + 0.811752D+00 else if ( q <= 10.0D+00 .and. kd == 3 ) then a0 = (( & 1.971096D-03 * q & - 5.482465D-02 ) * q & - 1.152218D+00 ) * q & + 1.10427D+00 else call cvql ( kd, m, q, a0 ) end if else if ( m == 2 ) then if ( q <= 1.0D+00 .and. kd == 1 ) then a0 = ((( & - 0.0036391D+00 * q2 & + 0.0125888D+00 ) * q2 & - 0.0551939D+00 ) * q2 & + 0.416667D+00 ) * q2 + 4.0D+00 else if ( q <= 1.0D+00 .and. kd == 4 ) then a0 = ( & 0.0003617D+00 * q2 & - 0.0833333D+00 ) * q2 + 4.0D+00 else if ( q <= 15.0D+00 .and. kd == 1 ) then a0 = ((( & 3.200972D-04 * q & - 8.667445D-03 ) * q & - 1.829032D-04 ) * q & + 0.9919999D+00 ) * q & + 3.3290504D+00 else if ( q <= 10.0D+00 .and. kd == 4 ) then a0 = (( & 2.38446D-03 * q & - 0.08725329D+00 ) * q & - 4.732542D-03 ) * q & + 4.00909D+00 else call cvql ( kd, m, q, a0 ) end if else if ( m == 3 ) then if ( q <= 1.0D+00 .and. kd == 2 ) then a0 = (( & 6.348D-04 * q & + 0.015625D+00 ) * q & + 0.0625 ) * q2 & + 9.0D+00 else if ( q <= 1.0D+00 .and. kd == 3 ) then a0 = (( & 6.348D-04 * q & - 0.015625D+00 ) * q & + 0.0625D+00 ) * q2 & + 9.0D+00 else if ( q <= 20.0D+00 .and. kd == 2 ) then a0 = ((( & 3.035731D-04 * q & - 1.453021D-02 ) * q & + 0.19069602D+00 ) * q & - 0.1039356D+00 ) * q & + 8.9449274D+00 else if ( q <= 15.0D+00 .and. kd == 3 ) then a0 = (( & 9.369364D-05 * q & - 0.03569325D+00 ) * q & + 0.2689874D+00 ) * q & + 8.771735D+00 else call cvql ( kd, m, q, a0 ) end if else if ( m == 4 ) then if ( q <= 1.0D+00 .and. kd == 1 ) then a0 = (( & - 2.1D-06 * q2 & + 5.012D-04 ) * q2 & + 0.0333333 ) * q2 & + 16.0D+00 else if ( q <= 1.0D+00 .and. kd == 4 ) then a0 = (( & 3.7D-06 * q2 & - 3.669D-04 ) * q2 & + 0.0333333D+00 ) * q2 & + 16.0D+00 else if ( q <= 25.0D+00 .and. kd == 1 ) then a0 = ((( & 1.076676D-04 * q & - 7.9684875D-03 ) * q & + 0.17344854D+00 ) * q & - 0.5924058D+00 ) * q & + 16.620847D+00 else if ( q <= 20.0D+00 .and. kd == 4 ) then a0 = (( & - 7.08719D-04 * q & + 3.8216144D-03 ) * q & + 0.1907493D+00 ) * q & + 15.744D+00 else call cvql ( kd, m, q, a0 ) end if else if ( m == 5 ) then if ( q <= 1.0D+00 .and. kd == 2 ) then a0 = (( & 6.8D-6 * q & + 1.42D-05 ) * q2 & + 0.0208333D+00 ) * q2 & + 25.0D+00 else if ( q <= 1.0D+00 .and. kd == 3 ) then a0 = (( & - 6.8D-06 * q & + 1.42D-05 ) * q2 & + 0.0208333D+00 ) * q2 & + 25.0D+00 else if ( q <= 35.0D+00 .and. kd == 2 ) then a0 = ((( & 2.238231D-05 * q & - 2.983416D-03 ) * q & + 0.10706975D+00 ) * q & - 0.600205D+00 ) * q & + 25.93515D+00 else if ( q <= 25.0D+00 .and. kd == 3 ) then a0 = (( & - 7.425364D-04 * q & + 2.18225D-02 ) * q & + 4.16399D-02 ) * q & + 24.897D+00 else call cvql ( kd, m, q, a0 ) end if else if ( m == 6 ) then if ( q <= 1.0D+00 ) then a0 = ( 0.4D-06 * q2 + 0.0142857 ) * q2 + 36.0D+00 else if ( q <= 40.0D+00 .and. kd == 1 ) then a0 = ((( & - 1.66846D-05 * q & + 4.80263D-04 ) * q & + 2.53998D-02 ) * q & - 0.181233D+00 ) * q & + 36.423D+00 else if ( q <= 35.0D+00 .and. kd == 4 ) then a0 = (( & - 4.57146D-04 * q & + 2.16609D-02 ) * q & - 2.349616D-02 ) * q & + 35.99251D+00 else call cvql ( kd, m, q, a0 ) end if else if ( m == 7 ) then if ( q <= 10.0D+00 ) then call cvqm ( m, q, a0 ) else if ( q <= 50.0D+00 .and. kd == 2 ) then a0 = ((( & - 1.411114D-05 * q & + 9.730514D-04 ) * q & - 3.097887D-03 ) * q & + 3.533597D-02 ) * q & + 49.0547D+00 else if ( q <= 40.0D+00 .and. kd == 3 ) then a0 = (( & - 3.043872D-04 * q & + 2.05511D-02 ) * q & - 9.16292D-02 ) * q & + 49.19035D+00 else call cvql ( kd, m, q, a0 ) end if else if ( 8 <= m ) then if ( q <= 3.0D+00 * m ) then call cvqm ( m, q, a0 ) else if ( m * m .lt. q ) then call cvql ( kd, m, q, a0 ) else if ( m == 8 .and. kd == 1 ) then a0 = ((( & 8.634308D-06 * q & - 2.100289D-03 ) * q & + 0.169072D+00 ) * q & - 4.64336D+00 ) * q & + 109.4211D+00 else if ( m == 8 .and. kd == 4 ) then a0 = (( & - 6.7842D-05 * q & + 2.2057D-03 ) * q & + 0.48296D+00 ) * q & + 56.59D+00 else if ( m == 9 .and. kd == 2 ) then a0 = ((( & 2.906435D-06 * q & - 1.019893D-03 ) * q & + 0.1101965D+00 ) * q & - 3.821851D+00 ) * q & + 127.6098D+00 else if ( m == 9 .and. kd == 3 ) then a0 = (( & - 9.577289D-05 * q & + 0.01043839D+00 ) * q & + 0.06588934D+00 ) * q & + 78.0198D+00 else if ( m == 10 .and. kd == 1 ) then a0 = ((( & 5.44927D-07 * q & - 3.926119D-04 ) * q & + 0.0612099D+00 ) * q & - 2.600805D+00 ) * q & + 138.1923D+00 else if ( m == 10 .and. kd == 4 ) then a0 = (( & - 7.660143D-05 * q & + 0.01132506D+00 ) * q & - 0.09746023D+00 ) * q & + 99.29494D+00 else if ( m == 11 .and. kd == 2 ) then a0 = ((( & - 5.67615D-07 * q & + 7.152722D-06 ) * q & + 0.01920291D+00 ) * q & - 1.081583D+00 ) * q & + 140.88D+00 else if ( m == 11 .and. kd == 3 ) then a0 = (( & - 6.310551D-05 * q & + 0.0119247D+00 ) * q & - 0.2681195D+00 ) * q & + 123.667D+00 else if ( m == 12 .and. kd == 1 ) then a0 = ((( & - 2.38351D-07 * q & - 2.90139D-05 ) * q & + 0.02023088D+00 ) * q & - 1.289D+00 ) * q & + 171.2723D+00 else if ( m == 12 .and. kd == 4 ) then a0 = ((( & 3.08902D-07 * q & - 1.577869D-04 ) * q & + 0.0247911D+00 ) * q & - 1.05454D+00 ) * q & + 161.471D+00 end if end if return end subroutine cv0