1 | include(dom.inc)
   2 | 
   3 |       SUBROUTINE BAND_INTEG(wkabs, kabs, kscat, Lb, Lo,                 &
   4 |      &                      WVNB_SI, DWVNB_SI,                          &
   5 |      &                      Gtot, Lbtot,  Htot, Q_rtot)
   6 | 
   7 |         USE mod_pmm
   8 |         USE mod_slave
   9 |         USE mod_inout
  10 | 
  11 |         IMPLICIT NONE
  12 | 
  13 |         include 'dom_constants.h'
  14 |         include 'pmm_constants.h'
  15 | 
  16 | !       IN
  17 |         DOM_REAL          :: Gmax, Gmin, Smax, Smin
  18 |         DOM_REAL          :: WVNB_SI, DWVNB_SI
  19 |         DOM_REAL,DIMENSION(is_nkabs)          :: wkabs
  20 |         DOM_REAL,DIMENSION(is_nkabs,is_ncells):: kabs
  21 | 
  22 |         DOM_REAL,DIMENSION(is_ncells)             :: Lb, kscat
  23 |         DOM_REAL,ALLOCATABLE, DIMENSION(:)        :: WSGG_W_cell
  24 |         DOM_REAL,ALLOCATABLE, DIMENSION(:)        :: WSGG_W_face
  25 |         DOM_REAL,DIMENSION(is_nbfaces)            :: Lo
  26 | 
  27 | !       OUT
  28 |         DOM_REAL,DIMENSION(3,is_ncells)          :: Q_rtot
  29 |         DOM_REAL,DIMENSION(is_ncells)            :: Gtot, Lbtot
  30 |         DOM_REAL,DIMENSION(is_nbfaces)           :: Htot
  31 | 
  32 | !       LOCAL
  33 |         DOM_INT           :: iquad, n_iter, i, j, maxiter, ibnd
  34 |         DOM_INT           :: mloc(1)
  35 |         DOM_REAL          :: error
  36 | 
  37 |         DOM_REAL,DIMENSION(is_ncells)              :: Gi, G
  38 |         DOM_REAL,DIMENSION(is_nfacesmax,is_ncells) :: Li
  39 |         DOM_REAL,DIMENSION(is_nbfaces)             :: H
  40 |         DOM_REAL,DIMENSION(3,is_ncells)            :: Qr
  41 | 
  42 | !       print*, " (",pmm_rank,") Doing band integration"
  43 | 
  44 | !       --------------------------!
  45 | !       Integrate the narrow band !
  46 | !       --------------------------!
  47 | 
  48 |         DO iquad=1,is_nkabs
  49 | 
  50 | !         --------------------------------------------------!
  51 | !         WSGG: Le traitement de chaque gaz se fait sur une !
  52 | !               fraction du corps noir A*Lb (avec sum(A)=1) !
  53 | !         --------------------------------------------------!
  54 |           IF (trim(mediumtype).eq.'WSGG') THEN
  55 | 
  56 |             IF (.not.ALLOCATED(WSGG_W_cell))                            &
  57 |      &          ALLOCATE(WSGG_W_cell(is_ncells))
  58 |             IF (.not.ALLOCATED(WSGG_W_face))                            &
  59 |      &          ALLOCATE(WSGG_W_face(is_nbfaces))
  60 | 
  61 |             CALL GATHER(s_WSGG_W(iquad,:), WSGG_W_cell, 1)
  62 |             Lb(:) = WSGG_W_cell(:) * Lb(:)
  63 | 
  64 |             CALL GATHER_FACES(s_WSGG_W(iquad,:), WSGG_W_face, 1)
  65 |             Lo(:) = WSGG_W_face(:) * Lo(:)
  66 | 
  67 |           ENDIF
  68 | 
  69 | !         -----------------------------------------------------!
  70 | !         Initialisation values before reflection subiteration !
  71 | !         -----------------------------------------------------!
  72 | 
  73 |           G  = 0.
  74 |           H  = 0.
  75 |           Li = 0.
  76 | 
  77 |           error   = 1.
  78 |           n_iter  = 0
  79 |           maxiter = 500
  80 | 
  81 |           DO WHILE ((error.gt.critconv).and.(n_iter.lt.maxiter))
  82 | 
  83 | !           -------------------------------------------!
  84 | !           Update values for reflection sub iteration !
  85 | !           -------------------------------------------!
  86 | 
  87 |             Gi = G
  88 | 
  89 | !           --------------------------------------------!
  90 | !           Update faces "in" lumminance for reflection !
  91 | !           --------------------------------------------!
  92 | 
  93 |             DO ibnd = 1,is_nbfaces
  94 | 
  95 |               j = is_bcell(ibnd)
  96 |               i = is_bface(ibnd)
  97 |               Li(i,j) = Lo(ibnd) + (1.-s_epsil(ibnd))*H(ibnd)/pi
  98 | 
  99 |             ENDDO
 100 | 
 101 | !           ---------------!
 102 | !           SPATIAL SCHEME !
 103 | !           ---------------!
 104 | 
 105 |             CALL SPATIAL_SCHEME(I_SCHEME, Lb, Gi, G, kabs(iquad,:),     &
 106 |      &                          kscat, Qr, Li, H, alpha)
 107 | 
 108 | !           ---------------------------------------!
 109 | !           Convergence test for reflexion problem !
 110 | !           ---------------------------------------!
 111 | 
 112 |             mloc  = MAXLOC(abs(G-Gi))
 113 |             IF ((n_iter.ne.0).and.((abs(Gi(mloc(1)))).ne.0.)) THEN
 114 |               error = (abs(G(mloc(1))-Gi(mloc(1))))/abs(Gi(mloc(1)))
 115 |             ELSE
 116 |               error = 1.
 117 |             ENDIF
 118 | 
 119 |             n_iter = n_iter + 1
 120 | 
 121 | !           IF (pmm_rank.eq.0) print*, "Nq=",iquad,"/",is_nkabs,"-",    &
 122 | !    &                         "n_iter=",n_iter,"- err=", error
 123 | 
 124 | !           IF (pmm_rank.eq.0) print*, "mloc=",mloc(1),"- G=",          &
 125 | !    &                         G(mloc(1)),"- Gi=",Gi(mloc(1))
 126 | 
 127 |           ENDDO
 128 | 
 129 | !         ----------------------------------------!
 130 | !         Sum over all slave directions and bands !
 131 | !         ----------------------------------------!
 132 | 
 133 |           G  = G  / pmm_n_p
 134 |           H  = H  / pmm_n_p
 135 |           Qr = Qr / pmm_n_p
 136 | 
 137 |           Gtot(:)  = Gtot(:)  + G(:)*kabs(iquad,:)*wkabs(iquad)*DWVNB_SI
 138 | 
 139 |           Lbtot(:) = Lbtot(:) + Lb(:)*4*pi*kabs(iquad,:)*wkabs(iquad)   &
 140 |      &               *DWVNB_SI
 141 | 
 142 |           Htot(:)  = Htot(:)  + H(:)*wkabs(iquad)*DWVNB_SI
 143 | 
 144 |           Q_rtot(1,:) = Q_rtot(1,:)+Qr(1,:)*wkabs(iquad)*DWVNB_SI
 145 |           Q_rtot(2,:) = Q_rtot(2,:)+Qr(2,:)*wkabs(iquad)*DWVNB_SI
 146 |           Q_rtot(3,:) = Q_rtot(3,:)+Qr(3,:)*wkabs(iquad)*DWVNB_SI
 147 | 
 148 |         ENDDO
 149 | 
 150 |       END SUBROUTINE BAND_INTEG