Changeset 413 for codes/icosagcm


Ignore:
Timestamp:
06/10/16 17:49:07 (8 years ago)
Author:
dubos
Message:

Improved output

Location:
codes/icosagcm/trunk/src
Files:
4 edited

Legend:

Unmodified
Added
Removed

  • codes/icosagcm/trunk/src/caldyn_gcm.f90

    r406 r413  
    3535    IF (is_master) PRINT *, 'caldyn_conserv=',def 
    3636 
     37    nqdyn=1 ! default value 
     38 
    3739    def='theta' 
    38     CALL getin('caldyn_thermo',def) 
     40    CALL getin('thermo',def) 
    3941    SELECT CASE(TRIM(def)) 
    4042    CASE('theta') 
    … …  
    5153       physics_thermo=thermo_fake_moist 
    5254    CASE('moist') 
    53        caldyn_thermo=thermo_moist 
     55       caldyn_thermo=thermo_moist_debug 
    5456       physics_thermo=thermo_moist 
     57       nqdyn = 2 
    5558    CASE DEFAULT 
    5659       IF (is_mpi_root) PRINT *,'Bad selector for variable caldyn_thermo : <', & 
    … …  
    5861       STOP 
    5962    END SELECT 
    60      
     63 
     64    IF(is_master) THEN 
     65       SELECT CASE(caldyn_thermo) 
     66       CASE(thermo_theta) 
     67          PRINT *, 'caldyn_thermo = thermo_theta' 
     68       CASE(thermo_entropy) 
     69          PRINT *, 'caldyn_thermo = thermo_entropy' 
     70       CASE(thermo_moist_debug) 
     71          PRINT *, 'caldyn_thermo = thermo_moist_debug' 
     72       CASE DEFAULT 
     73          STOP 
     74       END SELECT 
     75 
     76       SELECT CASE(physics_thermo) 
     77       CASE(thermo_dry) 
     78          PRINT *, 'physics_thermo = thermo_dry' 
     79       CASE(thermo_fake_moist) 
     80          PRINT *, 'physics_thermo = thermo_fake_moist' 
     81       CASE(thermo_moist) 
     82          PRINT *, 'physics_thermo = thermo_moist' 
     83       END SELECT 
     84 
     85       PRINT *, 'nqdyn =', nqdyn 
     86    END IF 
     87 
    6188    CALL allocate_caldyn 
    6289 

  • codes/icosagcm/trunk/src/earth_const.f90

    r406 r413  
    2020 
    2121  INTEGER, PARAMETER,PUBLIC :: thermo_theta=1, thermo_entropy=2, & 
    22        thermo_moist=3, thermo_dry=4, thermo_fake_moist=5  
     22       thermo_moist=3, thermo_moist_debug=10, thermo_dry=4, thermo_fake_moist=5  
    2323  INTEGER, PUBLIC :: caldyn_thermo, physics_thermo 
    2424  !$OMP THREADPRIVATE(caldyn_thermo)  

  • codes/icosagcm/trunk/src/observable.f90

    r404 r413  
    1616 
    1717  PUBLIC init_observable, write_output_fields_basic, f_theta 
    18   LOGICAL,SAVE :: first_output=.TRUE. 
     18 
    1919!$OMP THREADPRIVATE(first_output) 
    2020   
    … …  
    3636    CALL allocate_field(f_pmid,  field_t,type_real,llm,  name='pmid')       ! mid layer pressure 
    3737  END SUBROUTINE init_observable 
    38    
    39   SUBROUTINE write_output_fields_basic(f_phis, f_ps, f_mass, f_geopot, f_theta_rhodz, f_u, f_W, f_q) 
     38 
     39  SUBROUTINE write_output_fields_basic(init, f_phis, f_ps, f_mass, f_geopot, f_theta_rhodz, f_u, f_W, f_q) 
     40    USE xios 
     41    USE disvert_mod 
    4042    USE wind_mod 
    4143    USE output_field_mod 
    … …  
    4345    USE time_mod 
    4446    USE xios 
    45     USE disvert_mod 
    4647    USE earth_const 
    4748    USE pression_mod 
    4849    USE vertical_interp_mod 
    4950    USE theta2theta_rhodz_mod 
    50     USE wind_mod 
    5151    USE omega_mod 
    52      
    53     TYPE(t_field),POINTER :: f_phis(:), f_ps(:), f_mass(:), f_geopot(:), f_theta_rhodz(:), f_u(:), f_W(:), f_q(:) 
    54 !    IF (is_master) PRINT *,'CALL write_output_fields_basic' 
     52    LOGICAL, INTENT(IN) :: init 
     53    INTEGER :: l 
    5554    REAL :: scalar(1) 
    5655    REAL :: mid_ap(llm) 
    5756    REAL :: mid_bp(llm) 
    58     INTEGER :: l 
    59  
    60     IF (first_output) THEN 
    61       scalar(1)=dt 
    62       CALL xios_send_field("timestep", scalar) 
    63       scalar(1)=preff 
    64       CALL xios_send_field("preff", scalar) 
    65       CALL xios_send_field("ap",ap) 
    66       CALL xios_send_field("bp",bp) 
    67       DO l=1,llm 
    68         mid_ap(l)=(ap(l)+ap(l+1))/2 
    69         mid_bp(l)=(bp(l)+bp(l+1))/2 
    70       ENDDO 
    71       CALL xios_send_field("mid_ap",mid_ap) 
    72       CALL xios_send_field("mid_bp",mid_bp) 
    73  
    74       CALL output_field("phis",f_phis) 
    75          
    76       first_output=.FALSE. 
    77     ENDIF 
     57 
     58    TYPE(t_field),POINTER :: f_phis(:), f_ps(:), f_mass(:), f_geopot(:), f_theta_rhodz(:), f_u(:), f_W(:), f_q(:) 
     59!    IF (is_master) PRINT *,'CALL write_output_fields_basic' 
     60 
     61    IF(init) THEN 
     62       scalar(1)=dt 
     63       CALL xios_send_field("timestep", scalar) 
     64       scalar(1)=preff 
     65       CALL xios_send_field("preff", scalar) 
     66       CALL xios_send_field("ap",ap) 
     67       CALL xios_send_field("bp",bp) 
     68       DO l=1,llm 
     69          mid_ap(l)=(ap(l)+ap(l+1))/2 
     70          mid_bp(l)=(bp(l)+bp(l+1))/2 
     71       ENDDO 
     72       CALL xios_send_field("mid_ap",mid_ap) 
     73       CALL xios_send_field("mid_bp",mid_bp) 
     74 
     75       CALL output_field("phis",f_phis) 
     76       CALL output_field("Ai",geom%Ai)        
     77    END IF 
     78 
     79    CALL divide_by_mass(1, f_mass, f_theta_rhodz, f_buf_i) 
     80    IF(init) THEN 
     81       CALL output_field("theta_init",f_buf_i) 
     82    ELSE 
     83       CALL output_field("theta",f_buf_i) 
     84    END IF 
     85 
     86    IF(nqdyn>1) THEN 
     87       CALL divide_by_mass(2, f_mass, f_theta_rhodz, f_buf_i) 
     88       IF(init) THEN 
     89          CALL output_field("dyn_q_init",f_buf_i) 
     90       ELSE 
     91          CALL output_field("dyn_q",f_buf_i) 
     92       END IF 
     93    END IF 
     94 
     95    CALL theta_rhodz2temperature(f_ps,f_theta_rhodz,f_buf_i) ;  
     96    CALL Tv2T(f_buf_i,f_q,f_buf1_i) 
     97    IF(init) THEN 
     98       CALL output_field("temp_init",f_buf_i) 
     99    ELSE 
     100       CALL output_field("temp",f_buf_i) 
     101       CALL vertical_interp(f_ps,f_buf_i,f_buf_s,85000.) 
     102       CALL output_field("t850",f_buf_s) 
     103       CALL vertical_interp(f_ps,f_buf_i,f_buf_s,50000.) 
     104       CALL output_field("t500",f_buf_s) 
     105       CALL vertical_interp(f_ps,f_buf_i,f_buf_s,preff) 
     106       CALL output_field("SST",f_buf_s)        
     107    END IF 
    78108     
    79109    CALL progonostic_vel_to_horiz(f_geopot, f_ps, f_mass, f_u, f_W, f_buf_uh, f_buf_i) 
    80110    CALL transfert_request(f_buf_uh,req_e1_vect)  
    81     CALL output_field("uz",f_buf_i) 
    82     CALL vertical_interp(f_ps,f_buf_i,f_buf_s,85000.) 
    83     CALL output_field("w850",f_buf_s) 
    84     CALL vertical_interp(f_ps,f_buf_i,f_buf_s,50000.) 
    85     CALL output_field("w500",f_buf_s) 
    86      
    87      
    88111    CALL un2ulonlat(f_buf_uh, f_buf_ulon, f_buf_ulat) 
    89     CALL output_field("ulon",f_buf_ulon) 
    90     CALL output_field("ulat",f_buf_ulat) 
    91     CALL output_field("ps",f_ps) 
    92     CALL output_field("Ai",geom%Ai) 
    93  
    94     !       CALL output_field("dps",f_dps) 
    95     CALL output_field("mass",f_mass) 
    96     CALL output_field("geopot",f_geopot) 
    97     !       CALL output_field("dmass",f_dmass) 
    98     !       CALL output_field("vort",f_qv) 
    99      
    100      
    101     CALL theta_rhodz2temperature(f_ps,f_theta_rhodz,f_buf_i) ;  
    102     CALL Tv2T(f_buf_i,f_q,f_buf1_i)  
    103     CALL output_field("temp",f_buf_i) 
    104     CALL vertical_interp(f_ps,f_buf_i,f_buf_s,85000.) 
    105     CALL output_field("t850",f_buf_s) 
    106     CALL vertical_interp(f_ps,f_buf_i,f_buf_s,50000.) 
    107     CALL output_field("t500",f_buf_s) 
    108     CALL vertical_interp(f_ps,f_buf_i,f_buf_s,preff) 
    109     CALL output_field("SST",f_buf_s) 
    110  
    111  
    112     CALL extract_slice(f_theta_rhodz, f_buf_i,1) 
    113     CALL output_field("theta",f_buf_i) 
    114             
    115     !       CALL output_field("exner",f_pk) 
    116     !       CALL output_field("pv",f_qv) 
    117     CALL output_field("q",f_q) 
    118112    CALL pression_mid(f_ps, f_pmid) 
    119     CALL output_field("p",f_pmid) 
    120  
    121     CALL vertical_interp(f_ps,f_buf_ulon,f_buf_s,85000.) 
    122     CALL output_field("u850",f_buf_s) 
    123     CALL vertical_interp(f_ps,f_buf_ulon,f_buf_s,50000.) 
    124     CALL output_field("u500",f_buf_s) 
    125  
    126     CALL vertical_interp(f_ps,f_buf_ulat,f_buf_s,85000.) 
    127     CALL output_field("v850",f_buf_s) 
    128     CALL vertical_interp(f_ps,f_buf_ulat,f_buf_s,50000.) 
    129     CALL output_field("v500",f_buf_s) 
    130      
    131     CALL w_omega(f_ps, f_u, f_buf_i) 
    132     CALL output_field("omega",f_buf_i) 
    133     CALL vertical_interp(f_ps,f_buf_i,f_buf_s,85000.) 
    134     CALL output_field("omega850",f_buf_s) 
    135     CALL vertical_interp(f_ps,f_buf_i,f_buf_s,50000.) 
    136     CALL output_field("omega500",f_buf_s) 
    137  
    138  
    139      
     113    IF(init) THEN 
     114       CALL output_field("uz_init",f_buf_i) 
     115       CALL output_field("ulon_init",f_buf_ulon) 
     116       CALL output_field("ulat_init",f_buf_ulat) 
     117       CALL output_field("p_init",f_pmid) 
     118       CALL output_field("ps_init",f_ps) 
     119       CALL output_field("mass_init",f_mass) 
     120       CALL output_field("geopot_init",f_geopot) 
     121       CALL output_field("q_init",f_q) 
     122    ELSE 
     123       CALL output_field("uz",f_buf_i) 
     124       CALL output_field("ulon",f_buf_ulon) 
     125       CALL output_field("ulat",f_buf_ulat) 
     126       CALL output_field("p",f_pmid) 
     127       CALL output_field("ps",f_ps) 
     128       CALL output_field("mass",f_mass) 
     129       CALL output_field("geopot",f_geopot) 
     130       CALL output_field("q",f_q) 
     131 
     132       !       CALL output_field("exner",f_pk) 
     133       !       CALL output_field("pv",f_qv) 
     134        
     135       CALL vertical_interp(f_ps,f_buf_ulon,f_buf_s,85000.) 
     136       CALL output_field("u850",f_buf_s) 
     137       CALL vertical_interp(f_ps,f_buf_ulon,f_buf_s,50000.) 
     138       CALL output_field("u500",f_buf_s) 
     139        
     140       CALL vertical_interp(f_ps,f_buf_ulat,f_buf_s,85000.) 
     141       CALL output_field("v850",f_buf_s) 
     142       CALL vertical_interp(f_ps,f_buf_ulat,f_buf_s,50000.) 
     143       CALL output_field("v500",f_buf_s) 
     144 
     145       CALL vertical_interp(f_ps,f_buf_i,f_buf_s,85000.) 
     146       CALL output_field("w850",f_buf_s) 
     147       CALL vertical_interp(f_ps,f_buf_i,f_buf_s,50000.) 
     148       CALL output_field("w500",f_buf_s)     
     149 
     150       CALL w_omega(f_ps, f_u, f_buf_i) 
     151       CALL output_field("omega",f_buf_i) 
     152       CALL vertical_interp(f_ps,f_buf_i,f_buf_s,85000.) 
     153       CALL output_field("omega850",f_buf_s) 
     154       CALL vertical_interp(f_ps,f_buf_i,f_buf_s,50000.) 
     155       CALL output_field("omega500",f_buf_s) 
     156    END IF 
     157 
    140158  END SUBROUTINE write_output_fields_basic 
    141159   
    … …  
    344362     
    345363  END SUBROUTINE Tv2T 
    346    
     364 
     365  SUBROUTINE divide_by_mass(iq, f_mass, f_theta_rhodz, f_theta) 
     366    INTEGER, INTENT(IN) :: iq 
     367    TYPE(t_field), POINTER :: f_mass(:), f_theta_rhodz(:), f_theta(:) 
     368    REAL(rstd), POINTER :: mass(:,:), theta_rhodz(:,:,:), theta(:,:) 
     369    INTEGER :: ind 
     370    DO ind=1,ndomain 
     371       IF (.NOT. assigned_domain(ind)) CYCLE 
     372       CALL swap_dimensions(ind) 
     373       CALL swap_geometry(ind) 
     374       mass=f_mass(ind) 
     375       theta_rhodz=f_theta_rhodz(ind) 
     376       theta=f_theta(ind) 
     377       theta(:,:) = theta_rhodz(:,:,iq) / mass(:,:) 
     378    END DO 
     379  END SUBROUTINE divide_by_mass 
     380     
    347381END MODULE observable_mod 

  • codes/icosagcm/trunk/src/timeloop_gcm.f90

    r407 r413  
    3131 
    3232    CHARACTER(len=255) :: def 
    33      
    34 !    IF (xios_output) itau_out=1 
     33 
     34    CALL init_caldyn 
     35     
     36    IF (xios_output) itau_out=1 
    3537    IF (.NOT. enable_io) itau_out=HUGE(itau_out) 
    3638 
    … …  
    4143       STOP 
    4244    END IF 
    43  
    44     nqdyn = 1 ! one dynamical tracer = theta for the moment 
    4545 
    4646    def='ARK2.3' 
    … …  
    141141    CALL init_sponge 
    142142    CALL init_observable 
    143     CALL init_caldyn 
    144143    CALL init_guided 
    145144    CALL init_advect_tracer 
    … …  
    220219    CALL trace_on 
    221220 
     221    IF (xios_output) THEN ! we must call update_calendar before any XIOS output 
     222       CALL xios_update_calendar(1) 
     223    END IF 
     224    CALL write_output_fields_basic(.TRUE., f_phis, f_ps, f_mass, f_geopot, f_theta_rhodz, f_u, f_W, f_q) 
     225 
    222226    DO it=itau0+1,itau0+itaumax 
    223227 
    224228       IF (is_master) CALL print_iteration(it, itau0, itaumax, start_clock, rate_clock) 
     229 
    225230       IF (xios_output) THEN 
    226           CALL xios_update_calendar(it) 
     231          IF(it>itau0+1) CALL xios_update_calendar(it-itau0) 
    227232       ELSE 
    228233          CALL update_time_counter(dt*it) 
    … …  
    318323               f_ps,f_dps,f_u,f_theta_rhodz,f_phis) 
    319324          CALL check_conserve(f_ps,f_dps,f_u,f_theta_rhodz,f_phis,it)  
    320        ENDIF        
     325       ENDIF 
    321326 
    322327       IF (mod(it,itau_out)==0 ) THEN 
    323328          CALL transfert_request(f_u,req_e1_vect) 
    324           CALL write_output_fields_basic(f_phis, f_ps, f_mass, f_geopot, f_theta_rhodz, f_u, f_W, f_q) 
     329          CALL write_output_fields_basic(.FALSE.,f_phis, f_ps, f_mass, f_geopot, f_theta_rhodz, f_u, f_W, f_q) 
    325330       ENDIF 
    326331 
Note: See TracChangeset for help on using the changeset viewer.