source: trunk/LMDZ.COMMON/libf/evolution/evol_h2o_ice_s_mod_slope.F90 @ 2909

Last change on this file since 2909 was 2888, checked in by llange, 2 years ago

PEM

  • Following r-2886-5, debugging some issues ( conservation of H2O, water_reservoir not initialized, info_pem which was not working)
  • Cleaning of some redundant routines (e.g., stopping criterion of the PEM) and variables renamed (e.g., alpha_criterion now transofrmed in ice_criterion)
  • Ice table subroutine is now in a module and recalled "compute_ice_table_equilibrium" (v.s. dynamic ice table, efforts ongoing)
  • Abort_pem introduced to avoid just stopping the pem with raw "stop" in the codes
  • conf_pem has been improved so that values are not hard coded (e.g., geothermal flux, mass of water_reservoir)
  • Regolith thermal properties updated in a cleaner way

(Not atomic commit, sorry)
LL & RV
File size: 3.7 KB
Line 
1!
2! $Id $
3!
4SUBROUTINE evol_h2o_ice_s_slope(qsurf,tendencies_h2o_ice_phys,&
5                             iim_input,jjm_input,ngrid,cell_area,STOPPING,nslope)
6
7  USE temps_mod_evol, ONLY: dt_pem
8  use comslope_mod, ONLY: subslope_dist
9  use criterion_pem_stop_mod, only: criterion_waterice_stop
10      IMPLICIT NONE
11
12!=======================================================================
13!
14!  Routine that compute the evolution of the water ice
15!
16!=======================================================================
17
18!   arguments:
19!   ----------
20
21!   INPUT
22
23  INTEGER, intent(in) :: iim_input,jjm_input, ngrid,nslope   ! # of grid points along longitude/latitude/ total
24  REAL, intent(in) ::  cell_area(ngrid)
25
26!   OUTPUT
27  REAL, INTENT(INOUT) ::  qsurf(ngrid,nslope)                ! physical point field : Previous and actual density of water ice
28  LOGICAL, INTENT(INOUT) :: STOPPING
29  REAL, intent(inout) ::  tendencies_h2o_ice_phys(ngrid,nslope) ! physical point field : Evolution of perenial ice over one year
30
31!   local:
32!   ----
33
34  INTEGER :: i,j,ig0,islope                                  ! loop variable
35  REAL :: pos_tend, neg_tend, real_coefficient,negative_part
36  REAL ::  new_tendencies(ngrid,nslope)
37
38!=======================================================================
39
40  STOPPING=.false.
41
42  pos_tend=0.
43  neg_tend=0.
44
45  do i=1,ngrid
46     do islope=1,nslope
47     if (qsurf(i,islope).GT.0) then
48         if (tendencies_h2o_ice_phys(i,islope).GT.0) then
49            pos_tend=pos_tend+tendencies_h2o_ice_phys(i,islope)*cell_area(i)*subslope_dist(i,islope)
50         else
51            neg_tend=neg_tend-tendencies_h2o_ice_phys(i,islope)*cell_area(i)*subslope_dist(i,islope)
52         endif
53     endif
54     enddo
55  enddo
56
57  if(neg_tend.GT.pos_tend .and. pos_tend.GT.0) then
58     do i=1,ngrid
59       do islope=1,nslope
60       if(tendencies_h2o_ice_phys(i,islope).LT.0) then
61          new_tendencies(i,islope)=tendencies_h2o_ice_phys(i,islope)*(pos_tend/neg_tend)
62       else
63          new_tendencies(i,islope)=tendencies_h2o_ice_phys(i,islope)
64       endif
65       enddo
66     enddo
67  elseif(neg_tend.LT.pos_tend .and. neg_tend.GT.0) then
68     do i=1,ngrid
69       do islope=1,nslope
70       if(tendencies_h2o_ice_phys(i,islope).LT.0) then
71          new_tendencies(i,islope)=tendencies_h2o_ice_phys(i,islope)
72       else
73          new_tendencies(i,islope)=tendencies_h2o_ice_phys(i,islope)*(neg_tend/pos_tend)
74       endif
75       enddo
76     enddo
77  elseif(pos_tend.EQ.0 .OR. neg_tend.EQ.0) then
78    print *, "Reason of stopping : There is either no water ice sublimating or no water ice increasing !!"
79    print *, "Tendencies on ice sublimating=", neg_tend
80    print *, "Tendencies on ice increasing=", pos_tend
81    print *, "This can be due to the absence of water ice in the PCM run!!"
82      call criterion_waterice_stop(cell_area,1.,qsurf(:,:)*0.,STOPPING,ngrid,cell_area)
83      do i=1,ngrid
84         do islope=1,nslope
85          new_tendencies(i,islope)=0
86         enddo
87      enddo
88  endif
89
90  negative_part = 0.
91
92! Evolution of the water ice for each physical point
93  do i=1,ngrid
94    do islope=1, nslope
95      qsurf(i,islope)=qsurf(i,islope)+new_tendencies(i,islope)*dt_pem
96      if (qsurf(i,islope).lt.0) then
97        negative_part=negative_part-qsurf(i,islope)*cell_area(i)*subslope_dist(i,islope)
98        qsurf(i,islope)=0.
99        tendencies_h2o_ice_phys(i,islope)=0.
100      endif
101    enddo
102  enddo
103 
104 
105  if(pos_tend.eq.0) then
106   real_coefficient = 0.
107  else
108   real_coefficient = negative_part/pos_tend
109  endif
110  do i=1,ngrid
111    do islope=1, nslope
112     if(new_tendencies(i,islope).GT.0) then
113         qsurf(i,islope)=qsurf(i,islope)-new_tendencies(i,islope)*real_coefficient*dt_pem
114     endif
115    enddo
116  enddo
117
118
119
120END SUBROUTINE evol_h2o_ice_s_slope
Note: See TracBrowser for help on using the repository browser.