Changeset 424 for branches/GRISLIv3


Ignore:
Timestamp:
04/24/23 17:33:37 (15 months ago)
Author:
dumas
Message:

use only in subroutine read_recovery

Location:
branches/GRISLIv3/SOURCES
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • branches/GRISLIv3/SOURCES/3D-physique-gen_mod.f90

    r400 r424  
    347347  real,dimension(nx,ny) :: ramollo      !< pour ramollir les ice shelves 
    348348  real,dimension(nx,ny) :: Abar         !< coefficient de Glen integre 
    349   real,dimension(nx,ny) :: OLDU         !< pour cptr et dans diagnoshelf 
    350349  real,dimension(nx,ny) :: Uiter_centre !< pour iterations equation diagnostique 
    351350  real,dimension(nx,ny) :: tabtest      !< tableau de travail 
  • branches/GRISLIv3/SOURCES/out_cptr_mod.f90

    r260 r424  
    2727subroutine read_recovery(icpt)           !---- reprise d'un fichier compteur ----- 
    2828 
    29   use module3D_phy 
     29  use module3D_phy, only: reprcptr,num_forc,time,S,H,B,Bsoc,ibase,bmelt,hwater,uxbar,uybar,T,ux,uy,hmx,hmy,& 
     30       sdx,sdy,uxflgz,uyflgz,flot,debug_3D 
     31  use geography, only: dx,dy,nx,ny,nz,nzm 
     32  use runparam, only: itracebug,tbegin,num_tracebug 
    3033  use netcdf 
    31   use io_netcdf_grisli 
    32   use tracer_vars                  !aurel neem 
     34  use io_netcdf_grisli, only: ncdf_type,Read_ncdf_dim,Read_ncdf_var,lect_netcdf_type 
     35  use tracer_vars, only: xdep_out,ydep_out,tdep_out 
    3336  implicit none 
    3437 
     
    3639  integer                            :: nzzm 
    3740  integer                            :: icpt  ! icompteur 1->  tout (topo, T, U) 
    38                                               !           2 -> tout sauf topo  
    39                                               !           3 -> tout sauf topo et eau basale 
    40   
    41  
    42  real,dimension(nx,ny)              :: bidon         !< sert à sauter les lectures de Ss,H,B,... 
    43                                
     41  !           2 -> tout sauf topo  
     42  !           3 -> tout sauf topo et eau basale 
     43 
     44 
     45  real,dimension(nx,ny)              :: bidon         !< sert à sauter les lectures de Ss,H,B,... 
     46 
    4447  ! Variables specifiques pour les fichier .nc 
    45    
     48 
    4649  real*8, dimension(:,:),   pointer  :: tab => null()         !< tableau 2d real ecrit dans le fichier 
    4750  real*8, dimension(:,:,:), pointer  :: tab1  => null()       !< tableau 3d real 
     
    4952  real*8, dimension(:),     pointer  :: liste_time => null() 
    5053 
    51 ! rappel de f90 
    52 ! reprcptr est un string (ici le nom du fichier recovery) 
    53 ! index (reprcptr,'.cptr')    renvoie : 
    54 !                 0, si 'cptr' n'est pas dans reprcptr  
    55 !                 la position de la premiere occurence de cptr sinon 
     54  integer :: ni,nzz,nxx,nyy 
     55  integer :: i,j 
     56  real,dimension(nx,ny) :: oldu 
     57 
     58  ! rappel de f90 
     59  ! reprcptr est un string (ici le nom du fichier recovery) 
     60  ! index (reprcptr,'.cptr')    renvoie : 
     61  !                 0, si 'cptr' n'est pas dans reprcptr  
     62  !                 la position de la premiere occurence de cptr sinon 
    5663 
    5764  if (ncdf_type.eq.0) call lect_netcdf_type         !< pour lire la valeur de netcdf_type (machine dependant) 
    5865 
    5966 
    60 ascii_nc:  if ((index(reprcptr,'.cptr').ne.0) .or.(index(reprcptr,'.CPTR').ne.0))  then   !reprise fichier cptr 
    61    
     67  ascii_nc:  if ((index(reprcptr,'.cptr').ne.0) .or.(index(reprcptr,'.CPTR').ne.0))  then   !reprise fichier cptr 
     68 
    6269     if (itracebug.eq.1)  call tracebug(' Entree dans routine read_recovery : ascii') 
    6370 
     
    6673 
    6774 
    68 ! temps de la reprise 
    69 !----------------------- 
     75     ! temps de la reprise 
     76     !----------------------- 
    7077     if (tbegin.gt.1.d9) then       !    si tbegin est tres grand , on reprend le temps du cptr 
    7178        tbegin=time 
     
    8895     read(num_forc,*) 
    8996 
    90 lect_topo:     if (icpt.eq.1) then                     ! lecture topo seulement si icompteur = 1   
     97     lect_topo:     if (icpt.eq.1) then                     ! lecture topo seulement si icompteur = 1   
    9198        read(num_forc,*) S 
    9299        read(num_forc,*) 
     
    116123     end if lect_topo 
    117124 
    118                                                        ! lecture temperature dans tous les cas      
     125     ! lecture temperature dans tous les cas      
    119126     read(num_forc,*) Ibase 
    120127     read(num_forc,*) 
     
    123130     read(num_forc,*) 
    124131 
    125 lect_Hwater: if  (icpt.ne.3) then                      ! pas de lecture Hwater si icompteur = 3   
     132     lect_Hwater: if  (icpt.ne.3) then                      ! pas de lecture Hwater si icompteur = 3   
    126133        read(num_forc,*) Hwater 
    127134        read(num_forc,*) 
     
    133140     end if lect_Hwater 
    134141 
    135                                                        ! lecture vitesses dans tous les cas. 
    136         read(num_forc,*) Uxbar 
    137         read(num_forc,*) 
    138  
    139         read(num_forc,*) Uybar 
     142     ! lecture vitesses dans tous les cas. 
     143     read(num_forc,*) Uxbar 
     144     read(num_forc,*) 
     145 
     146     read(num_forc,*) Uybar 
    140147 
    141148 
     
    16016724   continue 
    161168 
    162 ! aurel pour faire des cptr avec les traceurs : 
    163       read(num_forc,*,end=25) 
    164       read(num_forc,*,end=25) 
    165       read(num_forc,*,end=25) Xdep_out  !xdepk ou xdep ? a voir 
    166 25    continue 
    167  
    168       read(num_forc,*,end=26) 
    169       read(num_forc,*,end=26) 
    170       read(num_forc,*,end=26) Ydep_out 
    171 26    continue 
    172  
    173       read(num_forc,*,end=27) 
    174       read(num_forc,*,end=27) 
    175       read(num_forc,*,end=27) tdep_out 
    176 27    continue 
    177       tdep_out(:,:,:)=tdep_out(:,:,:)+time 
     169     ! aurel pour faire des cptr avec les traceurs : 
     170     read(num_forc,*,end=25) 
     171     read(num_forc,*,end=25) 
     172     read(num_forc,*,end=25) Xdep_out  !xdepk ou xdep ? a voir 
     17325   continue 
     174 
     175     read(num_forc,*,end=26) 
     176     read(num_forc,*,end=26) 
     177     read(num_forc,*,end=26) Ydep_out 
     17826   continue 
     179 
     180     read(num_forc,*,end=27) 
     181     read(num_forc,*,end=27) 
     182     read(num_forc,*,end=27) tdep_out 
     18327   continue 
     184     tdep_out(:,:,:)=tdep_out(:,:,:)+time 
    178185 
    179186     close(num_forc) 
     
    184191 
    185192 
    186         if (.not.associated(liste_time)) then 
    187            allocate(liste_time(1))  
    188         end if 
    189         if (.not.associated(tab)) then 
    190            allocate(tab(nx,ny),tab1(nx,ny,nz),tab1T(nx,ny,nz+nzm)) 
    191         end if 
    192  
    193         call Read_ncdf_var('temps',TRIM(reprcptr),liste_time) 
    194  
    195 ! temps de la reprise 
    196 !----------------------- 
     193     if (.not.associated(liste_time)) then 
     194        allocate(liste_time(1))  
     195     end if 
     196     if (.not.associated(tab)) then 
     197        allocate(tab(nx,ny),tab1(nx,ny,nz),tab1T(nx,ny,nz+nzm)) 
     198     end if 
     199 
     200     call Read_ncdf_var('temps',TRIM(reprcptr),liste_time) 
     201 
     202     ! temps de la reprise 
     203     !----------------------- 
    197204     if (tbegin.gt.1.d9) then       !    si tbegin est tres grand , on reprend le temps du .nc 
    198205        tbegin=liste_time(1) 
     
    205212 
    206213 
    207 !        TBEGIN=liste_time(1)-liste_time(1) 
    208         
    209         call Read_ncdf_dim('x' ,TRIM(reprcptr),nxx) 
    210         call Read_ncdf_dim('y' ,TRIM(reprcptr),nyy) 
    211         call Read_ncdf_dim('z' ,TRIM(reprcptr),nzz) 
    212         call Read_ncdf_dim('zm',TRIM(reprcptr),nzzm) 
    213  
    214         if ((nxx.ne.nx).or.(nyy.ne.ny).or.(nzz.ne.nz)       & 
    215              .or.(nzzm-nzz.ne.nzm)) write(6,*) 'attention', & 
    216              'les tailles de tableaux ne sont pas compatibles' 
    217  
    218         !     lecture des tableaux 2D 
    219         ! ----------------------------- 
    220  
    221         if (icpt.eq.1) then                             ! lecture topo seulement si icompteur = 1 
    222            call Read_ncdf_var('S',TRIM(reprcptr),tab) 
    223            S(:,:)=tab(:,:) 
    224            call Read_ncdf_var('H',TRIM(reprcptr),tab) 
    225            H(:,:)=tab(:,:) 
    226            call Read_ncdf_var('B',TRIM(reprcptr),tab) 
    227            B(:,:)=tab(:,:) 
    228            call Read_ncdf_var('BSOC',TRIM(reprcptr),tab) 
    229            Bsoc(:,:)=tab(:,:) 
    230         end if 
    231  
    232                                                        ! lecture temperature dans tous les cas       
    233         call Read_ncdf_var('IBASE',TRIM(reprcptr),tab) 
    234         Ibase(:,:)=tab(:,:) 
    235         call Read_ncdf_var('BMELT',TRIM(reprcptr),tab) 
    236         Bmelt(:,:)=tab(:,:) 
    237  
    238         if (icpt.ne.3) then                            ! pas de lecture Hwater si icompteur = 3 
    239            call Read_ncdf_var('HWATER',TRIM(reprcptr),tab) 
    240            Hwater(:,:)=tab(:,:) 
    241         else 
    242            Hwater(:,:)=0. 
    243         end if 
    244  
    245                                                        ! lecture vitesses dans tous les cas. 
    246         call Read_ncdf_var('UXBAR',TRIM(reprcptr),tab) 
    247         UXBAR(:,:)=tab(:,:) 
    248         call Read_ncdf_var('UYBAR',TRIM(reprcptr),tab) 
    249         UYBAR(:,:)=tab(:,:) 
    250  
    251         do i=1,nx-1 
    252            do j=1,ny-1 
    253               oldu(i,j) = sqrt(((uxbar(i,j)+uxbar(i+1,j))/2)**2 & 
    254                    + ((uybar(i,j)+uybar(i,j+1))/2.)**2) 
    255            end do 
    256         end do 
    257  
    258         !     tableaux 3D 
    259         ! ----------------------------- 
    260  
    261         call Read_ncdf_var('T',TRIM(reprcptr),tab1T) 
    262         T(:,:,:)=tab1T(:,:,:) 
    263  
    264         call Read_ncdf_var('UX',TRIM(reprcptr),tab1) 
    265         Ux(:,:,:)=tab1(:,:,:) 
    266  
    267         call Read_ncdf_var('UY',TRIM(reprcptr),tab1) 
    268         Uy(:,:,:)=tab1(:,:,:) 
    269        
    270         if (itracebug.eq.1)  call tracebug('avant lectures variables traceur') 
    271  
    272  
    273         call Read_ncdf_var('XDEP',TRIM(reprcptr),tab1)  !aurel neem 
    274         Xdep_out(:,:,:)=tab1(:,:,:) 
    275  
    276         call Read_ncdf_var('YDEP',TRIM(reprcptr),tab1) 
    277         Ydep_out(:,:,:)=tab1(:,:,:) 
    278  
    279         call Read_ncdf_var('TDEP',TRIM(reprcptr),tab1) 
    280         Tdep_out(:,:,:)=tab1(:,:,:)+time 
    281   
    282      end if ascii_nc 
    283  
     214     !        TBEGIN=liste_time(1)-liste_time(1) 
     215 
     216     call Read_ncdf_dim('x' ,TRIM(reprcptr),nxx) 
     217     call Read_ncdf_dim('y' ,TRIM(reprcptr),nyy) 
     218     call Read_ncdf_dim('z' ,TRIM(reprcptr),nzz) 
     219     call Read_ncdf_dim('zm',TRIM(reprcptr),nzzm) 
     220 
     221     if ((nxx.ne.nx).or.(nyy.ne.ny).or.(nzz.ne.nz)       & 
     222          .or.(nzzm-nzz.ne.nzm)) write(6,*) 'attention', & 
     223          'les tailles de tableaux ne sont pas compatibles' 
     224 
     225     !     lecture des tableaux 2D 
     226     ! ----------------------------- 
     227 
     228     if (icpt.eq.1) then                             ! lecture topo seulement si icompteur = 1 
     229        call Read_ncdf_var('S',TRIM(reprcptr),tab) 
     230        S(:,:)=tab(:,:) 
     231        call Read_ncdf_var('H',TRIM(reprcptr),tab) 
     232        H(:,:)=tab(:,:) 
     233        call Read_ncdf_var('B',TRIM(reprcptr),tab) 
     234        B(:,:)=tab(:,:) 
     235        call Read_ncdf_var('BSOC',TRIM(reprcptr),tab) 
     236        Bsoc(:,:)=tab(:,:) 
     237     end if 
     238 
     239     ! lecture temperature dans tous les cas       
     240     call Read_ncdf_var('IBASE',TRIM(reprcptr),tab) 
     241     Ibase(:,:)=tab(:,:) 
     242     call Read_ncdf_var('BMELT',TRIM(reprcptr),tab) 
     243     Bmelt(:,:)=tab(:,:) 
     244 
     245     if (icpt.ne.3) then                            ! pas de lecture Hwater si icompteur = 3 
     246        call Read_ncdf_var('HWATER',TRIM(reprcptr),tab) 
     247        Hwater(:,:)=tab(:,:) 
     248     else 
     249        Hwater(:,:)=0. 
     250     end if 
     251 
     252     ! lecture vitesses dans tous les cas. 
     253     call Read_ncdf_var('UXBAR',TRIM(reprcptr),tab) 
     254     UXBAR(:,:)=tab(:,:) 
     255     call Read_ncdf_var('UYBAR',TRIM(reprcptr),tab) 
     256     UYBAR(:,:)=tab(:,:) 
    284257 
    285258     do i=1,nx-1 
    286259        do j=1,ny-1 
    287            oldu(i,j) = sqrt(((Uxbar(i,j)+Uxbar(i+1,j))/2)**2 & 
    288                 + ((Uybar(i,j)+Uybar(i,j+1))/2.)**2) 
     260           oldu(i,j) = sqrt(((uxbar(i,j)+uxbar(i+1,j))/2)**2 & 
     261                + ((uybar(i,j)+uybar(i,j+1))/2.)**2) 
    289262        end do 
    290263     end do 
    291264 
    292  
    293      !    moyenne de l'epaisseur et pentes 
    294      do j=2,ny 
    295         do i=2,nx 
    296            Hmy(i,j)=(H(i,j)+H(i,j-1))/2. 
    297            Hmx(i,j)=(H(i,j)+H(i-1,j))/2. 
    298            Sdx(i,j)=(S(i,j)-S(i-1,j))/dx 
    299            Sdy(i,j)=(S(i,j)-S(i,j-1))/dy 
    300         end do 
     265     !     tableaux 3D 
     266     ! ----------------------------- 
     267 
     268     call Read_ncdf_var('T',TRIM(reprcptr),tab1T) 
     269     T(:,:,:)=tab1T(:,:,:) 
     270 
     271     call Read_ncdf_var('UX',TRIM(reprcptr),tab1) 
     272     Ux(:,:,:)=tab1(:,:,:) 
     273 
     274     call Read_ncdf_var('UY',TRIM(reprcptr),tab1) 
     275     Uy(:,:,:)=tab1(:,:,:) 
     276 
     277     if (itracebug.eq.1)  call tracebug('avant lectures variables traceur') 
     278 
     279 
     280     call Read_ncdf_var('XDEP',TRIM(reprcptr),tab1)  !aurel neem 
     281     Xdep_out(:,:,:)=tab1(:,:,:) 
     282 
     283     call Read_ncdf_var('YDEP',TRIM(reprcptr),tab1) 
     284     Ydep_out(:,:,:)=tab1(:,:,:) 
     285 
     286     call Read_ncdf_var('TDEP',TRIM(reprcptr),tab1) 
     287     Tdep_out(:,:,:)=tab1(:,:,:)+time 
     288 
     289  end if ascii_nc 
     290 
     291 
     292  do i=1,nx-1 
     293     do j=1,ny-1 
     294        oldu(i,j) = sqrt(((Uxbar(i,j)+Uxbar(i+1,j))/2)**2 & 
     295             + ((Uybar(i,j)+Uybar(i,j+1))/2.)**2) 
    301296     end do 
    302  
     297  end do 
     298 
     299 
     300  !    moyenne de l'epaisseur et pentes 
     301  do j=2,ny 
    303302     do i=2,nx 
    304         Hmx(i,1)=(H(i,1)+H(i-1,1))/2. 
     303        Hmy(i,j)=(H(i,j)+H(i,j-1))/2. 
     304        Hmx(i,j)=(H(i,j)+H(i-1,j))/2. 
     305        Sdx(i,j)=(S(i,j)-S(i-1,j))/dx 
     306        Sdy(i,j)=(S(i,j)-S(i,j-1))/dy 
    305307     end do 
    306      do j=2,ny 
    307         Hmy(1,j)=(H(1,j)+H(1,j-1))/2. 
    308      end do 
    309  
    310      ! vitesse a la base 
    311      Uxflgz(:,:) = Ux(:,:,nz) 
    312      Uyflgz(:,:) = Uy(:,:,nz) 
    313       
    314      where (flot(:,:)) 
    315        debug_3D(:,:,122) = T(:,:,nz)+273.15 
    316      elsewhere 
    317        debug_3D(:,:,120) = T(:,:,nz)+273.15 
    318      endwhere 
     308  end do 
     309 
     310  do i=2,nx 
     311     Hmx(i,1)=(H(i,1)+H(i-1,1))/2. 
     312  end do 
     313  do j=2,ny 
     314     Hmy(1,j)=(H(1,j)+H(1,j-1))/2. 
     315  end do 
     316 
     317  ! vitesse a la base 
     318  Uxflgz(:,:) = Ux(:,:,nz) 
     319  Uyflgz(:,:) = Uy(:,:,nz) 
     320 
     321  where (flot(:,:)) 
     322     debug_3D(:,:,122) = T(:,:,nz)+273.15 
     323  elsewhere 
     324     debug_3D(:,:,120) = T(:,:,nz)+273.15 
     325  endwhere 
    319326 
    320327end subroutine read_recovery 
     
    391398  character(len=20),dimension(4)    :: dimnames3d      !< pour 3d troisieme dim est nz 
    392399  character(len=20),dimension(4)    :: dimnames3dT     !< pour 3d troisieme dim est nz+nzm 
    393 !  integer                           :: nrecs=1         !< compteur pour les enregistrements temps des variables   
    394 !  integer                           :: idef=0          !< pour savoir si la variable a ete definie ou non 
     400  !  integer                           :: nrecs=1         !< compteur pour les enregistrements temps des variables   
     401  !  integer                           :: idef=0          !< pour savoir si la variable a ete definie ou non 
    395402  integer                           :: isortie         !< pour le choix de type de fichier en sortie  
    396403 
     
    455462        write(num_file1,*) Uy(:,:,:) 
    456463 
    457 !aurel neem : traceurs dans le cptr 
     464        !aurel neem : traceurs dans le cptr 
    458465 
    459466        write(num_file1,*) 
     
    480487        filin=trim(filin)//'.nc' 
    481488 
    482 ! creation du fichier 
    483 !------------------------ 
     489        ! creation du fichier 
     490        !------------------------ 
    484491 
    485492        if (ncdf_type.eq.32) then 
     
    498505        call write_ncdf_var('temps','time',trim(filin),time_out,'double')     ! ecrit le temps du snapshot 
    499506 
    500 !       variables 2D 
     507        !       variables 2D 
    501508 
    502509        tab(:,:)= S(:,:) 
     
    527534        call write_ncdf_var('UYBAR',dimnames2d,trim(filin),tab,'double') 
    528535 
    529 !       variables 3D  
     536        !       variables 3D  
    530537 
    531538        tab1T(:,:,:)= T(:,:,:) 
     
    538545        call write_ncdf_var('UY',dimnames3d,trim(filin),tab1,'double') 
    539546 
    540 !       variables 3D  pour les traceurs (Aurel NEEM) 
    541 !       attention ce sont les tabelaux _out qui ont les bonnes dimensions. 
     547        !       variables 3D  pour les traceurs (Aurel NEEM) 
     548        !       attention ce sont les tabelaux _out qui ont les bonnes dimensions. 
    542549 
    543550        tab1(:,:,:)=XDEP_out(:,:,:) 
Note: See TracChangeset for help on using the changeset viewer.