[7541] | 1 | ! ================================================================================================================================= |
---|
| 2 | ! PROGRAM : testrouting |
---|
| 3 | ! |
---|
| 4 | ! CONTACT : jan.polcher _at_ lmd.jussieu.fr |
---|
| 5 | ! |
---|
| 6 | ! LICENCE : :) |
---|
| 7 | ! |
---|
| 8 | !>\BRIEF This program tests routing scheme (from routing.f90) which routes the water over the continents |
---|
| 9 | !! into the oceans and computes the water stored in floodplains or taken for irrigation. |
---|
| 10 | !! |
---|
| 11 | !!\n DESCRIPTION: None |
---|
| 12 | !! |
---|
| 13 | !! RECENT CHANGE(S): None |
---|
| 14 | !! |
---|
| 15 | !! REFERENCE(S) : |
---|
| 16 | !! |
---|
| 17 | !! SVN : |
---|
| 18 | !! $HeadURL : svn://forge.ipsl.jussieu.fr/orchidee/trunk/ORCHIDEE/[somewhere]/testrouting.f90 $ |
---|
| 19 | !! $Date: 2014-10-27 10:39:00 +0200 (Mon, 27 Oct 2014) $ |
---|
| 20 | !! $Revision: XXXX $ |
---|
| 21 | !! \n |
---|
| 22 | !_ ================================================================================================================================ |
---|
| 23 | PROGRAM testrouting |
---|
| 24 | ! |
---|
| 25 | USE ioipsl_para |
---|
| 26 | USE pft_parameters |
---|
| 27 | USE mod_orchidee_para |
---|
| 28 | USE control |
---|
| 29 | USE constantes_soil_var |
---|
| 30 | USE constantes_var |
---|
| 31 | USE constantes |
---|
| 32 | USE time |
---|
| 33 | USE routing |
---|
| 34 | USE timer |
---|
| 35 | USE grid |
---|
| 36 | ! |
---|
| 37 | USE getlandseamask |
---|
| 38 | ! |
---|
| 39 | IMPLICIT NONE |
---|
| 40 | ! |
---|
| 41 | INTEGER(i_std) :: nbseg !! |
---|
| 42 | ! |
---|
| 43 | INTEGER(i_std) :: iim !! Size in longitude of coarser grid |
---|
| 44 | INTEGER(i_std) :: jjm !! Size in latitude of coarser grid |
---|
| 45 | INTEGER(i_std) :: i, j !! Integer variable for loops |
---|
| 46 | INTEGER(i_std) :: ibegt, iendt |
---|
| 47 | INTEGER(i_std) :: ni !! For checking nbindex |
---|
| 48 | REAL(r_std) :: nbyears !! Lenght of simulation in years |
---|
| 49 | INTEGER(i_std) :: simlen !! Lenght of simulation: simlen = 365*48*nbyears |
---|
| 50 | REAL(r_std) :: dx, dy !! Lon/Lat resolution of coarser grid |
---|
| 51 | REAL(r_std), ALLOCATABLE, DIMENSION (:,:) :: lon, lat !! Lon/lat of coarser grid |
---|
| 52 | REAL(r_std), ALLOCATABLE, DIMENSION (:,:) :: orog !! New orography after interpolation |
---|
| 53 | REAL(r_std), ALLOCATABLE, DIMENSION (:) :: orog_land, orog_loc |
---|
| 54 | REAL(r_std), ALLOCATABLE, DIMENSION (:,:) :: lalo_land |
---|
| 55 | REAL(r_std), ALLOCATABLE, DIMENSION (:) :: contfrac_land, contfrac_loc |
---|
| 56 | REAL(r_std), ALLOCATABLE, DIMENSION (:,:) :: contfrac_2d |
---|
| 57 | ! |
---|
| 58 | REAL(r_std), DIMENSION (1) :: lev !! Number of level (requested by restini routine) (unitless) |
---|
| 59 | CHARACTER(LEN=80) :: histname !! Name of history file (can not find HISTNAME?) |
---|
| 60 | INTEGER(i_std) :: hori_id !! ID of the default horizontal longitude and latitude map. |
---|
| 61 | INTEGER(i_std) :: hist_id !! History file identification for ??? |
---|
| 62 | INTEGER(i_std) :: rest_id !! ID of the restart file |
---|
| 63 | ! |
---|
| 64 | REAL(r_std) :: date0 !! Initial date |
---|
| 65 | REAL(r_std) :: date0_rest !! Initial date from restart file |
---|
| 66 | REAL(r_std) :: date !! Current date |
---|
| 67 | REAL(r_std) :: dt !! Same as dtradia ??? |
---|
| 68 | REAL(r_std) :: dw !! 86400. ??? |
---|
| 69 | REAL(r_std) :: one_day_loc |
---|
| 70 | ! |
---|
| 71 | CHARACTER(LEN=40) :: flux_op !! Operations to be performed on fluxes |
---|
| 72 | CHARACTER(LEN=40) :: flux_scinsec !! Operation in seconds |
---|
| 73 | CHARACTER(LEN=40) :: avescatter, once_wrt !! The various operation to be performed |
---|
| 74 | ! |
---|
| 75 | ! Input for routing_main |
---|
| 76 | ! |
---|
| 77 | INTEGER(i_std) :: kjit !! Time step number |
---|
| 78 | INTEGER(i_std) :: nbindex !! Number of local continental points |
---|
| 79 | REAL(r_std) :: dtradia !! Timestep length |
---|
| 80 | ! |
---|
| 81 | INTEGER(i_std), ALLOCATABLE, DIMENSION (:) :: kindex_g !! Index of land point on 2D map (in local position) |
---|
| 82 | INTEGER(i_std), ALLOCATABLE, DIMENSION (:) :: kindex !! index of land point per proc |
---|
| 83 | REAL(r_std), ALLOCATABLE, DIMENSION (:) :: runoff !! Grid-point runoff (kg/m^2/dt) |
---|
| 84 | REAL(r_std), ALLOCATABLE, DIMENSION (:) :: drainage !! Grid-point drainage (kg/m^2/dt) |
---|
| 85 | REAL(r_std), ALLOCATABLE, DIMENSION (:) :: humrel !! Soil moisture stress, root extraction potential (unitless) |
---|
| 86 | REAL(r_std), ALLOCATABLE, DIMENSION (:,:) :: stempdiag !! Diagnostic soil temperature profile |
---|
| 87 | REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:) :: totfrac_nobio !! Total fraction of continental ice+lakes+cities+... |
---|
| 88 | REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:) :: veget_max !! Max. fraction of vegetation type (LAI -> infty, unitless) |
---|
| 89 | REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:) :: floodout !! Flow out of floodplains from hydrol |
---|
| 90 | REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:) :: transpot !! Potential Transpiration (needed for irrigation) |
---|
| 91 | REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:) :: precip_rain !! Rainfall (kg/m^2/dt) |
---|
| 92 | REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:) :: k_litt !! Averaged conductivity for saturated infiltration in the 'litter' layer (kg/m^2/dt) |
---|
| 93 | REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:) :: reinf_slope !! Coefficient which determines the reinfiltration ratio in the grid box due to flat areas (unitless;0-1) |
---|
| 94 | INTEGER(i_std) :: hist2_id !! Access to history file 2 (unitless) |
---|
| 95 | ! |
---|
| 96 | REAL(r_std), ALLOCATABLE, DIMENSION (:) :: evapot_corr !! Soil Potential Evaporation |
---|
| 97 | ! |
---|
| 98 | ! Output from routing_main |
---|
| 99 | ! |
---|
| 100 | REAL(r_std), ALLOCATABLE, DIMENSION (:) :: returnflow !! The water flow from lakes and swamps which returns to the grid box (kg/m^2/dt) |
---|
| 101 | REAL(r_std), ALLOCATABLE, DIMENSION (:) :: irrigation !! This is the water taken from the reservoirs and beeing put into the upper layers of the soil (kg/m^2/dt) |
---|
| 102 | REAL(r_std), ALLOCATABLE, DIMENSION (:) :: riverflow !! Outflow of the major rivers, will be located on the continental grid but this should be a coastal point (kg/dt) |
---|
| 103 | REAL(r_std), ALLOCATABLE, DIMENSION (:) :: coastalflow !! Outflow on coastal points by small basins, the water which flows in a disperse way into the ocean (kg/dt) |
---|
| 104 | REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:) :: reinfiltration !! Water flow from ponds and floodplains which returns to the grid box (kg/m^2/dt) |
---|
| 105 | REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:) :: flood_res !! Diagnostic of water amount in the floodplains reservoir (kg) |
---|
| 106 | REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:) :: flood_frac !! Flooded fraction of the grid box (unitless;0-1) |
---|
| 107 | ! |
---|
| 108 | ! |
---|
| 109 | !_ ================================================================================================================================ |
---|
| 110 | ! |
---|
| 111 | ! |
---|
| 112 | CALL Init_orchidee_para() |
---|
| 113 | ! |
---|
| 114 | CALL getlandseamask_init(iim, jjm, nbindex) |
---|
| 115 | ALLOCATE(lon(iim,jjm)) |
---|
| 116 | ALLOCATE(lat(iim,jjm)) |
---|
| 117 | ALLOCATE(orog(iim,jjm)) |
---|
| 118 | ALLOCATE(contfrac_2d(iim,jjm)) |
---|
| 119 | CALL getlandseamask_read(lon, lat, contfrac_2d, orog) |
---|
| 120 | ! |
---|
| 121 | ! ALLOCATE memory needed |
---|
| 122 | ! |
---|
| 123 | ALLOCATE(kindex_g(nbindex)) |
---|
| 124 | ALLOCATE(lalo_land(nbindex,2)) |
---|
| 125 | ALLOCATE(contfrac_land(nbindex)) |
---|
| 126 | ALLOCATE(orog_land(nbindex)) |
---|
| 127 | ! |
---|
| 128 | ! |
---|
| 129 | ! |
---|
| 130 | ni=0 |
---|
| 131 | DO j=1,jjm |
---|
| 132 | DO i=1,iim |
---|
| 133 | IF ( contfrac_2d(i,j) > 0.0 ) THEN |
---|
| 134 | ni = ni + 1 |
---|
| 135 | IF ( ni .GT. nbindex ) THEN |
---|
| 136 | WRITE(*,*) "We are expecting ", nbindex, "point." |
---|
| 137 | WRITE(*,*) "We are at : ", i, j, orog(i,j) |
---|
| 138 | STOP 'Too many continental points' |
---|
| 139 | ENDIF |
---|
| 140 | kindex_g(ni) = (j-1)*iim + i |
---|
| 141 | lalo_land(ni,1) = lat(i,j) |
---|
| 142 | lalo_land(ni,2) = lon(i,j) |
---|
| 143 | contfrac_land(ni) = contfrac_2d(i,j) |
---|
| 144 | orog_land(ni) = orog(i,j) |
---|
| 145 | ENDIF |
---|
| 146 | ENDDO |
---|
| 147 | ENDDO |
---|
| 148 | ! |
---|
| 149 | ! |
---|
| 150 | nbseg = 4 |
---|
| 151 | ! |
---|
| 152 | ! |
---|
| 153 | CALL grid_set_glo(iim, jjm, nbindex) |
---|
| 154 | CALL grid_allocate_glo(nbseg) |
---|
| 155 | ! |
---|
| 156 | CALL bcast(nbindex) |
---|
| 157 | ALLOCATE(index_g(nbindex)) |
---|
| 158 | IF ( is_root_prc ) index_g(:) = kindex_g(:) |
---|
| 159 | CALL bcast(index_g) |
---|
| 160 | ! |
---|
| 161 | WRITE(*,*) "GOING INTO Init_orchidee_data_para_driver", nbindex, index_g(1), SIZE(kindex_g) |
---|
| 162 | CALL Init_orchidee_data_para_driver(nbindex, index_g) |
---|
| 163 | WRITE(*,*) "OUT OF Init_orchidee_data_para_driver" |
---|
| 164 | CALL init_ioipsl_para |
---|
| 165 | ! |
---|
| 166 | WRITE(*,*) mpi_rank, "DIMENSIONS of grid on processor : iim, jjm, nbindex = ", iim, jjm, nbindex, nbp_loc |
---|
| 167 | ! |
---|
| 168 | CALL grid_init (nbp_loc, nbseg, regular_lonlat, "ForcingGrid") |
---|
| 169 | ! |
---|
| 170 | !========================================================================== |
---|
| 171 | ! |
---|
| 172 | ! Transfer the global grid variables to the root proc |
---|
| 173 | ! *_glo -> *_g |
---|
| 174 | ! Variables *_g were allocated with the CALL init_grid |
---|
| 175 | ! |
---|
| 176 | IF ( is_root_prc) THEN |
---|
| 177 | ! |
---|
| 178 | lalo_g(:,:) = lalo_land(:,:) |
---|
| 179 | lon_g(:,:) = lon(:,:) |
---|
| 180 | lat_g(:,:) = lat(:,:) |
---|
| 181 | contfrac_g(:) = contfrac_land(:) |
---|
| 182 | ! |
---|
| 183 | ENDIF |
---|
| 184 | ! |
---|
| 185 | CALL grid_stuff(nbindex, iim, jjm, lon_g, lat_g, index_g) |
---|
| 186 | ! |
---|
| 187 | ! |
---|
| 188 | ! Distribute the grid to all processors |
---|
| 189 | ! |
---|
| 190 | ! Redistribute the indeces on all procs (apple distribution of land points) |
---|
| 191 | ! |
---|
| 192 | ALLOCATE(kindex(nbp_loc)) |
---|
| 193 | ALLOCATE(orog_loc(nbp_loc), contfrac_loc(nbp_loc)) |
---|
| 194 | CALL bcast(lon_g) |
---|
| 195 | CALL bcast(lat_g) |
---|
| 196 | CALL scatter(index_g, kindex) |
---|
| 197 | CALL scatter(lalo_land, lalo) |
---|
| 198 | CALL scatter(orog_land, orog_loc) |
---|
| 199 | CALL scatter(contfrac_land, contfrac_loc) |
---|
| 200 | ! |
---|
| 201 | ! |
---|
| 202 | ! Apply the offset needed so that kindex refers to the index of the land point |
---|
| 203 | ! on the current region, i.e. the local lon lat domain. |
---|
| 204 | ! |
---|
| 205 | kindex(1:nbp_loc)=kindex(1:nbp_loc)-(jj_begin-1)*iim |
---|
| 206 | ! |
---|
| 207 | ! |
---|
| 208 | !========================================================================================== |
---|
| 209 | ! |
---|
| 210 | ! The grid is in place and we can start to prepare the time of integration. |
---|
| 211 | ! |
---|
| 212 | CALL ioconf_calendar("gregorian") |
---|
| 213 | ! |
---|
| 214 | ! |
---|
| 215 | ! Determine initial step or restart one |
---|
| 216 | ! |
---|
| 217 | !Config Key = RESTART_IN |
---|
| 218 | !Config Desc = Name of restart file to read at restart |
---|
| 219 | !Config If = [-] |
---|
| 220 | !Config Def = NONE |
---|
| 221 | !Config Help = This function allows to select a restart file with which |
---|
| 222 | ! the simulation will be initialized. |
---|
| 223 | !Config Units = [-] |
---|
| 224 | !- |
---|
| 225 | ! |
---|
| 226 | restname_in="NONE" |
---|
| 227 | CALL getin('RESTART_IN', restname_in) |
---|
| 228 | ! |
---|
| 229 | !Config Key = RESTART_OUT |
---|
| 230 | !Config Desc = Name of restart file to be written at the end of the simulation. |
---|
| 231 | !Config If = [-] |
---|
| 232 | !Config Def = NONE |
---|
| 233 | !Config Help = This function allows to select a restart file which will be written by the |
---|
| 234 | ! model and which can be used as input for a future restart. |
---|
| 235 | !Config Units = [-] |
---|
| 236 | !- |
---|
| 237 | ! |
---|
| 238 | restname_out="restart_out.nc" |
---|
| 239 | CALL getin('RESTART_OUT', restname_out) |
---|
| 240 | ! |
---|
| 241 | !Config Key = SIMULATION_LEN |
---|
| 242 | !Config Desc = Time step length in years for "testrouting" |
---|
| 243 | !Config If = [-] |
---|
| 244 | !Config Def = 1 |
---|
| 245 | !Config Help = This is time step length for testrouting |
---|
| 246 | !Config Units = [-] |
---|
| 247 | nbyears=1.0 |
---|
| 248 | CALL getin('SIMULATION_LEN', nbyears) |
---|
| 249 | ! |
---|
| 250 | !Config Key = DTRADIA |
---|
| 251 | !Config Desc = Time step length for "testrouting" |
---|
| 252 | !Config If = [-] |
---|
| 253 | !Config Def = 1800. |
---|
| 254 | !Config Help = This is time step length for testrouting |
---|
| 255 | !Config Units = [-] |
---|
| 256 | !- |
---|
| 257 | !DTRADIA = 1800. |
---|
| 258 | dtradia = 1800. |
---|
| 259 | CALL getin('DTRADIA', dtradia) |
---|
| 260 | ! |
---|
| 261 | !- Initial date |
---|
| 262 | CALL ymds2ju (2000,1,1,0.0, date0) |
---|
| 263 | date0_rest = date0 |
---|
| 264 | ! |
---|
| 265 | dt = dtradia |
---|
| 266 | dw = 86400. |
---|
| 267 | ! |
---|
| 268 | CALL control_initialize |
---|
| 269 | ! |
---|
| 270 | CALL ioget_calendar(one_year,one_day_loc) |
---|
| 271 | ! |
---|
| 272 | ! We have all we need and we can start to work |
---|
| 273 | ! |
---|
| 274 | ! |
---|
| 275 | IF (is_root_prc) THEN |
---|
| 276 | CALL restini(restname_in, iim, jjm, lon_g, lat_g, 1, lev, & |
---|
| 277 | & restname_out, ibegt, date0_rest, dtradia, rest_id, .FALSE.) |
---|
| 278 | ELSE |
---|
| 279 | rest_id=0 |
---|
| 280 | ENDIF |
---|
| 281 | CALL bcast (ibegt) |
---|
| 282 | CALL bcast (date0_rest) |
---|
| 283 | CALL bcast (dtradia) |
---|
| 284 | ! |
---|
| 285 | ! |
---|
| 286 | IF ( INDEX(restname_in, "NONE") > 0 ) THEN |
---|
| 287 | kjit = 1 |
---|
| 288 | ibegt = 1 |
---|
| 289 | ELSE |
---|
| 290 | kjit = ibegt |
---|
| 291 | date0 = date0_rest |
---|
| 292 | ENDIF |
---|
| 293 | WRITE(*,*) 'Out of restini : kjit=',kjit, " ibegt=", ibegt, " date0=", date0 |
---|
| 294 | ! |
---|
| 295 | !- time step length |
---|
| 296 | ! |
---|
| 297 | ! Set up the history file |
---|
| 298 | ! |
---|
| 299 | !Config Key = HISTNAME |
---|
| 300 | !Config Desc = Name of the history file |
---|
| 301 | !Config If = [-] |
---|
| 302 | !Config Def = out_testrouting |
---|
| 303 | !Config Help = The name of the file which will contain all the diagnostics of routing. |
---|
| 304 | !Config Units = [-] |
---|
| 305 | !- |
---|
| 306 | WRITE(flux_op,'("ave(scatter(X*",F8.1,"))")') one_day_loc/dt |
---|
| 307 | WRITE(flux_scinsec,'("ave(scatter(X*",F8.6,"))")') 1.0/dt |
---|
| 308 | avescatter = 'ave(scatter(X))' |
---|
| 309 | once_wrt = 'once(scatter(X))' |
---|
| 310 | ! |
---|
| 311 | histname="out_testrouting" |
---|
| 312 | CALL getin('HISTNAME', histname) |
---|
| 313 | ! |
---|
| 314 | CALL histbeg(histname, iim, lon, jjm, lat, 1, iim, 1, jjm, & |
---|
| 315 | & kjit-1, date0, dtradia, hori_id, hist_id, domain_id=orch_domain_id) |
---|
| 316 | ! |
---|
| 317 | CALL histdef(hist_id, 'Orog', 'Orography', ' ', & |
---|
| 318 | & iim, jjm, hori_id, 1,1,1, -99, 32, once_wrt, dt, dw) |
---|
| 319 | CALL histdef(hist_id, 'Contfrac', 'Fraction of continent', ' ', & |
---|
| 320 | & iim, jjm, hori_id, 1,1,1, -99, 32, once_wrt, dt, dw) |
---|
| 321 | CALL histdef(hist_id, 'Areas', 'Mesh areas', 'm2', & |
---|
| 322 | & iim,jjm, hori_id, 1,1,1, -99, 32, once_wrt, dt, dw) |
---|
| 323 | ! |
---|
| 324 | CALL histdef(hist_id, 'riversret', 'Return from endorheic rivers', 'mm/d', & |
---|
| 325 | & iim,jjm, hori_id, 1,1,1, -99, 32, flux_op, dt,dw) |
---|
| 326 | CALL histdef(hist_id, 'hydrographs', 'Hydrographs of gridbox outflow', 'm^3/s', & |
---|
| 327 | & iim,jjm, hori_id, 1,1,1, -99, 32, flux_scinsec, dt,dw) |
---|
| 328 | ! |
---|
| 329 | CALL histdef(hist_id, 'fastr', 'Fast flow reservoir', 'kg/m^2', & |
---|
| 330 | & iim,jjm, hori_id, 1,1,1, -99, 32, avescatter, dt,dw) |
---|
| 331 | CALL histdef(hist_id, 'slowr', 'Slow flow reservoir', 'kg/m^2', & |
---|
| 332 | & iim,jjm, hori_id, 1,1,1, -99, 32, avescatter, dt,dw) |
---|
| 333 | CALL histdef(hist_id, 'streamr', 'Stream flow reservoir', 'kg/m^2', & |
---|
| 334 | & iim,jjm, hori_id, 1,1,1, -99, 32, avescatter, dt,dw) |
---|
| 335 | CALL histdef(hist_id, 'pondr', 'Volume in pond reservoir', 'kg/m^2', & |
---|
| 336 | & iim,jjm, hori_id, 1,1,1, -99, 32, avescatter, dt,dw) |
---|
| 337 | CALL histdef(hist_id, 'lakevol', 'Volume in lake reservoir', 'kg/m^2', & |
---|
| 338 | & iim,jjm, hori_id, 1,1,1, -99, 32, avescatter, dt,dw) |
---|
| 339 | ! |
---|
| 340 | CALL histdef(hist_id, 'basinmap', 'Aproximate map of the river basins', ' ', & |
---|
| 341 | & iim,jjm, hori_id, 1,1,1, -99, 32, avescatter, dt,dw) |
---|
| 342 | CALL histdef(hist_id, 'nbrivers', 'Number or rivers in the outflow grid box', ' ', & |
---|
| 343 | & iim,jjm, hori_id, 1,1,1, -99, 32, avescatter, dt,dw) |
---|
| 344 | ! |
---|
| 345 | CALL histend(hist_id) |
---|
| 346 | ! |
---|
| 347 | ! Put a copy of the orography into the restart |
---|
| 348 | ! |
---|
| 349 | ! |
---|
| 350 | CALL histwrite_p(hist_id, 'Orog', kjit+1, orog_loc, nbp_loc, kindex) |
---|
| 351 | CALL histwrite_p(hist_id, 'Contfrac', kjit+1, contfrac_loc, nbp_loc, kindex) |
---|
| 352 | CALL histwrite_p(hist_id, 'Areas', kjit+1, area, nbp_loc, kindex) |
---|
| 353 | ! |
---|
| 354 | ! Override some settings as testrouting is not as flexible as the rest of ORCHIDEE. |
---|
| 355 | ! |
---|
| 356 | hist2_id=-1 |
---|
| 357 | almaoutput=.FALSE. |
---|
| 358 | !================================================================================ |
---|
| 359 | ! |
---|
| 360 | ! Set up the routing schemes |
---|
| 361 | ! |
---|
| 362 | ! |
---|
| 363 | ! Allocate all the physical variables |
---|
| 364 | ! |
---|
| 365 | ! Input variables |
---|
| 366 | ALLOCATE(runoff(nbp_loc), drainage(nbp_loc), humrel(nbp_loc), stempdiag(nbp_loc,nslm), transpot(nbp_loc,nvmc)) |
---|
| 367 | ALLOCATE(totfrac_nobio(nbp_loc), veget_max(nbp_loc,nvmc), floodout(nbp_loc)) |
---|
| 368 | ALLOCATE(precip_rain(nbp_loc), k_litt(nbp_loc),reinf_slope(nbp_loc)) |
---|
| 369 | ALLOCATE(evapot_corr(nbp_loc)) |
---|
| 370 | ! Output variables |
---|
| 371 | ALLOCATE(returnflow(nbp_loc), irrigation(nbp_loc), riverflow(nbp_loc), coastalflow(nbp_loc), reinfiltration(nbp_loc)) |
---|
| 372 | ALLOCATE(flood_frac(nbp_loc), flood_res(nbp_loc)) |
---|
| 373 | ! |
---|
| 374 | ! Get some fake value for input arrays |
---|
| 375 | ! |
---|
| 376 | runoff(:) = 1.0 |
---|
| 377 | drainage(:) = 1.0 |
---|
| 378 | humrel(:) = 0.75 |
---|
| 379 | stempdiag(:,:) = 273.5 |
---|
| 380 | transpot(:,:)=0.0 |
---|
| 381 | reinfiltration(:)=0.0 |
---|
| 382 | flood_frac(:)=0.0 |
---|
| 383 | flood_res(:)=0.0 |
---|
| 384 | totfrac_nobio = 0.1 |
---|
| 385 | veget_max = 0.2 |
---|
| 386 | floodout = 0.0 |
---|
| 387 | precip_rain = 0.0 |
---|
| 388 | k_litt = 0.0 |
---|
| 389 | reinf_slope = 0.1 |
---|
| 390 | evapot_corr(:) = 10.0 |
---|
| 391 | ! |
---|
| 392 | ! |
---|
| 393 | ! |
---|
| 394 | CALL routing_initialize( kjit, nbp_loc, kindex, & |
---|
| 395 | & rest_id, hist_id, hist2_id, lalo, & |
---|
| 396 | & neighbours, resolution, contfrac, stempdiag, & |
---|
| 397 | & returnflow, reinfiltration, irrigation, riverflow, & |
---|
| 398 | & coastalflow, flood_frac, flood_res) |
---|
| 399 | ! |
---|
| 400 | ! Do loop over a number of time-steps |
---|
| 401 | ! |
---|
| 402 | simlen = NINT(nbyears*365*one_day_loc/dtradia) |
---|
| 403 | ibegt=kjit |
---|
| 404 | iendt=kjit+simlen |
---|
| 405 | WRITE(*,*) "The simulation will go from ", ibegt, " to ", iendt |
---|
| 406 | ! |
---|
| 407 | DO kjit = ibegt,iendt |
---|
| 408 | ! |
---|
| 409 | date = date0 + (kjit-1)*(dtradia/one_day_loc) |
---|
| 410 | ! |
---|
| 411 | IF ( date < date0+1 ) THEN |
---|
| 412 | ! During one day one kg/m^2d divided up in runoff and drainage |
---|
| 413 | runoff(:) = 0.5/48. |
---|
| 414 | drainage(:) = 0.5/48. |
---|
| 415 | ELSE |
---|
| 416 | runoff(:) = 0.0 |
---|
| 417 | drainage(:) = 0.0 |
---|
| 418 | ENDIF |
---|
| 419 | ! |
---|
| 420 | CALL routing_main(kjit, nbp_loc, kindex, & |
---|
| 421 | & lalo, neighbours, resolution, contfrac, totfrac_nobio, veget_max, floodout, runoff, & |
---|
| 422 | & drainage, transpot, precip_rain, humrel, k_litt, flood_frac, flood_res, & |
---|
| 423 | & stempdiag, reinf_slope, returnflow, reinfiltration, irrigation, riverflow, coastalflow, rest_id, hist_id, hist2_id) |
---|
| 424 | ! |
---|
| 425 | WRITE(*,*) "Out of routing at time step = ",kjit,' Seconds since start', (kjit-ibegt)*dtradia |
---|
| 426 | ! |
---|
| 427 | ENDDO |
---|
| 428 | ! |
---|
| 429 | ! Shut everything down |
---|
| 430 | ! |
---|
| 431 | CALL routing_finalize(kjit, nbp_loc, rest_id, flood_frac, flood_res) |
---|
| 432 | ! |
---|
| 433 | ! |
---|
| 434 | CALL histclo |
---|
| 435 | IF ( is_root_prc) THEN |
---|
| 436 | CALL restclo |
---|
| 437 | ENDIF |
---|
| 438 | ! |
---|
| 439 | CALL Finalize_mpi |
---|
| 440 | ! |
---|
| 441 | END PROGRAM testrouting |
---|
| 442 | ! |
---|