<p><b>laura@ucar.edu</b> 2010-11-08 12:06:52 -0700 (Mon, 08 Nov 2010)</p><p>minor modifications to allow sourcecode to compile without warnings using g95<br>
</p><hr noshade><pre><font color="gray">Modified: branches/atmos_physics/src/core_physics/module_driver_convection_deep.F
===================================================================
--- branches/atmos_physics/src/core_physics/module_driver_convection_deep.F        2010-11-08 19:06:33 UTC (rev 601)
+++ branches/atmos_physics/src/core_physics/module_driver_convection_deep.F        2010-11-08 19:06:52 UTC (rev 602)
@@ -71,8 +71,8 @@
do j = jms,jme
do i = ims,ime
+ cu_act_flag(i,j) = .false.
pratec_p(i,j) = 0.
- cu_act_flag(i,j) = 0.
cubot_p(i,j) = 0.
cutop_p(i,j) = 0.
nca_p(i,j) = 0.
@@ -129,15 +129,15 @@
end subroutine convection_deep_deallocate
!=============================================================================================
- subroutine convection_deep_init(s)
+ subroutine convection_deep_init(state)
!=============================================================================================
- type (state_type), intent(in) :: s
+ type (state_type), intent(in) :: state
!local variables and arrays:
!---------------------------
logical:: allowed_to_read
- integer:: p_qi,p_qs,p_first_scalar
+ integer:: p_qi,p_qs,p_first
!---------------------------------------------------------------------------------------------
@@ -148,9 +148,9 @@
write(0,*) '--- begin kain-fritsch initialization:'
allowed_to_read = .false.
- p_first_scalar = s % moist_start + 1
- p_qi = s % index_qi
- p_qs = s % index_qs
+ p_first = state % moist_start + 1
+ p_qi = state % index_qi
+ p_qs = state % index_qs
f_qv = .false.
f_qc = .false.
@@ -158,15 +158,17 @@
f_qi = .true.
f_qs = .true.
- call kf_eta_init(rthcuten_p,rqvcuten_p, &
- rqccuten_p,rqrcuten_p, &
- rqicuten_p,rqscuten_p, &
- nca_p,w0avg_p,p_qi,p_qs, &
- svp1,svp2,svp3,svpt0, &
- p_first_scalar,restart,allowed_to_read, &
- ids,ide,jds,jde,kds,kde, &
- ims,ime,jms,jme,kms,kme, &
- its,ite,jts,jte,kts,kte)
+ call kf_eta_init( &
+ rthcuten = rthcuten_p , rqvcuten = rqvcuten_p , rqccuten = rqccuten_p , &
+ rqrcuten = rqrcuten_p , rqicuten = rqicuten_p , rqscuten = rqscuten_p , &
+ nca = nca_p , w0avg = w0avg_p , p_qi = p_qi , &
+ p_qs = p_qs , svp1 = svp1 , svp2 = svp2 , &
+ svp3 = svp3 , svpt0 = svpt0 , p_first_scalar = p_first , &
+ restart = restart , allowed_to_read = allowed_to_read , &
+ ids = ids , ide = ide , jds = jds , jde = jde , kds = kds , kde = kde , &
+ ims = ims , ime = ime , jms = jms , jme = jme , kms = kds , kme = kme , &
+ its = its , ite = ite , jts = jts , jte = jte , kts = kts , kte = kte &
+ )
write(0,*) '--- end kain-kritsch initialization:'
case default
@@ -177,14 +179,15 @@
!=============================================================================================
- subroutine convection_deep_driver(itimestep,mesh,state)
+ subroutine convection_deep_driver(itimestep,mesh,diag_physics,tend_physics)
!=============================================================================================
!input and output arguments:
!---------------------------
integer,intent(in):: itimestep
type(mesh_type),intent(in):: mesh
- type(state_type),intent(inout):: state
+ type(diag_physics_type),intent(inout):: diag_physics
+ type(tend_physics_type),intent(inout):: tend_physics
!local variables and arrays:
!---------------------------
@@ -203,7 +206,7 @@
!initialize instantaneous precipitation, and copy convective tendencies from the dynamics to
!the physics grid:
- call convection_from_MPAS(state)
+ call convection_from_MPAS(diag_physics,tend_physics)
!call convection scheme:
@@ -226,31 +229,6 @@
enddo
enddo
- !call to kain-fritsch-eta convection scheme:
- !write(0,*) '--- enter kf_eta_cps:'
- !call kf_eta_cps( &
- ! !wrf-like dimensions:
- ! ids,ide,jds,jde,kds,kde, &
- ! ims,ime,jms,jme,kms,kme, &
- ! its,itf,jts,jtf,kts,ktf, &
- ! dt_dyn,itimestep,dx,cudt_pass, &
- ! curr_secs_pass, &
- ! adapt_step_flag_pass, &
- ! rho_p,raincv_p,pratec_p, &
- ! nca_p, &
- ! u_p,v_p,th_p,t_p, &
- ! w_p,dz_p,pres_p,pi_p, &
- ! w0avg_p,xlv0,xlv1,xls0,xls1, &
- ! cp,r_d,g,ep_1,ep_2, &
- ! svp1,svp2,svp3,svpt0, &
- ! n_cu,cu_act_flag,warm_rain, &
- ! cutop_p,cubot_p,qv_p, &
- ! f_qv,f_qc,f_qr,f_qi,f_qs, &
- ! rthcuten_p,rqvcuten_p, &
- ! rqccuten_p,rqrcuten_p, &
- ! rqicuten_p,rqscuten_p &
- ! )
-
call kf_eta_cps ( &
dt = dt_dyn , ktau = itimestep , &
dx = dx , cudt = dt_cu , &
@@ -289,37 +267,34 @@
!copy instantaneous and accumulated precipitation, convective tendencies, and "other" arrays
!specific to convection parameterization back to the dynamics grid:
- call convection_to_MPAS(state)
+ call convection_to_MPAS(diag_physics,tend_physics)
write(0,*) '--- end subroutine convection_driver'
end subroutine convection_deep_driver
!=============================================================================================
- subroutine convection_from_MPAS(state)
+ subroutine convection_from_MPAS(diag_physics,tend_physics)
!=============================================================================================
!input arguments:
- type(state_type),intent(in):: state
+ type(diag_physics_type),intent(in):: diag_physics
+ type(tend_physics_type),intent(in):: tend_physics
!---------------------------------------------------------------------------------------------
do j = jts,jte
do i = its,ite
-! if(s % rainc % array(i) .GT. 0.) &
-! write(0,201) i,s%raincv%array(i),s%rainc%array(i), &
-! s%cubot%array(i),s%cutop%array(i)
-
raincv_p(i,j) = 0.
- rainc_p(i,j) = state % rainc % array(i)
+ rainc_p(i,j) = diag_physics % rainc % array(i)
do k = kts, ktf
- rthcuten_p(i,k,j) = state % rthcuten % array(k,i)
- rqvcuten_p(i,k,j) = state % rqvcuten % array(k,i)
- rqccuten_p(i,k,j) = state % rqccuten % array(k,i)
- rqrcuten_p(i,k,j) = state % rqrcuten % array(k,i)
- rqicuten_p(i,k,j) = state % rqicuten % array(k,i)
- rqscuten_p(i,k,j) = state % rqscuten % array(k,i)
+ rthcuten_p(i,k,j) = tend_physics % rthcuten % array(k,i)
+ rqvcuten_p(i,k,j) = tend_physics % rqvcuten % array(k,i)
+ rqccuten_p(i,k,j) = tend_physics % rqccuten % array(k,i)
+ rqrcuten_p(i,k,j) = tend_physics % rqrcuten % array(k,i)
+ rqicuten_p(i,k,j) = tend_physics % rqicuten % array(k,i)
+ rqscuten_p(i,k,j) = tend_physics % rqscuten % array(k,i)
enddo
enddo
@@ -332,12 +307,12 @@
do j = jts,jte
do i = its,ite
- nca_p(i,j) = state % nca % array(i)
- cubot_p(i,j) = state % cubot % array(i)
- cutop_p(i,j) = state % cutop % array(i)
+ nca_p(i,j) = diag_physics % nca % array(i)
+ cubot_p(i,j) = diag_physics % cubot % array(i)
+ cutop_p(i,j) = diag_physics % cutop % array(i)
do k = kts, ktf
- w0avg_p(i,k,j) = state % w0avg % array(k,i)
+ w0avg_p(i,k,j) = diag_physics % w0avg % array(k,i)
enddo
enddo
@@ -350,31 +325,36 @@
end subroutine convection_from_MPAS
!=============================================================================================
- subroutine convection_to_MPAS(state)
+ subroutine convection_to_MPAS(diag_physics,tend_physics)
!=============================================================================================
!inout arguments:
- type(state_type),intent(inout):: state
+ type(diag_physics_type),intent(inout):: diag_physics
+ type(tend_physics_type),intent(inout):: tend_physics
!---------------------------------------------------------------------------------------------
-!write(0,*) '--- enter convection_deep_end:'
+ write(0,*) '--- enter convection_deep_end:'
do j = jts,jte
do i = its,ite
- state % raincv % array(i) = raincv_p(i,j)
- state % rainc % array(i) = state % rainc % array(i) + state % raincv % array(i)
+ diag_physics % raincv % array(i) = raincv_p(i,j)
do k = kts, ktf
- state % rthcuten % array(k,i) = rthcuten_p(i,k,j)
- state % rqvcuten % array(k,i) = rqvcuten_p(i,k,j)
- state % rqccuten % array(k,i) = rqccuten_p(i,k,j)
- state % rqrcuten % array(k,i) = rqrcuten_p(i,k,j)
- state % rqicuten % array(k,i) = rqicuten_p(i,k,j)
- state % rqscuten % array(k,i) = rqscuten_p(i,k,j)
+ tend_physics % rthcuten % array(k,i) = rthcuten_p(i,k,j)
+ tend_physics % rqvcuten % array(k,i) = rqvcuten_p(i,k,j)
+ tend_physics % rqccuten % array(k,i) = rqccuten_p(i,k,j)
+ tend_physics % rqrcuten % array(k,i) = rqrcuten_p(i,k,j)
+ tend_physics % rqicuten % array(k,i) = rqicuten_p(i,k,j)
+ tend_physics % rqscuten % array(k,i) = rqscuten_p(i,k,j)
enddo
enddo
enddo
+
+ do i = its,ite
+ diag_physics % rainc % array(i) = diag_physics % rainc % array(i) &
+ + diag_physics % raincv % array(i)
+ enddo
convection_select: select case(conv_deep_scheme)
@@ -383,12 +363,12 @@
do j = jts,jte
do i = its,ite
- state % nca % array(i) = nca_p(i,j)
- state % cubot % array(i) = cubot_p(i,j)
- state % cutop % array(i) = cutop_p(i,j)
+ diag_physics % nca % array(i) = nca_p(i,j)
+ diag_physics % cubot % array(i) = cubot_p(i,j)
+ diag_physics % cutop % array(i) = cutop_p(i,j)
do k = kts, ktf
- state % w0avg % array(k,i) = w0avg_p(i,k,j)
+ diag_physics % w0avg % array(k,i) = w0avg_p(i,k,j)
enddo
enddo
</font>
</pre>