co2_sw.F90       coverage:  20.00 %func     1.36 %block


     1) module co2_sw_module
     2) 
     3)   ! module contains only 1 interface with other part, read as:
     4) 
     5)   ! co2_sw_prop(p,t,rho,dddt,dddp,fg,dfgdp,dfgdt,
     6)   !                 eng,ent,dhdt,dhdp,visc,dvdt,dvdp)
     7) 
     8)   use PFLOTRAN_Constants_module
     9) 
    10)       implicit none
    11) 
    12) #include "petsc/finclude/petscsys.h"
    13) 
    14)       private
    15)       save
    16)       PetscReal,  public :: co2_sw_c_t0_tab = -5.d0+273.15D0, co2_sw_c_p0_tab = 0.5d0
    17)       PetscReal,  public :: co2_sw_c_t1_tab = 5.d2+273.15D0, co2_sw_c_p1_tab = 250.d0
    18)       PetscReal,  public :: co2_sw_f_t0_tab = -5.d0+273.15D0, co2_sw_f_p0_tab = 0.5d0
    19)   !   PetscReal,  public :: co2_sw_f_t1_tab = 1.d2+273.15D0, co2_sw_f_p1_tab = 15.d0
    20)       PetscReal,  public :: co2_sw_f_t1_tab = 1.d2+273.15D0, co2_sw_f_p1_tab = 30.d0
    21) 
    22)       PetscInt :: ntab_t_c = 50, ntab_p_c = 100
    23)      ! PetscInt :: ntab_t_f = 105, ntab_p_f = 145
    24)       PetscInt :: ntab_t_f = 105, ntab_p_f = 295
    25)        
    26)       PetscReal :: dt_tab_c = 10.d0, dp_tab_c = 2.5d0
    27)       PetscReal :: dt_tab_f = 1.d0, dp_tab_f = 0.1d0
    28)   
    29)       PetscReal ,private, allocatable :: co2_prop_sw_c(:,:,:)
    30)       PetscReal ,private, allocatable :: co2_prop_sw_f(:,:,:)
    31)       PetscInt var_index(4)
    32)       logical ifinetable
    33) 
    34)       PetscReal, parameter, private ::  denc = 467.6d0, tc = 304.1282d0,&
    35)                                      rg = 0.1889241d0, pc = 7.3773d0
    36) 
    37)       public initialize_sw_interp, co2_sw_interp, init_span_wagner  
    38)  contains
    39) 
    40)      
    41) ! ************************************************************************** !
    42) 
    43) subroutine init_span_wagner(option)
    44)   ! 
    45)   ! init_span_wagner
    46)   ! 
    47)   ! Author: Chuan Lu
    48)   ! Date: 5/13/08
    49)   ! 
    50)   use Option_module
    51)   use co2_span_wagner_module
    52)   use co2_span_wagner_spline_module
    53) 
    54)   implicit none
    55)   type(option_type) :: option
    56)   PetscMPIInt :: myrank
    57) 
    58)   if (option%co2eos == EOS_SPAN_WAGNER) then
    59)     select case(option%itable)
    60)        case(0,1,2)
    61)          call initialize_span_wagner(option%itable, &
    62)                                      option%myrank, &
    63)                                      option)
    64)        case(4,5)
    65)          myrank = option%myrank
    66)          call initialize_span_wagner(ZERO_INTEGER,myrank, &
    67)                                      option)
    68)          call initialize_sw_interp(option%itable,myrank)
    69)        case(3)
    70)          call sw_spline_read
    71)        case default
    72)          print *, 'Wrong table option : STOP'
    73)          stop
    74)     end select
    75)   endif
    76) end subroutine init_span_wagner
    77) 
    78) ! ************************************************************************** !
    79) 
    80) subroutine initialize_sw_interp(itable,myrank)
    81)   ! 
    82)   ! prepare data table for interpolation
    83)   ! 
    84)       use co2_span_wagner_module, only: co2_span_wagner, vappr
    85)         
    86)         implicit none
    87)         PetscInt itable
    88)         PetscMPIInt myrank
    89)         
    90)         PetscInt i, j
    91)         PetscInt :: iflag = 1
    92)         
    93)         PetscReal p,t, ts,tmp2, tmp, dum1, dum2
    94) !       PetscReal dtmp,dddt,dddp
    95)         PetscReal dtemp,dpres
    96)         PetscReal rhodp,rhodt,fgdp,fgdt,engdp,engdt,entdp,entdt,vdp,vdt
    97) 
    98)         character*3 :: q
    99)         character*1 :: tab
   100)       
   101)       
   102)       
   103)       
   104) 
   105)         tab = char(9)
   106)          q = '","'
   107)         
   108)         allocate( co2_prop_sw_c(0:ntab_p_c,0:ntab_t_c,1:15)) 
   109)         allocate( co2_prop_sw_f(0:ntab_p_f,0:ntab_t_f,1:15))
   110)       
   111)         var_index(1) = 3; var_index(2) = 6; var_index(3) = 10; var_index(4) = 13;  
   112)         
   113)         
   114)         select case(itable)
   115)         case(5)
   116) ! prepare the coarse table  
   117)          tmp2=0D0
   118)          do j = 0, ntab_t_c
   119)             tmp=tmp2
   120)             t = co2_sw_c_t0_tab + dt_tab_c * real(j)
   121)            do i = 0, ntab_p_c
   122)              p = co2_sw_c_p0_tab + dp_tab_c * real(i)
   123)              co2_prop_sw_c(i, j, 1) = p
   124)              co2_prop_sw_c(i, j, 2) = t
   125)              co2_prop_sw_c(i,j,3) = tmp
   126)              call co2_span_wagner(p,t,co2_prop_sw_c(i,j,3),co2_prop_sw_c(i,j,4),&
   127)                 co2_prop_sw_c(i,j,5),co2_prop_sw_c(i,j,6),co2_prop_sw_c(i,j,7),&
   128)                 co2_prop_sw_c(i,j,8),co2_prop_sw_c(i,j,9),co2_prop_sw_c(i,j,10),&
   129)                 co2_prop_sw_c(i,j,11),co2_prop_sw_c(i,j,12),co2_prop_sw_c(i,j,13),&
   130)                 co2_prop_sw_c(i,j,14),co2_prop_sw_c(i,j,15),iflag)
   131)              
   132)                ts = 1000.D0
   133)                if (p .le.pc .and. t .le.tc) then
   134)                  call vappr (ts,p,dum1,dum2,TWELVE_INTEGER)
   135)                 endif
   136) 
   137)              
   138)              dpres = 1.e-3
   139)              if (t >ts) dpres = - dpres
   140)               call co2_span_wagner(p+dpres,t,rhodp,co2_prop_sw_c(i,j,4),&
   141)                       co2_prop_sw_c(i,j,5),fgdp,co2_prop_sw_c(i,j,7),&
   142)                       co2_prop_sw_c(i,j,8),engdp,entdp,&
   143)                       co2_prop_sw_c(i,j,11),co2_prop_sw_c(i,j,12),vdp,&
   144)                       co2_prop_sw_c(i,j,14),co2_prop_sw_c(i,j,15),iflag)
   145) 
   146)         
   147)                dtemp = 1.e-6
   148)                if (t  < ts) dtemp = -dtemp
   149)                call co2_span_wagner(p,t+dtemp,rhodt,co2_prop_sw_c(i,j,4),&
   150)                  co2_prop_sw_c(i,j,5),fgdt,co2_prop_sw_c(i,j,7),&
   151)                  co2_prop_sw_c(i,j,8),engdt,entdt,&
   152)                  co2_prop_sw_c(i,j,11),co2_prop_sw_c(i,j,12),vdt,&
   153)                  co2_prop_sw_c(i,j,14),co2_prop_sw_c(i,j,15),iflag)
   154) 
   155)                tmp = co2_prop_sw_c(i,j,3)
   156)                if (j==0) tmp2 = co2_prop_sw_c(i,j,3)
   157)       
   158)                 
   159)                    !dddt
   160)           co2_prop_sw_c(i,j,4) = (rhodt-co2_prop_sw_c(i,j,3))/dtemp
   161)           
   162)           !dfgdt
   163)           co2_prop_sw_c(i,j,8) = (fgdt-co2_prop_sw_c(i,j,6))/dtemp
   164)           
   165)           !dhdt
   166)           co2_prop_sw_c(i,j,11) = (entdt-co2_prop_sw_c(i,j,10))/dtemp
   167)           
   168)           !dvdt
   169)           co2_prop_sw_c(i,j,14) = (vdt-co2_prop_sw_c(i,j,13))/dtemp
   170) !       endif
   171) !       if (i>0) then
   172)           !dddp
   173)           co2_prop_sw_c(i,j,5) = (rhodp-co2_prop_sw_c(i,j,3))/dpres
   174)           
   175)           !dfgdp
   176)           co2_prop_sw_c(i,j,7) = (fgdp-co2_prop_sw_c(i,j,6))/dpres
   177)           
   178)           !dhdp
   179)           co2_prop_sw_c(i,j,12) = (entdp-co2_prop_sw_c(i,j,10))/dpres
   180)           
   181)           !dvdp
   182)           co2_prop_sw_c(i,j,15) = (vdp-co2_prop_sw_c(i,j,13))/dpres
   183)           
   184)             enddo
   185) !            print *,'co2_sw: ', i,j, p,t, co2_prop_sw_c(i,j,3)
   186)          enddo    
   187) 
   188) 
   189)  ! Prepare fine table
   190)          tmp2=0D0
   191)          do j = 0, ntab_t_f
   192)             tmp=tmp2
   193)             t = co2_sw_f_t0_tab + dt_tab_f * real(j)
   194)            do i = 0, ntab_p_f
   195)              p = co2_sw_f_p0_tab + dp_tab_f * real(i)
   196)              co2_prop_sw_f(i, j, 1) = p
   197)              co2_prop_sw_f(i, j, 2) = t
   198)              co2_prop_sw_f(i,j,3) = tmp
   199)              call co2_span_wagner(p,t,co2_prop_sw_f(i,j,3),co2_prop_sw_f(i,j,4),&
   200)                 co2_prop_sw_f(i,j,5),co2_prop_sw_f(i,j,6),co2_prop_sw_f(i,j,7),&
   201)                 co2_prop_sw_f(i,j,8),co2_prop_sw_f(i,j,9),co2_prop_sw_f(i,j,10),&
   202)                 co2_prop_sw_f(i,j,11),co2_prop_sw_f(i,j,12),co2_prop_sw_f(i,j,13),&
   203)                 co2_prop_sw_f(i,j,14),co2_prop_sw_f(i,j,15),iflag)
   204)              
   205)                ts = 1000.D0
   206)                if (p .le.pc .and. t .le.tc) then
   207)                  call vappr (ts,p,dum1,dum2,TWELVE_INTEGER)
   208)                 endif
   209) 
   210)              
   211)              dpres = 1.e-3
   212)              if (t >ts) dpres = - dpres
   213)               call co2_span_wagner(p+dpres,t,rhodp,co2_prop_sw_f(i,j,4),&
   214)                       co2_prop_sw_f(i,j,5),fgdp,co2_prop_sw_f(i,j,7),&
   215)                       co2_prop_sw_f(i,j,8),engdp,entdp,&
   216)                       co2_prop_sw_f(i,j,11),co2_prop_sw_f(i,j,12),vdp,&
   217)                       co2_prop_sw_f(i,j,14),co2_prop_sw_f(i,j,15),iflag)
   218) 
   219)         
   220)                dtemp = 1.e-6
   221)                if (t < ts) dtemp = -dtemp
   222)                call co2_span_wagner(p,t+dtemp,rhodt,co2_prop_sw_f(i,j,4),&
   223)                  co2_prop_sw_f(i,j,5),fgdt,co2_prop_sw_f(i,j,7),&
   224)                  co2_prop_sw_f(i,j,8),engdt,entdt,&
   225)                  co2_prop_sw_f(i,j,11),co2_prop_sw_f(i,j,12),vdt,&
   226)                  co2_prop_sw_f(i,j,14),co2_prop_sw_f(i,j,15),iflag)
   227) 
   228)                tmp = co2_prop_sw_f(i,j,3)
   229)                if (j==0) tmp2 = co2_prop_sw_f(i,j,3)
   230)       
   231)                 
   232)                    !dddt
   233)           co2_prop_sw_f(i,j,4) = (rhodt-co2_prop_sw_f(i,j,3))/dtemp
   234)           
   235)           !dfgdt
   236)           co2_prop_sw_f(i,j,8) = (fgdt-co2_prop_sw_f(i,j,6))/dtemp
   237)           
   238)           !dhdt
   239)           co2_prop_sw_f(i,j,11) = (entdt-co2_prop_sw_f(i,j,10))/dtemp
   240)           
   241)           !dvdt
   242)           co2_prop_sw_f(i,j,14) = (vdt-co2_prop_sw_f(i,j,13))/dtemp
   243) !       endif
   244) !       if (i>0) then
   245)           !dddp
   246)           co2_prop_sw_f(i,j,5) = (rhodp-co2_prop_sw_f(i,j,3))/dpres
   247)           
   248)           !dfgdp
   249)           co2_prop_sw_f(i,j,7) = (fgdp-co2_prop_sw_f(i,j,6))/dpres
   250)           
   251)           !dhdp
   252)           co2_prop_sw_f(i,j,12) = (entdp-co2_prop_sw_f(i,j,10))/dpres
   253)           
   254)           !dvdp
   255)           co2_prop_sw_f(i,j,15) = (vdp-co2_prop_sw_f(i,j,13))/dpres
   256)            
   257)            enddo
   258) !          print *,' co2_sw: ', i,j, p,t,co2_prop_sw_f(i,j,3)
   259)          enddo    
   260) 
   261) 
   262)       
   263)       if (myrank == 0) then
   264)         print *,'Writing Table lookup file: Coarse...'
   265)             if (myrank==0) print *,'--> open co2data_c.dat'
   266)            open(unit=122,file='co2data_c.dat',status='unknown')
   267)              write(122,'(''TITLE= "'',''co2data_c.dat'',''"'')')
   268)              write(122,'(''VARIABLES= "'',a6,100(a3,a6))') &
   269)           'p',q,'T',q,'d',q,'dddT',q,'dddp',q,'fg',q,'dfgdp',q,'dfgdT',q, &
   270)           'u',q,'h',q,'dhdT',q,'dhdp',q,'vis',q,'dvdT',q,'dvdp','"'
   271)              write(122,'(''ZONE T= "'',''",'','' I='',i4,'' , J='',i4)') ntab_t_c+1,ntab_p_c+1
   272)             do i = 0, ntab_p_c
   273)              
   274)              do j = 0, ntab_t_c
   275)                  write(122,'(1p15e14.6)') co2_prop_sw_c(i,j,1:15)
   276)              enddo
   277)           enddo
   278)         close (122)
   279)       
   280)          print *,'Writing Table lookup file: fine...'
   281)             if (myrank==0) print *,'--> open co2data_f.dat'
   282)            open(unit=122,file='co2data_f.dat',status='unknown')
   283)              write(122,'(''TITLE= "'',''co2data_c.dat'',''"'')')
   284)              write(122,'(''VARIABLES= "'',a6,100(a3,a6))') &
   285)           'p',q,'T',q,'d',q,'dddT',q,'dddp',q,'fg',q,'dfgdp',q,'dfgdT',q, &
   286)           'u',q,'h',q,'dhdT',q,'dhdp',q,'vis',q,'dvdT',q,'dvdp','"'
   287)              write(122,'(''ZONE T= "'',''",'','' I='',i4,'' , J='',i4)') ntab_t_f+1,ntab_p_f+1
   288)             do i = 0, ntab_p_f
   289)                 do j = 0, ntab_t_f
   290)                write(122,'(1p15e14.6)') co2_prop_sw_f(i,j,1:15)
   291)              enddo
   292)           enddo
   293)         close (122)
   294)      
   295)       endif
   296)  
   297)        
   298)       case(4)
   299)         if (myrank==0) print *,'Reading Table ...'
   300)         if (myrank==0) print *,'--> open co2data0.dat'
   301)         open(unit=122,file='co2data_c0.dat',status='unknown')
   302)         read(122,*)
   303)         read(122,*)
   304)         read(122,*)
   305)         do i = 0, ntab_p_c
   306)           do j = 0, ntab_t_c
   307)             read(122,'(1p15e14.6)') co2_prop_sw_c(i,j,1:15)
   308)           enddo
   309)         enddo
   310)         close (122)
   311) 
   312)         open(unit=122,file='co2data_f0.dat',status='unknown')
   313)         read(122,*)
   314)         read(122,*)
   315)         read(122,*)
   316)         do i = 0, ntab_p_f
   317)           do j = 0, ntab_t_f
   318)             read(122,'(1p15e14.6)') co2_prop_sw_f(i,j,1:15)
   319)           enddo
   320)         enddo
   321)         close (122)
   322)      end select
   323)       
   324) end subroutine initialize_sw_interp 
   325) 
   326) ! ************************************************************************** !
   327) 
   328) PetscReal function co2_prop_spwag(ip,it,iv)
   329)      implicit none 
   330)  !    PetscReal co2_prop_spwag
   331)      PetscInt ip,it,iv
   332)       
   333)       
   334)       if (ifinetable)then
   335)          co2_prop_spwag = co2_prop_sw_f(ip,it,iv) 
   336)        else
   337)           co2_prop_spwag = co2_prop_sw_c(ip,it,iv) 
   338)       endif 
   339)      
   340)   end function co2_prop_spwag
   341) 
   342) ! ************************************************************************** !
   343) 
   344) subroutine interp(x1,x2,y)
   345)   ! 
   346)   ! 2-d function interpolation
   347)   ! 
   348) 
   349)       use co2_span_wagner_module, only: vappr, co2_span_wagner
   350)       implicit none 
   351)   
   352) #include "petsc/finclude/petscsys.h"
   353)   
   354)       PetscReal :: x1,x2,y(15)   
   355)       
   356)       PetscReal factor(1:4), fac(2,2) , iindex, jindex
   357)       PetscReal funv(1:2,1:2,0:2), funi(1:2,1:2,0:2)
   358)              
   359)       PetscReal ps, tmp, tmp2, ntab_t, ntab_p, dt_tab, dp_tab, p0_tab, t0_tab 
   360)       PetscInt isucc, i1,i2,j1,j2, icross, i,j
   361)       PetscInt :: iflag = 1
   362)       
   363)       ifinetable = PETSC_FALSE
   364)       if (x2 <= co2_sw_f_t1_tab .and. x1<= co2_sw_f_p1_tab) then
   365)        ! within the fine grid table 
   366)         ifinetable = PETSC_TRUE
   367)       endif   
   368) 
   369) 
   370)     isucc= 0 
   371)     if (ifinetable)then
   372)       ntab_t = ntab_t_f
   373)       ntab_p = ntab_p_f
   374)       dt_tab = dt_tab_f 
   375)       dp_tab = dp_tab_f
   376)       p0_tab =  co2_sw_f_p0_tab
   377)       t0_tab =  co2_sw_f_t0_tab
   378) !     print *, 'using fine table', x1,x2
   379)     else
   380)       ntab_t = ntab_t_c
   381)       ntab_p = ntab_p_c
   382)       dt_tab = dt_tab_c 
   383)       dp_tab = dp_tab_c
   384)       p0_tab =  co2_sw_c_p0_tab
   385)       t0_tab =  co2_sw_c_t0_tab
   386)     endif
   387) 
   388) 
   389)     tmp = (x1 - p0_tab) / dp_tab; i1 = floor(tmp); i2 = i1+1; iindex=tmp 
   390)     tmp = (x2 - t0_tab) / dt_tab; j1 = floor(tmp); j2 = j1+1; jindex=tmp 
   391) 
   392)     isucc=1
   393) 
   394)  ! Check whether the table block covers the saturation line, missed special case.
   395)     icross =0
   396)     if (ifinetable) then
   397)       call vappr(co2_prop_spwag(i1,j1,TWO_INTEGER),ps,tmp,tmp2,ELEVEN_INTEGER)
   398)       if ((ps - co2_prop_spwag(i1,j1,ONE_INTEGER)) * (ps - co2_prop_spwag(i1,j2,ONE_INTEGER)) <0.D0)then
   399)         icross = 1; isucc=0
   400)       else
   401)         call vappr(co2_prop_spwag(i2,j1,TWO_INTEGER),ps,tmp,tmp,ELEVEN_INTEGER)
   402)         if ((ps - co2_prop_spwag(i2,j1,ONE_INTEGER)) * (ps - co2_prop_spwag(i2,j2,ONE_INTEGER)) <0.D0)then
   403)           icross = 1; isucc=0
   404)     !     else
   405)     !     call vappr(0.5D0 * (co2_prop_spwag(i1,j1,TWO_INTEGER)+co2_prop_spwag(i2,j1,TWO_INTEGER)),ps, tmp,tmp2,TWELVE_INTEGER)
   406)     !     if ((ts - co2_prop_spwag(i2,j1,TWO_INTEGER)) * (ts - co2_prop_spwag(i2,j2,TWO_INTEGER)) <0.D0)then
   407)     !       icross = 1; isucc=0
   408)     !     endif
   409)         endif
   410)       endif   
   411)     endif
   412) 
   413)     if (icross == 1) print *,'co2_sw: cross sat line'
   414) 
   415)     if (iindex > ntab_p .or. iindex < 0.d0 .or. jindex < 0.d0 .or. jindex > ntab_t) then
   416)       print  *,' Out of Table Bounds (interp): ', 'p=',x1,' t=',x2,' i=',iindex,' j=',jindex
   417)       isucc=0
   418)     endif
   419) 
   420) 
   421)     if (isucc>0 .and. icross ==0 )then
   422) 
   423)       factor(1)= (iindex-i2) * (jindex-j2)
   424)       factor(2)= -(iindex-i1) * (jindex-j2)
   425)       factor(3)= -(iindex-i2) * (jindex-j1)
   426)       factor(4)= (iindex-i1) * (jindex-j1)
   427) 
   428)       fac(1,1) = (iindex-i1 ); fac(2,1) = 1.D0 -  fac(1,1) 
   429)       fac(1,2) = (jindex-j1 ); fac(2,2) = 1.D0 -  fac(1,2)
   430)  
   431)        !  print *,  icross,isucc, factor, fac 
   432) 
   433)       i=1
   434)       y(1) = factor(1)*co2_prop_spwag(i1,j1,i) + factor(2)*co2_prop_spwag(i2,j1,i) &
   435)            + factor(3)*co2_prop_spwag(i1,j2,i) + factor(4)*co2_prop_spwag(i2,j2,i)
   436)       if (dabs(y(1)-x1)>1D-10 ) then
   437)         print *,' Error in intropolate::P',x1,iindex,factor;isucc=0
   438)       endif
   439)      !print *, 'Table: P ',iindex,jindex, factor,i  
   440)       i=i+1
   441)       y(2) = factor(1)*co2_prop_spwag(i1,j1,i) + factor(2)*co2_prop_spwag(i2,j1,i) &
   442)            + factor(3)*co2_prop_spwag(i1,j2,i) + factor(4)*co2_prop_spwag(i2,j2,i)
   443)       if (dabs(y(2)-x2)>1D-10 ) then
   444)         print *,' Error in intropolate:;T', x2,jindex,factor; isucc=0
   445)       endif
   446)     endif
   447) 
   448)   ! print *, 'Table: T',iindex,jindex,factor,i,isucc,itable
   449) 
   450) 
   451)     if (isucc==1)then
   452)        
   453)       do i =3,15
   454)       ! if (i==var_index(1) .or. i==var_index(2) .or.i==var_index(3) .or.i==var_index(4)) cycle
   455)         y(i) = factor(1)*co2_prop_spwag(i1,j1,i) + factor(2)*co2_prop_spwag(i2,j1,i) &
   456)            + factor(3)*co2_prop_spwag(i1,j2,i) + factor(4)*co2_prop_spwag(i2,j2,i)
   457)       enddo 
   458)        
   459) !#if 0      
   460)       do j=  1, 4
   461)        i= var_index(j)
   462)          funv(1,1,0) = co2_prop_spwag(i1,j1,i) 
   463)          funv(1,2,0) = co2_prop_spwag(i1,j2,i)
   464)          funv(2,1,0) = co2_prop_spwag(i2,j1,i)
   465)          funv(2,2,0) = co2_prop_spwag(i2,j2,i)
   466)           
   467)        select case(i)
   468)        case(3,10,13)     
   469)          funv(1,1,1) = co2_prop_spwag(i1,j1,i+2)   
   470)          funv(1,2,1) = co2_prop_spwag(i1,j2,i+2)
   471)          funv(2,1,1) = co2_prop_spwag(i2,j1,i+2)
   472)          funv(2,2,1) = co2_prop_spwag(i2,j2,i+2)
   473)          funv(1,1,2) = co2_prop_spwag(i1,j1,i+1)   
   474)          funv(1,2,2) = co2_prop_spwag(i1,j2,i+1)
   475)          funv(2,1,2) = co2_prop_spwag(i2,j1,i+1)
   476)          funv(2,2,2) = co2_prop_spwag(i2,j2,i+1)
   477)         case(6)
   478)          funv(1,1,1) = co2_prop_spwag(i1,j1,i+1)   
   479)          funv(1,2,1) = co2_prop_spwag(i1,j2,i+1)
   480)          funv(2,1,1) = co2_prop_spwag(i2,j1,i+1)
   481)          funv(2,2,1) = co2_prop_spwag(i2,j2,i+1)
   482)          funv(1,1,2) = co2_prop_spwag(i1,j1,i+2)   
   483)          funv(1,2,2) = co2_prop_spwag(i1,j2,i+2)
   484)          funv(2,1,2) = co2_prop_spwag(i2,j1,i+2)
   485)          funv(2,2,2) = co2_prop_spwag(i2,j2,i+2)
   486)        end select          
   487) 
   488)     
   489)         funi(:,:,0)= funv(:,:,0)   
   490)      
   491)         funi(1,1,1) = funv(1,1,1) - (funv(2,1,0) -funv(1,1,0)) /dp_tab
   492)         funi(1,1,2) = funv(1,1,2) - (funv(1,2,0) -funv(1,1,0)) /dt_tab
   493)         funi(2,1,1) = funv(2,1,1) - (funv(2,1,0) -funv(1,1,0)) /dp_tab
   494)         funi(2,1,2) = funv(2,1,2) - (funv(2,2,0) -funv(2,1,0)) /dt_tab
   495)         funi(1,2,1) = funv(1,2,1) - (funv(2,2,0) -funv(1,2,0)) /dp_tab
   496)         funi(1,2,2) = funv(1,2,2) - (funv(1,2,0) -funv(1,1,0)) /dt_tab
   497)         funi(2,2,1) = funv(2,2,1) - (funv(2,2,0) -funv(1,2,0)) /dp_tab
   498)         funi(2,2,2) = funv(2,2,2) - (funv(2,2,0) -funv(2,1,0)) /dt_tab
   499)   
   500)          y(i) = fac(2,2) *( fac(2,1)*(funi(1,1,0) + fac(1,1)* fac(2,1)* funi(1,1,1)+ fac(1,2)* fac(2,2)* funi(1,1,2))&    
   501)                            +fac(1,1)*(funi(2,1,0) - fac(1,1)* fac(2,1)* funi(2,1,1)- fac(1,2)* fac(2,2)* funi(2,1,2)))&  
   502)               + fac(1,2) *( fac(2,1)*(funi(1,2,0) + fac(1,1)* fac(2,1)* funi(1,2,1)+ fac(1,2)* fac(2,2)* funi(1,2,2))&    
   503)                            +fac(1,1)*(funi(2,2,0) - fac(1,1)* fac(2,1)* funi(2,2,1)- fac(1,2)* fac(2,2)* funi(2,2,2)))  
   504)         enddo
   505) !#endif
   506)        endif
   507)        
   508)       if ((icross ==1 .and. ifinetable).or. isucc < 1)then
   509)       ! print *, ' Exit table looking', icross, isucc, ifinetable
   510)        call co2_span_wagner(x1,x2,y(3),y(4),y(5),y(6),y(7),y(8), &
   511)         y(9),y(10),y(11),y(12),y(13),y(14),y(15),ZERO_INTEGER,iflag)
   512)        
   513)       endif  
   514) end subroutine 
   515) 
   516) ! ************************************************************************** !
   517) 
   518) subroutine co2_sw_interp(p,tc,rho,dddt,dddp,fg,dfgdp,dfgdt, &
   519)         eng,ent,dhdt,dhdp,visc,dvdt,dvdp,itable)
   520)          
   521)        implicit none
   522)        
   523)       PetscReal :: p,tc,rho,eng,ent,dhdt,dhdp,dddt,dddp,visc,dvdt,dvdp
   524)       PetscReal :: fg,dfgdp,dfgdt
   525)       PetscInt itable
   526)        
   527)       PetscReal :: t, prop(15)
   528)      
   529)   
   530)       t=tc+273.15D0
   531)       prop(1)=p
   532)       prop(2)=t
   533)        
   534)       call  interp(p,t,prop)       
   535)       
   536)       rho = prop(3)
   537)       dddt = prop(4)
   538)       dddp = prop(5)
   539)       fg = prop(6)
   540)       dfgdp = prop(7)
   541)       dfgdt = prop(8)
   542)       eng = prop(9)
   543)       ent = prop(10)
   544)       dhdt = prop(11)
   545)       dhdp = prop(12)
   546)       visc = prop(13)
   547)       dvdt = prop(14)
   548)       dvdp = prop(15)
   549)   
   550) !  contains
   551) end subroutine co2_sw_interp 
   552) 
   553)       
   554) end module co2_sw_module     

generated by
Intel(R) C++/Fortran Compiler code-coverage tool
Web-Page Owner: Nobody