Opened 8 years ago

Closed 8 years ago

#102 closed defect (fixed)

Negative drainage create problems in the routing

Reported by: jpolcher Owned by: somebody
Priority: major Milestone:
Component: Anthropogenic processes Version:
Keywords: Cc:

Description

In its current version, CWRR can produce negative drainage. This water has to be taken from the routing scheme ... if it is activated.

While we solve the problem in hydrol.f90, I propose to transform a fatal error in the routing when the slow reservoir become negative into only a warning. It has been verified that this causes no problems elsewhere in routing.f90.

So in some simulations we might have negative slow reservoirs in the routing. It will print warning but should not have other effects. In any case water will be conserved.

The proposed correction also fixes a small problem which occurs when irrigation is activated over a limited domain. The model was looking for extra water outside of the domain !

Attachments (1)

diff_routing.f90 (4.6 KB) - added by jpolcher 8 years ago.
A version of the diff easier to read

Download all attachments as: .zip

Change History (3)

comment:1 Changed 8 years ago by jpolcher

Proposed diff :

Index: routing.f90
===================================================================
--- routing.f90 (revision 1296)
+++ routing.f90 (working copy)
@@ -2011,7 +2011,7 @@

!> Water that flows into the pond_reservoir is withdrawn from the fast_reservoir.

!
error_reservoir=.FALSE.

  • loopresnbpt : DO ig=1,nbpt

+ DO ig=1,nbpt

DO ib=1,nbasmax

!
fast_reservoir(ig,ib) = fast_reservoir(ig,ib) + runoff(ig)*routing_area(ig,ib) - &

@@ -2061,19 +2061,18 @@

!
IF ( slow_reservoir(ig,ib) .LT. zero ) THEN

error_reservoir=.TRUE.

  • EXIT loopresnbpt

ENDIF
!

ENDDO

  • ENDDO loopresnbpt

+ ENDDO

IF (error_reservoir) THEN

  • WRITE(numout,*) 'There is a negative reservoir at :', ig, ib
  • WRITE(numout,*) 'slowr, slow_flow, drainage', &

+ WRITE(numout,*) 'WARNING : There is a negative reservoir at :', ig, ib
+ WRITE(numout,*) 'WARNING : slowr, slow_flow, drainage', &

& slow_reservoir(ig,ib), slow_flow(ig,ib), drainage(ig)

  • WRITE(numout,*) 'pondr, pond_inflow, pond_drainage', &

+ WRITE(numout,*) 'WARNING : pondr, pond_inflow, pond_drainage', &

& pond_reservoir(ig), pond_inflow(ig,ib), pond_drainage(ig,ib)

  • CALL ipslerr_p(3, 'routing_flow', 'ERROR, slow_reservoir can be safely be reset to zero.',,)

+ CALL ipslerr_p(2, 'routing_flow', 'ERROR, slow_reservoir cannot be safely be reset to zero.',,)

ENDIF


totflood(:) = zero

@@ -2258,21 +2257,29 @@

!
IF ( irrig_deficit_glo(ig,ib) > min_sechiba ) THEN

!

+ streams_around(:,:) = zero
+ !

DO in=1,8

ig2 = neighbours_g(ig,in)

  • streams_around(in,:) = stream_reservoir_glo(ig2,:)
  • igrd(in) = ig2

+ IF (ig2 .GT. 0 ) THEN
+ streams_around(in,:) = stream_reservoir_glo(ig2,:)
+ igrd(in) = ig2
+ ENDIF

ENDDO
!

  • ff=MAXLOC(streams_around)
  • ig2=igrd(ff(1))
  • ib2=ff(2)
  • !
  • IF ( routing_area_glo(ig2,ib2) .GT. 0 .AND. stream_reservoir_glo(ig2,ib2) > zero ) THEN
  • adduction = MIN(irrig_deficit_glo(ig,ib), stream_reservoir_glo(ig2,ib2))
  • stream_reservoir_glo(ig2,ib2) = stream_reservoir_glo(ig2,ib2) - adduction
  • irrig_deficit_glo(ig,ib) = irrig_deficit_glo(ig,ib) - adduction
  • irrig_adduct_glo(ig,ib) = irrig_adduct_glo(ig,ib) + adduction

+ IF ( MAXVAL(streams_around) .GT. zero ) THEN
+ !
+ ff=MAXLOC(streams_around)
+ ig2=igrd(ff(1))
+ ib2=ff(2)
+ !
+ IF ( routing_area_glo(ig2,ib2) .GT. 0 .AND. stream_reservoir_glo(ig2,ib2) > zero ) THEN
+ adduction = MIN(irrig_deficit_glo(ig,ib), stream_reservoir_glo(ig2,ib2))
+ stream_reservoir_glo(ig2,ib2) = stream_reservoir_glo(ig2,ib2) - adduction
+ irrig_deficit_glo(ig,ib) = irrig_deficit_glo(ig,ib) - adduction
+ irrig_adduct_glo(ig,ib) = irrig_adduct_glo(ig,ib) + adduction
+ ENDIF
+ !

ENDIF
!

ENDIF

@@ -4135,8 +4142,8 @@

callsign = "routing_basins"
ok_interpol = .FALSE.
!

  • nix=INT(MAXVAL(resolution_g(:,1))/MAXVAL(resol_lu(:,:,1)))+2
  • njx=INT(MAXVAL(resolution_g(:,2))/MAXVAL(resol_lu(:,:,2)))+2

+ nix=INT(MAXVAL(resolution_g(:,1))/MINVAL(resol_lu(:,:,1)))+2
+ njx=INT(MAXVAL(resolution_g(:,2))/MINVAL(resol_lu(:,:,2)))+2

nbvmax = nix*njx*2
!
! We are on the root processor here as this routine is not in parallel. So no need to broadcast.

@@ -7712,9 +7719,9 @@

! Some lmargin is taken.
!
IF (is_root_prc) THEN

  • nix=INT(MAXVAL(resolution_g(:,1))/MAXVAL(resol_lu(:,:,1)))+2
  • njx=INT(MAXVAL(resolution_g(:,2))/MAXVAL(resol_lu(:,:,2)))+2
  • nbpmax = nix*njx

+ nix=INT(MAXVAL(resolution_g(:,1))/MINVAL(resol_lu(:,:,1)))+2
+ njx=INT(MAXVAL(resolution_g(:,2))/MINVAL(resol_lu(:,:,2)))+2
+ nbpmax = nix*njx*2

ENDIF
CALL bcast(nbpmax)
!

Changed 8 years ago by jpolcher

A version of the diff easier to read

comment:2 Changed 8 years ago by jgipsl

  • Resolution set to fixed
  • Status changed from new to closed

Correction done by Jan Polcher commited on the trunk in revison 1323.

Note: See TracTickets for help on using tickets.