Changeset 97

Timestamp:

08/04/12 00:13:44 (12 years ago)

Author:

ymipsl

Message:

Put conforming IO for dcmip testcase

YM

Location:

codes/icosagcm/trunk/src

Files:

• caldyn_gcm.f90 (modified) (4 diffs)
• icosa_gcm.f90 (modified) (3 diffs)
• time.f90 (modified) (4 diffs)
• timeloop_gcm.f90 (modified) (1 diff)

codes/icosagcm/trunk/src/caldyn_gcm.f90

@@ -643 +643 @@
     USE omega_mod
     USE write_field
+    USE vertical_interp_mod
     TYPE(t_field),POINTER :: f_ps(:), f_phis(:), f_u(:), f_theta_rhodz(:), f_dps(:), &
          f_buf_i(:), f_buf_v(:), f_buf_i3(:), f_buf1_i(:), f_buf2_i(:), f_buf_s(:), f_buf_p(:)
+    
+    REAL(rstd) :: out_pression_lev
+    CHARACTER(LEN=255) :: str_pression
+    
+    out_pression_level=0
+    CALL getin("out_pression_level",out_pression_level) 
+    WRITE(str_pression,*) INT(out_pression_level/100)
+    str_pression=ADJUSTL(str_pression)
     
     CALL writefield("ps",f_ps)
@@ -654 +663 @@
     CALL w_omega(f_ps, f_u, f_buf_i)
     CALL writefield('omega', f_buf_i)
+    IF (out_pression_level<=preff) THEN
+      CALL vertical_interp(f_ps,f_buf_i,f_buf_s,out_pression_level)
+      CALL writefield("omega"//TRIM(str_pression),f_buf_s)
+    ENDIF
     
     ! Temperature
     CALL theta_rhodz2temperature(f_ps,f_theta_rhodz,f_buf_i) ;
     CALL writefield("T",f_buf_i)
+    
+    IF (out_pression_level<=preff) THEN
+      CALL vertical_interp(f_ps,f_buf_i,f_buf_s,out_pression_level)
+      CALL writefield("T"//TRIM(str_pression),f_buf_s)
+    ENDIF
     
     ! velocity components
@@ -663 +681 @@
     CALL writefield("ulon",f_buf1_i)
     CALL writefield("ulat",f_buf2_i)
+
+    IF (out_pression_level<=preff) THEN
+      CALL vertical_interp(f_ps,f_buf1_i,f_buf_s,out_pression_level)
+      CALL writefield("ulon"//TRIM(str_pression),f_buf_s)
+      CALL vertical_interp(f_ps,f_buf2_i,f_buf_s,out_pression_level)
+      CALL writefield("vlon"//TRIM(str_pression),f_buf_s)
+    ENDIF
     
     ! geopotential
@@ -670 +695 @@
     CALL writefield("theta",f_buf1_i) ! potential temperature
     CALL writefield("pk",f_buf2_i) ! Exner pressure
+  
   
   END SUBROUTINE write_output_fields

codes/icosagcm/trunk/src/icosa_gcm.f90

@@ -6 +6 @@
   USE wind_mod
   USE mpipara
+  USE vertical_interp_mod
   IMPLICIT NONE
   
@@ -48 +49 @@
   CALL compute_geometry
   CALL init_disvert 
-  
+  CALL init_vertical_interp
+
   CALL allocate_field(sum_ne,field_T,type_real)
   
@@ -79 +81 @@
 !  CALL WriteField("sum_ne",sum_ne)
   CALL write_apbp
-  CALL create_time_counter_header
+  CALL init_time
   CALL timeloop
 

codes/icosagcm/trunk/src/time.f90

@@ -1 +1 @@
 MODULE time_mod
-
+  use prec
   PRIVATE
 
@@ -7 +7 @@
   INTEGER,SAVE :: it
 
-  PUBLIC create_time_counter_header, update_time_counter, close_time_counter
+  REAL(rstd),SAVE :: dt
+  REAL(rstd),SAVE :: write_period
+  INTEGER,SAVE :: itau_out
+  INTEGER,SAVE :: itaumax
+  
+  PUBLIC create_time_counter_header, update_time_counter, close_time_counter, init_time
 
 
 
 CONTAINS
+  
+  SUBROUTINE init_time
+  USE earth_const
+  USE ioipsl
+  IMPLICIT NONE
+  REAL(rstd) :: run_length
+  
+    dt=90.
+    CALL getin('dt',dt)
+
+    itaumax=100
+    CALL getin('itaumax',itaumax)
+
+    run_length=dt*itaumax
+    CALL getin('run_length',run_length)
+    itaumax=run_length/dt
+    PRINT *,'itaumax=',itaumax
+    dt=dt/scale_factor
+
+    write_period=0
+    CALL getin('write_period',write_period)
+    write_period=write_period/scale_factor
+    itau_out=FLOOR(.5+write_period/dt)
+    PRINT *, 'Output frequency (scaled) set to ',write_period, ' : itau_out = ',itau_out
+    
+    CALL create_time_counter_header
+    
+  END SUBROUTINE init_time
 
   SUBROUTINE create_time_counter_header
@@ -21 +54 @@
   INTEGER :: timeid, dtid
   REAL(rstd) :: dt
+  CHARACTER(LEN=255) :: time_frequency
+  
     status = NF90_CREATE('time_counter.nc', NF90_CLOBBER, ncid)
     status = NF90_DEF_DIM(ncid,'time_counter',NF90_UNLIMITED,timeid)
@@ -28 +63 @@
     status = NF90_PUT_ATT(ncid,time_counter_id,"calendar","noleap")
     status = NF90_DEF_VAR(ncid,'mdt',NF90_DOUBLE,varid=dtid)
+
+    WRITE(time_frequency,*) write_period
+    PRINT*,TRIM(time_frequency)
+    status = NF90_PUT_ATT(ncid,NF90_GLOBAL,"time_frequency",TRIM(time_frequency)//"s")
+
     status = NF90_ENDDEF(ncid) 
 

codes/icosagcm/trunk/src/timeloop_gcm.f90

@@ -102 +102 @@
   
   CALL etat0(f_ps,f_phis,f_theta_rhodz,f_u, f_q)
+  CALL writefield("phis",f_phis,once=.TRUE.)
   CALL transfert_request(f_q,req_i1)