! module co2_module use params_module,only: nlevp1,nlonp4,nlat implicit none ! ! real,dimension(nlonp4,nlat) :: co2_ubc ! upper boundary ! real,dimension(nlonp4,3,nlat) :: co2_lbc ! lower boundary ! real,dimension(nlevp1,nlonp4,nlat) :: ! | co2_prod, ! production of co2 ! | co2_loss ! loss of co2 ! ! Boundary conditions and production and loss terms are allocated ! subdomains by sub alloc_co2 (called from allocdata.F). ! real,allocatable,dimension(:,:) :: co2_ubc ! upper boundary (i,j) real,allocatable,dimension(:,:,:) :: co2_lbc ! lower boundary (i,3,j) real,allocatable,dimension(:,:,:) :: | co2_prod, ! production of co2 (k,i,j) | co2_loss ! loss of co2 (k,i,j) real :: phi_co2(3) = (/1.199, 3.91, 1.201/) real,parameter :: alfa_co2 = 0. ! thermal diffusion coefficient ! contains !----------------------------------------------------------------------- subroutine comp_co2(tn,o2,o1,n4s,n2d,ne,o2p,op,o1d,co,oh, | barm,xnmbar,lev0,lev1,lon0,lon1,lat) ! ! Advance CO2 by one time step: ! use cons_module,only: rmass_co2,p0,boltz,expz, | expzmid_inv,rmassinv_o2,rmassinv_co,rmassinv_oh, | rmassinv_o1,rmassinv_n2 use chemrates_module,only: beta1,beta2,beta3,beta6,beta8,beta9, | beta9n,beta17,rk5,rk13,gam14,rkm42,gam15,co2mix use fields_module,only: tlbc ! lower boundary interface level of TN use qrj_module,only: | pdco2t ! total photodissociation of co2 (F(NJCO2T)) ! ! Args: integer,intent(in) :: lev0,lev1,lon0,lon1,lat real,dimension(lev0:lev1,lon0-2:lon1+2),intent(in) :: | tn, ! neutral temperature (deg K) | o2, o1, ! o2, o mass mixing ratios | n4s, n2d, ! n4s, n2d mass mixing ratios | ne, o2p, ! electron density, o2+ | op,o1d,co,oh, ! additional inputs | barm, ! mean molecular weight | xnmbar ! p0*e(-z)*barm/kT ! ! Local: integer :: k,i integer ::i0,i1,nk,nkm1 ! for addfsech real,dimension(lev0:lev1,lon0:lon1) :: | xn2 ! N2 (mmr) ! i0 = lon0 i1 = lon1 nk = lev1-lev0+1 nkm1 = nk-1 ! ! Diffusive equilibrium at upper boundary: do i=lon0,lon1 co2_ubc(i,lat) = 0. ! ! Lower boundary: ! (note tn(lev1,i) assumes bottom boundary of tn is in top slot): ! 10/3/05 btf: bottom boundary of TN is in tlbc: ! co2_lbc(i,1,lat) = 0. co2_lbc(i,2,lat) = 1. co2_lbc(i,3,lat) = -co2mix*rmass_co2*boltz*tlbc(i,lat)/ | (p0*expzmid_inv*expz(lev0)*barm(lev0,i)) enddo ! i=lon0,lon1 ! ! Sources: do i=lon0,lon1 do k=lev0,lev1-1 xn2(k,i) = 1.-o2(k,i)-o1(k,i) ! co2_prod(k,i,lat) = xnmbar(k,i)**2*(rkm42(k,i,lat)*co(k,i)* | rmassinv_co*oh(k,i)*rmassinv_oh+gam14(k,i,lat)*co(k,i)* | rmassinv_co*o1(k,i)*rmassinv_o1*xnmbar(k,i)/ | (0.5*(barm(k,i)+barm(k+1,i)))) ! ! co2_loss(k,i,lat) = -(rk13*op(k,i)+pdco2t(k,i,lat)+gam15* | o1d(k,i)*rmassinv_o1*xnmbar(k,i)) enddo ! k=lev0,lev1-1 enddo ! i=lon0,lon1 ! call addfsech('XNMBAR' ,' ',' ',xnmbar ,i0,i1,nk,nkm1,lat) ! call addfsech('CO2_PROD',' ',' ',co2_prod(:,i0:i1,lat), ! | i0,i1,nk,nkm1,lat) ! call addfsech('CO2_LOSS',' ',' ',co2_loss(:,i0:i1,lat), ! | i0,i1,nk,nkm1,lat) ! end subroutine comp_co2 !----------------------------------------------------------------------- subroutine minor_co2(tn,o2,o1,co2,co2_nm1,co2_out,co2_nm1_out, | lev0,lev1,lon0,lon1,lat0,lat1) use cons_module,only: rmass_co2 ! ! Input args at full task subdomains: integer,intent(in) :: lev0,lev1,lon0,lon1,lat0,lat1 real,dimension(lev0:lev1,lon0-2:lon1+2,lat0-2:lat1+2),intent(in):: | tn, ! neutral temperature (deg K) | o2, ! molecular oxygen (mmr) | o1, ! atomic oxygen (mmr) | co2, ! CO2 (mmr) | co2_nm1 ! CO2 at time n-1 ! ! Output args also at full task subdomains: real,dimension(lev0:lev1,lon0-2:lon1+2,lat0-2:lat1+2), | intent(out) :: | co2_out, ! CO2 output | co2_nm1_out ! CO2 output at time n-1 ! ! Local: integer :: lat integer ::nk,nkm1 ! for addfsech ! nk = lev1-lev0+1 nkm1 = nk-1 ! ! Minor returns co2_out and co2_nm1_out. Module data co2_prod, ! co2_loss, etc, were defined by comp_co2. ! ! subroutine minor(tn,o2,o1,fcomp,fcomp_tm1,fcomp_out,fcomp_tm1_out, ! | sloss,sprod,flbc,fubc,rmx,phix,alfax,lev0,lev1,lon0,lon1,lat0, ! | lat1,idebug) call minor(tn,o2,o1,co2,co2_nm1,co2_out,co2_nm1_out, | co2_loss,co2_prod,co2_lbc,co2_ubc,rmass_co2,phi_co2,alfa_co2, | lev0,lev1,lon0,lon1,lat0,lat1,0) ! do lat=lat0,lat1 ! call addfsech('CO2_OUT' ,' ',' ',co2_out(:,lon0:lon1,lat), ! | lon0,lon1,nk,nkm1,lat) ! call addfsech('CO2_TM1' ,' ',' ',co2_nm1_out(:,lon0:lon1,lat), ! | lon0,lon1,nk,nkm1,lat) ! enddo ! lat=lat0,lat1 end subroutine minor_co2 !----------------------------------------------------------------------- subroutine alloc_co2(lon0,lon1,lat0,lat1) ! ! Allocate subdomains (without ghost cells) to module data for boundary ! conditions and production and loss terms. This is called once per run ! from sub allocdata (allocdata.F). ! ! Args: integer,intent(in) :: lon0,lon1,lat0,lat1 ! ! Local: integer :: istat ! ! Allocate (i,j) subdomains to boundary conditions: allocate(co2_ubc(lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_no: error allocating', | ' co2_ubc: stat=',i3)") istat allocate(co2_lbc(lon0:lon1,3,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_no: error allocating', | ' co2_lbc: stat=',i3)") istat ! ! Allocate (k,i,j) subdomains to production and loss: allocate(co2_prod(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_no: error allocating', | ' co2_prod: stat=',i3)") istat allocate(co2_loss(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> alloc_no: error allocating', | ' co2_loss: stat=',i3)") istat ! write(6,"('alloc_co2: allocated module data')") ! end subroutine alloc_co2 !----------------------------------------------------------------------- end module co2_module