!--------1---------2---------3---------4---------5---------6---------7---------8
subroutine radius_check
!--------1---------2---------3---------4---------5---------6---------7---------8
  use global_constants, only : rkind, time0, time, year
  use material_parameters, only : iw, ir &
    & , ii, ii_old, ir, ir_old, iw, iw_old, radius_ice, radius_water, mn, frez_start
  use global_parameters, only : halfdr, radius, temp
  implicit none
  radius_ice = radius(ii)+halfdr
  ii_old = ii
  if(ii /= 0 .AND. radius_water /= 0._rkind) then
    call water_check
    if(frez_start) then
      ii = iw+1
      radius_ice = radius(ii)+halfdr
    endif
    if(iw /= 0) then
      if(.NOT.frez_start) call rock_check
    endif
  endif
  return
end subroutine radius_check
!--------1---------2---------3---------4---------5---------6---------7---------8
subroutine water_check
!--------1---------2---------3---------4---------5---------6---------7---------8
  use global_constants, only : rkind, time, year
  use global_parameters, only : radius, halfdr, irmax
  use material_parameters, only : ii, radius_water, iw, iw_old, mn &
    & , radius_rock, drock, ir
  implicit none
  integer :: i
  do i = 0, irmax
    if(radius_water < radius(i)+halfdr) then
      iw = i - 1
      exit
    endif
  enddo
  return
end subroutine water_check
!--------1---------2---------3---------4---------5---------6---------7---------8
subroutine rock_check
!--------1---------2---------3---------4---------5---------6---------7---------8
    use global_constants, only : rkind, volume_ratio_rock, time, year, density &
      & , W_ALTERATION, mass_ratio_rock_new
    use global_parameters, only : radius, temp, dr, halfdr
    use heat_parameters, only : thermal_conductivity, thermometric_diffusivity 
    use material_parameters, only : radius_rock, radius_water, ir, iw, ir_old &
      & , drock, mn, melt, dr_rock, iw_old, material_select&
      & , mass_ratio_rock, mass_ratio_ice, mass_ratio_water &
      & , mass_ratio_anhydrous_rock, mass_ratio_hydrous_rock, dehyd &
      & , mass_ratio_vapor
    use ice_parameters, only : temperature_melting_I
    use energy, only : energy_rock
    implicit none
    integer :: i
    if(ir>=1) then
      radius_rock = radius_rock + drock/(radius_rock**2._rkind)/volume_ratio_rock(15)
    elseif(ir==0) then
      radius_rock = radius_water*(((mass_ratio_hydrous_rock(0)/density(4) &
        & +mass_ratio_anhydrous_rock(0)/density(4))*density(6)/volume_ratio_rock(15))**(1._rkind/3._rkind))
      if(radius_rock < halfdr) return
    endif
    drock = 0._rkind
    if(radius_rock < halfdr) return
    do i = 0, iw
      if(radius_rock < radius(i)+halfdr) then
        dr_rock = radius_rock - radius(i-1)
        volume_ratio_rock(5) = (dr_rock-halfdr)/dr
        ir = i - 1
        mn(i) = 5
        call material_select(i)
        exit
      endif
    enddo
    if(ir > ir_old) energy_rock(1) = energy_rock(0)
    ir = max(ir, ir_old)
    do i = 0, ir
      mn(i) = 15
      if(mass_ratio_vapor(i)/=0._rkind) then 
        mn(i) = 17
        if(mass_ratio_water(i)==0._rkind) then
          mn(i) = 16
        endif
     endif
      mass_ratio_rock(i) = mass_ratio_rock_new
      if(W_ALTERATION) then
        if(dehyd(i)) cycle
        if(mn(i)==15) cycle
        mass_ratio_anhydrous_rock(i) = 0._rkind
        mass_ratio_hydrous_rock(i) = 1._rkind
      elseif(.NOT.W_ALTERATION) then
        mass_ratio_anhydrous_rock(i) = mass_ratio_rock_new
        mass_ratio_hydrous_rock(i) = 0._rkind
      endif
      if(mn(i)==15) mass_ratio_water(i) = 1._rkind - mass_ratio_rock_new
    enddo
    do i = ir+2, iw
      mn(i) = 1
      mass_ratio_rock(i) = 0._rkind
      mass_ratio_anhydrous_rock(i) = 0._rkind
      mass_ratio_hydrous_rock(i) = 0._rkind
      mass_ratio_water(i) = 1._rkind
    enddo
    melt(0:iw) = .TRUE.
    mass_ratio_ice(0:iw) = 0._rkind
    return
  end subroutine rock_check