00001 subroutine calc_Jome_int(typeR,omega,Jome,Jome_int)
00002 use phys_constant, only : long
00003 use def_matter_parameter
00004 implicit none
00005 real(long), intent(inout) :: omega, Jome
00006 real(long), intent(out):: Jome_int
00007
00008
00009 real(long) :: pmfac
00010 character(len=5) :: typeR
00011
00012
00013
00014
00015 select case (typeR)
00016 case('type0')
00017 if (index_DR.ne.2.0d0) then
00018 Jome_int = A2DR*omega**2*((ome/omega)**index_DR/(2.0d0-index_DR)-0.5d0) &
00019 & - index_DR*A2DR*ome**2/(4.0d0-2.0d0*index_DR)
00020 else
00021 Jome_int = - A2DR*ome**2*dlog(ome/omega) + 0.5d0*A2DR*(ome**2-omega**2)
00022 end if
00023
00024 case('type1','type2')
00025 if (typeR.eq.'type1') pmfac = -1.0d0
00026 if (typeR.eq.'type2') pmfac = 1.0d0
00027 Jome_int = A2DR*(omega-ome)*((1.0d0+index_DR)*omega+ome) &
00028 & *(pmfac*(omega/ome - 1.0d0))**index_DR &
00029 & /((1.0+index_DR)*(2.0+index_DR))
00030
00031 case('OJ2nd')
00032 Jome_int = &
00033 & A2DR*ome**2*((-0.7853981633974483*A2DR)/B2DR + &
00034 & (1.*A2DR**3*Jome*ome**3)/(Jome**4 + A2DR**4*ome**4) + &
00035 & (1.*A2DR**3*Jome**2*ome**2)/ &
00036 & (B2DR*Jome**4 + A2DR**4*B2DR*ome**4) + &
00037 & (0.3535533905932738 + (0.5*A2DR)/B2DR)* &
00038 & atan(1. - (1.4142135623730951*Jome)/(A2DR*ome)) + &
00039 & (-0.3535533905932738 + (0.5*A2DR)/B2DR)* &
00040 & atan(1. + (1.4142135623730951*Jome)/(A2DR*ome)) + &
00041 & 0.1767766952966369*Log(Jome**2 - &
00042 & 1.4142135623730951*A2DR*Jome*ome + A2DR**2*ome**2) - &
00043 & 0.1767766952966369*Log(Jome**2 + &
00044 & 1.4142135623730951*A2DR*Jome*ome + A2DR**2*ome**2))
00045
00046 case('OJ4th')
00047 Jome_int = &
00048 & (-0.1414213562373095*A2DR* &
00049 & (3.0575518754255118*A2DR**2 + 0.7217900873323894*B2DR**2)* &
00050 & ome**2)/B2DR**2 + 0.05*A2DR*ome**2* &
00051 & ((20.*A2DR**4*Jome**3*ome**2)/ &
00052 & (B2DR**2*(Jome**5 + A2DR**5*ome**5)) + &
00053 & (20.*A2DR**4*Jome*ome**4)/(Jome**5 + A2DR**5*ome**5) + &
00054 & (-7.608452130361228 + (4.702282018339785*A2DR**2)/B2DR**2)* &
00055 & atan((0.2628655560595668* &
00056 & (4.*Jome + 1.2360679774997898*A2DR*ome))/(A2DR*ome)) + &
00057 & (4.702282018339785 + (7.608452130361228*A2DR**2)/B2DR**2)* &
00058 & atan((0.42532540417601994* &
00059 & (-4.*Jome + 3.23606797749979*A2DR*ome))/(A2DR*ome)) - &
00060 & 4.*Log(Jome + A2DR*ome) - &
00061 & (4.*A2DR**2*Log(Jome + A2DR*ome))/B2DR**2 + &
00062 & 3.23606797749979*Log(Jome**2 - &
00063 & 1.618033988749895*A2DR*Jome*ome + A2DR**2*ome**2) - &
00064 & (1.2360679774997898*A2DR**2* &
00065 & Log(Jome**2 - 1.618033988749895*A2DR*Jome*ome + &
00066 & A2DR**2*ome**2))/B2DR**2 - &
00067 & 1.2360679774997898* &
00068 & Log(Jome**2 + 0.6180339887498949*A2DR*Jome*ome + &
00069 & A2DR**2*ome**2) + &
00070 & (3.23606797749979*A2DR**2* &
00071 & Log(Jome**2 + 0.6180339887498949*A2DR*Jome*ome + &
00072 & A2DR**2*ome**2))/B2DR**2)
00073
00074 case('OJjco')
00075 Jome = max(Jome,1.0d-30)
00076 Jome_int = &
00077 & (-0.5*Jome*(2.*(1. + index_DRp)**2*Jome**(1. + index_DRp) - &
00078 & 2.*A2DR*index_DRp*(2. + index_DRp)*Jome**index_DRp*ome + &
00079 & B2DR**index_DRp*(2. + 3.*index_DRp + index_DRp**2)*Jome* &
00080 & ome**index_DRp))/ &
00081 & (A2DR*B2DR**(1.*index_DRp)*(1. + index_DRp)*(2. + index_DRp)* &
00082 & ome**(1.*index_DRp))
00083 end select
00084
00085 end subroutine calc_Jome_int
