module resultados

! objetivo: calcular solução numérica e 
!           apresentar gráficos e listagem dos resultados

!-------------------------------------------------

use coeficientes
use solvers

implicit none

contains

!-------------------------------------------------

  subroutine solucao_numerica

    real*8 hmedio
    character*20 :: nome ! nome do arquivo iterativo

    nome = trim(caso)//'.itn'

    open(8,file=trim(nome))

    write(8,4) 
    4 format(//,t3,'Iteração', t14,'T máximo V.C.', t38,'T mínimo V.C.', &
               t62,'Norma L1 dos resíduos')

    call cpu_time ( tcpu1 )

    do it = 1, itmax    ! ciclo iterativo externo

       call coef_tipo_5
       call coef_tipo_2
       call coef_tipo_1
       call coef_tipo_3
       call coef_tipo_4
       call coef_tipo_7
       call coef_tipo_8
       call coef_tipo_6
       call coef_tipo_9

       if ( problema == 1 ) then
          call coef_tipo_16
          call coef_tipo_15
          call coef_tipo_17
       end if

       if ( problema == 2 ) then
          call coef_tipo_10
          call coef_tipo_11
          call coef_tipo_12
          call coef_tipo_13
          call coef_tipo_14
       end if

       call adi (m,n,nitm,aw,ae,as,an,ap,bp,T)

       norma = 0.0d0
       call norma_L1 (m,n,aw,ae,as,an,ap,bp,T,norma)
       if ( it == 1 ) norma1 = norma

	   if(it==1.or.it==itmax.or.mod(it,freq)==0) &
          write(8,2) it, maxval(T), minval(T), norma/norma1
          2 format(i8,2x,3(1pe24.15))

    end do

    close(8)

    call cpu_time ( tcpu2 )

    tcpu = tcpu2-tcpu1

    deallocate ( aw, ae, as, an, ap, bp )

    ! escreve campo de temperaturas

    if ( escreve == 1 ) then
       write(10,11)
       11 format(//,3x,'Campo de Temperaturas', &
                 //,5x,'Em cada coluna, do norte para o sul')
       do i = 1, m1
          write(10,12) i
          do j = n1, 1, -1
             write(10,13) T(i,j)
          end do
       end do
       do i = m1+1, m
          write(10,12) i
          do j = n, 1, -1
             write(10,13) T(i,j)
          end do
       end do
    end if

    12 format(/,'i = ',i5)
    13 format(1pe15.6)

    write(10,19) minval(T), maxval(T)
    19 format(/,5x,'Tmin: temperatura mínima da parede no centro dos volumes de controle (K) = ',1pe24.15,//, &
                5x,'Tmax: temperatura máxima da parede no centro dos volumes de controle (K) = ',1pe24.15)

    ! calcula e escreve as temperaturas dos contornos

    call face_tipo_4
    call face_tipo_6

    if ( problema == 1 ) then
       call face_tipo_16
       call face_tipo_17
    end if

    write(10,21)
    21 format(/,3x,'Ts: temperatura do contorno sul (do oeste para leste)')
    do i = 1, m
       write(10,13) Ts(i)
    end do

    write(10,22)
    22 format(/,3x,'Twi: temperatura do contorno oeste interna (do norte para o sul)')
    do j = n, n1+1, -1
       write(10,13) Tw(j)
    end do

    write(10,32)
    32 format(/,3x,'Twe: temperatura do contorno oeste externa (do norte para o sul)')
    do j = n1, 1, -1
       write(10,13) Tw(j)
    end do

    write(10,23)
    23 format(/,3x,'Te: temperatura do contorno leste (do norte para o sul)')
    do j = n, 1, -1
       write(10,13) Te(j)
    end do

    write(10,26)
    26 format(/,3x,'Tni: temperatura do contorno norte interno (do oeste para leste)')
    do i = 1, m1
       write(10,13) Tni(i)
    end do

    write(10,24)
    24 format(/,3x,'Tne: temperatura do contorno norte externo (do oeste para leste)')
    if ( problema == 1 ) then
       do i = m1+1, m
          write(10,13) Tn(i)
       end do
    else
       do i = 1, m
          write(10,13) Tn(i)
       end do
       write(10,25)
       25 format(/,3x,'Ts: temperatura do contorno sul interno (do oeste para leste)')
       do i = 1, m1
          write(10,13) Tsi(i)
       end do
    end if

    ! calcula a solução 1D

    call modelo_1D

    ! calcula a taxa de transferência de calor 2D

    q = 0.0d0
    j = 1
    do i = 1, m
       q = q + hg*(Taw-Ts(i))*dx(i,j)
    end do

    Tmax2D = maxval(Ts)

    write(10,20) Tmax2D, Tmax2D/Tmax1D, q, q/q1d
    20 format(/,5x,'Tmax-2D: Temperatura máxima da parede (K) = ',1pe24.15, &
             //,5x,'razão Tmax-2D/1D = ',1pe24.15, &
             //,5x,'q2D: taxa de transferência de calor (W) = ',1pe24.15, &
             //,5x,'razão q2D / q1D = ',1pe24.15)

    note_caso = 'notepad '//trim(nome)
    dos = system(note_caso) ! lista arquivo iterativo

    write(10,31) tcpu
    31 format(//,f14.3,' = tempo de CPU (s)')

    open(12,file=trim(adjustl(caso))//'.Richardson_3p0')

    hmedio = dsqrt ( (a+d)*(b+e) / ( m*n ) )

    write(12,55) hmedio,  Tmax2D/Tmax1D, q/q1d
 
    55 format(  2x, 1pe25.15, ' = h médio (m)',/, &
             ' 1', 1pe25.15, ' = razão Tmax-2D/1D',/, &
             ' 2', 1pe25.15, ' = razão q2D / q1D')

    close(12)

    ! lista arquivo de resultados para o programa Richardson
    dos = system('notepad '// trim(adjustl(caso)) // '.Richardson_3p0')

  end subroutine solucao_numerica

!-------------------------------------------------

  subroutine modelo_1D 

    real*8 :: Twh, Twc, eta, U,  P, At, &
              Rg,  Rw,  Rc,  RT, G, L

    L = a + d

    Rg = 1.0d0 / (hg * L)

    Rw = e / (k1 * L)

    P = 2.0d0

    At = 2.0d0 * d

    U = dsqrt(hc*P/(k1*At))

    eta = dtanh(U*b)/(U*b)

    G = a + eta*b

	Rc = 1.0d0 / (hc * G)

    RT = Rg + Rw + Rc

    q1d = (Taw - Tc) / RT

    Twh = Taw - q1d * Rg

    Twc = Twh - q1d * Rw

    Tmax1D = Twh

    write(10,1) Twc, Twh, q1d
    1 format(/,3x,'MODELO UNIDIMENSIONAL', &
             /,5x,'Twc: temperatura da parede fria   (K) = ',1pe24.15, &
             /,5x,'Twh: temperatura da parede quente (K) = ',1pe24.15, &
             /,5x,'q1D: taxa de transferência de calor (W) = ',1pe24.15, &
             /,7x,'Em q1D considera-se profundidade unitária em z')

  end subroutine modelo_1D

! -----------------------------------------------

end module resultados
