#include module heatnirco2_module use addfld_module,only: addfld implicit none ! integer,parameter :: ndpara=62,ncolgr=10 real,dimension(ndpara) :: xspara,co2stand,zppara real,dimension(ndpara,ncolgr) :: colmpara,corrnormpara ! ! Heatnirco2 and heatkday outputs are allocated by sub heatnirco2_init ! (see below) at subdomain dimensions: (nlevp1,lon0:lon1,lat0:lat1) ! real,allocatable,dimension(:,:,:) :: | heatnirco2, ! ergs | heatkday ! deg/day ! contains !----------------------------------------------------------------------- subroutine heatnearirco2(o2,o1,n2,co2,barm,vco2,scco2, | lev0,lev1,lon0,lon1,lat) use params_module,only: dz use cons_module,only: rmassinv_o2,rmassinv_o1,rmassinv_n2, | rmassinv_co2 use chapman_module,only: chi use qrj_module,only: qtotal ! (nlevp1,lon0:lon1,lat0:lat1) ! ! Args: integer,intent(in) :: lev0,lev1,lon0,lon1,lat real,dimension(lev0:lev1,lon0-2:lon1+2),intent(in) :: | o2,o1,n2,co2,barm,vco2,scco2 ! ! Local: integer :: i,k,kk,id,icolm,izp,kdpara real,dimension(ndpara,lon0:lon1) :: cvmu,colzp,scco2_dzp, | vmro2_dzp,vmro_dzp,vmrn2_dzp,vmrco2_dzp,heatnirco2_dzp real,parameter :: smallvalue=1.0E-20 real,parameter :: coeftr=8.31441E7/86400.0 real,dimension(lev0:lev1,lon0:lon1) :: vmro2,vmro,vmrn2,vmrco2 real :: reldcolm,avmol logical :: dblzp real :: chi_ki(lev0:lev1,lon0:lon1) ! for addfld ! dblzp = .false. ! normal vertical resolution if (dz==0.25) dblzp = .true. ! double vertical resolution ! ! Convert densities from mmr (f-array) to volume mixing ratio: do i=lon0,lon1 do k=lev0,lev1 vmro2(k,i) = o2(k,i)*barm(k,i)*rmassinv_o2 vmro (k,i) = o1(k,i)*barm(k,i)*rmassinv_o1 vmrn2(k,i) = n2(k,i)*barm(k,i)*rmassinv_n2 ! ! 12/03/04 btf bugfix: change rmassinv_n2 to rmassinv_co2 vmrco2(k,i) = co2(k,i)*barm(k,i)*rmassinv_co2 enddo enddo do i=lon0,lon1 chi_ki(:,i) = chi(i,lat) enddo ! ! call addfld('chi',' ',' ',chi_ki, ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('scco2' ,' ',' ',scco2 (:,lon0:lon1), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('vmro2' ,' ',' ',vmro2 (lev0:lev1,:) , ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('vmro' ,' ',' ',vmro (lev0:lev1,:) , ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('vmrn2' ,' ',' ',vmrn2 (lev0:lev1,:) , ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('vmrco2',' ',' ',vmrco2(lev0:lev1,:), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! ! Define densities at double resolution up to zp=-1.75 (k==ndpara). ! If normal resolution (2 grid points per scale height), interpolate ! to double resolution: ! if (dblzp) then ! double res kdpara = ndpara vmro2_dzp (:,:) = vmro2 (1:kdpara,:) vmro_dzp (:,:) = vmro (1:kdpara,:) vmrn2_dzp (:,:) = vmrn2 (1:kdpara,:) vmrco2_dzp(:,:) = vmrco2(1:kdpara,:) scco2_dzp (:,:) = scco2 (1:kdpara,lon0:lon1) else ! normal res (interpolate to double) do k=lev0,lev1 kk = k*2-1 if (kk <= ndpara) then vmro2_dzp (kk,:) = vmro2 (k,:) vmro_dzp (kk,:) = vmro (k,:) vmrn2_dzp (kk,:) = vmrn2 (k,:) vmrco2_dzp(kk,:) = vmrco2(k,:) scco2_dzp (kk,:) = scco2 (k,lon0:lon1) endif enddo ! k=lev0,lev1 do k=2,ndpara-1,2 vmro2_dzp (k,:) = (vmro2_dzp (k-1,:)+vmro2_dzp (k+1,:))*0.5 vmro_dzp (k,:) = (vmro_dzp (k-1,:)+vmro_dzp (k+1,:))*0.5 vmrn2_dzp (k,:) = (vmrn2_dzp (k-1,:)+vmrn2_dzp (k+1,:))*0.5 vmrco2_dzp(k,:) = (vmrco2_dzp(k-1,:)+vmrco2_dzp(k+1,:))*0.5 scco2_dzp (k,:) = (scco2_dzp (k-1,:)+scco2_dzp (k+1,:))*0.5 enddo ! ! 5/18/04 btf: define at ndpara: vmro2_dzp (ndpara,:) = vmro2_dzp (ndpara-1,:) vmro_dzp (ndpara,:) = vmro_dzp (ndpara-1,:) vmrn2_dzp (ndpara,:) = vmrn2_dzp (ndpara-1,:) vmrco2_dzp(ndpara,:) = vmrco2_dzp(ndpara-1,:) scco2_dzp (ndpara,:) = scco2_dzp (ndpara-1,:) endif ! ! real,dimension(ndpara,lon0:lon1) :: cvmu,colzp,scco2_dzp, !| vmro2_dzp,vmro_dzp,vmrn2_dzp,vmrco2_dzp,heatnirco2_dzp ! heatnirco2_dzp(:,:) = smallvalue colzp(:,:) = 0. cvmu(:,:) = 0. do i=lon0,lon1 do k=1,ndpara if (chi(i,lat) < 1.57) then ! write(6,"('heatnirco2: lat=',i3,' i=',i3,' k=',i3, ! | ' vmrn2_dzp=',e12.4,' vmro2_dzp=',e12.4,' vmro_dzp=', ! | e12.4)") lat,i,k,vmrn2_dzp(k,i),vmro2_dzp(k,i), ! | vmro_dzp(k,i) avmol=28.0*vmrn2_dzp(k,i)+32.0*vmro2_dzp(k,i)+ | 16.0*vmro_dzp(k,i) cvmu(k,i)=coeftr/avmol*(7.0/2.0*(1.0-vmro_dzp(k,i))+5.0/2.0* | vmro_dzp(k,i)) colzp(k,i)=scco2_dzp(k,i) endif ! chi ! ! colzp is 0 at j==36, i==57 (lon 90, 6 am) ! tgcm24 is 0 here only if init to 0 is added before the above loop. ! write(6,"('heatnirco2: k=',i3,' i=',i3,' j=',i3, ! | ' chi=',f15.7,' scco2_dzp=',e12.4,' colzp=',e12.4)") ! | k,i,lat,chi(i,lat),scco2_dzp(k,i),colzp(k,i) enddo ! k=1,ndpara enddo ! i=lon0,lon1 ! call addfld('colzp' ,' ',' ',colzp(lev0:lev1,:) , ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! do i=lon0,lon1 do izp=1,ndpara if(colzp(izp,i) > 0.0.and.vmrco2_dzp(izp,i) > 0.0) | colzp(izp,i)=alog(colzp(izp,i)) enddo enddo ! real,dimension(ndpara,lon0:lon1) :: cvmu,colzp,scco2_dzp, ! call addfld('cvmu' ,' ',' ',cvmu(lev0:lev1,:) , ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('colzp' ,' ',' ',colzp(lev0:lev1,:) , ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! do i=lon0,lon1 do izp=1,ndpara if (chi(i,lat) < 1.57) then ! ! Linear interpolation over column density for given altitude point ! if (colzp(izp,i) < colmpara(izp,1)) then heatnirco2_dzp(izp,i) = corrnormpara(izp,1) elseif (colzp(izp,i) >= colmpara(izp,ncolgr)) then heatnirco2_dzp(izp,i) = corrnormpara(izp,ncolgr) else loop1: do icolm=1,ncolgr-1 if (colzp(izp,i) >= colmpara(izp,icolm).and. | colzp(izp,i) < colmpara(izp,icolm+1)) then reldcolm=(colzp(izp,i)-colmpara(izp,icolm))/ | (colmpara(izp,icolm+1)-colmpara(izp,icolm)) heatnirco2_dzp(izp,i)=corrnormpara(izp,icolm)+ | (corrnormpara(izp,icolm+1)-corrnormpara(izp,icolm))* | reldcolm exit loop1 endif enddo loop1 endif ! ! From normalized value to the one corresponding to the given vmrco2 heatnirco2_dzp(izp,i)=heatnirco2_dzp(izp,i)* | vmrco2_dzp(izp,i)/co2stand(izp) endif ! chi enddo ! i enddo ! izp ! call addfld('qco2dzp' ,' ',' ',heatnirco2_dzp(lev0:lev1,:), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! ! Prepare outputs at model's vertical resolution: heatnirco2(:,:,:) = smallvalue heatkday(:,:,:) = smallvalue if (dblzp) then ! double res kdpara = ndpara heatkday(1:kdpara,:,lat) = heatnirco2_dzp(:,:) heatnirco2(1:kdpara,:,lat) = heatnirco2_dzp(:,:)*cvmu(:,:) ! k/day->ergs else ! normal res do k=lev0,lev1 kk = k*2-1 if (kk <= ndpara) then heatkday(k,:,lat) = heatnirco2_dzp(kk,:) heatnirco2(k,:,lat) = heatnirco2_dzp(kk,:)*cvmu(kk,:) ! k/day->ergs endif enddo ! k=lev0,lev1 endif ! ! Add co2 heating to total: do i=lon0,lon1 do k=lev0,lev1 qtotal(k,i,lat) = qtotal(k,i,lat)+heatnirco2(k,i,lat) enddo ! k=lev0,lev1 enddo ! i=lon0,lon1 ! ! call addfld('QTOT_CO2',' ',' ',qtotal(:,lon0:lon1,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('qco2ergs',' ',' ',heatnirco2(:,lon0:lon1,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! call addfld('qco2kday',' ',' ',heatkday(:,lon0:lon1,lat), ! | 'lev',lev0,lev1,'lon',lon0,lon1,lat) ! end subroutine heatnearirco2 !----------------------------------------------------------------------- subroutine heatnirco2_init use params_module,only: nlevp1 use mpi_module,only: lon0,lon1,lat0,lat1 ! ! Called once per run from init (init.F) to define module data. ! ! Allocate co2 heating: ! real,allocatable,dimension(:,:,:) :: heatnirco2,heatkday ! integer :: istat allocate(heatnirco2(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> heatnirco2_init: error allocating', | ' heatnirco2: stat=',i3)") istat allocate(heatkday(nlevp1,lon0:lon1,lat0:lat1),stat=istat) if (istat /= 0) write(6,"('>>> heatnirco2_init: error allocating', | ' heatkday: stat=',i3)") istat ! zppara = (/ | -17.00, -16.75, -16.50, -16.25, -16.00, | -15.75, -15.50, -15.25, -15.00, -14.75, | -14.50, -14.25, -14.00, -13.75, -13.50, | -13.25, -13.00, -12.75, -12.50, -12.25, | -12.00, -11.75, -11.50, -11.25, -11.00, | -10.75, -10.50, -10.25, -10.00, -9.75, | -9.50, -9.25, -9.00, -8.75, -8.50, | -8.25, -8.00, -7.75, -7.50, -7.25, | -7.00, -6.75, -6.50, -6.25, -6.00, | -5.75, -5.50, -5.25, -5.00, -4.75, | -4.50, -4.25, -4.00, -3.75, -3.50, | -3.25, -3.00, -2.75, -2.50, -2.25, | -2.00, -1.75/) xspara = (/ | 4.4164E+00, 4.6664E+00, 4.9164E+00, 5.1664E+00, 5.4164E+00, | 5.6664E+00, 5.9164E+00, 6.1664E+00, 6.4164E+00, 6.6664E+00, | 6.9164E+00, 7.1664E+00, 7.4164E+00, 7.6664E+00, 7.9164E+00, | 8.1664E+00, 8.4164E+00, 8.6664E+00, 8.9164E+00, 9.1664E+00, | 9.4164E+00, 9.6664E+00, 9.9164E+00, 1.0166E+01, 1.0416E+01, | 1.0666E+01, 1.0916E+01, 1.1166E+01, 1.1416E+01, 1.1666E+01, | 1.1916E+01, 1.2166E+01, 1.2416E+01, 1.2666E+01, 1.2916E+01, | 1.3166E+01, 1.3416E+01, 1.3666E+01, 1.3916E+01, 1.4166E+01, | 1.4416E+01, 1.4666E+01, 1.4916E+01, 1.5166E+01, 1.5416E+01, | 1.5666E+01, 1.5916E+01, 1.6166E+01, 1.6416E+01, 1.6666E+01, | 1.6916E+01, 1.7166E+01, 1.7416E+01, 1.7666E+01, 1.7916E+01, | 1.8166E+01, 1.8416E+01, 1.8666E+01, 1.8916E+01, 1.9166E+01, | 1.9416E+01, 1.9666E+01/) co2stand = (/ | 3.6000E-04, 3.6000E-04, 3.6000E-04, 3.6000E-04, 3.6000E-04, | 3.6000E-04, 3.6000E-04, 3.6000E-04, 3.6000E-04, 3.6000E-04, | 3.6000E-04, 3.6000E-04, 3.6000E-04, 3.6000E-04, 3.6000E-04, | 3.6000E-04, 3.6000E-04, 3.6000E-04, 3.6000E-04, 3.6000E-04, | 3.6000E-04, 3.6000E-04, 3.6000E-04, 3.6000E-04, 3.6000E-04, | 3.6000E-04, 3.6000E-04, 3.6000E-04, 3.6000E-04, 3.6000E-04, | 3.6000E-04, 3.6000E-04, 3.5867E-04, 3.5534E-04, 3.5134E-04, | 3.4401E-04, 3.3235E-04, 3.1668E-04, 2.9902E-04, 2.8102E-04, | 2.6235E-04, 2.4335E-04, 2.2435E-04, 2.0535E-04, 1.8635E-04, | 1.6735E-04, 1.4835E-04, 1.3002E-04, 1.1202E-04, 9.5350E-05, | 8.2013E-05, 7.1344E-05, 6.2010E-05, 5.3675E-05, 4.6341E-05, | 3.9341E-05, 3.2341E-05, 2.6006E-05, 2.0672E-05, 1.5672E-05, | 1.1338E-05, 8.0032E-06/) ! colmpara(1,:) = (/ | 4.5112E+01, 4.5647E+01, 4.6183E+01, 4.6719E+01, 4.7254E+01, | 4.7789E+01, 4.8325E+01, 4.8861E+01, 4.9396E+01, 4.9931E+01/) colmpara(2,:) = (/ | 4.4862E+01, 4.5398E+01, 4.5934E+01, 4.6469E+01, 4.7004E+01, | 4.7540E+01, 4.8076E+01, 4.8611E+01, 4.9147E+01, 4.9682E+01/) colmpara(3,:) = (/ | 4.4613E+01, 4.5148E+01, 4.5684E+01, 4.6220E+01, 4.6755E+01, | 4.7290E+01, 4.7826E+01, 4.8362E+01, 4.8897E+01, 4.9432E+01/) colmpara(4,:) = (/ | 4.4363E+01, 4.4899E+01, 4.5435E+01, 4.5970E+01, 4.6506E+01, | 4.7041E+01, 4.7577E+01, 4.8112E+01, 4.8647E+01, 4.9183E+01/) colmpara(5,:) = (/ | 4.4114E+01, 4.4650E+01, 4.5185E+01, 4.5721E+01, 4.6256E+01, | 4.6791E+01, 4.7327E+01, 4.7863E+01, 4.8398E+01, 4.8933E+01/) colmpara(6,:) = (/ | 4.3865E+01, 4.4401E+01, 4.4936E+01, 4.5471E+01, 4.6007E+01, | 4.6542E+01, 4.7078E+01, 4.7613E+01, 4.8149E+01, 4.8684E+01/) colmpara(7,:) = (/ | 4.3615E+01, 4.4151E+01, 4.4686E+01, 4.5222E+01, 4.5757E+01, | 4.6293E+01, 4.6828E+01, 4.7364E+01, 4.7899E+01, 4.8435E+01/) colmpara(8,:) = (/ | 4.3366E+01, 4.3901E+01, 4.4437E+01, 4.4972E+01, 4.5508E+01, | 4.6043E+01, 4.6579E+01, 4.7114E+01, 4.7650E+01, 4.8185E+01/) colmpara(9,:) = (/ | 4.3116E+01, 4.3652E+01, 4.4188E+01, 4.4723E+01, 4.5258E+01, | 4.5794E+01, 4.6329E+01, 4.6865E+01, 4.7400E+01, 4.7936E+01/) colmpara(10,:) = (/ | 4.2867E+01, 4.3402E+01, 4.3938E+01, 4.4473E+01, 4.5009E+01, | 4.5544E+01, 4.6080E+01, 4.6616E+01, 4.7151E+01, 4.7686E+01/) colmpara(11,:) = (/ | 4.2617E+01, 4.3153E+01, 4.3688E+01, 4.4224E+01, 4.4759E+01, | 4.5295E+01, 4.5830E+01, 4.6366E+01, 4.6901E+01, 4.7437E+01/) colmpara(12,:) = (/ | 4.2368E+01, 4.2903E+01, 4.3439E+01, 4.3974E+01, 4.4510E+01, | 4.5045E+01, 4.5581E+01, 4.6116E+01, 4.6651E+01, 4.7187E+01/) colmpara(13,:) = (/ | 4.2118E+01, 4.2654E+01, 4.3189E+01, 4.3725E+01, 4.4260E+01, | 4.4796E+01, 4.5331E+01, 4.5867E+01, 4.6402E+01, 4.6938E+01/) colmpara(14,:) = (/ | 4.1869E+01, 4.2404E+01, 4.2940E+01, 4.3475E+01, 4.4011E+01, | 4.4546E+01, 4.5082E+01, 4.5617E+01, 4.6153E+01, 4.6688E+01/) colmpara(15,:) = (/ | 4.1619E+01, 4.2155E+01, 4.2690E+01, 4.3225E+01, 4.3761E+01, | 4.4296E+01, 4.4832E+01, 4.5367E+01, 4.5903E+01, 4.6438E+01/) colmpara(16,:) = (/ | 4.1369E+01, 4.1905E+01, 4.2440E+01, 4.2976E+01, 4.3511E+01, | 4.4046E+01, 4.4582E+01, 4.5118E+01, 4.5653E+01, 4.6188E+01/) colmpara(17,:) = (/ | 4.1119E+01, 4.1655E+01, 4.2190E+01, 4.2726E+01, 4.3261E+01, | 4.3797E+01, 4.4332E+01, 4.4868E+01, 4.5403E+01, 4.5938E+01/) colmpara(18,:) = (/ | 4.0869E+01, 4.1405E+01, 4.1940E+01, 4.2476E+01, 4.3011E+01, | 4.3547E+01, 4.4082E+01, 4.4618E+01, 4.5153E+01, 4.5688E+01/) colmpara(19,:) = (/ | 4.0619E+01, 4.1154E+01, 4.1690E+01, 4.2226E+01, 4.2761E+01, | 4.3296E+01, 4.3832E+01, 4.4367E+01, 4.4903E+01, 4.5438E+01/) colmpara(20,:) = (/ | 4.0368E+01, 4.0904E+01, 4.1439E+01, 4.1975E+01, 4.2510E+01, | 4.3046E+01, 4.3581E+01, 4.4117E+01, 4.4652E+01, 4.5188E+01/) colmpara(21,:) = (/ | 4.0117E+01, 4.0653E+01, 4.1189E+01, 4.1724E+01, 4.2259E+01, | 4.2795E+01, 4.3330E+01, 4.3866E+01, 4.4402E+01, 4.4937E+01/) colmpara(22,:) = (/ | 3.9866E+01, 4.0402E+01, 4.0937E+01, 4.1473E+01, 4.2008E+01, | 4.2544E+01, 4.3080E+01, 4.3615E+01, 4.4150E+01, 4.4686E+01/) colmpara(23,:) = (/ | 3.9615E+01, 4.0150E+01, 4.0686E+01, 4.1221E+01, 4.1757E+01, | 4.2292E+01, 4.2828E+01, 4.3363E+01, 4.3899E+01, 4.4434E+01/) colmpara(24,:) = (/ | 3.9363E+01, 3.9898E+01, 4.0434E+01, 4.0969E+01, 4.1505E+01, | 4.2040E+01, 4.2576E+01, 4.3111E+01, 4.3647E+01, 4.4182E+01/) colmpara(25,:) = (/ | 3.9110E+01, 3.9645E+01, 4.0181E+01, 4.0716E+01, 4.1252E+01, | 4.1787E+01, 4.2323E+01, 4.2858E+01, 4.3394E+01, 4.3929E+01/) colmpara(26,:) = (/ | 3.8856E+01, 3.9391E+01, 3.9927E+01, 4.0462E+01, 4.0998E+01, | 4.1533E+01, 4.2069E+01, 4.2604E+01, 4.3140E+01, 4.3675E+01/) colmpara(27,:) = (/ | 3.8601E+01, 3.9136E+01, 3.9672E+01, 4.0207E+01, 4.0743E+01, | 4.1278E+01, 4.1814E+01, 4.2349E+01, 4.2885E+01, 4.3420E+01/) colmpara(28,:) = (/ | 3.8344E+01, 3.8879E+01, 3.9415E+01, 3.9950E+01, 4.0486E+01, | 4.1021E+01, 4.1557E+01, 4.2092E+01, 4.2628E+01, 4.3163E+01/) colmpara(29,:) = (/ | 3.8085E+01, 3.8620E+01, 3.9156E+01, 3.9692E+01, 4.0227E+01, | 4.0762E+01, 4.1298E+01, 4.1833E+01, 4.2369E+01, 4.2904E+01/) colmpara(30,:) = (/ | 3.7823E+01, 3.8359E+01, 3.8894E+01, 3.9430E+01, 3.9965E+01, | 4.0501E+01, 4.1036E+01, 4.1572E+01, 4.2107E+01, 4.2642E+01/) colmpara(31,:) = (/ | 3.7558E+01, 3.8093E+01, 3.8628E+01, 3.9164E+01, 3.9699E+01, | 4.0235E+01, 4.0771E+01, 4.1306E+01, 4.1842E+01, 4.2377E+01/) colmpara(32,:) = (/ | 3.7287E+01, 3.7823E+01, 3.8358E+01, 3.8894E+01, 3.9429E+01, | 3.9965E+01, 4.0500E+01, 4.1036E+01, 4.1571E+01, 4.2107E+01/) colmpara(33,:) = (/ | 3.7011E+01, 3.7547E+01, 3.8082E+01, 3.8617E+01, 3.9153E+01, | 3.9689E+01, 4.0224E+01, 4.0759E+01, 4.1295E+01, 4.1830E+01/) colmpara(34,:) = (/ | 3.6728E+01, 3.7263E+01, 3.7799E+01, 3.8334E+01, 3.8870E+01, | 3.9405E+01, 3.9941E+01, 4.0476E+01, 4.1012E+01, 4.1547E+01/) colmpara(35,:) = (/ | 3.6437E+01, 3.6972E+01, 3.7508E+01, 3.8043E+01, 3.8578E+01, | 3.9114E+01, 3.9650E+01, 4.0185E+01, 4.0720E+01, 4.1256E+01/) colmpara(36,:) = (/ | 3.6136E+01, 3.6672E+01, 3.7207E+01, 3.7743E+01, 3.8278E+01, | 3.8814E+01, 3.9349E+01, 3.9885E+01, 4.0420E+01, 4.0956E+01/) colmpara(37,:) = (/ | 3.5827E+01, 3.6363E+01, 3.6898E+01, 3.7434E+01, 3.7969E+01, | 3.8505E+01, 3.9040E+01, 3.9575E+01, 4.0111E+01, 4.0647E+01/) colmpara(38,:) = (/ | 3.5511E+01, 3.6046E+01, 3.6582E+01, 3.7117E+01, 3.7653E+01, | 3.8188E+01, 3.8724E+01, 3.9259E+01, 3.9795E+01, 4.0331E+01/) colmpara(39,:) = (/ | 3.5189E+01, 3.5724E+01, 3.6259E+01, 3.6795E+01, 3.7331E+01, | 3.7866E+01, 3.8402E+01, 3.8937E+01, 3.9472E+01, 4.0008E+01/) colmpara(40,:) = (/ | 3.4860E+01, 3.5395E+01, 3.5931E+01, 3.6466E+01, 3.7002E+01, | 3.7538E+01, 3.8073E+01, 3.8608E+01, 3.9144E+01, 3.9679E+01/) colmpara(41,:) = (/ | 3.4525E+01, 3.5061E+01, 3.5596E+01, 3.6132E+01, 3.6667E+01, | 3.7203E+01, 3.7738E+01, 3.8274E+01, 3.8809E+01, 3.9345E+01/) colmpara(42,:) = (/ | 3.4184E+01, 3.4719E+01, 3.5254E+01, 3.5790E+01, 3.6326E+01, | 3.6861E+01, 3.7396E+01, 3.7932E+01, 3.8467E+01, 3.9003E+01/) colmpara(43,:) = (/ | 3.3835E+01, 3.4370E+01, 3.4906E+01, 3.5441E+01, 3.5977E+01, | 3.6512E+01, 3.7048E+01, 3.7583E+01, 3.8119E+01, 3.8654E+01/) colmpara(44,:) = (/ | 3.3478E+01, 3.4014E+01, 3.4549E+01, 3.5085E+01, 3.5620E+01, | 3.6155E+01, 3.6691E+01, 3.7227E+01, 3.7762E+01, 3.8297E+01/) colmpara(45,:) = (/ | 3.3112E+01, 3.3648E+01, 3.4183E+01, 3.4719E+01, 3.5254E+01, | 3.5790E+01, 3.6325E+01, 3.6861E+01, 3.7396E+01, 3.7932E+01/) colmpara(46,:) = (/ | 3.2738E+01, 3.3273E+01, 3.3808E+01, 3.4344E+01, 3.4880E+01, | 3.5415E+01, 3.5950E+01, 3.6486E+01, 3.7021E+01, 3.7557E+01/) colmpara(47,:) = (/ | 3.2354E+01, 3.2889E+01, 3.3425E+01, 3.3960E+01, 3.4496E+01, | 3.5031E+01, 3.5567E+01, 3.6102E+01, 3.6638E+01, 3.7173E+01/) colmpara(48,:) = (/ | 3.1963E+01, 3.2498E+01, 3.3034E+01, 3.3569E+01, 3.4105E+01, | 3.4640E+01, 3.5176E+01, 3.5711E+01, 3.6247E+01, 3.6782E+01/) colmpara(49,:) = (/ | 3.1567E+01, 3.2103E+01, 3.2638E+01, 3.3174E+01, 3.3709E+01, | 3.4244E+01, 3.4780E+01, 3.5316E+01, 3.5851E+01, 3.6386E+01/) colmpara(50,:) = (/ | 3.1171E+01, 3.1707E+01, 3.2242E+01, 3.2778E+01, 3.3314E+01, | 3.3849E+01, 3.4384E+01, 3.4920E+01, 3.5456E+01, 3.5991E+01/) colmpara(51,:) = (/ | 3.0777E+01, 3.1313E+01, 3.1848E+01, 3.2384E+01, 3.2919E+01, | 3.3455E+01, 3.3990E+01, 3.4526E+01, 3.5062E+01, 3.5597E+01/) colmpara(52,:) = (/ | 3.0379E+01, 3.0914E+01, 3.1450E+01, 3.1985E+01, 3.2521E+01, | 3.3056E+01, 3.3592E+01, 3.4127E+01, 3.4663E+01, 3.5198E+01/) colmpara(53,:) = (/ | 2.9970E+01, 3.0506E+01, 3.1041E+01, 3.1577E+01, 3.2112E+01, | 3.2648E+01, 3.3183E+01, 3.3719E+01, 3.4254E+01, 3.4790E+01/) colmpara(54,:) = (/ | 2.9547E+01, 3.0083E+01, 3.0618E+01, 3.1154E+01, 3.1689E+01, | 3.2225E+01, 3.2760E+01, 3.3296E+01, 3.3831E+01, 3.4367E+01/) colmpara(55,:) = (/ | 2.9104E+01, 2.9640E+01, 3.0175E+01, 3.0711E+01, 3.1246E+01, | 3.1782E+01, 3.2317E+01, 3.2853E+01, 3.3388E+01, 3.3924E+01/) colmpara(56,:) = (/ | 2.8636E+01, 2.9171E+01, 2.9707E+01, 3.0242E+01, 3.0778E+01, | 3.1313E+01, 3.1849E+01, 3.2384E+01, 3.2920E+01, 3.3455E+01/) colmpara(57,:) = (/ | 2.8140E+01, 2.8675E+01, 2.9211E+01, 2.9746E+01, 3.0282E+01, | 3.0817E+01, 3.1353E+01, 3.1888E+01, 3.2424E+01, 3.2960E+01/) colmpara(58,:) = (/ | 2.7615E+01, 2.8151E+01, 2.8686E+01, 2.9221E+01, 2.9757E+01, | 3.0292E+01, 3.0828E+01, 3.1364E+01, 3.1899E+01, 3.2435E+01/) colmpara(59,:) = (/ | 2.7048E+01, 2.7583E+01, 2.8119E+01, 2.8654E+01, 2.9190E+01, | 2.9725E+01, 3.0261E+01, 3.0796E+01, 3.1332E+01, 3.1867E+01/) colmpara(60,:) = (/ | 2.6417E+01, 2.6952E+01, 2.7488E+01, 2.8023E+01, 2.8558E+01, | 2.9094E+01, 2.9630E+01, 3.0165E+01, 3.0701E+01, 3.1236E+01/) colmpara(61,:) = (/ | 2.5690E+01, 2.6226E+01, 2.6761E+01, 2.7297E+01, 2.7832E+01, | 2.8368E+01, 2.8903E+01, 2.9439E+01, 2.9974E+01, 3.0510E+01/) colmpara(62,:) = (/ | 2.4753E+01, 2.5288E+01, 2.5824E+01, 2.6359E+01, 2.6895E+01, | 2.7430E+01, 2.7966E+01, 2.8501E+01, 2.9036E+01, 2.9572E+01/) ! corrnormpara(1,:) = (/ | 1.0127E+00, 7.6998E-01, 5.8224E-01, 4.2806E-01, 2.9084E-01, | 1.9402E-01, 1.3023E-01, 9.1149E-02, 6.4704E-02, 4.7688E-02/) corrnormpara(2,:) = (/ | 1.0187E+00, 7.7869E-01, 5.9453E-01, 4.4852E-01, 3.1375E-01, | 2.1065E-01, 1.3653E-01, 9.1144E-02, 6.1994E-02, 4.4571E-02/) corrnormpara(3,:) = (/ | 1.0270E+00, 7.8203E-01, 5.9534E-01, 4.5240E-01, 3.2403E-01, | 2.2302E-01, 1.4233E-01, 8.9070E-02, 5.6653E-02, 3.8611E-02/) corrnormpara(4,:) = (/ | 1.0373E+00, 7.8500E-01, 5.9297E-01, 4.4837E-01, 3.2384E-01, | 2.2687E-01, 1.4450E-01, 8.4366E-02, 4.9511E-02, 3.1025E-02/) corrnormpara(5,:) = (/ | 1.0490E+00, 7.8929E-01, 5.9130E-01, 4.4231E-01, 3.1728E-01, | 2.2191E-01, 1.4024E-01, 7.6380E-02, 4.1139E-02, 2.2905E-02/) corrnormpara(6,:) = (/ | 1.0617E+00, 7.9456E-01, 5.9057E-01, 4.3581E-01, 3.0571E-01, | 2.0872E-01, 1.2728E-01, 6.2907E-02, 3.0269E-02, 1.4243E-02/) corrnormpara(7,:) = (/ | 1.0773E+00, 8.0212E-01, 5.9180E-01, 4.2987E-01, 2.8925E-01, | 1.8822E-01, 1.0372E-01, 3.4590E-02, 1.2018E-02, 3.8835E-03/) corrnormpara(8,:) = (/ | 1.1003E+00, 8.1612E-01, 5.9836E-01, 4.2760E-01, 2.8986E-01, | 1.8413E-01, 9.7462E-02, 4.0988E-02, 7.2898E-03, -9.4405E-03/) corrnormpara(9,:) = (/ | 1.1409E+00, 8.4433E-01, 6.1714E-01, 4.3559E-01, 2.9400E-01, | 1.8117E-01, 8.9445E-02, 3.9187E-02, -2.7993E-03, -2.3482E-02/) corrnormpara(10,:) = (/ | 1.2102E+00, 8.9617E-01, 6.5632E-01, 4.6226E-01, 3.0804E-01, | 1.8559E-01, 8.5822E-02, 3.2302E-02, -1.1882E-02, -3.1824E-02/) corrnormpara(11,:) = (/ | 1.3155E+00, 9.7740E-01, 7.2026E-01, 5.1267E-01, 3.4425E-01, | 2.1050E-01, 1.0137E-01, 4.2290E-02, -6.6371E-03, -2.9570E-02/) corrnormpara(12,:) = (/ | 1.4568E+00, 1.0872E+00, 8.0860E-01, 5.8668E-01, 3.7564E-01, | 2.2742E-01, 1.1522E-01, 6.0656E-02, 1.5509E-02, -1.3943E-02/) corrnormpara(13,:) = (/ | 1.6291E+00, 1.2209E+00, 9.1448E-01, 6.7723E-01, 4.4810E-01, | 2.8628E-01, 1.6087E-01, 8.2129E-02, 3.7453E-02, 1.2558E-02/) corrnormpara(14,:) = (/ | 1.8231E+00, 1.3678E+00, 1.0284E+00, 7.7316E-01, 5.4154E-01, | 3.6841E-01, 2.3079E-01, 1.3201E-01, 7.5752E-02, 4.5007E-02/) corrnormpara(15,:) = (/ | 2.0245E+00, 1.5178E+00, 1.1409E+00, 8.6451E-01, 6.2741E-01, | 4.4552E-01, 2.9960E-01, 1.9064E-01, 1.2103E-01, 7.7739E-02/) corrnormpara(16,:) = (/ | 2.2205E+00, 1.6645E+00, 1.2499E+00, 9.4698E-01, 7.0103E-01, | 5.1247E-01, 3.5917E-01, 2.4008E-01, 1.5973E-01, 1.0640E-01/) corrnormpara(17,:) = (/ | 2.4095E+00, 1.8035E+00, 1.3532E+00, 1.0211E+00, 7.6117E-01, | 5.6591E-01, 4.0665E-01, 2.7886E-01, 1.9017E-01, 1.2938E-01/) corrnormpara(18,:) = (/ | 2.5945E+00, 1.9315E+00, 1.4449E+00, 1.0869E+00, 8.0855E-01, | 6.0410E-01, 4.3994E-01, 3.0564E-01, 2.1118E-01, 1.4559E-01/) corrnormpara(19,:) = (/ | 2.7711E+00, 2.0436E+00, 1.5179E+00, 1.1386E+00, 8.4343E-01, | 6.2755E-01, 4.5825E-01, 3.2012E-01, 2.2249E-01, 1.5473E-01/) corrnormpara(20,:) = (/ | 2.9222E+00, 2.1322E+00, 1.5680E+00, 1.1720E+00, 8.6482E-01, | 6.3853E-01, 4.6350E-01, 3.2373E-01, 2.2527E-01, 1.5753E-01/) corrnormpara(21,:) = (/ | 3.0229E+00, 2.1936E+00, 1.5980E+00, 1.1857E+00, 8.7186E-01, | 6.3904E-01, 4.5877E-01, 3.1914E-01, 2.2164E-01, 1.5567E-01/) corrnormpara(22,:) = (/ | 3.0587E+00, 2.2270E+00, 1.6153E+00, 1.1860E+00, 8.6554E-01, | 6.3031E-01, 4.4655E-01, 3.0841E-01, 2.1341E-01, 1.5074E-01/) corrnormpara(23,:) = (/ | 3.0304E+00, 2.2330E+00, 1.6267E+00, 1.1810E+00, 8.4969E-01, | 6.1375E-01, 4.2865E-01, 2.9308E-01, 2.0168E-01, 1.4361E-01/) corrnormpara(24,:) = (/ | 2.9488E+00, 2.2160E+00, 1.6360E+00, 1.1770E+00, 8.2992E-01, | 5.9196E-01, 4.0652E-01, 2.7399E-01, 1.8695E-01, 1.3468E-01/) corrnormpara(25,:) = (/ | 2.8341E+00, 2.2280E+00, 1.6729E+00, 1.1866E+00, 8.3395E-01, | 5.9046E-01, 4.0163E-01, 2.7945E-01, 1.8343E-01, 1.2405E-01/) corrnormpara(26,:) = (/ | 2.6932E+00, 2.1887E+00, 1.6860E+00, 1.2020E+00, 8.3694E-01, | 5.8129E-01, 3.8726E-01, 2.7069E-01, 1.7119E-01, 1.1151E-01/) corrnormpara(27,:) = (/ | 2.5405E+00, 2.0861E+00, 1.6620E+00, 1.2163E+00, 8.0881E-01, | 5.3543E-01, 3.3890E-01, 2.1588E-01, 1.3605E-01, 9.7129E-02/) corrnormpara(28,:) = (/ | 2.3765E+00, 2.0097E+00, 1.6517E+00, 1.2423E+00, 8.4083E-01, | 5.4485E-01, 3.3131E-01, 2.1040E-01, 1.2250E-01, 8.1111E-02/) corrnormpara(29,:) = (/ | 2.2055E+00, 1.9224E+00, 1.6294E+00, 1.2636E+00, 8.7801E-01, | 5.6100E-01, 3.2565E-01, 2.0542E-01, 1.0852E-01, 6.4218E-02/) corrnormpara(30,:) = (/ | 2.0096E+00, 1.7961E+00, 1.5707E+00, 1.2597E+00, 8.8893E-01, | 5.5944E-01, 3.0422E-01, 1.8113E-01, 8.4620E-02, 4.6343E-02/) corrnormpara(31,:) = (/ | 1.8388E+00, 1.6814E+00, 1.5111E+00, 1.2497E+00, 8.9772E-01, | 5.6143E-01, 2.8445E-01, 1.5706E-01, 6.0365E-02, 2.9833E-02/) corrnormpara(32,:) = (/ | 1.6396E+00, 1.5305E+00, 1.4078E+00, 1.1951E+00, 8.7151E-01, | 5.4223E-01, 2.4998E-01, 1.2041E-01, 2.5746E-02, 7.6820E-03/) corrnormpara(33,:) = (/ | 1.5550E+00, 1.4667E+00, 1.3663E+00, 1.1810E+00, 8.7444E-01, | 5.4784E-01, 2.3483E-01, 9.4771E-02, -5.0855E-03, -1.1990E-02/) corrnormpara(34,:) = (/ | 1.5314E+00, 1.4440E+00, 1.3477E+00, 1.1757E+00, 8.8300E-01, | 5.6149E-01, 2.2941E-01, 7.1049E-02, -4.1493E-02, -3.9332E-02/) corrnormpara(35,:) = (/ | 1.7270E+00, 1.6064E+00, 1.4824E+00, 1.2918E+00, 9.9076E-01, | 6.5768E-01, 2.9160E-01, 9.1000E-02, -5.7683E-02, -6.3116E-02/) corrnormpara(36,:) = (/ | 2.0397E+00, 1.8621E+00, 1.6875E+00, 1.4588E+00, 1.1422E+00, | 7.9506E-01, 3.9605E-01, 1.3928E-01, -6.2587E-02, -9.5685E-02/) corrnormpara(37,:) = (/ | 2.2668E+00, 2.0382E+00, 1.8186E+00, 1.5593E+00, 1.2376E+00, | 8.9109E-01, 4.7986E-01, 1.7537E-01, -7.7256E-02, -1.5153E-01/) corrnormpara(38,:) = (/ | 2.3460E+00, 2.0887E+00, 1.8443E+00, 1.5690E+00, 1.2543E+00, | 9.1963E-01, 5.1685E-01, 1.8720E-01, -9.8736E-02, -2.1873E-01/) corrnormpara(39,:) = (/ | 2.3330E+00, 2.0637E+00, 1.8110E+00, 1.5281E+00, 1.2244E+00, | 9.0429E-01, 5.2074E-01, 1.9076E-01, -1.0517E-01, -2.6267E-01/) corrnormpara(40,:) = (/ | 2.2348E+00, 1.9681E+00, 1.7201E+00, 1.4351E+00, 1.1484E+00, | 8.4754E-01, 4.9071E-01, 1.8348E-01, -9.6347E-02, -2.6481E-01/) corrnormpara(41,:) = (/ | 2.1025E+00, 1.8445E+00, 1.6021E+00, 1.3071E+00, 1.0438E+00, | 7.6659E-01, 4.4062E-01, 1.6942E-01, -7.5640E-02, -2.2199E-01/) corrnormpara(42,:) = (/ | 1.9475E+00, 1.7005E+00, 1.4648E+00, 1.1508E+00, 9.1109E-01, | 6.6102E-01, 3.6639E-01, 1.3488E-01, -6.7300E-02, -1.6573E-01/) corrnormpara(43,:) = (/ | 1.7615E+00, 1.5319E+00, 1.3075E+00, 9.7985E-01, 7.5181E-01, | 5.1793E-01, 2.4758E-01, 5.6814E-02, -1.0011E-01, -1.3994E-01/) corrnormpara(44,:) = (/ | 1.5496E+00, 1.3439E+00, 1.1379E+00, 8.1610E-01, 5.7226E-01, | 3.1847E-01, 4.8379E-02, -9.4852E-02, -1.9907E-01, -1.7544E-01/) corrnormpara(45,:) = (/ | 1.3282E+00, 1.1519E+00, 9.7038E-01, 6.7744E-01, 3.8117E-01, | 4.9884E-02, -2.7429E-01, -3.5478E-01, -3.8500E-01, -2.7837E-01/) corrnormpara(46,:) = (/ | 1.1036E+00, 9.5974E-01, 8.0834E-01, 5.6022E-01, 1.8634E-01, | -2.7016E-01, -7.2208E-01, -7.3811E-01, -6.7652E-01, -4.4553E-01/) corrnormpara(47,:) = (/ | 8.8569E-01, 7.7516E-01, 6.5602E-01, 4.5814E-01, 2.0128E-02, | -5.5226E-01, -1.1772E+00, -1.1949E+00, -1.0850E+00, -6.9480E-01/) corrnormpara(48,:) = (/ | 6.9026E-01, 6.0924E-01, 5.1984E-01, 3.6822E-01, -7.9705E-02, | -6.9601E-01, -1.4791E+00, -1.6532E+00, -1.6410E+00, -1.1362E+00/) corrnormpara(49,:) = (/ | 5.1458E-01, 4.5782E-01, 3.9370E-01, 2.8169E-01, -1.1179E-01, | -6.7907E-01, -1.5384E+00, -2.0210E+00, -2.3049E+00, -1.9254E+00/) corrnormpara(50,:) = (/ | 3.6681E-01, 3.2819E-01, 2.8391E-01, 2.0345E-01, -1.0147E-01, | -5.6092E-01, -1.3819E+00, -2.1793E+00, -2.8408E+00, -2.9029E+00/) corrnormpara(51,:) = (/ | 2.4572E-01, 2.1994E-01, 1.8976E-01, 1.3317E-01, -8.4483E-02, | -4.2521E-01, -1.1176E+00, -2.0715E+00, -2.9853E+00, -3.5512E+00/) corrnormpara(52,:) = (/ | 1.3765E-01, 1.2118E-01, 1.0134E-01, 6.3601E-02, -8.1307E-02, | -3.1424E-01, -8.3137E-01, -1.7228E+00, -2.6429E+00, -3.4631E+00/) corrnormpara(53,:) = (/ | 3.9753E-02, 2.9513E-02, 1.7194E-02, -6.4057E-03, -9.6451E-02, | -2.4263E-01, -5.8604E-01, -1.2695E+00, -2.0044E+00, -2.7729E+00/) corrnormpara(54,:) = (/ | -4.9228E-02, -6.7635E-02, -8.4465E-02, -1.0228E-01, -1.5492E-01, | -2.4047E-01, -4.4807E-01, -9.0507E-01, -1.4082E+00, -1.9908E+00/) corrnormpara(55,:) = (/ | -1.2942E-01, -1.5302E-01, -1.7275E-01, -1.8669E-01, -2.1598E-01, | -2.6364E-01, -3.8122E-01, -6.5929E-01, -9.6948E-01, -1.3556E+00/) corrnormpara(56,:) = (/ | -1.9101E-01, -2.0825E-01, -2.2230E-01, -2.3063E-01, -2.4588E-01, | -2.7113E-01, -3.3375E-01, -4.8833E-01, -6.6192E-01, -8.8584E-01/) corrnormpara(57,:) = (/ | -2.2531E-01, -2.3179E-01, -2.3703E-01, -2.4020E-01, -2.4747E-01, | -2.5995E-01, -2.9111E-01, -3.6874E-01, -4.5632E-01, -5.6825E-01/) corrnormpara(58,:) = (/ | -2.3044E-01, -2.3034E-01, -2.3035E-01, -2.3068E-01, -2.3375E-01, | -2.3952E-01, -2.5404E-01, -2.8977E-01, -3.3031E-01, -3.7962E-01/) corrnormpara(59,:) = (/ | -2.1246E-01, -2.1209E-01, -2.1172E-01, -2.1162E-01, -2.1269E-01, | -2.1506E-01, -2.2143E-01, -2.3725E-01, -2.5530E-01, -2.7715E-01/) corrnormpara(60,:) = (/ | -1.7345E-01, -1.7352E-01, -1.7358E-01, -1.7338E-01, -1.7352E-01, | -1.7428E-01, -1.7682E-01, -1.8326E-01, -1.9064E-01, -1.9988E-01/) corrnormpara(61,:) = (/ | -1.2718E-01, -1.2718E-01, -1.2722E-01, -1.2698E-01, -1.2675E-01, | -1.2695E-01, -1.2775E-01, -1.2992E-01, -1.3249E-01, -1.3566E-01/) corrnormpara(62,:) = (/ | -8.6627E-02, -8.6640E-02, -8.6627E-02, -8.6447E-02, -8.6173E-02, | -8.6140E-02, -8.6281E-02, -8.6855E-02, -8.7563E-02, -8.8431E-02/) end subroutine heatnirco2_init end module heatnirco2_module