Changeset 146 for branches


Ignore:
Timestamp:
10/20/17 09:31:39 (7 years ago)
Author:
aquiquet
Message:

Grisli-iLoveclim branch: merged to trunk at revision 145

Location:
branches/iLoveclim
Files:
27 edited
6 copied

Legend:

Unmodified
Added
Removed

  • branches/iLoveclim

  • branches/iLoveclim/SOURCES/3D-physique-gen_mod.f90

    r126 r146  
    3636  integer ::  err                       !< pour l'allocation des tableaux 
    3737  integer ::  igrdline                  !< si 1 fixe la position en jouant sur la fusion shelf 
     38  integer ::  ibmelt_inv                !< si 1 inversion du bmelt (avec igrdline=1) 
    3839  integer ::  i_resolmeca               !< defini le type d'association SIA-L1 
    3940  integer ::  iter_beta                 !< pour la determination du frottement 
    … …  
    225226  integer,dimension(nx,ny) :: FRONTFACEY!< type de front sur les faces y 
    226227  integer,dimension(nx,ny) :: gr_line_schoof ! points ou on impose le flux de schoof (pour sorties) 
     228  integer,dimension(nx,ny) :: gr_line   !< points grounding line pour les sorties 
    227229  integer,dimension(nx,ny) :: ICE       !< presence de glace au point considere, 
    228230  !< seuil 0 si pose, 1 si flottant 
    … …  
    443445 
    444446  real,dimension(nx,ny) :: debug_2D 
    445   real,dimension(nx,ny,199) :: debug_3D 
     447  real,dimension(nx,ny,220) :: debug_3D 
    446448 
    447449  !  tableaux qui permettent de faire des sorties netcdf avec les variables qui ne sont pas globales 

  • branches/iLoveclim/SOURCES/Ant40_files/lect-anteis_mod.f90

    r123 r146  
    228228    do J=1,NY 
    229229       do I=1,NX 
    230           if (((BSOC(I,J)+H(I,J)*RO/ROW -SEALEVEL).LT.0.).and.(H(I,J).gt.1.E-3)) then 
     230          if ((BSOC(I,J)+H(I,J)*RO/ROW -SEALEVEL).LT.0.) then 
    231231             FLOT(I,J)=.TRUE. 
    232232          else 

  • branches/iLoveclim/SOURCES/Ant40_files/output_anta40_mod-0.4.f90

    r123 r146  
    1515USE module3D_phy 
    1616use bilan_eau_mod 
    17  
     17use netcdf 
     18use io_netcdf_grisli 
     19use bilan_flux_mod 
    1820 
    1921implicit none 
    2022       
    2123!real ::  vol ; !integer :: np 
     24integer :: nflot                      !< nbr de point flottant 
    2225real ::  bmean                        !< 
    2326real ::  accmean                      !< accumulation moyenne 
    … …  
    3538real ::  hdotmean                     !< moyenne derivee / temps de h 
    3639real ::  bdotmean                     !< moyenne bedrock derive / temps 
    37 real :: volf                          !< volume au dessus de la flottaison 
     40real ::  volf                         !< volume au dessus de la flottaison 
     41 
     42real :: lim                           !< Total ice mass 
     43real :: iareag                        !< surface posee 
     44real :: iareaf                        !< surface flottante 
     45real :: tendacabf                     !< Total SMB flux 
     46real :: tendlibmassbf                 !< Total Basal mass balance flux 
     47real :: tendlicalvf                   !< Total calving flux 
     48real :: tendligroundf                 !< Total grounding line flux 
     49 
     50real,dimension(nx,ny) :: corrsurf     !< facteur de correction de la surface 
     51!real,parameter :: ice_density=910.    !< densite de la glace pour conversion en masse 
     52 
     53! variables netcdf 
     54integer,parameter :: ncshortout=1 ! 1 sorties netcdf short initMIP 
     55integer,parameter :: nvar=8 ! nombre de variables dans le fichier de sortie temporel netcdf  
     56integer :: ncid 
     57integer :: status 
     58integer :: timeDimID 
     59integer :: timeVarID 
     60integer,dimension(nvar) :: varID 
     61integer :: nbtimeout  ! index time output 
     62 
     63real,dimension(nvar) :: var_shortoutput 
    3864 
    3965 
    … …  
    4167 
    4268subroutine init_outshort 
     69 
     70double precision,dimension(:,:),pointer      :: tab               !< tableau 2d real pointer 
     71character(len=100),dimension(nvar) :: namevar ! name, standard_name, long_name, units 
     72character(len=100),dimension(nvar) :: standard_name 
     73character(len=100),dimension(nvar) :: long_name 
     74character(len=100),dimension(nvar) :: units 
    4375 
    4476!ndisp sorite courte tous les ndisp 
    4577NDISP=100 
     78 
     79 
     80 
     81if (ncshortout.eq.1) then  ! ecriture netcdf 
     82 
     83! lecture du fichier avec les corrections de surface 
     84  call Read_Ncdf_var('z',trim(DIRNAMEINP)//'/corrsurf-initMIP-16km.grd',tab) 
     85  corrsurf(:,:)=tab(:,:) 
     86 
     87  open(568,file=trim(dirsource)//'/Fichiers-parametres/short-initMIPnc.dat',status='old') 
     88! lecture en-tete 
     89  read(568,*) 
     90  read(568,*) 
     91  read(568,*) 
     92  read(568,*) 
     93  read(568,*) 
     94! lecture des infos sur les variables : 
     95  do k=1,nvar 
     96    read(568,'(a100)') namevar(k) 
     97    read(568,'(a100)') standard_name(k) 
     98    read(568,'(a100)') long_name(k) 
     99    read(568,'(a100)') units(k) 
     100    read(568,*) 
     101  enddo 
     102  close(568) 
     103! Fichier Netcdf initMIP 
     104! creation du fichier Netcdf : 
     105  status=nf90_create(path = 'short'//runname//'.nc', cmode = nf90_clobber, ncid = ncid) 
     106  if (status /= nf90_noerr) call handle_err(status) 
     107 
     108! definition des dimension : 
     109  status=nf90_def_dim(ncid, name="time", len=NF90_UNLIMITED, dimid=timeDimID) 
     110  if (status /= nf90_noerr) call handle_err(status) 
     111  status=nf90_def_var(ncid, name="time", xtype=nf90_float, dimids=(/ timeDimID/), varid=timeVarID)  
     112  if (status /= nf90_noerr) call handle_err(status) 
     113  status=nf90_put_att(ncid, timeVarID, "standard_name", "time") 
     114  if (status /= nf90_noerr) call handle_err(status) 
     115  status=nf90_put_att(ncid, timeVarID,"units", "years since 2007-01-01 00:00:00") 
     116  if (status /= nf90_noerr) call handle_err(status) 
     117 
     118! definition des variables de sortie : 
     119  do k=1,nvar  ! boucle sur le nbr de variable a definir 
     120    status=nf90_def_var(ncid, name=trim(namevar(k)), xtype=nf90_float, dimids= & 
     121          (/ timeDimID /), varid=varID(k)) 
     122    if (status /= nf90_noerr) call handle_err(status) 
     123    status=nf90_put_att(ncid, varID(k), "standard_name", trim(standard_name(k))) 
     124    if (status /= nf90_noerr) call handle_err(status) 
     125    status=nf90_put_att(ncid, varID(k), "long_name", trim(long_name(k))) 
     126    if (status /= nf90_noerr) call handle_err(status) 
     127    status=nf90_put_att(ncid, varID(k), "units", trim(units(k))) 
     128    if (status /= nf90_noerr) call handle_err(status) 
     129  enddo 
     130 
     131! fin de la definition du fchier : 
     132  status=nf90_enddef(ncid) 
     133  if (status /= nf90_noerr) call handle_err(status) 
     134  nbtimeout = 0 ! initialisation compteur sorties axe time 
     135else ! pas de sortie netcdf et sans correction de surface 
     136  corrsurf(:,:)=1. 
     137endif 
     138 
    46139end subroutine init_outshort 
    47140 
    … …  
    54147!------------------ 
    55148real ::  smax 
     149 
    56150      vol=0.  
    57151      np=0 
     152      nflot=0 
    58153      hmax=0.  
    59154      smax=0. 
    … …  
    74169      bdotmean=0. 
    75170      volf=0. 
    76   
     171      lim=0. 
     172      tendacabf=0. 
     173      iareag=0. 
     174      iareaf=0. 
     175      tendlicalvf=0. 
     176      tendligroundf=0. 
     177      tendlibmassbf=0. 
    77178! 2_preparing outputs 
    78179!--------------------       
    79     do i=1,nx  
    80       do j=1,ny 
    81         if (.not.flot(i,j)) then  
    82 !       if (h(i,j).gt.1.) then  
    83           np=np+1 
    84           vol=vol+h(i,j)  
    85  
    86 !        calcul de la hauteur au dessus de la flottaison 
    87          if (sealevel-B(i,j).le.0.) then    ! socle au dessus du niveau des mers 
    88                volf=volf+h(i,j) 
    89          else 
    90             volf=volf+h(i,j)-row/ro*(sealevel-b(i,j)) 
    91          endif 
    92  
    93  
    94           if (h(i,j).gt.hmax) hmax=h(i,j)  
    95           if (s(i,j).gt.smax) smax=s(i,j)  
    96           bmean=bm(i,j)+bmean  
    97           accmean=acc(i,j)+accmean 
    98           tbmean=tbmean+t(i,j,nz) 
    99           tbdotmean=tbdotmean+tbdot(i,j) 
    100           vsmean=vsmean+sqrt(ux(i,j,1)**2+uy(i,j,1)**2) 
     180    do j=1,ny 
     181      do i=1,nx 
     182        if (ice(i,j).eq.1) then ! point englace 
     183          if (.not.flot(i,j)) then ! point pose  
     184            np=np+1 
     185            vol=vol+h(i,j) 
     186            iareag=iareag+1.*corrsurf(i,j)                            ! surface englacee posee 
     187 
     188!         calcul de la hauteur au dessus de la flottaison 
     189            if (sealevel-B(i,j).le.0.) then    ! socle au dessus du niveau des mers 
     190              volf=volf+h(i,j)*corrsurf(i,j)                          ! volume au-dessus de la flottaison 
     191            else 
     192              volf=volf+(h(i,j)-row/ro*(sealevel-b(i,j)))*corrsurf(i,j) ! volume au-dessus de la flottaison 
     193            endif 
     194 
     195            if (h(i,j).gt.hmax) hmax=h(i,j)  
     196            if (s(i,j).gt.smax) smax=s(i,j)  
     197            bmean=bm(i,j)+bmean  
     198            accmean=acc(i,j)+accmean 
     199            tbmean=tbmean+t(i,j,nz) 
     200            tbdotmean=tbdotmean+tbdot(i,j) 
     201            vsmean=vsmean+sqrt(ux(i,j,1)**2+uy(i,j,1)**2) 
    101202!          vsdotmean=vsdotmean+vsdot(i,j) 
    102203!          uzsmean=uzsmean+uz(i,j,1) 
    103           uzsdotmean=uzsdotmean+uzsdot(i,j) 
    104           uzkmean=uzkmean+uzk(i,j) 
    105           hdotmean=hdotmean+abs(hdot(i,j))  
    106           bdotmean=bdotmean+abs(bdot(i,j))  
    107           bmeltmean=bmeltmean+bmelt(i,j) 
     204            uzsdotmean=uzsdotmean+uzsdot(i,j) 
     205            uzkmean=uzkmean+uzk(i,j) 
     206            hdotmean=hdotmean+abs(hdot(i,j))  
     207            bdotmean=bdotmean+abs(bdot(i,j))  
     208            bmeltmean=bmeltmean+bmelt(i,j) 
     209          else ! point flottant 
     210            iareaf=iareaf+1.*corrsurf(i,j)                               ! surface flottante 
     211          endif 
     212          lim=lim+h(i,j)*corrsurf(i,j)                                   ! volume total de glace 
     213          tendacabf=tendacabf+bm_dtt(i,j)*corrsurf(i,j)/dtt              ! smb surface 
     214          tendlibmassbf=tendlibmassbf-bmelt_dtt(i,j)*corrsurf(i,j)/dtt   ! fonte basale      
    108215        endif 
    109         calvmean=calvmean+calv(i,j)   
    110         ablbordmean=ablbordmean+ablbord(i,j) 
     216        tendlicalvf=tendlicalvf+calv_dtt(i,j)*corrsurf(i,j)/dtt        ! calving 
     217        tendlibmassbf=tendlibmassbf-ablbord_dtt(i,j)*corrsurf(i,j)/dtt ! partie ablbord de la fonte basale 
     218        tendligroundf=tendligroundf+grline_dtt(i,j)*corrsurf(i,j)/dtt  ! flux a la grounding line 
    111219      end do 
    112       end do 
     220    end do 
    113221 
    114222 
    … …  
    116224        hmean=vol/np  
    117225        vol=vol*dx*dy  
    118         volf=volf*dx*dy 
     226        volf=volf*dx*dy*ice_density 
    119227        bmean=bmean/np  
    120228        accmean=accmean/np  
    … …  
    130238        uzsdotmean=uzsdotmean/np 
    131239        uzkmean=uzkmean/np 
    132         hdotmean=hdotmean/np  
    133       endif  
    134  
     240        hdotmean=hdotmean/np 
     241      endif 
     242      lim=lim*dx*dy*ice_density 
     243      iareag=iareag*dx*dy  
     244      iareaf=iareaf*dx*dy 
     245      tendacabf=tendacabf*dx*dy*ice_density/secyear 
     246      tendlibmassbf=tendlibmassbf*dx*dy*ice_density/secyear 
     247      tendlicalvf=tendlicalvf*dx*dy*ice_density/secyear 
     248      tendligroundf=tendligroundf*ice_density/secyear 
     249       
     250       
    135251      bdotmean=bdotmean/nx/ny   
    136252 
    … …  
    154270!940   format('%%%% ',a,'   time=',f8.0,' %%%%') 
    155271 
     272 
     273if (ncshortout.eq.1) then  ! ecriture netcdf 
     274  ! Total ice mass 
     275  var_shortoutput(1)=lim 
     276  ! Mass above floatation 
     277  var_shortoutput(2)=volf 
     278  ! Grounded ice area 
     279  var_shortoutput(3)=iareag 
     280  ! Floating ice area 
     281  var_shortoutput(4)=iareaf 
     282  ! Total SMB flux  
     283  var_shortoutput(5)=tendacabf 
     284  ! Total Basal mass balance flux 
     285  var_shortoutput(6)=tendlibmassbf 
     286  ! Total calving flux  
     287  var_shortoutput(7)=tendlicalvf 
     288  ! Total grounding line flux 
     289  var_shortoutput(8)=tendligroundf 
     290 
     291  nbtimeout=nbtimeout+1 
     292   
     293  status=nf90_put_var(ncid, timeVarID, time, start=(/nbtimeout/)) 
     294  if (status /= nf90_noerr) call handle_err(status) 
     295 
     296  do k=1,nvar  ! boucle sur le nbr de variable a ecrire 
     297    status=nf90_put_var(ncid, varID(k), var_shortoutput(k),start=(/nbtimeout/)) 
     298    if (status /= nf90_noerr) call handle_err(status) 
     299  enddo 
     300  status=nf90_sync(ncid) 
     301  if (status /= nf90_noerr) call handle_err(status) 
     302endif 
     303 
     304 
    156305end subroutine shortoutput 
     306 
     307subroutine handle_err(status) 
     308  integer, intent(in) :: status 
     309  if (status /= nf90_noerr) then 
     310     print*,trim(nf90_strerror(status)) 
     311     stop "stopped" 
     312  end if 
     313end subroutine handle_err 
     314 
    157315end module  output_antarcti_mod 

  • branches/iLoveclim/SOURCES/Fichiers-parametres/hemin40_param_list.dat

    r123 r146  
    4343!___________________________________________________________ 
    4444&topo_file 
    45  topo_ref        = "Greenland_Bamber-et-al_composite_hemin40.nc" !"etopo1_grisliN_GMT.g40"     
    46  topo_dep        = "Greenland_Bamber-et-al_composite_hemin40.nc" !"etopo1_grisliN_GMT.g40"     
     45 topo_ref        = "Greenland_Bamber-et-al_composite_hemin40.nc" !"SHB_GLAC1D_21k_hemin40.nc" !"ICE6G_21k.nc" !"Greenland_Bamber-et-al_composite_hemin40.nc" 
     46 topo_dep        = "Greenland_Bamber-et-al_composite_hemin40.nc" !"SHB_GLAC1D_21k_hemin40.nc" !"ICE6G_21k.nc" !"Greenland_Bamber-et-al_composite_hemin40.nc" 
    4747 grid_topo       = "coord-nord-40km.dat" 
    4848 ghf_fich        = "ijphi_hemin40.nc" 
    … …  
    7979 hmax_till      =    20.00000     
    8080 poro_till      =   0.5000000     
    81  kond0          =   1.E-6 !10.000000E-05 
     81 kond0          =   10.000000E-05 
    8282/ 
    8383! hmax_till (m) : epaisseur max du sediment  
    … …  
    128128&drag_neff_slope              ! nom du bloc dragging neff slope 
    129129 
    130 cf              =       2.e-5             ! 1.e-4 
     130cf              =       20.e-5             ! 1.e-4 
    131131betamax         =       1000. 
    132132betamin         =       10. 
    … …  
    169169  
    170170Hcoup_plateau   =       50     !   tres petit quand shelves fixes sinon 250 
    171 Hcoup_abysses   =       300     !   tres petit quand shelves fixes sinon 250 
    172 prof_plateau   =       50     !   tres petit quand shelves fixes sinon 250 
    173 prof_abysses   =       1000     !   tres petit quand shelves fixes sinon 250 
     171Hcoup_abysses   =       400     !   tres petit quand shelves fixes sinon 250 
     172prof_plateau   =       100     !   tres petit quand shelves fixes sinon 250 
     173prof_abysses   =       800     !   tres petit quand shelves fixes sinon 250 
    174174ifrange         =       4 
    175175meth_Hcoup      =       0 

  • branches/iLoveclim/SOURCES/Hemin40_files/output_hemin40_mod.f90

    r126 r146  
    3030! afq for CONSEAU REAL, dimension(nx,ny) :: old_H_dtt      ! Epaisseur de glace au pas de temps precedent 
    3131integer, dimension(nx,ny) :: write_mask 
     32real :: ilc_units                     !afq, unit change from Grisli to iLoveclim 
    3233 
    3334CONTAINS 
    … …  
    184185! afq for CONSEAU      REAL, dimension(nx,ny) :: delta_H_dtt 
    185186 
     187      ilc_units = dx*dy/DICE 
    186188 
    187189!     open(unit=4145,file='reg_output_nord.dat') 
    … …  
    330332 
    331333! afq -- iLOVECLIM water conservation:  
    332       trendWAC        = diff_H*dx*dy !or water_bilan 
    333       smbWAC(:,:)     = (ablbord_dtt(:,:)+ Bm_dtt(:,:))*dx*dy/dtt 
    334       bmeltWAC(:,:)   = Bmelt_dtt(:,:)*dx*dy/dtt 
    335       calvingWAC(:,:) = Calv_dtt(:,:)*dx*dy/dtt 
     334!        units: GRISLI is in m i.e. / yr 
     335!               iLOVECLIM expects m3 w.e. /yr 
     336      trendWAC        = diff_H*ilc_units !or water_bilan 
     337      smbWAC(:,:)     = (-ablbord_dtt(:,:)+ Bm_dtt(:,:))*ilc_units/dtt 
     338      bmeltWAC(:,:)   = -Bmelt_dtt(:,:)*ilc_units/dtt 
     339      calvingWAC(:,:) = Calv_dtt(:,:)*ilc_units/dtt 
    336340!check-> 
    337341      if (abs(1-diff_H/(water_bilan+1e-20)).gt.0.01) then 
    338          write(*,*) "Water is not conserved in GRISLI", diff_H, water_bilan, abs(1-diff_H/(water_bilan+1e-20)) 
     342         write(*,*) "Water not conserved in GRISLI", diff_H, water_bilan, abs(1-diff_H/(water_bilan+1e-20)) 
    339343      endif 
    340344 

  • branches/iLoveclim/SOURCES/Makefile.grisli.inc

    r123 r146  
    4343        climat-forcage-insolation_mod.o climat_GrIce2sea_years_mod.o \ 
    4444        climat_GrIce2sea_years_perturb_mod.o \ 
     45        climat_InitMIP_years_perturb_mod.o \ 
    4546        climat-perturb_mod-0.4.o \ 
    4647        ablation_mod.o no_ablation_mod.o \ 
    … …  
    7778        calving_frange.o no_calving.o no_lakes.o \ 
    7879        out_profile_mod.o printtable_mod.o mix-SIA-L1_mod.o \ 
    79         furst_schoof_mod.o \ 
     80        furst_schoof_mod.o bilan_flux_output_mod.o \ 
    8081        relaxation_water_diffusion.o \ 
    8182        prescribe-H-i2s_mod.o  \ 
    … …  
    9798        calving_frange.o no_calving.o no_lakes.o \ 
    9899        out_profile_mod.o printtable_mod.o mix-SIA-L1_mod.o \ 
    99         furst_schoof_mod.o \ 
     100        furst_schoof_mod.o bilan_flux_output_mod.o \ 
    100101        relaxation_water_diffusion.o \ 
    101102        prescribe-H-i2s_mod.o   \ 
    … …  
    187188 
    188189 
     190Liste_ANT16 = output_anta40_mod-0.4.o \ 
     191                dragging_prescr_beta_mod.o \ 
     192                dragging_prescr_beta_nolin_mod.o \ 
     193                lect-Ant_gen2010_dat.o \ 
     194                lect-anteis_mod.o \ 
     195                bmelt-ant-regions_mod.o \ 
     196                bmelt-ant-regions-initmip_mod.o \ 
     197                fake-routines-ant_mod.o \ 
     198                beta_iter_vitbil_mod.o \ 
     199                module_choix-ant16km.o \ 
     200                massb-ant_perturb_Tparam.o \ 
     201                track_ant40_mod.o 
     202 
    189203Liste_ANT15-LBq = output_anta_mod-0.4.o \ 
    190204                dragging_prescr_beta_mod.o \ 
    … …  
    197211                module_choix_antar15_LBq.o \ 
    198212                massb-ant_perturb_Tparam.o \ 
    199                 track_ant40_mod.o \ 
     213                track_ant40_mod.o 
    200214 
    201215Liste_hudson = sedim_declar_hudson_mod.o climat-hudson_mod.o \ 
    … …  
    296310 
    297311Dim_GrIce2sea = paradim-GrIce2sea-cut_Tamsin.o geography-GrIce2sea.o 
     312 
     313Dim_ANT16 = paradim-ant16_mod.o geography-ant16.o 
    298314 
    299315Dim_ANT15-LBq    = paradim-ant15_LBq_mod.o geography-Ant15_LBq.o    
    … …  
    382398%.o: Ant40_files/%.f90 
    383399        $(FT) $(NCDF_INC) -c Ant40_files/$*.f90 
     400 
     401# Ant16_files 
     402%.o: Ant16_files/%.f90 
     403        $(FT) $(NCDF_INC) -c Ant16_files/$*.f90 
    384404         
    385405# ANT15-LBq_files 
    … …  
    436456#--------------------------------------- 
    437457 
     458Ant-16 : $(Dim_ANT16) $(mod_dim_communs) \ 
     459        $(toy_recul) \ 
     460        $(mod_communs) \ 
     461        $(mod_clim_tof)  \ 
     462        $(mod_no_tracers) \ 
     463        $(mod_ell) $(Liste_ANT16) \ 
     464        $(diagnoshelf) \ 
     465        $(Liste_Netcdf) \ 
     466        $(routines_communes) steps_time_loop.o \ 
     467        $(routine_elliptiques) \ 
     468        $(Liste_BLAS) 
     469 
     470        $(LK) -o ../bin/Ant-16 \ 
     471        $(Dim_ANT16) $(mod_dim_communs) \ 
     472        $(toy_recul)  \ 
     473        $(mod_communs) \ 
     474        $(mod_clim_tof)  \ 
     475        $(mod_no_tracers) \ 
     476        $(mod_ell) $(Liste_ANT16) \ 
     477        $(diagnoshelf) \ 
     478        $(Liste_Netcdf) \ 
     479        $(routines_communes) steps_time_loop.o \ 
     480        $(routine_elliptiques) $(NCDF_LIB) $(MKL_LIB) $(Liste_BLAS) 
     481 
     482Ant-16_iterbeta : $(Dim_ANT16) $(mod_dim_communs) \ 
     483        $(toy_recul) \ 
     484        $(mod_communs) \ 
     485        $(mod_clim_tof)  \ 
     486        $(mod_no_tracers) \ 
     487        $(mod_ell) $(Liste_ANT16) \ 
     488        $(diagnoshelf) \ 
     489        $(Liste_Netcdf) \ 
     490        $(routines_communes) steps_time_loop_avec_iterbeta.o \ 
     491        $(routine_elliptiques) \ 
     492        $(Liste_BLAS) 
     493 
     494        $(LK) -o ../bin/Ant-16_iterbeta \ 
     495        $(Dim_ANT16) $(mod_dim_communs) \ 
     496        $(toy_recul)  \ 
     497        $(mod_communs) \ 
     498        $(mod_clim_tof)  \ 
     499        $(mod_no_tracers) \ 
     500        $(mod_ell) $(Liste_ANT16) \ 
     501        $(diagnoshelf) \ 
     502        $(Liste_Netcdf) \ 
     503        $(routines_communes) steps_time_loop_avec_iterbeta.o \ 
     504        $(routine_elliptiques) $(NCDF_LIB) $(MKL_LIB) $(Liste_BLAS) 
     505 
    438506Ant-15 : $(Dim_ANT15-LBq) $(mod_dim_communs) \ 
    439507        $(toy_recul) \ 

  • branches/iLoveclim/SOURCES/Netcdf-routines/Description_Variables.dat

    r123 r146  
    14161416"drag_centre"                      !~  description  
    14171417------------------------------------------------------ 
     1418"acabf"                            ! nom de la variable 
     1419206      acabf     5             ! number, name in grisli, class 
     1420"o"                                ! node type : "o" ,"^", ">" 
     1421"land_ice_surface_specific_mass_balance_flux"                              !~  long name 
     1422"land_ice_surface_specific_mass_balance_flux"                              !~  standard name 
     1423"kg/m2/s"                          !~  unit  
     1424"Surface mass balance flux"                            !~  description  
     1425------------------------------------------------------ 
     1426"libmassbf"                        ! nom de la variable 
     1427207      libmassbf     5          ! number, name in grisli, class 
     1428"o"                                ! node type : "o" ,"^", ">" 
     1429"land_ice_basal_specific_mass_balance_flux"                        !~  long name 
     1430"land_ice_basal_specific_mass_balance_flux"                        !~  standard name 
     1431"kg/m2/s"                          !~  unit  
     1432"Basal mass balance flux"                        !~  description  
     1433------------------------------------------------------ 
     1434"dlithkdt"                         ! nom de la variable 
     1435208      dlithkdt     5           ! number, name in grisli, class 
     1436"o"                                ! node type : "o" ,"^", ">" 
     1437"tendency_of_land_ice_thickness"                           !~  long name 
     1438"tendency_of_land_ice_thickness"                           !~  standard name 
     1439"m/s"                              !~  unit  
     1440"Ice thickness imbalance"                         !~  description  
     1441------------------------------------------------------ 
     1442"licalvf"                          ! nom de la variable 
     1443209      licalvf     5        ! number, name in grisli, class 
     1444"o"                                ! node type : "o" ,"^", ">" 
     1445"land_ice_specific_mass_flux_due_to_calving"                       !~  long name 
     1446"land_ice_specific_mass_flux_due_to_calving"                          !~  standard name 
     1447"kg/m2/s"                          !~  unit  
     1448"Calving flux"                          !~  description  
     1449------------------------------------------------------ 
     1450"ligroundf"                        ! nom de la variable 
     1451210      ligroundf     5          ! number, name in grisli, class 
     1452"o"                                ! node type : "o" ,"^", ">" 
     1453"land_ice_specific_mass_flux_due_at_grounding_line"                        !~  long name 
     1454"land_ice_specific_mass_flux_due_at_grounding_line"                        !~  standard name 
     1455"kg/m2/s"                          !~  unit  
     1456"Grounding line flux"                        !~  description  
     1457------------------------------------------------------ 
     1458"uvelsurf"                        ! nom de la variable 
     1459211      uvelsurf     5          ! number, name in grisli, class 
     1460"o"                                ! node type : "o" ,"^", ">" 
     1461"land_ice_surface_x_velocityuvelsurf"                              !~  long name 
     1462"land_ice_surface_x_velocityuvelsurf"                        !~  standard name 
     1463"m/s"                      !~  unit  
     1464"Surface velocity in x"                        !~  description  
     1465------------------------------------------------------ 
     1466"vvelsurf"                        ! nom de la variable 
     1467212      vvelsurf     5          ! number, name in grisli, class 
     1468"o"                                ! node type : "o" ,"^", ">" 
     1469"land_ice_surface_y_velocity"                              !~  long name 
     1470"land_ice_surface_y_velocity"                        !~  standard name 
     1471"m/s"                      !~  unit  
     1472"Surface velocity in y"                        !~  description  
     1473------------------------------------------------------ 
     1474"wvelsurf"                        ! nom de la variable 
     1475213      wvelsurf     5          ! number, name in grisli, class 
     1476"o"                                ! node type : "o" ,"^", ">" 
     1477"land_ice_surface_upward_velocity"                         !~  long name 
     1478"land_ice_surface_upward_velocity"                        !~  standard name 
     1479"m/s"                      !~  unit  
     1480"Surface velocity in z"                        !~  description  
     1481------------------------------------------------------ 
     1482"uvelbase"                        ! nom de la variable 
     1483214      uvelbase     5          ! number, name in grisli, class 
     1484"o"                                ! node type : "o" ,"^", ">" 
     1485"land_ice_basal_x_velocity"                        !~  long name 
     1486"land_ice_basal_x_velocity"                        !~  standard name 
     1487"m/s"                      !~  unit  
     1488"Basal velocity in x"                        !~  description  
     1489------------------------------------------------------ 
     1490"vvelbase"                        ! nom de la variable 
     1491215      vvelbase     5          ! number, name in grisli, class 
     1492"o"                                ! node type : "o" ,"^", ">" 
     1493"land_ice_basal_y_velocity"                        !~  long name 
     1494"land_ice_basal_y_velocity"                        !~  standard name 
     1495"m/s"                      !~  unit  
     1496"Basal velocity in y"                        !~  description  
     1497------------------------------------------------------ 
     1498"wvelbase"                        ! nom de la variable 
     1499216      wvelbase     5          ! number, name in grisli, class 
     1500"o"                                ! node type : "o" ,"^", ">" 
     1501"land_ice_basal_upward_velocity"                           !~  long name 
     1502"land_ice_basal_upward_velocity"                        !~  standard name 
     1503"m/s"                      !~  unit  
     1504"Basal velocity in z"                        !~  description  
     1505------------------------------------------------------ 
     1506"strbasemag"                        ! nom de la variable 
     1507217      strbasemag     5          ! number, name in grisli, class 
     1508"o"                                ! node type : "o" ,"^", ">" 
     1509"magnitude_of_land_ice_basal_drag"                         !~  long name 
     1510"magnitude_of_land_ice_basal_drag"                        !~  standard name 
     1511"Pa"                                !~  unit  
     1512"Basal drag"                        !~  description  
     1513------------------------------------------------------ 
     1514"sftgrf"                        ! nom de la variable 
     1515218      sftgrf     5          ! number, name in grisli, class 
     1516"o"                                ! node type : "o" ,"^", ">" 
     1517"grounded_ice_sheet_area_fraction"                         !~  long name 
     1518"grounded_ice_sheet_area_fraction"                        !~  standard name 
     1519"-"                                 !~  unit  
     1520"Grounded ice sheet area fraction"                        !~  description  
     1521------------------------------------------------------ 
     1522"sftflf"                        ! nom de la variable 
     1523219      sftflf     5          ! number, name in grisli, class 
     1524"o"                                ! node type : "o" ,"^", ">" 
     1525"floating_ice_sheet_area_fraction"                         !~  long name 
     1526"floating_ice_sheet_area_fraction"                        !~  standard name 
     1527"-"                                 !~  unit  
     1528"Floating ice sheet area fraction"                        !~  description  
     1529------------------------------------------------------ 
     1530"litempbot"                        ! nom de la variable 
     1531220      litempbot     5          ! number, name in grisli, class 
     1532"o"                                ! node type : "o" ,"^", ">" 
     1533"land_ice_basal_temperature"                       !~  long name 
     1534"land_ice_basal_temperature"                        !~  standard name 
     1535"K"                                 !~  unit  
     1536"Basal temperature"                        !~  description  
     1537------------------------------------------------------ 
     1538 

  • branches/iLoveclim/SOURCES/New-remplimat/diagno-L2_mod.f90

    r123 r146  
    273273  enddo 
    274274  !$OMP END DO 
    275  
     275  !$OMP BARRIER 
     276 
     277  ! afq -- initMIP outputs: 
     278  !$OMP DO 
     279  do j=2,ny-1  
     280     do i=2,nx-1 
     281        taub(i,j) = sqrt( ((tobmx(i,j)+tobmx(i+1,j))*0.5)**2 & 
     282             + ((tobmy(i,j)+tobmy(i,j+1))*0.5)**2 ) 
     283     end do 
     284  end do 
     285  !$OMP END DO 
     286  !$OMP WORKSHARE 
     287  debug_3d(:,:,117) = taub(:,:) 
     288  !$OMP END WORKSHARE 
     289 
     290   
    276291  ! Mise ne réserve des vitesses flgzmx et flgzmy 
    277292  !$OMP WORKSHARE 

  • branches/iLoveclim/SOURCES/Temperature-routines/icetemp_mod.f90

    r123 r146  
    305305!~            nx,ny,temps,t_cpu,norme 
    306306 
    307  
     307    debug_3D(:,:,120) = T(:,:,nz)+273.15 
     308     
    308309    If (Itracebug.Eq.1) Write(num_tracebug,*)' fin routine icetemp' 
    309310    return 

  • branches/iLoveclim/SOURCES/ablation_bord.f90

    r123 r146  
    5555ablbord_dtt(:,:) = ablbord_dtt(:,:) + ablbord(:,:) 
    5656 
    57 ablbord(:,:)=0. 
     57!cdc deplace dans bilan_flux_output_mod ablbord(:,:)=0. 
    5858 
    5959end subroutine ablation_bord 

  • branches/iLoveclim/SOURCES/bilan_eau_mod.f90

    r126 r146  
    3131real,dimension(nx,ny) :: Bmelt_dtt    !< basal melting on ice points accumulated during dtt 
    3232real,dimension(nx,ny) :: calv_dtt     !< calving sur dtt (pour calcul bilan d'eau) 
     33real,dimension(nx,ny) :: archimtab    ! point pose si > 0 
     34real,dimension(nx,ny) :: grline_dtt   ! grounding line flux during dtt 
    3335 
    3436real :: sum_H_old 
    3537real :: diff_H 
     38real :: alpha_flot                      ! ro/row  
    3639 
    3740real,dimension(nx,ny) :: bm_dt,bmelt_dt 
    … …  
    4245subroutine init_bilan_eau 
    4346 ! initialisation des variables 
    44         diff_H=0. 
    45         sum_H_old = sum(H(2:nx-1,2:ny-1),mask=ice(2:nx-1,2:ny-1)==1) 
    46         tot_water(:,:)=0. 
    47         bm_dt(:,:)=0. 
    48         bmelt_dt(:,:)=0. 
     47  diff_H=0. 
     48  sum_H_old = sum(H(2:nx-1,2:ny-1),mask=ice(2:nx-1,2:ny-1)==1) 
     49  tot_water(:,:)=0. 
     50  bm_dt(:,:)=0. 
     51  bmelt_dt(:,:)=0. 
     52  alpha_flot=ro/row 
    4953! iLOVECLIM initialisation of water conservation related variables 
    50         trendWAC=0. 
    51         smbWAC(:,:)=0. 
    52         bmeltWAC(:,:)=0. 
    53         calvingWAC(:,:)=0. 
     54  trendWAC=0. 
     55  smbWAC(:,:)=0. 
     56  bmeltWAC(:,:)=0. 
     57  calvingWAC(:,:)=0. 
    5458end subroutine init_bilan_eau 
    5559 
    56          
     60 
    5761 
    5862 
    … …  
    8892 
    8993where (ice(2:nx-1,2:ny-1).eq.1) 
    90         Bm_dtt(2:nx-1,2:ny-1) = Bm_dtt(2:nx-1,2:ny-1) + Bm_dt(2:nx-1,2:ny-1) !* dt ! somme Bm sur dt 
    91         bmelt_dtt(2:nx-1,2:ny-1) = bmelt_dtt(2:nx-1,2:ny-1) + bmelt_dt(2:nx-1,2:ny-1) ! * dt ! somme bmelt sur dt 
     94   Bm_dtt(2:nx-1,2:ny-1) = Bm_dtt(2:nx-1,2:ny-1) + Bm_dt(2:nx-1,2:ny-1) !* dt ! somme Bm sur dt 
     95   bmelt_dtt(2:nx-1,2:ny-1) = bmelt_dtt(2:nx-1,2:ny-1) + bmelt_dt(2:nx-1,2:ny-1) ! * dt ! somme bmelt sur dt 
    9296endwhere 
     97 
     98archimtab(:,:) = Bsoc(:,:)+H(:,:)*alpha_flot - sealevel 
     99gr_line(:,:)=0 
     100do j=1,ny 
     101  do i=1,nx 
     102     !afq    if ((H(i,j).gt.0.).and.(archimtab(i,j).GE.0.).and.(Bsoc(i,j).LE.sealevel)) then ! grounded with ice 
     103     if ((H(i,j).gt.0.).and.(archimtab(i,j).GE.0.)) then ! grounded with ice 
     104      if (archimtab(i-1,j).LT.0..and.Uxbar(i,j).LT.0..and..not.flot_marais(i-1,j)) gr_line(i,j)=1 
     105      if (archimtab(i+1,j).LT.0..and.Uxbar(i+1,j).GT.0..and..not.flot_marais(i+1,j)) gr_line(i,j)=1 
     106      if (archimtab(i,j-1).LT.0..and.Uybar(i,j).LT.0..and..not.flot_marais(i,j-1)) gr_line(i,j)=1 
     107      if (archimtab(i,j+1).LT.0..and.Uybar(i,j+1).GT.0..and..not.flot_marais(i,j+1)) gr_line(i,j)=1 
     108    endif 
     109  enddo 
     110enddo 
     111 
     112where (gr_line(:,:).eq.1) 
     113!~   grline_dtt(:,:)= (((uxbar(:,:)+eoshift(uxbar(:,:),shift=1,boundary=0.,dim=1))**2+ & 
     114!~                     (uybar(:,:)+eoshift(uybar(:,:),shift=1,boundary=0.,dim=2))**2)**0.5)*0.5 & 
     115!~                     *H(:,:) + grline_dtt(:,:)  
     116   grline_dtt(:,:)= - sqrt(                                                                        & 
     117                           ( (uxbar(:,:)+eoshift(uxbar(:,:),shift=1,boundary=0.,dim=1))*dy/2. )**2 +  & 
     118                           ( (uybar(:,:)+eoshift(uybar(:,:),shift=1,boundary=0.,dim=2))*dx/2. )**2 )  & 
     119                           * H(:,:) * dt + grline_dtt(:,:) 
     120endwhere 
     121 
    93122 
    94123if (isynchro.eq.1) then 
    … …  
    103132    
    104133! bilan d'eau sur la grille : 
    105         water_bilan=sum(tot_water(:,:)) 
    106         diff_H = diff_H/dtt 
     134  water_bilan=sum(tot_water(:,:)) 
     135  diff_H = diff_H/dtt 
    107136 
    108137!999 format(f0.2,1x,e15.8,1x,i10,8(1x,e15.8)) 
    109138!       write(6,999),time,sum_H,count(ice(:,:)==1),diff_H,water_bilan,sum(calv_dtt(:,:))/dtt,sum(ablbord_dtt(:,:))/dtt,sum(bmelt_dtt(:,:),mask=ice(:,:)==1)/dtt,sum(bm(:,:),mask=ice(:,:)==1),sum(Bm_dtt(:,:))/dtt,sum(bmelt_dtt(:,:))/dtt 
    110         diff_H_water_bilan(2:nx-1,2:ny-1)=tot_water(2:nx-1,2:ny-1)-diff_H_2D(2:nx-1,2:ny-1) 
     139  diff_H_water_bilan(2:nx-1,2:ny-1)=tot_water(2:nx-1,2:ny-1)-diff_H_2D(2:nx-1,2:ny-1) 
    111140 
    112141endif 

  • branches/iLoveclim/SOURCES/bmelt-seuil-profondeur_mod.f90

    r77 r146  
    1616module  bmelt_seuil_prof          
    1717 
    18 ! prametrise la fusion basale (ice shelves) 
    19 ! Pour l'actuel : 2 valeurs  pour 2 domaines de profondeur 
    20 ! + une valeur pour bmgrz 
    21 ! A choisir dans le module_choix 
     18  ! prametrise la fusion basale (ice shelves) 
     19  ! Pour l'actuel : 2 valeurs  pour 2 domaines de profondeur 
     20  ! + une valeur pour bmgrz 
     21  ! A choisir dans le module_choix 
    2222 
    23 use module3d_phy 
     23  use module3d_phy 
    2424 
    2525 
    26 implicit none 
     26  implicit none 
    2727 
    28 real :: bm_grz                        !< valeur prescrite a la grounding zone 
    29 real,dimension(nx,ny) ::  bmgrz       !< tabelau fusion basale a la grounding zone 
     28  real :: bm_grz                        !< valeur prescrite a la grounding zone 
     29  real,dimension(nx,ny) ::  bmgrz       !< tabelau fusion basale a la grounding zone 
    3030 
    31 real :: bmshelf_plateau               !< valeur prescrite sur le plateau cont. 
    32 real :: bmshelf_abysses               !< valeur prescrite au dessus de l'ocean profond 
    33 real :: depth_talus                   !< profondeur de transition  
    34 real,dimension(nx,ny) ::  bmshelf     !< tableau fusion basale sous shelf 
     31  real :: bmshelf_plateau               !< valeur prescrite sur le plateau cont. 
     32  real :: bmshelf_abysses               !< valeur prescrite au dessus de l'ocean profond 
     33  real :: depth_talus                   !< profondeur de transition  
     34  real,dimension(nx,ny) ::  bmshelf     !< tableau fusion basale sous shelf 
    3535 
    3636 
    3737contains 
    38 !------------------------------------------------------------------------------- 
    39 !> SUBROUTINE: init_bmelt 
    40 !! Cette routine fait l'initialisation pour la fusion basale. 
    41 !! @note Elle est appelée par inputfile-vec-0.5.f90 
    42 !! 
    43 !> 
    44 subroutine init_bmelt 
     38  !------------------------------------------------------------------------------- 
     39  !> SUBROUTINE: init_bmelt 
     40  !! Cette routine fait l'initialisation pour la fusion basale. 
     41  !! @note Elle est appelée par inputfile-vec-0.5.f90 
     42  !! 
     43  !> 
     44  subroutine init_bmelt 
    4545 
    4646 
    4747 
    48 ! Cette routine fait l'initialisation pour la fusion basale. 
    49 ! Elle est appelée par inputfile-vec-0.5.f90 
     48    ! Cette routine fait l'initialisation pour la fusion basale. 
     49    ! Elle est appelée par inputfile-vec-0.5.f90 
    5050 
    5151 
    52 namelist/bmelt_seuil/bm_grz,bmshelf_plateau,bmshelf_abysses,depth_talus 
     52    namelist/bmelt_seuil/bm_grz,bmshelf_plateau,bmshelf_abysses,depth_talus 
    5353 
    54 if (itracebug.eq.1)  call tracebug('entree dans init_bmelt de bmelt_seuil_prof') 
     54    if (itracebug.eq.1)  call tracebug('entree dans init_bmelt de bmelt_seuil_prof') 
    5555 
    56 rewind(num_param)        ! pour revenir au debut du fichier param_list.dat 
    57 read(num_param,bmelt_seuil) 
    58 write(num_rep_42,bmelt_seuil) 
     56    rewind(num_param)        ! pour revenir au debut du fichier param_list.dat 
     57    read(num_param,bmelt_seuil) 
     58    write(num_rep_42,bmelt_seuil) 
    5959 
    60 !    ecriture dans le fichier parametres 
     60    !    ecriture dans le fichier parametres 
    6161 
    6262428 format(A) 
    6363 
    64 write(num_rep_42,428)'!___________________________________________________________'  
    65 write(num_rep_42,428) '&bmelt_seuil                     ! module  bmelt_seuil_prof' 
    66 write(num_rep_42,*)   'bm_grz           =',bm_grz 
    67 write(num_rep_42,*)   'bmshelf_plateau  =',bmshelf_plateau 
    68 write(num_rep_42,*)   'bmshelf_abysses  =',bmshelf_abysses 
    69 write(num_rep_42,*)   'depth_talus      =',depth_talus 
    70 write(num_rep_42,*)'/'                       
    71 write(num_rep_42,428) '! Pour l actuel : bm_grz a la grounding line' 
    72 write(num_rep_42,428) '!                 bmshelf_plateau sur le plateau continental' 
    73 write(num_rep_42,428) '!                 bmshelf_abysses pour les grandes profondeurs' 
    74 write(num_rep_42,428) '!       depth_talus, negative, separation entre les 2 domaines' 
    75 write(num_rep_42,*) 
    76 write(num_rep_42,428)'!___________________________________________________________'  
     64    write(num_rep_42,428)'!___________________________________________________________'  
     65    write(num_rep_42,428) '&bmelt_seuil                     ! module  bmelt_seuil_prof' 
     66    write(num_rep_42,*)   'bm_grz           =',bm_grz 
     67    write(num_rep_42,*)   'bmshelf_plateau  =',bmshelf_plateau 
     68    write(num_rep_42,*)   'bmshelf_abysses  =',bmshelf_abysses 
     69    write(num_rep_42,*)   'depth_talus      =',depth_talus 
     70    write(num_rep_42,*)'/'                       
     71    write(num_rep_42,428) '! Pour l actuel : bm_grz a la grounding line' 
     72    write(num_rep_42,428) '!                 bmshelf_plateau sur le plateau continental' 
     73    write(num_rep_42,428) '!                 bmshelf_abysses pour les grandes profondeurs' 
     74    write(num_rep_42,428) '!       depth_talus, negative, separation entre les 2 domaines' 
     75    write(num_rep_42,*) 
     76    write(num_rep_42,428)'!___________________________________________________________'  
    7777 
    7878 
    7979 
    80 bmgrz(:,:) = bm_grz 
     80    bmgrz(:,:) = bm_grz 
    8181 
    82 where (Bsoc0(:,:).lt.depth_talus) 
    83    bmshelf(:,:)=bmshelf_abysses 
    84 elsewhere 
    85    bmshelf(:,:)=bmshelf_plateau 
    86 end where 
     82    where (Bsoc0(:,:).lt.depth_talus) 
     83        bmshelf(:,:)=bmshelf_abysses 
     84    elsewhere 
     85        bmshelf(:,:)=bmshelf_plateau 
     86    end where 
    8787 
    88 debug_3D(:,:,34)=bmshelf(:,:) 
    89 return 
    90 end subroutine init_bmelt 
     88    debug_3D(:,:,34)=bmshelf(:,:) 
     89    return 
     90  end subroutine init_bmelt 
    9191 
    92 !________________________________________________________________________________ 
     92  !________________________________________________________________________________ 
    9393 
    94 !> SUBROUTINE: bmeltshelf 
    95 !! Cette routine calcule la fusion basale proprement dite pour les points flottants 
    96 !! @note coefbmshelf a ete calcule par forclim et permet de faire varier en fonction du climat 
    97 !> 
     94  !> SUBROUTINE: bmeltshelf 
     95  !! Cette routine calcule la fusion basale proprement dite pour les points flottants 
     96  !! @note coefbmshelf a ete calcule par forclim et permet de faire varier en fonction du climat 
     97  !> 
    9898 
    99 subroutine bmeltshelf 
     99  subroutine bmeltshelf 
    100100 
    101101 
    102 ! cette routine calcule la fusion basale proprement dite pour les points flottants 
    103 ! coefbmshelf a ete calcule par forclim et permet de faire varier en fonction du climat 
     102    ! cette routine calcule la fusion basale proprement dite pour les points flottants 
     103    ! coefbmshelf a ete calcule par forclim et permet de faire varier en fonction du climat 
    104104 
    105 integer :: ngr           ! nombre de voisins flottants 
    106 real    :: coef_talus    ! pour ne pas changer la fusion au dessus de l'ocean profond 
     105    integer :: ngr           ! nombre de voisins flottants 
     106    real    :: coef_talus    ! pour ne pas changer la fusion au dessus de l'ocean profond 
    107107 
    108 if (itracebug.eq.1)  call tracebug('entree dans bmeltshelf de bmelt_seuil_prof') 
     108    if (itracebug.eq.1)  call tracebug('entree dans bmeltshelf de bmelt_seuil_prof') 
    109109 
    110 where (Bsoc0(:,:).lt.depth_talus) 
    111    bmshelf(:,:)=bmshelf_abysses 
    112 elsewhere 
    113    bmshelf(:,:)=bmshelf_plateau 
    114 end where 
     110    where (Bsoc0(:,:).lt.depth_talus) 
     111        bmshelf(:,:)=bmshelf_abysses 
     112    elsewhere 
     113        bmshelf(:,:)=bmshelf_plateau 
     114    end where 
    115115 
    116 do j=1,ny 
    117    do i=1,nx 
     116    do j=1,ny 
     117      do i=1,nx 
    118118 
    119 talus_nochange : if (Bsoc0(i,j).lt.depth_talus) then  
    120            coef_talus = 1         ! pas de changement au dessus de l'ocecan profond 
     119        talus_nochange: if (Bsoc0(i,j).lt.depth_talus) then  
     120            coef_talus = 1         ! pas de changement au dessus de l'ocecan profond 
    121121        else 
    122            coef_talus = coefbmshelf 
     122            coef_talus = coefbmshelf 
    123123        endif talus_nochange 
    124124 
    125125 
    126 shelf:    if (flot(i,j)) then    ! partie flottante 
     126        shelf: if (flot(i,j)) then    ! partie flottante 
    127127 
    128128            bmelt(i,j)=coef_talus*bmshelf(i,j) 
    129129 
    130 ! fbm est vrai si le point est flottant mais un des voisins est pose 
    131 ! calcule dans flottab 
     130            ! fbm est vrai si le point est flottant mais un des voisins est pose 
     131            ! calcule dans flottab 
    132132 
    133133            if (fbm(i,j)) then  
    134                bmelt(i,j)=coef_talus*bmgrz(i,j) 
    135             endif    
     134                bmelt(i,j)=coef_talus*bmgrz(i,j) 
     135            endif 
    136136 
    137137 
    138138 
    139 ! ATTENTION le bloc suivant sert a determiner la fusion basale d'equilibre 
    140 ! pour les shelves stationnaires 
     139            ! ATTENTION le bloc suivant sert a determiner la fusion basale d'equilibre 
     140            ! pour les shelves stationnaires 
    141141 
    142   if (igrdline.eq.1) then 
    143      corrbmelt(i,j)=hdot(i,j)+bmelt(i,j)   ! le bmelt d'equilibre 
    144      debug_3D(i,j,28)=corrbmelt(i,j) 
    145   endif 
     142            if ((igrdline.eq.1).and.(ibmelt_inv.eq.1)) then 
     143                !     corrbmelt(i,j)=hdot(i,j)+bmelt(i,j)   ! le bmelt d'equilibre 
     144                !     bmelt(i,j)=corrbmelt(i,j) 
     145                corrbmelt(i,j)=corrbmelt(i,j)+hdot(i,j)*0.85 
     146                bmelt(i,j)=bmelt(i,j)+corrbmelt(i,j)      
     147            endif 
    146148 
    147149 
    148150        else                   ! point posé, on compte le nombre de voisins flottants 
    149            ngr=0 
    150            if ((i.ne.1).and.(i.ne.nx).and.(j.ne.1).and.(j.ne.ny)) then  
    151               if (flot(i+1,j)) ngr=ngr+1 
    152               if (flot(i-1,j)) ngr=ngr+1 
    153               if (flot(i,j+1)) ngr=ngr+1 
    154               if (flot(i,j-1)) ngr=ngr+1 
    155            end if 
    156                
     151            ngr=0 
     152            if ((i.ne.1).and.(i.ne.nx).and.(j.ne.1).and.(j.ne.ny)) then  
     153                if (flot(i+1,j)) ngr=ngr+1 
     154                if (flot(i-1,j)) ngr=ngr+1 
     155                if (flot(i,j+1)) ngr=ngr+1 
     156                if (flot(i,j-1)) ngr=ngr+1 
     157            end if 
    157158 
    158 !   la fusion des points limites est une combinaison entre valeur posée et valeur flottante 
    159 !   en fonction du nombre de points flottants 
    160159 
    161            bmelt(i,j)= ngr/4.*bmgrz(i,j)*coef_talus+(1.-ngr/4.)*bmelt(i,j) 
    162           
     160            !   la fusion des points limites est une combinaison entre valeur posée et valeur flottante 
     161            !   en fonction du nombre de points flottants 
    163162 
     163            bmelt(i,j)= ngr/4.*bmgrz(i,j)*coef_talus+(1.-ngr/4.)*bmelt(i,j) 
     164 
     165            if ((igrdline.eq.1).and.(ibmelt_inv.eq.1)) then 
     166                corrbmelt(i,j)=corrbmelt(i,j)+hdot(i,j)  !*0.85 
     167                bmelt(i,j)=bmelt(i,j)+corrbmelt(i,j) 
     168            endif 
    164169        endif shelf 
    165  
    166       end do  
     170      end do 
    167171    end do 
    168172 
    169  
    170 if (igrdline.eq.1) then 
    171     bmelt(:,:)=0.      ! hdot donne alors le -bmelt 
     173if ((igrdline.eq.1).and.(ibmelt_inv.eq.0)) then 
     174    where (i_Hp(:,:).eq.1) 
     175        bmelt(:,:)=0.    ! hdot donne alors le -bmelt 
     176    endwhere 
    172177endif 
    173  
    174178 
    175179return 

  • branches/iLoveclim/SOURCES/calving_frange.f90

    r126 r146  
    122122  mass_total=0d0 
    123123!afq for CONSEAU  calvtest=0. 
    124   mass_total=0d0 
    125124 
    126125! calcul du dhdt lagrangien 
    … …  
    159158 
    160159!ifint:    if((front(i,j).gt.0).and.(front(i,j).lt.4)) then 
    161         ifint:    if(front(i,j).lt.4) then 
     160!cdc pb avec front        ifint:    if(front(i,j).lt.4) then 
     161        ifint:    if((H(i-1,j).lt.2.).or.(H(i+1,j).lt.2.).or.(H(i,j-1).lt.2.).or.(H(i,j+1).lt.2)) then ! si on est au bord avec test sur H 
    162162! ifint: le point doit avoir au - 1 voisin mais ne pas etre entouré de glace 
    163163! ce qui evite la formation des polynies dans les shelfs 
    … …  
    491491                     S(i,j)=H(i,j)*(1.-ro/row) + sealevel 
    492492                     B(i,j)=S(i,j) - H(i,j) 
     493          ! ATTENTION ne pas mettre ice=0 sinon degradation bilan d'eau (bm et bmelt non comptabilises dans ce cas) 
    493494                   endif 
    494495                   !H(i,j)=1              

  • branches/iLoveclim/SOURCES/conserv-mass-adv-diff_sept2009_mod.f90

    r123 r146  
    416416 
    417417!    if ((time.lt.time_gl(1)).or.(nt.lt.2)) then      ATTENTION test seulement si version prescribe-H_mod.f90 
    418        call prescribe_present_H_gl       ! prescribe present thickness at the grounding line 
     418   if (ibmelt_inv.eq.0) then 
     419      call prescribe_present_H_gl       ! prescribe present thickness at the grounding line 
     420   else if (ibmelt_inv.eq.1) then ! inversion bmelt 
     421      call prescribe_present_H_gl_bmelt ! prescribe only grounding line points 
     422   else 
     423      print*,'ibmelt_inv value must be 0 or 1 !!!' 
     424      stop 
     425   endif 
    419426!    else 
    420427!       call prescribe_retreat 
    … …  
    550557H(:,:)=max(H(:,:),0.)       ! pas d'epaisseur negative 
    551558 
    552 !~ ! calcul du masque "ice" 
    553 !~ where (flot(:,:))        ! points flottants, sera éventuellement réévalué dans flottab 
    554 !~    where(H(:,:).gt.max(Hmin,Hmin+BM(:,:)-Bmelt(:,:))) 
    555 !~       ice(:,:)=1 
    556 !~    elsewhere 
    557 !~       ice(:,:)=0 
    558 !~       H(:,:)=max(1.,min(0.,(sealevel - Bsoc(:,:))*row/ro-0.01)) 
    559 !~    endwhere 
    560 !~ !   H(:,:)=max(H(:,:),1.) ! dans la partie marine l'épaisseur mini est 1 m 
    561 !~ !   H(:,:)=max(1.,min(0.,(sealevel - Bsoc(:,:))*row/ro-0.01)) 
    562 !~ elsewhere                ! points posés 
    563 !~    where(H(:,:).gt.0.) 
    564 !~       ice(:,:)=1 
    565 !~    elsewhere 
    566 !~       ice(:,:)=0 
    567 !~    endwhere 
    568 !~ endwhere 
    569  
    570559! eventuellement retirer apres spinup ou avoir un cas serac 
    571560Hdot(:,:)=(H(:,:)-vieuxH(:,:))/dt 
    572561!$OMP END WORKSHARE 
    573562 
    574 !~ if (igrdline.ne.3) then 
    575 !~    !$OMP WORKSHARE 
    576 !~    Hdot(:,:)=min(Hdot(:,:),10.) 
    577 !~    Hdot(:,:)=max(Hdot(:,:),-10.) 
    578 !~    !$OMP END WORKSHARE 
    579 !~ endif 
     563if (ibmelt_inv.eq.1) then ! inversion du bmelt avec igrdline=1 
     564   !$OMP WORKSHARE 
     565   where (i_Hp(:,:).eq.1) 
     566      H(:,:)=Hp(:,:) 
     567   endwhere 
     568   !$OMP END WORKSHARE 
     569endif 
    580570 
    581571!$OMP WORKSHARE 
    582 where (i_hp(:,:).ne.1) 
    583    H(:,:)=vieuxH(:,:)+Hdot(:,:)*dt 
    584 end where 
    585  
    586  
    587  
    588572! si mk0=-1, on garde l'epaisseur precedente 
    589573! en attendant un masque plus propre 

  • branches/iLoveclim/SOURCES/dragging_param_beta_mod.f90

    r123 r146  
    3636real :: betamin   ! betamin : (Pa) frottement mini sous les streams 
    3737 
    38 real :: beta_slope       ! = A in: logB = A log Neff +K 
    39 real :: beta_intercept   ! = K in: logB = A log Neff +K 
     38real :: beta_slope        ! = A in:     B = A x Neff ** K 
     39real :: beta_expo         ! = K in:     B = A x Neff ** K 
    4040 
    4141real,dimension(nx,ny) :: neff   ! pression effective noeuds majeurs 
    … …  
    6161implicit none 
    6262 
    63 ! local variables, defining the rugosity-enhanced flow 
    64 real :: expo_slope 
    65 real :: pente_min, pente_max 
    66  
    67 namelist/drag_param_beta/beta_intercept,beta_slope,betamax,betamin,coef_ile 
     63namelist/drag_param_beta/beta_slope,beta_expo,betamax,betamin,coef_ile 
    6864 
    6965if (itracebug.eq.1)  call tracebug(' dragging avec neff et slope') 
    … …  
    8278write(num_rep_42,428) '&drag_param_beta          ! nom du bloc dragging param beta' 
    8379write(num_rep_42,*) 
    84 write(num_rep_42,*) 'beta_intercept        = ', beta_intercept 
    8580write(num_rep_42,*) 'beta_slope            = ', beta_slope 
     81write(num_rep_42,*) 'beta_expo             = ', beta_expo 
    8682write(num_rep_42,*) 'betamax               = ', betamax 
    8783write(num_rep_42,*) 'betamin               = ', betamin 
    8884write(num_rep_42,*)'/'                             
    89 write(num_rep_42,428) '! Intercept & slope of linear reg on loglog' 
     85write(num_rep_42,428) '! Slope & expo of beta = - slope x Neff ** expo' 
    9086 
    9187!------------------------------------------------------------------- 
    … …  
    181177!$OMP END DO 
    182178 
    183  
    184179!$OMP WORKSHARE 
    185180 
    186181! new parametrisation of beta on Neff: 
    187 !betamx(:,:)=exp(beta_intercept * log(10.) )*neffmx(:,:)**(beta_slope) 
    188 !betamy(:,:)=exp(beta_intercept * log(10.) )*neffmy(:,:)**(beta_slope) 
    189 betamx(:,:)= (10.**beta_intercept)*(neffmx(:,:)**beta_slope) 
    190 betamy(:,:)= (10.**beta_intercept)*(neffmy(:,:)**beta_slope) 
    191  
     182betamx(:,:)= beta_slope*(neffmx(:,:)**beta_expo) 
     183betamy(:,:)= beta_slope*(neffmy(:,:)**beta_expo) 
    192184 
    193185where (ilemx(:,:)) betamx(:,:) = betamx(:,:) * coef_ile 
    … …  
    203195where ( hwatmy(:,:) .lt. 0.5 ) betamy(:,:) = betamax 
    204196 
    205 ! aurel, we add the neff threshold: 
    206 !where ((ibase(:,:).ne.1).and.(.not.flot(:,:)).and.(H(:,:).gt.1.)) fleuve(:,:)=.true. 
    207 !where ( (betamx(:,:).lt.betamax) .or. (betamy(:,:).lt. betamax) ) fleuve(:,:)=.true. 
    208197!$OMP END WORKSHARE 
    209198 

  • branches/iLoveclim/SOURCES/flottab2-0.7.f90

    r123 r146  
    6868 integer ::  numtache 
    6969 integer ::  nb_pt 
     70 real    ::  petit_H=0.001 ! pour test ice sur zone flottante 
    7071contains 
    7172! ----------------------------------------------------------------------------------- 
    … …  
    593594    where (flot(:,:)) 
    594595!cdc 1m       where (H(:,:).gt.max(Hmin,Hmin+BM(:,:)-Bmelt(:,:))) 
    595                          where (H(:,:).gt.0.) 
     596      where (H(:,:).gt.max(BM(:,:)-Bmelt(:,:)+petit_H,0.)*dt) 
     597!cdc                     where (H(:,:).gt.0.) 
    596598          ice(:,:)=1 
    597599       elsewhere  
    … …  
    745747 
    746748    call determin_marais 
     749 
     750    ! pour sorties initMIP: 
     751    debug_3D(:,:,118) = ice(:,:)*(1-mk(:,:)) 
     752    debug_3D(:,:,119) = ice(:,:)*mk(:,:) 
     753 
    747754     
    748755  end subroutine flottab 
    … …  
    12701277end do 
    12711278!$OMP END DO  
    1272  
     1279!$OMP BARRIER 
    12731280! traitement des coins 
    12741281 
    … …  
    12791286 
    12801287! Les points à plus d'un point de grille du bord sont front=-1 
     1288!$OMP WORKSHARE 
    12811289nofront(:,:)=0 
     1290!$OMP END WORKSHARE 
     1291!$OMP BARRIER 
    12821292 
    12831293!$OMP DO 
    … …  
    12901300enddo 
    12911301!$OMP END DO 
     1302!$OMP BARRIER 
     1303 
     1304!$OMP WORKSHARE 
    12921305where (nofront(:,:).eq.-1) 
    12931306        front(:,:)=-1 
    12941307endwhere 
    1295  
     1308!$OMP END WORKSHARE 
    12961309!$OMP END PARALLEL               
    12971310 

  • branches/iLoveclim/SOURCES/initial-0.3.f90

    r123 r146  
    3030  use sorties_ncdf_grisli 
    3131  use util_recovery 
     32  use bilan_flux_mod 
    3233 
    3334  !------------------------------------------------------------------------------------- 
    … …  
    137138  call init_calving 
    138139 
     140  call init_bilan_flux 
    139141! 
    140142!------------------------------------------------------------------------------------- 

  • branches/iLoveclim/SOURCES/initial-phy-2.f90

    r123 r146  
    2727 
    2828  namelist/runpar/runname,icompteur,iout,reprcptr,itracebug,num_tracebug,comment_run 
    29   namelist/grdline/igrdline,Schoof 
     29  namelist/grdline/igrdline,Schoof,ibmelt_inv 
    3030  namelist/timesteps/dtmin,dtmax,dtt,testdiag,tbegin,tend 
    3131 
    … …  
    127127  write(num_rep_42,*) 'igrdline     = ',igrdline 
    128128  write(num_rep_42,*) 'Schoof       = ',Schoof 
     129  write(num_rep_42,*) 'ibmelt_inv   = ',ibmelt_inv 
    129130  write(num_rep_42,*)'/'                            
    130   write(num_rep_42,428)'! igrdline :  1 ligne d echouage fixée, sinon 0' 
    131   write(num_rep_42,428)'! Schoof   :  0 pas de Schoof, 1 flux de Schoof' 
     131  write(num_rep_42,428)'! igrdline   : 1 ligne d echouage fixée, sinon 0' 
     132  write(num_rep_42,428)'! Schoof     : 0 pas de Schoof, 1 flux de Schoof' 
     133  write(num_rep_42,428)'! ibmelt_inv : 0 cas std, 1 inversion du bmelt (avec igrdline=1)' 
    132134  write(num_rep_42,*) 
    133135 

  • branches/iLoveclim/SOURCES/main3D-0.4-40km.f90

    r126 r146  
    237237 
    238238 
    239  
    240   if (iter_beta.eq.0) then 
    241  
    242      if (itracebug.eq.1)  call tracebug(' Avant appel routine icethick3') 
    243      call icethick3 
    244      debug_3D(:,:,88) = S(:,:) 
    245      if (itracebug.eq.1)  call tracebug(' Apres appel routine icethick3') 
    246   end if 
     239!cdc supprime pour initialisation propre 
     240!~   if (iter_beta.eq.0) then 
     241 
     242!~      if (itracebug.eq.1)  call tracebug(' Avant appel routine icethick3') 
     243!~      call icethick3 
     244!~      debug_3D(:,:,88) = S(:,:) 
     245!~      if (itracebug.eq.1)  call tracebug(' Apres appel routine icethick3') 
     246!~   end if 
    247247 
    248248 
    … …  
    305305           end do 
    306306 
    307  
    308307           call Neffect() 
    309308           call flottab() 
    310  
    311  
     309           call calving 
     310           call ablation_bord 
     311           call flottab 
    312312           call Neffect() 
    313  
    314  
    315313           call diffusiv() 
    316314           call SIA_velocities() 
    317  
    318315        endif 
    319316 

  • branches/iLoveclim/SOURCES/mix-SIA-L1_mod.f90

    r77 r146  
    171171endif 
    172172 
     173 
     174! afq -- initMIP outputs: 
     175 
     176  !$OMP WORKSHARE 
     177  debug_3d(:,:,111) = ( (ux  (:,:,1)  + eoshift( ux(:,:,1) ,shift=1,boundary=0.,dim=1))/2.  )/secyear 
     178  debug_3d(:,:,112) = ( (uy  (:,:,1)  + eoshift( uy(:,:,1) ,shift=1,boundary=0.,dim=2))/2.  )/secyear 
     179  debug_3d(:,:,113) = ( uzr (:,:,1) * ice(:,:)                                              )/secyear 
     180  debug_3d(:,:,114) = ( (ux  (:,:,nz) + eoshift( ux(:,:,nz),shift=1,boundary=0.,dim=1))/2.  )/secyear 
     181  debug_3d(:,:,115) = ( (uy  (:,:,nz) + eoshift( uy(:,:,nz),shift=1,boundary=0.,dim=2))/2.  )/secyear 
     182  debug_3d(:,:,116) = ( uzr (:,:,nz) * ice(:,:)                                             )/secyear 
     183  !$OMP END WORKSHARE 
     184 
     185 
    173186end subroutine mix_SIA_L1 
    174187 

  • branches/iLoveclim/SOURCES/out_cptr_mod.f90

    r77 r146  
    312312     Uyflgz(:,:) = Uy(:,:,nz) 
    313313 
     314     debug_3D(:,:,120) = T(:,:,nz) + 273.15 
    314315 
    315316 

  • branches/iLoveclim/SOURCES/prescribe-H-i2s_mod.f90

    r4 r146  
    127127       hp(:,:)   = Hp0(:,:) 
    128128    end where 
    129 if (itracebug.eq.1)  call tracebug(' fin prescribe_present_H_gl ') 
     129if (itracebug.eq.1)  call tracebug(' fin prescribe_present_H_gl') 
    130130 
    131131 
    132132  end subroutine prescribe_present_H_gl 
     133   
     134   !______________________________________________________________________________________ 
     135  !>  function prescribe_present_H_gl_bmelt 
     136  !! calculate the initial grounding line position 
     137  !! only grounding line points are prescribed to present value 
     138  !! @return i_hp(:,:) and hp(:,:) 
     139 
     140  subroutine prescribe_present_H_gl_bmelt 
     141 
     142    implicit none 
     143 
     144    if (itracebug.eq.1)  call tracebug(' Entree dans routine  prescribe_present_H_gl_bmelt') 
     145 
     146!    where (MK_gl0(:,:) .eq. 1) ! grounding line only !cdc pour calcule fonte basale 
     147!       i_hp(:,:) = 1                                 ! thickness prescribed to present value 
     148!       hp(:,:)   = Hp0(:,:) 
     149!    end where 
     150    i_hp(:,:) = 0 
     151if (itracebug.eq.1)  call tracebug(' fin prescribe_present_H_gl_bmelt') 
     152 
     153 
     154  end subroutine prescribe_present_H_gl_bmelt 
    133155 
    134156  !______________________________________________________________________________________ 

  • branches/iLoveclim/SOURCES/resol_adv_diff_2D-sept2009.f90

    r77 r146  
    278278 
    279279relax_loop: do while(.not.stopp) 
    280 ntour=ntour+1 
    281 !$OMP PARALLEL 
    282 !$OMP DO PRIVATE(reste) 
    283 do j=2,ny-1 
    284    do i=2,nx-1 
    285  
    286       reste = (((arelax(i,j)*newH(i-1,j) +drelax(i,j)*newH(i,j-1)) & 
    287            + (brelax(i,j)*newH(i+1,j) + erelax(i,j)*newH(i,j+1))) & 
    288            + crelax(i,j)*newH(i,j))- frelax(i,j) 
    289  
    290 !if (ntour.eq.1)  debug_3D(i,j,49)=(((arelax(i,j)*newH(i-1,j) +drelax(i,j)*newH(i,j-1)) & 
    291 !           + (brelax(i,j)*newH(i+1,j) + erelax(i,j)*newH(i,j+1))) & 
    292 !           + crelax(i,j)*newH(i,j)) 
    293  
    294  
    295       deltaH(i,j) = reste/crelax(i,j) 
    296  
     280   ntour=ntour+1 
     281   !$OMP PARALLEL 
     282   !$OMP DO PRIVATE(reste) 
     283   do j=2,ny-1 
     284      do i=2,nx-1 
     285         reste = (((arelax(i,j)*newH(i-1,j) +drelax(i,j)*newH(i,j-1)) & 
     286               + (brelax(i,j)*newH(i+1,j) + erelax(i,j)*newH(i,j+1))) & 
     287               + crelax(i,j)*newH(i,j))- frelax(i,j) 
     288 
     289         deltaH(i,j) = reste/crelax(i,j) 
     290      end do 
    297291   end do 
    298 end do 
    299 !$OMP END DO 
    300  
    301 !debug_3D(:,:,50)=arelax(:,:) 
    302 !debug_3D(:,:,51)=brelax(:,:) 
    303 !debug_3D(:,:,52)=crelax(:,:) 
    304 !debug_3D(:,:,53)=drelax(:,:) 
    305 !debug_3D(:,:,54)=erelax(:,:) 
    306 !debug_3D(:,:,55)=frelax(:,:) 
    307  
    308  
    309 !$OMP WORKSHARE 
    310 newH(:,:)=newH(:,:)-deltaH(:,:) 
    311 !$OMP END WORKSHARE 
    312 !$OMP END PARALLEL 
     292   !$OMP END DO 
     293 
     294   !$OMP WORKSHARE 
     295   newH(:,:)=newH(:,:)-deltaH(:,:) 
     296   !$OMP END WORKSHARE 
     297   !$OMP END PARALLEL 
    313298 
    314299 
    315300! critere d'arret: 
    316301! ----------------          
    317  
    318 delh=0 
    319  
    320 !$OMP PARALLEL 
    321 !$OMP DO REDUCTION(+:delh) 
    322 do j=2,ny-1 
    323    do i=2,nx-1 
    324       delh=delh+deltaH(i,j)**2 
     302   delh=0 
     303 
     304   !$OMP PARALLEL 
     305   !$OMP DO REDUCTION(+:delh) 
     306   do j=2,ny-1 
     307      do i=2,nx-1 
     308         delh=delh+deltaH(i,j)**2 
     309      end do 
    325310   end do 
    326 end do 
    327 !$OMP END DO 
    328 !$OMP END PARALLEL 
    329  
    330 if (delh.gt.0.) then 
    331    testh=sqrt(delh)/((nx-2)*(ny-2)) 
    332 else 
    333    testh=0. 
    334 endif 
    335 stopp = (testh.lt.1.e-4).or.(ntour.gt.100) 
    336  
    337  
     311   !$OMP END DO 
     312   !$OMP END PARALLEL 
     313 
     314   if (delh.gt.0.) then 
     315      testh=sqrt(delh)/((nx-2)*(ny-2)) 
     316   else 
     317      testh=0. 
     318   endif 
     319   stopp = (testh.lt.1.e-4).or.(ntour.gt.100) 
    338320 
    339321end do relax_loop 
    340322 
    341  
    342 ! thickness at the upwind node 
    343 !do j = 1, ny 
    344 !   do i = 2, nx 
    345 !      debug_3D(i,j,92) = c_west(i,j) * newH(i-1,j) + c_east(i,j) * newH(i,j) 
    346 !   end do 
    347 !end do 
    348 !do j = 2, ny 
    349 !   do i = 1, nx 
    350 !      debug_3D(i,j,93) = c_south(i,j) * newH(i,j-1) + c_north(i,j) * newH(i,j) 
    351 !   end do 
    352 !end do 
    353  
    354323if (itracebug.eq.1)  call tracebug(' fin routine resolution_diffusion') 
    355  
    356324  
    357325return 

  • branches/iLoveclim/SOURCES/spinup_mod.f90

    r123 r146  
    102102case default 
    103103  print*,'type_vitbil valeur invalide dans spinup_vitbil' 
     104  stop 
    104105end select 
    105106 
    … …  
    257258    ddbx(:,:)=0. 
    258259    ddby(:,:)=0. 
    259  
     260    gzmx(:,:)=.true. 
     261    gzmy(:,:)=.true. 
     262    flgzmx(:,:)=.false. 
     263    flgzmy(:,:)=.false. 
     264    fleuvemx(:,:)=.false. 
     265    fleuvemy(:,:)=.false. 
     266     
     267    do j=2,ny 
     268      do i=2,nx 
     269          if (flot(i,j).and.(flot(i-1,j))) then 
     270            flotmx(i,j)=.true. 
     271            flgzmx(i,j)=.true. 
     272            gzmx(i,j)=.false. 
     273          end if 
     274          if (flot(i,j).and.(flot(i,j-1))) then 
     275            flotmy(i,j)=.true. 
     276            flgzmy(i,j)=.true. 
     277            gzmy(i,j)=.false. 
     278          end if 
     279      end do 
     280    end do    
     281     
     282    where (coefmxbmelt(:,:).le.0.) 
     283      gzmx(:,:)=.false. 
     284    endwhere 
     285    where (coefmybmelt(:,:).le.0.) 
     286      gzmy(:,:)=.false. 
     287    endwhere 
     288    flgzmx(:,:) = flgzmx(:,:) .or. gzmx(:,:) 
     289    flgzmy(:,:) = flgzmy(:,:) .or. gzmy(:,:) 
     290    fleuvemx(:,:)=gzmx(:,:) 
     291    fleuvemy(:,:)=gzmy(:,:) 
     292     
    260293    if (itracebug.eq.1)  call tracebug(' Subroutine  calc_coef_vitbil') 
    261294    call slope_surf 

  • branches/iLoveclim/SOURCES/steps_time_loop.f90

    r123 r146  
    2222  use icetempmod 
    2323  use diagno_mod 
    24   use bilan_eau_mod  
     24  use bilan_eau_mod 
     25  use bilan_flux_mod 
    2526  !  use track_debug  
    2627 
    … …  
    5657 
    5758     call icethick3 
    58  
    5959     call flottab 
    60  
    6160     call calving 
    62  
    6361     call ablation_bord 
    64  
    6562     call bilan_eau 
    66      if (isynchro.eq.1) then 
    67          call shortoutput 
    68          diff_H = 0. 
    69          Bm_dtt(:,:) = 0. 
    70          bmelt_dtt(:,:) = 0. 
    71          calv_dtt(:,:)=0. 
    72          ablbord_dtt(:,:)=0. 
    73          diff_H_2D(:,:)=0. 
    74       endif 
     63     call bilan_flux_output 
     64 
     65!~      if (isynchro.eq.1) then 
     66!~          call shortoutput 
     67!~          diff_H = 0. 
     68!~          Bm_dtt(:,:) = 0. 
     69!~          bmelt_dtt(:,:) = 0. 
     70!~          calv_dtt(:,:)=0. 
     71!~          ablbord_dtt(:,:)=0. 
     72!~          diff_H_2D(:,:)=0. 
     73!~          grline_dtt(:,:)=0. 
     74!~       endif 
    7575      
    7676     call flottab 
    … …  
    158158  call testsort_time_ncdf(time) 
    159159  if (iglob_ncdf .EQ. 1) call sortie_ncdf_cat 
    160  
     160  if (isynchro.eq.1) then 
     161    call shortoutput 
     162    diff_H = 0. 
     163    Bm_dtt(:,:) = 0. 
     164    bmelt_dtt(:,:) = 0. 
     165    calv_dtt(:,:)=0. 
     166    ablbord_dtt(:,:)=0. 
     167    diff_H_2D(:,:)=0. 
     168    grline_dtt(:,:)=0. 
     169  endif 
    161170 
    162171  !   sortie compteur tous les dtcpt ans  
    … …  
    449458        end if 
    450459  
    451         iter_visco= 10                                 ! warning, test sur l'impact iteration dragging et visco  
     460!cdc        iter_visco= 10                                 ! warning, test sur l'impact iteration dragging et visco  
    452461                                                       ! Cat 18 janv 2013 
    453462 

  • branches/iLoveclim/SOURCES/util_recovery.f90

    r32 r146  
    7676    character(len=*),intent(inout) :: filin 
    7777    character(len=10)              :: ntime 
     78    character(len=1)               :: signe,unite 
    7879 
    7980    if (itracebug.eq.1)  call tracebug(' Entree dans routine testout_recovery ') 
    … …  
    8889    ! un identifiant "grestart" pour identifier ces .nc restarts 
    8990    ! des autres 
    90      
     91   
     92    predef: do indice=1,size(tab_time) 
     93      difftime=abs(TIME-tab_time(indice)) 
     94      if ( difftime.lt.dtmin ) then 
     95        ipredef=1 
     96        time_out(1)=tab_time(indice) 
     97        exit predef 
     98      end if 
     99    end do predef 
     100 
    91101    if ((mod(abs(dble(TIME)),dble(DTOUT)).lt.dble(dtmin)).OR.   & 
    92102         (ABS(TIME+297000).LT.dtmin) .OR.   & 
    … …  
    94104         (ABS(TIME+21000) .LT.dtmin) .OR.   & 
    95105         (ABS(TIME+16000) .LT.dtmin) .OR.   &  
    96          (TIME .EQ. TEND)                  )    THEN 
    97  
    98        predef: do indice=1,size(tab_time) 
    99           difftime=abs(TIME-tab_time(indice)) 
    100           if ( difftime.lt.dtmin ) then 
    101              ipredef=1 
    102              time_out(1)=tab_time(indice) 
    103              exit predef 
    104           end if 
    105        end do predef 
     106         (TIME .EQ. TEND).OR.(ipredef.eq.1) )    THEN 
    106107 
    107108       if (ipredef .EQ. 1) then 
    108           write(ntime,'(i10)')  INT(abs(dble(tab_time(indice)))/dble(DTOUT))   
    109           if (tab_time(indice) .Lt. 0) then 
    110              filin=runname//'-grestart-k'//trim(ADJUSTL(ntime)) 
    111           else  
    112              filin=runname//'-grestart+k'//trim(ADJUSTL(ntime)) 
    113           end if 
     109         ! pour changer de signe entre le passe et le futur 
     110         if (tab_time(indice).GT.0.) then 
     111          signe= '+' 
     112         else 
     113          signe= '-' 
     114         endif 
    114115 
    115        else !ipredef==0  
    116           if (TIME .ge. Tend) then 
    117              time_out(1)=TIME 
    118              filin=runname//'-grestart_2' 
    119           else 
    120              time_out(1)=TIME 
    121              filin=runname//'-grestart_1' 
    122           end if 
    123        end if 
    124        logic_out=.TRUE. 
    125     endif 
     116        if (int(mod(abs(tab_time(indice)),1000.)).eq.0) then 
     117!     temps multiple de 1000 
     118          unite='k' 
     119          write(ntime,'(i10)') int(abs(tab_time(indice)/1000.)) 
     120        else if (int(mod(abs(tab_time(indice)),100.)).eq.0) then 
     121!     temps multiple de 100 
     122          unite='c' 
     123          write(ntime,'(i10)') int(abs(tab_time(indice)/100.)) 
     124        else if (int(mod(abs(tab_time(indice)),10.)).eq.0) then 
     125!     temps multiple de 10 
     126          unite='d' 
     127          write(ntime,'(i10)') int(abs(tab_time(indice)/10.)) 
     128        else 
     129!     temps en annees 
     130          unite='a' 
     131          write(ntime,'(i10)') int(abs(tab_time(indice)/1.)) 
     132        endif 
     133        filin=runname//'-grestart'//signe//unite//trim(ADJUSTL(ntime)) 
     134 
     135     else !ipredef==0  
     136        if (TIME .ge. Tend) then 
     137           time_out(1)=TIME 
     138           filin=runname//'-grestart_2' 
     139        else 
     140           time_out(1)=TIME 
     141           filin=runname//'-grestart_1' 
     142        end if 
     143     end if 
     144     logic_out=.TRUE. 
     145  endif 
    126146 
    127147  end subroutine testout_recovery 
Note: See TracChangeset for help on using the changeset viewer.