heat_solve_main.f90

Go to the documentation of this file.

!-------------------------------------------------------------------------------
! Copyright (c) 2019 FrontISTR Commons
! This software is released under the MIT License, see LICENSE.txt
!-------------------------------------------------------------------------------
module m_heat_solve_main
contains
 
  subroutine heat_solve_main(hecMESH, hecMAT, hecMESHmpc, hecMATmpc, &
     & fstrSOLID, fstrPARAM, fstrHEAT, ISTEP, iterALL, next_time, delta_time)
    use m_fstr
    use m_heat_mat_ass_conductivity
    use m_heat_mat_ass_capacity
    use m_heat_mat_ass_boundary
    use m_solve_lineq
    use m_fstr_elemact
    implicit none
    integer(kind=kint) :: i, iterALL, ISTEP, bup_n_dof
    real(kind=kreal)   :: delta_time, next_time
    type(hecmwst_local_mesh)  :: hecmesh
    type(hecmwst_matrix)      :: hecMAT
    type(fstr_solid)          :: fstrSOLID
    type(fstr_heat)           :: fstrHEAT
    type(fstr_param)          :: fstrPARAM
    type(hecmwst_local_mesh), pointer :: hecMESHmpc
    type(hecmwst_matrix), pointer :: hecMATmpc
    logical :: is_congerged
 
    if( fstrheat%elemact%ELEMACT_egrp_tot>0 ) &
      &  call fstr_update_elemact_heat( hecmesh, fstrheat%elemact, next_time, fstrsolid%elements )
 
    iterall = 0
    do
      iterall = iterall + 1
      hecmat%X = 0.0d0
 
      call heat_mat_ass_conductivity(hecmesh, hecmat, fstrsolid, fstrheat, fstrheat%beta)
      if(fstrheat%is_steady == 0) call heat_mat_ass_capacity(hecmesh, hecmat, fstrsolid, fstrheat, delta_time)
      call heat_mat_ass_boundary(hecmesh, hecmat, hecmeshmpc, hecmatmpc, &
       & fstrsolid, fstrheat, next_time, delta_time)
 
      hecmatmpc%Iarray(97) = 1 !Need numerical factorization
      bup_n_dof = hecmesh%n_dof
      hecmesh%n_dof = 1
      call solve_lineq(hecmeshmpc, hecmatmpc)
      hecmesh%n_dof = bup_n_dof
      call hecmw_mpc_tback_sol(hecmesh, hecmat, hecmatmpc)
 
      do i = 1, hecmesh%n_node
        fstrheat%TEMPC(i) = fstrheat%TEMP(i)
        fstrheat%TEMP (i) = hecmat%X(i)
      enddo
 
      if( fstrheat%elemact%ELEMACT_egrp_tot>0 ) &
        &  call fstr_updatedof_elemact( 1, hecmesh, fstrheat%elemact, fstrsolid%elements, fstrheat%TEMP0, fstrheat%TEMP )
 
      call heat_check_convergence(hecmesh, fstrheat, fstrparam, istep, iterall, is_congerged)
 
      if(is_congerged .or. fstrheat%is_iter_max_limit .or. fstrheat%beta == 0.0d0) exit
    enddo
  end subroutine heat_solve_main
 
  subroutine heat_check_convergence(hecMESH, fstrHEAT, fstrPARAM, ISTEP, iterALL, is_congerged)
    use m_fstr
    implicit none
    integer(kind=kint) :: i, iterALL, iter_max, ISTEP
    real(kind=kreal)   :: val, iter_eps
    type(hecmwst_local_mesh)  :: hecmesh
    type(fstr_heat)           :: fstrHEAT
    type(fstr_param)          :: fstrPARAM
    logical :: is_congerged
 
    is_congerged = .false.
    fstrheat%is_iter_max_limit = .false.
    iter_max = fstrparam%itmax(istep)
    iter_eps = fstrparam%eps(istep)
 
    val = 0.0d0
    do i = 1, hecmesh%nn_internal
      val = val + (fstrheat%TEMP(i) - fstrheat%TEMPC(i))**2
    enddo
    call hecmw_allreduce_r1(hecmesh, val, hecmw_sum)
    val = dsqrt(val)
 
    if(hecmesh%my_rank == 0)then
      write(*,'(i8,1pe16.6)') iterall, val
    endif
 
    if(val < iter_eps)then
      if(hecmesh%my_rank == 0) write(imsg,*) ' !!! CONVERGENCE ACHIEVED '
      is_congerged = .true.
    endif
 
    if(iter_max <= iterall)then
      if(hecmesh%my_rank == 0) write(*,*) ' !!! ITERATION COUNT OVER : MAX = ', iter_max
      fstrheat%is_iter_max_limit = .true.
    endif
  end subroutine heat_check_convergence
end module m_heat_solve_main