!> \file precribe-H_mod.f90 !! !! Determine mask and thickness for prescription of H !! !< !> \namespace prescribe_H !! !! Determine mask and thickness for prescription of H !! !! \author CatRitz !! \date june 2009 !! !! @note use module3D_phy !! @note use param_phy_mod !< module prescribe_H use module3D_phy use param_phy_mod use interface_input implicit none ! ces dimensionnements sont maintenant dans 3D ! real,dimension(nx,ny) :: Hp !< H value if prescribed ! real,dimension(nx,ny) :: Hp0 !< H value if prescribed (reference value) ! real,dimension(nx,ny) :: Delta_H !< Delta_H value if prescribed ! integer,dimension(nx,ny) :: i_delta_H !< 1 if Delta_H is prescribed on this node, else 0 ! integer,dimension(nx,ny) :: i_Hp !< 1 if H is prescribed on this node, else 0 ! integer,dimension(nx,ny) :: i_Hp0 !< i_hp mask reference value does not change with time ! integer,dimension(nx,ny) :: MK_gl0 !< mask grounding line initial ! integer,dimension(nx,ny) :: MK_flot0 !< mask float initial ! pour grounding line retrea, ice2sea integer,dimension(nx,ny) :: i_status !< 0 pas prescrit, 1 prescrit retrait, 2 prescrit shelf real :: time_onset ! moment ou la grounding line commence a reculer ! for prescribed grounding line retreat : ancienne version integer :: voisin !< work variable real (kind=kind(0.d0)) :: t_break !< temps : en double precision integer,parameter :: iretreat_max = 20 !< maximum number of retreat steps integer :: iretreat real (kind=kind(0.d0)),dimension(iretreat_max) :: time_gl !< time of the retreat steps integer, dimension(nx,ny) :: tag_gl !< tag of all the retreat points !< floating -> 10000 !< grounded not gl -> 1000 integer :: floating_tag = -10000 integer :: grounded_tag = 1000 integer, dimension(nx,ny) :: upwind_gl !< 1 if gl can still retreat, else 0 real, dimension(nx,ny) :: H_gl !< flotation thickness at the time real, dimension(nx,ny) :: dH_gl_dt !< to smoothly decrease ice thickness !< for next gl points ! nouvelle version recul grounding line pour ice2sea ! ------------------------------------------------------ ! real, dimension(nx,ny) :: H_gl !< flotation thickness at the time (deja declare autre version) real, dimension(nx,ny) :: H_beggin !< thickness when the point begins to be prescribed real, dimension(nx,ny) :: time_beggin !< time when the point begins to be prescribed real, dimension(nx,ny) :: time_just_gr !< time when the point is just grounded contains !______________________________________________________________________________________ !> initialise prescribe H !! calculate the initial grounding line position !! @return MK_flot0 and Mk_gl0 subroutine init_prescribe_H implicit none if (itracebug.eq.1) call tracebug(' Entree dans routine init_prescribe_H') ! keep the mask of floating points where (flot(:,:)) MK_flot0(:,:)= 1 elsewhere MK_flot0(:,:)= 0 end where ! determine the initial grounding line MK_gl0(:,:)=0 HP0(:,:)=H0(:,:) do j=2,ny-1 do i=2,nx-1 voisin = MK_flot0(i-1,j) + MK_flot0(i+1,j) + MK_flot0(i,j-1) + MK_flot0(i,j+1) if ((MK_flot0(i,j).eq.0).and.(voisin.gt.0)) then ! grounded and at least one neighbour floating MK_gl0(i,j) = 1 ! tagged mask grounding line end if end do end do call prescribe_grid_boundary_0 !call init_prescribe_retreat time_gl(1)= 1000000 ! pour un run sans recul t_break = 1000000. ! pour un run avec ice shelves prescrits ! Dans le cas du Groenland ice2sea, mk_init = 3 -> noeuds qui ne bougent pas if (geoplace(1:4).eq.'GI2S') then do j = 1,ny do i = 1,nx if (mk_init(i,j).eq.3) then i_hp0 (i,j) = 1 Hp0 (i,j) = H0(i,j) endif end do end do end if if (igrdline.eq.3) then call init_prescr_ice2sea_retreat end if end subroutine init_prescribe_H !______________________________________________________________________________________ !> function prescribe_present_H_gl !! calculate the initial grounding line position !! @return i_hp(:,:) and hp(:,:) subroutine prescribe_present_H_gl implicit none if (itracebug.eq.1) call tracebug(' Entree dans routine prescribe_present_H_gl') where ((MK_flot0(:,:).eq. 1).or.(MK_gl0(:,:) .eq. 1)) ! floating or grounding line i_hp(:,:) = 1 ! thickness prescribed to present value hp(:,:) = Hp0(:,:) end where if (itracebug.eq.1) call tracebug(' fin prescribe_present_H_gl ') end subroutine prescribe_present_H_gl !______________________________________________________________________________________ !> function prescribe_fixed_points !! fixes thickness for some points !! @return i_hp(:,:) and hp(:,:) subroutine prescribe_fixed_points if (itracebug.eq.1) call tracebug(' Entree dans routine prescribe_fixed_points') where (i_hp0(:,:).eq.1) ! les points i_hp0 le sont pour tout le run i_hp(:,:) = i_hp0(:,:) Hp(:,:) = Hp0(:,:) end where if (itracebug.eq.1) call tracebug(' fin prescribe_fixed_points') end subroutine prescribe_fixed_points !______________________________________________________________________________________ !> function prescribe_paleo_gl_shelf !! calculate the grounding line thickness !! no ice shelves !! @return i_hp(:,:) and hp(:,:) subroutine prescribe_paleo_gl_shelf implicit none if (itracebug.eq.1) call tracebug(' Entree dans routine prescribe_paleo_gl_shelf') ! noeuds qui doivent être imposés where ((MK_flot0(:,:).eq. 1).or.(MK_gl0(:,:) .eq. 1)) ! floating or grounding line i_hp(:,:) = 1 ! thickness prescribed end where ! valeur imposee where (MK_flot0(:,:).eq. 1) ! paleo shelf epaisseur a 1 hp(:,:) = 1. end where where (MK_gl0(:,:) .eq. 1) hp(:,:) = (sealevel-Bsoc(:,:))*row/ro +1. end where end subroutine prescribe_paleo_gl_shelf !______________________________________________________________________________________ !> function prescribe_grid_boundary !! prescribe the grid boundary thickness to 0 !! @return i_hp(:,:) and hp(:,:) subroutine prescribe_grid_boundary_0 implicit none if (itracebug.eq.1) call tracebug(' Entree dans routine prescribe_grid_boundary_0') ! prescribe thickness on the grid boundaries i_hp(:,:) = 0 i_hp(1,:) = 1 i_hp(nx,:)= 1 hp(1,:) = 0 hp(nx,:) = 0 i_hp(:,1) = 1 i_hp(:,ny)= 1 hp(:,1) = 0 hp(:,ny) = 0 ! valeurs de reference i_hp0(1,:) = 1 i_hp0(nx,:) = 1 hp0(1,:) = 0 hp0(nx,:) = 0 i_hp0(:,1) = 1 i_hp0(:,ny) = 1 hp0(:,1) = 0 hp0(:,ny) = 0 end subroutine prescribe_grid_boundary_0 !______________________________________________________________________________________ !______________________________________________________________________________________ !> lecture of a prescribe field of delta H !! !! @return Delta_H, i_delta_H subroutine lect_prescribed_deltaH implicit none character(len=100) :: delta_H_file ! variation d'epaisseur imposee if (itracebug.eq.1) call tracebug(' Entree dans routine lect_prescribed_deltaH') if (itracebug.eq.1) call tracebug(' Attention ne dvrait pas passer par la') ! lecture delta_H_file = 'sensibxg-25km.dat' delta_H_file = trim(dirnameinp)//trim(delta_H_file) call lect_input(3,'Delta_H',1,Delta_H,delta_H_file,trim(dirnameinp)//trim(runname)//'.nc') !call lect_datfile(nx,ny,Delta_H,1,delta_H_file) do j=1,ny do i=1,nx if (Delta_H(i,j).ne.0.) then i_delta_H(i,j)=1 i_Hp(i,j)=1 ! i_delta_H might be not completely included in i_Hp Hp0(i,j) = H(i,j)-Delta_H(i,j) ! warning Delta_H is substracted if (MK_flot0(i,j).eq. 1) then ! si flotte, pas de changement d'epaisseur Hp0(i,j)=H0(i,j) else if (Hp0(i,j).lt. -Bsoc(i,j)*row/ro) then ! devient flottant Hp0(i,j) = -Bsoc(i,j)*row/ro-10. ! floating thickness - 10 MK_flot0(i,j)= 1 ! est tagge masque flottant end if Hp(i,j) = Hp0(i,j) else i_delta_H(i,j)=0 end if end do end do end subroutine lect_prescribed_deltaH !______________________________________________________________________________________ !> substract Delta_H to the current thickness !! !! @return subroutine prescribe_deltaH implicit none if (itracebug.eq.1) call tracebug(' Entree dans routine prescribe_deltaH') where (i_delta_H(:,:).eq.1) Hp(:,:) = Hp0(:,:) i_Hp(:,:)=1 end where end subroutine prescribe_deltaH !______________________________________________________________________________________ !> break ice shelves !! !! subroutine break_all_ice_shelves implicit none if (itracebug.eq.1) call tracebug(' Entree dans routine break_all_ice_shelves ') debug_3D(:,:,89)=0. where (flot(:,:)) ! floating from the beginning and eventually from grounding line retreat i_hp(:,:) = 1 ! thickness prescribed to 1 m hp (:,:) = 1. hp0 (:,:) = 1. H(:,:) =1. debug_3D(:,:,89)=1. end where where (MK_flot0(:,:).eq. 1) i_hp(:,:) = 1 ! thickness prescribed to 1 m hp (:,:) = 1. hp0 (:,:) = 1. H(:,:) =1. end where end subroutine break_all_ice_shelves ! !______________________________________________________________________________________ !> break ice shelves !! !! subroutine melt_ice_shelves implicit none integer :: nbvoisins,nbdeglac,iv,jv real :: hmoy,hmin,hmax if (itracebug.eq.1) call tracebug(' Entree dans routine melt_ice_shelves') do j = 2,ny-1 do i=2, nx-1 if (flot(i,j)) then ! floating i_hp(i,j) = 1 ! thickness prescribed hmoy=0. hmin=0. hmax=0. nbvoisins=0 nbdeglac=0 do jv=-1,1 do iv=1,1 if (flot(i+iv,j+jv)) then nbvoisins=nbvoisins+1 hmoy=hmoy+h(i+iv,j+jv) hmin=min(hmin,h(i+iv,j+jv)) hmax=max(hmax,h(i+iv,j+jv)) if (h(i+jv,j+jv).lt.10.) nbdeglac=nbdeglac+1 end if end do end do hmoy=hmoy/(max(nbvoisins,1)) if(H(i,j).gt.hmoy) then ! melt only if smooth hp(i,j)=H(i,j)-1.*dt else if (hmax-hmin.lt.10.) then hp(i,j)=H(i,j)-1.*dt else if ((H(i,j).lt.50.).and.(Hmin.lt.30.)) then hp(i,j)=1. else if ((H(i,j).lt.50.).and.(nbdeglac.ge.4)) then hp(i,j)=1. endif hp(i,j) = max(hp(i,j),1.) end if end do end do end subroutine melt_ice_shelves ! !______________________________________________________________________________________ !> prescribe grounding line retreat !! !! time_gl(iretreat_max) : times of the retreat steps (in) !! tag_gl(nx,ny) : tag that gives at which time a node is prescribed !! as a grounding line node !! @return Hp !! @return i_hp subroutine prescribe_retreat implicit none integer, dimension(nx,ny) :: new_tag ! working array if (itracebug.eq.1) call tracebug(' Entree dans routine prescribe_retreat') if (itracebug.eq.1) write(num_tracebug,*) 'time_gl(1)=', time_gl(1) new_tag(:,:) = tag_gl(:,:) ! initialize new_tag start_time: if (time.lt.time_gl(1)) then ! before the perturbation if (itracebug.eq.1) call tracebug(' Avant appel prescribe_present_H_gl') else if (itracebug.eq.1) call tracebug('avant entree step_gl') iretreat_loop: do iretreat = 1,iretreat_max-1 if (itracebug.eq.1) write(num_tracebug,*) 'iretreat',iretreat find_step: if (abs(time-time_gl(iretreat)).lt.dtmin/2.) then ! exactly on a retreat step ! write(6,*) time, 'sur un step',iretreat where (tag_gl(:,:).eq.iretreat) H_gl(:,:) = -Bsoc(:,:)*row/ro ! flotation thickness Hp(:,:) = H_gl(:,:)+5. ! +5 to be sure to be grounded ! even with small isostatic changes i_Hp(:,:) = 1 elsewhere (tag_gl(:,:).eq.iretreat+1) ! points that will be gl at next step H_gl(:,:) = -Bsoc(:,:)*row/ro ! flotation thickness dH_gl_dt(:,:) = (H(:,:)-H_gl(:,:))-5. ! thickness to be reduced to reach gl dH_gl_dt(:,:) = dH_gl_dt(:,:) / (time_gl(iretreat+1)-time_gl(iretreat)) ! rate of decrease elsewhere ((tag_gl(:,:).eq.iretreat-1).and. & ! points that become floating (upwind_gl(:,:).eq.1)) H_gl(:,:) = -Bsoc(:,:)*row/ro ! flotation thickness Hp(:,:) = H_gl(:,:)-10. ! to be sure to be floating i_Hp(:,:) = 1 new_tag(:,:) = floating_tag end where ! pour une raison mysterieuse la boucle suivante marche et celle en where ne marche pas ! (memory falut) do j=1,ny do i=1,nx if ((tag_gl(i,j).eq.iretreat-1).and. & ! points that stay gl (upwind_gl(i,j).eq.0)) then H_gl(i,j) = -Bsoc(i,j)*row/ro ! flotation thickness Hp(i,j) = H_gl(i,j)+5. ! to be sure to be grounded i_Hp(i,j) = 1 new_tag(i,j) = iretreat ! tagged with the new step number ! write(6,*) i,j,tag_gl(i,j), upwind_gl(i,j) endif end do end do !!$ where ((tag_gl(:,:).eq.iretreat-1).and. & ! points that stay gl !!$ (upwind_gl(:,:).eq.0)) !!$ !!$ H_gl(:,:) = -Bsoc(:,:)*row/ro ! flotation thickness !!$ Hp(:,:) = H_gl(:,:)+5. ! to be sure to be grounded !!$ i_Hp(:,:) = 1 !!$ new_tag(:,:) = iretreat ! tagged with the new step number !!$ end where else if ((time.gt.time_gl(iretreat)+dtmin/2.).and.(time.lt.time_gl(iretreat+1))) then ! write(6,*) time, 'entre',time_gl(iretreat),time_gl(iretreat+1) where (tag_gl(:,:).eq.iretreat) H_gl(:,:) = -Bsoc(:,:)*row/ro ! flotation thickness Hp(:,:) = H_gl(:,:)+5. ! +5 to be sure to be grounded ! even with small isostatic changes i_Hp(:,:) = 1 elsewhere (tag_gl(:,:).eq.iretreat+1) ! points that will be gl at next step Hp(:,:) = H(:,:)-dH_gl_dt(:,:)*(time-time_gl(iretreat)) i_Hp(:,:) = 1 elsewhere (tag_gl(:,:).eq.floating_tag) Hp(:,:) = H(:,:) i_Hp(:,:) = 1 end where else if (time.gt.time_gl(iretreat_max)) then ! write(6,*) time, 'apres iretreat_max',iretreat,time_gl(iretreat_max) where (tag_gl(:,:).eq.iretreat_max) H_gl(:,:) = -Bsoc(:,:)*row/ro ! flotation thickness Hp(:,:) = H_gl(:,:)+5. ! +5 to be sure to be grounded ! even with small isostatic changes i_Hp(:,:) = 1 elsewhere (tag_gl(:,:).eq.floating_tag) Hp(:,:) = H(:,:) i_Hp(:,:) = 1 end where !faire une sortie de la boucle step_gl end if find_step end do iretreat_loop end if start_time tag_gl(:,:) = new_tag(:,:) debug_3D(:,:,84) = tag_gl(:,:) debug_3D(:,:,85) = i_Hp(:,:) debug_3D(:,:,86) = Hp(:,:) end subroutine prescribe_retreat ! !______________________________________________________________________________________ !> initialise prescribe retreat !! !! @return tag_gl !! @return upwind_gl !! @return time_gl subroutine init_prescribe_retreat implicit none real, dimension(nx,ny) :: lect_tab if (itracebug.eq.1) call tracebug(' Entree dans routine init_prescribe_retreat') ! read the prescribed gl file call lect_input(3,'lecttab',1,lect_tab,trim(dirnameinp)//'masque_gl_15km_final.dat',trim(dirnameinp)//trim(runname)//'.nc') !call lect_datfile(nx,ny,lect_tab,1,trim(dirnameinp)//'masque_gl_15km_final.dat') tag_gl(:,:)=nint(lect_tab(:,:)) ! look if neighbors could become grounding line later upwind_gl(:,:) = 0 do j=2,ny-1 do i=2,nx-1 if ((tag_gl(i,j).ne.floating_tag).and.(tag_gl(i,j).ne.grounded_tag)) then if (any(tag_gl(i-1:i+1,j-1:j+1).eq.tag_gl(i,j)+1)) then upwind_gl(i,j) = 1 end if end if end do end do ! time of the retreat steps do i=1,iretreat_max time_gl(i)=(i-1)*30.+100000. end do time_gl(iretreat)=200000. write(num_rep_42,*) '! Grounding line retreat prescribed' write(num_rep_42,*) '! time_gl begin et suivant : ',time_gl(1),time_gl(2), & ' time_glmax=',time_gl(iretreat) ! beginning of breaking or melting ice shelves t_break = 1000000.d00 write(num_rep_42,*) '! t_break=', t_break end subroutine init_prescribe_retreat !______________________________________________________________________________________ !______________________________________________________________________________________ !> function prescribe_present_H_gl !! calculate the initial grounding line position !! @return i_hp(:,:) and hp(:,:) subroutine prescribe_present_H_gl_copy implicit none if (itracebug.eq.1) call tracebug(' Entree dans routine prescribe_present_H_gl-copy') do j=1,ny do i=1,nx if ((MK_flot0(i,j).eq. 1).or.(MK_gl0(i,j) .eq. 1)) then i_hp(i,j) = 1 ! thickness prescribed to present value hp(i,j) = Hp0(i,j) end if end do end do !!$where ((MK_flot0(:,:).eq. 1).or.(MK_gl0(:,:) .eq. 1)) ! floating or grounding line !!$ i_hp(:,:) = 1 ! thickness prescribed to present value !!$! hp(:,:) = 1 !Hp0(:,:) !!$end where end subroutine prescribe_present_H_gl_copy !______________________________________________________________________________________ !> initialise prescr_ice2sea_retreat !! !! @return tag_gl !! @return upwind_gl !! @return time_gl ! nouvelle routine de recul impose de grounding line pour ice2sea : decembre 2011 ! cette routine transforme les informations distance en information temps ou un noeud sera just grounded" subroutine init_prescr_ice2sea_retreat implicit none character(len=100) :: dist_file !< nom fichier contenant les distances character(len=100) :: file_ncdf !< nom bidon pour lecture netcdf real, dimension(nx,ny) :: dist_gl ! distance à la grounding line ( a lire) real :: vretrait ! vitesse de retrait real :: delta_t real :: time_max ! pour calcul time_max integer :: i1 integer :: j1 namelist/retreat_ice2sea/dist_file,vretrait,time_onset if (itracebug.eq.1) call tracebug(' Entree dans routine init_prescr_ice2sea_retreat') ! lecture des noms de fichiers et parametres !-------------------------------------------- 428 format(A) rewind(num_param) ! pour revenir au debut du fichier param_list.dat read(num_param,retreat_ice2sea) write(num_rep_42,428)'!___________________________________________________________' write(num_rep_42,428)'! prescr_ice2sea_retreat ' write(num_rep_42,retreat_ice2sea) write(num_rep_42,428)'!___________________________________________________________' ! lecture de la carte de distance dites !------------------------------------- ! distance file : distances en km dist_file = trim(dirnameinp)//trim(dist_file) call lect_input(1,'dist_gl',1,dist_gl,dist_file,file_ncdf) ! les zones interdites sont taggees -1 -> les transformer en +100000 where ((dist_gl(:,:).eq.-1.).or.(Bsoc0(:,:).gt.0.)) dist_gl(:,:) = 100000. end where ! calcul de l'epaisseur H_gl, correction des distances pour mieux prendre en compte le premier recul where(Mk_flot0(:,:).eq.0) ! points poses y compris GL h_gl(:,:) = (sealevel-Bsoc(:,:))*row/ro ! epaisseur flottaison dist_gl(:,:) = dist_gl(:,:)+2.+2.5 ! 2 Pour avoir au moins 0 sur la GL ! 2.5 = demi maille 5 km time_just_gr(:,:) = dist_gl(:,:)/Vretrait + time_onset ! aussi valable pour les points GL H_beggin(:,:) = H0(:,:) ! pour initialiser, strictement vrai pour la GL elsewhere ! points flottants toujours prescrit time_just_gr(:,:) = time_onset-100000. time_beggin(:,:) = time_onset-110000 i_status(:,:) = 2 H_beggin(:,:) = H0(:,:) end where ! calcul specifique a la grounding line where(Mk_gl0(:,:).eq.1) time_beggin(:,:) = time_onset H_beggin(:,:) = H0(:,:) end where ! recherche de Time_beggin : pour les points grounded mais pas gl delta_t = dx/1000./Vretrait -0.1 ! temps pour traverser une maille do j=2,ny-1 do i=2,nx-1 if((Mk_flot0(i,j).eq.0).and.(Mk_gl0(i,j).eq.0)) then ! grounded mais pas GL time_max = time_onset do j1 = j-1,j+1 do i1 = i-1,i+1 if (Time_just_gr(i1,j1).lt.Time_just_gr(i,j)-delta_t) then ! 0.1 par securite time_max = max(time_max,Time_just_gr(i1,j1)) endif end do end do if (Time_just_gr(i,j)-time_max.lt.2.*delta_t) then time_beggin(i,j) = time_max else time_beggin(i,j) = Time_just_gr(i,j)-2.*delta_t end if end if end do end do ! impression pour verifier ( a enlever apres) !do j=1,ny ! do i=1,nx ! write(666,'(i3,1x,i3,3(f8.0,1x))') i,j,Time_just_gr(i,j),time_beggin(i,j),Time_just_gr(i,j)-time_beggin(i,j) ! end do !end do end subroutine init_prescr_ice2sea_retreat !______________________________________________________________________________________ subroutine prescr_ice2sea_retreat implicit none real :: t_local ! pour les tests car time est en double precision if (itracebug.eq.1) call tracebug(' Entree dans routine prescr_ice2sea_retreat') !temps dans la boucle t_local = time ! avant le demarrage la grounding line est fixee a sa valeur actuelle if (t_local.lt.time_onset) then where ((MK_flot0(:,:).eq. 1).or.(MK_gl0(:,:) .eq. 1)) ! floating or grounding line i_hp(:,:) = 1 ! thickness prescribed to present value hp(:,:) = Hp0(:,:) end where else ! calcul du H prescrit avec recul do j=2,ny-1 do i= 2,nx-1 ! point pas encore impose time_test: if (t_local.lt.time_beggin(i,j)) then i_Hp(i,j) = 0 i_status(i,j) = 0 ! point en periode de retrait else if ( (t_local.ge.time_beggin(i,j)).and.(t_local.le.time_just_gr(i,j)) ) then retrait: if (i_status(i,j).eq.0) then ! on vient de rentrer dans la periode retrait i_status(i,j) = 1 H_beggin(i,j) = H(i,j) i_Hp(i,j) = 0 ! on n'impose pas H le premier coup else if (i_status(i,j) .eq. 1) then i_Hp(i,j) = 1 ! on impose ce point h_gl(i,j) = (sealevel-Bsoc(i,j))*row/ro ! epaisseur flottaison ! decroissance lineaire entre time_beggin et time_just_gr Hp(i,j) = (H_beggin(i,j) - H_gl(i,j)) * ( (time_just_gr(i,j) - t_local) & / (time_just_gr(i,j) - time_beggin(i,j))) + H_gl(i,j) end if retrait ! point en periode flottantte else if (t_local.gt.time_just_gr(i,j)) then h_gl(i,j) = (sealevel-Bsoc(i,j))*row/ro ! epaisseur flottaison i_status(i,j) = 2 i_Hp(i,j) = 1 ! on impose ce point Hp(i,j) = min(H_gl(i,j)-10., H0(i,j)) ! minimum flottaison ou epaisseur initiale end if time_test end do end do end if end subroutine prescr_ice2sea_retreat end module prescribe_H