Changeset 140 for LMDZ.3.3

Timestamp:

Oct 16, 2000, 5:18:52 PM (24 years ago)

Author:

lmdzadmin

Message:

Version debugge coupleur
LF

File:

• LMDZ.3.3/branches/rel-LF/libf/phylmd/interface_surf.F90 (modified) (19 diffs)

LMDZ.3.3/branches/rel-LF/libf/phylmd/interface_surf.F90

@@ -364 +364 @@
 
       tsurf_temp = tsurf_new
-      cal = 0.
       dif_grnd = 0.
-      beta = 1.
+      beta = 1.0
 
 !    else if (ocean == 'slab  ') then
@@ -376 +375 @@
      &  tsurf_new, pctsrf_new)
 
-      cal = calice
-      where (snow > 0.0) cal = calsno
+      tsurf_temp = tsurf
+      dif_grnd = 1.0 / tau_gl
       beta = 1.0
-      dif_grnd = 1.0 / tau_gl
-      tsurf_temp = tsurf
-    endif
+    endif
+
+    cal = calice
+    where (snow > 0.0) cal = calsno
 
     call calcul_fluxs( knon, nisurf, dtime, &
@@ -437 +437 @@
     abort_message = 'Index surface non valable'
     call abort_gcm(modname,abort_message,1)
-  endif
-
-  if (check) then
-    write(*,*)'In ',modname
-    do ii = 1, nbsrf
-      write(*,*) 'surface, pctsrf_new',ii,pctsrf_new(:,ii)
-    enddo
   endif
 
@@ -735 +728 @@
   real, dimension(knon), intent(IN) :: precip_rain, precip_snow
   real, dimension(knon), intent(IN) :: tsurf, fder, albsol, taux, tauy
-  INTEGER              :: nexca, npas
+  INTEGER              :: nexca, npas, kstep
   real, dimension(klon), intent(IN) :: zmasq
 
@@ -745 +738 @@
 
 ! Variables locales
-  integer                    :: j, error, sum_error, ig
+  integer                    :: j, error, sum_error, ig, cpl_index,i
   character (len = 20) :: modname = 'interfoce_cpl'
   character (len = 80) :: abort_message
@@ -772 +765 @@
 ! variable tampon
   real, dimension(klon)       :: tamp
-  real, dimension(knon)       :: tamp_sic
+  real, dimension(klon)       :: tamp_sic
+! sauvegarde des fractions de surface d'un pas de temps a l'autre apres 
+! l'avoir lu
+  real, allocatable,dimension(:,:),save :: pctsrf_sav
   real, dimension(iim, jjm+1, 2) :: tamp_srf
   integer, allocatable, dimension(:), save :: tamp_ind
@@ -778 +774 @@
   real, dimension(iim, jjm+1) :: deno
   integer                     :: idtime
+  integer, allocatable, dimension(:,:) :: isst
+  integer, allocatable,save,dimension(:) :: unity
 ! 
   logical, save    :: first_appel = .true.
+  logical          :: print
 
 !
@@ -787 +786 @@
  
   if (first_appel) then
+    error = 0
+    allocate(unity(klon), stat = error)
+    if ( error  /=0) then
+      abort_message='Pb allocation variable unity'
+      call abort_gcm(modname,abort_message,1)
+    endif
+    allocate(pctsrf_sav(klon,2), stat = error)
+    if ( error  /=0) then
+      abort_message='Pb allocation variable pctsrf_sav'
+      call abort_gcm(modname,abort_message,1)
+    endif
+
+    do ig = 1, klon
+      unity(ig) = ig
+    enddo
     sum_error = 0
     allocate(cpl_sols(knon,2), stat = error); sum_error = sum_error + error
@@ -825 +839 @@
     idtime = int(dtime)
     call inicma(npas , nexca, idtime,(jjm+1)*iim)
-!
-! 1ere lecture champs ocean
-!
-!    if (nisurf == is_oce) then
-!      call fromcpl(itime - 1,(jjm+1)*iim,                                  &
-!     &        read_sst, read_sic, read_sit, read_alb_sic)
-!
-! je voulais utiliser des where mais ca ne voulait pas compiler dans un 
-! if construct sur sun
-!
-!      do j = 1, jjm + 1
-!        do ig = 1, iim
-!          if (abs(1. - read_sic(ig,j)) < 0.00001) then
-!            read_sst(ig,j) = RTT - 1.8
-!            read_sit(ig,j) = read_sit(ig,j) / read_sic(ig,j)
-!            read_alb_sic(ig,j) = read_alb_sic(ig,j) / read_sic(ig,j)
-!          else if (abs(read_sic(ig,j)) < 0.00001) then
-!            read_sst(ig,j) = read_sst(ig,j) / (1. - read_sic(ig,j))
-!            read_sit(ig,j) = read_sst(ig,j)
-!            read_alb_sic(ig,j) =  0.6
-!          else
-!            read_sst(ig,j) = read_sst(ig,j) / (1. - read_sic(ig,j))
-!            read_sit(ig,j) = read_sit(ig,j) / read_sic(ig,j)
-!            read_alb_sic(ig,j) = read_alb_sic(ig,j) / read_sic(ig,j)
-!          endif
-!        enddo
-!      enddo
-!    endif
 
     first_appel = .false.
@@ -861 +847 @@
 ! calcul des fluxs a passer
 
-  cpl_sols(:,nisurf) = cpl_sols(:,nisurf) + swdown      / FLOAT(nexca)
-  cpl_nsol(:,nisurf) = cpl_nsol(:,nisurf) + lwdown      / FLOAT(nexca)
-  cpl_rain(:,nisurf) = cpl_rain(:,nisurf) + precip_rain / FLOAT(nexca)
-  cpl_snow(:,nisurf) = cpl_snow(:,nisurf) + precip_snow / FLOAT(nexca)
-  cpl_evap(:,nisurf) = cpl_evap(:,nisurf) + evap        / FLOAT(nexca)
-  cpl_tsol(:,nisurf) = cpl_tsol(:,nisurf) + tsurf       / FLOAT(nexca)
-  cpl_fder(:,nisurf) = cpl_fder(:,nisurf) + fder        / FLOAT(nexca)
-  cpl_albe(:,nisurf) = cpl_albe(:,nisurf) + albsol      / FLOAT(nexca)
-  cpl_taux(:,nisurf) = cpl_taux(:,nisurf) + taux        / FLOAT(nexca)
-  cpl_tauy(:,nisurf) = cpl_tauy(:,nisurf) + tauy        / FLOAT(nexca)
-  cpl_rriv(:,nisurf) = cpl_rriv(:,nisurf) + run_off     / FLOAT(nexca)/dtime
-  cpl_rcoa(:,nisurf) = cpl_rcoa(:,nisurf) + run_off     / FLOAT(nexca)/dtime
+  cpl_index = 1
+  if (nisurf == is_sic) cpl_index = 2
+  do ig = 1, knon
+    cpl_sols(ig,cpl_index) = cpl_sols(ig,cpl_index) &
+     &                               + swdown(ig)      / FLOAT(nexca)
+    cpl_nsol(ig,cpl_index) = cpl_nsol(ig,cpl_index) &
+     &                               + lwdown(ig)      / FLOAT(nexca)
+    cpl_rain(ig,cpl_index) = cpl_rain(ig,cpl_index) &
+     &                               + precip_rain(ig) / FLOAT(nexca)
+    cpl_snow(ig,cpl_index) = cpl_snow(ig,cpl_index) &
+     &                               + precip_snow(ig) / FLOAT(nexca)
+    cpl_evap(ig,cpl_index) = cpl_evap(ig,cpl_index) &
+     &                               + evap(ig)        / FLOAT(nexca)
+    cpl_tsol(ig,cpl_index) = cpl_tsol(ig,cpl_index) &
+     &                               + tsurf(ig)       / FLOAT(nexca)
+    cpl_fder(ig,cpl_index) = cpl_fder(ig,cpl_index) &
+     &                               + fder(ig)        / FLOAT(nexca)
+    cpl_albe(ig,cpl_index) = cpl_albe(ig,cpl_index) &
+     &                               + albsol(ig)      / FLOAT(nexca)
+    cpl_taux(ig,cpl_index) = cpl_taux(ig,cpl_index) &
+     &                               + taux(ig)        / FLOAT(nexca)
+    cpl_tauy(ig,cpl_index) = cpl_tauy(ig,cpl_index) &
+     &                               + tauy(ig)        / FLOAT(nexca)
+    cpl_rriv(ig,cpl_index) = cpl_rriv(ig,cpl_index) &
+     &                               + 0.     / FLOAT(nexca)/dtime
+    cpl_rcoa(ig,cpl_index) = cpl_rcoa(ig,cpl_index) &
+     &                               + 0.     / FLOAT(nexca)/dtime
+  enddo
 
   if (mod(itime, nexca) == 1) then
@@ -900 +902 @@
         enddo
       enddo
+!
+! transformer read_sic en pctsrf_sav
+!
+    call cpl2gath(read_sic, tamp_sic , klon, klon,iim,jjm, unity)
+    do ig = 1, klon
+      IF (pctsrf(ig,is_oce) > epsfra .OR.            &
+     &             pctsrf(ig,is_sic) > epsfra) THEN
+            pctsrf_sav(ig,is_sic) = tamp_sic(ig) * pctsrf(ig,is_sic)
+            pctsrf_sav(ig,is_oce) = pctsrf(ig,is_oce)    &
+     &                        - (pctsrf_sav(ig,is_sic)-pctsrf(ig,is_sic))
+      endif
+    enddo
+    if (minval(pctsrf_new(:,is_oce)) < 0.) then
+      write(*,*)'Pb fraction ocean inferieure a 0'
+      write(*,*)'au point ',minloc(pctsrf_new(:,is_oce))
+      write(*,*)'valeur = ',minval(pctsrf_new(:,is_oce)) 
+      abort_message = 'voir ci-dessus'
+      call abort_gcm(modname,abort_message,1)
+    endif
+    if (minval(pctsrf_new(:,is_sic)) < 0.) then
+      write(*,*)'Pb fraction glace inferieure a 0'
+      write(*,*)'au point ',minloc(pctsrf_new(:,is_sic))
+      write(*,*)'valeur = ',minval(pctsrf_new(:,is_sic)) 
+      abort_message = 'voir ci-dessus'
+      call abort_gcm(modname,abort_message,1)
+    endif
+       if (check) then
+         write(47,*)'Sortie fromcpl apres bidouille'
+         write(47,*)' read_sst = '
+         write(47,'(72f8.3)')read_sst
+         call flush(47)
+!        allocate(isst(iim, jjm+1), stat = error)
+!        isst = 0
+!        where (read_sst >0.) isst = 1
+!        write(46,'(72i1)')isst
+       endif
     endif 
   endif                         ! fin mod(itime, nexca) == 1
@@ -1008 +1046 @@
   if (nisurf == is_oce) then
     call cpl2gath(read_sst, tsurf_new, klon, knon,iim,jjm, knindex)
-    call cpl2gath(read_sic, tamp_sic , klon, knon,iim,jjm, knindex)
-!
-! transformer tamp_sic en pctsrf_new
-!
-    do ig = 1, klon
-      IF (pctsrf(ig,is_oce) > epsfra .OR.            &
-     &             pctsrf(ig,is_sic) > epsfra) THEN
-            pctsrf_new(ig,is_oce) = pctsrf(ig,is_oce)    &
-     &                        - (tamp_sic(ig)-pctsrf(ig,is_sic))
-            pctsrf_new(ig,is_sic) = tamp_sic(ig)
-      endif
-    enddo
-    if (check) write(*,*)'In ',modname
-    if (check) write(*,*)' surface, pctsrf_new',is_oce,pctsrf_new(:,is_oce)
-    if (check) write(*,*)' surface, pctsrf_new',is_sic,pctsrf_new(:,is_sic)
   else if (nisurf == is_sic) then
-      call cpl2gath(read_sit, tsurf_new, klon, knon,iim,jjm, knindex)
-      call cpl2gath(read_alb_sic, alb_new, klon, knon,iim,jjm, knindex)
+    call cpl2gath(read_sit, tsurf_new, klon, knon,iim,jjm, knindex)
+    call cpl2gath(read_alb_sic, alb_new, klon, knon,iim,jjm, knindex)
   endif
+  pctsrf_new(:,nisurf) = pctsrf_sav(:,nisurf)
   
 !  if (lafin) call quitcpl
@@ -1132 +1156 @@
   
   if (check) write(*,*)modname,' :: jour_lu, deja_lu', jour_lu, deja_lu 
+  if (check) write(*,*)modname,' :: itime, lmt_pas ', itime, lmt_pas,dtime
 
 ! Tester d'abord si c'est le moment de lire le fichier
@@ -1339 +1364 @@
 ! Variables locales
   integer     :: ii
-  integer     :: lmt_pas     ! frequence de lecture des conditions limites 
+  integer,save :: lmt_pas     ! frequence de lecture des conditions limites 
                              ! (en pas de physique)
   logical,save :: deja_lu_sur! pour indiquer que le jour a lire a deja
@@ -1372 +1397 @@
   
   if (check) write(*,*)modname,':: jour_lu_sur, deja_lu_sur', jour_lu_sur, deja_lu_sur 
+  if (check) write(*,*)modname,':: itime, lmt_pas', itime, lmt_pas
+  call flush(6)
 
 ! Tester d'abord si c'est le moment de lire le fichier
@@ -1517 +1544 @@
   real                  :: bilan_f, fq_fonte
   real, parameter :: t_grnd = 271.35, t_coup = 273.15
-  logical         :: check = .true.
+  logical         :: check = .false.
   character (len = 20)  :: modname = 'calcul_fluxs'
   logical         :: fonte_neige = .false.
@@ -1523 +1550 @@
   character (len = 80) :: abort_message 
 
-  if (check) write(*,*)'Entree ', modname
+  if (check) write(*,*)'Entree ', modname,' surface = ',nisurf
 
   if (size(run_off) /= knon .AND. nisurf == is_ter) then
@@ -1685 +1712 @@
     tamp(ig) = champ_in(i)
   enddo    
-  champ_out(:,1) = tamp(1)
+  ig = 1
+  champ_out(:,1) = tamp(ig)
   do j = 2, jjm
     do i = 1, iim
-      champ_out(i,j) = tamp((j-2)*jjm + i + 1)
+      ig = ig + 1
+      champ_out(i,j) = tamp(ig)
     enddo
   enddo
-  champ_out(:,jjm+1) = tamp(klon)
+  ig = ig + 1
+  champ_out(:,jjm+1) = tamp(ig)
 
   END SUBROUTINE gath2cpl
@@ -1722 +1752 @@
   integer                   :: i, ig, j
   real, dimension(klon)     :: tamp
-
-  tamp(1) = champ_in(1,1)
+  logical                   :: check = .false.
+
+  ig = 1
+  tamp(ig) = champ_in(1,1)
   do j = 2, jjm
     do i = 1, iim
-      tamp((j-2)*jjm + i + 1) = champ_in(i,j)
+      ig = ig + 1
+      tamp(ig) = champ_in(i,j)
     enddo
   enddo
-  tamp(klon) = champ_in(1,jjm+1)
+  ig = ig + 1
+  tamp(ig) = champ_in(1,jjm+1)
 
   do i = 1, knon